Mamba-3 and the exponential-trapezoidal integrator

10.1 Two upgrades to the selective SSM

Chapter 9 built the selective SSM: input-dependent $(\stepsize_t, \inputmat_t, \outputmat_t)$ make the dynamics linear time-varying, the convolution view collapses, and the selective scan plus the semiseparable-matrix view (SSD) carry the computation. Mamba-1 and Mamba-2 both discretize with zero-order hold (ZOH) and keep a real diagonal state matrix $\statemat$.

Mamba-3 Lahoti et al. (2026) changes two things of dynamical-systems substance and leaves the rest of the SSD scaffold (the semiseparable structure of Theorem 9.4, the scan, the attention duality of Theorem 9.5) intact:

1. A second-order integrator. ZOH is first-order accurate; Mamba-3 adopts a second-order exponential-trapezoidal scheme — the exponential-integrator idea of Chapter 4 §4.5 in its trapezoidal-quadrature form (§10.2) — while keeping ZOH’s exact treatment of the homogeneous dynamics (§10.2-10.3).
2. A complex state. A real diagonal mode $\statemat < 0$ can only decay; a complex mode $\statemat = \log\rho + i\theta$ decays and rotates, so a single head becomes a damped oscillator (§10.4).

Both changes are dynamical-systems moves, and both have been waiting since the foundations chapters: the integrator is the Chapter 4-5 story applied to a real architecture, and the complex state is the spiral phase portrait of Chapters 1-2. The C1 pilot’s empirical question — whether a second-order, structure-preserving integrator actually buys accuracy on a selective SSM without costing stability — is the through-line.

10.2 The exponential-trapezoidal scheme

Recall the continuous mode $\dot{\statevec} = \statemat \statevec + u(t)$, whose exact one-step map is $\statevec_k = e^{\statemat\stepsize}\statevec_{k-1} + \int_0^{\stepsize} e^{\statemat(\stepsize - s)}\, u(t_{k-1}+s)\,ds$. Every exponential integrator keeps the transition $e^{\statemat\stepsize}$ exact (the philosophy of Chapter 4 §4.5) and approximates only the forcing integral. ZOH holds the input constant at the left endpoint; the exponential-trapezoidal scheme applies the trapezoidal quadrature rule to the forcing integral — averaging the integrand at the two endpoints. The integrand is $e^{\statemat\stepsize} u_{k-1}$ at $s = 0$ and $u_k$ at $s = \stepsize$, so the rule gives $\tfrac{\stepsize}{2}(e^{\statemat\stepsize} u_{k-1} + u_k)$. Weighting the two endpoints by $\lambda \in [0,1]$ gives the three-term update

The contrast with ZOH is the third term: where ZOH reads only the left endpoint $u_{k-1}$, the trapezoid reads both endpoints. At $\lambda = \tfrac12$ this is the symmetric trapezoidal rule (globally second-order); at $\lambda = 1$ it collapses to a shifted ZOH (first order). Mamba-3 makes $\lambda$ itself input-dependent (a sigmoid of a learned projection), one more selective knob.

This single fact — that $\alpha$ is identical for the two schemes — has a sharp consequence for how the order can be measured.

The proof is Exercise 10.4 (§10.9). The practical reading is a trap worth naming: because the homogeneous transition is exact, a convergence study run on an autonomous system shows ZOH and exp-trapezoidal as identical lines at roundoff — the second order is simply not visible. The companion confirms this directly: on the autonomous mode both schemes sit at $\sim 2 \times 10^{-15}$, while on the forced mode the fitted convergence slopes separate cleanly into $1.00$ (ZOH) and $2.00$ (exp-trapezoidal). This is the same lesson Chapter 4 flagged in Exercise 4.3 (“setting $u \equiv 0$ makes ZOH exact; the slope is meaningless”), now load-bearing for verifying a real architecture’s headline claim.

10.3 Stability and accuracy decouple

Why adopt an exponential integrator rather than a higher-order Runge-Kutta method? Because the exponential transition is exact, the discrete amplification factor is $|\alpha| = e^{\operatorname{Re}(\statemat)\stepsize}$, which is $< 1$ for any stable mode ($\operatorname{Re}(\statemat) < 0$) at any step size. ZOH and exp-trapezoidal are therefore unconditionally stable on the homogeneous part — A-stable over the entire left half-plane — and, crucially, stability is set by $\statemat$ alone, independent of the order of the forcing quadrature. For a Runge-Kutta method a single stability polynomial controls both accuracy and the stability region; the exponential integrator decouples them.

The proof is Exercise 10.6 (§10.9). The companion makes the contrast quantitative: at a stiff $\statemat\stepsize = -30$, the exponential schemes give $|\alpha| \approx 9.4 \times 10^{-14}$ (fully damped) while bilinear gives $|\alpha| = 0.875$ (still ringing). This is why an exponential integrator is the right tool when the input-dependent $\statemat(u_t)$ of a selective SSM produces stiff dynamics, as Chapter 4 §4.5 anticipated for exactly this chapter. One caveat matters for §10.4: the theorem bounds the magnitude $|\alpha|$ only. For a complex mode it says nothing about how faithfully the discrete rotation tracks the continuous phase $\theta\stepsize$ — phase fidelity is a separate question, and the open one the frontier paragraph below raises.

The geometric-integrator frontier (where the pilot goes next). Exp-trapezoidal is a second-order exponential integrator; it is not symplectic. For the oscillatory complex modes introduced in §10.4, the question the symplectic-integrator pilot asks is whether a geometric integrator (Chapter 6) — one that preserves a phase-space invariant rather than merely a damping bound — controls long-horizon phase error better than exp-trapezoidal does. That is an open empirical question, not a settled result, and it is the heart of the “discretization stability atlas” the pilot pursues; Chapter 17 returns to it. This chapter establishes the settled part: exp-trapezoidal is second-order, A-stable, and exact on the homogeneous part.

10.4 Complex state: decay and oscillation in one head

A real diagonal mode $\statemat < 0$ gives $\alpha = e^{\statemat\stepsize} \in (0,1)$: pure geometric decay, a single forgetting timescale. Mamba-3 lets the mode be complex, $\statemat = \log\rho + i\theta$, so

The magnitude decays at rate $\log\rho$ per unit step and the phase advances by $\theta\stepsize$: a single head is now a damped oscillator, tracing the logarithmic spiral of Chapters 1-2 rather than a monotone decay. This lets one head represent an oscillatory feature — a frequency — instead of only a timescale.

In production Mamba-3 never stores a complex dtype. It keeps a real state and applies rotary position embeddings (RoPE) — 2-D rotations — to the projections $\inputmat, \outputmat$. The justification is the elementary isomorphism between complex multiplication and $2\times 2$ rotation-scaling.

The proof is Exercise 10.5 (§10.9). The companion pins the two recurrences equal to $\sim 3 \times 10^{-15}$ under a random complex drive, and confirms the spiral’s decay rate equals $\log\rho$ exactly — the empirical content of “decay and oscillation in one head.”

10.5 The trapezoidal SSD pass

The second-order stencil must coexist with the SSD machinery of Chapter 9 — the semiseparable matrix that lets the layer run as a parallel matmul or a linear-time scan. It does, because only the injection stencil changes. Unrolling the trapezoidal recurrence with the state-transition operator $\Phi(k,j) = \prod_{i=j+1}^{k} \alpha_i = \exp(\statemat \sum_{i=j+1}^{k}\stepsize_i)$ (the same $\Phi$ as Definition 9.2) splits the output into two semiseparable streams that share the one decay $\Phi$: a right-endpoint stream weighting source $j$ by $\lambda\stepsize_j$, and a left-endpoint stream weighting source $j$ by $(1-\lambda)\stepsize_{j+1}$ across the strictly-lower triangle. Both reuse the Chapter 9 segment-sum decay; the order-2 upgrade is two cheap masked matmuls in place of one.

The companion verifies that the dense two-stream matmul reproduces the sequential trapezoidal recurrence to $\sim 2 \times 10^{-14}$ (for both real and complex modes), and that the end-to-end output is genuinely second-order — halving the step size cuts the error by a factor of four. The SSD duality of §9.5 survives the integrator upgrade intact: same matrix, two schedules, now order 2.

10.6 MIMO and the production block

Two further Mamba-3 changes sit outside the integrator-and-state story this chapter develops, and we treat them briefly. MIMO rank-$R$ state mixing (the input and output projections each gain a rank-$R$ dimension) raises decode-time arithmetic intensity at fixed state size and, Mamba-3 reports, also lifts model quality through the added expressivity (a $\sim 1.2$-point average downstream gain over the SISO variant) — but neither effect changes the discretization or stability story. The full production block — SwiGLU gating, QK-normalization, the $\inputmat\outputmat$ bias — wraps the layer studied here; those too are orthogonal to the integrator-and-state question, and we point to the reference implementation rather than reproduce them Lahoti et al. (2026) .

What this chapter isolates is the dynamical-systems substance: a selective SSM (Chapter 9) run with a second-order exponential integrator (Chapters 4-5) and a complex state (Chapters 1-2, 8). The open program — whether geometric integrators (Chapter 6) preserve long-horizon phase structure better than the exponential scheme over the learned eigenvalue spectrum — is the discretization stability atlas the C1 pilot builds, and Chapter 17 returns to it.

10.7 What’s next

Mamba-3 kept Chapter 9’s selectivity and changed two things: a second-order exponential-trapezoidal integrator (accuracy in the forcing, stability from the exact exponential) and a complex state (decay and oscillation in one head), realized via RoPE and computed through the trapezoidal SSD pass. This closes the core SSM line — HiPPO (Chapter 7) to LTI S4/S4D/S5 (Chapter 8) to selective Mamba (Chapter 9) to Mamba-3 (here).

Chapter 11 opens the “beyond SSM” part by taking the SSD duality of §9.6 the other way: linear attention and the Hyena long-convolution lineage enter from the attention side of Theorem 9.5, where the structured contractions of these four chapters reappear wearing transformer clothing.

10.8 Exercises

Six problems: three short (inline solutions), three longer (full solutions in §10.9).

Exercise 10.1 (short)

For a complex mode $\statemat = \log\rho + i\theta$ with $\rho = 0.9$, $\theta = \pi/6$, and step $\stepsize = 1$, compute $|\alpha|$ and the phase advance of $\alpha = e^{\statemat\stepsize}$ over two steps. Which is the decay timescale and which is the oscillation frequency?

Exercise 10.2 (short, code)

Run companions/ch10/jax/discretization.py. Confirm the fitted convergence slopes are $\approx 1$ for ZOH and $\approx 2$ for exp-trapezoidal on the forced system. Then explain, in one sentence, why the same script reports both schemes at $\sim 10^{-15}$ on the homogeneous system.

Exercise 10.3 (short)

Write the exp-trapezoidal coefficients $(\alpha, \beta, \gamma)$ for $\lambda = 0$, $\lambda = \tfrac12$, and $\lambda = 1$. Show that $\lambda = 1$ reduces to a first-order shifted ZOH (no current-input term carrying the decay).

Exercise 10.4 (theory) — solution in §10.9

Prove Theorem 10.1: for $\lambda = \tfrac12$ the exponential-trapezoidal scheme has local truncation error $O(\stepsize^3)$ on a forced system, and show that on the homogeneous system both ZOH and exp-trapezoidal are exact.

Exercise 10.5 (theory) — solution in §10.9

Prove Theorem 10.3: that a 2-D rotation by $\theta$ on $[\operatorname{Re} z, \operatorname{Im} z]$ equals multiplication of $z$ by $e^{i\theta}$, and conclude the complex and real-RoPE recurrences coincide.

Exercise 10.6 (theory) — solution in §10.9

Prove Theorem 10.2: that the exponential schemes are A-stable with $|\alpha| \to 0$ on stiff modes, and contrast the bilinear Padé factor, showing $|\alpha| \to 1$ as $\statemat\stepsize \to -\infty$.

10.9 Full solutions to theory exercises

Solution to Exercise 10.4

Write the exact one-step map by variation of constants over $[t_{k-1}, t_k]$:

The transition $e^{\statemat\stepsize}$ is reproduced exactly by both schemes, so all error is in the forcing integral. Exp-trapezoidal approximates $\int_0^{\stepsize} g(s)\,ds$, with integrand $g(s) = e^{\statemat(\stepsize-s)}\,u(t_{k-1}+s)$, by the trapezoidal quadrature rule $\tfrac{\stepsize}{2}\big(g(0) + g(\stepsize)\big)$ — the $\lambda = \tfrac12$ endpoint weights of §10.2. The trapezoidal rule on an interval of width $\stepsize$ has error $-\tfrac{\stepsize^3}{12}\,g''(\xi)$ for some interior $\xi$; since $g$ is $C^2$ whenever $u$ is $C^2$, the local truncation error is $O(\stepsize^3)$ per step, and summing $O(\stepsize^{-1})$ steps over a fixed horizon yields global order 2. For the homogeneous system $u \equiv 0$ the forcing integral vanishes identically and $\statevec_k = e^{\statemat t_k}\statevec_0$ exactly for both schemes, regardless of order — there is no forcing term to approximate, so the order distinction disappears. $\blacksquare$

Solution to Exercise 10.5

Let $z = a + bi$ and write $[z] = [a, b]^\top$. Then

whose real-vector form is $[a\cos\theta - b\sin\theta,\; a\sin\theta + b\cos\theta]^\top = R(\theta)[z]$ with $R(\theta) = \left(\begin{smallmatrix}\cos\theta & -\sin\theta\\ \sin\theta & \cos\theta\end{smallmatrix}\right)$. Scaling by $\rho$ commutes with the real isomorphism, so multiplication by $\rho e^{i\theta}$ is $\rho R(\theta)$ on $[z]$. Applying this at every step of the affine recurrence (the additive drive $b_k$ maps to $[b_k]$ componentwise, which the isomorphism also preserves) shows the complex recurrence $x_k = (\rho e^{i\theta})x_{k-1} + b_k$ and the real recurrence $s_k = \rho R(\theta) s_{k-1} + [b_k]$ produce identical trajectories under $[\cdot]$. $\blacksquare$

Solution to Exercise 10.6

For the exponential schemes $\alpha = e^{\statemat\stepsize}$, so $|\alpha| = e^{\operatorname{Re}(\statemat)\stepsize}$. If $\operatorname{Re}(\statemat) < 0$ and $\stepsize > 0$ then the exponent is negative and $|\alpha| < 1$ for every $\stepsize$ (A-stability over the whole left half-plane), and as the mode stiffens ($\operatorname{Re}(\statemat)\stepsize \to -\infty$) we have $|\alpha| \to 0$. For the bilinear scheme $\alpha = (1 + z/2)/(1 - z/2)$ with $z = \statemat\stepsize$; writing $z = x + iy$ with $x < 0$, $|1 + z/2| < |1 - z/2|$ (the left point is closer to $-2$, the right to $+2$), so $|\alpha| < 1$ — also A-stable. But as $x \to -\infty$ the ratio $(1 + z/2)/(1 - z/2) \to -1$, so $|\alpha| \to 1$: stiff modes are mapped to the unit circle and ring rather than decaying. This is the qualitative gap the exponential transition closes. $\blacksquare$

10.10 Companion code

The companions live in companions/ch10/{jax,julia,torch}/. The JAX module covers the integrator, the complex state, and the trapezoidal SSD pass; the PyTorch module mirrors the integrator and complex state (cross-framework parity); the Julia module is the discretization atlas only — the numerical-analysis core — and joins the Chapter 4-6 Julia family that the C1 pilot’s stability atlas builds on.

The discretization study (order slopes, homogeneous-blindness, stability regions):

PYTHONPATH=. python companions/ch10/jax/discretization.py
julia --project=companions/ch10/julia companions/ch10/julia/discretization.jl

The complex state and RoPE equivalence (the spiral figure):

PYTHONPATH=. python companions/ch10/jax/complex_state.py

The trapezoidal SSD pass (matmul equals the sequential oracle):

PYTHONPATH=. python companions/ch10/jax/trapezoidal_ssd.py

The PyTorch mirror (integrator and complex state, parity-checked against JAX):

PYTHONPATH=. python companions/ch10/torch/discretization.py
PYTHONPATH=. python companions/ch10/torch/complex_state.py

Run the suites — JAX and Julia in the fast loop, torch separately:

make companion-jax-tests
make companion-julia-tests
make companion-torch-tests