The electromagnetic fine-structure constant in primordial nucleosynthesis revisited

We study the dependence of the primordial nuclear abundances as a function of the electromagnetic fine-structure constant α\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\alpha $$\end{document}, keeping all other fundamental constants fixed. We update the leading nuclear reaction rates, in particular the electromagnetic contribution to the neutron-proton mass difference pertinent to β\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\beta $$\end{document}-decays, and go beyond certain approximations made in the literature. In particular, we include the temperature-dependence of the leading nuclear reactions rates and assess the systematic uncertainties by using four different publicly available codes for Big Bang nucleosynthesis. Disregarding the unsolved so-called lithium-problem, we find that the current values for the observationally based 2H\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${}^{2}\text {H}$$\end{document} and 4He\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$${}^{4}\text {He}$$\end{document} abundances restrict the fractional change in the fine-structure constant to less than 2%\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$2\%$$\end{document} , which is a tighter bound than found in earlier works on the subject.


I. INTRODUCTION
Since the early work of Dirac [1], variations of the fundamental constants of physics have been considered in a variety of scenarios, see [2] for a recent status report on possible spatial and temporal variations of the electromagnetic fine-structure constant α, the gravitational constant G and the proton magnetic moment µ p .See also the recent review [3].
As is well known, primordial or Big Bang nucleosynthesis (BBN) is a fine laboratory to test our understanding of the fundamental physics describing the generation of the light elements.In particular, it sets bounds on the possible variation of the parameters of the Standard Model of particle physics as well as the Standard Model of cosmology (ΛCDM).For recent reviews, see e.g.Refs.[4][5][6].Here, we are interested in bounds on the electromagnetic fine-structure constant α derived from the element abundances in primordial nucleosynthesis.For earlier work on this topic, see e.g.[7][8][9][10] and references therein.This work is part of a larger program that tries to map out the habitable universes in the sense that the pertinent nuclei needed to generate life as we know it are produced in the Big Bang and in stars in a sufficient amount, see e.g.[11,12] for reviews.
Here, we focus largely on the nuclear and particle physics underlying the element generation in primordial nucleosynthesis.In particular, we reassess the dependence of the nuclear reactions rates on the fine-structure constant, overcoming on one side certain approximations made in the literature and on the other side providing new and improved parameterizations for the most important reactions in the reaction network, using modern determinations of the ingredients whenever possible, such as the Effective Field Theory (EFT) description of the leading nuclear reaction n + p → d + γ 1 and the calculation of the nuclear Coulomb energies based on Nuclear Lattice Effective Field Theory.For β-decays, we also use up-to-date information on the neutron-proton mass difference based on dispersion relations (Cottingham sum rule).Most importantly, as already done in Ref. [19], we utilize four different publicly available codes for BBN [20][21][22][23][24][25][26] to address the systematic uncertainties related to the modeling of the BBN network.In particular these codes differ in the number of nuclei and reactions taken into account as well as in the specific parameterization of the nuclear rates entering the coupled rate equations for the BBN network.Moreover, in determining the sensitivity of primordial abundances on nuclear parameters, we account for the temperature dependence of the variation of some rates on the value of the fine-structure constant α.To our knowledge, such a comparative study where this temperature dependence was explicitly considered has not been published before.We further note that we keep all other constants, like e.g. the light quark masses m u , m d fixed at their physical values.
The paper is organized as follows: In Sect.II we collect the basic formulas needed for discussing the finestructure constant dependence in BBN.In this section we discuss the various dependences of the reaction rates 1 There are some ab initio calculations of other reactions in the BBN network such as [13][14][15][16][17][18], mainly concerned with radiative capture reactions.The calculations in the framework of so-called "halo-EFT" potentially offer the possibility to study the α dependence of the cross sections analytically, but the implementation is numerically rather involved and we thus refrained from doing so in the present context.
on the value of the electromagnetic fine-structure constant α.The actual calculation of the reaction rates is treated in Sect.III.The BBN response matrix is introduced in Sect.IV.The numerical results of this study are presented in Sect.V and discussed in Sect.VI.We also present a detailed comparison to results obtained in earlier works.In Appendix A we give the novel parameterizations of 18 leading nuclear reactions in the BBN network.

II. BASIC FORMALISM
As discussed in Ref. [19], the basic quantities to be determined in BBN are the nuclear abundances Y ni , where n i denotes some nuclide.The evolution of the nuclear abundance Y n1 is then generically given by Ẏn1 = n 2 , . . ., np m 1 , . . ., mq N n1 Γ m1,...,mq→n1,...np Y where the dot denotes the time derivative in a comoving frame, and N ni is the stochiometric coefficient of species n i in the reaction.Further, for a two-particle reaction a + b → c + d , Γ ab→cd = n B γ ab→cd is the reaction rate with n B the baryon volume density.This can readily be generalised to reactions involving more (or less) particles, see [26].These equations are coupled via the corresponding energy densities to the standard Friedmann equation describing the cosmological expansion in the early universe, for details and basic assumptions, see also [22,25,26].In what follows, we discuss the various types of reactions in the BBN network and their dependence on the electromagnetic fine-structure constant.

A. Reaction rates
The average reaction rate γ ab→cd = N A ⟨σ ab→cd v⟩ for a two-particle reaction a+b → c+d is obtained by folding the cross section σ ab→cd (E) with the Maxwell-Boltzmann velocity distribution in thermal equilibrium (2) conventionally multiplied by Avogadro's number N A , where µ ij is the reduced mass of the nuclide pair ij, µ ij = m i m j /(m i + m j ), E is the kinetic energy in the center-of-mass system (CMS), T is the temperature and k the Boltzmann constant.Defining y = E/(kT ) this can be written in the form γ ab→cd (T ) = N A 8 kT π µ ab ∞ 0 dy σ ab→cd (kT y) y e −y .
(3) This is suited for numerical computation e.g. with a Gauß-Laguerre integrator.In fact, in order to deal with cases with singular cross sections for E → 0 it is even better to split the integral and write and evaluate the first integral with a Gauß-Legendre and the second with a Gauß-Laguerre integrator for some suitable value of y.Note that in the first term the substitution x = √ y was performed.
With the detailed balance relation where are the CMS kinetic energies in the entrance and exit channels, respectively, and g i is the spin multiplicity of particle i, energy conservation implies (7) in terms of the Q-value for the forward reaction.In thermal equilibrium the inverse reaction rate is then related to the forward rate as This brings us to the central question of this paper, namely how the value of the electromagnetic finestructure constant influences the reaction rates?This clearly depends on the reaction type.With the exception of the leading n + p → d + γ nuclear reaction, to be discussed in some detail below, no ab initio expressions for most of the reaction cross sections is available and accordingly one has to rely on model assumptions concerning the fine-structure constant dependence of the cross sections and thus of the reaction rates (see also the discussion in Sect.VI on this issue).These shall be discussed in the following subsections separately for direct reactions of the type radiative capture reactions of the type and β-decays, We shall start with a brief discussion of the Coulomb penetration factor for charged particles, relevant for what follows.

Coulomb-penetration factor
The Coulomb-penetration factor for an l-wave is given by, see e.g.[27,28] , where F ℓ , G ℓ are the regular and irregular Coulomb functions, respectively, that are the linearly independent solutions of the radial Schrödinger equation where we defined for the Coulomb-scattering of charges Z a e , Z b e with masses m a , m b and the reduced mass at the energy E of the relative motion, subject to the Coulomb-potential Approximate parameterizations of v ℓ (η, ρ) have been extensively discussed in the literature, see e.g.[29] in particular for the dependence on the nuclear distance r where this is to be evaluated for a specific reaction.As argued in [27,28], this distance is not well defined and the cross section should not depend on such an unobservable parameter.Accordingly, if one takes, as in [27,28], lim ρ → 0 the penetration factor for an ℓ-wave then reads Therefore, we shall use as Coulomb penetration factor the expression for s-waves: Note that the corrections due to a variation of α in ε 2 ℓ for ℓ > 0 according to Eq. ( 18) are of higher order in α and thus small anyway .Here we defined in terms of the so-called Gamow energy for a two-particle reaction channel ij and the CMS energy E or E + Q for the entrance and exit channel, respectively.
For a direct reaction of this type the Q-value is given by where the nuclear mass of each nuclide i with mass number A i and charge number Z i reads with B i the nuclear binding energy.Thus, because of baryon number and charge conservation ( Now the binding energy can be written as where B N i denotes the strong contribution to the binding energy and is the expectation value of the Coulomb contribution proportional to the value of the electromagnetic finestructure constant.Considering its variation in the form α = α 0 (1 + δ α ) , where is the current experimental value from Ref. [30], the Qvalue varies as One therefore needs an estimate of the Coulomb contribution to the nuclear masses, this we shall discuss in Sect.II G.We shall assume that the cross section for a direct reaction a + b → c + d depends on α as where f is some function independent of α and P i (x i ), P f (x f ) are the penetration factors given by Eq. ( 19) reflecting the Coulomb repulsion in the entrance and in the exit channel, respectively.The first factor in Eq. ( 29) accounts for the exit channel momentum dependence of the cross section of the direct reaction a + b → c + d.Here, are the arguments of the penetration factors, with the Gamow-energies in the entrance and the exit channel, respectively, and µ ij = m i m j /(m i + m j ) the corresponding reduced masses.Although in order to calculate the linear response of the abundances one could proceed by calculating first order partial derivatives etc. , we prefer not to presume linearity and rather calculate a variation of the cross section with a variation α = α 0 (1 + δ α ) through the expression where specifically the first factor reads and the remaining factors are given by We note that these factors are energy-dependent and therefore the change in the rate i.e. the factor depends on the temperature T and as it stands requires a numerical evaluation of Eq. (38) .
Similar considerations hold for radiative capture reactions.The cross section of a reaction a + b → c + γ is assumed to depend on α as with f some α-independent function and P i (x i ) the penetration factor, see Eqs. (19,30,32), for the entrance channel.The first factor accounts for the fact that in the amplitude for a radiative capture reaction the photon coupling is proportional to e, leading to a factor proportional to α = e 2 /(ℏc) in the cross section.The second factor reflects the final momentum dependence assuming dipole dominance of the radiation 2 .We thus calculate a variation of the cross section for radiative capture with a variation α = α 0 (1 + δ α ) via 2 Note, however, that this is not always the case, exceptions with appreciable E2 (electric quadrupole) contributions are e.g. the reactions: where the first factor is the same as in Eq. (36) .Again note that both factors are energy-dependent and therefore a change in the rate, see Eq. (38), is temperaturedependent.

D. Approximate treatment of α-dependent factors
As mentioned twice, the variation of the cross sections with a variation of α induces energy-dependent factors, that in turn lead to temperature-dependent variations in the corresponding reaction rates, that can be fully accounted for only via a numerical integration of Eq. (38).In fact this is what was done in the present work for the most important 18 nuclear reactions in the BBN network, listed in Sect.V.For the remaining reactions we relied on the following approximations, that turned out to be effective.For neutron induced reactions where for a non-resonant reaction R(E) is a weakly dependent function of the CMS kinetic energy E, see e.g.[31] .If we make the extreme approximation that R(E) ≈ const.the maximum of the remaining energy dependent factors in the integrand of Eq. ( 38) is reached at the energy Likewise, assuming that for the astrophysical S-factor for charged particle induced reactions holds, one finds that the maximum of the remaining energy-dependent factors in the rate is reached at Substituting E → E n , E c in the expressions in Eqs.(35,41) then leads to temperature-dependent factors, that can be taken in front of the integral in Eq. ( 38) and thus merely multiply the corresponding rates.The quality of this approximation may be inferred from Fig. 1, where we compare the results of the numerical calculation of the rates according to Eq. (38) (yellow areas for a variation δ α ∈ [−0.1, 0.1]) with the approximation discussed in this subsection, represented by blue lines for δ α = −0.1 and 0.1 .Note that in the present treatment we preferred to account for the Coulomb suppression in an entrance or exit channel with charged particles by the penetration factor of Eq. ( 19) and do not rely on a simple Gamowfactor ∝ e −x = e − √ E G (α)/E , with E being the CMS energy of the relevant channel with charged particles.We found that doing so would lead to overestimating the αdependence in the rates by a factor ≈ 1.5, while the temperature dependence would still roughly follow the same trends as in Fig. 1.

E. Coulomb-effects in β-decays
Next, we consider the various β-decays in the BBN network.The rate for β-decays a → c + e ± + (−) ν can be written as, see e.g.[32], where G is Fermi's weak coupling constant, M ac the nuclear matrix element and see Eq. (2.158) of Ref. [33] , where we defined q = Q/m e = (m a − m c )/m e .Further, F (±Z, E) is the socalled Fermi-function with the definitions where Z is the atomic number of the daughter nucleus c and R its radius, see Eqs. (2.121)-(2.125),(2.131) of Ref. [33] .The upper/lower sign holds for β − /β + decays, respectively.For Z α ≪ 1 we can approximate or, with such that lim Z→0 F 0 (±Z, E) = 1.Accordingly, setting a = 2π Z α then  Defining also p = q 2 − 1 (i.e. the maximal momentum in β-decay divided by the electron mass) and with the substitution y = √ x 2 − 1/p we can rewrite the expression for f (±Z, q) as which is slightly better suited for a numerical implementation, e.g with a Gauß-Legendre-integrator.Note that both a and p depend on α.

Electromagnetic contribution to the proton-neutron mass difference
The neutron-proton mass difference plays an important role in BBN, see e.g.[41].According to Refs.[34,35] the proton-neutron mass difference is given by where the nominal electromagnetic contribution is somewhat smaller than the value ∆m QED = 0.7 ± 0.3 MeV given earlier in Ref. [36].
We also note that the splitting in strong and electromagnetic contributions is convention-dependent, for a pedagogical discussion see [37].For a comparison of these results with lattice QCD and other phenomenological determinations of the electromagnetic contribution to the neutron-proton mass difference, we refer to [34].
The neutron-proton mass difference is a crucial parameter both in the various n ↔ p weak interactions in the early phase of BBN and in all β-decays.We shall start with a discussion of the latter.

Implications for β-decays
Writing the Q-value for the β-decay a → c + e ∓ + (−) ν e depends on a variation of α = α 0 (1 + δ α ) as where, as in Sect.II B, B N i is the nuclear (strong) contribution to the binding energy of nuclide i and V C i the expectation value of the Coulomb-interaction to the binding energy of nuclide i .
One thus finds for the variation of the β-decay rate with a variation of α, where and the factor determining the variation of the β-decay rate with a variation of α , see Eq. ( 57), is determined by evaluating f (ã, p) and f (a, p) via Eq.( 53) .We note that ≥ m e implies an upper limit for δ α : As for other cases where during a variation of α the Qvalue of a reaction becomes negative, we have put the corresponding rate to zero.We also note that for the neutron decay n → p+e − +ν the variation of the rate with a variation of α merely implies a variation of the neutron lifetime τ n ∝ 1/λ n→p .

Implications for the weak n ↔ p reactions
As detailed in Ref. [26] the six reactions determine the evolution of the neutron abundance in the early phase of BBN and hence are crucial for all other primordial nuclear abundances.Assuming local thermodynamical equilibrium in terms of a temperature T and a distinct neutrino temperature T ν in the so-called infinite nucleon mass approximation the n → p (angular averaged) reaction rate can be written, see e.g.[26] for details, as with the Fermi-Dirac distribution function.The ratio T ν /T follows from the cosmological evolution, see the black 2. Variation of the rate ratios Rn→p (Eq.( 63), blue hatched area) and Rp→n (Eq.( 64), red hatched area) with decreasing temperature (T9 = T /[10 9 K]) for δα in the range δα = −0.1 (lower curves) up to δα = 0.1 (upper curves).Also shown is the ratio Tν /T (black curve).curve in Fig. 2. The constant K is fixed by requiring that Γ n→p (∆m; 0) = 1/τ n , with τ n the neutron lifetime.The p → n rate is simply given by substituting ∆m → −∆m in Eq. ( 61) above.In this case Γ p→n (∆m; 0) = 0 .As discussed in [26] and [31] there are a number of corrections to the n → p and p → n rates as given above, viz. the Coulomb correction (as discussed above in section II E 2), electromagnetic radiative corrections, finite nucleon mass corrections, plasma corrections and non-instantaneous neutrino decoupling effects.Some of these involve the fine-structure constant α, but since these effects are small corrections anyway, the most relevant effect when varying α is through the change ∆m(α) = ∆m(α 0 ) − ∆m QED δ α .This effect is illustrated in Fig. 2, where the double rate ratios and obtained by a numerical integration according to Eq. ( 61) (with a method similar to that of Eq. ( 4)) are plotted as a function of T 9 = T /[10 9 K]).This double ratio was chosen such that for T → 0 the n → p curves tend to unity and the p → n curves to zero; the α-dependence of the n → p rate in this low-temperature limit is then given by the expressions in the preceding section II E 2.
As is evident from this figure the variation of the rates with varying α is non-linear and strongly temperature dependent.
F. The n + p → d + γ reaction Fortunately, for the n + p → d + γ reaction an accurate treatment within the framework of pionless EFT [38,39] is available.Accordingly, for this leading nuclear reaction in the BBN network it is possible to study dependences of the cross section and hence of the reaction rate on various nuclear parameters, such as the binding energy of the deuteron, np scattering lengths, effective ranges etc. as was done in [19].Here we shall focus on the α-dependence.
The cross section was given in [39] as where p is the relative momentum, m N = (m p + m n )/2 denotes the nucleon mass, γ = √ B d m N is the so-called binding momentum, with B d = 2.225 MeV the binding energy of the deuteron, and , are the dimensionless amplitudes for isovector electric dipole, isovector magnetic dipole, isoscalar magnetic dipole and isoscalar electric quadrupole contributions, respectively.For the energies relevant in BBN only the isovector contributions are significant and these were calculated at N4LO and N2LO for the electric and magnetic parts, respectively.The overall theoretical uncertainty is claimed to be better than 1% for CMS energies E ≤ 1 MeV .The expression of Eq. ( 65) with all terms included was used to calculate the cross section for this reaction throughout the present investigation.
Concerning the variation of this cross section when varying α = α 0 (1 + δ α ) it is evident that the dominant effect is simply Note that there is no Coulomb-contribution to the binding energy of the deuteron, while the expectation value ⟨v EM ⟩ of the electromagnetic interaction, mainly due to the magnetic dipole-dipole interaction moment term, see [40] for a treatment based on the Argonne v 18 nucleon-nucleon potential, is very small, ⟨v EM ⟩ = 0.018 MeV.Hence the effects of a change of the Q-value of the reaction with varying α, as discussed in the previous subsections, are considered to be negligible in the present context.Moreover, in the expression of Eq.( 65), as well as in the expressions for the amplitudes χ of Ref. [39] the nucleon mass m N = (m p + m n )/2 occurs at various instances.Although a moderately accurate value for the electromagnetic contribution to the neutronproton mass difference is available (and was, in fact, used in our discussion of the β-decays in Sect.II E 1), only rough estimates are available for the electromagnetic contribution to the neutron and proton mass separately.In Eq. ( 12.3) of Ref. [42] the estimates m Born p ≈ 0.63 MeV, m Born n ≈ −0.13 MeV (with an estimated accuracy of ≈ 0.3 MeV) for the total electromagnetic self-energy can be found, which, via m Born N ≈ 0.25 MeV , would imply that m N varies with α (putting m N ≈ 1 GeV for this estimate) as For |δ α | < 0.1 , as considered here, this would lead to effects well below the theoretical accuracy quoted above and therefore this effect was neglected and the variation of the n + p → d + γ cross section with α is assumed to be entirely given by Eq. ( 66) .

G. Coulomb energies
A variation of the value of the fine-structure constant α implies a variation of the nuclear binding energies and hence a variation of the Q-values of the reactions, which in turn leads to a variation of the cross sections and the corresponding rates.Therefore the present study requires an estimate of the electromagnetic contribution to the nuclear masses or equivalently to the nuclear binding energies.A rough estimate is provided by the Coulomb term in the Bethe-Weiszsäcker formula (for a recent determination, see e.g.[43] and references therein): approximately accounting for the Coulomb repulsion by the protons in a nucleus.However, this formula is not very precise when applied to the light nuclei relevant here.We therefore prefer to use the expectation values of the Coulomb interaction as determined from a recent ab initio calculation of light nuclear masses in the framework of Nuclear Lattice Effective Field Theory (NLEFT) [44] , listed in Table I .We also compare the calculated binding energies to the experimental data as used here in order to give an impression of the quality of the calculation.

III. CALCULATION OF THE REACTION RATES
For the 18   a from [44].The errors quoted in parentheses include all the statistic and systematic uncertainties.In case of 2 H, the error is entirely given by the variation of the np phase shifts at N3LO within their uncertainties.b from [45], as used in the present work.c extrapolated from a least-squares fit to the other data with the rates and their variations with α are calculated by a numerical integration of Eq. ( 38) and tabulated for 60 temperatures in the range 0.001 ≤ T 9 = T /[10 9 K] ≤ 10.0.These values were then used via a cubic spline interpolation in the four publicly available BBN codes as outlined in Sect.IV.The resulting rates and their variations with α = α 0 (1 + δ α ) in the range δ α ∈ [−0.1, 0.1] are displayed in Fig. 1 in Sect.II D. To this end, we made new fits to the cross sections (or equivalently of the corresponding astrophysical S-factors) for the reactions listed above.The parameterizations can be found in Appendix A. In addition in Fig. 3 the resulting reaction rates for α = α 0 are compared to the rates implemented in the original versions of the four programmes considered here.In Fig. 3 we also display the rates obtained with the NACRE II database, see [46], which served as a further check on our calculated reaction rates at α = α 0 .
The rates of all other reactions were taken as in the original implementation of the codes and the variation of the rates with α was calculated as discussed in Sect.II D.
The variation of the β-decay rates according to Eq. ( 57) was implemented directly in the various codes.In Fig. 4 it is shown how the β-rates at low temperature (i.e.T ≪ T 9 ) change by a variation of α = α 0 (1 + δ α ).In particular the rates of the tritium decay and the 14 Cdecay strongly depend on the value of δ α , the effect of (relatively large) changes in the (relatively small) Qvalues due to changes in the Coulomb contribution to the binding energies being dominant.
As already touched upon in section II E 3 the variation of the weak n ↔ p rates with α is dominated by the variation of proton-neutron mass difference with α and is strongly temperature dependent in the early phase of BBN.In the default version of the Kawano code NUC123 [20] this temperature dependence is parameterized as outlined in Appendix F of Ref. [20], but a numerical integration along Eq. ( 61) can be enforced and was in fact used to implement the α-dependence of these rates.The PArthENoPE code [23][24][25] contains a slightly more sophisticated parameterization, see e.g.Appendix C of Ref. [31], accounting also for some higher order corrections.Here we used the α-dependence of the n ↔ p rates as illustrated in Fig. 2 as a factor multiplying the parameterized rate.In the AlterBBN code [21,22] the temperature dependence of the weak n ↔ p rates was already determined numerically as in Eq. ( 61) and the α dependence can be accounted for by an appropriate variation of ∆m .In this code also the Coulomb correction, see Eq. ( 57), was included in the integrand of Eq. ( 61), but this was found to have no significant impact on the final abundances to be discussed below in section V.The PRIMAT [26] implementation offers the possibility to study the α-dependence of the weak n ↔ p reactions in all detail including all the higher order electromagnetic corrections mentioned in section II E 3 .In fact this code was used to verify that the variation of the rates through the variation of ∆m with α as discussed in section II E 3 is indeed the dominant effect.Indeed, ignoring the α dependence in the higher order corrections implemented in PRIMAT led to response coefficients that differ at most by 0.5% from the values listed in Table III below.Accordingly, in spite of the fact that the n ↔ p reactions are treated at various levels of sophistication, the resulting primordial abundances and their variation with α, to be discussed in section V, were found to be rather consistent.

IV. THE BBN RESPONSE MATRIX
We estimated the linear dependence of the primordial abundances Y n on small changes in the value of the finestructure constant α = α 0 (1 + δ α ) by calculating the abundance of the nuclide n, with n ∈ { 2 H , 3 H + 3 He , 4   C) as well as an implementation as a mathematica-notebook, PRIMAT [26] .To this end we performed least-squares fits of a quadratic polynomial to the abundances: such that will be called an element of the linear nuclear BBN response matrix.It represents the dimensionless fractional change in the primordial abundance ratio Y n /Y H due to a fractional change α in linear approximation.Deviations from a linear response are then given by the coefficient c 2 .

V. RESULTS AND DISCUSSION
In most of what follows we shall use η = 6.14 • 10 −10 from Ref. [30] as the nominal baryon-to-photon density ratio while varying α.The programs were modified as indicated in Sect. 4 of Ref. [19] and the rates for the most relevant reactions listed in Sect.III, resulting from the new fits of the cross sections presented in Appendix A, were used in all programmes.
The resulting nominal (i.e. at α = α 0 ) abundances at the end of the BBN epoch in terms of the number ratios Y2 H /Y H , Y3 H+ 3 He /Y H , Y6 Li /Y H , Y7 Li+ 7 Be /Y H , and the mass ratio for 4 He are compared to the values quoted in Ref. [19] and experimental data in Table II .Although the mass ratio for 4 He and, to a minor extend, the number ratio for deuterium did not change significantly with respect to the values obtained in Ref [19], the 3 H + 3 He number ratio increased by approximately 10% , the 6 Li number ratio was found to be larger by about 70 − 80% , while the 7     3. Reaction rates γ(T9) for 18 leading nuclear reactions in the BBN network, where T9 = T /[10 9 K] .The rates resulting from the new parameterizations of the S-factors in Appendix A are represented by solid red curves (color online).The rates in the original version of the programmes are given by green curves for NUC123 [20], magenta curves for PArthENoPE [25], blue curves for AlterBBN [22] and cyan curves for the PRIMAT [26] code.Also shown as a thin black curve is the result from the NACRE II database, see [46] .a factor of three, a phenomenon known as the lithiumproblem, which is thus unsolved even using the updated cross sections used here.As stated previously in [19], in spite of this unresolved issue in BBN the consistency of the cosmic microwave background observations with the determined abundances of deuterium and helium is considered to be a non-trivial success.Accordingly, we think that this issue is no obstacle for the study presented here.
The elements of the response matrix were then determined by a polynomial fit, as explained above in Sect.IV for the abundances relative to the hydrogen abundance, namely Y2 H /Y H , Y3 H+ 3 He /Y H , Y6 Li /Y H , Y7 Li+ 7 Be /Y H , and the mass ratio for 4 He .
The dependence of these ratios on the value of the fine structure constant α = α 0 (1 + δ α ) is displayed in Fig. 5  abundance ratios is found to be very similar for all four publicly available codes, in spite of the fact that these codes differ in details, such as the number of reactions in the BBN network or the manner in which the rate equations are solved numerically.Note, however, that in the present study the rates calculated for the major reactions listed in Sect.III and their variation with α are the same.
Of course this then also applies to the values for the resulting response matrix elements.The response matrix elements ∂ log (Y n /Y H )/∂ log α = c 1 and the coefficients of the quadratic term in Eq. ( 73) are given and compared to some results from the literature in Table III.Note that with the exception of 6 Li, we have |c 2 | ≃ |c 1 |, so that due to the smallness of α, the second order contribution to the response is of minor importance.All programs were run with the full network implemented in the original version codes.We checked that if we run the programs with a smaller network the results listed in Tables II,III change only in the last digit and therefore conclude that the approximation, see Sect.II D, we made for rate changes in the reactions beyond the reactions listed in Eqs.(69)(70)(71) are without any effect for the present investigation.
Apart from the values of c 1 for 2 H(≈ 3.6) and for 6 Li(≈ 6.8) the values obtained in the present study, although consistent among each other, differ appreciably from the values obtained in Refs.[7,8,10].In particular in the present calculations the linear response for 3 H + 3 He is much larger while the linear response for 7 Li + 7 Be is appreciably smaller in magnitude, although there seems to be at least a consensus concerning the sign.
In order to clarify this issue, we shall discuss in some detail the relevance of the various factors that reflect the α-dependence of the nuclear rates: • First of all we list in Table IV the linear response of the BBN abundances to a variation of α in the β-decay rates only.
• In Table V we display the linear response of the BBN abundances to a variation of the nuclear reaction rates.The relevance of the variation of the binding energies with α may be appreciated by the linear response due to changes in α accounting for the effects due to the Coulomb penetration factors only, i.e. without accounting for changes in the binding energies, listed in Table VI .Here we also compared our results to the results presented in Table I of Ref. [10] for the dependence of the abundances on the nuclear rate variation with α, that thus differ from our results significantly for c 1 ( 3 H + 3 He) and c 1 ( 7 Li + 7 Be) , our results being larger in magnitude for the former and smaller for the latter.
Indeed, if we substitute our values, as well as the results we obtained in [19] for the linear response of the abundances on binding energies and the neutron life-time τ n for the values of the response matrix C of Table I in [10] and furthermore account for the smaller value ∂ log τ n /∂ log α ≈ 2.90, obtained via Eq.( 57) (instead of 3.86 in [10]) and the smaller value ∂ log Q N /∂ log α ≈ −0.45 (instead of −0.59 in [10]), due to the smaller new value for ∆m QED , and use our values for the response of the binding energies ∂ log B i /∂ log α that are smaller by about 10% in Table IV of [10] we find approximately for the linear responses , AlterBBN [22], PArthENoPE [25], PRIMAT [26] .Here, we use η = 6.14 • 10 −10 and τn = 879.4s .Also shown are the solid curves obtained by the fits according to Eq. ( 73) with the parameters listed in Table III.The experimental values cited in PDG [30] (thick red lines) are indicated by yellow-highlighted regions (color online) representing the 1σ limits by red lines.
rather close to our values for c 1 given in Table III.Most of the effects listed above are, although significant, of minor importance only, and accordingly the difference can be traced back to the fact that our results for the variation of the rates with a variation of α when ignoring the effects based on Q-value changes, as listed in Table VI differ appreciably from those of [10] .Unfortunately in the latter reference no results on the α-dependence of the rates are explicitly given.In appendix A.2 of [10] it is mentioned that parameterizations of the S-factors were used, the parameters determined by fitting the NETGEN rates as closely as possible.In order to check our parameterizations of the nuclear rates we compared our rates with results generated by the NETGEN-tool [46] in Fig. 3 and found that these are indeed compatible for all reactions, except for the reaction 7 Be + n → 4 He + 4 He, where the NETGEN-tool merely uses the THALYS nuclear reaction model [47].We instead used data, see also Fig. 26 for our fit of the S-factor.Therefore the difference must be due to the different way the Coulomb penetration effects are  73) at η = 6.14 • 10 −10 and τn = 879.4s .Yn/YH are the number ratios of the abundances relative to hydrogen; Yp is conventionally the 4 He/H mass ratio.The results obtained with the four BBN codes NUC123 [20], PArthENoPE [25], AlterBBN [22], PRIMAT [26] are given in four subsequent rows and compared to earlier results from Refs.[7,8,10].treated.Note that, as emphasized in Sect.II A 1, we did not rely on temperature-independent penetration factors taken as a Gamow-factor, but rather accounted for the penetration dependences in the cross-section, which then leads to temperature-dependent changes in the rates.
Our results also differ from the results in Refs.[8] and [7] published even earlier.Concerning the treatment in [8], it is noted that, although the authors present a detailed discussion of the α-dependence in the penetration factors, even accounting for additional αdependent effects due to the peripheral nature of some radiative capture reactions such as e.g. the 3 He + 4 He → 7 Be + γ , an effect taken into account also in the present treatment.Nevertheless, in contrast to our treatment, α-dependent effects seem to be treated merely by temperature-independent factors multiplying the rates.In Ref. [7], the changes in the reaction rates due to changes in α were treated through approximate expressions based on expansions of the S-factors, whereas we preferred to make no further approximations beyond the modeling of the penetration factors discussed in Sect.II A 1. Note that a comparison with the work of [9] is not possible, since there any variation of the fine-structure constant is tied to the variation of certain Yukawa couplings.
All in all our results indicate that the BBN abundance for 7 Li + 7 Be is less sensitive and the abundance of 3 H+ 3 He is more sensitive to variations of the value of the electromagnetic fine-structure constant α than what was determined earlier.Note that such a reduced sensitivity on nuclear quantities, such as binding energies etc. was also observed in [19] .There it was also found that this is mainly due to inclusion of the temperature-dependent changes in the rates.Unfortunately, the primordial abundance of 3 H + 3 He is not known precisely enough to lead to any implications and the nominal prediction for the 7 Li + 7 Be abundance is too large anyway.
If we focus on the deuterium and 4 He abundance ratios alone we can extract from the observationally based data bounds on the variation δ α of the value of the finestructure constant as listed in Table VII for the four programs considered here, showing that one can allow for a variation of the fine-structure constant α by less than 2% on the basis of the results obtained with all programs considered here using the current value for the baryonto-photon ratio η = 6.14•10 −10 , given in [30] .The values ) limits for the variation δα of the fine-structure constant α = α0 (1 + δα) determined such that the resulting abundance lies within the error bounds of the observationally based abundance ratios for 2 H and 4 He given in [30].for the 4 He mass ratio Y p obtained with all four programs are rather consistent and the range [−0.018, 0.006] which is more restrictive than the rough estimate |δ α | < 0.1 quoted in [7,8] and the limit |δ α | ≤ 0.019 mentioned in [10].The values found on the basis of the deuterium number ratio show a larger spread, mainly because the nominal values, see Table II, vary more strongly for the four programs.In spite of this we can determine the range [−0.007, 0.008] , also still more restrictive than the (1σ) range [−0.04, 0.10] of [8].Our new restrictions on the variation of α are also stronger than found earlier in the NLEFT analysis of the triple-alpha process in hot, old stars [48,49].From a comparison of Tables IV-VI we also see that the linear response for Y p due to variations in the βdecay is the dominant effect.Indeed, as argued in [7], the variation of Y p with α mainly depends on the variation of the proton-neutron mass difference with α , i.e. on δm QED that enters the n → p weak decay.
As was done previously in Ref. [8] we also studied to what extend the results presently obtained vary with variations of the baryon-to-photon ratio η and found that our results for the linear response coefficients c 1 do not change significantly if η is varied within the error range quoted in [30], η 10 = η • 10 10 = 6.143 ± 0.190 .With the values of the primordial abundance ratios for d, 4 He and 7 Li+ 7 Be mentioned in PDG [30] we can derive parameter ranges for restricting δ α and η as presented in Figs.6-8 .Note that we here allowed for a variation of η well beyond the currently accepted limits quoted in [30] .The results are similar to those obtained in Ref. [8] although the regions of possible values for the δ α -and η-values are narrower here due to the newer, more precise observational data quoted in [30] .The comparison of these results again show that the value of the Li/Be abundance is incompatible with the other data and we therefore refrain from any conclusions concerning possible variations of α on the basis of the 7 Li observation.

VI. SUMMARY
In the present paper we investigated the impact of variations in the value of the fine-structure constant α on the abundances of the light elements, viz.n + p → d + γ reaction in the framework of pionless EFT.For all other reactions of this kind we rely on modifications of experimentally determined reaction cross sections, trying to account for electromagnetic effects, such as penetration factors, modeling the suppression due to the Coulomb barrier in channels involving charged particles as well as changes in the binding energies of nuclides due to the Coulomb repulsion of the protons and hence the Q-values of the nuclear reactions where these are involved.To this end we made new parameterizations of the cross sections of the 18 leading nuclear reactions in the BBN network using current experimental data compiled by EXFOR.We made an assumption about the α dependence of the penetration factors which differs from the Gamow-factor form that was used in previous investigations and used novel estimates for the Coulomb contribution to nuclear binding energies based on a recent ab initio calculation in the framework of NLEFT in order to determine the α-dependence of the nuclear binding energies and the corresponding Q-values.A further new ingredient for studying the α-dependence of the weak β-decays in the BBN network is a novel value for the electromagnetic contribution to the neutron-proton mass difference, which is slightly smaller than what has been used before.All these new inputs were then used to determine the variation of the reaction rates with varying α.Here, we found in particular that the variation of the reaction rates depends on the temperature, a feature that seems to have been ignored in previous investigations.We found consistent results with all four codes mentioned above and hence conclude that the model-dependence concerning the specific treatment of the BBN network is of minor importance for the α-dependence of the primordial abundances studied here.The results for the linear response do, however, deviate significantly from older results, in particular for the α-dependence of the abun-dances of 3 H + 3 He and 7 Li + 7 Be, the former being much larger and the latter much smaller than found previously.Unfortunately, in the standard Big Bang scenario used here, the nominal abundance ratio Y7 Li+ 7 Be /Y H exceeds the current observationally based determination by a factor of three, a feature known as the lithium-problem that is not solved in the present treatment.This then also impedes a determination of consistent bounds on the value of the fine-structure constant from all available primordial abundance data.Using the observations for 2 H and 4 He alone, we can nevertheless state that these data would limit a possible variation of α to |δ α | < 0.02.This is a stronger bound than found earlier in comparable investigations.
An investigation of the kind presented here heavily relies on the modeling of electromagnetic effects in the cross section data (or, equivalently the astrophysical S-factors) of the relevant nuclear reactions in the BBN-network.Here we opted for a specific form of Coulomb penetration factors that differ from Gamowfactors used before and stressed the relevance of the temperature dependence of the variation of α in the reaction rates that resulted from numerically integrating γ(α; T ) ∝ ∞ 0 dE E σ(α; E) exp(E/kT ).It seems that further progress with the purpose of using primordial nucleosynthesis as a laboratory for exploring our understanding of fundamental physics, apart from astrophysical or cosmological aspects will be feasible only if ab initio theories describing the relevant nuclear reactions including electromagnetic effects become available.NLEFT appears to be a promising framework for doing just that, see e.g.Refs.[50,51].
The cross section for the leading nuclear reaction of BBN, namely n + p → d + γ, was calculated according to the formulas , viz.Eqs.(3.3)- (3.16) given in [39] with the parameters quoted there, also see Sect.II F. In Fig. 9 this description is compared to the existing data as compiled in [54] .

Other radiative capture reactions
The parameters found by a fit of the parameters in Eqs.(A2,A3) to the data are displayed in Table VIII for most radiative capture reactions treated here.
The parameterizations are compared to experimental data compiled by EXFOR [54]  FIG. 10.
Fit (red curve, color online) of the S-factor for the d + p → 3 He + γ reaction compared to data compiled by EXFOR [54] .
The only exception is the reaction d + 4 He → 6 Li + γ, where a resonance appears.In this case the S-factor is given by the sum of a cubic polynomial in E and a relativistic Breit-Wigner function: with the parameters listed in Tables VIII and XI.This parameterization is compared to experimental data compiled by EXFOR [54] in Fig. 15 .FIG. 13.Fit (red curve, color online) of the S-factor for the 3 He + 4 He → 7 Be + γ reaction compared to data compiled by EXFOR [54] .

S(E) [MeV mb]
own fit [91] [92] [93] FIG.14. Fit (red curve, color online) of the S-factor for the 6 Li + p → 7 Be + γ reaction compared to data compiled by EXFOR [54] .TABLE IX.Fit parameters of the S-factor, see Eq.(A1), according to Eqs. (A2,A3,A7,A8) and (A9) for charged particle reactions.S0 is given in MeV mb; a k and q k in units of MeV −k .For these reactions q3 = 0 .For this reaction the fit of the S-factor is compared to experimental data compiled by EXFOR [54] in Fig. 24 .

Reaction
For the reaction 7 Be + n → p + 7 Li the following parameterization in terms of a non-relativistic Breit-Wigner function and a polynomial in √ E was used: where again the CMS energy E is given in MeV and the parameters of the Breit-Wigner function can be found in Table XI.For this reaction the S-factor fit is compared to experimental data compiled by EXFOR [54] in Fig. 25 .
Finally, the form of the parameterization for the reaction 7 Be + n → 4 He + 4 He reads The S-factor is compared to experimental data compiled by EXFOR [54] in Fig. 26 .

TABLE I .
Binding energies B (calculated (cal) and experimental (exp) values) and expectation values for the Coulomb interaction V C of light nuclei.

TABLE III .
BBN response matrix c1 = ∂ log (Yn/YH )/∂ log α and the coefficients c2 of the quadratic term in Eq. (

TABLE IV .
BBN response matrix c1 = ∂ log(Yn/YH )/∂ log α accounting for the variation of the β-decay rates only.See also the caption of TableIII.

TABLE V .
BBN response matrix c1 = ∂ log(Yn/YH )/∂ log α accounting for the variation of the nuclear rates only, but also including the variation of the binding energies and thus of the Q-values of the reactions.See also the caption of TableIII.

TABLE VI .
BBN response matrix c1 = ∂ log(Yn/YH )/∂ log α accounting for the variation of the nuclear rates only, but excluding the variation of the binding energies.Also see caption to TableIII.

TABLE VIII .
Fit parameters of the S-factor, see Eq.(A1), according to Eqs. (A2,A3) and Eq.(A6) for radiative capture reactions.S0 is given in MeV mb; a k and q k in units of MeV −k .

TABLE X .
Fit parameters of the function S[n], see Eq.(A10), according to Eqs. (A12,A3) neutron-induced reactions.S is given in MeV 1/2 mb, and a k and q k in units of MeV −k .
a here S[n] 0 in units of MeV mb; to be divided by √ E, E in MeV. in order to yield S[n].
[54]18.Fit (red curve, color online) of the S-factor for the 3 H + d → n + 4 He reaction compared to data compiled by EXFOR[54].