Time-Dependent Hamiltonian Notebook

This notebook guide makes time ordering numerically unavoidable. A rotating-field two-level Hamiltonian has nonzero commutators at unequal times, yet it also admits an exact rotating-frame solution. It is therefore an unusually clean benchmark for showing why

$$\exp\left[ - \frac{i}{\hbar} \int_0^T H(t)\,dt \right]$$

is generally not the evolution operator.

The formal distinction belongs to Time-Dependent Hamiltonians and Time Ordering. This page owns the reproducible numerical comparison among an exact solution, a deliberately naive exponential, chronologically ordered short-step products, and an adaptive ODE solver.

The existing source notebook notebooks/wave-mechanics-canonical-systems/two-level-systems/two-level-system-dynamics.ipynb supplies compatible Pauli-matrix and state-validation patterns. The present calculation adds explicit drive time dependence and unequal-time noncommutativity.

Purpose

The notebook should demonstrate that:

• explicit time dependence alone does not invalidate an ordinary exponential; unequal-time noncommutativity does;
• the naive exponential can converge accurately to the wrong operator;
• chronological products of exact short-step exponentials converge to the time-ordered evolution;
• left-endpoint and midpoint products have different global convergence orders;
• an adaptive ODE solver provides an independent numerical route but is not automatically unitary;
• norm conservation is necessary but insufficient for validating a quantum propagator.

The comparison should be operator-level whenever possible. Testing one initial state can miss an error on the orthogonal state or conceal a global phase.

Rotating-Field Model

Use

$$H(t) = \frac{\hbar}{2} \mathbf b(t)\cdot\boldsymbol\sigma,$$

with

$$\mathbf b(t) = \left( \Omega\cos\nu t,\, \Omega\sin\nu t,\, \Delta \right).$$

Here $\Omega$ is the transverse coupling, $\Delta$ is the static longitudinal angular frequency, and $\nu$ is the angular frequency of the rotating transverse field. All three have units of inverse time.

The Hamiltonian is Hermitian and traceless. Its instantaneous eigenvalues are constant,

$$E_\pm = \pm \frac{\hbar}{2} \sqrt{\Omega^2+\Delta^2},$$

even though its eigenvectors rotate in time. Constant instantaneous eigenvalues do not make the dynamics equivalent to a time-independent Hamiltonian.

Unequal-Time Commutator

For Pauli vectors,

$$\left[ \mathbf a\cdot\boldsymbol\sigma,\, \mathbf c\cdot\boldsymbol\sigma \right] = 2i \left( \mathbf a\times\mathbf c \right) \cdot\boldsymbol\sigma.$$

Therefore

$$\begin{aligned} [H(t_1),H(t_2)] &= \frac{i\hbar^2}{2} \left[ \mathbf b(t_1)\times \mathbf b(t_2) \right] \cdot\boldsymbol\sigma. \end{aligned}$$

This is generically nonzero. For example, at

$$t_1=0, \qquad t_2=\frac{\pi}{2\nu},$$

the cross product is

$$\mathbf b(t_1)\times\mathbf b(t_2) = \left( -\Delta\Omega,\, -\Delta\Omega,\, \Omega^2 \right).$$

The notebook should evaluate the commutator norm numerically before attempting time evolution. This prevents the key assumption from remaining implicit.

Exact Rotating-Frame Solution

Define the rotation

$$R(t) = \exp\left( - \frac{i\nu t}{2}\sigma_z \right).$$

It obeys

$$R(t)\sigma_xR^\dagger(t) = \cos\nu t\,\sigma_x + \sin\nu t\,\sigma_y.$$

The laboratory Hamiltonian can be written

$$H(t) = R(t)H_0R^\dagger(t),$$

where

$$H_0 = \frac{\hbar}{2} \left( \Omega\sigma_x+\Delta\sigma_z \right).$$

Set

$$\lvert\psi(t)\rangle = R(t)\lvert\chi(t)\rangle.$$

The rotating-frame state satisfies

$$i\hbar \frac{d}{dt} \lvert\chi(t)\rangle = H_{\rm rot} \lvert\chi(t)\rangle,$$

with the constant Hamiltonian

$$H_{\rm rot} = \frac{\hbar}{2} \left[ \Omega\sigma_x + (\Delta-\nu)\sigma_z \right].$$

Thus the exact propagator from $0$ to $T$ is

$$U_{\rm ex}(T,0) = R(T) \exp\left( - \frac{i}{\hbar} H_{\rm rot}T \right).$$

This formula is the primary benchmark. Verify it by checking $U_{\rm ex}(0,0)=I$ and the operator differential equation

$$i\hbar \frac{\partial U_{\rm ex}}{\partial T} = H(T)U_{\rm ex}.$$

Closed Pauli Exponentials

For any real vector $\mathbf a$,

$$\begin{aligned} \exp\left( - \frac{i}{2} \mathbf a\cdot\boldsymbol\sigma \right) &= \cos\left( \frac{\lvert\mathbf a\rvert}{2} \right)I \\ &\quad- i\sin\left( \frac{\lvert\mathbf a\rvert}{2} \right) \frac{ \mathbf a\cdot\boldsymbol\sigma }{ \lvert\mathbf a\rvert }. \end{aligned}$$

Use the continuous limit $I$ when $\lvert\mathbf a\rvert=0$. This identity evaluates every exact and short-step exponential without depending on a general matrix-exponential routine. A library routine remains useful as an independent cross-check.

The Naive Exponential

The time-integrated effective vector is

$$\begin{aligned} \mathbf B(T) &= \int_0^Tdt\, \mathbf b(t) \\ &= \left( \frac{\Omega}{\nu}\sin\nu T,\, \frac{\Omega}{\nu} \left[1-\cos\nu T\right],\, \Delta T \right). \end{aligned}$$

The naive approximation is

$$U_{\rm naive}(T) = \exp\left[ - \frac{i}{2} \mathbf B(T)\cdot\boldsymbol\sigma \right].$$

This expression uses the exact time integral. Its error is not quadrature error. It discards all information about the order in which differently oriented infinitesimal rotations act.

At one drive period,

$$T_d=\frac{2\pi}{\nu},$$

the transverse components of $\mathbf B$ vanish:

$$\mathbf B(T_d) = \left( 0,0,\Delta T_d \right).$$

On resonance, $\Delta=\nu$, the naive result is

$$U_{\rm naive}(T_d)=-I.$$

It predicts no transition between $\sigma_z$ basis states.

The exact result on resonance is

$$U_{\rm ex}(T_d,0) = - \exp\left( - \frac{i\Omega T_d}{2}\sigma_x \right).$$

For an initial $\lvert0\rangle$ with $\sigma_z\lvert0\rangle=\lvert0\rangle$,

$$P_{0\to1}^{\rm ex}(T_d) = \sin^2\left( \frac{\pi\Omega}{\nu} \right),$$

while

$$P_{0\to1}^{\rm naive}(T_d)=0.$$

This is a direct physical consequence of time ordering, not a small numerical correction.

Chronologically Ordered Products

Divide $[0,T]$ into $N$ equal steps,

$$\Delta t=\frac{T}{N}.$$

For sampling times $t_j^*$, define the exact frozen-Hamiltonian step

$$S_j = \exp\left[ - \frac{i}{\hbar} H(t_j^*)\Delta t \right].$$

The chronological product is

$$U_N(T,0) = S_{N-1}S_{N-2}\cdots S_1S_0.$$

The earliest step acts first and appears on the right. In an implementation that loops forward in $j$, update by left multiplication:

$$U\leftarrow S_jU.$$

Two baseline choices are:

$$t_j^*=j\Delta t$$

for the left-endpoint product, and

$$t_j^* = \left(j+\frac12\right)\Delta t$$

for the exponential midpoint product.

For smooth $H(t)$, the expected global operator errors are

$$E_{\rm left} = O(\Delta t),$$

and

$$E_{\rm mid} = O(\Delta t^2).$$

Each individual step and the full product are exactly unitary in exact arithmetic. Unitarity does not prove accuracy: reversing the product order also multiplies unitary matrices but approaches the wrong ordered evolution.

A Discretized Naive Control

To separate quadrature convergence from time-ordering convergence, also compute

$$\begin{aligned} U_{{\rm sum},N} &= \exp\left[ - \frac{i\Delta t}{\hbar} \sum_{j=0}^{N-1} H(t_j^*) \right]. \end{aligned}$$

As $N$ grows, the sum converges to $\int H(t)\,dt$. For a noncommuting model, however, $U_{{\rm sum},N}$ approaches $U_{\rm naive}$ rather than $U_{\rm ex}$. Its error therefore plateaus at a nonzero value even while the integral quadrature converges.

This is an important diagnostic pattern: numerical convergence of an intermediate quantity does not establish that the mathematical formula built from it is correct.

Adaptive ODE Route

As an independent route, integrate

$$\frac{dU}{dt} = - \frac{i}{\hbar} H(t)U(t), \qquad U(0)=I,$$

with a high-order adaptive ODE solver. Flatten the complex $2\times2$ matrix into a length-four vector and reshape it inside the derivative function. Use a complex initial array so the solver stays in the complex domain.

A high-precision explicit method such as DOP853 is appropriate for this small smooth benchmark. Record relative tolerance, absolute tolerance, accepted step count, function evaluations, and any maximum-step constraint. Tighten tolerances and reduce the maximum step until the operator error stabilizes.

Adaptive local-error control does not enforce exact unitarity. Report the unitarity residual rather than renormalizing columns after each step. The exact rotating-frame result remains the standard against which the adaptive solver is judged.

Commuting Control Case

Run the same workflow for

$$H_c(t) = \frac{\hbar}{2} \left( \Delta+a\cos\nu t \right)\sigma_z.$$

All unequal-time Hamiltonians commute:

$$[H_c(t_1),H_c(t_2)]=0.$$

The exact propagator is the ordinary exponential

$$\begin{aligned} U_c(T,0) &= \exp\left\{ - \frac{i\sigma_z}{2} \left[ \Delta T + \frac{a}{\nu}\sin\nu T \right] \right\}. \end{aligned}$$

The naive, ordered-product, and adaptive routes should all converge to this answer. This control establishes that the failure in the rotating-field model comes from noncommutativity rather than from explicit time dependence by itself.

Baseline Parameters

Set

$$\hbar=1, \qquad \Omega=0.7, \qquad \Delta=\nu=1.3,$$

and evolve for one drive period,

$$T=T_d=\frac{2\pi}{\nu}.$$

Use the initial state

$$\lvert\psi_0\rangle = \begin{pmatrix} 1\\ 0 \end{pmatrix}.$$

The exact transition probability is

$$P_{0\to1}^{\rm ex} = \sin^2\left( \frac{0.7\pi}{1.3} \right) \approx 0.98547.$$

The naive exponential predicts zero. For the convergence study, begin with

$$N=8,16,32,64,128,256.$$

With double precision and the closed Pauli step, these values should expose the first- and second-order regimes before roundoff becomes relevant.

Notebook Workflow

Organize the notebook into reproducible stages:

1. define Pauli matrices and verify their commutators;
2. define $\mathbf b(t)$ and $H(t)$;
3. check Hermiticity and a nonzero unequal-time commutator;
4. construct $U_{\rm ex}$ from the rotating-frame formula;
5. verify its initial condition and differential equation;
6. compute $U_{\rm naive}$ from the analytic integral;
7. implement left-endpoint and midpoint chronological products;
8. implement the exponentiated-sum control and the deliberately reversed product;
9. integrate the operator ODE adaptively;
10. compare operator, state, transition-probability, and Bloch-vector results;
11. estimate observed convergence orders;
12. repeat the calculation for the commuting control Hamiltonian;
13. run assertions before making plots.

Record parameter units, transform conventions, multiplication order, solver tolerances, step counts, and error norms.

Operator and State Errors

Use the phase-sensitive Frobenius error

$$E_F = \left\lVert U_{\rm num}-U_{\rm ex} \right\rVert_F.$$

Also report the global-phase-insensitive gate error

$$E_{\rm gate} = 1- \frac{ \left\lvert \operatorname{Tr} \left( U_{\rm ex}^\dagger U_{\rm num} \right) \right\rvert^2 }{4}.$$

For a chosen state, use

$$E_{\rm state} = 1- \left\lvert \langle\psi_{\rm ex} \vert \psi_{\rm num}\rangle \right\rvert^2.$$

The unitarity residual is

$$E_U = \left\lVert U_{\rm num}^\dagger U_{\rm num}-I \right\rVert_F.$$

A wrong ordered product can have $E_U$ near roundoff while $E_F$ is large. Conversely, an adaptive Runge–Kutta solution can have a small physical error and a small but nonzero unitarity residual.

Observed Convergence Order

For an error sequence $E_N$, estimate

$$p_N = \frac{ \log(E_N/E_{2N}) }{ \log2 }.$$

In the asymptotic regime,

$$p_N\to1$$

for left-endpoint ordering and

$$p_N\to2$$

for midpoint ordering.

Plot error against $\Delta t$ on logarithmic axes and report the numerical slopes. A slope inferred from only two coarse points is not convincing; include enough refinements to identify a stable regime.

Bloch and Energy Diagnostics

For a pure state, compute

$$r_j(t) = \langle\psi(t)\vert \sigma_j \lvert\psi(t)\rangle.$$

The Bloch-vector norm should remain

$$\lvert\mathbf r(t)\rvert=1.$$

Plot the exact and numerical trajectories on equal scales or compare their Cartesian components versus time. The rotating laboratory field and rotating-frame effective field should not be conflated.

For the closed driven system,

$$\frac{d}{dt} \langle H(t)\rangle = \left\langle \frac{\partial H(t)}{\partial t} \right\rangle.$$

An optional finite-difference check of this identity tests both state evolution and explicit Hamiltonian time dependence. The interpretation of energy exchange with the external drive continues in Driven Closed Quantum Systems.

Validation Checks

The minimum validation suite is:

Check	Target
Pauli algebra	$[\sigma_i,\sigma_j]=2i\epsilon_{ijk}\sigma_k$
Hamiltonian structure	Hermitian and traceless at sampled times
Noncommutativity	selected unequal-time commutator norm is nonzero
Exact benchmark	initial condition and operator Schrödinger residual vanish
Ordered products	left and midpoint routes converge to $U_{\rm ex}$
Observed order	slopes approach one and two, respectively
Exponentiated sum	integral converges but operator error plateaus
Reversed product	remains unitary but does not converge to forward evolution
Adaptive ODE	error decreases under tolerance and maximum-step refinement
Transition probability	approaches $0.98547$ for the baseline
Unitarity	ordered exponential products remain unitary to roundoff
Composition	$U(T,0)=U(T,t_*)U(t_*,0)$ with absolute-time sampling
Commuting control	naive and ordered routes converge to the same answer
Bloch norm	remains one for exact pure-state evolution

For the composition test, build both subinterval propagators using the Hamiltonian at their actual laboratory times. Restarting the drive phase at $t_*$ changes the problem.

Expected Results

For the baseline model:

• $U_{\rm naive}$ predicts zero transition after one period;
• the exact transition probability is approximately $0.98547$;
• the exponentiated-sum result approaches the naive answer as quadrature is refined;
• the left chronological product shows first-order global convergence;
• the midpoint chronological product shows second-order global convergence;
• both ordered products remain unitary to roundoff;
• the adaptive ODE result approaches the exact operator as tolerances tighten;
• a reversed product can be exactly unitary and physically wrong;
• all correct routes agree for the commuting control Hamiltonian.

These outcomes distinguish algebraic correctness, numerical accuracy, and structure preservation.

Common Numerical Failures

• Exponentiating the integrated Hamiltonian without checking commutators. Accurate quadrature cannot restore missing time ordering.
• Multiplying step operators on the wrong side. The earliest factor belongs on the right for column-state evolution.
• Calling each frozen step exact and therefore calling the full product exact. Time sampling still creates global error.
• Using only norm conservation. Any unitary matrix preserves norm, including the wrong reversed product.
• Renormalizing an adaptive solution after every output point. This hides solver drift and changes the numerical method.
• Testing one initial state only. Operator errors can act entirely in another state direction.
• Ignoring global phase in every metric. Transition probabilities may be insensitive to a phase that matters in composition or controlled evolution.
• Forgetting the spinor sign $R(T_d)=-I$. A $2\pi$ rotation is not the identity on spinors.
• Restarting the drive phase on each subinterval. Nonautonomous propagators depend on both endpoint times.
• Tightening solver tolerances without limiting step size. A rapidly varying drive may still be undersampled.
• Comparing methods at different Hamiltonian sample times. This confounds method order with a convention change.
• Treating energy nonconservation as numerical failure. The external drive can do work on the closed system.

Extensions

• Move off resonance and compare the exact rotating-frame transition formula.
• Replace the circularly rotating field with an elliptically polarized drive, removing the simple exact benchmark.
• Compare the midpoint product with a commutator-corrected Magnus approximation on short intervals.
• Add piecewise control pulses and verify chronological factor ordering.
• Continue to Floquet Operators by diagonalizing the one-period propagator.
• Compare the time-dependent benchmark with the operator-splitting errors in the Trotter Evolution Notebook.

Cross-Links

• Time-Dependent Hamiltonians
• Time Ordering
• Dyson Expansion as Formal Evolution
• Rotating Frames
• Driven Closed Quantum Systems
• Two-State Hamiltonians
• Time-Stepping Methods
• Matrix Exponentials Numerically
• Convergence Tests
• Validation Tests

References

• J. J. Sakurai and J. Napolitano, Modern Quantum Mechanics, 3rd ed., Cambridge University Press, 2020.
• D. J. Tannor, Introduction to Quantum Mechanics: A Time-Dependent Perspective, University Science Books, 2007.
• S. Blanes, 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.
• E. Hairer, C. Lubich, and G. Wanner, Geometric Numerical Integration, 2nd ed., Springer, 2006.
• SciPy Developers, solve_ivp documentation, consulted for complex-valued adaptive integration and DOP853 solver controls.

Exercises

1. Evaluate an unequal-time commutator

Derive the commutator at $t_1=0$ and $t_2=\pi/(2\nu)$.

Solution

The two effective vectors are

$$\mathbf b(0) = (\Omega,0,\Delta),$$

and

$$\mathbf b\left( \frac{\pi}{2\nu} \right) = (0,\Omega,\Delta).$$

Their cross product is

$$\mathbf b(0)\times \mathbf b\left( \frac{\pi}{2\nu} \right) = (-\Delta\Omega,-\Delta\Omega,\Omega^2).$$

Hence

$$\begin{aligned} [H(0),H(\pi/(2\nu))] &= \frac{i\hbar^2}{2} \left( -\Delta\Omega\sigma_x -\Delta\Omega\sigma_y +\Omega^2\sigma_z \right), \end{aligned}$$

which is nonzero when $\Omega\ne0$.

2. Derive the rotating-frame Hamiltonian

Starting from $\lvert\psi\rangle=R\lvert\chi\rangle$, derive $H_{\rm rot}$.

Solution

Differentiate the state:

$$i\hbar \left( \dot R\lvert\chi\rangle + R\lvert\dot\chi\rangle \right) = RH_0\lvert\chi\rangle.$$

Multiply by $R^\dagger$:

$$i\hbar\lvert\dot\chi\rangle = \left( H_0-i\hbar R^\dagger\dot R \right) \lvert\chi\rangle.$$

Since

$$R^\dagger\dot R = -\frac{i\nu}{2}\sigma_z,$$

one obtains

$$\begin{aligned} H_{\rm rot} &= H_0-\frac{\hbar\nu}{2}\sigma_z \\ &= \frac{\hbar}{2} \left[ \Omega\sigma_x + (\Delta-\nu)\sigma_z \right]. \end{aligned}$$

3. Explain the resonant one-period failure

At $\Delta=\nu$, show why the naive result predicts no transition while the exact result generally does.

Solution

Over one period the transverse integrals vanish, and

$$\mathbf B(T_d) = (0,0,2\pi).$$

Therefore

$$U_{\rm naive}(T_d) = e^{-i\pi\sigma_z} = -I.$$

This changes only the global phase of a $\sigma_z$ basis state. In the rotating frame, resonance gives

$$H_{\rm rot} = \frac{\hbar\Omega}{2}\sigma_x.$$

Since $R(T_d)=-I$,

$$U_{\rm ex}(T_d) = - e^{-i\Omega T_d\sigma_x/2}.$$

The $\sigma_x$ rotation produces

$$P_{0\to1} = \sin^2\left( \frac{\Omega T_d}{2} \right) = \sin^2\left( \frac{\pi\Omega}{\nu} \right).$$

4. Prove stepwise unitarity

Show that a chronological product of frozen-Hamiltonian exponentials is unitary, regardless of the sampling rule.

Solution

For Hermitian $H(t_j^*)$,

$$S_j = e^{-iH(t_j^*)\Delta t/\hbar}$$

obeys

$$S_j^\dagger S_j=I.$$

For

$$U_N=S_{N-1}\cdots S_0,$$

one has

$$\begin{aligned} U_N^\dagger U_N &= S_0^\dagger\cdots S_{N-1}^\dagger S_{N-1}\cdots S_0 \\ &= I. \end{aligned}$$

The cancellations prove exact unitarity in exact arithmetic. They do not prove that the sampled product approximates the correct time-ordered operator accurately.

5. Derive the commuting control propagator

Solve the evolution for $H_c(t)$ and explain why time ordering is unnecessary.

Solution

Every Hamiltonian is proportional to $\sigma_z$, so all pairs commute. Therefore

$$U_c(T,0) = \exp\left[ - \frac{i}{\hbar} \int_0^T H_c(t)\,dt \right].$$

The scalar integral is

$$\int_0^T \left( \Delta+a\cos\nu t \right)dt = \Delta T + \frac{a}{\nu}\sin\nu T.$$

Substitution gives

$$U_c(T,0) = \exp\left\{ - \frac{i\sigma_z}{2} \left[ \Delta T + \frac{a}{\nu}\sin\nu T \right] \right\}.$$

6. Distinguish accuracy from unitarity

Give an example from this notebook of a unitary approximation that is inaccurate and a potentially accurate approximation that is not exactly unitary.

Solution

The reversed product of exact frozen-Hamiltonian exponentials is unitary because it is a product of unitary matrices, but it uses the wrong chronological order and does not converge to the forward propagator. A high-order adaptive Runge–Kutta solution can approximate the exact propagator very accurately, but its finite-step update is not constrained to be exactly unitary. Operator error and unitarity residual must therefore be reported separately.