Magnus Expansion Numerics

Magnus integrators approximate the time-ordered propagator by the exponential of a matrix built from Hamiltonians and their commutators. For Hermitian Hamiltonians, the truncated exponent remains anti-Hermitian, so its exact exponential is unitary. This structural property coexists with time-discretization error and with whatever error is introduced when applying the exponential.

Helpful background. Matrix exponentials supplies exponential evaluation; the time-evolution operator explains time ordering.

Commutators encode the order of the drive

Write $\dot U=A(t)U$, with $A(t)=-iH(t)/\hbar$. On one interval, write $U=e^\Omega$ and expand

$$\begin{aligned} \Omega_1&=\int_{t_0}^{t_0+h}A(t_1)\,dt_1,\\ \Omega_2&=\frac12\int_{t_0}^{t_0+h}dt_1 \int_{t_0}^{t_1}dt_2\,[A(t_1),A(t_2)]. \end{aligned}$$

If Hamiltonians commute at all pairs of times, all commutator terms vanish and the ordinary exponential of the integral is exact. Otherwise $\Omega_2$ accounts for the noncommuting order of the drive.

Commutators of anti-Hermitian matrices are anti-Hermitian. Hence truncating the Magnus series preserves that algebraic structure. A common sufficient convergence condition for finite matrices is $\int\|A(t)\|_2dt<\pi$ on the interval; it is not necessary and does not by itself set an accuracy tolerance. Long evolution is handled through short steps, not by assuming one truncated exponent works globally.

Practical second- and fourth-order steps

The exponential midpoint rule uses

$$\Omega^{[2]}=hA(t_0+h/2).$$

For sufficiently smooth $A$, it has local error $O(h^3)$ and global order two. To obtain order four, evaluate at two Gauss nodes,

$$A_{1,2}=A\!\left(t_0+h\left[\frac12\mp\frac{\sqrt3}{6}\right]\right),$$

and use

$$\Omega^{[4]} =\frac h2(A_1+A_2) -\frac{\sqrt3h^2}{12}[A_1,A_2].$$

The sign follows the evolution convention $\dot U=A U$ and the chronological ordering of the nodes. Reversing a commutator without reversing its coefficient changes the method. Apply $e^{\Omega^{[4]}}$ with a controlled dense or Krylov exponential action.

Explicit commutators can destroy sparsity. Their action can sometimes be evaluated as $A_1(A_2v)-A_2(A_1v)$ without forming a matrix, although this increases each action’s cost. Commutator-free exponential compositions are another option; their coefficients and order conditions should come from a verified scheme.

A linear-generator sign check

On $[0,h]$, let $A(t)=A_0+tB$. Then

$$[A(t_1),A(t_2)]=(t_2-t_1)[A_0,B],$$

so

$$\Omega_1=hA_0+\frac{h^2}{2}B,\qquad \Omega_2=-\frac{h^3}{12}[A_0,B].$$

The two-node formula reproduces both terms because $[A_1,A_2]=(h/\sqrt3)[A_0,B]$. This symbolic check catches a frequent sign error before any large propagation. The midpoint exponential includes $\Omega_1$ exactly for this example but omits the leading ordering correction.

A constant Hamiltonian gives $A_1=A_2$, so the formula reduces to the exact constant-generator exponential. This tests the limiting case but does not measure fourth-order convergence.

Refinement and nonsmooth controls

For a smooth noncommuting drive, compare final states at $h,h/2,h/4$ and verify the expected factor of about $16$ in differences for the fourth-order method. Keep exponential-action error below this scale. Check observables and phase-sensitive distances as well as unitarity.

Split integration at pulse jumps or derivative discontinuities. High-order quadrature across an unresolved discontinuity loses its smooth-problem order. A narrow pulse can require small steps even when its integrated area is moderate. If the Hamiltonian depends on the evolving state, a linear Magnus step with frozen state dependence needs an additional nonlinear treatment.

Exercises

If $A(t)=f(t)A_0$ for a fixed matrix $A_0$, what error remains after all commutator terms vanish?

Solution

The exact exponent is $A_0\int f(t)dt$. A numerical method can still approximate that scalar integral inaccurately, and the exponential action can have its own error. Commutativity removes time-ordering error, not quadrature or matrix-function error.

References

• Blanes, S., F. Casas, J. A. Oteo, and J. Ros, “The Magnus expansion and some of its applications,” Physics Reports 470, 151–238 (2009), doi:10.1016/j.physrep.2008.11.001 — convergence, quadrature formulas, and numerical integrators.
• Saad, Y., “Analysis of some Krylov subspace approximations to the matrix exponential operator,” SIAM Journal on Numerical Analysis 29, 209–228 (1992), doi:10.1137/0729014 — controlled exponential actions.