Nuclear foundations, primordial nucleosynthesis and LiBeB

What a nucleus is, what it means to fuse one

An atom is made of an extremely dense nucleus at the center and a cloud of electrons around it. The sizes involved are not schoolbook sizes: if the atom were as large as a football pitch, the nucleus would be a pinhead at the center of the center circle. And yet that pinhead contains practically all the mass: the density of nuclear matter is of order $2{.}3 \times 10^{17}\,\mathrm{kg/m^3}$, and a teaspoon of nuclear matter weighs as much as a terrestrial mountain. The nucleus is composed of protons, positively charged, and neutrons, electrically neutral; together we call them nucleons. The number of protons $Z$ establishes which element it is — $Z=1$ is hydrogen, $Z=6$ carbon, $Z=26$ iron, $Z=79$ gold — and the number of neutrons $N$ establishes which isotope of that element one is facing. The sum $A = Z + N$ is the mass number, and each combination $(Z, N)$ defines a nuclide: a specific nuclear species. Carbon-12 and carbon-14 are two different nuclides of the same element.

To fuse two nuclei means bringing them to nuclear distances, of order one femtometer, so that the strong nuclear force can win over the Coulomb repulsion and form a compound state or a bound product. This is difficult because two positively charged nuclei repel each other with increasing energy as the distance decreases: to make fusion probable one needs temperatures of millions or billions of kelvin, or densities and times large enough to exploit the quantum tail of the energy distribution. These conditions are found in stellar centers, in thermonuclear and core-collapse explosions, on the surfaces of accreting compact objects and, before all that, in the first cosmological minutes. The fact that we observe the products of this physics — the atoms we are made of — is the sign that, somewhere and at some time, those conditions were reached. The history of nucleosynthesis is the history of where, when and with which rates.

The chart of nuclides

All known nuclides, stable and unstable, can be represented on a graph with the number of neutrons $N$ on the abscissa and the number of protons $Z$ on the ordinate. This chart of nuclides — the nuclear equivalent of the chemical periodic table, but much more extended — is the geographic map of nucleosynthesis, and we shall meet it constantly: every process (s, r, p, rp) traces a characteristic path on this chart. The chart today contains more than 3,300 experimentally known nuclides, of which only about 250 are stable or quasi-stable on cosmological timescales; the others decay with timescales from microseconds to billions of years [Brookhaven National Laboratory] .

The stable nuclides cluster along a curved line called the valley of stability, which for low $Z$ coincides with $N \approx Z$ and for high $Z$ bends toward $N > Z$, reflecting the fact that heavy nuclei need an excess of neutrons to neutralize the growing Coulomb repulsion among the protons. Far from the valley, in the direction of increasing $N$, one reaches the neutron drip line: beyond this line the addition of a neutron does not produce a bound state, because the nucleus spontaneously emits a neutron. The proton drip line, symmetric on the other side, is reached more quickly because of the Coulomb repulsion. Between the two drip lines lies the region of bound nuclides, and it is within it that all nucleosynthesis unfolds. The r-process runs along the neutron-rich region of the chart, the rp-process along the proton-rich one; the s-process flows right next to the valley of stability.

Binding energy

A nucleus is a bound state of $Z$ protons and $N$ neutrons: like every bound state in quantum mechanics, it has a mass lower than the sum of the masses of its separated constituents. The difference, converted into energy according to $E = mc^2$, is the binding energy

$$B(Z,N) = \left[ Z m_p + N m_n - m(Z,N) \right] c^2.$$

This formula uses nuclear masses. In practice, however, the precision tables such as the AME almost always report atomic masses, that is, masses of the nucleus plus the bound electrons of the neutral atom. In the Q-values of many nuclear reactions the atomic electrons cancel automatically if the total number of neutral electrons is balanced on the two sides; in $\beta^+$ decays, in electron capture and in reactions with a net change of charge one must instead use atomic masses and electronic thresholds with care. This bookkeeping detail is a common source of numerical errors in first network calculations.

For $^{4}\mathrm{He}$, for example, $B \approx 28{.}3$ MeV — the energetic “weld” of the nucleus, equal to about 7 MeV per nucleon. This figure, the binding energy per nucleon $B/A$, is the key quantity of nucleosynthesis: it indicates how stable a nucleus is with respect to its neighbors. If we plot it as a function of $A$ we obtain a curve universally known as the Aston curve: it starts from zero at $A=1$ (hydrogen, which has a single nucleon bound to nothing), rises rapidly to a peak of about 8.8 MeV per nucleon around $A = 56$ ($^{56}\mathrm{Fe}$ is almost the maximum, $^{62}\mathrm{Ni}$ holds the record by a few keV), and then descends slowly to $\sim 7{.}5$ MeV/nucleon for the heaviest elements such as $^{238}\mathrm{U}$.

The consequence is immediate and dominates the entire discipline: the fusion of light nuclei toward the iron group is in general exoergic, it releases energy; the fission of very heavy nuclei is also exoergic; but at the iron group the thermonuclear “mine” runs out. A reaction that starts from nuclei around iron and produces a single, heavier one cannot be a net source of stellar energy: it must be powered from outside. This is why massive stars, after building a core dominated by iron-group nuclei, stop, and stop catastrophically — it is precisely the inability to extract further energy by fusion that contributes to the collapse of the core and the subsequent supernova explosion.

The semi-empirical Bethe-Weizsäcker formula

An approximate description of $B(Z,N)$, sufficient for qualitative reasoning over most of the chart of nuclides, is the semi-empirical mass formula proposed by Carl von Weizsäcker in 1935 and refined by Hans Bethe in the following years. The idea is to treat the nucleus as a charged liquid drop and to sum distinct physical contributions:

$$B(Z,N) = a_V A - a_S A^{2/3} - a_C \frac{Z(Z-1)}{A^{1/3}} - a_A \frac{(N-Z)^2}{A} + \delta(A,Z)$$

with typical coefficients $a_V \approx 15{.}8$ MeV, $a_S \approx 17{.}8$ MeV, $a_C \approx 0{.}71$ MeV, $a_A \approx 23{.}7$ MeV. Each term tells a part of the story. The first, the volume term, says that each added nucleon contributes the same energy “gain” because it interacts only with its immediate neighbors (the nuclear force is short-range). The second, the surface term, is a negative correction: the nucleons at the surface of the nucleus have fewer neighbors, and are therefore “less bound”; it weighs relatively less for large nuclei. The third, the Coulomb term, is the electric repulsion among the $Z$ protons: it grows quickly with $Z^2$, and is the reason why heavy nuclei need more neutrons to compensate it. The fourth, the asymmetry term, is a purely quantum-statistical effect — the Pauli principle penalizes too large an excess of $N$ over $Z$ — and identifies for each $A$ an optimal composition $(N, Z)$ which is precisely the valley of stability. The last, the pairing term, is positive for nuclei with $N$ and $Z$ both even, roughly zero for nuclei with odd $A$, negative for odd-odd nuclei: it makes particularly bound nuclei such as $^{4}\mathrm{He}$, $^{12}\mathrm{C}$, $^{16}\mathrm{O}$, $^{20}\mathrm{Ne}$, $^{24}\mathrm{Mg}$ — the “alpha-multiples” that dominate the nucleosynthesis of the advanced phases.

The Bethe-Weizsäcker formula is not exact — it has typical errors of a few MeV on individual nuclides — but it captures the essential, and from its derivatives many properties can be read off directly: the valley of stability is found by imposing $\partial B / \partial Z |_A = 0$, the spontaneous-fission barrier is estimated from the balance between the Coulomb and surface terms, and so on. For quantitative nucleosynthesis, however, one needs masses measured nuclide by nuclide: the reference compilation is the Atomic Mass Evaluation [Wang et al. 2021] , updated periodically with typical uncertainties of a few keV for stable nuclei and up to hundreds of keV for nuclei far from the valley.

Magic numbers and the shell model

There is a fact that the liquid-drop formula does not explain: some nuclei are much more bound and stable than the formula would predict, at specific values of $Z$ or $N$. These values are the magic numbers:

$$2,\ 8,\ 20,\ 28,\ 50,\ 82,\ 126.$$

Nuclei with magic $Z$, or magic $N$, or both (the “doubly magic” ones such as $^{4}\mathrm{He}$, $^{16}\mathrm{O}$, $^{40}\mathrm{Ca}$, $^{48}\mathrm{Ca}$, $^{208}\mathrm{Pb}$) are particularly bound, have anomalously small neutron-capture cross sections, and accumulate in the nucleosynthesis chains. The explanation comes from the shell model proposed independently by Maria Goeppert Mayer and Hans Jensen at the end of the 1940s (Nobel Prize in physics, shared with Wigner, in 1963): the nucleons occupy quantized energy levels very similar to those of atomic electrons, and the magic numbers correspond to closed shells — configurations in which a shell is completely filled and the next one is separated by an energy gap.

The consequences for the cosmic abundance curve are direct and visible to the naked eye. The two abundance peaks at $A \approx 80$ and $A \approx 130$ correspond to the r-process passing through the magic nuclei at $N = 50$ and $N = 82$ respectively; the two peaks at $A \approx 90$ and $A \approx 138$, slightly shifted to the right, correspond to the s-process passing through the same magic numbers but with different $Z$, because the s-process flows along the valley of stability while the r runs far from it, and the $\beta$ decays following neutron capture bring the isotope back toward the valley at higher $Z$. The peak at $A \approx 208$ for lead corresponds to $N = 126$, the last experimentally known magic number. The shell structure of the nucleus is written directly in the periodic table of cosmic abundances, and B²FH were the first to recognize it precisely through this signature.

Precision masses for exotic nuclei, far from the valley, are critical for two specific groups: the neutron-rich nuclei along the r-process path, where the masses set the Q-values of the $\beta$ decays and hence the position of the abundance peaks, and the proton-rich nuclei along the rp-process, where the masses set the Q-values of proton capture and the characteristic “waiting points”. Precision experiments on exotic nuclei are carried out with Penning traps (ISOLTRAP at CERN-ISOLDE, JYFLTRAP at Jyväskylä, LEBIT at NSCL/FRIB) and with time-of-flight mass spectrometry in storage rings; the reference databases are NNDC [Brookhaven National Laboratory] , ENSDF, and the Nuclear Data Section of the IAEA [International Atomic Energy Agency] .

Decays

When a nucleus is not the most stable among its neighbors on the chart of nuclides — it is not at the “optimal position” for that combination of $A$ — it decays: it spontaneously transforms into another, more stable nuclide, emitting particles or radiation. Decays are just as important as fusion reactions for nucleosynthesis, because they determine where the produced nuclei end up, on which timescale, and with which sign of transformation $Z \to Z' \neq Z$.

$\beta^-$ decay

A nucleus with an excess of neutrons can convert a neutron into a proton, emitting an electron and an antineutrino:

$$n \to p + e^{-} + \bar{\nu}_e \qquad \Longrightarrow \qquad (Z, N) \to (Z+1, N-1).$$

The nuclide moves on the chart one cell up and one to the left, and $A$ remains unchanged. The $\beta^-$ decay is not the inverse reaction of neutron capture — the microscopic inverse of $(n,\gamma)$ is $(\gamma,n)$ — but it is the mechanism that converts neutrons into protons at nearly constant $A$ after a sequence of captures. After a nuclide has captured one or more neutrons, moving far from the valley of stability, a sequence of $\beta^-$ decays brings it back toward the valley. This is precisely the dynamics of the r-process and of the s-process. The $\beta^-$ half-lives vary enormously across the chart of nuclides: from milliseconds for extremely neutron-rich nuclei up to billions of years for nuclei close to the valley ($^{40}\mathrm{K}$ has $t_{1/2} \approx 1{.}25 \times 10^{9}$ years).

A subtlety that becomes crucial in stellar environments: laboratory half-lives are measured for neutral atoms, with their retinue of atomic electrons, while in a hot plasma the atoms can be partially or completely ionized and the nuclei can find themselves in thermally populated excited states. For heavy nuclei this can change the weak rate significantly: ionization can suppress bound-electron-capture channels, open or amplify $\beta^-$ decay into bound electronic states, and the electron density of the plasma can make capture from free electrons dominant. The rates “in stellar plasma” are compiled by Takahashi-Yokoi for evolved stars [Takahashi & Yokoi 1987] and by Langanke-Martínez-Pinedo for supernovae [Langanke & Martínez-Pinedo 2003] , and must be used instead of the laboratory rates when temperature and density require it.

$\beta^+$ decay and electron capture

Specularly, a nucleus with an excess of protons can convert a bound proton into a neutron in two alternative ways. It can emit a positron and a neutrino:

$$(Z,N) \to (Z-1,N+1) + e^{+} + \nu_e.$$

Or it can capture an atomic electron (typically from an inner orbital):

$$(Z,N) + e^{-} \to (Z-1,N+1) + \nu_e.$$

The free proton cannot decay in this way: the transition is possible only inside a nucleus if the nuclear mass difference supplies the necessary energy. Electron capture competes with $\beta^+$ emission, and becomes the only channel when the atomic mass difference is below $2 m_e c^2 = 1{.}022$ MeV — the threshold below which $\beta^+$ emission is kinematically impossible. In stellar environments with high electron densities (such as the core of a massive star in the final phases, where the matter is degenerate), electron capture is drastically amplified: nuclei stable in the laboratory can become unstable in the plasma. This phenomenon — neutronization — is central to the collapse of the iron core preceding a supernova: massive electron capture on protons and on iron-group nuclei reduces $Y_e$, removes electron degeneracy pressure and accelerates the collapse [Langanke & Martínez-Pinedo 2003] .

$\alpha$ decay

The heaviest nuclei, beyond $A \approx 150$, can decay by emitting a $^{4}\mathrm{He}$ nucleus — an alpha particle:

$$(Z, N) \to (Z-2, N-2) + {}^{4}\mathrm{He}.$$

The mechanism is once again the tunnel effect, applied this time not to a collision but to the escape of a preformed $\alpha$ cluster from the nuclear well. The $\alpha$ half-lives range from a fraction of a second (for the most extreme isotopes) to billions of years (for $^{238}\mathrm{U}$, $t_{1/2} \approx 4{.}47 \times 10^{9}$ years) and depend extraordinarily sensitively on the available energy — the Geiger-Nuttall law predicts an exponential dependence between $\log t_{1/2}$ and $1/\sqrt{Q_\alpha}$. For nucleosynthesis, $\alpha$ decay and spontaneous fission contribute to setting the practical upper limit of the natural periodic table: beyond lead and bismuth there are no stable nuclides, and beyond uranium the half-lives are in general so short that the transuranic elements appear in nature only in traces produced by decays, neutron captures or extremely rare primordial residues. Spontaneous fission, active for the heaviest nuclei, completes the picture: even if the r-process can in principle produce nuclei beyond $A \approx 270$, these dissolve by fission on short timescales, re-emitting intermediate-mass fragments and contributing through fission cycling loops to a part of the r abundance.

Q-values and energy budgets

In every reaction or decay, the quantity that fixes its energetics is the Q-value, defined as

$$Q = \left( \sum_{i\,\text{initial}} m_i - \sum_{f\,\text{final}} m_f \right) c^2 = E_{\text{final kinetic}} - E_{\text{initial kinetic}}.$$

A reaction with $Q > 0$ is exoergic: it releases energy. A reaction with $Q < 0$ is endoergic: it absorbs energy, and therefore has a minimum threshold in the kinetic energy of the initial system in order to occur. Most of the fusions that carry light nuclei toward the iron group have $Q > 0$; photodisintegrations instead have a negative Q-value with respect to the direct capture, because the photon must supply the separation energy. The Q-value is the quantity most frequently used in nucleosynthesis calculations, and must be computed from the AME masses — not from the Bethe-Weizsäcker formula, except in qualitative cases. In $\alpha$ decays the Q-value fixes almost all the total kinetic energy of the products; in $\beta$ decays it is instead shared among the electron or positron, the neutrino and the nuclear recoil, and what is measured in the electron spectrum is an endpoint, not a monoenergetic line. This distinction provides an immediate consistency check between models and measurements.

Coulomb barrier and quantum tunneling

When two nuclei approach each other, a “hill” of electric repulsion rises between them. To fuse them classically one would need enough energy to clear the top. But stars are lukewarm compared with that summit: at 15 million kelvin at the center of the Sun, practically no pair of protons has enough energy to climb over it. Classical physics said: stellar fusion cannot work. It is a historically beautiful case: quantum mechanics comes to the rescue, and at the microscopic level a particle can pass through a barrier even without the energy to climb over it. It is the tunnel effect, discovered by Gamow in 1928 (see chapter 1). The tunneling probability is extremely low but not zero, and grows quickly with energy. Without this effect stars could not burn; and since stars do burn, we know that quantum mechanics must be right.

Quantitatively, the Coulomb barrier between two nuclei of charge $Z_1$ and $Z_2$ has height

$$E_C = \frac{Z_1 Z_2 e^2}{R}, \qquad R \approx 1{.}2 \, (A_1^{1/3} + A_2^{1/3}) \,\mathrm{fm}$$

where $R$ is the sum of the nuclear radii, expressed through the well-known relation $R = r_0 A^{1/3}$ with $r_0 \approx 1{.}2$ fm, and $e^2$ is written in the usual nuclear units, that is, it includes the Coulomb factor and equals $e^2 \simeq 1{.}44\,\mathrm{MeV\,fm}$. For p+p, $E_C \approx 0{.}55$ MeV; at the center of the Sun the thermal energy $kT \approx 1{.}3$ keV is about four hundred times lower. The tunneling probability through a Coulomb barrier is approximated by the Gamow factor:

$$P_G(E) = \exp\!\left[ - 2 \pi \eta(E) \right], \qquad \eta(E) = \frac{Z_1 Z_2 e^2}{\hbar v} = \alpha Z_1 Z_2 \sqrt{\frac{\mu c^2}{2E}}$$

where $\eta$ is the Sommerfeld parameter, $\mu = m_1 m_2 / (m_1 + m_2)$ the reduced mass, $v$ the relative velocity, $\alpha \approx 1/137$ the fine-structure constant. The analytic form of the Gamow factor is a WKB (Wentzel-Kramers-Brillouin) solution of the pure Coulomb scattering problem, valid in the low-energy regime $\eta \gg 1$ — exactly the regime of astrophysical reactions. For higher energies the transmission approaches unity, but higher energies are rarely relevant for the star, because the thermal distribution of the particles exponentially penalizes the high-energy tails.

Corrections to pure tunneling

The analytic form just written is valid for a two-body system in vacuum, in $s$-wave, with a pure Coulomb barrier. In a real environment there are corrections that must be taken into account if one wants percent precision or better. The first correction is the centrifugal barrier: for partial waves with angular momentum $\ell > 0$, a term $\ell(\ell+1) \hbar^2 / (2 \mu R^2)$ adds to the Coulomb repulsion, raising the effective barrier and suppressing the tunneling. For most stellar fusion reactions the dominant contribution is the $s$-wave one ($\ell = 0$), but for some reactions — such as the radiative capture $^{12}\mathrm{C}(\alpha,\gamma)^{16}\mathrm{O}$ — the $\ell > 0$ contributions are significant.

The second important correction is electron screening. In a plasma or in a laboratory target, the interacting nuclei are not isolated: the surrounding electrons envelop them and partially screen their charge, reducing the effective height of the Coulomb barrier and increasing the reaction rate with respect to the vacuum value. In the stellar environment the effect is quantified by the Salpeter formula (weakly coupled regime) for non-degenerate plasmas, and by the corrections of Mitler, Itoh and others for the intermediate and strong coupling regimes typical of degenerate stars or of the dense cores of pre-Ia white dwarfs. The rate amplification from stellar screening is generally 5-15%, but can reach 50% in extreme conditions; for low-charge reaction rates like the pp it is negligible, for high-charge reactions like $^{12}\mathrm{C} + {}^{12}\mathrm{C}$ it is relevant and must be incorporated in the codes. In the laboratory, screening manifests differently and in general more strongly — the electrons of the target amplify the low-energy cross section by factors that at $E \lesssim 10$ keV can exceed 2 — and must be carefully subtracted to extract the “bare” $S(E)$ factor relevant for astrophysics. It is one of the reasons why the LUNA measurements at Gran Sasso are conducted with gas targets and carefully characterized screening-subtraction procedures.

Cross section and astrophysical factor

The quantity that measures how probable it is for two nuclei to react is called the cross section $\sigma$ and has units of area: it conceptually represents the “size” of the target seen by the projectile. If you fire $N_p$ projectiles per second at a thin target containing $n_t$ nuclei per unit area, the number of reactions per second is $N_p \, n_t \, \sigma$. The traditional unit of the nuclear cross section is the barn = $10^{-24}\,\mathrm{cm^2}$; astrophysical reactions at low energies typically have cross sections of nanobarns ($10^{-33}\,\mathrm{cm^2}$) or even less.

For stellar reactions, $\sigma$ is a strongly energy-dependent function: at very low energy it is dominated by the penetration of the Coulomb barrier, at higher energies the nuclear structure and the dynamical effects of the reaction come into play. Nuclear physicists “remove” the Gamow part to obtain a smoother, more easily extrapolated function, called the astrophysical $S(E)$ factor:

$$\sigma(E) = \frac{S(E)}{E} \exp\!\left[ -2\pi\eta(E) \right].$$

The $S(E)$ factor has units of energy times area and can vary by many orders of magnitude from one reaction to another; its usefulness is not in being “large” or “small”, but in varying much more slowly than $\sigma(E)$ when the dominant behavior is Coulomb penetration. For non-resonant reactions, $S(E)$ is therefore extrapolatable from the energies measured in the laboratory (tens or hundreds of keV) to the astrophysical Gamow peak (a few keV or tens of keV) with relatively small uncertainties if sufficient data are available [Adelberger et al. 2011] . For resonant reactions, $S(E)$ has peaks at the positions of the resonances and the extrapolation requires more sophisticated tools (see the next section).

The thermal rate and the Gamow peak

A star does not see a single pair of nuclei at a definite energy: it sees a plasma with a thermal Maxwell-Boltzmann distribution of velocities. The mean reaction rate per pair, as a function of temperature, is

$$\langle \sigma v \rangle = \left( \frac{8}{\pi \mu} \right)^{1/2} (kT)^{-3/2} \int_0^\infty E \, \sigma(E) \, e^{-E/kT} \, dE.$$

Substituting the expression for $\sigma(E)$ through the astrophysical factor, the integrand contains the product

$$e^{-E/kT} \cdot e^{-2\pi\eta(E)}$$

— a thermal Maxwell-Boltzmann term multiplied by the exponential tunneling factor. The thermal term decreases toward high energies; the tunneling factor is almost zero at low energies and grows rapidly to the right. Their product therefore has a sharp maximum at an energy $E_0$ well above $kT$ but well below $E_C$ — the Gamow peak:

$$E_0 = \left[ \frac{(\pi \alpha Z_1 Z_2)^2 \, \mu c^2 \, (kT)^2}{2} \right]^{1/3} = \left( \frac{E_G (kT)^2}{4} \right)^{1/3}, \qquad E_G = 2\mu c^2(\pi\alpha Z_1Z_2)^2.$$

For p+p fusion in the Sun at $T = 1{.}5 \times 10^{7}$ K, the Gamow peak is at $E_0 \approx 6$ keV; for the $\alpha+\alpha$ entrance feeding the $3\alpha$ in red giants at $T = 10^{8}$ K, it is at about 80 keV; for $^{12}\mathrm{C}(\alpha,\gamma)^{16}\mathrm{O}$ in helium burning at $T \sim 2 \times 10^{8}$ K, it falls around 300 keV; for silicon burning in pre-supernovae at $T \approx 3 \times 10^{9}$ K, $E_0$ is of order one MeV. It is the energy around which most reactions actually occur: the plasma selects by itself its nuclear working band. The width of the peak is $\Delta = 4 \sqrt{E_0 \, kT / 3}$, generally narrower than the central energy, and in most practical cases the convolution is evaluated in the saddle-point approximation, giving

$$\langle \sigma v \rangle \approx \left( \frac{2}{\mu} \right)^{1/2} \frac{\Delta}{(kT)^{3/2}} S(E_0) \, e^{-3 E_0/kT}.$$

The exponent $\propto E_0/kT \propto T^{-1/3}$ produces the very strong temperature dependence of the rates typical of thermonuclear nucleosynthesis. Around the relevant stellar conditions, local approximations such as $\epsilon_{\mathrm{pp}} \propto T^{4}$, $\epsilon_{\mathrm{CNO}} \propto T^{16-20}$, $\epsilon_{3\alpha} \propto T^{30-40}$ are often used: they are not universal laws, but they describe the sensitivity well in the temperature interval in which those processes operate. These extreme dependences are what makes stars extremely sensitive “thermostats”: a small variation in temperature produces an enormous variation in energy production, and this is the deep mechanism of the nuclear stability of stars in hydrostatic equilibrium.

Standard compilations of $\langle \sigma v \rangle(T)$ in polynomial form (REACLIB parameterizations with 7 coefficients) are provided by NACRE-II [Xu et al. 2013] for the reactions of the main chains and by the JINA REACLIB library [Cyburt et al. 2010] for the large network of explosive nucleosynthesis. The discussion of the solar reactions, in particular $S(0)$ and the extrapolation to the solar Gamow peak, is crystallized in the Solar Fusion Cross Sections II [Adelberger et al. 2011] , the standard reference for the rates of the pp chain and the CNO cycle at the level currently reachable by laboratory data.

Resonances and non-resonant reactions

Not all nuclear reactions have a smooth, gently curved $S(E)$ factor. Some occur with much higher probability at a specific energy: one says they have a resonance at that energy. It is a phenomenon analogous to the glass that shatters when a singer hits the right note, or to a swing pushed at the exact moment of its proper period. A nuclear resonance corresponds to the formation of an excited state of the compound nucleus during the reaction: the two colliding nuclei are temporarily “welded” into an intermediate nucleus in a quasi-bound state, which lives for a time equal to $\hbar/\Gamma$ — where $\Gamma$ is the width of the state — before disintegrating into one of the open channels.

Resonances can radically change the picture. We have already seen the paradigmatic example in the historical chapter: the $3\alpha$, which without the resonance of the Hoyle state at 7.65 MeV in $^{12}\mathrm{C}$ would be astrophysically irrelevant, becomes with it the engine of carbon production in the universe. Without that single resonance, carbon nucleosynthesis would be lower by some ten orders of magnitude.

The Breit-Wigner form

An isolated resonance of width $\Gamma$ and energy $E_R$ contributes to the cross section according to the Breit-Wigner formula:

$$\sigma_{BW}(E) = \pi \lambdabar^{2} \, \omega \, \frac{\Gamma_a \, \Gamma_b}{(E - E_R)^2 + (\Gamma/2)^2}$$

with $\lambdabar = \hbar / \sqrt{2 \mu E}$ the reduced de Broglie wavelength, $\Gamma = \Gamma_a + \Gamma_b + \dots$ the sum of the partial widths of all open channels, and $\omega = (2J+1)/[(2J_1+1)(2J_2+1)]$ a statistical spin factor (with $J$ the spin of the resonance and $J_{1,2}$ the spins of the reactants). The contribution of the resonance to the thermal rate is dominated by the region of energies around $E_R$, and after integration over the Maxwellian distribution it is written in terms of the resonance strength:

$$\omega \gamma \equiv \omega \, \frac{\Gamma_a \Gamma_b}{\Gamma}$$

which concentrates in a single laboratory-measurable number all the astrophysical information of the resonance. The thermal rate is then

$$\langle \sigma v \rangle_R = \left( \frac{2 \pi}{\mu kT} \right)^{3/2} \hbar^2 \, (\omega\gamma) \, e^{-E_R/kT}$$

for an isolated, narrow resonance ($\Gamma \ll E_R, kT$). A resonance at energy $E_R \sim E_0$ — at the Gamow peak — can completely dominate the rate, suppressing the non-resonant contribution by many orders of magnitude. Key reactions in which resonances, subthreshold tails or interferences modify the rate include $^{14}\mathrm{N}(p,\gamma)^{15}\mathrm{O}$ (the slow step of the CNO cycle, measured by LUNA), $^{13}\mathrm{C}(\alpha,n)^{16}\mathrm{O}$ and $^{22}\mathrm{Ne}(\alpha,n)^{25}\mathrm{Mg}$ (neutron sources for the s-process), $^{12}\mathrm{C}(\alpha,\gamma)^{16}\mathrm{O}$ (the C/O ratio after helium burning) and many others. Each subsequent chapter of this book will present new ones.

Subthreshold resonances and R-matrix

There is a subtlety that deserves to be pointed out, because it is at the origin of many of the residual uncertainties in astrophysical rates: the subthreshold resonances. They are excited states of the compound nucleus that lie below the threshold energy of a given entrance channel, but which still have a “tail” extending above the threshold by virtue of their intrinsic width. A subthreshold resonance does not correspond to a true resonance in the laboratory — the cross section shows no sharp peak — but contributes significantly to $S(E)$ at low energies, deforming its extrapolated shape. The $^{12}\mathrm{C}(\alpha,\gamma)^{16}\mathrm{O}$ has two critical subthreshold resonances, one with $J^\pi = 1^-$ at $-45$ keV below threshold and one with $2^+$ at $-245$ keV, whose tails dominate $S(300\,\mathrm{keV})$ — the value at the Gamow peak of helium burning in red giants — and which remain the main source of uncertainty in the C/O ratio at the end of helium burning [deBoer et al. 2017] .

To treat consistently several resonances, possible interferences among them, and non-resonant contributions, the standard tool is the R-matrix formalism of Wigner and Eisenbud (1947). The R-matrix parameterizes the scattering matrix $S$ in terms of a small number of poles (the resonances) and channel parameters, and allows the simultaneous fitting of cross-section data from different channels — $(p,\gamma)$, $(p,\alpha)$, $(p,p')$ on the same compound nucleus, for example — coherently constraining the partial widths. Public programs such as AZURE2 and SAMMY are today standard tools for precision R-matrix analyses [Azuma et al. 2010] , and are the way in which the uncertainties on the rates of the key reactions are progressively reduced by combining experiments covering different energy ranges and different channels. For the “holy grail” reaction $^{12}\mathrm{C}(\alpha,\gamma)^{16}\mathrm{O}$, the R-matrix analysis of all existing data today returns $S(300\,\mathrm{keV}) = 140 \pm 21$ keV·barn — a 15% uncertainty, still insufficient for the most stringent applications, and one that justifies the new experiments under way at LUNA-MV and at JLab [deBoer et al. 2017] .

Neutron capture and photodisintegration

There is another large family of nuclear reactions, fundamental for the production of the heavy metals: neutron capture. Neutrons, being electrically neutral, see no Coulomb barrier and can be absorbed by nuclei even at very low thermal energies, of order the kT of the plasma. It is the mechanism by which more than half of the elements heavier than iron are produced: gold, platinum, the lanthanides, lead, uranium. A chain of $(n,\gamma)$ neutron captures alternating with $\beta^-$ decays makes it possible to cross the entire trans-iron region of the periodic table, and depending on the relative speed of n capture versus $\beta$ decay one has the s-process (slow capture, along the valley of stability) or the r-process (rapid capture, far from the valley).

For neutron capture $(n,\gamma)$ at low energies, the cross section often follows the so-called $1/v$ law, because a neutron with velocity $v$ “sees” the nucleus for a time $\propto 1/v$. The natural quantity to tabulate is the Maxwellian-Averaged Cross Section (MACS), obtained by integrating $\sigma$ over the thermal distribution of the neutrons:

$$\mathrm{MACS}(kT) = \frac{2}{\sqrt{\pi}} \int_0^\infty \sigma(E) \, \frac{E}{kT} \, e^{-E/kT} \, d\left( \frac{E}{kT} \right).$$

For the s-process in AGB stars one works at $kT \approx 8$-$30$ keV, and the standard compilation of the MACS in this interval is the KADoNiS database maintained at Karlsruhe [Karlsruhe Institute of Technology] ; the review by Käppeler and collaborators [Käppeler et al. 2011] is the methodological reference for the whole discipline.

Photodisintegration and detailed balance

In very hot environments, where the temperature of the plasma corresponds to a thermal radiation background with photons of energy comparable to the nuclear Q-values, the inverse process of capture also operates: the photons “strip” a neutron, a proton or an alpha particle from a nucleus. This is photodisintegration, written as $(\gamma, n)$, $(\gamma, p)$, $(\gamma, \alpha)$. It becomes relevant when the temperature exceeds $T \gtrsim 2 \times 10^{9}$ K (the blackbody radiation then contains photons with MeV energies in the Wien tail), and dominates at still higher temperatures. It is the mechanism underlying the final stages of silicon burning, the p-process (the synthesis of proton-rich nuclei by photodisintegration of s/r nuclei), and the freeze-out of statistical equilibrium.

The relation between capture and photodisintegration is fixed by detailed balance: in strict thermodynamic equilibrium, the rate of the direct reaction and that of the inverse one must coincide, and this determines the ratio between the thermal capture rate and the inverse photonuclear rate. For the generic reaction $A + a \rightleftharpoons B + \gamma$ one writes, in simplified form and without excited partition factors,

$$\lambda_{\gamma}(B \to A+a) = \frac{(2J_A+1)(2J_a+1)}{(2J_B+1)} \left( \frac{\mu kT}{2\pi\hbar^2} \right)^{3/2} \langle\sigma v\rangle_{A+a\to B+\gamma}\, e^{-Q/kT}.$$

Here $\lambda_\gamma$ is a rate per nucleus, not a cross section, and $Q$ is the Q-value of the direct capture. In complete calculations the nuclear partition functions and the thermally populated excited states are also included. The practical consequence is important: for every capture rate measured or calculated reliably, detailed balance provides the rate of the inverse reaction, with no need to measure it directly. It is an essential experimental economy, because many photodisintegration reactions on nuclei far from stability are simply impossible to measure with current techniques.

The Hauser-Feshbach statistical model

For the tens of thousands of reactions that enter a complete explosive nucleosynthesis network, only a minority fraction has been measured in the laboratory. For all the others — in particular those involving nuclei far from stability along the r and rp paths — the rate must be computed from a model. The standard statistical model is that of Hauser-Feshbach, based on the hypothesis that the reaction proceeds through a compound nucleus with a large number of overlapping levels, and that the cross section can be written as an incoherent sum over the open channels weighted by transmission coefficients $T_a$:

$$\sigma_{HF}(a \to b) = \frac{\pi}{k^2} \frac{1}{(2J_a+1)(2J_A+1)} \sum_{J,\pi} (2J+1) \frac{T_a \, T_b}{\sum_c T_c}.$$

The model is applicable when the level density of the compound nucleus is high — typically $A \gtrsim 40$ and at excitation energies above a few MeV — and is the basis of the nucleosynthesis codes for “unresolved” reactions. Public programs such as TALYS, EMPIRE and NON-SMOKER are used to generate Hauser-Feshbach rates systematically; a reference compilation of rates for the complete network is in [Rauscher et al. 2013] . The uncertainty of a Hauser-Feshbach rate is often a factor of a few close to stability, but can grow considerably for nuclei far from the valley, where $\gamma$-strength functions, level densities, optical potentials, masses and thresholds are not measured. For the r-process these uncertainties add to those on masses, $\beta$ half-lives, fission and astrophysical trajectories; reducing them is one of the main motivations of the new experiments at FRIB and FAIR (see the final section).

Nuclear statistical equilibrium

There is a regime of nucleosynthesis that deserves a section of its own because it is not a chain of reactions in the ordinary sense: it is a state of equilibrium. When the temperature of a nuclear plasma becomes high enough to make fast all the fusion, photodisintegration and capture reactions — in both the direct and the inverse directions — the composition of the plasma ceases to depend on the kinetic details and becomes a function only of three thermodynamic parameters: temperature $T$, density $\rho$, and the electron fraction

$$Y_e = \frac{n_p}{n_p + n_n} = \sum_i Z_i Y_i$$

where $n_p$ and $n_n$ denote total protons and neutrons, free or bound in nuclei; the sum runs over all nuclides and $Y_i = X_i / A_i$ is the number fraction of nuclide $i$. In this regime — nuclear statistical equilibrium or NSE — the relative abundances are given by the Saha equation applied to every equilibrium of the type $Z \cdot p + N \cdot n \rightleftharpoons (Z,N)$:

$$Y(Z,N) = G(Z,N) \, A^{3/2} \, 2^{-A} \left( \frac{\rho N_A}{\theta} \right)^{A-1} e^{B(Z,N)/kT} \, Y_p^{Z} \, Y_n^{N}$$

where $G$ is the partition function of the nuclide, $\theta = (m_u kT / 2\pi \hbar^2)^{3/2}$, and $Y_p$, $Y_n$ are the free fractions of protons and neutrons, determined by the conservation laws $\sum_i A_i Y_i = 1$ and $\sum_i Z_i Y_i = Y_e$. The exponential dependence on the binding energy $B(Z,N)$ has a spectacular consequence: at fixed $Y_e$, and when density and entropy allow heavy nuclei to dominate over free nucleons and $\alpha$ particles, the equilibrium composition concentrates near the most bound nuclides compatible with that $Y_e$.

In conditions of neutron/proton symmetry ($Y_e = 0{.}5$), the favored nuclide in the iron group is $^{56}\mathrm{Ni}$ — and indeed NSE at $Y_e$ close to 0.5 produces large quantities of nickel-56, which will subsequently decay into $^{56}\mathrm{Co}$ and then into $^{56}\mathrm{Fe}$. This is exactly what happens in the thermonuclear Type Ia supernovae: the explosive burning of the carbon-oxygen material of the white dwarf brings part of the ejecta into NSE at $Y_e \approx 0{.}5$, and the ejecta is dominated by $^{56}\mathrm{Ni}$ and nearby iron-group isotopes. The radioactive chain $^{56}\mathrm{Ni} \to {}^{56}\mathrm{Co} \to {}^{56}\mathrm{Fe}$ powers the light curve of the supernova in the following months. The quantity of $^{56}\mathrm{Ni}$ produced is one of the physical parameters that make SNe Ia standardizable candles for cosmology.

If instead the plasma is enriched in neutrons ($Y_e < 0{.}5$, the typical condition of the neutronized regions of the collapsing core of a core-collapse supernova), the most bound nuclide at fixed $Y_e$ shifts toward more neutron-rich isotopes of the iron group, and $^{58}\mathrm{Ni}$, $^{60}\mathrm{Ni}$, and then nuclei of the region $A \approx 80$-$90$ are produced. It is the nuclear basis of the neutronized contribution to the iron peak of the cosmic abundances, and it is the reason why the isotopic compositions of solar iron and nickel carry the imprint of the mean $Y_e$ of the progenitors that produced them.

Freeze-out and the $\alpha$-rich process

NSE is an instantaneous regime: it does not survive cooling. When the temperature drops below $T \lesssim 5 \times 10^{9}$ K, the triple and quadruple reactions necessary to maintain the equilibrium are no longer fast enough, and the composition “freezes” — it is the NSE freeze-out. If the expansion is fast enough (as in core-collapse supernovae), a fraction of free alpha particles remains trapped outside the final equilibrium, and this asymptotic regime is called the $\alpha$-rich process: it produces an excess of $\alpha$-multiple nuclei such as $^{40}\mathrm{Ca}$, $^{44}\mathrm{Ti}$, $^{48}\mathrm{Cr}$, $^{52}\mathrm{Fe}$. $^{44}\mathrm{Ti}$ in particular, with a half-life of about 60 years, is directly observable in the X-ray and gamma decay lines and has been mapped by NuSTAR in Cas A [Grefenstette et al. 2014] and detected in SN 1987A [Boggs et al. 2015] — a direct observational confirmation that the $\alpha$-rich freeze-out leaves measurable radioactive tracers in core-collapse supernovae.

For the r-process, finally, NSE in the extreme neutron-rich regime ($Y_e \lesssim 0{.}3$, $n_n \gtrsim 10^{24}\,\mathrm{cm^{-3}}$) generates an equilibrium distribution dominated by the competition neutron capture/photodisintegration on isotopes very far from the valley of stability, near the neutron drip line. When the explosion expands and cools, the equilibrium breaks, the capture chains stop, and the subsequent sequences of $\beta^-$ decays bring the nuclides back toward the valley of stability producing the final distribution of r-process abundances. All of this will be treated in detail in the dedicated chapters (silicon burning, supernovae, neutron-star mergers); but the conceptual structure — NSE, freeze-out, decay asymptotics — is the unifying framework of nucleosynthesis at extreme temperatures.

Measuring nucleosynthesis

Everything we have discussed so far — cross sections, thermal rates, astrophysical factors, nuclear masses — is in the end a number, and every number that enters a nucleosynthesis calculation comes from an experimental measurement or from a model tuned on measurements. The discipline that provides these numbers is low-energy experimental nuclear physics, and it deserves a concise overview, because without its silent work nucleosynthesis would be an exercise in blind extrapolation.

Accelerators at the surface and underground

The cross sections of fusion reactions between charged particles at low energies are so small — picobarns or less at the Gamow peak — that their direct measurement is limited by the cosmic-ray background. At the surface, every square meter of detector registers of order hundreds of muons per second from cosmic rays, generating a background that at low nuclear energies can swamp any signal. The solution, adopted since 1992 by the LUNA laboratory at the Gran Sasso National Laboratories, is to put the accelerator under the rock: the 1,400 meters of limestone and dolomite above the laboratory reduce the muon flux by about six orders of magnitude, cutting the background to a level that allows cross sections to be measured down to energies close to the Gamow peak of many astrophysical reactions [Broggini et al. 2010] . The sequence of reactions measured by LUNA over the last thirty years includes $^{3}\mathrm{He}(^{3}\mathrm{He}, 2p)^{4}\mathrm{He}$ (pp chain), $^{14}\mathrm{N}(p,\gamma)^{15}\mathrm{O}$ (the slow step of the CNO cycle), $^{3}\mathrm{He}(\alpha,\gamma)^{7}\mathrm{Be}$ (critical for solar neutrinos and BBN), $^{22}\mathrm{Ne}(p,\gamma)^{23}\mathrm{Na}$ (NeNa cycles in AGB), and many others; the new LUNA-MV accelerator (3.5 MV) has extended the program to reactions at higher energy. JUNA, at the Chinese Jinping laboratory, and CASPAR, at the Sanford Underground Research Facility in the United States, follow the same model with complementary programs.

Radioactive beams

A significant part of nucleosynthesis occurs through radioactive nuclei that live milliseconds or seconds, and that cannot be used as fixed targets: the reactions must be measured “in reverse”, with a beam of the radioactive nucleus bombarding a hydrogen or helium target. This is the domain of radioactive beams, a technique that has transformed the discipline over the last twenty years. The main facilities include FRIB (Facility for Rare Isotope Beams) at Michigan State University, operational for users since 2022 [U.S. Department of Energy Office of Science] ; FAIR (Facility for Antiproton and Ion Research) at Darmstadt, a European infrastructure under construction and progressive commissioning [FAIR GmbH] ; HIE-ISOLDE at CERN; RIBF at RIKEN in Japan; and TRIUMF-ISAC in Canada. These facilities measure half-lives, masses, indirect cross sections, breakup and decays of nuclei along the r and rp paths. Over the next decade they will reduce many systematic uncertainties, without closing them entirely: the astrophysical rates will remain tied also to fission, to masses not yet reachable, to excited states and to the hydrodynamic conditions of the site.

Precision traps

For the masses of exotic nuclei, the instrument of absolute precision is the Penning trap: the ions are confined in a magnetic-electrostatic field, their cyclotron frequency is measured with parts-per-billion precision, and the mass is extracted with uncertainties of a few keV even for isotopes that live milliseconds. ISOLTRAP at CERN, JYFLTRAP at Jyväskylä, LEBIT at FRIB are the main laboratories of the sector, and their data flow into the periodic update of the AME [Wang et al. 2021] . A complementary technique, Storage Ring Mass Spectrometry, is employed at GSI in Germany and at CSRe in China for the masses of extremely exotic nuclei: the ions circulate in a storage ring for milliseconds, and the mass is measured from the revolution frequency.

Databases and uncertainty propagation

All this experimental work flows into databases that feed the nucleosynthesis codes. The three most important have already been named, but it is worth recapping them with their role:

• AME (Atomic Mass Evaluation) [Wang et al. 2021] : the authoritative compilation of nuclear masses, updated every ~4 years.
• NNDC ENSDF [Brookhaven National Laboratory] : the energy levels, the lifetimes, the decay schemes.
• IAEA Nuclear Data Section [International Atomic Energy Agency] : a broader infrastructure that includes evaluated cross sections, reactor data, astrophysical data.
• JINA REACLIB [Cyburt et al. 2010] : thermonuclear rates in parameterized form for extended networks.
• KADoNiS [Karlsruhe Institute of Technology] : MACS and neutron-capture data for the s-process.

The astrophysical rates tabulated as functions of temperature for use in nucleosynthesis codes are compiled primarily in the JINA REACLIB library, with tens of thousands of reactions in standardized parameterizations. The uncertainties on these rates are propagated onto the nucleosynthesis results through Monte Carlo techniques: each rate is sampled within its uncertainty distribution (typically lognormal), the nucleosynthesis calculation is repeated thousands of times, and confidence bands on the final abundances are produced [Longland et al. 2010] . This methodology, today standard in codes and libraries such as PRIMAT for BBN and SkyNet for explosive nucleosynthesis [Lippuner & Roberts 2017] , makes it possible to answer the question how much do we really know about the nucleosynthesis of the various regions of the periodic table, and to direct the experimental measurement campaigns toward the reactions that dominate the global uncertainty (the so-called sensitivity studies).

With this conceptual and experimental apparatus in hand we are ready to go back to the first twenty minutes after the big bang and apply all of this to primordial nucleosynthesis — the first of the concrete astrophysical sites we shall meet in the book.

The first three minutes

The universe, in its very first instants, is an extremely hot and dense soup of elementary particles: quarks and gluons at one microsecond, and then — past the confinement transition of chromodynamics — protons, neutrons, electrons, positrons, neutrinos and photons in thermal equilibrium. The expansion cools it according to a well-defined law, $T \propto 1/a(t)$, which links the temperature to the scale factor of the universe, and $a(t) \propto t^{1/2}$ in the radiation-dominated phase. The sequence of nuclear events unfolds entirely within an extremely narrow window — between the first second and the twentieth minute of cosmic life — during which the temperature passes from some ten billion kelvin to a few hundred million. It is, in every sense, a windowed nucleosynthesis: before one second it is too hot, because any nucleus that managed to form would be immediately dissociated by the very-high-energy photons; after twenty minutes it is too dilute, because the nucleons no longer meet frequently enough to fuse before the expansion separates them definitively.

Between these two temporal extremes unfolds primordial nucleosynthesis or Big Bang Nucleosynthesis (BBN), the first of the concrete astrophysical sites that nucleosynthesis treats. The final balance is simple and quantitatively precise: about 75% of the baryonic mass remains as hydrogen, about 25% ends up in $^{4}\mathrm{He}$, and a minute trace — one atom in a hundred thousand — in deuterium; still smaller amounts in $^{3}\mathrm{He}$ and $^{7}\mathrm{Li}$, and practically nothing of all the rest. No carbon, no oxygen, no iron: the elements heavier than helium require environments that BBN cannot provide, and to produce them one must wait for the ignition of the first stars a few hundred million years later. The first stellar generations are therefore born in a chemically almost virgin gas, and must cook the whole periodic table from scratch starting from primordial H and He. This is why BBN is not just a chapter of cosmology, but the initial condition of all the subsequent chemistry of the cosmos.

The precise prediction of the primordial abundances is one of the most solid quantitative successes of the standard cosmological model. With a single free parameter — the baryon-to-photon ratio $\eta = n_B/n_\gamma$, which fixes the baryon density today — the theory predicts $^{4}\mathrm{He}$ with sub-percent uncertainty, while the observational comparison is today limited by systematics at the percent level; for $D/H$ theory and observations both reach a precision of order one percent. The only persistent tension concerns the $^{7}\mathrm{Li}$ measured in the oldest stars. The same $\eta$ is today independently determined from the anisotropies of the cosmic microwave background (CMB) measured by Planck with 0.5% precision [Collaboration 2020] ; that the two determinations — one based on nuclear physics at $T \sim$ MeV in the first minute, the other on photon-baryon plasma physics at $T \sim$ eV after 380,000 years — coincide within the uncertainties is one of the most stringent confirmations of the whole of standard cosmology.

The thermodynamic picture

The prelude to BBN is the history of the neutron/proton ratio. At temperatures $T \gtrsim 1$ MeV ($t \lesssim 1$ s) the rates of the weak reactions

$$n + e^{+} \rightleftharpoons p + \bar\nu_e, \qquad n + \nu_e \rightleftharpoons p + e^{-}, \qquad n \rightleftharpoons p + e^{-} + \bar\nu_e$$

are much faster than the expansion rate $H(T)$ of the universe, and keep $n$ and $p$ in statistical equilibrium. The equilibrium ratio is fixed by the mass difference $Q = (m_n - m_p) c^2 \approx 1{.}293$ MeV:

$$\left( \frac{n_n}{n_p} \right)_{\mathrm{eq}} = \exp\!\left( - \frac{Q}{kT} \right).$$

At $T = 10$ MeV the ratio is almost 1; at $T = 1$ MeV it is already $\sim 0{.}25$. As $T$ falls, the exponential imposes $n_n \ll n_p$, and if the equilibrium were maintained down to the cold there would be no neutrons left, and hence neither deuterium nor helium. What saves BBN is that the equilibrium is not maintained: at $T \approx 0{.}7$ MeV, the dilution of the universe makes the weak rates too slow compared with $H(T)$, and the $n/p$ ratio freezes at its equilibrium value at that moment, $n_n/n_p \approx 1/6$. This freezing — the neutron freeze-out — is the first event that marks BBN, and the physics that determines it is purely nuclear-weak (the mass of the neutron, the Fermi constant, the $\tau_n$ of $\beta$ decay) combined with the Friedmann expansion rate.

Immediately after the freeze-out, the residual neutrons decay freely with a mean lifetime of order $\tau_n \simeq 880$ s, reducing the ratio from $1/6$ to $\sim 1/7$ by the time the actual nucleosynthesis begins. The datum $\tau_n$ is one of the few numbers of basic nuclear physics whose uncertainty still propagates appreciably onto the prediction of $Y_p$ (the mass fraction of $^{4}\mathrm{He}$); between the two measurement methods — beam (beams of cold neutrons crossing a detection region, with counting of the protons produced) and bottle (ultracold neutrons confined in a magnetic or material bottle, with counting of the survivors) — there is a historical discrepancy not yet fully resolved, on which dedicated experiments are under way [Particle Data Group 2024] .

At this point, one would naturally expect the nucleons to begin fusing in pairs to form deuterium:

$$p + n \to d + \gamma, \qquad Q = 2{.}224 \, \mathrm{MeV}.$$

The cross section of this reaction is large (it is an M1 capture without Coulomb barrier, the neutron being neutral), and naively one would expect deuterium to begin forming as soon as $kT \lesssim Q$, that is, around 0.3 MeV ($t \sim 10$ s). It is not so: the deuterium is continuously photodissociated by the high-energy tail of the Planck distribution of the photons, because the photon-to-baryon ratio is enormous ($n_\gamma / n_B = 1/\eta \approx 1{.}6 \times 10^{9}$), and a few photons in the tail suffice to keep the deuterium dissociated down to much lower temperatures. This is the famous deuterium bottleneck: nucleosynthesis cannot proceed beyond deuterium until the photon tail has sufficiently emptied, and this happens only around $T \approx 80$ keV ($t \approx 3$ min). The quantitative criterion is written as $\eta \, e^{-Q/kT} \sim 1$, whence $kT_{\mathrm{BBN}} \sim Q / |\ln \eta| \approx 0{.}07$ MeV — well below $Q$ but fixed by the logarithm of $\eta$, a merely logarithmic dependence that explains why the opening temperature of the bottleneck is little sensitive to the cosmological details.

Once the bottleneck is dissolved, the reactions chain together extremely rapidly:

$$d + d \to {}^{3}\mathrm{He} + n, \qquad d + d \to t + p, \qquad d + p \to {}^{3}\mathrm{He} + \gamma$$ $${}^{3}\mathrm{He} + d \to {}^{4}\mathrm{He} + p, \qquad t + d \to {}^{4}\mathrm{He} + n.$$

Practically all the remaining neutrons end up in $^{4}\mathrm{He}$, the most bound nuclide among those accessible in these conditions ($B/A \approx 7$ MeV), and the mass fraction of helium produced is simply fixed by the $n/p$ ratio at the beginning of the chain:

$$Y_p \equiv \frac{4 n_{^{4}\mathrm{He}}}{n_B} \approx \frac{2 \, n_n/n_p}{1 + n_n/n_p} \approx 0{.}247.$$

It is a prediction of surprising simplicity: $Y_p$ depends almost only on weak physics and on $\tau_n$, and is insensitive at first approximation to the precise value of $\eta$. The reaction chain then quickly exhausts itself as $T$ falls further: Coulomb capture between $^{4}\mathrm{He}$ and other nuclei requires temperatures higher than those available, and the absence of stable nuclides at $A = 5$ and $A = 8$ (the $^{8}\mathrm{Be}$ decays into two α particles in $10^{-16}$ s) creates a second bottleneck that BBN, unlike the stars, fails to overcome. There remain in traces $^{3}\mathrm{He}$ (residues of the incomplete burning of deuterium) and $^{7}\mathrm{Li}$ (produced both by $t(\alpha,\gamma)^{7}\mathrm{Li}$ at low $\eta$ and by $^{3}\mathrm{He}(\alpha,\gamma)^{7}\mathrm{Be}$ with subsequent decay at high $\eta$); the abundances of $^{6}\mathrm{Li}$, $^{9}\mathrm{Be}$, $^{10,11}\mathrm{B}$ predicted by standard BBN are many orders of magnitude below their observed value, a sign that these rare nuclides are produced elsewhere — in the cosmic rays interacting with the interstellar medium, as we shall see in the final part of this chapter, devoted to the origin of LiBeB.

The physics hidden behind the formulas

It is worth noting how almost every element of BBN is a combination of physics known at the laboratory level. The weak reaction constants that fix $n/p$ at freeze-out are parameterized by the Fermi constant $G_F$ measured in muon decays and in $\beta$ decay; the neutron lifetime is measured directly in the laboratory; the cross sections of the nuclear reactions $d(p,\gamma){}^{3}\mathrm{He}$, $d(d,n){}^{3}\mathrm{He}$, $d(d,p)t$, $^{3}\mathrm{He}(d,p){}^{4}\mathrm{He}$, $t(d,n){}^{4}\mathrm{He}$ have been measured for decades with percent precision at laboratory energies overlapping the BBN Gamow peak. Gravity enters via $H(T) = \sqrt{8 \pi G \rho_{\mathrm{rad}} / 3}$, with $\rho_{\mathrm{rad}}$ fixed by the relativistic species present (photons, neutrinos, $e^{\pm}$ before annihilation). And thermal statistics enters via the Maxwell-Boltzmann distributions for nuclei and the Bose-Einstein/Fermi-Dirac ones for photons, $e^{\pm}$ and neutrinos. There are no exotic ingredients, no ad hoc parameters, no speculative assumptions: BBN is — from the point of view of the physics used — trivial, in the noblest sense of the word. And yet it produces quantitative predictions about an epoch 13.8 billion years in the past.

The primordial yields

The final BBN abundances, predicted by the modern codes with $\eta = (6{.}10 \pm 0{.}04) \times 10^{-10}$ from Planck [Collaboration 2020] , are compared with the observational measurements in the following table:

Each of these rows deserves a few words, because the physics behind them is qualitatively different.

Helium-4 is the most abundant species produced by BBN, and its mass fraction $Y_p$ is the most robustly predicted quantity. As anticipated, it depends almost only on $\tau_n$ and on weak physics, and is practically insensitive to $\eta$: in the interval $\eta \in [10^{-10}, 10^{-9}]$, $Y_p$ varies by less than 5%. Its observational measurement is carried out in metal-poor HII regions — clouds of ionized gas in dwarf galaxies poor in heavy elements — where the emission lines of He and H allow the helium fraction to be estimated and extrapolated to zero metallicity ($O/H \to 0$) to obtain the “primordial” $Y_p$, purged of the contribution of subsequent stellar nucleosynthesis. The current uncertainty, $\delta Y_p \sim 0{.}003$, is dominated by systematics of collisional emission, electron temperature and plasma density, not by statistics; recent surveys such as EMPRESS (Subaru) have added a hundred galaxies at even lower metallicity and refined the estimate, but the level of precision needed to discriminate sub-percent models on $N_{\nu}^{\mathrm{eff}}$ (see the final section) still requires a systematic leap.

Deuterium, at the opposite end, is the species most sensitive to $\eta$: $D/H \propto \eta^{-1{.}6}$ approximately, because at higher $\eta$ the burning of deuterium into $^{3}\mathrm{He}$ and $^{4}\mathrm{He}$ is more efficient and less of it remains. This makes it an excellent baryometer — an instrument to measure $\eta$ — and its measurement is today the most precise BBN-based determination of $\Omega_b h^2$. The standard technique is absorption spectroscopy of Lyman-α gas clouds at high redshift ($z = 2$-$4$) along the line of sight of luminous quasars: at these epochs the gas is still chemically almost virgin, stellar enrichment is minimal, and the measured $D/H$ is a good approximation of the primordial value. The measurements of Pettini, Cooke and collaborators on a sample of seven Damped Lyman-α systems have reached a precision of about 1% [Cooke et al. 2018] : the result $D/H = (2{.}547 \pm 0{.}033) \times 10^{-5}$ is in full agreement with the BBN prediction based on $\eta_{\mathrm{Planck}}$, and provides — combined with the CMB — one of the most precise tests of the effective number of neutrino species $N_{\nu}^{\mathrm{eff}}$.

The theoretical predictions of $D/H$, however, required considerable nuclear work to reach that precision. The rate of the reaction $d(p,\gamma){}^{3}\mathrm{He}$ — the first step of deuterium burning — was for a long time the dominant uncertainty; the LUNA measurement at Gran Sasso reduced this uncertainty to less than 3% at the BBN Gamow peak ($E \sim 100$ keV) [Mossa et al. 2020] , and the corresponding predicted value of $D/H$ shifted slightly downward, bringing it into even better agreement with the observations. It is a concrete example of the way BBN has been transforming, over the last twenty years, from an approximate cosmological test into a precision cosmological measurement.

Helium-3 is produced in small quantity by BBN and in small quantity also by stars, and its subsequent chemistry is complicated: low-mass stars produce it in the H-burning branch and return it in part to the interstellar medium, but high-mass stars destroy it. Its measurement, conducted mostly in Galactic HII regions through the hyperfine line of $^{3}\mathrm{He}^{+}$ at 3.46 cm (Bania, Rood, Balser), is difficult and provides only order-of-magnitude limits; it is not today a stringent test of BBN, but it is one of the few direct measurements of active nucleosynthesis in the Galaxy.

Lithium-7 is the problematic species. The abundances measured in the most metal-poor stars of our Galactic halo — the famous Spite plateau discovered in 1982 by François and Monique Spite [Spite & Spite 1982] — are remarkably constant over a wide interval of metallicity, around $A(\mathrm{Li}) = \log_{10}(n_{\mathrm{Li}}/n_H) + 12 \approx 2{.}2$, corresponding to $^{7}\mathrm{Li}/H \approx 1{.}6 \times 10^{-10}$. The BBN prediction, with $\eta_{\mathrm{Planck}}$, is $(4{.}7 \pm 0{.}7) \times 10^{-10}$: a factor $\sim 3$ more than the observed value, a discrepancy of about $3\sigma$ that has now lasted two decades and has taken the name of the cosmological lithium problem.

The cosmological lithium problem

On the origins of the lithium problem three types of explanation have been hypothesized, in decreasing order of current consensus:

1. Stellar depletion. $^{7}\mathrm{Li}$ is a nuclearly fragile species: it burns by proton capture already at temperatures of $\sim 2{.}5 \times 10^{6}$ K, well below the central temperatures of a solar-mass main-sequence star. If in the atmospheres of the old Pop II stars — those in which the Spite plateau is measured — any mixing mechanism carries the lithium from the surface to a depth sufficient for burning, the measured surface value turns out lower than the primordial value inherited from the protostellar cloud. Standard stellar models do not predict efficient mixing in these stars, but the combination of atomic diffusion (gravitational settling) — which makes the heavy elements migrate downward — and of empirically calibrated turbulence allows reductions of $0{.}2$-$0{.}4$ dex, sufficient to bridge the discrepancy. The direction of consensus over the last fifteen years has been consolidating on this scenario, supported by abundance measurements in stars of globular clusters such as NGC 6397 (Korn et al. 2007) [Korn et al. 2007] and by more recent observations in metal-poor field stars [Mucciarelli et al. 2022] .

2. Nuclear uncertainties. The reactions that produce $^{7}\mathrm{Li}$ — directly via $t(\alpha,\gamma)^{7}\mathrm{Li}$ at low $\eta$, and indirectly via $^{3}\mathrm{He}(\alpha,\gamma)^{7}\mathrm{Be}$ with subsequent decay $^{7}\mathrm{Be}(e^{-},\nu)^{7}\mathrm{Li}$ at high $\eta$ — and those that destroy it — $^{7}\mathrm{Be}(n,p)^{7}\mathrm{Li}$ and $^{7}\mathrm{Be}(n,\alpha)^{4}\mathrm{He}$ — have been the object of intensive measurement campaigns over the last fifteen years at LUNA, n_TOF (CERN), JYFLTRAP, and other laboratories. The cumulative uncertainties on the rates have been reduced below 5%, and the predicted value of $^{7}\mathrm{Li}/H$ is consequently well constrained: there is no residual nuclear margin to explain the discrepancy.

3. New physics. Models in which massive particles of intermediate lifetime ($\tau \sim 10^{3}$-$10^{7}$ s) decay during or after BBN can in principle selectively destroy $^{7}\mathrm{Be}$ (and hence reduce the final $^{7}\mathrm{Li}$) without altering $Y_p$ or $D/H$. Scenarios of this type include gravitinos, light neutralinos, sneutrinos, and sterile-neutrino portals. All these proposals are however strongly constrained by the ensemble of the other BBN+CMB observables, and none is favored today; most require fine tuning of the parameters not to break the agreement on deuterium and helium.

The emerging consensus, reflected in the most recent reviews [Fields et al. 2020] [Particle Data Group 2024] , is that the lithium problem is stellar — a problem of transport and mixing in Pop II atmospheres — and not cosmological. The situation is one of slow erosion: no single experimental breakthrough, but an accumulation of evidence (reduction of the nuclear uncertainties, new abundance measurements, more sophisticated diffusion + turbulence simulations) that has shifted the weight from “new physics” to “poorly understood stellar physics”. It is an interesting case also from the methodological point of view, because it illustrates how an apparently “fundamental” discrepancy can over time turn out to be a masked measurement systematic.

BBN as a laboratory of fundamental physics

The precision with which the BBN predictions match the observations — lithium aside — has transformed the discipline into a quantitative probe of physics beyond the Standard Model in the conditions of the primordial universe. Practically every ingredient of standard cosmology can be varied and its effects constrained by comparing the new set of BBN abundances with the measurements.

The most important constraint concerns the number of relativistic species at the BBN epoch, parameterized by the parameter $N_{\nu}^{\mathrm{eff}}$, defined so that in the Standard Model it equals about $3{.}045$ (the three neutrinos plus small corrections from the non-instantaneity of decoupling, finite-temperature QED, oscillations; in much of the literature the value is rounded to $3{.}046$). A value of $N_{\nu}^{\mathrm{eff}}$ above the standard one would correspond to the presence of an additional relativistic species — a thermalized sterile neutrino, a light axion, a new boson — at the time of BBN; it would increase $H(T)$, anticipate the neutron freeze-out (leaving a higher $n/p$ ratio), and hence increase $Y_p$. The combination of the observed value of $Y_p$ with the deuterium today provides $N_{\nu}^{\mathrm{eff}} = 2{.}88 \pm 0{.}27$ [Pitrou et al. 2018] , perfectly compatible with the Standard Model, and excludes with good robustness a fourth fully thermalized active species. The BBN + CMB combination tightens the constraint further, and the next-generation CMB analyses aim to bring $\sigma(N_{\nu}^{\mathrm{eff}})$ to a few hundredths, a precision sufficient to test many models of dark radiation.

A second class of constraints concerns the variation of the fundamental constants between the BBN epoch and today. A variation of the fine-structure constant $\alpha$ of one part in $10^{4}$ would produce a percent-level change in the rates of the charged reactions and hence a visible deviation in $D/H$; similarly for the ratio $m_e/m_p$, or for the scale $\Lambda_{\mathrm{QCD}}$, which enters the masses of $n$ and $p$. The gravitational constants are constrained at the level $|\delta G/G| \lesssim 0{.}1$ at $t \sim 100$ s, and the constraints on the lepton asymmetry are $|\xi_e| = |\mu_{\nu_e}/T| \lesssim 0{.}02$. These are limits that, for such remote epochs, no direct laboratory experiment can reach.

There is a third family of constraints, more speculative but very binding: non-thermal nucleosynthesis induced by the decays of relic particles with lifetimes in the interval $10^{3}$-$10^{12}$ s. Particles of this type would decay after the end of standard BBN, producing high-energy photons or hadrons that would fragment the primordial nuclei (photodisintegration of $D$ and $^{4}\mathrm{He}$, secondary formation of $^{6}\mathrm{Li}$, $^{3}\mathrm{He}$, $^{7}\mathrm{Be}$). The BBN limits on the abundance and lifetime of these particles constrain scenarios of gravitinos, KK-modes, light neutralinos, and extended supersymmetry setups, with a sensitivity that often exceeds that of direct experiments.

BBN, in other words, is not just a test of the Standard Model: it is a probe of the state of the universe before recombination, and it is — together with the CMB — the oldest observational window we have. Its temporal geometry and its information density make it irreplaceable for any discussion of physics beyond the Standard Model in the first cosmological phases.

The modern BBN calculation networks

The modern treatment of BBN solves a coupled system of equations that includes: the Friedmann equations for the thermal evolution of the background (with $g_{*} \simeq 10{.}75$ when photons, electrons, positrons and neutrinos all contribute as relativistic species, and $g_{*} \simeq 3{.}36$ after the $e^{\pm}$ annihilation, when photons and neutrinos cooled with respect to the photons remain), a Boltzmann equation for the $n/p$ ratio with the complete weak rates and the electromagnetic corrections, finite-temperature QED, and one-loop contributions, and a network of nuclear reactions with about 25 dominant reactions (plus dozens more of secondary importance). The propagation of the uncertainties on the nuclear rates and on $\tau_n$ is carried out via Monte Carlo, sampling each rate within its uncertainty distribution (typically lognormal) and propagating the calculation thousands of times to produce confidence bands on the final abundances.

The public reference codes are PArthENoPE [Pisanti et al. 2021] , AlterBBN, PRIMAT [Pitrou et al. 2018] and the more recent PRyMordial (2023, in Python). The predictions of the four codes agree today at the sub-percent level on the main yields; the residual differences derive mainly from the choice of the cross-section compilations used (NACRE-II vs. more recent evaluations) and from the treatment of the sub-leading finite-temperature QED contributions. The current precision on the predicted primordial yields is $0{.}1\%$ on $Y_p$, $1\%$ on $D/H$, and about $10\%$ on $^{7}\mathrm{Li}/H$ (dominated by the residual experimental uncertainties on the Be reactions).

Dominant reactions and their uncertainties

The BBN network is small enough to allow an extremely detailed sensitivity analysis: each rate can be varied individually within its uncertainty and the effect on the final yields observed. The sensitivity is not uniform, and most of the final uncertainties are concentrated in a few key reactions:

• $\tau_n$ (neutron lifetime): the current uncertainty $\sim 0{.}07\%$ contributes $\delta Y_p / Y_p \sim 0{.}04\%$, and the beam/bottle discrepancy introduces an additional unresolved systematic;
• $p(n,\gamma)d$: measured in the laboratory with 1% precision, contributes $\delta D/D \sim 0{.}5\%$;
• $d(p,\gamma){}^{3}\mathrm{He}$: the LUNA measurement [Mossa et al. 2020] reduced the uncertainty to 3%, leaving a contribution of $\delta D/D \sim 1\%$;
• $d(d,n){}^{3}\mathrm{He}$ and $d(d,p)t$: uncertainties $\sim 1\%$ each, similar contributions to $D/H$ and $^{3}\mathrm{He}/H$;
• $^{3}\mathrm{He}(\alpha,\gamma)^{7}\mathrm{Be}$: measured by LUNA and by Notre Dame at 5%, the main source of uncertainty on $^{7}\mathrm{Li}/H$;
• $^{7}\mathrm{Be}(n,p)^{7}\mathrm{Li}$ and $^{7}\mathrm{Be}(n,\alpha)^{4}\mathrm{He}$: measured at n_TOF at CERN, they have reduced the uncertainty on $^{7}\mathrm{Li}/H$ to about $10\%$.

The direction of the coming years is clear: even more precise measurements of $^{3}\mathrm{He}(\alpha,\gamma)^{7}\mathrm{Be}$ and $d(p,\gamma){}^{3}\mathrm{He}$ in underground laboratories and at low energy, definitive resolution of the $\tau_n$ discrepancy, and new observations of $D/H$ with very-high-resolution spectroscopy on ELT-class telescopes. With these improvements, BBN could enter the joint sub-percent precision band with the CMB and provide constraints on the $\eta$-$N_{\nu}^{\mathrm{eff}}$ plane at a level competitive with the great cosmological surveys of the coming decades.

With BBN we have fixed the chemical initial condition of the universe. Everything we shall see in the following chapters — hydrogen burning in main-sequence stars, helium in red giants, neutron capture in AGB stars and in neutron-star mergers — is a process of progressive enrichment starting from this initial composition: 75% H, 25% He, and a trace of deuterium. The first stars to be born burn this primordial mixture, with very little metal cooling and with a structure determined precisely by the initial absence of carbon, oxygen, silicon and iron. The next step is therefore natural: to study hydrogen burning, the stable engine that transforms the primordial chemical simplicity into the first lasting stellar nucleosynthesis.

The gap at Z = 3, 4, 5

If we line up the elements of the periodic table by cosmic abundance, after hydrogen and helium there should be lithium, beryllium and boron — the next three elements in the sequence of atomic numbers. The observed abundances show instead a deep collapse: $\mathrm{Li}$, $\mathrm{Be}$ and $\mathrm{B}$ are rare by a factor of a thousand or more compared with the neighboring elements, forming a characteristic gap in the cosmic abundance curve between He ($Z = 2$) and C ($Z = 6$). They are not rare by chance, however: they are the most fragile nuclei of the periodic table, with low binding energies per nucleon and destruction thresholds at temperatures accessible in the atmospheres and sub-surface shells of stars.

The reason for the collapse is twofold. First: stars, when they encounter LiBeB in their sub-surface layers, destroy them via $(p,\alpha)$ reactions at temperatures of a few million kelvin — well below the temperatures of the standard quiescent burnings. Second: their production does not occur mainly in stars, but in interstellar space through a qualitatively different mechanism, the spallation by cosmic rays (GCR), in which a proton or an $\alpha$ nucleus accelerated to GeV energies strikes a C, N or O atom of the interstellar medium, fragmenting it into lighter nuclei among which $^{6,7}\mathrm{Li}$, $^{9}\mathrm{Be}$, $^{10,11}\mathrm{B}$. It is a channel of Galactic (not stellar) nucleosynthesis, and it contributes dominantly to the cosmic budget of Be and B, and significantly to Li.

Quantitatively, the proton-destruction thresholds of the LiBeB are low and well characterized: $^{7}\mathrm{Li}(p,\alpha)^{4}\mathrm{He}$ at $T \gtrsim 2{.}5 \times 10^{6}$ K, $^{9}\mathrm{Be}(p,\alpha)^{6}\mathrm{Li}$ at $T \gtrsim 3{.}5 \times 10^{6}$ K, $^{10,11}\mathrm{B}(p,\alpha)$ at $T \gtrsim 5 \times 10^{6}$ K. In a main-sequence star, the convective envelope transports the surface material to depths where $T$ exceeds these thresholds, and the lithium observed at the surface decreases over time — it is the classic Li depletion, well documented in open clusters of different ages (Pleiades, Hyades, Praesepe, M67). To maintain the cosmic LiBeB budget in the face of this continuous destruction in stars, primary sources external to the stellar cycle are needed. The known production channels are five: the spallation of Galactic cosmic rays (GCR) on interstellar CNO; the inverse spallation (GCR-CNO striking ambient H, He); the BBN contribution to the primordial $^{7}\mathrm{Li}$ (previous part of this chapter); the $\nu$-process reaction in core-collapse SNe for $^{7}\mathrm{Li}$ and $^{11}\mathrm{B}$; the contribution of novae to $^{7}\mathrm{Li}$ via $^{3}\mathrm{He}(\alpha,\gamma)^{7}\mathrm{Be}$ and fast convective transport (chapter 5).

The idea that cosmic rays were responsible for the production of LiBeB was introduced by Reeves, Fowler and Hoyle in 1970 [Reeves et al. 1970] , replacing the original x-process of B²FH (1957), which hypothesized an unidentified specific stellar site. The proposal passed two crucial tests: the relative abundance $^{6}\mathrm{Li}/^{7}\mathrm{Li}/^{9}\mathrm{Be}/^{10}\mathrm{B}/^{11}\mathrm{B}$ predicted by standard spallation is in good agreement with the solar values, and the ratio $^{11}\mathrm{B}/^{10}\mathrm{B} \approx 4$ requires a contribution from low-energy GCR (or from the $\nu$-process), because the standard high-energy cross-section ratio predicts $^{11}\mathrm{B}/^{10}\mathrm{B} \approx 2{.}5$ — the first observational hint of the need for an explosive neutrino contribution to the LiBeB budget. Later updates (Ramaty [Ramaty et al. 1997] , Prantzos [Prantzos 2012] ) refined the decomposition between low-energy GCR (LECR) and high-energy GCR, and quantitatively introduced the contributions of the $\nu$-process and of novae to the total budget.

Galactic spallation: mechanics and budget

Galactic cosmic rays are particles accelerated to relativistic velocities, composed for $90\%$ of protons, for $9\%$ of He nuclei, and for $1\%$ of heavier nuclei, with a spectrum extending from a few MeV to $\sim 10^{20}$ eV. The sub-PeV component is accelerated mainly in the shock waves of core-collapse supernovae via Diffusive Shock Acceleration (Krymsky, Bell, Blandford-Ostriker), and crosses the Galaxy confined by the Galactic magnetic field with a residence time of $\sim 10^{7}$-$10^{8}$ years. When one of these GCR particles strikes a nucleus of the interstellar gas — an atom of C, N or O — it can “chip” its structure through a highly inelastic collision, tearing off nuclear fragments: thus are born $^{6,7}\mathrm{Li}$, $^{9}\mathrm{Be}$, $^{10,11}\mathrm{B}$ and other minor fragments. The GCR flux is low but constant in time, and multiplied by the billions of years of the Galaxy’s life it produces a cumulative LiBeB budget that reasonably matches what is observed in the Sun and in Pop I stars — one of the success stories of post-B²FH nuclear astrophysics.

In quantitative terms, the production rate of a nuclide $i$ by spallation is given by the convolution of the GCR flux with the fragmentation cross sections:

$$\dot N_i = \int dE\,\phi(E)\,\sum_j n_j\,\sigma_{ji}(E),$$

with $\phi(E)$ the differential GCR flux per energy per nucleon, $n_j$ the density of the target atoms in the interstellar medium, $\sigma_{ji}(E)$ the fragmentation cross section of target $j$ into product $i$. The analytic structure of the problem allows two complementary kinematic regimes to be distinguished. The “secondary” production corresponds to the case of light GCR (p, $\alpha$) striking CNO nuclei of the interstellar medium: the production rate scales as $\dot N \propto Z_{\mathrm{ISM}}$ — it grows linearly with the metallicity of the surrounding gas, because the density of the targets increases. The “primary” production corresponds to the inverse case: heavy GCR (C, N, O accelerated in SNe) striking the ambient H, He of the interstellar medium; the rate scales as $\dot N \propto Z_{\mathrm{GCR}}$, substantially constant in time because the GCR come from contemporary local SNe, independent of the historical ambient metallicity. At low metallicity ($[\mathrm{Fe/H}] \lesssim -1$) the primary contribution dominates; at high metallicity the secondary dominates. The standard prediction is that $\log(\mathrm{Be/H})$ vs $[\mathrm{Fe/H}]$ has a slope close to $1$ in the purely secondary regime, and the observations of Smiljanic et al. (2009) and later ones show a slope compatible with $1$ within the uncertainties over a wide interval of metallicity — strong confirmation of the mechanism.

The fragmentation cross sections are measured experimentally in accelerator experiments (Read & Viola 1984, Webber et al., Tatischeff 2007 and later), with nuclear-data compilations updated periodically. The databases in use for the transport of GCR in the Galaxy are GALPROP and DRAGON, coupled to an integrated LiBeB nucleosynthesis network (the formalism of Prantzos [Prantzos 2012] ). A delicate point is the low-energy cosmic ray spectrum (LECR) below $\sim 1$ GeV/nucleon: the GCR observed from Earth have a spectrum deformed by the modulation of the solar wind, and the true flux below the modulation is uncertain until direct measurements outside the heliopause are made. The production of $^{6}\mathrm{Li}$ is particularly sensitive to LECR because the channel $\alpha + \alpha \to {}^{6,7}\mathrm{Li}$ has a low threshold ($\sim 10$ MeV/nucleon). Direct evidence of LECR on the ISM comes from the $\gamma$ observations at $4{.}4$ MeV (the emission line of de-excited $^{12}\mathrm{C}^\ast$) and $1{.}3$ MeV ($^{14}\mathrm{N}^\ast$) by INTEGRAL/SPI: the observed flux implies a LECR density in the Galactic disk of $\sim 5 \times 10^{-3}$ erg/cm$^{3}$, consistent with the standard transport models. The direct measurements of Voyager 1 and 2 (beyond the heliopause, post-2012 and 2018 respectively) have provided the first unmodulated GCR spectra at $\sim 10$ MeV/nucleon, and suggest LECR slightly more intense than expected — with a potential upward revision of the $^{6}\mathrm{Li}$ production at low metallicities.

The contributions to $^{7}\mathrm{Li}$

The case of lithium is particularly complex among the three LiBeB, because it has several simultaneous sources contributing in different proportions at different cosmological epochs. Primordial lithium is produced in the first minutes after the Big Bang by BBN (previous part of this chapter); cosmic lithium is produced by GCR spallation in the interstellar medium over the billions of years of Galactic evolution; stellar lithium is produced in small quantities by some classes of stars (intermediate-mass AGB with active HBB, novae); explosive lithium is produced by the $\nu$-process in core-collapse supernovae through the neutrino spallation of $\alpha$ nuclei. Distinguishing who contributes how much requires looking at $A(\mathrm{Li}) \equiv \log\epsilon(\mathrm{Li}) = \log(N_{\mathrm{Li}}/N_{\mathrm{H}}) + 12$ as a function of the metallicity $[\mathrm{Fe/H}]$ in stars of different ages and environments.

The key observable in stellar spectroscopy is therefore $A(\mathrm{Li})$ vs $[\mathrm{Fe/H}]$, with three distinct regimes. For $[\mathrm{Fe/H}] \lesssim -1{.}5$ in warm MS stars ($T_{\mathrm{eff}} > 5800$ K) one observes the Spite plateau: $A(\mathrm{Li}) \approx 2{.}2$ independent of $[\mathrm{Fe/H}]$, with intrinsic scatter $< 0{.}1$ dex. The standard interpretation is residual primordial $^{7}\mathrm{Li}$ after a small diffusive stellar depletion (see the previous part of this chapter and chapter 8 for the cosmological lithium problem). For $-1{.}5 \lesssim [\mathrm{Fe/H}] \lesssim 0$ one observes an increase of $A(\mathrm{Li})$ with $[\mathrm{Fe/H}]$ up to the meteoritic plateau $A(\mathrm{Li}) \approx 3{.}3$, which corresponds to the Solar System value measured in CI chondrites. The interpretation is the cumulative contribution of GCR + stars + novae above the primordial value. For individual super-Li-rich stars (red giants with $A(\mathrm{Li}) > 4$) the origin is contaminated: HBB in intermediate-mass AGB with efficient dredge-up, or accretion from a binary companion. The BBN prediction combined with Planck gives $A(\mathrm{Li})_{\mathrm{BBN}} \approx 2{.}7$, while the observed plateau sits at $\approx 2{.}2$: the $\sim 0{.}5$ dex difference is the cosmological lithium problem, significantly reduced but not eliminated by stellar depletion (chapter 8).

The contribution of novae to Galactic $^{7}\mathrm{Li}$ was long conjectured theoretically and finally verified observationally in a direct way. In 2015, Tajitsu et al. detected the signature of $^{7}\mathrm{Be}$ (which decays into $^{7}\mathrm{Li}$ via electron capture with $\tau_{1/2} = 53$ d) in the spectra of Nova Cygni 2013, confirming the beryllium transport mechanism of Cameron-Fowler: in a nova explosion, the $^{3}\mathrm{He}$ accreted onto the surface of the white dwarf produces $^{7}\mathrm{Be}$ via $^{3}\mathrm{He}(\alpha,\gamma)^{7}\mathrm{Be}$ at temperatures of a few $10^{8}$ K, and the rapid convection in the explosive envelope transports the $^{7}\mathrm{Be}$ toward the cold surface before it can be destroyed. The later detections in Nova Sgr 2015, Nova Aql 2021 and in other classical novae have systematically confirmed the presence of $^{7}\mathrm{Be}$ in the ejecta, with a mass of $^{7}\mathrm{Li}$ produced per event of $\sim 10^{-8}$-$10^{-7}\,M_\odot$. With a Galactic nova rate of $\sim 30$ per year, one arrives at a significant contribution to the $^{7}\mathrm{Li}$ budget of the Solar System, of order $15\%$.

The $\nu$-process in core-collapse SNe contributes to $^{7}\mathrm{Li}$ and $^{11}\mathrm{B}$ through the neutrino spallation of $\alpha$ nuclei and of CNO nuclei during the passage of the shock through the outer layers of the massive star. The key reactions are $^{4}\mathrm{He}(\nu,\nu' p)^{3}\mathrm{H}$ followed by $^{3}\mathrm{H}(\alpha,\gamma)^{7}\mathrm{Li}$ for lithium, and the direct $^{12}\mathrm{C}(\nu,\nu' p)^{11}\mathrm{B}$ for boron, with neutrino-nucleus cross sections of order $10^{-42}$ cm$^{2}$ at $\langle E_\nu \rangle \sim 10$-$20$ MeV (Woosley & Haxton 1988, Heger et al. 2005 with updated shell model). The contribution to Galactic $^{11}\mathrm{B}$ is estimated at around $50\%$ of the total, and it is the main candidate to explain quantitatively the observed ratio $^{11}\mathrm{B}/^{10}\mathrm{B} \approx 4$ that standard GCR spallation under-predicts (kinematic ratio $\sim 2{.}5$). For $^{7}\mathrm{Li}$ the $\nu$-process contribution is of order $5\%$ of the Galactic budget.

Residual uncertainties and prospects

The LiBeB story is today reasonably closed in its main lines: GCR spallation + BBN + minor contributions from novae and the $\nu$-process quantitatively explain the cosmic abundance curve of the three elements with residual uncertainties of a factor of 2 at most. The Galactic LiBeB budget estimated by Prantzos [Prantzos 2012] is summarized in the following table, with $\pm 30\%$ uncertainties on the individual values:

Two main quantitative problems remain open. The first is the cosmological lithium problem: BBN combined with the Planck parameters predicts about $0{.}5$ dex more $^{7}\mathrm{Li}$ than observed in the Spite-plateau stars of the Galactic halo — a discrepancy that has persisted for twenty years. The emerging consensus solution is the slow stellar depletion in Pop II stars through turbulence-regulated diffusion (Korn et al. 2007 [Korn et al. 2007] in NGC 6397, Richard et al. 2005, Mucciarelli et al. 2022 [Mucciarelli et al. 2022] ), which reduces the photospheric abundance by $\sim 0{.}4$ dex on Gyr timescales — enough to reconcile the BBN prediction with the plateau within the residual systematic uncertainty of $\sim 0{.}1$-$0{.}2$ dex. The discrepancy has shrunk significantly but is not yet considered definitively resolved (chapter 8). The second problem is the precise quantification of the contribution of novae to the Galactic $^{7}\mathrm{Li}$ budget, which depends critically on the nova rate as a function of metallicity — a quantity measured today with $30$-$50\%$ uncertainties, being reduced thanks to the new-generation time-domain surveys (Vera Rubin LSST).

The open frontiers of the field include several converging directions. The direct measurement of the GCR flux below the solar modulation with experiments such as Voyager 1/2 (beyond the heliopause) has already provided the first unmodulated spectra, and suggests LECR slightly more intense than expected, with a potential upward revision of the $^{6}\mathrm{Li}$ production at low metallicity. The identification of the $^{6}\mathrm{Li}/^{7}\mathrm{Li}$ ratio in Spite-plateau stars has been the object of controversy: some early-2000s studies (Asplund 2006) found it high ($\sim 0{.}05$), suggesting a mysterious pre-Galactic $^{6}\mathrm{Li}$ population, but later 3D NLTE analyses (Lind et al. 2013) significantly reduced the values, showing that the $^{6}\mathrm{Li}$ lines were artifacts of the 1D treatment of turbulent mixing and eliminating the so-called “second lithium problem”. The contribution of dark matter decay to the primordial LiBeB budget, explored in classes of SUSY models with a metastable neutralino that predicted destruction/production of Li-Be-B in BBN, has been stringently constrained and excluded for most parameter values by Pospelov & Pradler (2010) and later works. On the experimental front, the LECR measurement campaigns on fragmentation cross sections at low energies (HIE-ISOLDE, NICA) and the laboratory measurements of the rate $^{7}\mathrm{Be}(n, \alpha)^{4}\mathrm{He}$ at n_TOF at CERN (relevant for BBN) will further reduce the uncertainties on the primordial and Galactic yields over the next five years. The state-of-the-art synthesis of modern BBN incorporating these new nuclear data is Cyburt, Fields, Olive and Yeh (2016) [Cyburt et al. 2016] .

The origin of Li, Be and B closes the picture of the nucleosynthesis of the light elements: together with primordial BBN, it covers the region $Z = 1$-$10$ of the periodic table with a mix of cosmological sites, Galactic sites external to the stars (GCR), and residual stellar ones (novae, $\nu$-process). What remains to be told in detail is the physics of stellar evolution and of the quiescent burnings from hydrogen to silicon, the subject of the next chapter.