Nuclear Physics

Hartree-Fock-Bogoliubov — exercises

Solve the BCS / HFB gap-and-number equations self-consistently for (Δ, λ) given a single-particle spectrum {e_k, Ω_k}, pairing strength G, and target particle number N; read v_k² as the orbital occupation and 2Δ as the excitation gap; identify the critical G below which pairing collapses.

1 worked example · 7 practice problems · 2 check problems

Worked example: two-level BCS at half-filling

Problem. A model nucleus has two single-particle levels at energies $e_{1} = - 1$ and $e_{2} = + 1$ (MeV), each with pair-degeneracy $Ω_{1} = Ω_{2} = 1$ . With pairing strength $G = 2$ MeV and target particle number $N = 1$ (half-filling), find the chemical potential $λ$ , gap $Δ$ , and the orbital occupations $v_{k}^{2}$ . Verify the gap and number equations are simultaneously satisfied.

Diagnosis. Two equations, two unknowns. The level structure has a reflection symmetry ( $e_{k} \to - e_{k}$ , which swaps levels 1 and 2) and the target $N = 1$ is exactly half of the total pair degeneracy $Ω_{1} + Ω_{2} = 2$ . The symmetry forces $λ = 0$ . That leaves just the gap equation, in one unknown $Δ$ .

Predict before reading on: why does the reflection symmetry pin $λ = 0$ ? Think about what the number equation does as $λ$ moves away from zero.

Symmetry argument for $λ = 0$ . Under the reflection $e_{k} \to - e_{k}$ , $λ \to - λ$ , the quasi-particle energies $E_{k} = (e_{k} - λ)^{2} + Δ^{2}$ are invariant when levels are swapped, but the occupation $v_{k}^{2} = \frac{1}{2} [1 - (e_{k} - λ) / E_{k}]$ picks up the substitution $(e_{k} - λ) \to - (e_{k} - λ)$ , which means $v_{1}^{2} \leftrightarrow v_{2}^{2}$ . Total $N = \sum_{k} Ω_{k} v_{k}^{2}$ is invariant only if it's already at the symmetric value $\sum_{k} Ω_{k} /2 = 1$ . So $N = 1$ is the symmetric point, and self-consistency forces $λ = 0$ .

Gap equation. With $λ = 0$ , $E_{1} = E_{2} = 1 + Δ^{2}$ . The gap equation

1 = \frac{G}{2} k \sum \frac{Ω _{k}}{E _{k}} = \frac{G}{2} \cdot \frac{2}{1 + Δ ^{2}} = \frac{G}{1 + Δ ^{2}}

collapses to $1 + Δ^{2} = G$ , so $Δ^{2} = G^{2} - 1 = 1$ and $Δ = 1$ MeV.

Predict before reading on: what would happen if $G < 1$ ? Look at the equation $1 + Δ^{2} = G$ and ask: can it be satisfied?

Occupations. With $Δ = 1, λ = 0$ , $E_{k} = 2$ :

v_{1}^{2} = \frac{1}{2} [1 - \frac{- 1 - 0}{2}] = \frac{1}{2} [1 + \frac{1}{2}] = \frac{1}{2} + \frac{1}{2 2} \approx 0.854

v_{2}^{2} = \frac{1}{2} [1 - \frac{1}{2}] \approx 0.146

Check both equations.

Number: $Ω_{1} v_{1}^{2} + Ω_{2} v_{2}^{2} = 0.854 + 0.146 = 1.000$ ✓
Gap: $(G /2) \sum_{k} Ω_{k} / E_{k} = (2 /2) (2/ 2) = 1$ ✓

Physical reading. The HF answer for $N = 1$ on these two levels would put one particle in level 1 (fully occupied) and zero in level 2 (empty) — a sharp step. HFB at $Δ = 1$ (comparable to the level spacing $2$ ) smears that to $(0.85, 0.15)$ . The 15% of probability that level 2 picks up is the pairing correlation: a fraction of a Cooper pair has "leaked" up to the empty level, lowering the total energy.

Articulate: state in one sentence what $Δ$ physically represents and what units it must have. What experimental signature in nuclear-mass data corresponds to a non-zero $Δ$ ?

Practice problems

Seven problems. The first three vary parameters on the same two-level model to test mechanics and the phase transition. The next two test interpretation (gap, HF limit). The last two extend to richer settings (even-odd staggering, sub-shell closure).

P.1 phase transition at critical G

Same two-level system as the worked example ( $e = \pm 1, Ω_{k} = 1, N = 1$ ). What is the critical pairing strength $G_{c}$ below which $Δ = 0$ is the only solution? Plot or describe how $Δ (G)$ behaves for $G \in [0, 2]$ .

Find the analogue: the worked example solved the gap equation $1 + Δ^{2} = G$ at $G = 2$ . What happens to this equation as $G$ decreases? At what $G$ does it stop having a real solution?

show answer

$1 + Δ^{2} = G \Rightarrow Δ = G^{2} - 1$ when $G \geq 1$ .

For $G < 1$ : no real $Δ > 0$ solution. The only fixed point is $Δ = 0$ — pairing collapses. $G_{c} = 1$ .

Shape of $Δ (G)$ : identically zero on $[0, 1]$ ; for $G > 1$ the gap turns on continuously, with $Δ (G) \to 0^{+}$ as $G \to 1^{+}$ . Square-root singularity: $Δ \approx 2 (G - 1)$ near the critical point.

This is the textbook second-order phase transition. The order parameter $Δ$ plays the role of the magnetization in Ising, with $G - G_{c}$ the analogue of $T_{c} - T$ , and the mean-field exponent $β_{MF} = 1/2$ matches the BCS mean-field result.

P.2 particle-number tuning via λ

Same two-level system, now with target $N = 0.5$ instead of 1 (quarter-filling). What sign does $λ$ need to have, qualitatively? Setting up the equations carefully, find $λ$ and $Δ$ numerically at $G = 2$ . (Hint: $λ$ is no longer pinned to zero; solve the coupled system.)

Find the analogue: the worked example used reflection symmetry to fix $λ = 0$ . Now you've broken that symmetry by shifting $N$ . The number equation becomes a non-trivial constraint that picks $λ$ .

show answer

Sign of $λ$ . $N = 0.5 < 1$ means we want fewer particles than half-filled; $λ$ should sit below the symmetric point so fewer states are occupied. $λ < 0$ , somewhere near or below $e_{1} = - 1$ .

Coupled equations. Define $ξ_{k} = e_{k} - λ$ :

Number: $\frac{1}{2} [1 - ξ_{1} / E_{1}] + \frac{1}{2} [1 - ξ_{2} / E_{2}] = 0.5$ with $E_{k} = ξ_{k}^{2} + Δ^{2}$
Gap: $1 = (G /2) [1/ E_{1} + 1/ E_{2}]$

Numerical solution (e.g., scipy fsolve or nested bisection like the page's solver):

import numpy as np
from scipy.optimize import fsolve

e = np.array([-1.0, 1.0])
G = np.sqrt(2)
N_target = 0.5

def eqs(p):
    lam, Delta = p
    xi = e - lam
    E  = np.sqrt(xi**2 + Delta**2)
    return [0.5*((1 - xi/E).sum()) - N_target,
            0.5*G*(1/E).sum() - 1]

lam, Delta = fsolve(eqs, [-1.0, 0.5])
print(f"λ = {lam:.4f}, Δ = {Delta:.4f}")

Result: $λ \approx - 1.155$ , $Δ \approx 0.667$ . The chemical potential drops just below $e_{1}$ , and the gap is smaller than in the half-filled case because we're now near the bottom of the spectrum where pairing has less "headroom."

P.3 three-level symmetric extension

Extend the worked example to three levels: $e = (- 2, 0, + 2)$ MeV, all with $Ω_{k} = 1$ . Target $N = 1.5$ (half-filling, since total degeneracy is 3). At $G = 0.8$ MeV, find $λ, Δ$ , and the three occupations $v_{k}^{2}$ .

Find the analogue: the same reflection symmetry pins $λ = 0$ at half-filling. The middle level $e_{2} = 0$ sits exactly at $λ$ , so $E_{2} = Δ$ directly — the gap shows up as a "naked" eigenvalue.

show answer

Symmetry pins $λ = 0$ . Gap equation:

$1 = (G /2) [1/ 4 + Δ^{2} + 1/Δ + 1/ 4 + Δ^{2}] = (G /2) [2/ 4 + Δ^{2} + 1/Δ]$

This is transcendental — solve numerically (bisection, scipy, hand iteration). At $G = 0.8$ :

$Δ \approx 0.646$ MeV, $λ = 0$ .

Quasi-particle energies: $E_{1} = E_{3} = 4 + 0.418 \approx 2.102$ , $E_{2} = Δ = 0.646$ .

Occupations: $v_{1}^{2} = 0.976, v_{2}^{2} = 0.500, v_{3}^{2} = 0.024$ . Sum: $1.500$ ✓.

Reading: level 1 (well below $λ$ ) is almost fully occupied; level 3 (well above) is almost empty; level 2 (right at $λ$ ) is exactly half-occupied. The smearing concentrates on levels within $\sim Δ$ of the chemical potential — the page's "step softens into a sigmoid of width $\sim Δ$ ."

P.4 excitation gap from quasi-particle spectrum

For each system below, compute the lowest quasi-particle energy $E_{m i n}$ , which equals the energy cost to create one quasi-particle excitation.

(a) Worked example: $e = \pm 1, Ω = (1, 1), Δ = 1, λ = 0$ .

(b) Worked example with $G = 2$ (so $Δ = 3$ ): compute $E_{m i n}$ .

State the general rule: when does $E_{m i n} = Δ$ exactly, and when is it larger?

Find the analogue: $E_{k} = (e_{k} - λ)^{2} + Δ^{2}$ is the cost of breaking up a Cooper pair to put a quasi-particle in level $k$ . The minimum over $k$ gives the excitation gap.

show answer

(a) $E_{1} = E_{2} = 1 + 1 = 2 \approx 1.414$ . $E_{m i n} = 2$ , not $Δ = 1$ .

(b) $E_{k} = 1 + 3 = 2$ for both. $E_{m i n} = 2$ , larger than $Δ = 3 \approx 1.732$ .

General rule. $E_{m i n} = Δ$ exactly when there exists a single-particle level sitting at the chemical potential, $e_{k} = λ$ . Otherwise $E_{m i n} = (e_{closest} - λ)^{2} + Δ^{2} > Δ$ . So the textbook statement "the BCS gap is $Δ$ " only applies to (effectively) continuous spectra where some level is always near $λ$ . In finite nuclei with discrete spectra, the observed gap can be noticeably larger than $Δ$ , especially for closed-shell-adjacent configurations where no level sits at $λ$ .

The two-quasi-particle excitation (breaking a pair) has energy $2 E_{m i n}$ , which is what the empirical mass formula $Δ \approx 12/ A$ MeV is fitting to. In the page's notation, the "pairing gap from mass-staggering" and the " $Δ$ from the BCS equations" agree only when a level sits at $λ$ .

P.5 HF limit recovery

Verify that as $G \to 0$ (no pairing), the HFB occupations collapse to the HF answer: each level is either fully occupied ( $v_{k}^{2} = 1$ ) or empty ( $v_{k}^{2} = 0$ ), with the boundary set by $λ$ .

(a) Start from $v_{k}^{2} = \frac{1}{2} [1 - (e_{k} - λ) / E_{k}]$ with $E_{k} = (e_{k} - λ)^{2} + Δ^{2}$ . Take $Δ \to 0$ and compute the limit.

(b) Explain what happens at a level exactly equal to $λ$ (e.g., the middle level in P.3 at $G \to 0^{+}$ ). What does HFB say there?

Find the analogue: the page's viz showed $Δ = 0$ recovering the HF step; this problem confirms it algebraically.

show answer

(a) With $Δ = 0$ : $E_{k} = ∣ e_{k} - λ ∣$ . The occupation becomes

$v_{k}^{2} = \frac{1}{2} [1 - (e_{k} - λ) /∣ e_{k} - λ ∣] = \frac{1}{2} [1 - sgn (e_{k} - λ)]$

For $e_{k} < λ$ : $sgn = - 1$ , so $v_{k}^{2} = 1$ (occupied). For $e_{k} > λ$ : $sgn = + 1$ , so $v_{k}^{2} = 0$ (empty). The Heaviside step from the HFB occupation page. ✓

(b) For $e_{k} = λ$ : $(e_{k} - λ) / E_{k} = 0/0$ — indeterminate! Resolving via $lim_{Δ \to 0^{+}} (e_{k} - λ) / (e_{k} - λ)^{2} + Δ^{2}$ with $e_{k} = λ$ gives 0, so $v_{k}^{2} = 1/2$ . HFB "splits the difference" at the chemical potential.

This is exactly the partial-filling situation that plain HF cannot handle: with discrete levels and a target $N$ falling between two configurations, HF has no solution; HFB regularizes it with a half-filled level. Pairing is precisely what fixes the partial-filling singularity of HF, even at infinitesimal $G$ .

P.6 even-odd mass staggering

The empirical pairing gap is extracted from nuclear masses via the three-point mass formula:

$Δ_{n}^{(3)} (N) = \frac{1}{2} [B (N + 1) - 2 B (N) + B (N - 1)]$

where $B (N)$ is the binding energy and the formula is evaluated at an odd $N$ . In an HFB-like model, removing one nucleon from an even- $N$ ground state costs an extra $Δ$ compared to the smooth single-particle change, because the unpaired nucleon can't pair with anyone.

For the worked example's two-level system at $N = 0, 1, 2$ :

(a) Compute the HFB ground-state energy at each $N$ in the constant- $G$ approximation (use $E_{HFB} = \sum_{k} Ω_{k} e_{k} v_{k}^{2} - Δ^{2} / G$ ).

(b) Compute the two-nucleon separation energy $S_{2 n} (2) = E (0) - E (2)$ .

Find the analogue: this is the "S_2n" formula from the HFB concept page applied to the worked example's mini-system. The HFB code from the page does this automatically; here you do it by hand on a 2-level model.

show answer

For $N = 1$ (worked example): $Δ = 1, λ = 0, v_{1}^{2} = 0.854, v_{2}^{2} = 0.146$ . Energy:

$E (1) = e_{1} v_{1}^{2} + e_{2} v_{2}^{2} - Δ^{2} / G = - 0.854 + 0.146 - 1/ 2 = - 0.708 - 0.707 = - 1.415$ MeV.

For $N = 0$ (empty): both $v_{k}^{2} = 0$ , $Δ = 0$ (no particles, no pairing). $E (0) = 0$ .

For $N = 2$ (filled): both $v_{k}^{2} = 1$ , $Δ = 0$ . $E (2) = e_{1} + e_{2} = - 1 + 1 = 0$ .

(b) $S_{2 n} (2) = E (0) - E (2) = 0 - 0 = 0$ MeV — removing both particles from the fully-filled state costs nothing in this symmetric toy model (the single-particle energies cancel).

Hmm — so in this toy model the most bound state is the half-filled one (with pairing). Adding or removing particles costs +1.415 MeV. That's the "pairing energy" the system pays per pair to organize at half-filling. In the full sd-shell numbers in the HFB code, this becomes $\sim 1.6$ MeV at $G = 0.4$ , matching the empirical $Δ \approx 12/ A \sim 2$ MeV for $A \sim 30$ .

P.7 sd-shell sub-shell closure

Look at the HFB concept page's sd-shell results. At $N = 8$ (full $0 d_{5/2}$ + $1 s_{1/2}$ ), the table shows $Δ = 0$ all the way up to $G = 0.40$ MeV, and only at $G = 0.60$ does pairing turn on. Explain in your own words:

(a) Why does pairing collapse at this particular $N$ ?

(b) What does the critical $G$ for sub-shell-closure pairing depend on?

Find the analogue: P.1 derived the critical $G$ for a symmetric two-level system in closed form. The sd-shell is a richer setting, but the same logic applies: $G_{c}$ is set by the level spacing across the chemical potential.

show answer

(a) At $N = 8$ , the levels $0 d_{5/2}$ (e = −3.93 MeV) and $1 s_{1/2}$ (e = −3.16) are fully occupied, and the next available shell $0 d_{3/2}$ sits at +2.11 — a 5.3 MeV gap to the next empty state. For pairing to turn on, $G$ has to overcome this large single-particle gap.

(b) The critical $G_{c}$ for a sub-shell closure depends on the energy spacing between the highest filled level and the lowest empty level. In rough terms, $G_{c} \sim$ (spacing across the gap)/(typical degeneracy near $λ$ ). A small single-particle gap makes pairing easy to turn on (small $G_{c}$ ); a large gap suppresses pairing (large $G_{c}$ ). This is why magic numbers resist pairing — they correspond to single-particle gaps of order 5-10 MeV, well above any realistic $G$ .

(c) At $N = 4$ , the $0 d_{5/2}$ shell is partially filled (4 of 6 m-states). There's no single-particle gap to overcome — the chemical potential sits in the middle of the partially-occupied shell, with empty states immediately available at essentially the same energy. Any positive $G$ produces a non-zero $Δ$ because the gap equation $1 = (G /2) \sum Ω_{k} / E_{k}$ has a divergent term (one $E_{k} \to 0$ ) at $Δ = 0$ , so it's always satisfiable for any $G > 0$ . This is the same phenomenon as the partial-filling discussion in P.5(b): pairing regularizes what would otherwise be a singularity.

The general picture: closed-shell nuclei (magic numbers) don't need HFB — HF is enough, pairing collapses. Open-shell nuclei need HFB — pairing is essential. This is the empirical rule "pair the open shells" formalized.

Check problems

Two problems that don't pattern-match against the practice set. The first derives a clean closed-form result; the second tests understanding of the deeper structural cost HFB pays.

Check 1 derivation

Derive a closed-form expression for the critical pairing strength $G_{c}$ of a symmetric two-level system with arbitrary level spacing $2 a$ (levels at $\pm a$ ) and arbitrary pair-degeneracies $Ω_{1} = Ω_{2} = Ω$ , at half-filling $N = Ω$ .

(a) Set up the gap equation explicitly using symmetry to fix $λ$ .

(b) Take the $Δ \to 0^{+}$ limit and solve for the critical $G_{c}$ .

(c) Sanity-check against the worked example ( $a = 1, Ω = 1$ gives $G_{c} = 1$ ) and against the limit $a \to 0$ (degenerate levels: pairing turns on at any $G > 0$ ).

show solution sketch

(a) Symmetry pins $λ = 0$ at half-filling. Both levels have $E_{k} = a^{2} + Δ^{2}$ . Gap equation:

$1 = (G /2) \sum_{k} Ω_{k} / E_{k} = (G /2) \cdot 2Ω/ a^{2} + Δ^{2} = G Ω/ a^{2} + Δ^{2}$

(b) Take $Δ \to 0^{+}$ :

$1 = G_{c} Ω/ a \Rightarrow G_{c} = a /Ω$

For $G > G_{c}$ : the gap equation has a real $Δ > 0$ solution $Δ = (G Ω)^{2} - a^{2}$ . For $G < G_{c}$ : only the trivial $Δ = 0$ solution exists.

Check 2: $a \to 0$ (degenerate levels): $G_{c} \to 0$ . Any positive $G$ creates a gap when the single-particle levels are degenerate — there's no spacing to overcome. ✓ This is why pairing is so robust in nearly-degenerate j-shells: $G_{c} \to 0$ as the level degeneracy increases, so pairing turns on for any non-zero interaction.

The result $G_{c} = a /Ω$ reads: larger spacings need more interaction to overcome; richer degeneracies need less. Both effects are real in nuclear physics — high- $j$ shells (large $Ω$ ) are easy to pair; large single-particle gaps (magic numbers) are hard to pair. The formula quantifies the trade-off.

Check 2 articulation

The HFB ground state is not an eigenstate of particle number $\hat{N}$ . The expected value $⟨ \hat{N} ⟩$ is pinned to $A$ by the Lagrange multiplier $λ$ , but $⟨ \hat{N}^{2} ⟩ - ⟨ \hat{N} ⟩^{2} > 0$ .

In 150-250 words, explain in concrete terms:

What does it mean for the HFB ground state to have a variance in particle number?
What goes physically wrong if you try to compute, say, a beta-decay matrix element with this state without first projecting onto definite $N$ ?
Why does the HFB community accept this cost, given that exact particle-number conservation is "just" a symmetry?

show solution sketch

An HFB ground state is a coherent superposition of components with different particle numbers — schematically, $∣Φ ⟩ = \sum_{N^{'}} c_{N^{'}} ∣ Ψ_{N^{'}} ⟩$ where $∣ c_{N^{'}} ∣^{2}$ peaks at the target $N$ but is non-zero for $N \pm 2, N \pm 4, \dots$ — Cooper pairs being added and removed in the same wavefunction. The variance $⟨ \hat{N}^{2} ⟩ - ⟨ \hat{N} ⟩^{2} > 0$ measures the width of this distribution. For a typical paired open-shell nucleus, the std deviation in $N$ is a few — small compared to $A \sim 100$ , but not zero.

When you compute a one-body matrix element like $⟨ Φ^{'} ∣ \hat{O} ∣Φ ⟩$ between HFB states for the initial and final nucleus in beta decay, the matrix element picks up cross terms between $N^{'} \neq = N$ components — spurious transitions that don't correspond to a real beta decay. The computed matrix element is contaminated by amplitudes for adding-or-removing pairs, which physically can't accompany a single-nucleon beta-decay process. Without particle-number projection, your half-life estimate will be systematically wrong, sometimes by orders of magnitude.

The HFB community accepts the cost because (i) for energetics and densities, the variance is small enough to be a quantitative rather than qualitative error — masses, radii, and binding energies come out within experimental precision; (ii) particle-number projection after variation (PAV) is a one-shot fix, cheap to apply when you need it; and (iii) full variation after projection (VAP) is sometimes done for matrix elements where the bias matters. The trade — exact $N$ for the pairing field — is the right one for the vast majority of HFB observables, and the projections are tools you reach for only when an exact- $N$ answer is required.