format: text

text:
# 1 Introduction

## 1 Introduction

The relativistic field theory of hadrons known as quantum hadrodynamics (QHD) has become a very useful tool for describing bulk and single-particle properties of nuclear matter and finite nuclei in the mean field approximation . Compared with the non-relativistic approach to the nuclear many-body problem, the relativistic model explicitly includes the mesonic degrees of freedom and treats the nucleons as Dirac particles. At the mean field (Hartree) level, nucleons interact in a relativistic covariant way by exchanging virtual mesons: an isoscalar-vector $`\omega `$ meson, an isoscalar-scalar $`\sigma `$ meson and an isovector-vector $`\rho `$ meson. With these ingredients the mean field treatment of QHD automatically takes into account the spin–orbit force, the finite range and the density dependence of the nuclear force. Adjusting some coupling constants and meson masses from the properties of a small number of finite nuclei, the relativistic mean field (RMF) model produces excellent results for binding energies, root-mean-square radii, quadrupole and hexadecapole deformations and other properties of spherical and deformed nuclei .

The original linear $`\sigma \omega `$ model of Walecka  was complemented with cubic and quartic non-linearities of the $`\sigma `$ meson  (non-linear $`\sigma \omega `$ model) to improve the results for the incompressibility and for finite nuclei. Since these models were proposed to be renormalizable, the scalar self-interactions were limited to a quartic polynomial and scalar-vector or vector-vector interactions were not allowed. Recently, and inspired by effective field theory (EFT), Furnstahl, Serot and Tang  abandoned the idea of renormalizability and extended the RMF theory by including other non-linear scalar-vector and vector-vector self-interactions as well as tensor couplings .

The EFT Lagrangian has an infinite number of terms since it contains all the non-renormalizable couplings consistent with the underlying QCD symmetries. Therefore, it is mandatory to develop a suitable scheme of expansion and truncation. At normal nuclear densities the scalar ($`\mathrm{\Phi }`$) and vector ($`W`$) meson fields are small compared with the nucleon mass ($`M`$), and they vary slowly with position in finite nuclei. This indicates that the ratios $`\mathrm{\Phi }/M`$, $`W/M`$, $`|\mathbf{}\mathrm{\Phi }|/M^2`$ and $`|\mathbf{}W|/M^2`$ can be used as the expansion parameters. With the help of the concept of naturalness, it is then possible to compute the contributions of the different terms in the expansion and to truncate the effective Lagrangian at a given level of accuracy . None of the couplings should be arbitrarily dropped out to the given order without a symmetry argument.

References  have shown that it suffices to go to fourth order in the expansion. At this level one recovers the standard non-linear $`\sigma \omega `$ model plus a few additional couplings, with thirteen free parameters in all. These parameters have been fitted (parameter sets G1 and G2) to reproduce some observables of magic nuclei . The fits display naturalness (i.e., all coupling constants are of the order of unity when written in appropriate dimensionless form), and the results are not dominated by the last terms retained. This evidence confirms the utility of the EFT concepts and justifies the truncation of the effective Lagrangian at the first lower orders.

Recent applications of the models based on EFT include studies of pion-nucleus scattering  and of the nuclear spin-orbit force , as well as calculations of asymmetric nuclear matter at finite temperature with the G1 and G2 sets . In a previous work  we have analyzed the impact of each one of the new couplings introduced in the EFT models on the nuclear matter saturation properties and on the nuclear surface properties. In Ref.  we have looked for constraints on the new parameters by demanding consistency with DBHF calculations and the properties of finite nuclei. During the last years a large amount of work has been devoted to measuring masses of nuclei far from stability . This body of experimental data has been used as a benchmark to test the predictions of the currently existent (relativistic and non-relativistic) nuclear effective forces . This fact motivates us to investigate in the present work the behaviour of the parameter sets G1 and G2 derived from EFT in regions far from the stability line.

To study ground-state properties of spherical open-shell nuclei one has to take into account the pairing correlations. Relativistic mean field calculations near the $`\beta `$-stability line have usually included pairing in a constant gap BCS approximation , with the gaps fitted to empirical odd-even mass differences. This approach works properly when the main effect of the pairing correlations is a smearing of the Fermi surface. Since the BCS pairing energy diverges for large momenta, a cut-off has to be introduced in the pairing channel to simulate phenomenologically the finite range of the particle-particle force. The limitations of this simple BCS method appear when one deals with nuclei far from the $`\beta `$-stability line. Close to the drip lines the Fermi level falls near the particle continuum and it is known that the BCS model does not provide a correct description of the coupling between bound and continuum states . In the non-relativistic framework this difficulty was overcome by the unified description of the mean field and the pairing correlations provided by the Hartree–Fock–Bogoliubov (HFB) theory , with Skyrme  or Gogny forces .

The same unified treatment was developed by Kucharek and Ring  in the relativistic framework. However, a quantitative description of the pairing correlations in nuclei cannot be achieved with relativistic mean field parametrizations because the meson exchange forces are not properly adapted to large momentum transfer . Later, Ring and coworkers  have used the RMF interaction for the particle-hole channel plus the pairing part of the Gogny force  (with the D1S parameters ) for the particle-particle channel, in relativistic Hartree-plus-Bogoliubov (RHB) calculations. Other authors have employed a density-dependent zero-range pairing force  instead of the Gogny pairing force .

Recent calculations with non-relativistic Skyrme forces and a zero-range force in the particle-particle channel have shown that a BCS approach is able to provide a good qualitative estimate of the drip lines if some quasibound states due to their centrifugal barrier (plus the Coulomb barrier for protons) are included in the calculation . In this work we will use a similar BCS approach with quasibound states to approximately take into account the effects of the continuum contributions near the drip lines. We will employ a constant pairing strength which can be considered as a simplification of the zero-range pairing force and which gives similar results to those obtained with a delta force for spherical nuclei .

The paper is organized as follows. We summarize the mean field approximation to the EFT nuclear model in the second section. In the third section we describe our modified BCS approach with quasibound states, and perform some calculations to test its possibilities and limitations by comparing with Bogoliubov results available from the literature. The fourth section is devoted to the detailed study with the EFT parametrizations G1 and G2 of properties such as separation energies, particle densities and spin–orbit potentials of nuclei belonging to chains of isotopes (isotones) with magic proton (neutron) number. Our conclusions are laid in the last section.

## 2 Relativistic mean field approach from effective field theory

The effective field theory approach to QHD has been developed in the recent years. The theory and the equations for nuclear matter and finite nuclei can be found in the literature  and here we shall only outline the formalism. We start from Ref.  where the field equations were derived from an energy density functional containing Dirac baryons and classical scalar and vector mesons. This functional can be obtained from the effective Lagrangian in the Hartree approximation, but it can also be considered as an expansion in terms of the ratios of the meson fields and their gradients to the nucleon mass of a general energy density functional that contains the contributions of correlations within the spirit of density functional theory .

According to Refs.  the energy density for finite nuclei can be written as

$`(𝐫)`$ $`=`$ $`{\displaystyle \underset{\alpha }{}}\phi _\alpha ^{}\{i𝜶\mathbf{}+\beta (M\mathrm{\Phi })+W+{\displaystyle \frac{1}{2}}\tau _3R+{\displaystyle \frac{1+\tau _3}{2}}A`$ (1)
$`{\displaystyle \frac{i}{2M}}\beta 𝜶(f_v\mathbf{}W+{\displaystyle \frac{1}{2}}f_\rho \tau _3\mathbf{}R+\lambda \mathbf{}A)+{\displaystyle \frac{1}{2M^2}}(\beta _s+\beta _v\tau _3)\mathrm{\Delta }A\}\phi _\alpha `$
$`+\left({\displaystyle \frac{1}{2}}+{\displaystyle \frac{\kappa _3}{3!}}{\displaystyle \frac{\mathrm{\Phi }}{M}}+{\displaystyle \frac{\kappa _4}{4!}}{\displaystyle \frac{\mathrm{\Phi }^2}{M^2}}\right){\displaystyle \frac{m_s^2}{g_s^2}}\mathrm{\Phi }^2{\displaystyle \frac{\zeta _0}{4!}}{\displaystyle \frac{1}{g_v^2}}W^4`$
$`+{\displaystyle \frac{1}{2g_s^2}}\left(1+\alpha _1{\displaystyle \frac{\mathrm{\Phi }}{M}}\right)\left(\mathbf{}\mathrm{\Phi }\right)^2{\displaystyle \frac{1}{2g_v^2}}\left(1+\alpha _2{\displaystyle \frac{\mathrm{\Phi }}{M}}\right)\left(\mathbf{}W\right)^2`$
$`{\displaystyle \frac{1}{2}}\left(1+\eta _1{\displaystyle \frac{\mathrm{\Phi }}{M}}+{\displaystyle \frac{\eta _2}{2}}{\displaystyle \frac{\mathrm{\Phi }^2}{M^2}}\right){\displaystyle \frac{m_{v}^{}{}_{}{}^{2}}{g_{v}^{}{}_{}{}^{2}}}W^2{\displaystyle \frac{1}{2g_\rho ^2}}\left(\mathbf{}R\right)^2{\displaystyle \frac{1}{2}}\left(1+\eta _\rho {\displaystyle \frac{\mathrm{\Phi }}{M}}\right){\displaystyle \frac{m_\rho ^2}{g_\rho ^2}}R^2`$
$`{\displaystyle \frac{1}{2e^2}}\left(\mathbf{}A\right)^2+{\displaystyle \frac{1}{3g_\gamma g_v}}A\mathrm{\Delta }W+{\displaystyle \frac{1}{g_\gamma g_\rho }}A\mathrm{\Delta }R,`$
where the index $`\alpha `$ runs over all occupied states $`\phi _\alpha (𝐫)`$ of the positive energy spectrum, $`\mathrm{\Phi }g_s\varphi _0(𝐫)`$, $`Wg_vV_0(𝐫)`$, $`Rg_\rho b_0(𝐫)`$ and $`AeA_0(𝐫)`$. Variation of the energy density (1) with respect to $`\phi _\alpha ^{}`$ and the meson fields gives the Dirac equation fulfilled by the nucleons and the meson field equations, which are solved self-consistently by numerical iteration. We refer the reader to Ref.  for the expressions of the variational equations.

The terms with $`g_\gamma `$, $`\lambda `$, $`\beta _s`$ and $`\beta _v`$ take care of effects related with the electromagnetic structure of the pion and the nucleon (see Ref. ). Specifically, the constant $`g_\gamma `$ concerns the coupling of the photon to the pions and the nucleons through the exchange of neutral vector mesons. The experimental value is $`g_\gamma ^2/4\pi =2.0`$. The constant $`\lambda `$ is needed to reproduce the magnetic moments of the nucleons. It is defined by

$`\lambda ={\displaystyle \frac{1}{2}}\lambda _p(1+\tau _3)+{\displaystyle \frac{1}{2}}\lambda _n(1\tau _3),`$ (2)
with $`\lambda _p=1.793`$ and $`\lambda _n=1.913`$ the anomalous magnetic moments of the proton and the neutron, respectively. The terms with $`\beta _s`$ and $`\beta _v`$ contribute to the charge radii of the nucleon .

In this work we will employ the EFT parameter sets G1 and G2 of Refs. . The masses of the nucleon and the $`\omega `$ and $`\rho `$ mesons take their experimental values: $`M=939`$ MeV, $`m_v=782`$ MeV and $`m_\rho =770`$ MeV. The thirteen parameters $`m_s`$, $`g_s`$, $`g_v`$, $`g_\rho `$, $`\eta _1`$, $`\eta _2`$, $`\eta _\rho `$, $`\kappa _3`$, $`\kappa _4`$, $`\zeta _0`$, $`f_v`$, $`\alpha _1`$ and $`\alpha _2`$ were fitted by a least-squares optimization procedure to twenty-nine observables (binding energies, charge form factors and spin–orbit splittings near the Fermi surface) of the nuclei 16O, 40Ca, 48Ca, 88Sr and 208Pb, as described in Ref. . The constants $`\beta _s`$, $`\beta _v`$ and $`f_\rho `$ were then chosen to reproduce the experimental charge radii of the nucleon. The fits yielded two best, distinct parameter sets (G1 and G2) with essentially the same $`\chi ^2`$ value .

We report in Table 1 the values of the parameters and the saturation properties of G1 and G2. One observes that the fitted parameters differ significantly between both interactions. For example, G2 presents a positive value of $`\kappa _4`$, as opposite to G1 and to many of the most successful RMF parametrizations, such as the NL3 parameter set . Formally a negative value of $`\kappa _4`$ is not acceptable because the energy spectrum then has no lower bound . Furthermore, the wrong sign in the $`\mathrm{\Phi }^4`$ coupling constant may cause troubles in obtaining stable solutions in light nuclei like 12C. We note that the value of the effective mass at saturation $`M_{\mathrm{}}^{}/M`$ in the EFT sets ($`0.65`$) is somewhat larger than the usual value in the RMF parameter sets ($`0.60`$). This fact is related with the presence of the tensor coupling $`f_v`$ of the $`\omega `$ meson to the nucleon, which has an important bearing on the spin–orbit force .

One should mention that the EFT perspective also has been helpful to elucidate the empirical success of the usual non-linear $`\sigma \omega `$ models that incorporate less couplings (just up to cubic and quartic self-interactions of the scalar field): the EFT approach accounts for the success of these RMF models and provides an expansion scheme at the mean field level and for going beyond it . In practice it has been seen that the mean field phenomenology of bulk and single-particle nuclear observables does not constrain all of the new parameters of the EFT model unambiguosly. That is, the constants of the EFT model are underdetermined by the observables currently included in the fits and different parameter sets with low $`\chi ^2`$ (comparable to G1 and G2) can be found . However, the extra couplings could prove to be very useful for the description of further observables. Indeed, for densities above the normal saturation density, and owing to the additional non-linear couplings, the EFT models are able  to give an equation of state and nuclear matter scalar and vector self-energies in much better agreement with the microscopic Dirac–Brueckner–Hartree–Fock predictions than the standard non-linear $`\sigma \omega `$ parametrizations (the latter completely fail in following the DBHF trends as the nuclear density grows ).

The sets G1 and G2 were fitted including centre-of-mass corrections in both the binding energy and the charge radius. Therefore, we will utilize the same prescription of Ref.  in our calculations with G1 and G2. Namely, a correction

$`E_{CM}={\displaystyle \frac{17.2}{A^{1/5}}}\mathrm{MeV}`$ (3)
to the binding energy and a correction

$`{\displaystyle \frac{3}{4}}{\displaystyle \frac{1}{(2MAE_{CM})}}\mathrm{fm}^2`$ (4)
to the mean-square charge radius.

## 3 The pairing calculation

It is well known that pairing correlations have to be included in any realistic calculation of medium and heavy nuclei. In principle the microscopic HFB theory should be used for this purpose. However, for pairing calculations of a broad range of nuclei not too far from the $`\beta `$-stability line, a simpler procedure is usually considered in which a seniority potential acts between time-reversed orbitals. In this section we want to discuss and test a straightforward improvement of this simple approximation to be able to describe in addition nuclei near the drip lines, at least on a qualitative level. Without the complications intrinsic to a full Bogoliubov calculation, our faster approximation will allow us later on to perform extensive calculations of chains of isotopes and isotones with the relativistic parameter sets.

The pairing correlation will be considered in the BCS approach . One assumes that the pairing interaction $`v_{\mathrm{pair}}`$ has non-zero matrix elements only between pairs of nucleons invariant under time reversal:

$`\alpha _2\stackrel{~}{\alpha }_2|v_{\mathrm{pair}}|\alpha _1\stackrel{~}{\alpha }_1=G,`$ (5)
where $`|\alpha =|nljm`$ and $`|\stackrel{~}{\alpha }=|nljm`$ (with $`G>0`$ and $`m>0`$). Most often the BCS calculations in the RMF model have been performed using a constant gap approach . Instead, here we choose a seniority-type interaction with a constant value of $`G`$ for pairs belonging to the active pairing shells.

The contribution of the pairing interaction to the total energy, for each kind of nucleon (neutrons or protons), is

$`E_{\mathrm{pair}}=G\left\{{\displaystyle \underset{\alpha >0}{}}\left[n_\alpha (1n_\alpha )\right]^{1/2}\right\}^2G{\displaystyle \underset{\alpha >0}{}}n_\alpha ^2,`$ (6)
where $`n_\alpha `$ is the occupation probability of a state with quantum numbers $`\alpha \{nljm\}`$ and the sum is restricted to positive values of $`m`$. One has

$`n_\alpha ={\displaystyle \frac{1}{2}}\left[1{\displaystyle \frac{\epsilon _\alpha \mu }{\sqrt{(\epsilon _\alpha \mu )^2+\mathrm{\Delta }^2}}}\right].`$ (7)
The Lagrange multiplier $`\mu `$ is called the chemical potential and the gap $`\mathrm{\Delta }`$ is defined by

$`\mathrm{\Delta }=G{\displaystyle \underset{\alpha }{}}\left[n_\alpha (1n_\alpha )\right]^{1/2}.`$ (8)
As usual the last term in Eq. (6) will be neglected. It is not a very important contribution and its only effect is a renormalization of the pairing energies .

Assuming constant pairing matrix elements (5) in the vicinity of the Fermi level one gets

$`{\displaystyle \frac{G}{2}}{\displaystyle \underset{\alpha >0}{}}{\displaystyle \frac{1}{\sqrt{(\epsilon _\alpha \mu )^2+\mathrm{\Delta }^2}}}`$ $`=`$ $`1,`$ (9)
$`{\displaystyle \frac{1}{2}}\left[1{\displaystyle \frac{\epsilon _\alpha \mu }{\sqrt{(\epsilon _\alpha \mu )^2+\mathrm{\Delta }^2}}}\right]`$ $`=`$ $`A,`$ (10)
where $`A`$ is the number of neutrons or protons involved in the pairing correlation. The solution of these two coupled equations allows one to find $`\mu `$ and $`\mathrm{\Delta }`$. Using Eqs. (7) and (8) the pairing energy for each kind of nucleon can be written as

$`E_{\mathrm{pair}}={\displaystyle \frac{\mathrm{\Delta }^2}{G}}.`$ (11)
This simple approach breaks down for nuclei far from the stability line. The reason is that in this case the number of neutrons (for isotopes) or protons (for isotones) increases, the corresponding Fermi level approaches zero and the number of available levels above it is clearly reduced. Moreover, in this situation the particle-hole and pair excitations reach the continuum. Ref.  showed that if one performs a BCS calculation using the same quasiparticle states as in a HFB calculation, then the BCS binding energies are close to the HFB ones but the r.m.s. radii (i.e., the single-particle wave functions) dramatically depend on the size of the box where the calculation is performed. This is due to the fact that there are neutrons (protons) that occupy continuum states for which the wave functions are not localized in a region, thus giving rise to an unphysical neutron (proton) gas surrounding the nucleus.

Recent non-relativistic calculations near the drip lines with Skyrme forces  have shown that the above problem of the BCS approach can be corrected, in an approximate manner, by taking into account continuum effects by means of the so-called quasibound states, namely, states bound because of their own centrifugal barrier (centrifugal-plus-Coulomb barrier for protons). When the quasibound states are included in the BCS calculation (from now on a qb-BCS calculation), it is necessary to prevent the unrealistic pairing of highly excited states and to confine the region of influence of the pairing potential to the vicinity of the Fermi level. Instead of using a cutoff factor as in Ref. , in our calculations we will restrict the available space to one harmonic oscillator shell above and below the Fermi level.

In order to check this approach we have performed with the G1 parameter set ($`G_n=21/A`$ MeV, see next section) calculations of the binding energy and r.m.s. radius of the 120Sn and 160Sn nuclei in boxes of sizes between 15 and 25 fm (as in the non-relativistic calculations of Ref. ). The results taking into account the quasibound levels $`1h_{9/2}`$, $`2f_{5/2}`$ and $`1i_{13/2}`$ for 120Sn, and $`1i_{11/2}`$ and $`1j_{13/2}`$ for 160Sn, are compared in Figure 1 with the output of a standard BCS calculation with only bound levels. It turns out that in the qb-BCS case the results are essentially independent of the size of the box where the calculations are carried out. When the quasibound levels are included the binding energies are larger than when only the bound levels are taken into account, due to the damping of the pairing correlation caused by disregarding the continuum states in the standard BCS calculation . We also show in Figure 1 the results of a BCS calculation using all bound and unbound levels (i.e., without restricting ourselves to quasibound levels) in the considered range. It is obvious that in this case the results are box dependent, as the binding energy and neutron r.m.s. radius of 160Sn evidence.

Another test of the qb-BCS approach concerns the asymptotic behaviour of the particle densities . In Figure 2 we display the radial dependence of the neutron density of 150Sn (as in Ref. ) calculated with the G1 parameter set in boxes of radii between 15 and 25 fm. For large enough distances the density decreases smoothly when the size of the box increases (except very near of the edge, where the density suddenly drops to zero because of the $`\phi _\alpha =0`$ boundary condition). This means that no neutron gas surrounding the nucleus has appeared. In a Bogoliubov calculation the asymptotic behaviour of the particle density is governed by the square of the lower component of the single-quasiparticle wave function corresponding to the lowest quasiparticle energy . This asymptotic behaviour is displayed by the (almost straight) dotted line in Figure 2. It can be seen that the density obtained with our approach decreases more slowly than the RHB density, i.e., asymptotically the qb-BCS density is not able to follow the RHB behaviour. This coincides with the conclusion of Ref.  (see Figure 19 of that work) where non-relativistic HFB densities are compared for large distances with the densities obtained in the qb-BCS approach with a state-dependent pairing .

Although the qb-BCS densities do not display the right asymptotic behaviour, it was conjectured in Ref.  that such an approach could allow one to compute properties of nuclei much closer to the drip lines than in a standard BCS calculation. Very recently, RHB calculations up to the drip lines of the two-neutron separation energy $`S_{2n}`$ for nickel isotopes  and of the charge and neutron r.m.s. radii for tin isotopes  have been carried out using the NL-SH parameter set  plus a density-dependent zero-range pairing force. We have repeated these calculations with our qb-BCS method for both isotopic chains (with a pairing interaction strength $`G_n=22.5/A`$ in the case of NL-SH).

We display the values of the $`S_{2n}`$ separation energies for the Ni chain in Figure 3a. The RHB calculation predicts the drip line at the isotope 100Ni and shows shell effects at $`N=28`$ and 50 (and to a minor extent at $`N=70`$). These features are well predicted by our simpler qb-BCS calculation. The differences between these qb-BCS and RHB results also come in part from the different pairing forces used in the calculations. To investigate this point we show in Figure 3b the neutron pairing energy obtained in our approach \[Eq. (11)\] for the isotopes of the Ni chain. It vanishes at $`N=28`$, 50 and 70, in agreement with the shell structure shown in Figure 3a by the $`S_{2n}`$ separation energies. The largest pairing energies are found in the middle of two closed shells and they are enhanced by increasing $`N`$. Figure 3b can be compared with the RHB values displayed in Figure 2 of Ref. . The tendencies are the same, though the qb-BCS pairing energies are slightly larger than in the RHB calculation. In Figure 4 we draw our results for the radii of the Sn isotopes, and compare them with the RHB values. In the case of the charge radii the agreement is excellent. The neutron radii obtained in our method closely follow the behaviour of the RHB neutron radii and the kink at $`N=132`$ is qualitatively reproduced.

We have furthermore computed the binding energies of nuclei of the $`N=20`$ isotonic chain for which RHB results exist with the NL3 parameter set . We present the extracted two-proton separation energies $`S_{2p}`$ in Table 2. The agreement between the qb-BCS and RHB approaches again is very good. In both models the last stable nucleus is 46Fe, as in experiment. Notice that in the present case the first levels with positive energy correspond to those of the $`pf`$ shell. Due to the Coulomb barrier all these levels become quasibound in our approach, and it is expected that they will lie close to the canonical levels. This explains the goodness of the qb-BCS energies for this isotonic chain.

¿From the previous comparisons we see that the simple qb-BCS calculation is able to reasonably follow the main trends of the more fundamental RHB pairing calculation. One can also conclude that the consideration of quasibound states in the BCS approach is, actually, a key ingredient to eliminate the spurious nucleon gas arising near the drip lines.

## 4 Results for EFT parameter sets

We want to analyze the ability of the G1 and G2 parameter sets based on effective field theory  to describe nuclear properties far from the stability line, i.e., far from the region where the parameters were fitted. To our knowledge such calculations have not been explored so far. We will contrast the results with experiment and with those predicted by the NL3 set, that we take as one of the best representatives of the usual RMF model with only scalar self-interactions.

As indicated, we shall use a schematic pairing with a state-independent matrix element $`G_\tau =C_\tau /A`$, where $`C_\tau `$ is a constant and $`\tau =n,p`$ for neutrons or protons, respectively. We fix the constant $`C_n`$ for neutrons by looking for the best agreement of our calculation with the known experimental binding energies of Ni and Sn isotopes. Similarly, we determine $`C_p`$ for protons from the experimental binding energies of the isotones of $`N=28`$ and $`N=82`$. The values obtained from this fit are $`C_n=21`$ MeV and $`C_p=22.5`$ MeV for the G1 set, $`C_n=19`$ MeV and $`C_p=25`$ MeV for G2, and finally $`C_n=20.5`$ MeV and $`C_p=23`$ MeV for NL3. Figure 5 shows that the neutron and proton state-independent gaps ($`\mathrm{\Delta }_n`$ and $`\mathrm{\Delta }_p`$) predicted by our calculation with G1 are scattered around the empirical average curve $`12/\sqrt{A}`$ . A similar picture is found with the parameter sets G2 and NL3.

### 4.1 Two-particle separation energies

In Figure 6a we present the two-neutron separation energies $`S_{2n}`$ for the chain of Ni isotopes. Clear shell effects arise at $`N=28`$ and 50. The three relativistic interactions (G1, G2 and NL3) slightly overestimate the shell effect at $`N=28`$ as compared with the experimental value, which also happens in more sophisticated RHB calculations with NL3 . In our qb-BCS approach some disagreement with experiment is found for the $`N=38`$ and $`N=40`$ isotopes. Again, this also occurs in the RHB calculations of Refs.  with NL3. However, if we compare Figure 6a with the results that we have shown in Figure 3a for the NL-SH parameter set, we see that NL-SH achieves a better agreement with experiment for these $`N=38`$ and $`N=40`$ isotopes.

We stop our calculation towards the neutron drip line when the two-neutron separation energy vanishes or when the neutron chemical potential becomes positive. The fact that $`S_{2n}`$ is not always zero at the drip line is connected with the quenching of the shell structure with $`N`$, which is a force-dependent property . This effect is illustrated in Figure 25 of Ref.  for HFB calculations with different non-relativistic forces. We find similar situations with the considered relativistic sets in our qb-BCS calculations of separation energies. In the case of the Ni isotopes we reach the drip line at $`N=66`$ with the G1 and NL3 sets and at $`N=68`$ with the G2 set. This agrees nicely with the value $`N=66`$ obtained in HFB calculations with the Skyrme forces SIII  and SkP . For NL-SH our qb-BCS scheme predicts the drip line at $`N=72`$ (see Figure 3a), the same value found in the RHB calculations of Ref. .

In Figure 7a we display our qb-BCS results for the two-neutron separation energies of the Sn isotopic chain. In Ref.  it was claimed that pure BCS calculations in the constant gap approach (with NL3) are not suitable for the Sn isotopes. We observe in Figure 7a that below $`N=60`$, as one moves towards $`N=50`$, some discrepancies with the experimental values appear, which also arise in the RHB calculations . The three forces slightly overestimate the shell effect at $`N=82`$ (as the RHB results of Refs.  for NL3). We have computed Sn isotopes up to $`A=176`$, when $`S_{2n}`$ vanishes for NL3 (in good agreement with RHB results for NL-SH  and HFB results for the Skyrme force SkP ). For G1 and G2 we find that $`S_{2n}`$ does not yet vanish at $`N=126`$, and it is not possible to increase the neutron number due to the shell closure at $`N=126`$ (the neutron chemical potential becomes positive for the $`N=128`$ isotope). This means that the quenching of the shell effect at $`N=126`$ for NL3 (and NL-SH) is larger than for the G1 and G2 parameter sets.

Our calculated $`S_{2n}`$ energies for Pb isotopes are shown in Figure 8. The experimental shell effect at $`N=126`$ is reasonably well reproduced by G1, G2 and NL3. The drip line is found at $`N=184`$ with a non-vanishing two-neutron separation energy, as in the calculations performed with the extended Skyrme force SLy4 in Ref. , where a similar approach (quasibound states and a state-dependent gap) to ours was used. The relatively large shell effect found at $`N=184`$ means that there is no quenching for this magic number in our qb-BCS approximation for the studied parameter sets. Indeed, to verify this point a full RHB calculation should be performed.

To analyze the proton pairing we have studied the two-proton separation energies in chains of isotones of $`N=28`$ (Figure 9a) and $`N=82`$ (Figure 10a). In the case of $`N=28`$, shell gaps appear at $`Z=20`$ and $`Z=28`$. For $`Z=20`$ the predicted gap is larger than in experiment. For $`Z=28`$, G1 and G2 agree better than NL3 with experiment. $`S_{2p}`$ vanishes at $`Z=30`$ for NL3, whereas it vanishes at $`Z=32`$ for G1 and G2. The isotones of $`N=82`$ display a clear shell effect at $`Z=50`$, in agreement with the non-relativistic calculation of Ref. . It is slightly larger for G2 than for G1 and NL3. Experimental information for this shell effect is not available. NL3 would predict another shell effect at $`Z=58`$, which does not appear experimentally. The effect is less pronounced in G1 and it does not show up in G2. The three forces indicate that the proton drip line is reached after the 156W isotope, in agreement with experimental information .

Figures 11a and 11b show, respectively, the calculated $`S_{2p}`$ separation energies for the $`N=50`$ and $`N=126`$ isotone chains. Note that we did not use any information about these nuclei in our fit of the $`G_p`$ pairing strength. For $`N=50`$ the set G2 follows the experimental data very well, specially for the larger $`Z`$. The trend of G1 and NL3 is only a little worse. The proton drip line is located at 100Sn in the three parametrizations, in good accordance with experiment. The quenching of the shell effect at $`Z=50`$ is larger for G2 than for G1 and NL3. The available data for two-proton separation energies of $`N=126`$ isotones are reasonably well estimated by the relativistic sets. However, the trend of NL3 is worse than that of G1 and G2. It would then be very interesting to perform RHB calculations of this chain to confirm the behaviour of NL3. The last nucleus of the chain stable against two-proton emission is 218U according to G1 and NL3, and 220Pu according to G2. The three sets predict a shell effect at $`Z=92`$, though it is relatively quenched for G2.

### 4.2 One-particle separation energies

We have computed one-neutron (one-proton) separation energies for Ni and Sn isotopes (for $`N=28`$ and $`N=82`$ isotones). The results are displayed in Figures 6b and 7b (9b and 10b), respectively. To deal with odd mass number nuclei we have used a spherical blocking approximation. One pair of conjugate states $`|\alpha `$ and $`|\stackrel{~}{\alpha }`$ is blocked, i.e., taken out of the pairing scheme . In the spherical approximation one replaces the blocked single-particle state by an average over the degenerate states in its $`j`$ shell. This way the rotational and time-reversal invariance of the many-body system is restored in the intrinsic frame . In this approach the contribution of the $`j`$ shell that contains the blocked state to the number of active particles and the pairing energy is

$`A_j`$ $`=`$ $`(2j1)n_j+1,`$ (12)
$`E_{\mathrm{pair},j}`$ $`=`$ $`G(2j1)[n_j(1n_j)]^{1/2},`$ (13)
respectively. The remaining active shells contribute in the usual manner \[Eqs. (6) and (10)\]. Due to rearrangement effects, blocking the single-particle state with smallest quasiparticle energy $`E_\alpha =\sqrt{(\epsilon _\alpha \mu )^2+\mathrm{\Delta }^2}`$ in the even $`A1`$ nucleus, does not necessarily lead to the largest binding energy of the odd $`A`$ nucleus. Therefore, in some cases one has to repeat the calculation blocking in turn the different single-particle states that lie around the Fermi level to find the configuration of largest binding energy .

The one-neutron (one-proton) separation energies lie over two different curves for even and odd neutron (proton) number. For Ni isotopes (Figure 6b) the three parameter sets G1, G2 and NL3 reproduce reasonably the experimental values. The shell effect at $`N=28`$ is, again, overestimated by the three forces. The heaviest Ni isotope stable against one-neutron emission is found at $`N=55`$ with NL3 and G1 and at $`N=57`$ with G2. For Sn isotopes (Figure 7b) the shell effect at $`N=82`$ is slightly overestimated by the studied forces. The predictions of the three parametrizations are roughly similar up to $`N=110`$, where the behaviour of NL3 starts to depart from G1 and G2 due to the large quenching of the $`N=126`$ shell effect shown by NL3 as compared with G1 and G2. In our calculations, the odd Sn isotopes become unstable against the emission of one neutron around $`N=110`$. This value is larger than the values found with the non-relativistic SkP interaction in HFB ($`N=103`$) or HF+BCS ($`N=101`$) calculations . The origin of this discrepancy lies in the fact that the shell distribution in tin isotopes around $`N=126`$ for SkP is rather different from that of the relativistic sets .

The one-proton separation energies for the isotones of $`N=28`$ (Figure 9b) show an overall good agreement with the experimental data. The shell effects at $`Z=20`$ and $`Z=28`$ are rather well reproduced by the forces analyzed here. G1 and G2 predict the heaviest nucleus stable against one-proton emission to be 57Cu, as in experiment, while it is unstable in the NL3 calculation. For the isotones of $`N=82`$ (Figure 10b) the shell effect predicted by NL3 and G1 at $`Z=50`$ is similar. Again, as for $`S_{2p}`$ (Figure 10a), NL3 predicts a shell effect at $`Z=58`$ which is not found experimentally, whereas for G1 this effect is clearly smaller and it does not appear for G2. The last stable nucleus against one-proton emission is 151Tm according to the three parameter sets.

### 4.3 One-body densities and potentials

The nuclear densities for chains of isotopes of light and medium size nuclei have recently been studied in the RHB approximation . As $`N`$ grows the neutron and mass densities extend outwards and the r.m.s. radii and the surface thickness increase. Special attention has been payed to the isospin dependence of the spin–orbit interaction. The magnitude of the spin–orbit potential is reduced when one approaches the neutron drip line and, as a consequence, there is a reduction of the energy splittings between spin–orbit partner levels . To our knowledge, for isotones such an study has only been carried out in the $`N=20`$ chain . It is to be remarked that the EFT parametrizations G1 and G2 contain a tensor coupling of the $`\omega `$ meson to the nucleon which plays a very important role in the spin–orbit force because there exists a trade-off between the size of this coupling and the size of the scalar field .

In Figures 12a and 12b we display, respectively, the neutron and proton densities of some $`N=28`$ isotones from $`Z=16`$ to $`Z=32`$ as predicted by the G2 set in our qb-BCS approach. Figures 13a and 13b show the results for some $`N=82`$ isotones from $`Z=40`$ to $`Z=70`$. Since $`N`$ is fixed in each isotonic chain, the spatial extension of the neutron densities is very similar for the different nuclei of the chain. In any case, as one goes from the lightest to the heaviest isotone of the chain, the neutron densities tend to be depressed in the interior region and their surface thickness (90%–10% fall-off distance) shows a decreasing tendency. The proton densities of the isotones exhibit a strong dependence on $`Z`$: by adding more protons they are raised at the interior and their surface is pushed outwards. For $`N=28`$ the surface thickness of $`\rho _p`$ remains roughly constant up to $`Z=28`$ and increases for heavier isotones as a consequence of the growing occupation of the $`1f_{7/2}`$ shell. At the origin the proton densities show a bump when $`Z20`$ because the $`2s_{1/2}`$ level is occupied. The $`Z=16`$ isotone shows a dip at the centre due to, precisely, the emptiness of this $`2s_{1/2}`$ level. For the considered nuclei of $`N=82`$, the proton densities have an approximately constant surface thickness and present a hole at the centre owing to the Coulomb repulsion. In Figure 14 we display the neutron and proton r.m.s. radii of the $`N=82`$ isotonic chain obtained with G2 and NL3. It turns out that the predictions of both sets are very similar. The proton radii increase uniformly with $`Z`$, similarly to the behaviour found for $`N=20`$ isotones with NL3 and the RHB scheme in Ref. . The neutron radii remain roughly constant with $`Z`$. They just show a slight decrease with increasing $`Z`$ till $`Z50`$ and slightly increase afterwards. This behaviour may be related with the shell effect for protons at $`Z=50`$.

The spin–orbit interaction is automatically included in the RMF approximation. It appears explicitly when the lower spinor of the relativistic wave function is eliminated in favour of the upper spinor. This way one obtains a Schrödinger-like equation with a term $`V_{SO}(r)`$ that has the structure of the single-particle spin–orbit potential. Including the contribution of the tensor coupling of the $`\omega `$ meson, the spin–orbit term reads

$`H_{SO}`$ $`=`$ $`{\displaystyle \frac{1}{2M^2}}V_{SO}(r)𝐋𝐒,`$ (14)
$`V_{SO}(r)`$ $`=`$ $`{\displaystyle \frac{M^2}{\overline{M}^2}}{\displaystyle \frac{1}{r}}\left({\displaystyle \frac{d\mathrm{\Phi }}{dr}}+{\displaystyle \frac{dW}{dr}}\right)+2f_v{\displaystyle \frac{M}{\overline{M}}}{\displaystyle \frac{1}{r}}{\displaystyle \frac{dW}{dr}},`$ (15)
where $`\overline{M}=M\frac{1}{2}(\mathrm{\Phi }+W)`$. We have checked numerically that the contribution to the spin–orbit potential of the $`f_\rho `$ tensor coupling of the $`\rho `$ meson is very small, even when one approaches the drip lines. Hence we have not written this contribution in Eq. (15).

The spin–orbit potential (15) for some lead isotopes computed with G2 and NL3 is displayed in Figures 15a and 15b, respectively. As a general trend, for both G2 and NL3, when the number of neutrons is increased the depth of the spin–orbit potential decreases gradually and the position of the bottom of the well is shifted outwards, which implies a significant weakening of the spin–orbit interaction. The same effect arises in other isotopic chains in RHB calculations . Comparing the spin–orbit potentials obtained with the G2 and NL3 sets, one sees that they have a similar strength for all the isotopes analyzed and that the minima of the wells are located at similar positions (slightly shifted to larger values of $`r`$ in G2). The higher effective mass of G2 at saturation ($`M_{\mathrm{}}^{}/M=0.664`$) with respect to NL3 ($`M_{\mathrm{}}^{}/M=0.595`$) is compensated by the tensor coupling included in G2 ($`f_v=0.692`$). To ascertain the relative importance of the tensor coupling we have drawn in the insert of Figure 15a, for 228Pb, the full potential (15) and the contribution resulting from setting $`f_v=0`$ in Eq. (15). We see that the full $`V_{SO}(r)`$ is much deeper and wider. The maximum depth of $`V_{SO}(r)`$ changes from $`68`$ MeV fm-2 (right scale of the insert) to $`44`$ MeV fm-2 when $`f_v=0`$. That is, the tensor coupling accounts for roughly one third of the total spin–orbit strength in the G2 parameter set.

One expects that the weakening of the spin–orbit potential in going to the neutron drip line will bring about a reduction of the spin–orbit splittings

$`\mathrm{\Delta }\epsilon =\epsilon _{nl,j=l1/2}\epsilon _{nl,j=l+1/2}`$ (16)
of the neutron levels . Figure 16 displays the energy splittings of some spin–orbit partner levels of neutrons for lead isotopes, obtained with the G2 and NL3 parameter sets. The splittings predicted by G2 and NL3 are very close as a consequence of the similarity of the corresponding spin–orbit potentials. Partner levels with high angular momentum undergo some reduction in the splitting along the Pb isotopic chain, but partners with small angular momentum show an almost constant splitting. By comparison of their RHB results for Ni and Sn, the authors of Ref.  pointed out that the weakening of the spin–orbit interaction should be less important for heavier isotopic chains. Our calculations for Pb would confirm this statement. All the single-particle levels involved in Figure 16 are bound. Of course, one should not expect the results for $`\mathrm{\Delta }\epsilon `$ to be so reliable in our qb-BCS approach if one, or both, of the partner levels lies at positive energy. The reason is that the single-particle energies of the quasibound levels do not exactly reproduce the energies of the corresponding canonical states of a RHB calculation.

In Figures 17a and 17b we show the spin–orbit potential for isotones of $`N=82`$ from $`Z=40`$ to $`Z=70`$, for the G2 and NL3 parametrizations. Similarly to what is found for isotopes, the results obtained from G2 and from NL3 are comparable and the spin–orbit potential well $`V_{SO}(r)`$ moves outwards with the addition of protons, following the tendency of the proton density. However, for isotones we find that the behaviour of the depth of the spin–orbit potential well is not so monotonous: it increases when one goes from the neutron drip line up to the $`\beta `$-stable region, while it decreases afterwards as more protons are added.

## 5 Summary and conclusion

We have analyzed the pairing properties of some chains of isotopes and isotones with magic $`Z`$ and $`N`$ numbers in the relativistic mean field approach. The study has been performed for the G1 and G2 parametrizations that were obtained in Ref.  from the modern effective field theory approach to relativistic nuclear phenomenology. We have compared the results with those obtained with the NL3 parameter set which is considered to be very successful for dealing with nuclei beyond the stability line.

For accurate calculations of pairing far from the valley of $`\beta `$-stability in the relativistic models, the relativistic Hartree–Bogoliubov approach should be applied. However, we have presented a simpler modified BCS approach which allows one to obtain pairing properties near the drip lines fast and confidently. The method has been used previously in non-relativistic calculations with Skyrme forces . The key ingredient is to take into account the continuum contributions through quasibound levels due to their centrifugal barrier. To further simplify the calculations we have assumed pairing matrix elements of the type $`G=C/A`$ instead of, e.g., a state-dependent pairing with a zero-range force.

The considered quasibound levels are mainly localized in the classically allowed region and decrease exponentially outside it. This eliminates the unphysical nucleon gas which, near the drip lines, surrounds the nucleus when all available positive energy levels are included in the usual BCS approach. Normally, the quasibound levels have high angular momentum and lie close in energy to the corresponding RHB canonical levels. One of the limitations of the qb-BCS approach employed here is the fact that the nuclear density does not follow the asymptotic fall-off of the densities computed with the relativistic Hartree-plus-Bogoliubov theory. In spite of this shortcoming, we have shown by comparison with available RHB results that the qb-BCS approach is able to predict the position of the drip lines, or the behaviour of the neutron and charge radii for nuclei far from the stability line, in a reasonable way. Also, the obtained pairing gaps are nicely scattered around the empirical average $`12/\sqrt{A}`$.

We have applied the qb-BCS approach to the Ni, Sn and Pb isotopic chains, and to the $`N=28`$, 50, 82 and 126 isotonic chains. The two-neutron (two-proton) and one-neutron (one-proton) separation energies, as well as the resulting shell gaps, are similar for the three studied relativistic parametrizations (G1, G2 and NL3) and in general they reproduce the available experimental data, at least qualitatively. The neutron and proton drip lines are usually reached at the same place with the three forces, though one may find a shift of one or two units of $`A`$ among them. We have paid some attention to the quenching of the shell structure near the drip lines. For example, the quenching of the shell effect at $`N=126`$ for Sn isotopes is larger in NL3 than in G1 and G2, while for Pb isotopes none of the three sets exhibits a quenching of the shell effect at $`N=184`$ in our qb-BCS calculation. The EFT parametrizations G1 and G2 contain tensor couplings that are not present in the RMF parametrizations like NL3 and have a larger effective mass at saturation. However, the predicted spin–orbit potentials along the isotopic and isotonic chains do not differ much from those obtained with NL3.

Our analysis shows that the parameter sets based on EFT are able to describe nuclei far from the $`\beta `$-stability line, after adding a phenomenological pairing residual interaction. Only experimental information about some magic nuclei was utilized in the fit of the constants of the G1 and G2 sets and, thus, the results for nuclei near the drip lines are veritable predictions of the model. In spite of the fact that the EFT sets include more couplings and parameters than the conventional RMF sets like NL3, both models reproduce the experimental systematics with a similar quality. In fact, the studied properties away from the valley of $`\beta `$-stability do not seem to provide further constraints on the EFT parameters, not even in the isovector sector. In conclusion, extended sets like G1 and G2 will serve almost the same purposes for normal systems as the conventional parameter sets. However, some of the extra parameters of the general EFT functional may be used to better describe regions of the equation of state at higher density or temperature  without spoiling the systematics for finite nuclei.

## 6 Acknowledgments

The authors would like to acknowledge support from the DGICYT (Spain) under grant PB98-1247 and from DGR (Catalonia) under grant 1998SGR-00011. S.K.P. thanks the Spanish Education Ministry grant SB97-OL174874 for financial support and the Departament d’Estructura i Constituents de la Matèria of the University of Barcelona for kind hospitality.

## Figure captions

Dependence of the binding energy (left) and neutron r.m.s. radius (right) of the nuclei 120Sn and 160Sn on the size of the box used in the calculations. Solid, dashed and dotted lines correspond to a BCS approach including quasi-bound levels, only bound levels and all available levels, respectively. The results are for the G1 parameter set.

Neutron density of 150Sn for different sizes of the box used in the qb-BCS calculations (for the set G1). The dotted line denotes the asymptotic behaviour expected from a Bogoliubov calculation .

Part (a): the two-neutron separation energies for Ni isotopes calculated in the qb-BCS approach are compared with the RHB results of Ref.  and with experiment. Part (b): the neutron pairing energy obtained in the qb-BCS approach. The results are for the set NL-SH.

Charge and neutron r.m.s. radii of Sn isotopes in qb-BCS and RHB  calculations performed with the NL-SH set.

The state-independent pairing gaps predicted by our qb-BCS approach for Ni, Sn and Pb isotopes (top) and for $`N=28`$, 50, 82 and 126 isotones (bottom). The G1 set has been used. The empirical average curve $`12/\sqrt{A}`$  is depicted by a solid line.

Two-neutron (a) and one-neutron (b) separation energies for Ni isotopes computed with the qb-BCS approach for the parametrizations G1, G2 and NL3, in comparison with the experimental data.

Same as Figure 6 for Sn isotopes.

Two-neutron separation energies for Pb isotopes computed with the qb-BCS approach, in comparison with experiment.

Two-proton (a) and one-proton (b) separation energies for $`N=28`$ isotones computed with the qb-BCS approach, in comparison with experiment.

Same as Figure 9 for $`N=82`$ isotones.

Two-proton separation energies for $`N=50`$ (a) and $`N=126`$ (b) isotones computed with the qb-BCS approach, in comparison with experiment.

Radial dependence of the neutron (a) and proton (b) densities of some $`N=28`$ isotones obtained with the G2 set.

Same as Figure 12 for some $`N=82`$ isotones.

Neutron and proton r.m.s. radii of $`N=82`$ isotones obtained with the G2 and NL3 sets.

Spin–orbit potential for some Pb isotopes obtained with the G2 set (a) and with the NL3 set (b).

Energy splitting of some spin–orbit partner levels of neutrons in Pb isotopes, calculated in the qb-BCS approach for the G2 and NL3 sets.

Spin–orbit potential for some $`N=82`$ isotones obtained with the G2 set (a) and with the NL3 set (b).

## Tables


text:
# Inhomogeneous Big Bang Nucleosynthesis: Upper Limit on Ω_𝑏 and Production of Lithium, Beryllium, and Boron

## I Introduction

The possibility that cosmic baryon number fluctuations may have existed on small scales in the early universe has received considerable attention between the late eighties and mid nineties . Such fluctuations in baryon number would have impact on the production of light elements during BBN provided the baryonic mass of individual lumps exceeds $`M_b10^{21}M_{}`$. It was speculated that production of inhomogeneities could result during a first-order QCD phase transition or even possibly during a scenario of electroweak baryogenesis. Initially it was hoped for that such scenarios could make BBN consistent with $`\mathrm{\Omega }_b=1`$ and therefore eliminate the need for “exotic” non-baryonic dark matter. Detailed calculations revealed that inhomogeneous BBN (hereafter, IBBN) scenarios may not be consistent with a universe closed in baryons due to considerable overproduction of 7Li and/or 4He. At present, $`\mathrm{\Omega }_b=1`$ seems also hardly desirable because of a variety of other arguments, such as the cosmological baryon budget , and the success of a structure formation scenario employing cold dark matter, among others.

Recently, observations of the cosmic microwave background radiation (hereafter; CMBR) on intermediate angular scales by the BOOMERANG and MAXIMA balloon missions have achieved unprecedented accuracy . These missions have allowed for a first stab at an estimate of a number of cosmological parameters such as the total cosmic density parameter $`\mathrm{\Omega }_{tot}`$ and $`\mathrm{\Omega }_b`$, among others. Common to both studies is the observation of a relatively suppressed CMBR power spectrum on scales where a secondary peak in the spectrum is anticipated compared to the power at the location of the first peak. Though preliminary, the conclusion of a number of authors is that, within the parameters commonly allowed to be varied, an increased $`\mathrm{\Omega }_b`$ could most easily account for such a suppression . (for alternative explanations cf. to Ref. ). This has lead to the preliminary claim that $`\mathrm{\Omega }_b`$ as inferred from CMBR anisotropy observations may be in conflict with the best estimate $`\mathrm{\Omega }_bh^20.02\pm 0.002`$ from SBBN , in particular, the CMBR data would prefer $`\mathrm{\Omega }_bh^20.03`$  ($`h`$ is the Hubble constant in units of 100 km s-1Mpc-1). Note that even though, at first glance the deviation between these two values seems relatively small, it is clear that a baryonic density parameter of $`\mathrm{\Omega }_b0.03h^2`$ can not be achieved within a SBBN scenario. For such large $`\mathrm{\Omega }_b`$, SBBN production of deuterium can neither account for the deuterium as observed in quasar absorption systems , nor for the inferred D abundance in the presolar nebula and only barely for the D as observed in the local interstellar medium.

Newly developed high-resolution spectrographs (such as UVES on the VLT) allow for a significant increase in the number of stars with claimed detections for the elements 6Li and 9Be. Whereas for a long time there had been only two claimed 6Li/7Li  detections in low-metallicity PopII halo stars, this number is/will rapidly increase in the immediate future. Moreover, 6Li/7Li detections have now also been claimed for disk stars at relatively high metallicities . The preliminary picture which emerges is that 6Li/H abundances in stars are remarkably similar over a wide range in metallicities, though interpretation of the data has to account for the possibility of stellar 6Li astration. Recently, there has been an interesting 9Be/H detection within the atmosphere of a very low-metallicity star . The 9Be/H abundance in this star is higher than expected from extrapolation of the approximately linear 9Be/H versus \[Fe/H\] relation such that this observation may represent tentative evidence for a flattening of the 9Be/H versus \[Fe/H\] slope at metallicities below \[Fe/H\] $`<3`$.

In light of the above, it seems worth reinvestigating BBN with an inhomogeneous baryon distribution. Whereas production of 6Li and 9Be in standard BBN is essentially negligible, it is known that production of these isotopes in IBBN may be significantly enhanced . Furthermore, it should be of interest not only to find out of how much the upper limit on $`\mathrm{\Omega }_b`$ in IBBN may be relaxed compared to that from standard BBN, but also in how much of the parameter space, characterizing the inhomogeneities, IBBN abundance yields may agree with observational constraints.

## II Inhomogeneous Big Bang Nucleosynthesis Calculations

We have performed detailed numerical computations of IBBN by employing the IBBN code described in Ref . This code treats all the relevant baryon diffusion of neutrons, protons, and lighter nuclei. In the Appendix we summarize the employed baryon diffusion constants for protons and neutrons, which includes a reevaluation of some diffusion constants and a correction for mistakes in the literature. The employed code is still the only existing code with a detailed treatment of the effects of photon diffusion and hydrodynamic expansion on the evolution of high-density regions . It is known that these dissipative processes operating at lower temperatures $`T30`$keV may affect the predicted abundances of 6Li, 7Li, 9Be, and 11B in some part of the parameter space (in particular, for compact high-density regions). For example, 7Li produced in form of 7Be may be prematurely destroyed by the reaction sequence 7Be $`(n,p)`$7Li $`(p,\alpha )\alpha `$ when enough neutrons may be delivered to the high-density regions where most of the 7Be is produced. The magnitude of this process depends on the efficiency of hydrodynamic expansion which increases the surface area of high-density regions but also on the correct neutron- and proton- diffusion constants at low temperatures. Note that the distribution (and diffusion) of protons affects the diffusion of neutrons through neutron-proton nuclear scattering .

We have updated the nuclear reaction rates employed in the IBBN code from those based on the compilation by Caughlan & Fowler, as described in Smith et al. , to include the improved charged nuclei induced reactions as compiled by the NACRE collaboration . Note that the modifications in predicted abundances, when the central values of the improved nuclear reaction rates of the NACRE compilation are employed, are fairly small for 2H and 4He but can be in the $``$ 20% - 30% range for 7Li, 6Li, 9Be, and 11B (cf. to Ref. ). Additional uncertainties in the computed abundances arise from appreciable error bars quoted in the NACRE compilations for a few reactions, such as D(p,$`\gamma `$)3He, 3He($`\alpha ,\gamma `$)7Be, and D($`\alpha ,\gamma `$)6Li. In the context of SBBN, these additional uncertainties are of similar magnitude to those quoted above, with the exception of the 6Li abundance which is subject to very large uncertainties of $``$ factor 3-4 in either direction . Though it is beyond the scope of the present work to present a detailed systematic analysis of uncertainties in the prediction of abundances in IBBN scenarios due to reaction rate uncertainties, we will comment below if such uncertainties could impact our main conclusions. Note that all 4He mass fractions $`Y_p`$ presented below are corrected by $`\mathrm{\Delta }Y_p=+0.0049`$ to account for a variety of physical effects as detailed in Ref. .

From the multitude of conceivable initial conditions for the baryon inhomogeneities (including stochastic ones as treated in Ref. ) we chose a regular lattice of spherical symmetric domains, approximating the possible outcome of baryon fluctuations generated during a first-order (e.g. QCD) phase transition around the shrinking bubbles of high-temperature phase. The spherical computation domain is then characterized by its physical length, $`l_{100}`$, specified at temperature of $`T=100`$MeV (specifically, 1 m at $`T=100`$MeV is to be understood as a length of $`5.96\times 10^{11}`$m at the present epoch ). Within this domain we assume a region of high baryon density with baryon-to-photon ratio $`\eta _h`$ occupying volume fraction $`f_V`$ and a low density region at $`\eta _l=\eta _h/R`$ occupying the remainder of the volume, with an initial discontinuity at the boundary of both regions (which softens after some baryon diffusion). This yields an average baryon-to-photon ratio

$$\eta =f_VR\eta _l+(1f_V)\eta _l.$$
(1)
Given these initial conditions there exist still four parameters to be specified, namely $`\eta `$, $`l_{100}`$, $`f_V`$, and $`R`$. The initial parameter space is reduced by assuming $`f_VR=200`$. Physically $`f_VR1`$ corresponds to essentially all baryons residing in the high-density region and none in the low-density region, such that for $`f_VR10`$ one obtains results essentially independent of the exact value of this parameter combination. Though the opposite limit, $`f_VR1`$, may be interesting for the production of significant amounts of isotopes with nucleon number $`A12`$ , in much of the parameter space it yields only minor changes in the D, 4He, and 7Li as compared to a SBBN scenario at the same $`\eta `$. Our calculations employ two different initial “geometries”: (a) spherical condensed - where the high-density region resides at the center of the spherical domain and (b) spherical shell - where the high-density region occupies a shell at the outer edge of the computational domain. These spherical domains are finite-differenced into 24 zones. By increasing the number of zones we estimate that the relative error in predicted abundances does not exceed $`0.5`$% for 4He, $`3`$% for D, and $``$ 10% for the shown isotopes with $`A7`$.

## III Results and Discussion

Figures 1 - 4 show computed abundance yields in IBBN scenarios for the isotopes of D, 4He, 7Li, as well as 6Li, 9Be, and 11B as a function of the length scale of the domains (approximately corresponding to the mean separation between fluctuations) for differing $`\mathrm{\Omega }_b`$ and a wide range of parameters describing the baryon inhomogeneities. The choice of the parameter space for which abundance yields are shown is supposed to bracket most potentially interesting $`\mathrm{\Omega }_b`$ (taking values of $`\mathrm{\Omega }_bh^2=0.012`$, $`0.025`$, $`0.038`$ and $`0.051`$ shown by the solid, dotted, short-dashed, and long-dashed lines in each figure, respectively), as well as to illustrate the general trends of changing “geometry” (Figures 1 and 2 are for spherical condensed fluctuations, whereas Figure 3 and 4 are for spherical shells) and changing volume fractions $`f_V`$ of the high-density regions. We have deliberately not indicated observationally inferred limits on the primordial abundances in these figures, as these are likely to change over the course of time, and since we are more interested in a qualitative understanding of IBBN abundance yields and their potential agreement/disagreement with observationally inferred abundance limits.

It is well known that there exists an “optimum” length scale where the 4He yield may be lower than in a SBBN scenario at the same $`\mathrm{\Omega }_b`$ due to the fact that neutrons may diffuse out of the high-density regions, subsequently decaying in the low-density region before they may be incorporated into 4He. This effect is more pronounced when $`f_V`$ is small (Figures 1 and 3) since back-diffusion of neutrons into the high-density regions is less efficient. Similarly, one finds an “optimum” $`l_{100}`$ where 7Li production is minimized, generally somewhat smaller than the length scale for minimum 4He production. Nevertheless, even at this “optimum” distance the 7Li yields are often higher, at best somewhat lower, than the 7Li yields in a SBBN scenario at the same $`\mathrm{\Omega }_b`$ and they increase with decreasing $`f_V`$ . These trends are due to the two different production mechanisms for 7Li (direct production in the low-density region and production of 7Be in the high-density region) and the relative efficiency of these mechanisms at either lower or higher $`\eta `$ than the approximate $`\eta `$ inferred from SBBN.

There are currently two mutually inconsistent observationally inferred values for the 4He mass fraction, i.e. a high value $`Y_p0.244`$  and a low value $`Y_p0.234`$ . It also becomes more and more appreciated that the inference of $`Y_p`$ from observations of HII regions is subject to systematic errors of possibly considerable magnitude. Much progress has been made in the determination of primordial D abundances in quasar absorption line systems (QASs). There are now several QASs seemingly indicating low D/H $`2.54\times 10^5`$  and only one which favors high D/H $`2\times 10^4`$ . In light of this, one should probably demand from a sucessful BBN scenario to have $`Y_p<0.25`$ (conservative) and low D/H $`25\times 10^5`$. (Note that even within the context of a SBBN scenario, a D/H $`3\times 10^5`$ abundance implies seemingly uncomfortably high $`Y_p0.247`$ . The figures illustrate that one may find considerable IBBN parameter space where these requirements on the primordial 4He and D/H abundances may be met, even for $`\mathrm{\Omega }_bh^2`$ as large as $`0.05`$. This is due to IBBN scenarios often yielding less 4He and more D production than a SBBN scenario at the same $`\mathrm{\Omega }_b`$.

Nevertheless, these considerations disregard observational limits on the 7Li/H abundance. Typical 7Li yields in the IBBN parameter space which agree with observational limits on $`Y_p`$ and D/H strongly depend on $`\mathrm{\Omega }_b`$ (as well as on geometry and $`f_V`$), ranging between about $`3\times 10^{10}`$ and $`10^8`$ for $`\mathrm{\Omega }_bh^2=0.025`$ and $`3\times 10^9`$ to $`3\times 10^8`$ for $`\mathrm{\Omega }_bh^2=0.051`$, This is typically well in excess of the claimed primordial 7Li/H $`1.7\times 10^{10}`$ as derived from observations of lithium abundances in metal-poor halo stars belonging to the Spite plateau. Though 7Li may in principle be depleted in these stars, there are strong arguments against this possibility, such as the claimed absence of intrinsic dispersion of 7Li abundances in stars belonging to the Spite plateau. Further, the more fragile 6Li isotope which should have been astrated as well, is by now observed in a few of these stars . Recently, a primordial 7Li/H even as low as $`1.2\times 10^{10}`$ has been claimed , which results from correcting 7Li abundances for galactic cosmic ray production of this isotope. Nevertheless, even an SBBN scenario (at $`\mathrm{\Omega }_bh^2=0.02`$) yields 7Li/H $`3.8\times 10^{10}`$, in excess of the primordial 7Li/H determination, when one assumes D/H $`3\times 10^5`$ as favored by the QAS data. One thus would have to resort to a small amount of 7Li depletion and/or systematic errors due to, for example, the use of inappropriate stellar atmospheric models, to reconcile these values. If one demands the IBBN yield of 7Li to be at, or below, the quoted SBBN reference value, one finds that the $`\mathrm{\Omega }_bh^2`$ should be below $`0.025`$, precluding a substantial increase of $`\mathrm{\Omega }_b`$ over that inferred from SBBN. Only if one were to relax the 7Li/H limit to about (seemingly unreasonable) $`10^9`$ could $`\mathrm{\Omega }_bh^2`$ in IBBN scenarios be consistent with the currently preferred $`\mathrm{\Omega }_bh^2`$ from CMBR anisotropy measurements $`0.03`$.

We have tested if these conclusions could be changed due to existing reaction rate uncertainties quoted by the NACRE collaboration . We have changed the following reactions to their quoted limits: D(p,$`\gamma `$)3He (lower limit), 3H($`\alpha ,\gamma `$)7Li (lower limit), 3He($`\alpha ,\gamma `$)7Be (lower limit), 7Li(p,$`\alpha `$)4He (upper limit), and D($`\alpha ,\gamma `$)6Li (upper limit). These changes have been “designed” to minimize 7Li production and maximize 6Li production. With these modified rates we have performed two calculations (a) spherical shell, $`\mathrm{\Omega }_bh^2=0.038`$, $`f_V^{1/3}=0.025`$, $`l_{100}=724`$m, and (b) spherical condensed, $`\mathrm{\Omega }_bh^2=0.038`$, $`f_V^{1/3}=0.0125`$, $`l_{100}=32`$m, where the length scales have been chosen close to the “optimum” distance for minimum 7Li production. This has lead to a 7Li/H yield of $`1.45\times 10^9`$ (compared to $`2\times 10^9`$ when the central values of the reaction rates are used) in case (a) and $`1.02\times 10^8`$ (compared to $`1.4\times 10^8`$) in case (b), illustrating that the uncertainty in the predicted 7Li remains within bounds. Nevertheless, a large uncertainty exists in the prediction for the 6Li abundance: we have found a factor $`3`$ and 4 increase in case (a) and (b), respectively.

It should be interesting to explore if the abundance yields of 6Li, 9Be, and 11B in IBBN scenarios may be as large as the abundances of these isotopes observed in the most metal-poor stars to date. Such a comparison could yield, in principle, independent confirmation/rejection of IBBN scenarios. In the parameter space where IBBN yields are consistent with observationally inferred limits on 4He and D/H, we find production of 6Li/H $`10^{14}7\times 10^{13}`$ (including the large reaction rate uncertainty in the D($`\alpha ,\gamma )`$6Li rate), implying a typical “maximum” enhancement factor for this isotope of about 10-30 compared to a SBBN scenario at the same $`\mathrm{\Omega }_b`$. For the 9Be and 11B isotopes one finds ranges of 9Be/H $`10^{18}\mathrm{a}\mathrm{few}\times 10^{15}`$, and 11B/H $`10^{16}10^{13}`$. Typical IBBN yields of these isotopes seem therefore still much below the observed 6Li/H $`7\times 10^{12}`$ , 9Be/H $`5\times 10^{14}`$ , and 11B/H $`10^{12}`$  in the lowest metallicity stars to date where such observations have been performed. These observations are thus inconclusive with regards to a validation of IBBN scenarios. Though it is not easy to completely rule out the possibility that there indeed exist very specific initial conditions for the baryon inhomogeneities which yield primordial production of 6Li, 9Be, and 11B in abundance as high as currently observed in the lowest metallicity stars, it seems clear that this is not the typical case.

In summary, we have performed numerical simulations of BBN in the presence of an inhomogeneous baryon distribution for wide ranges of the parameters describing the inhomogeneities and for a few representive baryon-to-photon ratios. Our choice of initial conditions is limited to scenarios where essentially all baryons are within overdense pockets and the remainder of the volume is initially void of baryons. We found that such scenarios may be consistent with observational limits on the primordial 4He and D abundance for $`\mathrm{\Omega }_bh^2`$ as large as $`0.05`$, however, they result in significantly overabundant production of 7Li with respect to the 7Li/H ratio as observed in stars belonging to the Spite plateau. We note here that similar conclusions have been recently drawn by Ref. . Typical production of 6Li, 9Be, and 11B in such scenarios are found to be still below the abundances of these isotopes observed in the most metal-poor stars to date. Unless 7Li in stars on the Spite plateau has been significantly astrated, which seems unlikely, IBBN scenarios thus do not allow for a significant increase of $`\mathrm{\Omega }_b`$ over that inferred from a SBBN scenario.

We wish to acknowledge several useful discussions with In-Saeng Suh and Naoki Yoshida.

## A Reevaluation of Baryon Diffusion Constants

In this appendix we summarize the baryon diffusion constants which we used in our inhomogeneous Big Bang nucleosynthesis calculations. Some of these diffusion constants have been reevaluated. Such a reevaluation seemed necessary not only since prior work on the subject  yielded partially conflicting results (e.g., the proton diffusion constant due to proton-electron scattering as computed by Applegate, Hogan, & Scherrer and Banerjee & Chitre), but also due to improvement on approximations, such as an energy independent neutron-proton cross section. Furthermore, we correct for the electron diffusion constant due to electron-photon scattering as given in . Rather than going over the partially lengthy details of the calculations we performed, we will state our results, outline by what procedure we obtained them, and highlight the differences to prior evaluations.

Banerjee & Chitre  (hereafter; BC) computed diffusion constants by using the first-order Chapman-Enskog approximation for arbitrarily relativistic particles as thoroughly discussed in the monograph by de Groot, van Leeuwen, & van Weert  (hereafter; GLW). The master equation given in BC for the computation of diffusion constants (i.e., Eqs. (1), (3), and (4) in BC) are not directly evident from GLW but involve a fairly detailed computation. We have therefore redone the calculation of this master equation and arrive at the same result  as BC. Note that the first-order Chapman-Enskog approximation is typically accurate to within 20-30 % .

(a) neutron-electron scattering

At higher temperatures ($`T50100`$keV), and when the local baryon-to-photon ratio ($`\eta `$) is not too large, the diffusion of neutrons is limited by magnetic moment scattering off electrons and positrons. Using the master equations of BC, under the assumption of an energy-independent cross section, and to lowest non-trivial order in the small quantities $`m_e/m_N`$ and $`T/m_N`$, where $`m_e`$, $`m_N`$, and $`T`$ are electron mass, nucleon mass, and temperature, respectively, but for arbitrary $`T/m_e`$, we find

$$D_{ne}=\frac{3}{8}\sqrt{\frac{\pi }{2}}\frac{1}{n_{e^\pm }\sigma _{ne}^t}\frac{1}{z_e^{1/2}}\frac{K_2(z_e)}{K_{5/2}(z_e)},$$
(A1)
in agreement with BC. In this expression $`z_e=m_e/T`$, the quantity $`n_{e^\pm }`$ is the total number density of electrons and positrons, the transport cross section is

$$\sigma _{ne}^t=𝑑\mathrm{\Omega }\frac{d\sigma _{ne}}{d\mathrm{\Omega }}(1cos\theta ),$$
(A2)
to be evaluated in the center-of-mass system, and $`K_n`$ are modified Bessel functions of the second kind and of n’th order, i.e.

$$K_2(z)=\frac{1}{z^2}_0^{\mathrm{}}k^2\mathrm{exp}(\sqrt{k^2+z^2})dk,$$
(A3)
and

$$K_{5/2}(z)=\sqrt{\frac{\pi }{2z}}e^z\left(1+\frac{3}{z}+\frac{3}{z^2}\right).$$
(A4)
The expression Eq. (A1) does agree with that derived by Applegate, Hogan, & Scherrer  (hereafter; AHS) via considering the drag force exerted by $`e^\pm `$ on neutrons and using the Einstein relation. As noted by Ref. , both derivations of the diffusion constant approximate the fermionic occupation number by a relativistic Maxwellian, i.e. $`f=(\mathrm{exp}(E/T)+1)^1\mathrm{exp}(E/T)`$, which nevertheless should only result in a small error. Using $`\sigma _{ne}^t=3\pi (\alpha \kappa /m_N)^2`$, with $`\alpha `$ the fine structure constant and $`\kappa =1.91`$, we may give a numerical value for the diffusion constant

$$D_{ne}=1.87\times 10^4\frac{\mathrm{m}^2}{\mathrm{s}}\left(\frac{T}{\mathrm{MeV}}\right)^{1/2}\frac{1}{(n_{e^\pm }/\mathrm{MeV}^3)}\frac{K_2(z_e)}{K_{5/2}(z_e)},$$
(A5)
where $`n_{e^\pm }`$ is given in natural units, i.e. $`\mathrm{}=c=1`$. The correct density $`n_{e^\pm }`$ is easily obtained from the BBN code.

Recently Suh & Mathews  have considered finite-temperature effects on the neutron diffusion constant. They find a transport cross section due to neutron-electron (positron) magnetic moment scattering which significantly increases over $`\sigma _{ne}^t=3\pi (\alpha \kappa /m_N)^2`$ at low $`T0.5`$MeV. Since finite-temperature effects should vanish in the limit $`T0`$ their result is fairly surprising. We have therefore reevaluated the neutron-electron cross section by using the result given in Ref.  and confirm that neutron-electron scattering is independent of energy for electron energies much below the nucleon mass. Similarly, the authors claim a significant increase of $`D_{ne}`$ at high $`T`$ . Within the context of their analysis, such an increase could only occur due to a change in the electron mass and/or the transport cross section. Nevertheless, according to their own analysis both quantities don’t seem to deviate much from their zero-temperature limits for temperatures below $`T35`$MeV. In light of these inconsistencies we therefore prefer to use the standard AHS and BC results in our calculations.

(b) neutron-proton scattering

Neutron-proton nuclear scattering limits diffusion of neutrons at lower temperatures and/or high $`\eta `$. At the low energies relevant for BBN the scattering cross section is dominated by scattering-angle independent s-wave scattering (zero angular momentum) resulting in a transport cross section which equals the total cross section,

$$\sigma _{np}=\frac{\pi a_s^2}{(a_sk)^2+(1\frac{1}{2}r_sa_sk^2)^2}+\frac{3\pi a_t^2}{(a_tk)^2+(1\frac{1}{2}r_ta_tk^2)^2}.$$
(A6)
In this expression $`k`$ is the nucleon wave vector in the center-of-mass (!) system , and the parameters $`a_s=23.71`$fm, $`a_t=5.432`$fm, $`r_s=2.73`$fm, and $`r_t=1.749`$fm. Only at very low temperatures the above cross section is approximately independent of energy. In that case, one may derive

$$D_{np}=\frac{3\sqrt{\pi }}{8}\frac{1}{n_p\sigma _{np}}\left(\frac{T}{m_N}\right)^{1/2},$$
(A7)
for the diffusion constant, which is a factor of two smaller than the result given in BC. Nevertheless, the energy dependence of the cross section becomes large at $`T50`$keV and should be taken into account. This may be done properly by evaluating the diffusion constant via the master equations in BC with the appropriate cross section Eq.(A6). Following this procedure we derived the cross section

$`D_{np}=2.82\times 10^5{\displaystyle \frac{\mathrm{m}^2}{\mathrm{s}}}\left({\displaystyle \frac{T}{\mathrm{MeV}}}\right)^{1/2}{\displaystyle \frac{1}{(n_p/\mathrm{MeV}^3)}}\times `$ (A8)
$`\times {\displaystyle \frac{1}{I(a_1,b_1)+0.16I(a_2,b_2)}},`$ (A9)
where

$$I(a,b)=\frac{1}{2}_0^{\mathrm{}}dx\frac{x^2e^x}{ax+(1bx/2)^2},$$
(A10)
are integrals to be evaluated numerically. Here the parameters $`a`$ and $`b`$ are given by $`a_1=13.59(T/\mathrm{MeV})`$, $`b_1=1.56(T/\mathrm{MeV})`$, $`a_2=0.71(T/\mathrm{MeV})`$, and $`b_2=0.23(T/\mathrm{MeV})`$. We note that the integral $`I(a,b)`$ converges only slowly against its limiting value $`I(0,0)=1`$ as the temperature is decreased. In that limit, Eq. (A9) converges against Eq. (A7) with $`\sigma _{np}=\pi a_s^2+3\pi a_t^2`$. (c) proton-electron scattering The diffusion of protons in IBBN scenarios is only significant at lower temperatures. As long as Debeye screening of proton charge in the plasma is effective, protons may diffuse independently of the additional (net) electrons required by charge neutrality (cf. electron-photon scattering). In that case, proton diffusion is limited by Coulomb scattering on $`e^\pm `$. Both, AHS and BC compute the proton diffusion constant due to Coulomb scattering $`D_{pe}`$. However, their results differ by as much as a factor of eight at low temperatures. We have therefore recomputed $`D_{pe}`$ by using the master equations given in BC and employing the Mott scattering cross section. To lowest order in $`m_e/m_N`$, and accurate to first order in $`T/m_e`$, we obtain

$$D_{pe}=\frac{3}{4\sqrt{2\pi }}\frac{T^2}{\alpha ^2n_{e^\pm }}\left(\frac{T}{m_e}\right)^{1/2}\frac{1+\frac{15}{8}\frac{T}{m_e}}{\mathrm{\Lambda }+2\frac{T}{m_e}(\mathrm{\Lambda }1)},$$
(A11)
where $`\mathrm{\Lambda }`$ is the well-known Coulomb logarithm with $`\mathrm{\Lambda }\mathrm{ln}(T^2m_e/2\pi \alpha n_{e^\pm })^{1/2}`$. In the limit $`T/m_e0`$ Eq. (A11) reproduces the result of AHS. We conclude that $`D_{pe}`$ as calculated by BC is a factor of eight to small at low $`T`$. Deviations between Eq. (A11) and the result of AHS to first order in $`T/m_e`$ are due to AHS approximating the electron energy in the Mott scattering cross section by $`m_e`$. Eq. (A11) yields for the numerical value of the proton diffusion constant

$`D_{pe}=9.29\times 10^2{\displaystyle \frac{\mathrm{m}^2}{\mathrm{s}}}\left({\displaystyle \frac{T}{\mathrm{MeV}}}\right)^{5/2}{\displaystyle \frac{1}{(n_{e^\pm }/\mathrm{MeV}^3)}}\times `$ (A12)
$`\times {\displaystyle \frac{1+\frac{15}{8}\frac{1}{z_e}}{(\mathrm{\Lambda }/5)+\frac{2}{z_e}((\mathrm{\Lambda }/5)1/5)}}.`$ (A13)
(d) electron-photon scattering With the decrease of temperature thermally produced electron-positron pairs become rare and Debeye screening of nuclear charge becomes inefficient. In this limit, electric forces which would rapidly be build up if the proton- and (net) electron- distributions differed, prevent the independent diffusion of these two species . One may show, by evaluating the electric fields which would be present due to differing proton- and electron- distributions in the presence of Debeye screening, and by comparison of the resulting proton flux due to the electric fields with the flux of protons due to diffusion, that protons may only diffuse independently, when $`n_{e^\pm }n_e^{}n_{e^+}`$. When this is not the case, electrons and protons diffuse together by ambipolar diffusion, with the effective diffusion constant given by twice that of the larger of electron- and proton- diffusion constants . Electron diffusion is rendered fairly inefficient due to Thomson scattering of electrons on the cosmic background photons. The diffusion constant may be computed by considering the drag force on an electron due to a photon blackbody

$$f_{\mathrm{drag}}=\frac{4}{3}\frac{\pi ^2}{15}\sigma _{\mathrm{Th}}T^4v,$$
(A14)
with $`\sigma _{\mathrm{Th}}6.65\times 10^{29}\mathrm{m}^2`$ the Thomson cross section and $`v`$ the velocity of the electron in the cosmic background photon rest frame. Note that Eq. (A14) is given in natural units. Using Eq. (A14) together with the Einstein relation $`D=Tb`$, where the mobility $`b`$ is defined as the proportionality constant between the terminal velocity $`v`$ which a particle reaches in a plasma when an external force $`f`$ is applied, i.e. $`v=bf`$, one finds for the effective proton diffusion constant

$$D_p^{\mathrm{eff}}=2D_e=7.86\times 10^2\frac{\mathrm{m}^2}{\mathrm{s}}\left(\frac{T}{\mathrm{MeV}}\right)^3,$$
(A15)
applicable when $`n_{e^+}<n_e^{}`$. We note here that Eq. (A15) is different from the simple estimate given in Ref. . We stress that the neglect of $`D_p^{\mathrm{eff}}`$ due to Thomson scattering, typically important at low $`T40`$keV, may lead to errors in the calculated 7Li abundances by more than an order of magnitude.