Higher-Order Trotter

The first-order Trotter formula accumulates a global error of $O(t^2/r)$ after $r$ steps. For many practical simulations that rate of convergence is too slow — the circuit depth required to reach a target accuracy $\varepsilon$ grows faster than necessary. Higher-order Trotter formulas, also called Suzuki–Trotter decompositions, arrange the same exponentials in a cleverer order so that leading error terms cancel, reducing the global error to $O(t^{2k+1}/r^{2k})$ for an order- $2k$ formula.

Recap: first-order Trotter

For a Hamiltonian split as $H = A + B$ , the Lie–Trotter product formula gives

e^{-i(A+B)\tau} = e^{-iA\tau}\,e^{-iB\tau} + O(\tau^2),

where $\tau = t/r$ is one step size. The local error per step is $O(\tau^2)$ ; after $r$ steps the global error is $r \cdot O(\tau^2) = O(t^2/r)$ .

The symmetric (second-order) formula

Suzuki observed in 1976 that interleaving the exponentials symmetrically causes the leading $O(\tau^2)$ term to cancel exactly. The second-order Suzuki–Trotter formula is

S_2(\tau) = e^{-iA\tau/2}\,e^{-iB\tau}\,e^{-iA\tau/2},

so that

e^{-i(A+B)\tau} = S_2(\tau) + O(\tau^3).

The local error is now $O(\tau^3)$ , and after $r$ steps the global error is $r \cdot O(\tau^3) = O(t^3/r^2)$ . For a fixed total time $t$ and a desired accuracy $\varepsilon$ , the required number of steps scales as $r = O(t^{3/2}/\varepsilon^{1/2})$ for the second-order formula, compared with $r = O(t^2/\varepsilon)$ for first order — a genuine improvement in circuit depth.

Why does symmetry kill the leading error?

The Baker–Campbell–Hausdorff (BCH) expansion of $e^X e^Y$ generates a series of nested commutators:

e^X e^Y = e^{X+Y + \frac{1}{2}[X,Y] + \frac{1}{12}([X,[X,Y]] - [Y,[Y,X]]) + \cdots}.

For $X = -iA\tau$ and $Y = -iB\tau$ the leading error term is $\frac{1}{2}[X,Y] = -\frac{\tau^2}{2}[A,B]$ . Now consider the symmetric product:

e^{-iA\tau/2}\,e^{-iB\tau}\,e^{-iA\tau/2}.

Applying BCH twice, the $[A,B]$ commutator appears once with coefficient $+\tfrac{1}{2}$ (from the left half) and once with coefficient $-\tfrac{1}{2}$ (from the right half), so the two contributions cancel. The next nonzero term is $O(\tau^3)$ , coming from double-commutators $[A,[A,B]]$ and $[B,[B,A]]$ .

Multi-term Hamiltonians

For $H = H_1 + H_2 + \cdots + H_k$ with $k$ terms, the second-order formula generalises to a symmetric sweep:

S_2(\tau) = e^{-iH_1\tau/2}\,e^{-iH_2\tau/2}\cdots e^{-iH_k\tau/2}\cdot e^{-iH_k\tau/2}\cdots e^{-iH_2\tau/2}\,e^{-iH_1\tau/2}.

This "forward then reverse" structure is sometimes called the Strang splitting in numerical analysis. Each step uses $2k - 1$ exponentials instead of $k$ , but the halved exponents at each boundary allow adjacent steps to be merged: the $e^{-iH_j\tau/2}$ ending step $r$ and the $e^{-iH_j\tau/2}$ starting step $r+1$ combine into a single $e^{-iH_j\tau}$ , keeping the total cost at roughly $k$ exponentials per effective step.

Fourth-order and fractal constructions

Suzuki showed that a $2k$ -th-order formula can be built recursively from any $(2k-2)$ -th-order formula $S_{2k-2}$ . The fourth-order formula is

S_4(\tau) = S_2(p\tau)^2\,S_2((1-4p)\tau)\,S_2(p\tau)^2,

where $p = (4 - 4^{1/3})^{-1}$ . This recursion yields fractal-like circuits that look complicated but give error $O(\tau^5)$ per step, so $r = O(t^{5/4}/\varepsilon^{1/4})$ steps suffice.

The table below summarises the scaling for order $2k$ :

| Order | Local error | Global error | Steps for accuracy $\varepsilon$ | |-------|-------------|--------------|----------------------------------| | 1 | $O(\tau^2)$ | $O(t^2/r)$ | $O(t^2/\varepsilon)$ | | 2 | $O(\tau^3)$ | $O(t^3/r^2)$ | $O(t^{3/2}/\varepsilon^{1/2})$ | | 4 | $O(\tau^5)$ | $O(t^5/r^4)$ | $O(t^{5/4}/\varepsilon^{1/4})$ | | $2k$ | $O(\tau^{2k+1})$ | $O(t^{2k+1}/r^{2k})$ | $O(t^{(2k+1)/(2k)}/\varepsilon^{1/(2k)})$ |

All scalings are in the limit of large $r$ with fixed total time $t$ .

Practical considerations

Gate count vs. error. Higher-order formulas use more gates per step but need fewer steps. For small $\varepsilon$ and large $t$ , order 4 or 6 is usually cheaper in total gate count than order 2. For large $\varepsilon$ (noisy-device regime) the overhead of extra exponentials may outweigh the benefit, making order 2 preferable.

Negative time steps. Fourth-order and higher formulas require some steps with negative $p$ coefficients, i.e., evolution backward in time. This is perfectly valid unitarily but can cause practical issues on noisy hardware because the backward gates may not cancel noise accumulated in the forward gates.

Randomised Trotter. Recent work (e.g., the qDRIFT algorithm by Campbell 2019) achieves a favourable system-size scaling by randomly sampling individual Hamiltonian terms with probability proportional to their weight. The gate count depends on $\sum_j \|H_j\|$ rather than on the spectral norm of the full Hamiltonian, which can be a large improvement for Hamiltonians with many weak terms — but the error per step still scales as $O(\tau^2)$ , the same exponent as first-order Trotter.

Summary

The symmetric Suzuki–Trotter formula $S_2(\tau) = e^{-iA\tau/2}\,e^{-iB\tau}\,e^{-iA\tau/2}$ achieves a global error of $O(t^3/r^2)$ by exploiting the symmetry of the BCH expansion to cancel the leading commutator term. Higher-order formulas are built recursively and achieve error $O(t^{2k+1}/r^{2k})$ , dramatically reducing the circuit depth needed for a fixed accuracy. The second-order formula is the workhorse of near-term Hamiltonian simulation and is the natural starting point whenever first-order Trotter proves insufficiently accurate.