University Physics V · The Schrödinger Equation · 4.3

Unitary Evolution & the Propagator

The formal solution of the last topic's law, and the object every quantum calculation actually runs on. Learn to read exp(−iĤt/ħ) as a rotation of the state vector: what an operator exponential means, how you really compute it, and why a unitary flow can still pump energy into an atom.

Open the lab ↓Practice & worked examples ↓

Build the model

Connect the measurement to the mechanism.

Read iħ d|Ψ⟩/dt = Ĥ|Ψ⟩ as an ordinary differential equation whose unknown is a vector and whose coefficient is an operator: for a Ĥ with no time dependence it integrates exactly as the scalar equation does, to |Ψ(t)⟩ = Û(t)|Ψ(0)⟩ with Û(t) = exp(−iĤt/ħ), defined by its power series and unambiguous because Ĥ commutes with itself at every pair of times. You never build it by summing operators: diagonalise Ĥ once and Û(t) = Σₙ e^{(−iEₙt/ħ)}|n⟩⟨n|, one phase clock per eigenvector, which is why scipy.linalg.eigh beats calling expm at every t. Self-adjointness is what makes Û unitary, and unitarity conserves far more than the norm — inner products, orthogonality, distinguishability, and reversibility, since Û(−t) undoes Û(t) exactly.

Stone's theorem makes the correspondence exact in both directions, so a merely symmetric Ĥ generates no dynamics at all. Two costs follow. Unitary is not energy-conserving: ⟨Ĥ⟩ stands still only because a static Û is a function of Ĥ and commutes with it, and a driven Ĥ(t) evolves the state unitarily while feeding it energy.

And once Ĥ moves, exp(−(i/ħ)∫Ĥ dt′) is simply the wrong answer unless the Ĥ(t) commute at different times; the true Û is a time-ordered Dyson series, which a solver approximates with steps short enough that the commutator does not matter.

Simple definition: The propagator Û(t) = exp(−iĤt/ħ) is the unitary operator that carries a closed system's state forward by a time t, |Ψ(t)⟩ = Û(t)|Ψ(0)⟩, generated by its time-independent Hamiltonian.
Example: In the energy basis Û is diagonal: for the two lowest levels of an electron in a 1.0 nm box, Û(t) = e^{(−i(0.376 eV)}t/ħ)|1⟩⟨1| + e^{(−i(1.504 eV)}t/ħ)|2⟩⟨2|, so a 50:50 superposition revives after h/(E₂ − E₁) = 3.67 fs.

Formal solution of the law of motion|Ψ(t)⟩ = Û(t)|Ψ(0)⟩, Û(t) = exp(−iĤt/ħ)

One object replaces the differential equation: build it once, apply it at any t, and run it backwards with Û(−t).

Time-independent Ĥ only. Ĥt/ħ is dimensionless: (J)(s)/(J s). Û is an operator, not a number.

Spectral form: what you actually computeÛ(t) = Σₙ e^{(−iEₙt/ħ)}|n⟩⟨n|

eigh once, then any t is a vector of phases. Calling expm at every step repeats an O(N³) job you already own.

Discrete spectrum; a continuum adds ∫dE |E⟩⟨E|. Eₙ in J, t in s, so Eₙt/ħ is an angle in rad.

Group law and unitarityÛ†Û = 𝟙, Û(t₁)Û(t₂) = Û(t₁ + t₂), Û(t)⁻¹ = Û(−t)

Norms and all overlaps survive, so orthogonal states stay orthogonal, probability stays at 1, and the flow is exactly reversible.

Û†Û = 𝟙 holds exactly when Ĥ† = Ĥ on its domain; Û(0) = 𝟙.

The generator — Stone's theoremĤ = iħ (dÛ/dt) at t = 0, Û(δt) = 𝟙 − iĤδt/ħ + O(δt²)

Ĥ is the infinitesimal generator of time translation. Symmetric is not enough: with no self-adjoint extension there is no Û at all.

A strongly continuous one-parameter unitary group has exactly one self-adjoint generator, and conversely.

Why ⟨Ĥ⟩ is frozen, and when it is not⟨Ĥ⟩(t) = ⟨Ψ₀|Û†ĤÛ|Ψ₀⟩ = ⟨Ĥ⟩(0), since [Ĥ, Û] = 0

Unitarity pins the norm; commuting pins the energy. Different facts — a driven system keeps the first while losing the second.

[Ĥ, Û] = 0 only because Û is a function of that same static Ĥ; it fails the moment Ĥ depends on t.

Time-dependent Ĥ: the ordered exponentialÛ(t) = 𝒯 exp(−(i/ħ)∫₀ᵗ Ĥ(t′)dt′)

Û stays unitary but is no longer a function of any one Ĥ. In code: short windows, Ĥ frozen at each midpoint, multiplied in time order.

𝒯 puts later times to the left. Dropping it is legal only if [Ĥ(t₁), Ĥ(t₂)] = 0 for all t₁, t₂.

Integrating an operator equation: what the exponential means

Write the law as d|Ψ⟩/dt = −(i/ħ)Ĥ|Ψ⟩ and it is dz/dt = λz with a vector for z and an operator for λ. The scalar solution is e^(λt)z(0); the operator solution is |Ψ(t)⟩ = exp(−iĤt/ħ)|Ψ(0)⟩, with the exponential defined by the same series, exp(Â) = Σₖ Âᵏ/k!, applied to Â = −iĤt/ħ. Two checks make that legitimate. Differentiate the series term by term and you get −(i/ħ)Ĥ multiplying the same series, so Û really does solve the equation, while Û(0) = 𝟙 supplies the initial condition. And there is no ordering ambiguity, because every term is a power of the single operator Ĥ, which commutes with itself — the innocent-looking fact that fails as soon as Ĥ depends on time. The exponent must be dimensionless, and it is: Ĥt/ħ carries J × s ÷ (J s). For an unbounded Ĥ the series does not converge on every vector in the space, which is why the honest definition of Û is the spectral one below. The series is the motivation; the spectral theorem is the machinery.

Diagonalise, do not exponentiate

Any function of a self-adjoint operator is defined by its action on the eigenvectors: if Ĥ|n⟩ = Eₙ|n⟩ then f(Ĥ) = Σₙ f(Eₙ)|n⟩⟨n|, so Û(t) = Σₙ e^{(−iEₙt/ħ)}|n⟩⟨n|. In the energy basis Û is a diagonal matrix of unit-modulus numbers and applying it costs one multiplication per level. That is the whole practical recipe: run scipy.linalg.eigh on the truncated Ĥ once, an O(N³) job, after which any t costs O(N²) — whereas scipy.linalg.expm(−1j*H*t/ℏ) pays O(N³) again for every t you ask for. In two dimensions you can do it in closed form. Because (n̂⋅σ̂)² = 𝟙, the series splits into even and odd parts and gives exp(−iθ n̂⋅σ̂) = cos θ 𝟙 − i sin θ n̂⋅σ̂. So the ammonia Hamiltonian Ĥ = −Δσ̂ₓ has Û(t) = cos(Δt/ħ)𝟙 + i sin(Δt/ħ)σ̂ₓ, and reading off ⟨2|Û|1⟩ = i sin(Δt/ħ) returns P₂ = sin²(Δt/ħ) with no differential equation solved at all. That is what a propagator buys: the dynamics of every initial state, in one object.

Unitarity conserves more than the norm

Differentiate Û†Û. From dÛ/dt = −(i/ħ)ĤÛ and its dagger dÛ†/dt = +(i/ħ)Û†Ĥ†, the product obeys d(Û†Û)/dt = (i/ħ)Û†(Ĥ† − Ĥ)Û, which vanishes when Ĥ† = Ĥ; since Û†Û = 𝟙 at t = 0, it is 𝟙 forever. Now push on what that buys. Any two states carried by the same Û keep their overlap: ⟨Φ(t)|Ψ(t)⟩ = ⟨Φ₀|Û†Û|Ψ₀⟩ = ⟨Φ₀|Ψ₀⟩. Norm conservation is only the special case Φ = Ψ. Orthogonality survives too, so two perfectly distinguishable states never blur together, and two states you cannot tell apart now you never will. That one line drives the no-cloning theorem: a copier would have to turn the overlap ⟨Φ|Ψ⟩ into ⟨Φ|Ψ⟩², and unitarity forbids any change at all unless the overlap is already 0 or 1. Reversibility comes free from the group law, Û(−t)Û(t) = Û(0) = 𝟙: the propagator has an exact inverse, so nothing is lost between measurements, and any irreversibility in physics is imported from outside this law.

Stone's theorem: the group and its generator are one object

The group law Û(t₁)Û(t₂) = Û(t₁ + t₂) with Û(0) = 𝟙 makes the family (Û(t)) a one-parameter group, and it is strongly continuous because Û(t)|Ψ⟩ → |Ψ⟩ as t → 0. Stone's theorem says every such group is exp(−iĤt/ħ) for exactly one self-adjoint Ĥ, and conversely that every self-adjoint Ĥ generates one. Two consequences matter here. First, expanding for small δt gives Û(δt) = 𝟙 − iĤδt/ħ + O(δt²), so Ĥ = iħ(dÛ/dt) at t = 0: the Hamiltonian is the infinitesimal generator of time translation, and calling its eigenvalues energies is a naming convention laid on top of that. Second, the theorem is strict about self-adjointness and merely symmetric will not do. Take p̂ = −iħ d/dx on the half-line x > 0: symmetric on smooth functions vanishing at the origin, but its deficiency indices are (1,0), so it has no self-adjoint extension. Physically, the group it would generate is a rigid translation, and sliding a function towards the origin pushes probability off the end of the line. Naming a domain is not pedantry — it is the difference between having dynamics and not.

Unitary is not the same as energy-conserving

It is tempting to read Û†Û = 𝟙 as a conservation law for energy. It is not. Norm conservation follows from Ĥ† = Ĥ alone. The constancy of ⟨Ĥ⟩ needs a separate and stronger fact: Û is built out of Ĥ, so [Ĥ, Û] = 0, and ⟨Ĥ⟩(t) = ⟨Ψ₀|Û†ĤÛ|Ψ₀⟩ = ⟨Ψ₀|ĤÛ†Û|Ψ₀⟩ = ⟨Ĥ⟩(0). Break that premise and the conclusion goes with it. Drive a two-level atom on resonance with Ĥ(t) = (ħω₀/2)σ̂_z + ħΩ_R cos(ω₀t)σ̂ₓ: the operator is self-adjoint at every instant, so the evolution is unitary and ⟨Ψ|Ψ⟩ never leaves 1 — yet a π-pulse carries the atom from the lower level to the upper one and ⟨Ĥ₀⟩ climbs by ħω₀. For an electron spin in 0.500 T that is 5.79 × 10⁻⁵ eV, delivered at 14.0 GHz, and it came from the drive rather than from nowhere. The general statement sits in the same commutator: d⟨Â⟩/dt = (i/ħ)⟨[Ĥ, Â]⟩ + ⟨∂Â/∂t⟩, so any observable commuting with a static Ĥ is conserved, and Ĥ itself is conserved only while it stands still.

When Ĥ moves, and what a solver must not do

For a time-dependent Ĥ the naive answer exp(−(i/ħ)∫₀ᵗ Ĥ(t′)dt′) averages the Hamiltonian and then exponentiates the average — legal only if [Ĥ(t₁), Ĥ(t₂)] = 0 for every pair of times. Baker–Campbell–Hausdorff prices the mistake: e^Âe^B̂ = exp(Â + B̂ + ½[Â, B̂] + …), so what you dropped is the commutator. The correct object is the time-ordered exponential, whose expansion is the Dyson series with later times always to the left. In code, chop the interval into windows short enough that Ĥ barely turns, exponentiate each with Ĥ frozen at its midpoint, and multiply them in time order; the error then falls as the square of the window length. What you must never do is truncate the exponential itself. The Euler step 𝟙 − iĤΔt/ħ multiplies an eigen-component by 1 − iα with α = EΔt/ħ, and |1 − iα| = √(1 + α²) > 1: at α = 0.1, a hundred steps inflate the norm by 64%. The Cayley form (𝟙 − iĤΔt/2ħ)(𝟙 + iĤΔt/2ħ)⁻¹ has modulus exactly 1 by construction, with a phase error of only α − 2 arctan(α/2) ≈ α³/12 per step — which is why Crank–Nicolson, not Euler, is the default in wavepacket codes.

Change one variable at a time

Make the relationship visible.

Interactive model

Phase per step α = EΔt/ħ0.20 rad

Number of steps n6

Push α up and the filled marker leaves the unit circle: a truncated exponential is not unitary, and the norm it invents grows as (1 + α²)^(n/2). Pull α back to 0.05 and add steps — same total turn, ring back on the circle. Small steps, or diagonalise Ĥ and exponentiate once.

EULER NORM |c| AFTER n1.125

EXACT TURN αn1.20 rad

EULER PHASE ERROR0.016 rad

CAYLEY PHASE ERROR0.0040 rad

Live interpretationEULER NORM |c| AFTER n: 1.125. EXACT TURN αn: 1.20 rad. EULER PHASE ERROR: 0.016 rad. CAYLEY PHASE ERROR: 0.0040 rad

Catch the common trap

Explain before calculating.

Misconception check

Two operator facts get lost because the exponential looks so much like the scalar case. First, for a time-dependent Ĥ the propagator is not exp(−(i/ħ)∫₀ᵗ Ĥ(t′)dt′). That expression averages Ĥ over the interval and exponentiates the average, silently assuming the Ĥ at different instants commute; Baker–Campbell–Hausdorff shows exactly what is dropped, the term ½[Â, B̂]. The right object is the time-ordered Dyson series, and the fix in code is shorter windows, not a bigger integral. Second, unitary does not mean energy-conserving. Û†Û = 𝟙 follows from Ĥ† = Ĥ alone and holds whatever Ĥ does; ⟨Ĥ⟩ stands still for the separate reason that a static Û commutes with its own Ĥ. Shine a laser on an atom and the norm never leaves 1 while the atom absorbs ħω₀ — unitary throughout, and energy manifestly not conserved.

A two-level atom is driven on resonance, so Ĥ(t) = (ħω₀/2)σ̂_z + ħΩ_R cos(ω₀t)σ̂ₓ — self-adjoint at every instant, but with [Ĥ(t₁), Ĥ(t₂)] ≠ 0. Which statement about the evolution from 0 to t is correct?

Choose an answer to test the model.

Practice & worked examples

Reason from the model, then test the result.

EasyAn electron in a 1.0 nm box starts in |Ψ(0)⟩ = (|1⟩ + |2⟩)/√2, with E₁ = 0.376 eV and E₂ = 1.504 eV. Write Û(t) in the energy basis, evolve the state, and find the first time it returns to |Ψ(0)⟩ up to a global phase.

Û is a function of Ĥ, so it is diagonal in the energy basis: Û(t) = e^{(−iE₁t/ħ)}|1⟩⟨1| + e^{(−iE₂t/ħ)}|2⟩⟨2|. No matrix exponential is needed — a function of a diagonal matrix acts entry by entry.
|Ψ(t)⟩ = Û(t)|Ψ(0)⟩ = (1/√2)(e^{(−iE₁t/ħ)}|1⟩ + e^{(−iE₂t/ħ)}|2⟩). Each coefficient keeps modulus 1/√2, so |c₁|² = |c₂|² = 0.5 for all t: Û never moves probability between the eigenlevels of its own Ĥ.
Pull out the common phase: |Ψ(t)⟩ = e^{(−iE₁t/ħ)}(1/√2)(|1⟩ + e^{(−i(E₂−E₁)}t/ħ)|2⟩). Only the relative phase is observable, so the state repeats when (E₂ − E₁)t/ħ = 2π.
E₂ − E₁ = 1.504 − 0.376 = 1.128 eV, so t = h/(E₂ − E₁) = 4.136 × 10⁻¹⁵ eV s ÷ 1.128 eV = 3.67 × 10⁻¹⁵ s.

AnswerÛ(t) = e^{(−iE₁t/ħ)}|1⟩⟨1| + e^{(−iE₂t/ħ)}|2⟩⟨2|; the level populations never move, and the state revives up to a global phase every h/(E₂ − E₁) = 3.67 fs.

MediumA spin-½ sits in a static field along z, so Ĥ = (ħω₀/2)σ̂_z with ω₀ = 8.80 × 10¹⁰ rad s⁻¹ (an electron in 0.500 T). The state starts as |+x⟩ = (|↑⟩ + |↓⟩)/√2. Build Û(t) from exp(−iθ n̂⋅σ̂) = cos θ 𝟙 − i sin θ n̂⋅σ̂, evolve the state, and contrast the period of the spinor with the period of the observable.

Ĥt/ħ = (ω₀t/2)σ̂_z, so with θ = ω₀t/2 and n̂ = ẑ: Û(t) = cos(ω₀t/2)𝟙 − i sin(ω₀t/2)σ̂_z = diag(e^{(−iω₀t/2)}, e^(+iω₀t/2)). The identity holds because σ̂_z² = 𝟙 splits the series into cosine and sine parts.
Unitarity check: Û†Û = (cos θ 𝟙 + i sin θ σ̂_z)(cos θ 𝟙 − i sin θ σ̂_z) = cos²θ 𝟙 + sin²θ σ̂_z² = 𝟙, the cross terms cancelling.
|Ψ(t)⟩ = (1/√2)(e^{(−iω₀t/2)}|↑⟩ + e^(+iω₀t/2)|↓⟩). With a = e^{(−iω₀t/2)}/√2 and b = e^(+iω₀t/2)/√2, a*b = ½e^(+iω₀t), so ⟨σ̂ₓ⟩ = 2Re(a*b) = cos(ω₀t) and ⟨σ̂_y⟩ = 2Im(a*b) = sin(ω₀t).
The spinor's phases turn at ω₀/2, so |Ψ⟩ needs t = 4π/ω₀ = 1.43 × 10⁻¹⁰ s to return to itself, while the observable turns at ω₀ and repeats in half that, T = 2π/ω₀ = 7.14 × 10⁻¹¹ s. The minus sign after one 2π is real but unobservable — a global phase.
⟨σ̂_z⟩ = |a|² − |b|² = 0 at every t, as [Ĥ, σ̂_z] = 0 requires, and ⟨Ĥ⟩ = (ħω₀/2)⟨σ̂_z⟩ = 0 throughout because Û is built from this same Ĥ.

AnswerÛ(t) = diag(e^{(−iω₀t/2)}, e^(+iω₀t/2)); ⟨σ̂ₓ⟩ = cos ω₀t and ⟨σ̂_y⟩ = sin ω₀t. The Bloch vector comes back in 71.4 ps, the spinor needs 143 ps, and ⟨Ĥ⟩ stays 0 because [Ĥ, Û] = 0.

HardA qubit is pulsed twice: Ĥ = (ħΩ/2)σ̂ₓ for 0 < t < τ, then Ĥ = (ħΩ/2)σ̂_z for τ < t < 2τ, with Ωτ = π/2. Starting from |↑⟩, find P(↓) at t = 2τ from the correctly ordered product, then from the naive exp(−(i/ħ)∫Ĥ dt′), and name the commutator responsible for the gap.

Each piece is constant, so each has an exact propagator: Û₁ = exp(−i(Ωτ/2)σ̂ₓ) = cos(π/4)𝟙 − i sin(π/4)σ̂ₓ = (1/√2)(𝟙 − iσ̂ₓ), and likewise Û₂ = (1/√2)(𝟙 − iσ̂_z).
Order matters — the later pulse acts last: Û(2τ) = Û₂Û₁ = ½(𝟙 − iσ̂_z)(𝟙 − iσ̂ₓ) = ½(𝟙 − iσ̂ₓ − iσ̂_y − iσ̂_z), using σ̂_zσ̂ₓ = iσ̂_y. Unitarity check: (½)² + 3 × (½)² = 1.
In matrices, σ̂ₓ + σ̂_y + σ̂_z = [[1, 1−i], [1+i, −1]], so Û = ½[[1−i, −(1+i)], [1−i, 1+i]] and Û|↑⟩ = ½(1−i, 1−i)ᵀ. Hence P(↓) = |1−i|²/4 = 2/4 = 0.500.
Now the naive route: ∫₀²τ Ĥ dt′ = (ħΩτ/2)(σ̂ₓ + σ̂_z), so the exponent is −i(π/4)(σ̂ₓ + σ̂_z) = −iθ n̂⋅σ̂ with n̂ = (x̂ + ẑ)/√2 and θ = π√2/4 = 1.1107 rad.
⟨↓|n̂⋅σ̂|↑⟩ = 1/√2, so ⟨↓|Û′|↑⟩ = −i sin θ/√2 and P(↓) = sin²(1.1107)/2 = 0.8029/2 = 0.401 — low by 0.099, a 20% error on the true value.
The culprit is [σ̂ₓ, σ̂_z] = −2iσ̂_y ≠ 0: Baker–Campbell–Hausdorff gives e^Âe^B̂ = exp(Â + B̂ + ½[Â, B̂] + …), and the naive exponential is exactly the version with ½[Â, B̂] thrown away.

AnswerOrdered product Û₂Û₁ gives P(↓) = 0.500; the naive single exponential gives sin²(1.1107)/2 = 0.401, a 20% error. The culprit is [σ̂ₓ, σ̂_z] = −2iσ̂_y ≠ 0, precisely the ½[Â, B̂] term Baker–Campbell–Hausdorff adds.