The shell model and spin–orbit splitting

The shell model trades the A-body Hamiltonian for a one-body one, Ĥ = p̂²/2m_N + V_WS(r) + Vₗₛ(r) l̂⋅ŝ/ħ², with V_WS(r) = −V₀/[1 + exp((r − R)/a)], V₀ ≈ 50 MeV, R = 1.25 A¹⁄³ fm and a ≈ 0.55 fm. Because Ĥ is one-body its eigenkets are single-nucleon orbitals, and the nucleus is built by filling them. Work in the coupled basis |n l j mⱼ⟩ = |n l⟩ ⊗ [Yₗ ⊗ χ_½]^j, and expand the radial factor in three-dimensional harmonic-oscillator functions Rₙₗ(r; b) with b = √(ħ/m_N ω). That basis choice is also the boundary condition: every Rₙₗ vanishes as r^l at the origin and dies as a Gaussian at infinity, so u(0) = 0 and u(r→∞) = 0 are imposed by construction rather than searched for by shooting. With ħ²/2m_N = 20.7 MeV fm², about fifteen radial functions per block converge the bound spectrum to tens of keV — and the truncation, not the physics, is the first thing to check.

Square ĵ = l̂ + ŝ: ĵ² = l̂² + ŝ² + 2 l̂⋅ŝ, so l̂⋅ŝ = ½(ĵ² − l̂² − ŝ²), and in |n l j mⱼ⟩ it is diagonal with eigenvalue (ħ²/2)[j(j+1) − l(l+1) − 3/4]. For j = l + ½ that is +lħ²/2; for j = l − ½ it is −(l+1)ħ²/2. Nothing else in Ĥ touches spin, so the matrix never connects different (l, j): it is block diagonal, one real symmetric radial block per (l, j), each perhaps 15 × 15, and numpy.linalg.eigh returns the lot in milliseconds. No Clebsch–Gordan coefficient is ever needed, because the coupled basis has already done the coupling, and the 2j + 1 values of mⱼ come out degenerate for free since Ĥ commutes with ĵ. Concretely 1f7/2 carries ⟨l⋅s⟩ = +1.5ħ² and 1f5/2 carries −2ħ²; those two numbers are the same in every nucleus, and only the coefficient Vₗₛ in front of them moves with A.

Subtract the two eigenvalues: ΔE = E(j = l − ½) − E(j = l + ½) = (2l+1)|Vₗₛ|/2. The gap is linear in l, so a p orbital splits by 3|Vₗₛ|/2 while an i orbital splits by 13|Vₗₛ|/2 — more than four times as much for the same coefficient. Where the two levels sit is then fixed by a sum rule: weight by 2j + 1 and note that (2l+2)(l) − (2l)(l+1) = 0, so l̂⋅ŝ leaves the multiplet's centre of gravity exactly where it was. The j = l + ½ level therefore falls by l/(2l+1) of ΔE and the j = l − ½ level rises by (l+1)/(2l+1) of it; the fatter level moves less. For 1g with |Vₗₛ| = 1.4 MeV, ΔE = 9 × 1.4/2 = 6.3 MeV, 1g9/2 drops 4/9 × 6.3 = 2.8 MeV, 1g7/2 climbs 5/9 × 6.3 = 3.5 MeV, and 10(−2.8) + 8(+3.5) = 0 as it must. Use that as an arithmetic check on every level scheme you draw.

The radial shape is Thomas-like: Vₗₛ(r) = −κ (1/r) dV_WS/dr. Woods–Saxon is flat inside, so its derivative is a bell of width a centred on r = R with peak V₀/4a = 50/(4 × 0.55) = 23 MeV fm⁻¹. At A = 90, R = 5.60 fm, and with κ ≈ 0.5 fm² the peak coefficient is 0.5 × 23/5.60 ≈ 2.0 MeV; averaging over a radial wavefunction that is not all in the skin brings it down to roughly 1.4 MeV. Two effects then point the same way — the 2l + 1 factor, and the centrifugal term l(l+1)ħ²/2m_N r² pushing high-l orbitals out into exactly the shell where the force acts. The sign is the part to memorise. Atomic fine structure has a positive coefficient and puts j = l + ½ higher; the nuclear term is negative and puts it lower. It is also twenty to thirty times too strong to be the electron's Thomas precession scaled up: it comes from the nucleon–nucleon interaction, and in a non-relativistic mean field both κ and its sign are fitted, not derived.

Count the plain oscillator first. Shell N = 2(n − 1) + l holds (N + 1)(N + 2) nucleons once spin is included, so the closures are 2, 8, 20, 40, 70, 112: three right, the rest wrong. Flattening the well adds a depression roughly proportional to l(l+1), which spreads each N shell and, in heavy nuclei, sinks its largest-l multiplet toward the shell below, so the gaps drift toward 58, 92 and 138, still not the observed ones. Spin–orbit is what regroups them: it tears the j = l + ½ member of the largest l in shell N away from the rest of that shell and opens a new gap above it. For 1f7/2 that is all that happens, and it stands alone between 20 and the rest of N = 3; from 1g9/2 on, the depression and l⋅s together carry that member down toward shell N − 1, through most or all of the gap ħω ≈ 41 A^(−1/3) MeV, until it sits on top of that shell or, for 1i13/2, among its levels. After that it is only addition. Shell N = 3 loses 1f7/2, which holds 8: 20 + 8 = 28. The remaining 12 of N = 3 take you to 40, and 1g9/2's 10 give 50. The remaining 20 of N = 4 reach 70, and 1h11/2's 12 give 82. The remaining 30 of N = 5 reach 112, and 1i13/2's 14 give 126. Every gap above 20 is made by one split-off orbital; 2, 8 and 20 need none, because the well alone makes them.

Even-even nuclei are 0⁺ in their ground state without exception, and the mean field does not explain that: it is the residual pairing interaction, left out of Ĥ by construction, coupling identical nucleons pairwise to J = 0. Grant it, and odd-A becomes a one-line prediction — J^π is the j of the last odd nucleon with parity (−1)^l for its orbital. ¹⁷O puts a neutron in 1d5/2 and is 5/2⁺; ⁴¹Ca puts one in 1f7/2 and is 7/2⁻; ⁹¹Zr's 51st neutron sits in 2d5/2 and gives 5/2⁺; ²⁰⁷Pb is a 3p1/2 hole in the N = 126 core and is 1/2⁻, a hole carrying the particle's j and parity. Then the failures. Odd-odd nuclei need Nordheim's rules to couple two odd nucleons. Deformed nuclei break spherical symmetry so j stops being a good quantum number: ¹⁹F should be 5/2⁺ and is 1/2⁺, ²³Na should be 5/2⁺ and is 3/2⁺, and both want Nilsson orbitals labelled Ω[N n_z Λ]. And a level scheme fixes which states exist, never how fast anything decays.