Primordial Nucleosynthesis with Non-Extensive Statistics
C. A. Bertulani Shubhchintak


Article Title]Primordial Nucleosynthesis with Non-Extensive Statistics










[1]C. A. Bertulanicarlos.bertulani@tamuc.edu 

2] Shubhchintakshubhchintak_phy@pbi.ac.in



*[1]Department of Physics and Astronomy, Texas A&M University-Commerce, Commerce, 75429, Texas, USA

[2]Department of Physics, Punjabi University, Patiala, Patiala, 147002, Punjab, India





The conventional Big Bang model successfully anticipates the initial abundances of ^2H(D), ^3He, and ^4He, aligning remarkably well with observational data. However, a persistent challenge arises in the case of ^7Li, where the predicted abundance exceeds observations by a factor of approximately three. Despite numerous efforts employing traditional nuclear physics to address this incongruity over the years, the enigma surrounding the lithium anomaly endures. In this context, we embark on an exploration of Big Bang nucleosynthesis (BBN) of light element abundances  with the application of Tsallis non-extensive statistics. A comparison is made between the outcomes obtained by varying the non-extensive parameter q away from its unity value and both observational data and abundance predictions derived from the conventional big bang model. A good agreement is found for the abundances of ^4He, ^3He and ^7Li, implying that the lithium abundance puzzle might be due to a subtle  fine-tuning of the physics ingredients used to determine the BBN. However, the deuterium abundance deviates from observations.



















[
    [
    24 April 2024
=================





§ BIG BANG NUCLEOSYNTHESIS


While the cosmological big bang model aligns with numerous observations crucial for our comprehension of the universe, the comparison between model-based calculations and actual observations is not straightforward. This complexity arises due to poorly understood time-sensitive effects and systematic errors affecting the data. Nevertheless, the big bang model stands as the primary tool for investigating the physics of the early universe within the time frame of 3-20 minutes. Beyond this point, the decreasing temperature and density prevent nuclear fusion, hindering the formation of elements heavier than beryllium <cit.>.

The model corresponds with the observed cosmic microwave background (CMB) radiation temperature of 2.275 K <cit.>, offering insights into various scientific realms, including nuclear and particle physics. Additionally, calculations based on the big bang model remain consistent with the number of light neutrino families, N_ν =3. The literature on the subject suggests a range for N_ν between 1.8 and 3.9 (see, for instance, Ref. <cit.>). Notably, measurements from experiments at CERN indicate N_ν = 2.9840 ± 0.0082 based on the Z_0 width <cit.>.


In the framework of the big bang model, nearly all neutrons ultimately combine to form ^4He, introducing a dependence of the abundance of ^4He on N_ν and the neutron lifetime τ_n. The dependence on τ_n influences Big Bang Nucleosynthesis (BBN) in two distinct ways. Firstly τ_n impacts weak reaction rates due to its connection to the weak coupling constant. A longer (or shorter)  τ_n implies that the reaction rates remain lower (or higher)  than the Hubble expansion rate until a higher (or lower)  freeze-out temperature, significantly affecting the neutron-to-proton equilibrium ratio at freeze-out.

The ratio of neutrons to protons in thermal equilibrium, denoted as n/p, is approximately characterized by n/p = exp[-Δ m/k T] ≈ 1/6, where  T stands for the temperature at the weak freeze-out, and Δ m∼ 1.3 MeV signifies the mass difference of neutron and proton. Another significant effect of τ_n arises from neutron decay occurring between the weak freeze-out (t∼ 1 s) and the onset of nucleosynthesis (t∼ 200 s), causing a decrease in the n/p ratio to approximately 1/7. A longer τ_n corresponds to a higher predicted BBN helium abundance. Here, we adhere to the most recent experimental determination of τ_n, specifically τ_n=877.75 ± 0.28_stat + 0.22/-0.16_syst s, as documented in recent experiments <cit.> (for an extensive review of neutron lifetime, see Ref. <cit.>). The ramifications of variations in τ_n on BBN predictions have been recently explored in Ref. <cit.>.

The baryonic density in the universe, derived from astronomical observations anisotropies in the CMB radiation, places constraints on the value of the baryons to photon ratio, denoted as η. This ratio remains constant as the universe expands. Calculations within the BBN model align with the value determined from WMAP observations, indicating η = 6.104 ± 0.055 × 10^-10 <cit.>.


Our focus in this brief review centers on the abundances of light elements. During the initial stages of the universe, specifically within the first 20 minutes, when conditions were sufficiently dense and hot to facilitate nuclear reactions, the temperature of the primordial plasma underwent a substantial decrease from a few million electronvolts (MeV) to approximately 10 kiloelectronvolts (keV). This environment gave rise to the production of light nuclides such as ^2H, ^3He, ^4He, and, to a lesser extent, ^7Li, through a complex network of nuclear reactions. The final abundances for these nuclides can be determined through various observational techniques, spanning diverse astrophysical environments.

However, discrepancies in lithium (Li) abundances observed in metal-poor stars, compared to inferences from WMAP data, have sparked extensive inquiries into both BBN and the stellar mixing processes that potentially influence Li abundance. The findings of these investigations remain inconclusive, emphasizing the need for further studies on the abundances of light elements in low-metallicity stars and in extragalactic HII regions, i.e., denser, collapsed regions in giant molecular clouds where stars are forming. Additionally, refining estimates using BBN becomes crucial to address this issue comprehensively. High-resolution spectroscopic analyses of stellar and interstellar matter need to be integrated with  nuclear physics experiments and nucleosynthesis models and nuclear theories to gain a deeper understanding of these intricacies <cit.>.

The calculation of nuclear reaction rates during BBN relies on fundamental inputs, one of which is the Maxwell-Boltzmann distribution describing the kinetic energy of particles in a plasma. This distribution is grounded on basic assumptions of the Boltzmann-Gibbs statistics, namely,  (a) the particle collision time is significantly shorter than the average time between collisions, (b) interactions are localized, (c) velocities of two particles are uncorrelated, and (d) there is a minimal energy transfer  to and from collective variables and therefore energy of the colliding particles is locally conserved. When conditions (a) and (b) are not satisfied, the effective two-body interaction becomes non-local, dependent on the energy and momentum of the ions. Assuming that the one-particle energy distribution follows the Maxwellian pattern, additional assumptions on particle correlations are needed to prove that the relative-velocity distribution is still Maxwellian.

While the Boltzmann-Gibbs (BG) model of statistical mechanics holds well in numerous situations, recent theoretical efforts have increasingly explored alternative approaches, treating BG statistics as a limit within more general theories <cit.> (see also, <cit.>). These alternative theories aim to describe systems with long-range interactions and memory effects, departing from the assumptions of ergodicity. One prominent alternative to BG statistics, as suggested by C. Tsallis <cit.>. This approach, known as non-extensive statistics, has gained popularity in recent years (for in-depth details, refer to comprehensive reviews <cit.>). While traditional statistical mechanics assumes energy to be an “extensive" variable, proportional to the system size, and entropy to be extensive, this justification may be valid for short-range interactions binding matter. However, when dealing with long-range interactions, particularly in the context of gravity, entropy is found to be non-extensive <cit.>.


In traditional statistics, determining the mean values of system quantities like energy, the number of molecules, and volume involves seeking the probability distribution that maximizes entropy while satisfying constraints for the correct average values of those quantities. Cognizant of this, Tsallis introduced a novel, non-extensive entropy, commonly known as the Tsallis entropy, and advocated for its maximization under the same constraints, replacing the conventional (BG) entropy. A family of Tsallis entropies are characterized by a real-valued parameter q that quantifies the departure from extensivity, reverting to the usual definition of entropy when q = 1. One has demonstrated in various scenarios that traditional results from statistical mechanics can be translated into this new framework <cit.>. These entropies,  result in probabilities following power laws instead of the typical exponential laws found in the standard statistics <cit.>. In many instances where the Tsallis formalism is employed, such as in Ref. <cit.>, the non-extensive parameter is held constant q ≃ 1, i.e., close to the value yielding conventional statistical mechanics. Some studies have explored substantial deviations of q from unity to elucidate various phenomena across diverse scientific domains <cit.>.

We will demonstrate that the foundational Maxwell-Boltzmann (MB) distribution, crucial to big bang and stellar evolution nucleosynthesis, undergoes significant modification under non-extensive statistics when q deviates substantially from unity. Consequently, this alteration exerts a profound impact on the predictions of BBN. A robust prior justification for a substantial deviation of q from unity in the BBN is lacking. Specifically, since matter is presumed to be in equilibrium with radiation during BBN, any modification to the MB distribution of velocities would influence the Planck distribution of photons. Analyses of temperature fluctuations in the CMB have demonstrated that a Planck distribution modified by the Tsallis statistics accurately describes temperature fluctuations in the CMB as observed with WMAP and q = 1.045 ± 0.005, a value close to unity but not precisely equal <cit.>. Notably, temperature distributions based on Boltzmann-Gibbs statistics,  inadequately describe CMB temperature fluctuations <cit.>. These fluctuations, possibly allow for even larger variations in q and might manifest during the BBN, consequently influencing the exponentially small tail of the velocity distribution.


Building on the accomplishments of the BBN model, it is reasonable to trust that it can impose stringent constraints on the parameter q of the non-extensive statistics. Within the scientific literature, endeavors to address the lithium problem often involve exploring various facets of “new physics" <cit.>. This study contributes to the array of novel approaches, with implications extending far beyond the expected scope for resolving the lithium problem within BBN <cit.>. If Tsallis statistics effectively characterizes deviations in the tails of statistical distributions, then BBN becomes a valuable tool for probing such deviations.  Indeed, the potential for a departure from the Maxwellian distribution and the repercussions of modifying the MB distribution's tail for nuclear fusion in stars have been explored in previous studies <cit.> and in recent studies it has been shown that a small departure of the parameter q from the unity is enough to solve the lithium puzzle <cit.>. 




§ EXTENSIVE VERSUS NON-EXTENSIVE STATISTICS




 §.§ Extensive Statistics

Consider two levels separated by an energy E with n_0 particles occupying the lowest energy level and n_1 the higher energy one. For identical particles, the number of ways that one can arrange the system is Ω = N!/n_0! n_1!, where N=n_0+n_1. The number Ω attains very large values for macroscopic systems and a better measure of the number of possibilities that the system can rearrange is  S=k lnΩ. This is the usual definition of entropy and k is identified as the Boltzmann constant. For the aforementioned two-level problem we get S=k[ln N! - ln n_0! - ln n_1!]. If an additional amount E of energy is added to the system, one particle is promoted to the higher level and the entropy increases to S' with the replacements n_0 → n_0 - 1 and  n_1 → n_1 + 1. The additional entropy gain is Δ S = S'-S=kln (n_0/n_1), for n_1 ≫ 1. If we compare to the thermodynamics relation for the change in entropy, we have Δ S = Δ E/T=E/T, which implies that n_1=n_0 exp(-E/kT), the famous Boltzmann law of statistics. 


Extending this logic to multiple energy levels, one defines occupation probabilities proportional to particle numbers in each level i. Normalizing these probabilities, yields

    P(E_i) = 1/ Zexp(- E_i/kT),      where    Z = ∑_i P(E_i),

for the occupation probability of the state i with energy E_i. Z is know as the partition function. 

The Boltzmann law is a cornerstone of thermodynamics and statistics. For non-interacting particles in thermal equilibrium, one can define f(v)dv as the probability of encountering a particle with velocity v. This leads to the Maxwell-Boltzmann (MB) distribution:

    f(v) = 4π v^2 ( m/2π kT)^3/2exp( -mv^2/2kT) ,

which we focus on in this review.

The property of extensivity can be easily understood if the Boltzmann-Gibbs entropy is expressed as

    S_BG=-k ∑_i p_i ln p_i ,

where p_i represents the probability of finding the system in the i-th microstate, and ∑_i p_i = 1. Considering two independent systems labeled A and B, the combined probability of finding a state i+j in a  microstate i of A and j  of B, is:

    p^A+B_i+j = p^A_i · p^B_j.

Inserting this into Eq. <ref>, implies that the Boltzmann-Gibbs entropy adheres to the relationship

    S_A+B=S_A +S_B.

Hence, entropy derived from Boltzmann-Gibbs statistics is recognized as an extensive quantity since, akin to volumes of thermodynamic systems, it adds extensively. 



 §.§ Non-Extensive Statistics

Boltzmann-Gibbs statistics provides a sound framework for treating ergodic systems, where given sufficient time, any state has a non-zero probability p_i of occurrence. However, particle correlations and the breakdown of ergodicity often push systems into out-of-equilibrium states, rendering Boltzmann-Gibbs statistics inadequate. Furthermore, some assumptions inherent to Boltzmann-Gibbs statistics, particularly when applied to interacting systems, prove overly restrictive. They are not universally valid, even in systems of particles at thermodynamic equilibrium. Indeed, alternatives to Boltzmann-Gibbs statistics are available and frequently utilized in the literature <cit.>.

Within the realm of non-extensive statistics <cit.>, a departure from tradition emerges as the conventional entropy yields to an alternative formulation, due to Tsallis <cit.>:

    S_q=k 1-∑_i p^q_i/q-1,

where q denotes a real number. When q=1, S_q=  S_BG, illustrating Tsallis statistics as a smooth extension of the Boltzmann-Gibbs entropy.

Following that, it can be easily shown that

    S_q(A+B)= S_q(A) + S_q(B)+(1-q)/k S_q(A) S_q(B).

Here, the parameter q acts as a gauge for non-extensivity, a concept pioneered by Tsallis in the development of statistical mechanics grounded in this generalized entropy <cit.>.

A significant implication of this non-extensive framework is the deviation of the distribution function maximizing S_q from the Maxwellian <cit.>. With q=1, the MB distribution is faithfully replicated. However, for q<1, high-energy velocities become notably more probable compared to the extensive distribution. Conversely, when q>1, high-energy velocities become less probable, with a subsequent cutoff beyond which no velocities are possible.



 §.§ Nuclear Reaction Rates with Extensive Statistics


Within stellar environments, the rate of thermonuclear reactions, weighted with a MD in the plasma, is given by <cit.>:

    R_ij   =   N_iN_j/1+δ_ij⟨σ v ⟩
       =   N_iN_j/1+δ_ij(8/πμ)^1/2(1/kT)^3/2∫_0^∞dE S(E) exp[-(E/kT+2πη(E))].

Here,  T denotes the temperature, S(E) stands for the astrophysical S-factor, σ denotes the nuclear fusion cross section, v signifies the relative velocity between the ij-pair, E=μ v^2/2 denotes their relative motion energy, μ stands for their reduced mass, N_i represents the number of nuclear species i,   and η=Z_iZ_je^2/ħ v is the Sommerfeld factor, where Z_i denotes the charge of the i-th particle. The factor 1+δ_ij in the denominator avoids double-counting when i=j. 
A modified expression for the reaction rates, namely, N_A ⟨σ v ⟩,  can be expressed in units of cm^3 mol^-1 s^-1, where N_A denotes Avogadro's number, and ⟨σ v ⟩ is given in  Eq. (<ref>) with the substitution of the Maxwell distribution f(E) by Eq. (<ref>) defined below. 

To factor out the steep energy dependence of reaction cross sections on energy, one usually defines:

    σ(E) = S(E)/Eexp[-2πη(E)].

We write 2πη=b/√(E), with

    b=0.9898 Z_iZ_j√(A)  MeV^1/2,

where A represents the reduced mass in atomic mass units. The S-factor exhibits a weaker dependence on E, unless it approaches  a resonance, where it manifests a strong peak. For a weak S(E) energy dependence, the reaction rate integral in Eq. <ref> is peaked at the “Gamow energy"  

    E_0=0.122(Z_i^2Z^2_jA)^1/3T_9^2/3  MeV,

and the “Gamow width" of the peak is 

    Δ E=0.2368(Z_i^2Z^2_jA)^1/6T_9^5/6  MeV,

where T_9 it the temperature in units of 10^9 K, and A denotes the reduced mass in atomic mass units (amu). Regarding charged particle reactions in the BBN, a high precision (within 0.1%) is attained by employing integration limits ranging from E_0 - 5Δ E to E_0+5Δ E, where Δ E is determined by Eq. (<ref>). 

Neutron-induced reactions are a subset of paramount importance for BBN, particularly including the p(n,γ)d, ^3He(n,p)t, and ^7Be(n,p)^7Li reactions. The cross section at low energies typically scales proportionally to 1/v, where the neutron velocity is represented by v=√(2mE)/ħ. Consequently, it's sometimes advantageous to reformulate Eq. (<ref>) as

    σ(E) = S(E)/E= R(E)/√(E),

where R(E) is a smooth function of energy akin to an S-factor. 





§ NUCLEAR REACTION RATES




 §.§ Nuclear Reaction Rates with Non-Extensive Statistics

The non-extensive characterization of the Maxwell-Boltzmann distribution is realized through the replacement f(E) → f_q(E), as detailed by <cit.>:

    f_q(E)=[1-q-1/kTE]^1/q-1   q→1⟶  exp(-E/kT),
       0<E<∞.

If q-1<0, this equation yields real values for any value of E≥ 0. On the other hand, if q-1>0, it is clear that f(E) is real only if the expression within square brackets remains positive. In other words,

    0 ≤ E ≤kT/q-1,   for  q≥ 1,      and    0≤ E,   for  q≤ 1.

Consequently, for 0<q<1, the range is 0<E<∞, while if 1<q<∞, it becomes 0<E<E_max=kT/(q-1).

Using this alternate statistics, the nuclear reaction rates transform into:

    R_ij= N_iN_j/1+δ_ij∫_0^E_maxdE S(E) M_q(E,T),

with the “modified" Gamow energy distribution given by

    M_q(E,T)   =    A(q,T)(1-q-1/kTE)^1/q-1e^-b/√(E)
       =    A(q,T) (1-q-1/0.08617T_9E)^1/q-1exp[-0.9898 Z_iZ_j√(A/E)],

representing the non-extensive counterpart of the nuclear astrophysics input in the BBN, where E_max=∞ if 0<q<1, and E_max=kT/(1-q) if 1<q<∞, and with E given in MeV units. The function A (q,T) is a normalization constant dependent on the temperature T and on the non-extensive parameter q.

In the case of neutron induced reactions, the modified Gamow distribution function in Eq. (<ref>) can be expressed as

    M_q(E,T)= A (q,T)f_q(E)= A (q,T)(1-q-1/kTE)^1/q-1.

The omission of the tunneling factor in this equation decreases the reliance of the reaction rates on the non-extensive parameter q.


The quest for a unified statistical framework encompassing Fermi-Dirac, Bose-Einstein, and Tsallis statistics for particles in a plasma has been a focal point of research, as discussed in Ref. <cit.>. In this pursuit, a non-extensive statistical distribution for both fermions and bosons was proposed, represented by the expression:

    n_q^± (E)=  1/[1-(q-1)(E-μ)/kT]^1/q-1± 1,

where μ denotes the chemical potential. This formulation seamlessly transitions to the Fermi distribution (n^+) as q approaches 1 and to the Bose-Einstein distribution (n^-) for photons when μ=0 and q tends to 1.
Furthermore, by incorporating the appropriate modifications, Planck's law governing the distribution of radiation is derived by augmenting n^- in the equation above (with μ=0) with ħω^2/(4π^2 c^2), where E=ħω.
For electrons, the number density can be derived from n^+ in the equation above, with considerations for relativistic or non-relativistic regimes, necessitating the inclusion of proper phase factors. Additionally, normalization factors A^± (q,T) are introduced to ensure consistency, as previously established in the literature.

During the early universe, electron density experiences significant variation with temperature, with electron number densities reaching approximately 10^32/cm^3 at T_9=10, vastly surpassing the electron density at the sun's core, approximately 10^26/cm^3. This heightened electron density stems from the copious production of e^+e^- pairs by abundant photons during BBN. However, despite the substantial electron densities, nuclear reactions remain largely unaffected, with the enhancement of nuclear reaction rates due to electron screening proving negligible, as demonstrated by reference <cit.>.







The electron Fermi energy at these densities remains much smaller than kT for energies relevant to BBN, permitting the utilization of μ=0 in Eq. (<ref>) for n^+. In Figure <ref>, the relative difference ratio (n_q^±-n^±)/n^± between non-extensive (n_q) and extensive (n=n_q→ 1) statistics is illustrated. Solid curves represent Fermi-Dirac (FD) statistics (n^+), while dashed curves depict Bose-Einstein (BE) statistics (n^-), with results shown for q=2 and q=0.5, at T_9=10. Notably, non-extensive distributions are amplified for q>1 and suppressed for q<1 compared to respective FD and BE quantum distributions. These deviations grow larger with energy, with the non-extensive electron distribution roughly twice as large as the usual FD distribution at E_e = 1 MeV, at T_9=10.

Although numerical results with modified FD and BE distributions are not presented here, it is anticipated that these generalizations will profoundly impact the freezout temperature and the neutron to proton (n/p) ratio. The freezout temperature, marking when the weak reaction rate, Γ∼<σ v>, for the weak reaction ν_e + n → p+e^- becomes slower than the universe's expansion rate, remains unaffected by non-extensive statistics due to low particle densities compared to kT, allowing μ to be set to zero. The adoption of non-extensive quantum distributions, as in Eq. (<ref>), yields the same temperature powers as those predicted by the FD and BE distributions, ensuring consistency in reaction rates.




 §.§ Non-Maxwellian Distribution for Relative Velocities



It's important to highlight that while a Maxwellian distribution governs individual particle energies, it doesn't inherently ensure a Maxwellian distribution of relative velocities. To achieve this, additional assumptions regarding particle correlations are necessary to establish the distribution of relative velocities, which is pivotal for accurate rate calculations. This matter has been extensively studied in references <cit.>, where it's demonstrated that non-Maxwellian distributions arise from non-extensive statistics.

In this context, we illustrate that if the particle velocity distribution deviates from Maxwellian behavior as described by Eq. (<ref>), then an adjustment can be made to the two-particle relative velocity to accommodate the recoil of the center-of-mass. Denoting the kinetic energy of a particle as E_i, this distribution is represented by:

    f_q^(i) = (1 - q - 1/kTE_i)^1/q - 1   q→1⟶exp[-(E_i/kT)].


Let's consider the equation:

    f_q^(i) = (1 - q - 1/kTE_i)^1/q - 1

Expressing the two-particle energy distribution as f_q^(12)=f_q^(1)f_q^(2), and exponentiating the Tsallis distribution, we obtain:

    f_q^(i) = exp{1/q - 1[ln(1 - q - 1/kT E_i)]}.

This leads to the simplification of the product f_q^(1)f_q^(2) as:

    f_q^(12) = exp{1/q - 1[ln(1 - q - 1/kT E_1)ln(1 - q - 1/kTE_2)]}.


Considering E_i = m_i v_i^2/2, where E_1 + E_2 = μ v^2/2 + M V^2/2, with μ representing the reduced mass of the two particles, M = m_1 + m_2, v being the relative velocity, and V as the center-of-mass velocity, the product inside the natural logarithm simplifies to:

    1 - 1 - q/kT(μ v^2/2 + M V^2/2) + (1 - q/kT)^2 μ v^2/2M V^2/2
          =(1 - 1 - q/kTμ v^2/2)(1 - 1 - q/kTM V^2/2)


Therefore, the two-body distribution can be factorized into a product of relative and center-of-mass components:

    f_q^(12)(v,V,T) = f_q^(rel)(v,T) f_q^(cm)(V,T),

where

    f_q^(rel)(v,T)   =    A_rel(q,T)(1 - 1 - q/kTμ v^2/2)^1/q - 1,
    
    f_q^(cm)(V,T)   =    A_cm(q,T)(1 - 1 - q/kTM V^2/2)^1/q - 1.


The normalization constants are determined by the condition

    ∫ d^3v d^3V f_q^(12)(v,V,T) = 1.

Given the factorization of the distribution, unit normalization can be obtained for the relative and center-of-mass distributions separately. Hence, the distribution required in the reaction rate formula is:

    f_q(v,T) = ∫ d^3V f_q^(12)(v,V,T)= f_q^(rel)(v,T),

which adopts the same form as the absolute distribution.

In the limit q → 1, the two-particle distribution reduces to a Gaussian, where the last term in the left-hand side of Eq. (<ref>) drops out:

    f_q^(12)=  A(q,T) exp[-(μ v^2/2 + MV^2/2/kT)],

as anticipated.




§ FIXING THE BBN LITHIUM PROBLEM




 §.§ Non-Extensive Energy Distributions

Referring to abundant literature on non-extensive statistics (see, for example, <cit.>), it is commonly expected that the non-extensive parameter q remains close to unity. Nonetheless, to explore the ramifications of non-Maxwellian distributions on BBN, we investigate scenarios where q significantly deviates from unity. This investigation aims to deepen our understanding of the underlying physics departing from Boltzmann-Gibbs (BG) statistics. Figure <ref> illustrates Gamow energy distributions of deuterons of importance for the reaction ^2H(d,p)^3H at T_9=1.  The solid line, representing q=1, refers to the use of a  MB distribution. Also depicted are the outcomes obtained if one employs non-extensive MB distributions. For the sake of understanding the deviation from the usual MB distribution, we exploit very large deviations for q from the unity, namely, q=0.5 (dotted line) and q=2 (dashed line).   It is noteworthy that  higher kinetic energies become more accessible when q<1 compared to the extensive case q=1, whereas high energies are less probable for q>1, with a distinct kinetic energy cutoff. Here, for q=2, the cutoff occurs at 0.086 MeV, equivalent to 86 keV.






We illustrate in Figure <ref> the reaction rates for ^2H(d,p)^3H as a function of T_9, considering q=0.5, 1, 2.  The solid curve represents the conventional Maxwell-Boltzmann distribution (q=1). On the other hand, the dashed and dotted curves correspond to q=2 and q=0.5 respectively. In both instances, deviations from the Maxwellian rate are evident. For q>1, the deviations are notably significant, leading to an overall suppression of the reaction rates, particularly at lower temperatures. This deviation stems from the non-Maxwellian energy cutoff, occurring at 0.086 MeV for this reaction, thereby hindering numerous reactions from occurring at higher energies. For q<1, the similarity in results to the Maxwell-Boltzmann distribution arises from a delicate balance between the suppression of reaction rates at low energies and their augmentation at high energies.






This clarifies why, at T_9=1, the relevant energy range for computing the reaction rate is delineated by the solid curve in Figure <ref>. For q<1, the Gamow window Δ E expands, contributing both to the suppression of reaction rates at low energies, akin to the Maxwell-Boltzmann distribution, and to a corresponding enhancement at higher energies. Consequently, this elucidates the nearly identical outcomes depicted in Figure <ref> for q=1 and q<1. This observation holds true for all charged particle reaction rates, except in cases where the S-factor exhibits a pronounced energy dependency at, or around, E=E_0. However, such behavior is absent in the most crucial charged-induced reactions during BBN. 





In Figure <ref>, we depict the distributions of the kinetic energy of nucleons relevant to the p(n,γ)d reaction at T_9=0.1 (upper panel) and T_9=10 (lower panel). q=1 is represented by the solid line,  with the dotted line for q=0.5 (dotted line) dashed line for q=2. Akin to the case of charged particles,  higher kinetic energies are more probable with q<1, whereas high energies are less accessible with q>1 and a cutoff is reached. A notable difference is the absence of the Coulomb barrier, leading to a lack of suppression of the reaction rates at low energies. When T_9 increases, the MB and non-MB distributions approach each other. 









 §.§ Backward Reactions

One usually defines the 1 + 2 → 3 + 4 reaction with positive Q-value as the forward reaction and  3 + 4 → 1 + 2  with a negative Q value as the reverse one. The ratio between the reverse and forward rates is  proportional to exp(-Q/kT). The non-extensive velocity distribution for reverse reactions is modified to

    M_q(E,T)   =   c× A(q,T)[1-(q-1)E+Q/kT]^1/q-1e^-b/√(E) ,

where, from the detailed balance theorem,

    c = m_1m_2 (2J_1+1)(2J_2+1)(1+δ_34)/m_3m_4 (2J_3+1)(2J_4+1)(1+δ_12).

It was first observed in Ref. <cit.> that the inclusion of Q-values in the reverse reactions has a large impact on the element abundance calculated in the BBN. 








 §.§ Solving the lithium puzzle, creating a deuterium problem


The primordial abundances are calculated using a standard BBN code and adopting the  cosmological baryon-to-photon ratio η = (6.104 ± 0.055) × 10^-10 <cit.>, and the neutron lifetime of τ = (878 ± 0.3) <cit.>.  The thermonuclear (forward and reverse) rates for those 17 principal reactions have been determined employing non-extensive statistics, including 11 reactions of primary importance and 6 secondary reactions in the primordial light-element nucleosynthesis. Standard MB rates have been adopted for the remaining reactions, as they play only a minor role during BBN.  Some reaction rates have been updated by using recent investigations <cit.>.

Table <ref>  shows that the predicted and observed abundances of D, ^4He, and ^7Li are in close agreement  when a non-extensive velocity distribution with 1.069 < q <1.082 is adopted.  The ^4He abundance is rather resilient to the variation of q and is in agreement with the data from Ref.  <cit.>.  The predicted ^3He abundance for the above range of q agrees at the 1.8σ level with the observed abundance <cit.>. The prediction of ^7Li is greatly improved and agrees with the observation <cit.>.  Thus,  the cosmological lithium problem can be solved with the non-extensive statistics. However, the deuterium abundance now deviates appreciably from the value one would obtain with the standard BG statistics (q=1).

The agreement between the ^7Li abundance and observational data, when considering non-extensive statistics, stems from the reduced production of ^7Li and radioactive ^7Be (which decays into ^7Li) when q>1. The production of these elements primarily occurs through radiative capture reactions, namely ^3H(α,γ)^7Li and ^3He(α,γ)^7Be. However, when q>1, the alpha-capture rates of these reactions decrease due to the decreased availability of high-energy nucleons compared to the standard MB (q=1) distribution. Conversely, the reverse photodisintegration rates remain unaffected by q due to the adoption of Planck’s radiation law for photon energy density. Consequently, the overall production of ^7Li and ^7Be decreases, resulting in a better match between predicted and observed primordial abundances. Figure <ref>, illustrates the temporal and thermal evolution of primordial abundances during BBN calculated using both MB and non-extensive distributions (with an average value of q=1.0755). It demonstrates a significant reduction in the predicted abundance of ^7Be (which eventually decays into ^7Li) compared to q=1, ultimately offering a solution to the ^7Li problem within this model.

Deuterium is a very fragile isotope that is surely destroyed after BBN throughout stellar evolution. The tension with the predicted abundances using q=1.0755 is explained because as seen from Figure <ref>, the reaction rates for ^2H(d,p)^3H tend to decrease for any value of q>1. Therefore, there is more deuterium left over from BBN. The fact that we find excellent agreement between observed and  predicted primordial abundances of ^7Li  in the interval 1.069<q<1.082 is intriguing. It shows that the solution of the so-called lithium puzzle might be due to a small change in the physics during the BBN. But, it seems that the solution proposed in Ref. <cit.>, while solving the lithium puzzle leads to a discrepancy with the observed value of 2.527 ± 0.030  (× 10^-5) for the  D/H abundance  <cit.>. With the standard MB statistics, the BBN deuterium abundance compares much better with the observations <cit.>.



§ CONCLUSIONS

 The concentrations of the four light elements, namely D, ^3He, ^4He, and ^7Li, serve as robust constraints for primordial abundance assessment. Across all these elements, non-extensive statistics with the non-extensive parameter q slightly above the unity manifest greater agreement with experimental findings. A chi-square fitting against observed elemental abundances suggests that the non-extensive parameter  q = 1.0755 <cit.> best describe observations. 

Presently, ^3He measurements are only viable within our Galaxy's interstellar medium, limiting assessments at low metallicity crucial for fair primordial element comparisons. Consequently, determining the primordial ^3He abundance remains unreliable.  Our analysis has excluded potential alterations in n/p conversion rates induced by non-extensive statistics, alongside expected electron distribution changes. These adjustments would influence freezeout temperatures, consequently impacting ^4He abundance.

Using a non-extensive MB statistics to replicate observed light element abundances without compromising other successful BBN predictions appears unfeasible if one considers its discrepancy with the deuterium observations. Should confirmation of a non-Maxwellian distribution stemming from Tsallis non-extensive statistics emerge, significant revisions to our comprehension of cosmic universe evolution would be imperative.


Acknowledgements

C.A.B. acknowledges support by the U.S. DOE Grant DE-FG02-08ER41533 and the Helmholtz Research Academy Hesse. Shubhchintak acknowledges SERB, DST, India for a Ramanujan Fellowship (RJF/2021/000176).

4
ABG48 R. Alpher, H. Bethe, and G. Gamow, Phys. Rev. 73, 803 (1948).
Alp50 R. Alpher and R.C. Herman, Rev. Mod. Phys. 22, 153 (1950).
Alp53 R. Alpher, J. Follin and R. Herman, Phys. Rev. 92, 1347 (1953).
WFH67 Robert Wagoner, William A. Fowler, and F. Hoyle, Astrophys. J.  148, 3 (1967).
Zel64Ya.B. Zel’dovich, Sov. Phys. Usp. 6, 475 (1964).
Peeb09 P.J.E. Peebles, L.A. Page and R.B. Partridge, Finding the Big Bang (Cambridge Univ. Press, 2009).
Peeb19 P.J.E. Peebles, Cosmology’s Century (Princeton University Press, 2019).
Not11 P. Noterdaeme, P. Petitjean,  R. Srianand, C.  Ledoux, and S.  Lopez,  Astron. and Astrophys. 526, L7 (2011).
Oli02 Keith A. Olive, Lectures at the First NCTS Workshop on Astroparticle Physics, December 6 - 9, 2001, Kenting, Taiwan, arXiv:astro-ph/0202486 (2002).

PDG20 P. A. Zyla et al. (Particle Data Group), Prog. Theor. Exp. Phys. 2020, 083C01 (2020).

Gonz21F. M. Gonzalez et al., (UCN Collaboration), Phys. Rev. Lett. 127, 162501 (2021).
FG11 Fred E. Wietfeldt and Geoffrey L. Greene, Rev. Mod. Phys.  83, 1173 (2011).

yeh23 T. H. Yeh, K. A. Olive, B. D. Fields, Universe 9, 183 (2023).

Agh2020 N. Aghanim, et al., Astron. Astrophys. 641, A6 (2020).
Fie11 Brian D. Fields, Ann. Rev. Nucl. Part. Sci. 61, 47 (2011).
cyburt16 Richard H. Cyburt, Brian D. Fields, Keith A. Olive, and Tsung-Han Yeh
Rev. Mod. Phys. 88, 015004 (2016).
Ren60 A. Rényi, Proc. of the 4th Berkeley Symposium on Mathematics, Statistics and Probability,  pp. 547561 (1960).
Ts88 C. Tsallis,  J. Stat. Phys.  52, 479 (1988).
GT04 M. Gell-Mann and C. Tsallis, eds., “Nonextensive Entropy-Interdisciplinary Applications, (Oxford University Press, New York, 2004). For a regularly updated bibliography see URL: http://tsallis.cat.cbpf.br/biblio.htm
TGS05 C. Tsallis, M. Gell-Mann and Y. Sato, Proc. Natl. Acad. Sc. USA 102, 15377 (2005).
Ts09 “Introduction to Nonextensive Statistical Mechanics", C. Tsallis, Springer Verlag (2009), ISBN 978-0-387-85358-1. 


FL01 K. Sau Fa and E.K. Lenzi, J. of Math. Phys. 42, 1148 (2001).
Lim02 J.A.S. Lima, R. Silva, and J. Santos, Astron. and Astrophys. 396, 309 (2002).
TS03 A. Taruya and M. Sakagami, Phys. Rev. Lett. 90, 181101 (2003).
TS04 A. Taruya and M. Sakagami, Cont.Mech. and Therm. 16, 279 (2004).
CS05 P.H. Chavanis and C. Sire, Physica A 356, (2005) 419.
Pes01 M.E. Pessah, D.F. Torres and H. Vucetich,  Physica A 297, 164 (2001).
Ber07 A. Bernui, C. Tsallis and T. Villela, Europhys. Lett. 78, 19001 ( 2007).
BFH13 C.A. Bertulani, J. Fuqua and M.S. Hussein, Astrophys. J. 767, 67 (2013).
Hou17 S.Q. Hou, J.J. He, A. Parikh, D. Kahl, C.A. Bertulani, T. Kajino, G.J. Mathews and G. Zhao, Astrop. J. 834, 165 (2017).
Oliv16 O. Oliveira,  C. A. Bertulani,  M. S. Hussein, W. de Paula, T. Frederico, Braz. J.  Phys. 46, 721 (2016). 
MQ05 M. Lissia and P. Quarati, Europhysics News 36, 211 (2005).
HK08 J. Haubold and D. Kumar, Astroparticle Physics 29 ,70 (2008).
Deg98 S. Degl'Innocenti et al, Phys. Lett. B 441, 291 (1998).
Cor99 M. Coraddu et al, Braz. J. Phys. 29, 153 (1999).
Sil98 R. Silva Jr., A.R. Plastino,  and J. Lima, Phys. Lett.  A249, 401 (1998).
Lim00 J.A. S. Lima, R. Silva Jr.,  and J. Santos,  Phys. Rev. E 61, 3260 (2000).
Mu06 V. Muñoz, Nonlin. Processes Geophys. 13, 237 (2006).
Fow67 W.A. Fowler, G.R.  Caughlan, and B.A. Zimmerman, Ann. Rev. Astron. Astrophys.  5, 525 (1967).
Buy00 F. Buyukkiliç and D. Demirhan,  Eur. Phys. J. B 14, 705 (2000).
WBB11 Biao Wang, C.A. Bertulani and A.B. Balantekin, Phys. Rev. C 83, 018801 (2011).
Kan05 G. Kaniadakis et al, Phys. Rev. E 71, 046128 (2005).
Hart18 M. Hartos, C. A. Bertulani, Shubhchintak, A. M. Mukhamedzhanov, and S. Hou, Astrophys. J. 862, 62 (2018).
BSM18 C.A. Bertulani, Shubhchintak, A.M. Mukhamedzhanov, EPJ Web of Conferences 184, 01002 (2018).
Mukh16 A. M. Mukhamedzhanov, Shubhchintak, and C.A. Bertulani, Phys. Rev. C 93, 045805 (2016).
BMS16 C.A. Bertulani, A. Mukhamedzhanov,  and Shubhchintak,  AIP Conf. Proc. 1753, 040001 (2016).
Lami17 L. Lamia, C. Spitaleri, C. A. Bertulani, S.Q. Hou, M. La Cognata, R.G. Pizzone, S. Romano, M.L. Sergi, A. Tumino, Astrophys. J. 850, 175 (2017).
Spar17 R. Sparta, C. Bertulani, C. Spitaleri, M. La Cognata, L. Lamia, A. Mukhamedzhanov, A. Tumino, J. Phys.: Conf. Series 940, 012017 (2017).
Aver21 E. Aver, D.A. Berg, K.A. Olive, R.W. Pogge, J.J. Salzer, and E.D. Skillman,  J. Cosmology Astropart. Phys. 03, 027 (2021).
BRB02T.M. Bania, R.T. Rood, D.S. Balser,  Nature 415, 54 (2002).
Sb10 L. Sbordone, et al., Astron. Astrop. 522, 26 (2010).
Coo18 R.J. Cooke, M. Pettini and C.C. Steidel,  Astrophys. J. 855, 102 (2018).
Pitr21 Cyril Pitrou,   Alain Coc,   Jean-Philippe Uzan,  and  Elisabeth Vangioni, MNRA  502,  2474 (2021).




