Papers

We model Herman Kahn's escalation ladder as an impartial combinatorial game. Reindexing each rung by its distance to the nuclear threshold turns the ladder into a subtraction game, the most tractable class in combinatorial game theory, and the doctrinal fact that no side wishes to fire first selects the misere convention. We prove that single-ladder stability is governed by a congruence (Theorem 4.1) and derive a ladder-design corollary that makes the burden of first escalation a function of ladder length and escalation granularity (Corollary 4.2). For simultaneous theaters we show, under normal play, that joint stability is the Nim-sum of the theater-wise escalation distances (Theorem 5.2), a condition that is neither additive nor dominated by the most dangerous theater. We then show the Nim-sum reduction fails under misere play, introduce the misere quotient as its replacement, and prove by exhaustive backward induction that for two-step escalation the quotient is the order-six monoid with generators a, b satisfying a^2=1 and b^3=b, with loss set {a,b^2} (Theorem 6.3). To our knowledge, impartial combinatorial game theory has not previously been applied to nuclear escalation ladders; the existing game-theoretic literature on escalation is classical and payoff-based.

We analyze via Evolutionary Gamma-convergence a stratified composite structure consisting of a thin adhesive layer with vanishing thickness and undergoing rate-independent damage, as well as two adjacent elastic adherents. As the width of the intermediate layer tends to zero, we prevent complete degradation of the material by assuming that the damage variable scales minimally like the thickness of the adhesive layer. As a result, we identify a limiting model that combines both a brittle constraint and an adhesive-type energy contribution featuring the jump of the admissible displacements.

Numerical simulation of quantum computing hardware and open quantum systems governed by the Lindblad equation is challenging due to the high dimensionality of the density matrix and the need to preserve fundamental physical properties. In our previous work, we developed an arbitrary-order, low-rank, completely positive and trace preserving (CPTP) method for the Lindblad equation with time-dependent Hamiltonians by nested Picard iteration (NPI). In this work, we develop Gregory NPI schemes, which are CPTP schemes constructed by Gregory-type quadrature on equispaced nodes. The methods, which are of order up to nine, substantially reduce the computational cost compared to our previously proposed NPI schemes with Gaussian quadrature rules, while retaining high-order accuracy and structure preservation. We analyze the stability of the resulting scheme for a physics-based test equation. Numerical experiments verify the convergence of the method and demonstrate the effectiveness of the low-rank approximation. We study the performance of a previously constructed CNOT gate for both closed and open quantum systems.

This paper presents a two-dimensional method-of-moments (MoM) solver for electromagnetic scattering from infinitely long perfectly electrically conducting (PEC) cylinders. Both TMz and TEz polarizations are considered. Starting from the scalar Helmholtz equation, the electric field integral equation (EFIE) is derived for TMz scattering and the magnetic field integral equation (MFIE) is derived for TEz scattering. The induced surface current on the PEC boundary is expanded using pulse basis functions, and the boundary integral equations are discretized using point matching at the segment centers. Circular cylinders with radii $R = {\lambda}$ and $R = 2{\lambda}$ are used as validation cases because analytical series solutions are available. The MoM-computed surface currents, total near fields, scattered near fields, and field-error distributions are compared against the analytical solutions. After validation, the same solver is applied to a square PEC cylinder, for which no simple closed-form analytical solution is used. The results show strong agreement between the MoM and analytical circular-cylinder solutions and demonstrate the geometry-dependent scattering behavior of the square cylinder.

Let $A$ be an abelian variety of dimension $g$ over a finite field $\mathbf{F}_q$. We show that if $q$ is sufficiently large relative to $g$, the $g$ point counts $\#A(\mathbf{F}_{q^i})$ for $1 \leq i \leq g$ determine the zeta function of $A$, equivalently the characteristic polynomial of its Frobenius endomorphism, and hence the isogeny class of $A$. This count is best possible for $g=2$ and $g=4$, but not in general: for $g=3$ two point counts already determine the zeta function, whereas a single count never does. The proof combines the functional equation of the $L$-polynomial with Newton's identities and an inductive error analysis that controls the power sums of the inverse Frobenius eigenvalues with enough precision to recover them, as integers, by rounding.

Let $k\ge 2$ be fixed. We study the distribution modulo one of the $n^k$ sums \begin{equation*} \sqrt{a_1} + \cdots + \sqrt{a_k}, \qquad 1\le a_1, \dots, a_k \le n, \end{equation*} counted with multiplicity. For \begin{equation*} S(h,n) = \sum_{n/2\le a\le n} \mathbf{e}(h\sqrt{a}), \qquad \mathbf{e}(x) = \exp(2\pi i x), \end{equation*} we prove second- and fourth-moment estimates matching the diagonal scale up to a factor $n^\varepsilon$. More precisely, \begin{equation*} \sum_{H/2\le h\le H} \left| S(h,n) \right|^2 \ll_{\varepsilon,\delta} Hn^{1+\varepsilon} \end{equation*} uniformly for $H\ge n^{1/2+\delta}$, and \begin{equation*} \sum_{H/2\le h\le H} \left| S(h,n) \right|^4 \ll_{\varepsilon,\delta} Hn^{2+\varepsilon} \end{equation*} uniformly for $n^{1/2+\delta} \le H \le n^{2/3}$, where $0<\delta<1/6$ in the fourth-moment estimate. Combining the second-moment bound with pointwise exponential-sum estimates and the Erdős–Turán inequality, we obtain \begin{equation*} D_k(n) \le n^{-\rho_k+o(1)}, \qquad \rho_k = \frac{71k+26}{26k+116}, \end{equation*} as $n\to\infty$, where $D_k(n)$ denotes the discrepancy with respect to arbitrary subintervals of $[0,1)$.

We propose a model of higher homotopy theory of $L_{\infty}[1]$-morphisms as a natural generalization of the $A_{\infty}$-homotopies defined by Fukaya-Oh-Ohta-Ono \cite{FOOO1}. Within this framework, we show that a filling condition holds for simplices whose vertices are assigned quasi-isomorphisms.

We give a single compositional setting in which gradient-based learning and Hamiltonian-style mechanics appear as functorial semantics. The syntax is an operad Arr whose objects are input-output interfaces (pairs of manifolds) and whose morphisms are *smooth adaptive arrangements*, which consist of a responsive parameter space, a lens given by smooth output and input maps, and a real-valued potential. The main technical result of the paper is what we call *lens internalization*, a lax symmetric monoidal functor Lens(C) $\to$ C associated to any symmetric monoidal closed category C. Using it, we provide two functors $\Phi_\text{phase}$, $\Phi_\text{conf}$: Arr $\to$ PC into the 2-category of polynomial coalgebras – input-output discrete dynamical systems – which we take as the semantics category. $\Phi_\text{phase}$ stores both position and momentum, whereas $\Phi_\text{conf}$ stores only position. When applied to a parameterized function, $\Phi_\text{conf}$ recovers the gradient descent training algorithm, with backpropagation as the lens' backward pass. When applied to harmonic particles wired together – in series, or according to any finite directed graph – one diagram yields two different regimes, both of which are governed by the graph Laplacian: $\Phi_\text{phase}$ gives the discrete wave equation, which is conservative and second-order, and $\Phi_\text{conf}$ gives the discrete heat equation, which is dissipative and first-order. They are two semantics of one adaptive arrangement, e.g. with the same potential in each case. And because Arr is an operad, such diagrams nest – larger systems wired from smaller ones – and each semantics assembles a system's dynamics functorially from its parts. These dynamics are moreover executable: a parameterized neural network and a graph of particles both compile, by the same construction, to explicit state machines one can run.

We study the static electromagnetic response of two spherical topological insulators embedded in a dielectric medium and subjected to a uniform external electric field. The gapped surface states are described by a piecewise constant axion field, which induces a topological magnetoelectric coupling localized at the spherical interfaces. {More generally, the same formalism applies to isotropic magnetoelectric media characterized by an effective scalar magnetoelectric response.} The electrostatic problem is solved at zeroth order using bispherical coordinates, allowing for an exact treatment of both parallel and perpendicular orientations of the external field relative to the center-to-center axis. The resulting mode expansions are determined by three-term recurrence relations, which are solved perturbatively for nonoverlapping spheres. The { magnetoelectric}-induced response is then computed to leading order in the fine-structure constant {(or, more generally, in the effective coupling strength)}. The induced sources are purely interfacial and generate distinct magnetostatic field configurations in the parallel and perpendicular geometries. Closed-form series representations for the induced vector potential and magnetic field are obtained in terms of the zeroth-order electrostatic coefficients. These results provide an analytically controlled description of {interaction-induced magnetostatics in coupled spherical magnetoelectric systems}.

We say that a square matrix $A$ is \emph{$1$-principled} if all its principal minors are equal to $1$. We show that over any well-behaved ring, any power $A^m$ of a $1$-principled matrix $A$ is again $1$-principled. Well-behaved rings include all reduced rings as well as all quotients of commutative rings modulo integrally closed ideals; in particular, all fields and all quotients of $\mathbb{Z}$ are well-behaved. We note that $m$ can be any integer, positive or negative. This generalizes Problem B5 of the 2021 Putnam contest in multiple directions. Over arbitrary commutative rings, we identify a stronger property that is always inherited by powers: We say that a matrix $A = \left(a_{i,j}\right)_{i,j\in\left[n\right]}$ is \emph{$1$-nullcyclic} if all its diagonal entries are $1$ and if all the cyclic products $a_{i_1, i_2} a_{i_2, i_3} \cdots a_{i_k, i_1}$ with $k>1$ and distinct $i_1,i_2,\ldots,i_k$ vanish. We show that if $A$ is $1$-nullcyclic, then so is $A^m$ for any integer $m$. Furthermore, every $1$-nullcyclic matrix is $1$-principled over any commutative ring, while the converse holds if the ring is well-behaved. Along the way, we prove analogous results that don't require the diagonal entries to be $1$. These are concerned with principled matrices (those whose principal minors equal the respective products of diagonal entries) and nullcyclic matrices (those whose cyclic products $a_{i_1, i_2} a_{i_2, i_3} \cdots a_{i_k, i_1}$ with $k>1$ and distinct $i_1,i_2,\ldots,i_k$ vanish); their diagonal entries can be arbitrary. A crucial auxiliary result, which holds for any $n\times n$-matrix $A$, is that the cyclic products $a_{i_1, i_2} a_{i_2, i_3} \cdots a_{i_k, i_1}$ (with $k>1$ and distinct $i_1,i_2,\ldots,i_k$) are integral over the ideal generated by the principal minors of $A$ minus the corresponding products of diagonal entries of $A$.

We propose a cohomological approach to computing the conformal anomaly. Using the example of current-current deformations of two-dimensional conformal field theories, we reproduce the well-known Cardy formula for the leading contribution to the perturbative beta function as the coefficient of the cocycle that realizes the obstruction to deforming the conformal algebra module structure on the state space. In addition to offering a novel conceptual perspective on the conformal anomaly, the proposed approach is anticipated to provide an efficient tool for computing higher-order coefficients of perturbative beta functions.

When datasets contain outliers, robust regression is a well-established alternative to Ordinary Least Squares. A commonly employed robust estimator is Least Trimmed Squares (LTS), which computes the regression coefficients from a subset of observations. Determining the exact solution corresponds to a combinatorial problem with prohibitive computational costs, even for instances of moderate dimension. Thus, the most prevalent approach in practice remains a heuristic known as Fast-LTS. Although the heuristic often performs effectively, certain elements of the approach remain open to improvement. In particular, its core procedure provides robust results only when initialized with a large number of starting points. To address the heuristic's limitations, this paper reformulates the LTS problem as a concave minimization problem subject to a capped simplex constraint, and proposes the successive Boosted Difference of Convex Functions Algorithm (sBDCA) as a solution method. Theoretically, we establish via the Łojasiewicz property that sBDCA converges to a local solution with a linear rate in the fastest case. To ensure robustness from a single initialization in practice, we derive and integrate a problem-specific preconditioning matrix into the algorithmic setup. Building on this theoretical foundation, we conduct numerical studies on various synthetic and real-world datasets to demonstrate the effectiveness of sBDCA with preconditioning. Specifically, we show that our approach is up to 3.25 times faster than Fast-LTS and achieves up to 90% lower objective function values, particularly in high-dimensional settings. As all code is openly available, this paper further provides a practical guide to robust regression in Python.

We prove that the Sinkhorn algorithm converges at the rate of $O(1/k)$ in $\ell_1$-norm marginal error and in joint relative entropy, which is known to be sharp in the asymptotically scalable case. The proof is based on examining the bipartite graph associated to the entropy-regularized optimal transport problem, and treating differently the edges that are assigned a positive mass in the optimal transport plan vs. those that are not. This yields a local convergence bound with the sharp rate, which is bootstrapped into a global bound using the author's previous result in arXiv:2604.26265 where we showed an almost-sharp rate up to a logarithmic factor.

An analytic solution has been recently developed for the Maxwell's equation in a medium with zero currents such as vacuum. The solution is attractive in the sense that it is formulated based on the Fourier expansion of the initial value. It has been used to study the properties of solutions like certain conservative laws and construct electromagnetic waves with certain features. In this paper, we study Maxwell's equation in a medium with non-zero currents. The structure of solutions in this setting turns out to be much more complicated than what has been achieved without currents, and a clean structure of analytic solutions as with zero current is no longer available in general. Nevertheless, we can still develop an algorithm to construct the solution effectively. Our efforts in seeking analytic solution focus on two special cases. First, we develop analytic solution under the assumption that Ohm's law is satisfied, i.e. the current density is proportional to electronic density; secondly, we add skew symmetric components under generalized Ohm's law, which is also refereed as Hall effect in literature, and study the properties of solutions. In addition, we consider the case where an independent local electromagnetic field is included and derive the analytical solution accordingly. As an application, we provide an example to use the analytic solution to construct parallel electronic and magnetic waves.

An explicit cubic Ramanujan–Sato formula for $1/\pi$ on $\Gamma_0(2)^+$ at $D=-163$ is presented. The construction produces a very small cubic CM parameter, giving about $15.01$ decimal digits of geometric contraction per term. This is slightly sharper than the classical Chudnovsky contraction, but the coefficients lie in a cubic algebraic field.

We study second-order fixed-ray Dini sensitivity of the efficient solution map \(S\) for the parametric vector problem \(\min_C f(p,x)\) subject to \(x\in H(p)\). Under a value-to-decision error bound (VDB) – an inverse stability estimate – the second-order semi-derivative of the marginal map \(\Phi\) transfers to an exact second-order Dini formula for \(S\). The result separates value sensitivity from the inverse recovery of efficient decisions. This mechanism allows a nonsingleton efficient-value fiber and requires neither injectivity nor strict monotonicity of \(\nabla_x f\). For structured feasible maps \(H(p)=\{x\in\Omega: g(p,x)\in D\}\), we derive primitive-data formulas for \(D^2H\), \(D^2\Phi\), and \(D^2S\) under Robinson metric regularity, second-order regularity of \(\Omega\) and \(D\), objective-aware recession conditions, and directional second-order semi-derivability of the data. Polyhedral inequality/equality systems give explicit specializations. A parametric multi-objective portfolio family verifies (VDB) with an explicit constant, obtained directly from the model's feasible directions, while a two-branch example shows how the formula operates on a nonsingleton decision fiber.

This first part of the series builds the analytic layer of noncommutative anisotropic diffusion in a separable Hilbert space. Let $\mu_0=\mathcal{N}(0,Q)$ be the reference Gaussian measure, with $Q\in L^1(\mathcal{H})$, and let $D(x)$ be a positive, state-dependent anisotropy. We do not assume that $[D(x),Q]=0$. Consequently, for the forward SDE with $\sigma(x)=D(x)^{1/2}Q^{1/2}$, the correct energy form is given not by the expression $\langle D\nabla u,\nabla v\rangle$ but by the consistent form $\Gamma_A(u,v)= \langle Q^{1/2}D(x)^{1/2}\nabla u, Q^{1/2}D(x)^{1/2}\nabla v\rangle$. We prove closability of the form, well-posedness of the forward dynamics, Galerkin convergence, stability of the $A$-LSI under a Mosco limit, the chain rule for relative entropy, and a general weak-bridge theorem. The main result of Part~I is a functional-analytic theorem: if $A$-consistency, a uniform $A$-LSI, and representability of the right-hand side of the backward weak form in the negative energy space all hold, then a backward weak drift $v=\mathsf{A}\nabla\Phi$ exists and the basic entropy dissipation estimate holds. In addition, we single out a three-dimensional tensor class of anisotropies, formulate a condition for the absence of diffusion degeneracy, and obtain a rate estimate for the homogenization limit, first on cylindrical subspaces and then on compact-tail classes, which yields strong resolvent convergence and convergence of the forward SDEs. The statistical closure, an independent isotropic benchmark, and an approximation theorem for $A$-adapted networks are treated in Part~II.

We establish Liouville-type theorems for the stationary fractional Navier-Stokes equations in $\mathbb{R}^n$ under suitable integrability conditions on the velocity field $u$ and a large-scale Morrey-type bound on the fractional energy. As a corollary, these assumptions are automatically satisfied if $u \in \dot{H}^{\frac{\alpha}{2}}(\mathbb{R}^n)$, yielding Liouville-type results under the finite fractional energy condition for $\frac{n}{3} \le \alpha < \frac{n+2}{3}$, where $\alpha$ denotes the order of the fractional Laplacian $(-\Delta)^{\frac{\alpha}{2}}$. This range reflects a scaling-critical correspondence between Liouville-type theorems in the finite-energy setting and the threshold arising in partial regularity theory. The proof relies on direct kernel estimates for the commutator of the fractional Laplacian, based on a dyadic decomposition of the tail term, which remain valid in the hyper-dissipative case. The argument also uses a bootstrap argument that propagates integrability from near the scaling-invariant exponent down to lower exponents, including the Sobolev embedding exponent.

In this paper, we provide an integration by parts formula for plurisubharmonic functions on a hyperconvex domain that are bounded outside a compact set. This extends a previous result of Urban Cegrell.

In the paper we describe hyperbolic surfaces filled by their systoles, where the total number of systoles is in $O(\frac{g}{\ln \,g})$, that is equivalent to the lower bound of Anderson, Parlier and Pittet \cite{APP}. Various papers \cite{SS}\cite{FB20}\cite{Sanki}\cite{ IM}\cite{ Mathieu} have investigated the same question, and the best previously known upper bounds where in $o(\frac{g}{{\sqrt{\ln \,g}}})$. Surprizingly the present approach is, in our opinion, much simpler than the methods of earlier papers.

This paper focuses on the regularity of viscosity solutions to normalized $p$-Laplacian equations with variable-exponent double phase type degeneracy/singularity and Hamiltonian terms. Based on a new improved oscillation-type estimate combined with a localized analysis, we establish sharp interior $C^{1,\alpha}$ regularity estimates in a unified way. In addition, in the degenerate case, we obtain improved gradient Hölder regularity results at points where the Hamiltonian coefficient and source term vanish, and establish a Schauder-type estimate at local extrema. Notably, our results are still novel even restricted to single power-type singularity or degeneracy law.

We solve large-$N$ Yang–Mills theory on $\mathbb{Z}^d$, for every $d\geq2$, at strong coupling, for structure group $\mathrm U(N)$ and for the heat-kernel action. More precisely, we prove that normalized Wilson loop expectations have infinite-volume large-$N$ limits, factorize at leading order, and admit an all-order $1/N$-expansion with exponentially local coefficients, whose leading order characterizes the master field. We also prove an area-law upper bound for the heat-kernel master field, with a stronger coefficientwise version. The proof is based on a rooted heat-kernel master loop equation. Unlike the Wilson-action equation or the two-dimensional Makeenko–Migdal equation, this equation does not close on Wilson loop observables alone; it closes on an extended space of loop observables coupled to compactly supported plaquette decorations. We prove a strong-coupling, order-truncated rooted trajectory expansion and then identify its leading term with the master field. The main inputs are the universal finite-$N$ duality formulas developed in the companion paper \cite{Lem26a} and large-$N$ heat-kernel estimates from \cite{LemMai25,LM2}.

In this paper, we study a quantitative Runge-type global approximation theorem for the linearized magnetohydrodynamic (MHD) system in bounded domains with arbitrary topology. In the context of magnetic relaxation, the interplay between the domain topology and magnetic field structure plays a crucial role. Recent studies illustrate a sharp contrast in the dynamics: while Enciso–Peralta-Salas (2025) highlights that the geometric complexity of magnetic fields acts as an obstruction to relaxation in non-resistive regimes, Kozono-Shimizu-Yanagisawa (2025) proves that in resistive regimes, the flow stably relaxes towards a harmonic equilibrium. Focusing on this resistive scenario, we adopt a control-theoretic viewpoint to quantitatively approximate the relaxation trajectory generated by the linearized initial-boundary value problem. Specifically, after decomposing the bounded-domain solution into the time-evolving part and the stationary part, we approximate it by a global solution on $\mathbb{R}^3$ under a remote forcing. An explicit dependence of the forcing cost on the approximation error is provided.

Urban vehicular networks (VNs) demand seamless connectivity and situational awareness within road-constrained environments, motivating the deployment of unmanned aerial vehicles (UAVs) platforms capable of simultaneously sensing vehicles and establishing communication with them. In this paper, we present a sensing-assisted UAV network that provides connectivity to the vehicles in an urban area. The road network of the urban area is modeled as Manhattan Poisson line process (MPLP), and the random location of vehicles on each road is modeled as one dimensional Poisson point processes (PPPs). UAVs are distributed in the urban area at a fixed altitude and provide connectivity after sensing the vehicles. Their locations are modeled as a two-dimensional homogeneous PPP. Combined with the fixed altitude, this results in a three-dimensional spatial configuration. We incorporate an elevation dependent blockage model and define the sensing radius based on detection probability (DP), showing that it is jointly limited by signal strength and blockage effects. We derive the DP and characterize the typical UAV's sensing region within the reliability requirements. We also derive the Laplace transform (LT) of aggregate interference accounting for directional patterns and sensing-driven activity, and analyze the resulting coverage probability (CP). Finally, we obtain the rate coverage (RC) of sensed vehicles falling within the UAV's sensing zone. Numerical results shows that increasing altitude degrades sensing and coverage performance, whereas RC exhibits a non-monotonic trend, first decreasing and then increasing with altitude.

The hydrodynamic limit to the barotropic Euler equations, including power-law pressure $P(\rho)=\rho^\gamma$, for a kinetic nonlinear Fokker–Planck equation with degenerate diffusion is established. This extends the well-known result of the derivation of isothermal Euler equations via Fokker–Planck equation with linear diffusion. We establish the asymptotic analysis using the relative entropy method by quantifying error estimates for pressures and employing the generalized Log-Sobolev inequality for degenerate diffusion.

In-context learning (ICL) is an emerging paradigm that employs the semantic information inherent in large language models (LLMs) for generating answers to user queries. While the remarkable performance of ICL has been widely known, a general modeling and a rigorous theoretical analysis of this paradigm are still lacking. This work presents a probabilistic model for ICL and derives the performance of ICL for both general parametric distributions and exponential families. Based on the derived results, the work explains the impact of multiple factors such as the number of demonstrations, the sensitivity of the probabilistic model to the variation of its parameters, as well as the similarity between the demonstrations and the query on the performance of ICL.

We propose and analyze a space-time discontinuous Galerkin method for the incompressible Stokes equations on moving domains within the arbitrary Lagrangian-Eulerian setting. We use a contravariant Piola map in the definition of the discrete velocity space to preserve the pointwise divergence-free property on the discrete level. We show that the method is inf-sup stable, with no constraints on the spatial mesh or the time partition. We also establish a priori error estimates in the energy norm for arbitrary degrees of approximation in space and time. For piecewise-constant and piecewise-linear approximations in time, we show that the method is also robust at low viscosity regimes, and provide numerical evidence suggesting that this property extends to high-order cases as well. We present several numerical experiments to validate our theoretical findings.

In this paper, we consider distributed parameter estimation with binary observations under measurement-side tampering: each node observes a thresholded output whose label may be flipped and exchanges information over a communication graph. We develop a distributed recursive projection algorithm based on the diffusion strategy. Without imposing independence, stationarity, or Gaussian assumptions, we establish almost sure upper bounds of both the accumulated regrets of the adaptive predictor and the distributed estimation error. Under a mild cooperative excitation condition, all nodes' estimate are consistent, even when each node is individually non-exciting. Simulations on a jointly exciting network corroborate the theory and show that the proposed algorithm converges, whereas non-cooperative and tampering-unaware baselines do not.

This paper establishes a pointwise gradient potential estimate for solutions to linearized Monge-Ampère equations and derives a modulus of continuity estimate for the gradient in terms of the associated section adapted potential. As applications, we obtain a sufficient criterion for the continuity of the gradient, local BMO and VMO estimates for the gradient with respect to Monge-Ampère sections, and a local $C^{1,\gamma}$ regularity criterion.

We study doppelsemigroups, i.e., algebraic structures equip\-ped with two associative binary operations satisfying a specified system of axioms. We investigate duality and isomorphisms of doppelsemigroups and examine the relationships between commutative, abelian, strong, and rectangular doppelsemigroups. Several examples are constructed, including nontrivial iso-opposite doppelsemigroups, noncommutative iso-dual doppelsemigroups, nonabelian iso-cross-dual doppelsemigroups, and nonstrong rectangular iso-opposite doppelsemigroups. Furthermore, we refine the complete classification of nonisomorphic doppelsemigroups of order~3. Finally, we present computer-assisted calculations yielding the numbers of all pairwise nonisomorphic doppelsemigroups and strong doppelsemigroups of orders up to~$5$, as well as all pairwise nonisomorphic commutative, abelian, and rectangular doppelsemigroups of orders up to~$6$, obtained using GAP, Python, and C++.