University Physics V · Finite Potential Wells and Quantum Tunneling · 6.3
Transfer Matrices for Piecewise-Constant Potentials
Any potential you can draw as a staircase is already solved, once you can solve one wall and one slab. This lesson turns each into a 2×2 matrix, multiplies them in the order the electron meets them, and then shows you exactly where that innocent product destroys itself in floating point.
Build the model
Connect the measurement to the mechanism.
A piecewise-constant potential leaves no differential equation to solve. Inside each slab the equation is ψ″ = −k²ψ with k constant, so the solution there is pinned by exactly two complex numbers: the right-going and left-going amplitudes (A, B). All the physics sits in how that pair changes from slab to slab, and because ψ and ψ′ are continuous across a finite jump, the change is linear — one 2×2 matrix per wall, one diagonal 2×2 per slab.
Chain them in the order the wave meets them and a staircase of any complexity collapses to a single M, with the whole scattering problem readable from two of its entries. What M is not is unitary. It conserves the σ₃-weighted probability current, M†σ₃M = (kL/kR)σ₃, and that indefinite metric leaves its entries unbounded: in a forbidden slab the propagation matrix is diag(e^−κd, e^+κd), stretching one direction by eκd and squeezing the other by the same factor.
That is the price. The exact staircase solution arrives in a representation whose condition number is e²κL, and by κL ≈ 18 the decaying solution has been rounded out of a double-precision number entirely. The transfer matrix is the right object to think with and, past a certain thickness, the wrong object to compute with.
- Simple definition
- A transfer matrix is the 2×2 matrix that carries the pair of plane-wave amplitudes (A, B) from one constant-potential region across to the next, so an entire staircase potential is solved by one ordered product of such matrices.
- Example
- At a single step from E = 2.00 eV onto a 1.00 eV plateau, s = kL/kR = √2 and M = ½[[1+√2, 1−√2], [1−√2, 1+√2]]; det M = √2 = kL/kR, and r = −M₂₁/M₂₂ = 0.1716, so R = 2.94%.
Two complex numbers replace the whole function in a slab, which is why a 2×2 matrix is enough to step from one slab to the next.
kⱼ = √(2m(E − Vⱼ))/ħ in m⁻¹; for E < Vⱼ write kⱼ = iκⱼ with κⱼ real. Aⱼ, Bⱼ are complex amplitudes.
One matrix per wall, fixed entirely by the ratio of the wavenumbers either side — the wall's position never enters it.
From continuity of ψ and ψ′ across a finite jump. det I = s, so it is unimodular only when the two wavenumbers agree.
Its condition number is e(2κⱼdⱼ): every forbidden slab crossed costs 2κd/ln 10 significant decimal digits.
dⱼ = xⱼ₊₁ − xⱼ is the slab width in m; det Pⱼ = 1, so slabs never touch the determinant.
The rightmost factor acts first, so the matrices are written in reverse spatial order, and they do not commute.
det M = kL/kR exactly: the interface determinants telescope and the slab determinants are 1.
The conserved quantity carries an indefinite metric, so M is not unitary and none of its entries is bounded by 1.
j = (ħk/m)(|A|² − |B|²) is the current in a lead, in s⁻¹. Entrywise: |M₁₁|² − |M₂₁|² = kL/kR.
T is a ratio of currents, not of amplitudes: |t|² can exceed 1 whenever the right lead is slower than the left.
From (t, 0)ᵀ = M (1, r)ᵀ. T collapses to (kL/kR)/|M₂₂|², and to 1/|M₂₂|² for matched leads.
Two amplitudes are the entire state inside a slab
Set V(x) = Vⱼ on the interval [xⱼ, xⱼ₊₁]. There H|ψ⟩ = E|ψ⟩ reduces to ψ″ = −kⱼ²ψ with kⱼ = √(2m(E − Vⱼ))/ħ constant, whose general solution is Aⱼe(ikⱼ(x−xⱼ)) + Bⱼe(−ikⱼ(x−xⱼ)). Two constants, so the state inside that slab is a vector (Aⱼ, Bⱼ) in ℂ² — the two-fold degeneracy of the continuum, written as a basis. Anchoring the exponentials at the slab's own left edge rather than at a global origin is a bookkeeping choice, and it is the one that keeps the matrices free of stray phases. When E < Vⱼ nothing in the algebra changes: kⱼ = iκⱼ with κⱼ = √(2m(Vⱼ − E))/ħ real, and the right-goer becomes the decaying exponential while the left-goer becomes the growing one. For an electron ħ²/2mₑ = 3.810 × 10⁻² eV nm², so k = 5.123 √(E/eV) nm⁻¹: one electronvolt buys about 5.1 inverse nanometres.
One matrix per wall, one per slab
At the wall xⱼ₊₁ the potential jumps but stays finite, so ψ and ψ′ are both continuous there. Let a = Aⱼe(ikⱼdⱼ) and b = Bⱼe(−ikⱼdⱼ) be the amplitudes carried across the slab of width dⱼ = xⱼ₊₁ − xⱼ. Continuity of ψ gives Aⱼ₊₁ + Bⱼ₊₁ = a + b; continuity of ψ′ gives kⱼ₊₁(Aⱼ₊₁ − Bⱼ₊₁) = kⱼ(a − b). Add and subtract, and with s = kⱼ/kⱼ₊₁ you get Aⱼ₊₁ = ½[(1+s)a + (1−s)b] and Bⱼ₊₁ = ½[(1−s)a + (1+s)b]. That is the interface matrix I = ½[[1+s, 1−s], [1−s, 1+s]], with det I = ¼[(1+s)² − (1−s)²] = s. The carrying step is the diagonal Pⱼ = diag(e(ikⱼdⱼ), e(−ikⱼdⱼ)), determinant 1, which in a classically forbidden slab reads diag(e(−κⱼdⱼ), e(+κⱼdⱼ)). Every piecewise-constant problem in one dimension is built from those two matrices and nothing else.
The staircase is a single ordered product
Compose them: M = I_(N−1→N) P_(N−1) ⋯ P₁ I_(0→1), so that (AR, BR)ᵀ = M (AL, BL)ᵀ. Because a matrix acts on what stands to its right, the leftmost interface is the rightmost factor — the product is written in reverse spatial order, and getting it backwards silently builds the mirror-image stack. Flux forces that mirror image to have the same T and the same R, but a different reflection phase and a completely different interior field, so the mistake hides from the two numbers you are most likely to check. Determinants are the free diagnostic: slab matrices contribute 1, interface matrices contribute kⱼ/kⱼ₊₁, and the chain telescopes to det M = k₀/kN = kL/kR. Code that returns anything else has an ordering bug, a wrong wavenumber, or a lost phase, and it should be caught there rather than in a transmission curve.
Flux, not norm: what M really conserves
The probability current in a lead is j = (ħkⱼ/m)(|Aⱼ|² − |Bⱼ|²) = (ħkⱼ/m) v†σ₃v with v = (Aⱼ, Bⱼ) and σ₃ = diag(1, −1). Equating the current in the two leads for every incoming state gives M†σ₃M = (kL/kR)σ₃: M preserves an indefinite form, not a norm. Taking determinants of that identity recovers |det M| = kL/kR; entry by entry it gives |M₁₁|² − |M₂₁|² = kL/kR, |M₂₂|² − |M₁₂|² = kL/kR, and M₁₁*M₁₂ = M₂₁*M₂₂. Because the metric is σ₃ and not the identity, the flux-rescaled matrix √(kR/kL) M sits in SU(1,1), a group of hyperbolic boosts whose entries grow without bound. Time reversal is a separate symmetry: for real V and propagating leads ψ* also solves the equation, which swaps (A, B) into (B*, A*) and forces σ₁Mσ₁ = M*, so M₂₂ = M₁₁* and M₂₁ = M₁₂*. Then M = [[a, b], [b*, a*]] with |a|² − |b|² = kL/kR — unimodular only when the leads match. Give V an absorptive imaginary part and the conjugate pairing goes with it.
Reading r, t, R and T off two entries
Scattering fixes the boundary condition: unit amplitude in from the left, nothing in from the right, so vL = (1, r) and vR = (t, 0). The second row of (t, 0)ᵀ = M(1, r)ᵀ gives 0 = M₂₁ + M₂₂ r, hence r = −M₂₁/M₂₂; the first row, tidied with the determinant, gives t = det M/M₂₂ = (kL/kR)/M₂₂. Probabilities are current ratios: R = |r|² and T = (kR/kL)|t|² = (kL/kR)/|M₂₂|², and R + T = 1 then follows from the pseudo-unitarity identity instead of being imposed by hand. Run it on a rectangular barrier with matched leads and the three factors multiply out to M₂₂ = cosh(κL) − (i/2)(k/κ − κ/k) sinh(κL). Squaring, |M₂₂|² = 1 + sinh²(κL)(k² + κ²)²/(4k²κ²) = 1 + V₀² sinh²(κL)/(4E(V₀ − E)), so T = 1/|M₂₂|² is the textbook barrier formula — reached by multiplying three matrices rather than solving eight matching equations.
Where the product dies, and what to use instead
Each forbidden slab contributes P = diag(e(−κd), e(+κd)), singular values e(∓κd) and condition number e(2κd). Those numbers multiply, so a total forbidden path L costs e(2κL) whether you cross it in one slab or five hundred: slicing finer never helps, because the exponents merely add. Double precision carries ε = 2.22 × 10⁻¹⁶, so significance is exhausted at 2κL = ln(1/ε) = 36.0, which is 15.65 decimal digits; for an electron 4.00 eV below a barrier top, κ = 10.25 nm⁻¹ and that width is only 1.76 nm. What survives is whatever lives in the growing direction: M₂₂ is dominated by e(+κL) and keeps full relative accuracy, so T = 1/|M₂₂|² is reliable even at 10⁻¹⁶. What dies is whatever needs the decaying direction — the field inside the barrier, the resonance peaks of a double-barrier stack where terms of size e(2κL) must cancel down to order 1, and the root M₂₂ = 0 that locates a bound state. The cure is to compose scattering matrices: S is unitary, every entry has modulus at most 1, and the star product t = t₂(1 − r₁′r₂)⁻¹t₁ holds no growing exponential.
Change one variable at a time
Make the relationship visible.
Take L to 1.6 nm with the barrier 6 eV above E: 2 κ L passes 36 and both dots cross the dashed lines, where the decaying solution is gone. Then raise the slab count - cutting the same barrier finer moves the boundaries but never the dots, because the exponents just add.
DECAY CONSTANT κ7.25 nm⁻¹
2κL AT THE EXIT7.25
DIGITS LOST3.1 of 15.65
DIGITS PER SLAB3.15
Live interpretationDECAY CONSTANT κ: 7.25 nm⁻¹. 2κL AT THE EXIT: 7.25. DIGITS LOST: 3.1 of 15.65. DIGITS PER SLAB: 3.15
Catch the common trap
Explain before calculating.
A transfer-matrix product for a barrier with matched leads (kL = kR) returns M₂₂ = 18.73 + 6.61i. What is the transmission probability T?
Choose an answer to test the model.
Practice & worked examples
Reason from the model, then test the result.
EasyAn electron with E = 2.00 eV crosses a single potential step at x = 0, from V = 0 up to V = 1.00 eV. Build the transfer matrix, then read off r, t, R and T, and check that flux balances.
- Wavenumbers first. With ħ²/2mₑ = 3.810 × 10⁻² eV nm², k = 5.123 √(E/eV) nm⁻¹. Left lead: kL = 5.123√2.00 = 7.245 nm⁻¹. To the right the electron keeps E − V = 1.00 eV, so kR = 5.123√1.00 = 5.123 nm⁻¹.
- One wall and no slab, so M is a single interface matrix with s = kL/kR = 7.245/5.123 = 1.4142: M = ½[[1+s, 1−s], [1−s, 1+s]] = [[1.2071, −0.2071], [−0.2071, 1.2071]].
- Check the determinant before trusting anything: det M = 1.2071² − 0.2071² = 1.4571 − 0.0429 = 1.4142, which is kL/kR exactly, as pseudo-unitarity demands.
- Impose no wave arriving from the right: (t, 0)ᵀ = M (1, r)ᵀ. The second row gives 0 = M₂₁ + M₂₂ r, so r = −M₂₁/M₂₂ = 0.2071/1.2071 = 0.1716. The first row then gives t = det M/M₂₂ = 1.4142/1.2071 = 1.1716.
- Convert to probabilities as current ratios. R = |r|² = 0.0294. Note |t|² = 1.3726 exceeds 1 and is therefore not a probability; T = (kR/kL)|t|² = 1.3726/1.4142 = 0.9706, and R + T = 0.0294 + 0.9706 = 1.0000.
Answerr = 0.1716, t = 1.1716, R = 0.0294, T = 0.9706, and R + T = 1.0000. The transmitted amplitude is larger than the incident one; the transmitted current is not, because the right lead is slower.
MediumAn electron of E = 1.00 eV meets a rectangular barrier of height V₀ = 3.00 eV and width L = 0.500 nm, with V = 0 on both sides. Form M as the ordered product I⋅P⋅I, extract M₂₂, and get T. Check it against the closed-form barrier result and against the opaque limit.
- Leads: k = 5.123√1.00 = 5.123 nm⁻¹ on both sides, so they match. Inside, E < V₀, so k → iκ with κ = 5.123√(3.00 − 1.00) = 7.245 nm⁻¹, giving κL = 7.245 × 0.500 = 3.6226.
- Three factors, applied right to left in the order the electron meets them: M = I(κ→k) P(κ) I(k→κ), with interface ratios s₁ = k/(iκ) = −0.7071i and s₂ = iκ/k = 1.4142i, and P = diag(e(−κL), e(+κL)) = diag(0.02671, 37.434).
- Only M₂₂ is needed. Multiplying out, the real parts collect into cosh and the imaginary parts into sinh: M₂₂ = cosh(κL) − (i/2)(k/κ − κ/k) sinh(κL) = 18.7303 − (i/2)(0.7071 − 1.4142)(18.7036) = 18.7303 + 6.6127i.
- |M₂₂|² = 350.82 + 43.73 = 394.55. Matched leads give det M = 1, so t = 1/M₂₂ and T = 1/|M₂₂|² = 2.535 × 10⁻³.
- Closed-form check: T = [1 + V₀² sinh²(κL)/(4E(V₀ − E))]⁻¹ = [1 + 9.00 × 349.82/8.00]⁻¹ = 1/394.55, identical.
- Opaque-limit check: 16E(V₀ − E)/V₀² × e(−2κL) = (32.0/9.00) × 7.137 × 10⁻⁴ = 2.537 × 10⁻³, high by 0.11%, because κL = 3.62 is only moderately large.
AnswerT = 2.535 × 10⁻³. The whole rectangular-barrier formula is one matrix entry, M₂₂ = cosh κL − (i/2)(k/κ − κ/k) sinh κL, squared and inverted — no eight-equation matching required.
HardAn electron sits 4.00 eV below the top of a barrier. Show that the transfer-matrix product's condition number is e(2κL), find the width at which IEEE double precision has no significant digits left, evaluate T there, and say which computed quantities survive and which do not.
- κ = 5.123√4.00 = 10.246 nm⁻¹. A forbidden slab of width d contributes P = diag(e(−κd), e(+κd)), whose singular values are e(−κd) and e(+κd), so cond₂(P) = e(2κd). Interface matrices have bounded entries, so a product spanning a total forbidden path L inherits cond ≈ e(2κL) — and cutting L into more slabs cannot help, since the exponents simply add.
- Double precision has ε = 2.22 × 10⁻¹⁶, so the growing and decaying directions become indistinguishable when e(2κL) = 1/ε, that is 2κL = ln(1/ε) = 36.04, or log₁₀(1/ε) = 15.65 decimal digits.
- Solve for the width: L* = 36.04/(2 × 10.246 nm⁻¹) = 1.759 nm. That is an ordinary tunnel gap rather than a contrived one; an STM works at roughly a third of it.
- Transmission there: 4.00 eV above a 1.00 eV electron means V₀ = 5.00 eV, so 16E(V₀ − E)/V₀² = 16(1.00)(4.00)/25.0 = 2.56 and T ≈ 2.56 e(−36.04) = 2.56 × 2.22 × 10⁻¹⁶ = 5.7 × 10⁻¹⁶.
- What survives: M₂₂ is dominated by the growing exponential, keeps full relative accuracy, and so T = 1/|M₂₂|² is trustworthy even at 10⁻¹⁶. What dies is anything in the decaying direction — the wavefunction inside the barrier, the resonance peaks of a double-barrier stack (a cancellation of e(2κL)-sized terms down to order 1), and the root M₂₂ = 0 that locates a bound state.
- The fix is to compose scattering matrices instead. Every entry of S has modulus at most 1, and the Redheffer star product t = t₂(1 − r₁′r₂)⁻¹t₁, r = r₁ + t₁′r₂(1 − r₁′r₂)⁻¹t₁ contains no growing exponential; |r₁′r₂| < 1 for any transmitting stack, so the inverse never blows up.
Answercond ≈ e(2κL); all 15.65 digits are gone by L* = 1.76 nm, where T ≈ 5.7 × 10⁻¹⁶. T itself survives, living in the growing entry M₂₂; interior fields, resonance peaks and bound-state roots do not, and those need the S-matrix recursion.