The scattering matrix and R + T = 1

Start from ρ = |ψ|² and iℏ ∂ψ/∂t = −(ℏ²/2m)ψ″ + Vψ. Then ∂ρ/∂t = ψ*∂ψ/∂t + ψ∂ψ*/∂t, and substituting the equation splits it in two. The kinetic piece reorganises into (iℏ/2m)(ψ*ψ″ − ψψ*″), which is exactly −∂/∂x of j = (ℏ/2mi)(ψ*ψ′ − ψψ*′) = (ℏ/m) Im(ψ*ψ′). The potential piece is (2/ℏ) Im(V)|ψ|², and it vanishes only because V is real: that is the one and only place hermiticity of H is used. So ∂ρ/∂t + ∂j/∂x = 0. Now impose stationarity. A solution of Hψ = Eψ, which for scattering is a δ-normalised generalised eigenfunction rather than a vector in the Hilbert space, evolves as ψ(x) exp(−iEt/ℏ), whose density does not move, so ∂ρ/∂t = 0 and dj/dx = 0. In one dimension that is far stronger than it looks: it is not a statement about a flow field, it says j is a single constant for the entire real line. Compute it at x = −∞ and you have computed it at x = +∞, through whatever sits in between.

S exists only where the potential flattens out. Require V → V_L as x → −∞ and V → V_R as x → +∞, and set k_L = √(2m(E−V_L))/ℏ, k_R = √(2m(E−V_R))/ℏ. For E above both floors the asymptotic solution is ψ → A exp(ik_L x) + B exp(−ik_L x) on the left and C exp(ik_R x) + D exp(−ik_R x) on the right. Four complex numbers, but the equation is second order, so only two are free. The split into incoming and outgoing is made by whether the group velocity points toward the scatterer, not by its sign or by position alone: A runs right and lives on the left, D runs left and lives on the right, so A and D are what you prepare and B and C are what comes back. Define the amplitude ratios r = B/A and t = C/A when D = 0, and r′ = C/D and t′ = B/D when A = 0. S maps the incoming pair onto the outgoing pair once each amplitude carries its flux weight √k, so its entries are r, r′, √(k_R/k_L) t and √(k_L/k_R) t′. Nothing in that definition mentions the shape of V in between. S is a boundary-value object, and two very different potentials sharing one S are indistinguishable to any scattering experiment.

Evaluate j on the left. With ψ = A exp(ikx) + B exp(−ikx), the product ψ*ψ′ = ik|A|² − ik|B|² + ik(A B* exp(2ikx) − A* B exp(−2ikx)). The bracket has the form z − z*, which is purely imaginary, so ik times it is purely real and drops out of Im(ψ*ψ′) exactly. The standing-wave ripple in |ψ|², the most conspicuous thing in any plot of the left region, therefore carries no current at all. What survives is j_L = (ℏk_L/m)(|A|² − |B|²), and identically j_R = (ℏk_R/m)(|C|² − |D|²). Set D = 0 for incidence from the left and impose dj/dx = 0: k_L|A|² = k_L|B|² + k_R|C|². Divide through by k_L|A|² and you have 1 = |r|² + (k_R/k_L)|t|², which is R + T = 1 with R and T read off as flux ratios. The factor k_R/k_L is not cosmetic: at the step in the worked examples it turns |t|² = 1.78 into T = 0.889.

Rescale so the conserved quantity is an ordinary norm: set a = √(k_L)A, b = √(k_L)B, c = √(k_R)C and d = √(k_R)D. Then j_L = j_R reads |a|² + |d|² = |b|² + |c|², so the map S taking the pair (a, d) to the pair (b, c) preserves length, and a length-preserving linear map on ℂ² is unitary. Hence S†S = 1. With S₁₁ = r, S₂₂ = r′, S₂₁ = √(k_R/k_L) t and S₁₂ = √(k_L/k_R) t′, its diagonal entries are R + T = 1 and R′ + T′ = 1; its off-diagonal entry is S₁₁*S₁₂ + S₂₁*S₂₂ = 0. That last relation is not decorative. Take moduli, then use the two diagonal relations to eliminate |S₂₁| and |S₁₂|, and you get |r|² = |r′|², hence |S₂₁| = |S₁₂|, which is T = T′: the transmission probability through a barrier is the same from either side, however lopsided the barrier. Numerically this is also the best diagnostic you have. Assemble S, form S†S − 1, and report its largest entry; on a healthy calculation it sits near 10⁻¹⁵.

M is defined by (C, D)ᵀ = M(A, B)ᵀ, and its virtue is that it composes: a staircase potential is the ordered product of one matrix per interface and per slab. Its vice is that it is not unitary. Flux conservation reads M†σ₃M = (k_L/k_R)σ₃ instead, which fixes |det M| = k_L/k_R and leaves M with one growing and one decaying direction, so its condition number climbs like exp(2κL). The two matrices are related by rearrangement, not by similarity. Solving (C, D)ᵀ = M(A, B)ᵀ for the outgoing pair gives r = −M₂₁/M₂₂, t′ = 1/M₂₂, t = det M/M₂₂ and r′ = M₁₂/M₂₂, and giving the two transmission amplitudes their flux weights turns that into S = (1/M₂₂)[[−M₂₁, √(k_L/k_R)],[√(k_R/k_L) det M, M₁₂]]. So the working rule is: multiply in M, then convert to S to read physics and to test unitarity. When a thick-barrier product has quietly lost its decaying component to round-off, the symptom that shows up is R + T ≠ 1.

Keep the derivation but let V = U − iW with U real and W positive. The potential term no longer cancels: ∂ρ/∂t + ∂j/∂x = −(2W/ℏ)|ψ|², so in a stationary state dj/dx = −(2W/ℏ)|ψ|², strictly negative wherever W acts, and the current steps down as it crosses the absorbing region. Then S†S ≤ 1 as a matrix inequality: its diagonal entries R + T and R′ + T′ are the survival probabilities for incidence from either side, its eigenvalues are the smallest and largest survival probability over all incoming superpositions, and P_abs = 1 − R − T is the flux removed. This is the optical model of nuclear physics: a nucleus is given a complex potential precisely because elastic scattering is not the only thing a neutron can do, and P_abs folds capture, inelastic excitation and fission into the one number the elastic model refuses to resolve. The figure below makes the trade concrete. A delta of strength α − iW, with W now in eV⋅nm like α, gives P_abs = 2w/((1+w)² + β²) with w = mW/(ℏ²k), which peaks at exactly 0.5: a purely absorbing point scatterer can never eat more than half the beam, because to absorb it must first fail to reflect.