Matrix Exponentials

For a time-independent finite Hamiltonian, the propagator is $U(t)=e^{-iHt/\hbar}$. The matrix exponential is defined by a convergent power series, but evaluating that series naively is rarely the most reliable algorithm. The computational choice is also different when only $U(t)v$ is needed rather than the complete dense matrix.

Helpful background. Matrix diagonalization supplies Hermitian eigensystems; the numerical Schrödinger equation supplies trajectory checks.

Spectral evaluation for Hermitian matrices

If $H=Q\Lambda Q^\dagger$ with orthonormal $Q$, then

$$U(t)=Q\,\operatorname{diag}(e^{-iE_jt/\hbar})\,Q^\dagger.$$

For one state, first compute $a=Q^\dagger v$, multiply its components by the phases, and transform back. Once diagonalization is complete, many times or initial vectors can reuse it. Degenerate eigenvectors may rotate, but the propagator remains invariant because equal eigenvalues receive equal phases.

The decomposition costs dense cubic work and quadratic storage. It is valuable as a reference for moderate matrices, but storing $U$ explicitly adds no benefit when an action on a few vectors suffices. For a non-Hermitian or highly non-normal generator, an eigenvector decomposition can be ill-conditioned or fail at a defective matrix; use an appropriate general matrix-function method.

Scaling, rational approximation, and actions

Scaling and squaring uses

$$e^A=\left(e^{A/2^s}\right)^{2^s}$$

and approximates the small-norm exponential, commonly by a Padé rational function. Linear solves evaluate the rational function; repeated squaring restores the scale. The scaling and approximation degree must balance truncation and roundoff. Large unnecessary scaling can itself worsen accuracy, so use a tested implementation instead of a fixed arbitrary number of squarings.

For large sparse problems, exponential-action algorithms approximate $e^Av$ without forming $e^A$, which is usually dense. Krylov propagation and scaled truncated-Taylor actions are two approaches. A sparse input does not imply a sparse exponential.

For $A=-iHt/\hbar$, exact $U$ is unitary. Numerical algorithms approximate it to a tolerance; check $\|U^\dagger U-I\|$ when the full matrix is available and norm changes for representative vectors otherwise. Neither check alone establishes phase accuracy.

Shifts and a scalar invariance test

A scalar shift gives

$$e^{-i(H-cI)t/\hbar} =e^{ict/\hbar}e^{-iHt/\hbar}.$$

Subtracting a convenient energy center can reduce the spectral scale of the generator. For one isolated state trajectory, this changes only a global phase. Preserve that phase when comparing amplitudes against a reference with a fixed energy origin.

For $H=\hbar\Omega\sigma_x/2$, $H^2=(\hbar\Omega/2)^2I$, so the power series separates into even and odd terms:

$$U(t)=\cos(\Omega t/2)I-i\sin(\Omega t/2)\sigma_x.$$

This exact expression checks the full matrix, not merely its norm. Also compare the semigroup identity $U(t+s)=U(t)U(s)$ for a fixed $H$. Floating-point discrepancies should decrease with numerical tolerance until roundoff dominates.

Time dependence and conditioning

For $H(t)$, exponentiating $-i\int H(t)dt/\hbar$ is generally wrong because Hamiltonians at distinct times need not commute. Products of short-time exponentials approximate time ordering; Magnus methods incorporate commutator corrections systematically.

A long propagation time magnifies an energy error into a phase error: $\delta\phi_j=t\,\delta E_j/\hbar$. Therefore highly accurate matrix exponentiation cannot compensate for inaccurate Hamiltonian eigenvalues from a coarse spatial basis. Bound the time horizon for the desired phase accuracy and refine the representation accordingly.

Exercises

Suppose two implementations use Hamiltonians differing by $cI$. Which comparisons should agree, and which need a phase adjustment?

Solution

Probabilities and expectation values of the same operators agree for the same initial state. State vectors and propagator matrices differ by the global phase $e^{-ict/\hbar}$. Their raw vector or matrix difference is therefore not a valid error measure until the energy-origin phase is accounted for.

References

• Al-Mohy, A. H., and N. J. Higham, “Computing the action of the matrix exponential, with an application to exponential integrators,” SIAM Journal on Scientific Computing 33, 488–511 (2011), author-hosted record — action algorithms and backward-error control.
• Higham, N. J., Functions of Matrices: Theory and Computation, SIAM, 2008 — scaling and squaring, conditioning, and general matrix functions.