Isospin-rich nuclei in neutron star matter

Stability of nuclei beyond the drip lines in the presence of an enveloping gas of nucleons and electrons, as prevailing in the inner crust of a neutron star, is studied in the temperature-dependent Thomas-Fermi framework. A limiting asymmetry in the isospin space beyond which nuclei cannot exist emerges from the calculations. The ambient conditions like temperature, baryon density and neutrino concentration under which these exotic nuclear systems can be formed are studied in some detail.


I. INTRODUCTION
In the outer region of the core of a neutron star, at densities ∼ 1.5 × 10 14 g/cm 3 or more, nucleons are distributed uniformly forming a homogeneous system [1][2][3][4].
At lower densities, in the inner crust of the star, inhomogeneities appear and then the positive charges get concentrated in individual nuclei of charge Z which are embedded in a sea of neutrons, electrons and possibly neutrinos, the whole system being charge neutral. The nuclei arrange themselves in a periodic lattice [5] due to the electrostatic interaction. In the typical conditions prevailing in the neutron star interior, the nuclei are in complete thermodynamic equilibrium with the environment and are assumed to be in β-equilibrium. They may well be beyond the neutron drip line [1,6] known for laboratory nuclei with (N − Z)/A ≤ 0.35. The excess neutrons of the very neutron rich nuclei present in the inner neutron star crust, which would otherwise decay under terrestrial conditions, are held together in stable equilibrium by the pressure exerted by the surrounding neutron sea.
The properties of the isospin-rich nuclei may be quite different from those of the terrestrial nuclei. For example, with increasing density, the nuclei may pass through different exotic shapes, namely from spheres to cylinders, slabs, cylindrical holes and spherical bubbles [7,8]. The presence of nonspherical nuclei could affect significantly the pinning of vortices and neutrino emission from neutron stars [9]. The existence of these exotic-shaped nuclei is, however, model-dependent [3,10]. The external gas surrounding the nuclei may also influence the density distribution of these nuclei and thus may modify their size [1].
The nuclear equation of state (EOS) plays the pivotal role in determining the macroscopic properties of a neutron star, such as its mass, radius and moment of inertia [11]. In this context, the presence of different nuclear species in the outer and inner crustal regions provides the starting point for the consideration of some important aspects, such as its superfluid and elastic properties. On the other hand, for nuclei immersed in a nucleonic gas there are two basic concerns from a microscopic viewpoint. One, a thermodynamically consistent treatment of the coexistence of the two phases of nuclear matter (namely, the nuclear liquid and the surrounding gas) [12][13][14], the other, a plausible description of the interface between the liquid and the gas [15][16][17]. It has been shown that this problem can be taken care of by solving the density profile in the subtraction procedure of Bonche, Levit and Vautherin (BLV) [18,19]. In this method, the density profile of the liquidplus-gas system and that of the gas are solved in a self-consistent procedure and then the extensive observables referring to the liquid are obtained as the difference between those of the two solutions. The influence of the surrounding gas on the surface energy is then automatically taken into account. The same concept can be applied to situations at zero temperature where drip nucleons occur and the system coexists with an outer nucleonic phase, as shown for very asymmetric cold semi-infinite nuclear matter in Ref. [17].
In this paper we focus on the effect of the external gas on the structure and stability of finite nuclei with large neutron excess, to have a broader understanding of the conditions under which these nuclei may exist in neutron star matter.
A preliminary study for isolated nuclei immersed in a nucleonic gas at zero temperature has been done recently [20]; in the present work we extend these ideas in the appropriate astrophysical context. We also extend the calculations to non-zero temperature which is relevant at the formation stage of the neutron stars.
We assume that in the nuclear matter at sub-nuclear density the nuclei are located in a lattice. In order to simplify our calculation, the Wigner-Seitz approximation [5] is applied, where each lattice volume is replaced by a spherical cell of radius R c , the nucleus being located at its centre. The matter in each cell is taken to be charge neutral, ie, the number of electrons is equal to the number of protons in the cell. The neutrino density in this cell is determined from the β-equilibrium condition. For a given average density of nuclear matter with a certain proton concentration, the density distribution is solved self-consistently. The calculations are performed in the finite temperature Thomas-Fermi formalism. The matter in the cell in which the protons and neutrons coexist does not define the nucleus itself, the nucleus is identified after subtraction of the gas part generated self-consistently, as in previous investigations of excited nuclei or asymmetric semi-infinite matter [17][18][19].
In Section II, the model employed in the calculations is introduceed. The results and discussions are presented in Section III. Section IV contains the concluding remarks.

II. MODEL
In the inner region of the crust of a neutron star, we consider a mixture of neutrons, protons (some of which may exist as bound nuclear clusters), electrons and neutrinos in thermodynamic equilibrium. We contemplate both cold matter as well as matter at a finite temperature. We ignore the contributions from alpha particles and also from photons. We further ignore the plasma effects. All these effects are known to be rather small [21]. The nuclear clusters are assumed to be arranged in a body-centered cubic lattice which we approximate by Wigner-Seitz (WS) cells defined as spheres with radius R c . Each cell is assumed to be electrically neutral and interactions among the cells are neglected. At the densities that we are interested in, the Fermi momenta of the electrons are much larger than their rest mass, the electrons are then extremely relativistic and can be assumed to be uniformly distributed in the cell. We assume the matter to be in β-equilibrium, ie, the chemical potentials of the considered species fulfil At a fixed temperature, under the condition of charge neutrality, the variables in the problem are the average baryon density ρ , the proton fraction Y p = Y e (the electron fraction), the neutrino fraction Y ν and the lattice radius R c . With the βequilibrium condition given by Eq. (1), if µ n , µ p and µ e are known, then µ ν can be determined which in turn yields ρ ν when the neutrinos are taken to be a degenerate Fermi gas. Then only three of the variables are left independent. The condition of β-equilibrium of the baryons with the electrons and neutrinos guarantees that the system has the minimum free energy. If one imposes further the constraint that the lattice sites contain a particular nuclear species with a given charge Z cl and baryon number A cl (immersed in a nucleonic gas or not), then there is only one independent variable.
To obtain the thermodynamic properties of the system of baryons, electrons and neutrinos, we minimize the free energy of the system in the WS cell with the constraint of number conservation of the individual species. Under the conditions of charge neutrality and β-equilibrium, the relevant grand potential is given by where F is the free energy and the index q = (n, p) refers to neutrons and protons.
The free energy has the following expression: In Eq. For the nuclear energy density, we use the Skyrme energy density functional. It is written as where ρ = ρ n + ρ p and the effective mass of the nucleons m * q is defined through In the numerical calculations we shall employ the SKM* interaction, whose parameters can be found in Ref. [22]. In the Thomas-Fermi approximation the kinetic energy density τ q is given by The fugacity η q is obtained as where V q is the single-particle potential experienced by nucleons (including the Coulomb part for the protons). The nucleonic density ρ q is related to η q by The functions J k (η q ) in Eqs. (7) and (9) are the standard Fermi integrals.
For the entropy density of the nucleons one has The direct part E d c of the Coulomb energy density E c of the charged particles is while the exchange part is computed in the Slater approximation: For f e and f ν we use the standard expressions [23]. For the minimisation procedure, we take recourse to the finite temperature Thomas-Fermi approximation which yields Eq. profile is found to be practically uniform which is identified as the low density gas.
To delineate the nucleus from the embedding low-density environment, we follow the BLV procedure which has been used earlier in the Hartee-Fock (HF) framework [18] as well as in a Thomas-Fermi (TF) scheme [19] for a hot nucleus in coexistence with the vapour surrounding it. The method is based on the existence of two solutions to the HF or TF equations, one corresponding to the liquid phase with the surrounding gas (LG) and the other corresponding to the gas (G) alone.
For an isolated hot nucleus in equilibrium with the gas surrounding it, in the absence of Coulomb forces, the densities ρ q LG and the gas density ρ q G can be obtained from the variational equations where Ω N LG and Ω N G are the grand potentials of the respective systems with the same chemical potential. In our procedure, the solution for ρ L = ρ LG − ρ G , which may be called the liquid profile, is independent of the box volume in which the calculation is done. In the presence of the Coulomb force, however, the Coulomb repulsion increases with the box size, driving the protons to the edge of the box and leading ultimately to a divergence problem. In order to work out a convergent prescription, Bonche et al. [18] calculated the Coulomb energy E c from the subtracted proton density ρ p LG −ρ p G and obtained the density profile from the variation of the subtracted grand potentialΩ = Ω N LG − Ω N G + E c with respect to both ρ q LG and ρ q G . In the astrophysical context, the situation is however different due to the presence of electrons. The Coulomb energies for the LG and G phases are The direct part of the single-particle Coulomb potentials (δE c LG /δρ p LG and δE c G /δρ p G ) for the LG are the G solutions is the same. It is given by Since the system is charge neutral, the divergence problem does not arise. The solutions to the density profiles can be directly obtained from the variations of the total (with Coulomb) grand potentials Ω LG for the liquid-plus-gas phase and Ω G for the gas phase with respect to ρ q LG and ρ q G , respectively. The resulting coupled equations are At zero temperature they reduce to In Eqs. (18)-(21), V q LG and V q G refer to the nuclear part of the single-particle potential corresponding to the liquid-plus-gas and the gas solutions. V c LG (V c G ) is the sum of the direct part of the Coulomb single-particle potential given by Eq. (17) and the exchange term −e 2 (3/π) 1/3 (ρ p LG(G) ) 1/3 for the protons in the LG (G) phase. At zero temperature, if the cluster is within the drip lines, the cell does not contain any nucleonic gas. Equations (20)-(21) then automatically lead to a gas solution which is zero throughout.
The nucleonic chemical potential is determined from the constraint of nucleon number conservation. In the case of a given nuclear cluster with neutron number N cl and proton number Z cl embedded in a nucleonic gas, the conservation of the number of nucleons in the cluster requires [from Eqs. (18) and (19)] Here A q refers to N cl or Z cl . An equation similar to (22) follows from Eqs. (20) and (21) at zero temperature.
On the other hand, when ρ is given in a WS cell with a given proton concentration, the total number of neutrons N and protons Z in the cell is defined. Then, the chemical potential is obtained as where now A q refers to N or Z.

III. RESULTS AND DISCUSSION
In an earlier paper by some of the present authors [20], calculations at zero temperature on some structural properties of isolated nuclei much beyond the drip line were reported. In that case, the excess pressure exerted by the enveloping gas stabilises the nucleus even beyond the nominal drip lines (which are defined by µ n = 0 or µ p = 0 in the Thomas-Fermi approach). We now place these calculations in a broader context related to the environment existing in the inner crust of a neutron star. We take cognizance of the presence of electrons and neutrinos existing in β-equilibrium with the neutrons and protons. We extend the calculations also at finite temperature. As mentioned earlier, the calculations are performed in the finite temperature Thomas-Fermi framework in a Wigner-Seitz lattice employing the SKM* interaction. For completeness, we also discuss the situation for the protonrich nuclei though they may not exist in the crustal matter of the neutron stars.
In the calculations on asymmetric infinite and semi-infinite nuclear matter in equilibrium with a drip phase, it was observed [12,16,17] that the neutron-proton asymmetry could be increased arbitrarily till the two phases merge into a uniform system when the densities and the asymmetries of both phases become equal. On the contrary, in the self-consistent TF calculation for isolated finite nuclei, it was found that one cannot add or remove neutrons from nuclei arbitrarily [20]. There exists a limiting neutron-proton asymmetry I = (N cl − Z cl )/A cl beyond which the system becomes thermodynamically unstable; establishing chemical equilibrium between the nuclear phase and the gas phase further becomes impossible.
In the realistic conditions considered in the present investigation, including the lattice effects, we calculate the stability limits at T = 0 as shown in the upper panel of Fig. 1. We find that both the limiting asymmetry (full line)  With increasing temperature, the central density is depleted with the appearance of a thickening tail. The gas densities ρ q G are shown in the bottom panels and the liquid densities ρ q L obtained as the difference between ρ q LG and ρ q G are shown in the central panels.
The density profile ρ q L of the isolated cluster 330 Pb is found to be independent of the size of the box. However, in a lattice, because of the modification of the Coulomb force in the presence of the electrons, the nucleonic densities in the cluster are boxsize dependent, though weakly. The liquid density vanishes at a distance r < R c even at T = 6 MeV. The neutron gas density is found to be practically constant throughout the box; at finite temperature, the gas density is larger as expected.
For this neutron-rich system, at T = 0, the proton gas does not exist; however, at finite temperature (T = 6 MeV), a few protons are present in the gaseous state. In finite temperature calculations for isolated nuclei, the proton gas density profile is strongly polarised [18] due to the repulsion from the nuclear core. In the presence of the lattice electrons, the density polarisation exists only in the vicinity of the core, at further distances the proton gas density is found to attain a nearly constant The proton fraction Y p as a function of the mass number is shown in the middle panel of Fig. 7. A kink is seen to occur at A cl = 143, the proton drip point below which Y p rises faster because of the presence of the proton gas. For A cl > 143, the charge Z cl is fixed. Hence, with increase of mass number, Y p decreases linearly up to the neutron drip point beyond which the slope of the curve changes as neutrons appear in the drip phase. With temperature, because of evaporation of nucleons taking place, the kinks smoothen out, but the overall qualitative features are the same.
The lower values of Y p at finite temperature reflect the predominance of neutrons in the gas for isotopes that are not very proton-rich. The variation of the neutrino fraction Y ν (= ρ ν / ρ ) with the isotope mass is displayed in the bottom panel. With increasing mass number, it is found to decrease. The effect of temperature on Y ν is very weak. Inspection of Fig. 7 reveals that for a particular ρ , A cl is generally double-valued. If the isotope is proton-rich, Y p is necessarily large and Y ν also comes out to be large. For the same value of ρ , the neutron-rich isotope is formed when Y ν is small.
So far our calculations have been done by fixing a given cluster (N cl , Z cl ), which yields various sets of values for the parameters of the problem. A common practice is, however, to fix the parameters like ρ and Y p and study the change in the internal structure of neutron star matter. When β-equilibrium is assumed, for a given value of ρ and Y p , the free energy density (equivalent to free energy per baryon in the present case) in the system of nucleons, electrons and neutrinos is a minimum in any given cell size R c . The cell size sets the periodicity in the neutron star matter; it is not a priori known and therefore is taken as a parameter. Even for given values of ρ and Y p , the neutrino concentration Y ν (under β-equilibrium conditions) would be different for different values of R c . One may look for a global minimum in the free energy density by varying Y ν (through variation of R c ). This is displayed in  Fig. 9. Here Y p is fixed at 0.3. Since these calculations are carried out at zero temperature, the total energy per baryon e tot = F/A is given by The different contributions e N , e lat , e e and e ν are the nuclear energy including the Coulomb interaction among the protons, the lattice energy (ie, the Coulomb interaction energy due to the presence of the electrons), the electron kinetic energy and the neutrino energy per baryon, respectively. It is seen that e N has a minimum at a lattice radius R c somewhat smaller compared to R min c , the value of R c at which F/A is minimum in Fig. 8 at the corresponding ρ . This may be understood from the variation of the cluster size A cl with R c as shown in the middle panels of the figure. The cluster size increases monotonically with R c and for these neutronrich clusters the minimum binding energy per nucleon occurs at around A cl ∼ 100 (for nuclei on the β-stability line, this occurs for 56 Fe). The lattice energy e lat decreases monotonically with R c . This is because for given Y p , the proton number and hence the electron number increases with increasing lattice size keeping the electron density constant, as a result of which e lat ∼ −A 2/3 [5]. The electron kinetic energy per baryon e e is a constant for a given density and Y p and therefore is not to decrease [26] with an increase in the so-called freeze-out volume or a lowered average density ρ . Figure 10 also tells us about the appearance of different nuclear clusters as the neutron star evolves in the formation stage. Initially, the temperature T may be as high as ∼ 10 MeV and then all the nucleons will likely be in the gas phase. As time flows by, the system cools down and seeds of nuclei with increasing size start appearing from the gas.

IV. CONCLUSIONS
We have investigated the structural properties of nuclei and nuclear matter as can be found in the environment of the inner crust of a neutron star. These nuclei are dipped in a sea of electrons, neutrinos and of low-density nucleons. Because of the pressure exerted by the surrounding nucleonic gas, nuclei may exist even far beyond the nominal nuclear drip line. Following the BLV subtraction procedure, we give a prescription in a Thomas-Fermi framework to properly isolate the nucleus from its environment. This is extremely helpful in further exploring the limits of stability.
From our calculations, a limiting asymmetry in the neutron-protron concentration emerges beyond which a nucleus even within the gaseous environment ceases to exist. A delicate balance between the Coulomb force and the diluted surface tension with increasing asymmetry and increasing density of the environment possibly plays a pivotal role here. We also investigate the ambient conditions, namely the average density, electron or proton fraction and the neutrino concentration, under which nuclei of a particular species can be formed at different temperatures in the stellar matter. We furthermore study the thermal evolution of the nuclear clusters at different densities, which may serve as a guide to understand the nucleation process of different species of nuclear clusters from the nucleonic sea as the neutron star cools down in the earlier stages of its formation.
We have left our calculations on a simple pedestal. That is, we have not included the extension of the Thomas-Fermi framework, we have not taken the shell-effects into account, possible plasma effects have been ignored and, similarly, the presence of α-particles at low densities has not been taken into consideration. We have worked with the SKM* interaction; at large asymmetry, its validity is not unquestionable.
Sophistications in the approach or the use of more suitable interactions may change the results somewhat, but the qualitative features that emerge from our calculations, we believe, will remain mostly unchanged.