Skyrme-Hartree-Fock calculations of nuclear properties in the drip-point region of neutron star crust
aa r X i v : . [ nu c l - t h ] D ec Skyrme-Hartree-Fock calculations of nuclear properties in the drip-pointregion of neutron star crust
Uwe Heinzmann, Igor N. Mishustin,
1, 2 and Stefan Schramm†
1, 31
Frankfurt Institute for Advanced Studies,Goethe University, D-60438 Frankfurt am Main, Germany National Research Center Kurchatov Institute, Moscow 123182, Russia Institute for Theoretical Physics, Goethe University,D-60438 Frankfurt am Main, Germany (Dated: January 1, 2021)
Abstract
In the present paper we explore the neutron-drip region of cold non-rotating isolated neutron stars. We haveperformed extended nuclear-structure calculations for nuclei embedded in the electron gas. For modeling theouter crust we use a set of Wigner-Seitz cells, where every cell contains one nucleus surrounded by a cloud ofrelativistic electrons. Above the drip point a non-relativistic neutron gas occurs in the cell. These calculationsare carried out within the Hartree-Fock approach in combination with Skyrme effective interactions. For everybaryon density we have determined the configuration with a minimal total energy. The drip elements andcorresponding drip densities have been determined for about 240 different parametrizations of Skyrme forcesused in the literature. We demonstrate that the calculated drip-point densities depend essentially on the Skyrmeparametrization used. Even the drip elements and the occupied shells in the drip region differ for differentparametrizations. We have found that the number of neutrons building the neutron gas at the drip point alsodepends essentially on the Skyrme force chosen. Nevertheless, the number density of the neutron gas in thedrip-region is more or less the same ( ∼ − neutronsf m ). The drip densities obtained within our approach aregenerally lower than predicted earlier. † deceased . INTRODUCTION After detecting gravitational waves from neutron-star mergers, compact objects are now again in thefocus of interest. From these catastrophic events one can obtain important information about the prop-erties of strongly-interacting matter in a very broad domain of baryon densities and temperatures. Inthe present paper we focus on cold isolated non-rotating neutron stars (NS), which have mass typ-ically of the order of a solar mass ( M ⊙ = , · kg) and radii of about 10 kilometers, see e.g.review [1]. They have only a small atmosphere and a thin solid crust forming the outer layers. Thethickness of the crust is about 1 km and the mass contained in the crust is about one percent of thetotal mass. But nevertheless, its structure is relevant for the interpretation of many observational data.For example glitches, sudden changes of the pulsar periods, are thought to be caused by the break-down of the crust due to the slowing down of the rotation. Even more violent destruction of the crustis expected in NS-NS merger events, see [2]. Also the cooling rate of magnetars is rather sensitive tothe composition and the thickness of the crust, see e.g. [3]. Generally, the whole information comingfrom inner layers of the neutron star is filtered by the crust material, see more in refs. [4], [5], [6].One usually divides the crust into two regions: the outer crust and the inner crust. The outer part of thecrust has a crystalline structure consisting of more or less spherical nuclei surrounded by electrons.The inner crust is defined as the transition layer between the neutron drip density ρ drip and the homo-geneous nucleonic matter (outer core) at about half of the nuclear saturation density, ρ ≈ .
16 fm − .At density ρ drip the neutron chemical potential becomes equal to m n c and neutrons drip out from thenuclei. At higher densities the nuclei are surrounded by a neutron gas whose pressure together withelectrons acts against gravity.In this paper we use fully microscopic Hartree-Fock approach to determine corresponding drip nucleiand their neighbours. These calculations are done for a large variety of Skyrme effective interactionsused in the literature. More specifically, we have considered about 240 different Skyrme parametriza-tions, which have been selected and analyzed in ref. [14]. Most previous calculations were done forone specific Skyrme force without analysing the differences. According to our results, the neutrondrip point is more likely between 3 , · gcm and 4 , · gcm than between 4 , · gcm and5 , · gcm , as found previously. Assuming a neutron drip at a lower density would cause a thinnerouter crust and, consequently, a thicker inner crust of the neutron star. Obviously, this may lead to sig-nificant phenomenological consequences, such as the tidal deformability of neutron stars in mergingevents, see e.g. refs. [15], [16]. 2 I. PREVIOUS ESTIMATES OF DRIP-POINT DENSITYA. Early estimates
First estimates of the drip density and drip elements have been done in ref. [7]. The authors de-scribed the outer crust within a thermodynamic approach using a phenomenological equation of stateobtained by extrapolating of known nuclear mass data. By minimizing the total energy of the systemof the nucleus, the electrons and their interaction (lattice energy), they estimated the drip density tobe 4 , gcm with the drip nucleus Kr .Later Negele and Vautherin [8] described the properties of the inner crust of a neutron star usingenergy-density functional method combined with the Wigner-Seitz cells [9, 10]. Within this approachthe crystal structure of the crust is approximated with a set of independent spherical cells, each con-taining one nucleus (N,Z) and Z electrons. Electrons form a degenerate Fermi gas, which becomesfully relativistic at a density of ρ ∼ gcm . Hence, at densities around the neutron drip point thecrystal consists of more or less spherical nuclei immersed in an uniform gas of ultra-relativistic elec-trons. B. Estimates using Liquid-drop model
Below we ignore small corrections due to deviations from spherical symmetry of electron distribu-tions. The volume per nucleus n N is defined using a sphere of radius r C : n N = π r C , (1)where r C is in the order of magnitude of the nuclear spacing. The electron density n e is defined in asimilar way as n e = π r e , (2)defining a corresponding length r e . Because of electrical neutrality the equation n e = Z · n N (3)must hold. Therefore the two lengths are correlated as r c = Z · r e . (4)3o get a rough estimate of the neutron drip density, we use a simple liquid drop model, as was firstdone in ref. [11]. The total energy of a nucleus with mass number A and proton number Z = x · A canbe expressed as E nucltot ( A , Z ) = ( A − Z ) · m n + Z · m p − a V · A + + a sym · ( − x ) · A + a S · A / + a C · x · A / , (5)where A = Z + N , m n and m p are the rest masses of neutron and proton, a V , a sym , a S and a C are thecorresponding empirical coefficients. We have used units ¯ h = c = Z electrons can be writtenas E tot A = x · ε e + E nucltot ( A , Z ) A , (6)where ε e is the average electron energy per particle. To find the optimum nucleus (A*,Z*) one shouldminimize this expression, first with respect to A , and then, with respect to x . The first minimizationyields a s A / = · a C x A / or A ∗ = a s a C x , (7)meaning that the optimal nucleus has surface energy per nucleon equal twice its Coulomb energy.Here one can see that with increasing nuclear mass, i.e. A , the proton fraction x = ZA of the optimalnucleus decreases. The next step is to find the optimal proton fraction. After substituting the expres-sion for A ∗ into E tot / A , one can differentiate it with respect to x at fixed N . Rewriting the sum ofsurface and Coulomb energies per nucleon as ( ) · (cid:0) a s a C x (cid:1) / , one finally gets µ e + m p − m n − a sym · ( − x ) + (cid:0) a C a s · x (cid:1) / = , (8)were µ e = ∂ E e ∂ N e is the electron chemical potential. This equation is another representation of the ther-modynamic equilibrium condition µ e + µ p = µ n , where µ p and µ n are the proton and neutron chem-ical potentials including the corresponding rest masses. This simply means that under the thermody-namic equilibrium, electron capture and neutron decay are in chemical equilibrium. Now this equationrelates the electron chemical potential to the proton fraction. At densities about 10 gcm , where theneutron drip starts, the electrons can be treated as an ultra-relativistic Fermi gas, e.g. µ e = p Fe , where p Fe is the electron Fermi momentum, related to the electron density as n e = π ( µ e ) = x · n , where n = AV is the nucleon density. Substituting this in (8) gives x as a function of n or alternatively of themass density ρ = m · n , where m = ( m n + m p ) is the average nucleon mass. The physical condition4or dripping out of neutrons is that the neutron chemical potential equals the neutron’s rest mass, i.e. µ n = m n . The corresponding chemical potential is obtained by differentiating the nuclear energy den-sity with respect to the total neutron number density, n n = ( − x ) · n . The result of this differentiationis µ n = ∂ E nucltot ( A , Z ) ∂ N ! Z = m n − a V + a sym · ( − x ) + ( a C a s ) / · x / . (9)Therefore, the condition for neutron drip reads a V = a sym · ( − x drip ) + ( · a C a s ) / x / drip . (10)Neglecting the surface term and Coulomb corrections (second term) yields x drip = ≃ (cid:16) − a v a sym (cid:17) .Using the values a V = . a sym =
23 MeV from ref. [12], one gets x drip = . A drip ≃
144 and thus Z drip ≃
41, which is the
Niobium isotope Nb . Pethick and Ravenhall [11] used somewhat different parameters a V ≃
16 MeV and a sym ≃ x drip = .
32. Therefore, they got A drip ≃
122 and Z drip ≃
39 which isthe
Yttrium isotope Y .By using µ e from the exact equation (8), after inserting a sym =
23 MeV, a C = .
715 MeV and a s = . µ e = .
544 MeV, which is ∼ µ e =
25 MeV. Using the parameters of the liquid-drop modelfrom ref. [12] a V =
16 MeV and a sym =
24 MeV one gets n e = . · − fm and r e = .
81 fm,that is about 20% lower than the value quoted by Pethick and Ravenhall, r e = .
14 fm.The corresponding mass density is calculated as ρ = ( A · m ) · n e Z = m · AZ · π r e . (11)Finally we get for the mass density at neutron drip point ρ drip = . or 6 , · gcm ,which is significantly higher than the result reported in [11], 3 , · gcm . Comparing these resultsone can conclude that even the parameters of the liquid-drop model do have astrophysical relevance.The situation is even more striking for calculations using multiple versions of Skyrme effective inter-actions, as will be demonstrated in Sect. III. 5 . Symmetry energy based estimates A first rough estimation of the neutron drip density using the symmetry energy S has been done inref. [13]. The authors use a simplified version of the nuclear mass formula neglecting Coulomb,surface and all other finite-size terms. Keeping only quadratic term in δ = ( N − Z ) / A one can writethe energy per nucleon in nucleus (A,Z) as E N ( A , Z ) A ≃ ε + S · δ , (12)where ε is the bulk energy per nucleon with respect to the rest mass and S is the symmetry energycoefficient for symmetric nuclear matter, both calculated at the saturation density. The neutron-protonmass difference is neglected. Following this approximation, the neutron chemical potential withrespect to the neutron mass is µ ′ n = ε + ( δ + δ ) S and the corresponding proton chemical potentialis µ ′ p = ε + ( − δ + δ ) · S . Now the value of δ D corresponding to the neutron-drip density ρ D canbe calculated from the condition µ ′ n = δ D = q − ε S −
1. For example, using values for theSkyrme force SKM* from ref. [14] ε = − .
77 MeV and S = .
03 MeV, we get δ D = . β -equilibrium condition we have µ e = µ n − µ p ≃ · S · δ . On the other hand, for relativisticelectrons µ e ≈ p Fe = (cid:0) π n b · ZA (cid:1) / . Using this identity and converting the number density into massdensity ρ one gets for µ e : µ e = . · ( ρ · ZA ) , (13)where ρ ≡ ρ , see details in ref. [13]. Computing x D by inserting δ D in x D = ZA = ( − δ D ) = . ρ D = (cid:18) · S · δ D . (cid:19) · x D . (14)Putting in the values for SKM*, one obtains ρ drip ≃ , · gcm , which is surprisingly close tothe drip density calculated in ref. [7] (BPS), 4 , · gcm , being now a commonly accepted value.With the nuclear saturation density n = .
16 fm − equivalent to ρ = . · gcm even BBP ( seeref. [17]) obtained in their calculations for the drip density ρ drip = , · gcm .More recently a systematic analysis of the outer crust structure in the vicinity of the drip pointhas been carried out in ref. [18]. A set of different nuclear models for the nuclear equation of statehas been used within a thermodynamic approach described in ref. [17]. The predicted densities forthe neutron drip point were found around the value ρ drip = · gcm . In the drip region the authorshave found a number of elements with proton numbers from Z=34 (Selenium) to
Z=38 (Strontium).6
II. REALISTIC SKYRME-HARTREE-FOCK APPROACHA. General remarks
The neutron drip area of neutron star crust is in the focus of our present research. For the description ofnuclei we use Hartree-Fock approach in combination with Skyrme effective interactions. We considerspherical Wigner-Seitz cells containing a nucleus immersed in relativistic electron gas, and performthousands of Skyrme-Hartree-Fock calculations running through all possible combinations of protonsand neutrons. After that we determine the ground state energy of the cell including the energy ofelectron gas as a function of ρ . The most stable nuclei are determined in a broad density range up tothe neutron drip point. We have found significant differences in predicted drip elements for differentSkyrme parametrizations. These uncertainties should be taken into account when calculating crustproperties and cooling rates. B. Skyrme energy-density functionals
Our goal is to model the neutron star crust within the framework of a Skyrme-Hartree-Fock (SHF)approach with BCS-pairing. We assume that only ground-state nuclei are present in the crust ofneutron stars, and consider below only even-even nuclei. There are different possibilities to builda Skyrme energy density functional which are discussed in ref. [19]. In our calculations we adopta parametrization previously used in refs. [20], [21]. A general energy-density functional for aninteracting finite system of neutrons (n) and protons (p) can be represented as E tot = E kin + E Skyrme + E Coulomb + E pair − E corr . (15)Here E kin is the kinetic energy of nucleons, calculated as E kin = R d r τ ( r ) , where τ = τ p + τ n is thesum of the kinetic energy densities of protons and neutrons. The particle densities are defined as ρ q = ∑ k ∈ q υ k | ψ k ( r ) | , (16)where we have introduced index q , which denotes proton (p) and neutrons (n). Here υ k are thevariational parameters in the wave function for even-even nuclei represented in the BCS mode as q ,see details in ref. [19]: | BCS i = ∏ k > ( u k + υ k ˆ a † k ˆ a †¯ k | i , (17)7here u k + υ k =
1. This equation follows from the normalization condition of the BCS state h BCS | BCS i = N . (18)The A-particle wave function | ψ i can be chosen as a superposition of all A-particle Slater determi-nants which can be written as | ψ i = ∑ i , i ,... i A c i , . . . , c i A · ˆ a † i · · · ˆ a † i A | i , (19)where the sets { i , i , . . . i A } represent the subspace of a complete basis of one-particle states.The kinetic energy densities can be calculated as τ q ( r ) = ∇ · ∇ ′ ρ q ( r , r ′ ) (cid:12)(cid:12) r = r ′ = ∑ k ∈ q υ k | ∇ψ k ( r ) | . (20)The contribution of the nuclear mean field is represented by the general Skyrme functional E Skyrme = Z d r (cid:26) b ρ − b ′ ∑ q ρ q + b ( ρτ − j ) − b ′ ∑ q ( ρ q τ q − j q ) − b ρ∆ρ + b ′ ∑ q ρ q ∆ρ q + b ρ α + − b ′ ρ α ∑ q ρ q (cid:27) + E LS , (21)where j q are the current densities defined as j q ( r ) = − i2 (cid:0) ∇ − ∇ ′ (cid:1) ρ q ( r , r ′ ) (cid:12)(cid:12) r = r ′ = − i2 ∑ k ∈ q v k { ψ † k ( r ) ∇ψ k ( r ) − [ ∇ψ † k ( r )] ψ k ( r ) } , (22)and the local single particle density is defined as ρ q ( r ) = ρ q ( r , r ) = ∑ k ∈ q υ k | ψ k ( r ) | . (23) E LS is the spin-orbit interaction term (see below). The nuclear Coulomb energy is calculated in thelocal density approximation including the exchange term, E Coulomb = e ZZ d r d r ′ ρ p ( r ) ρ p ( r ′ ) | r − r ′ | − e (cid:18) π (cid:19) / Z d r [ ρ p ( r )] / , (24)where the naked proton density is used, for details see ref. [22].8ollowing ref. [21], for the spin-orbit interaction we consider three possibilities: E stdLS = Z d r (cid:26) − b [ ρ∇ · J + s · ∇ × j + ∑ q ( ρ q ∇ · J q + s q · ∇ × j q )] (cid:27) , (25a) E ( J ) LS = E stdLS − Z d r (cid:26) ( t x + t x ) (cid:16) J − s · τ (cid:17) + ( t − t ) ∑ q ( J q − s q · τ q ) (cid:27) , (25b) E extLS = Z d r (cid:26) − b ( ρ∇ · J + s · ∇ × j ) − b ′ ∑ q ( ρ q ∇ · J q + s q · ∇ × j q ) (cid:27) . (25c)Finally E pair is the pairing energy density and the term E corr contains the center-of-mass correction.Higher order corrections to the exchange term do not play a role, as was shown in ref. [23].The variation of E Skyrme with respect to the density ρ q yields the self-consistent nuclear potential U q ( r ) = δ E δρ q ( r ) = b ρ ( r ) − b ′ ρ q ( r ) + b τ ( r ) − b ′ τ q ( r ) − b ∆ρ ( r ) + b ′ ∆ρ q ( r )+ b α + ρ α + ( r ) − b ′ ρ α q ( r ) − b ′ α ρ α − ( r ) ∑ q ′ ρ q ′ ( r ) − b ∇ · J ( r ) − b ′ ∇ · J q ( r ) . (26)For the protons one adds also the Coulomb term U c = e Z d r ′ ρ p ( r ′ ) | r − r ′ | − e (cid:18) π (cid:19) / ρ / p ( r ) . (27)The corresponding spin-orbit potential is calculated as W q ( r ) = δ E δ J q ( r ′ ) − ∇ δ E δ ( ∇ · J q ( r )) . (28)Depending on the choice made in eq.(20), one gets three possibilities: W stdLS = b · (cid:0) ∇ρ ( r ) + ∇ρ q ( r ) (cid:1) , (29a) W ( J ) LS = · ( t − t ) J q − · ( t x + t x ) · J + b · (cid:0) ∇ρ ( r ) + ∇ρ q ( r ) (cid:1) , (29b) W extLS = b ∇ρ ( r ) + b ′ ∇ρ q ( r ) . (29c)For E ext LS with b ′ = W rmf LS = ¯ h [ m − C eff ρ ( r )] C eff ∇ ( ρ ( r ) . (30)9he Skyrme parametrization using this choice is called SKI3 force. As shown in ref. [24], thisforce describes correctly the energy shifts in the lead nucleus. On the other hand, one can also vary b ′ while fitting. Hence one has an additional degree of freedom to describe the isoscalar and theisovector channels in the effective potential, like in the other terms of the Skyrme functional. In the SKI4 force this isospin dependence in the spin-orbit potential in was introduced, see ref. [24]. Theparameters b and b ′ have been adjusted to the spin-orbit splitting of O and the isotope shifts of Pb.The fitting has given b ′ ≈ − b .In the calculation of spin saturated systems, i.e. even-even -nuclei, the time odd currents are identicallyzero. In our approach we use a Density-Dependent-Delta-Interaction DDDI) for the usual pairingpotential υ pair ( r − r ′ ) = V (cid:20) − ( ρ ( r ) ρ ) γ · δ ( r − r ′ ) (cid:21) , (31)where γ = , ρ = .
16 fm − .Within our approach we have found that Skyrme forces like SkP or SLy7 do not predict the neutrondripping at all, i.e. the neutron chemical potential remains below m n c even at rather high densities.These Skyrme forces use E J LS from (25b) as a functional form of the spin-orbit interaction. The mostpopular Skyrme parametrizations including SkI3, SkI4, SkM*, SIII, SLy4, SLy6 and SLy230a neglectquadratic terms like J or J q and are consistent with the neutron drip-point. C. Implementation of electrons
In dense stellar matter a crucial role is played by electrons. The interaction between the nucleusand the electron background is very important and must be taken into account. In our Hartree-Fockcode this is done by using the method described in ref. [25]. The effects of inhomogeneity in theelectron distribution have been studied in ref. [26]. As follows from this analysis, the approximationof uniform density is good enough for the medium-size nuclei considered in this paper. But they maybecome significant for very big nuclei predicted in the inner crust, see ref. [26].The calculations below are done by using the WS method, i.e. the clusterized system is dividedinto WS cells, each cell containing a nucleus with charge number Z and Z electrons. Therefore the10equired electrical neutrality of the cell is automatically fulfilled. Using this in the calculation of theneutron star crust means a segmentation of the system into an ensemble of WS cells each containinga cluster of nuclear matter surrounded by an uniform background of electrons. In contrast to refs.[9, 10] we assume that the size of the WS cell depends only on the baryon density which is simulated.For a fixed baryon density n B and for each nuclear cluster containing Z protons and N neutrons theradius of the Wigner-Seitz cell R W S is calculated from the charge neutrality condition4 π R ws n e = Z , (32)where Z is the nuclear charge and n e = k F π is the constant background electron density. This gives k F = (cid:0) π Z (cid:1) · R W S ≈ . · Z / R W S , (33)where Z is the proton number in the nucleus. The value of k F calculated for every size of the WS cellat fixed baryon density n B determines the Coulomb interaction between the nucleus and the electronsand the energy of the relativistic electron gas. The radius of the WS cell depends on the mass numberof the nucleus A and the fixed baryon density n B , R W S = (cid:0) · A π · n B (cid:1) = ( π ) · n − B · A . (34)Inserting this in equation (33) gives a simple formula for the electron Fermi momentum, k F = ( π ) · n B · (cid:18) ZA (cid:19) . (35)Instead of baryon density n B , below we often use the nucleon mass density ρ B = m N × n B , where m N = ( m n + m p ) =
936 MeV / c .The interaction between the electrons and the nucleus in the WS cell is calculated in a self-consistentway, namely, by solving Poisson equation for the electrostatic potential created by both, protons andelectrons △ φ = − e · n ch ≡ − e · ( n p − n e ) . (36)To get a smooth charge distribution the electron density is parametrized with a smooth step-like func-tion n e ( r ) = n e + exp ( ( r − Rws ) a , where n e is the constant background electron density which ensures chargeneutrality in the cell. The diffuseness parameter a is taken to be 0.45 fm, as in ref. [25]. The electron11ensity is computed for each individual cell size depending on n B and ( A , Z ) . Finally, after the itera-tive procedure the electrostatic energy R ( ∇φ ) d r is added to the total energy of the cell.The calculational procedure is organized as follows:First, we choose the Skyrme force and start the calculation with a fixed baryon density dens in unitsof [ gcm ] . The loop over proton number runs from Z = Z =
130 and the loop over neutronsruns from N = Z − N = · Z , where we take only even Z and N. Then, the Wigner-Seitz cellparameters for every pair (Z, N) are calculated. These are the radius of the cell R W S , the Fermi mo-mentum k F , the electron chemical potential µ e and kinetic energy of the ultra-relativistic electron gas,which is E kine = · k F · Z . (37)Finally, the input files containing N , Z , k F and R W S are created and the Hartree-Fock program isexecuted with these input files. The numerical criterion of convergence of iterations to the optimumwithin the Hartree-Fock approach is chosen to be 10 − .For Z fixed, the nucleus (N,Z) with the minimal total energy per baryon E mintot / A is selected. Anexample of such calculations is presented in Fig. 1 for the SkM* force. It shows E tot / A for a setof nuclei (N, Z) with neutron numbers which minimize the total energy calculated at relatively lowdensities ( . − . ) · gcm . The nuclei with optimal proton number Z for several baryon densitieslie on the black line. These are ground states predicted at these baryon densities. IV. DESCRIPTION OF THE NUMERICAL CODE USED FOR NUCLEAR-STRUCTURE CALCU-LATIONS
We have performed HFB calculations of nuclear ground states with different Skyrme forces in a broadrange of baryon densities. We use the BCS model with delta-like pairing interaction for protons andneutrons. For the pairing strength of neutrons we choose -275,8 MeV and for the pairing strength ofprotons -291,7 MeV. The number of maximal iterations within the Hartree-Fock procedure is chosento be 1000 and the relative shift of energy level dE at the end of the iterations has to be less than 10 − ,as mentioned before. If dE < − is reached, the calculation is considered to be converged. The HFBcalculations are done on a grid with the maximal radius equal to the radius of the Wigner-Seitz cell R W S . The spacing of the grid points is chosen to be 0,1 fm.As a first step the grid functions defining the size and the fineness of the grid and the corresponding12rial wave functions on the grid are constructed. Grid size and grid spacing are taken from the inputdata. Then the choice of the representation of the baryons and their interaction are done. The nu-clear potentials are constructed by using the Skyrme parameters stored in a special file containing allparameter sets. The Skyrme forces used in the code are characterized by the following parameters: (cid:26) t , t , t , t , x , x , x , x , b , b ′ , α , ¯ h m , so_Curr, coul_Ex, cm_Corr, dens_Dep (cid:27) . (38)Using the baryon data from the input file the wave functions are calculated again with the help of theseconstructed potentials. With the new wave functions the potentials are calculated self-consistentlyagain.At the end of every iteration loop total energy of the nucleus is calculated as mentioned before. Thetotal kinetic energy is calculated by summing up the kinetic energy and adding the pairing energy foreach nucleon represented in a many-body wave function. Then the contributions of the Skyrme meanfields are added. The single-particle wave-functions are weighted with their occupation probabilitiesand degeneracy factors. The center of mass correction is done at the very last iteration step andnormalized with 1 / A .The Hartree-Fock equations during are solved iteratively using the damped gradient step method, seedetails in refs. [21, 27]. The true wave function is found from the iterative equation | φ n + α ) = O n | φ ( n ) α ) − ˆ D (cid:2) ˆ h − ( φ α | ˆ h | φ α ) (cid:3) | φ ( n ) α ) o , (39)where ˆ h is the Hamiltonian and only the projections on the diagonal elements using Π diag α = −| φ α ih φ α | are considered.The damping operator ˆ D is defined as ˆ D = x + ˆ tE (40)with the step parameter x , the kinetic energy ˆ t and the constant energy scale E which controls thedamping. When the iterative process does not conserve the orthonormality of the states , a new set ofwave functions has to be generated. The symbol O { ... } stands for this orthonormalization procedureafter the gradient step.The HFB code was initially developed in Fortran77 and was ported later to
C++ by Bender andRutz [21, 28]. Then the code was modified as described in ref. [25], where the uniform electronbackground was implemented. The code was further developed in ref. [26], where the neutron gasbuild by the dripping neutrons was simulated with a special choice of boundary conditions.13 . PREDICTED DRIP POINT DENSITIES AND DRIP POINT ELEMENTS
The drip points are defined as the configurations with a global minimum of the total energy E tot / A under the condition µ n = m n · c . It turned out that the results depend essentially on the Skyrme forceused in the calculation. For example, using the SKM* parametrization we get the drip point at abaryon density of 4 . · gcm with Ti as the drip nucleus. Figure 2 shows the total energy perparticle E tot / A as a function of proton number calculated for this density. We have also explored thechemical composition in a region slightly above the neutron drip region. Here free neutrons appear inthe cell. Figure 3 shows the evolution of the nuclear composition with increasing baryon density. Weobserve that the degree of neutronization increases, i.e. x = Z / A decreases, with increasing baryondensity. As mentioned before, we have to exclude some Skyrme forces, like SLy7 or SKP, whichdo not provide the dripping of neutrons at all. The total neutron number in the cell is determinedfrom the minimization of the total energy E tot / A . By performing these calculations one obtainsthe chemical composition of the neutron drip region as a function of baryon density. Due to theappearance of the neutron gas after the neutron drip in the WS cell one can expect a change of theslope in the EoS of cold nuclear matter, because the lighter neutrons generate an additional pressure(for details see ref. [29]).In the neutron-drip region the strong forces are attractive until the density 1 , × gcm is reached,where the transition to the outer core is expected. Hence the thickness of the inner crust dependsstrongly on the neutron-drip density. Below we demonstrate precisely how drip-point elements anddrip-point densities vary for different Skyrme forces. Comparing the drip-element with the pre-dripelement we observe that in every calculation the elements remain the same but the number of drippedneutrons varies from 2 to 12 depending on the Skyrme force used.Finally we present the results of our extended calculations carried out for a large variety of Skyrmeforces used in the literature. The most popular 7 parametrizations are discussed in some details.Predictions for δ , x = ZA and ρ ND of other more than 230 parametrizations are presented in theappendix. The Fermi momentum of the electron gas k F at the drip point where the slope of the EoSchanges can be calculated using equation (35). With this method only x = Z / A can be obtained, thespecial drip element can not be determined. All the values for x with parameters taken from [14] arein the same order: x ∼ .
3. 14n this section we present our calculations for the 7 most popular Skyrme forces: SkM*, SkI3,SkI4,SIII, SLy4, SLy6 and the newer SLy230a.
A. Skyrme force SKM ∗ The previously used Skyrme force SkM has been extensively studied for both spherical and deformednuclei through Hartree-Fock plus BCS calculations [30]. Ground-state radii and multipole momentsare found in excellent agreement with experimental data. But nevertheless binding energies weresystematically too high and fission barriers were significantly too low. The modified Skyrme forceSkM ∗ was the first Skyrme force with reasonable incompressibility as well as better fission properties[31]. Neutron and proton single-particle energies in the Zr , Pb and Pu were computed andcompared with experimental energies taken from the compilation of ref. [32]. With SkM* it waspossible to calculate the fission barrier of Pu rather accurately. This Skyrme force uses E s tdLS forthe spin-orbit interaction. Calculations with SKM* yield Ti as the pre-drip element at a density of4 , · gcm . The neutron drip occurs at a density of 4 , · gcm starting with the Titanium isotope Ti , which corresponds to x = . , · gcm obtained in [17]. The evolution of the ground state in the drip region is illustratedin Fig. 3. One can see that during the neutron drip the element remains the same, but with increasingdensity the neutronization increases, i.e. x = Z / A decreases. With the estimate based on the symmetryenergy described before the drip occurs at a density of 4 , · gcm with x = .
383 (see chapter II).Skyrme force ρ pre − drip [ gcm ] pre-drip nuclei drip element ρ drip [ gcm ] E tot A [MeV] x = ZA µ n [MeV]SkM ∗ , ·
11 8022 Ti Ti , · -1,3540 0.268 2 , · − . Skyrme forces SkI3 and SkI4 The SkIx Skyrme forces are based on calculations performed by P.-G. Reinhard and H. Flocard in1995 [24]. They used least-square fit of nuclear ground-state properties [33] and took experimentaldata of exotic nuclei into account. The Skyrme force SkI1 has the standard spin-orbit coupling E s tdLS and b ′ = b . This Skyrme force does not provide the dripping of neutrons and therefore is notconsidered here. The SkI3 force has a generalized spin-orbit coupling which is the non-relativisticlimit of a relativistic mean-field model ( E e xtLS and b ′ = E e xtLS . The pre-drip element calculated with SkI3 is Se witha magic N=82 shell at a density of 3 , · gcm . The drip element is the Selenium isotope Se at adensity of 3 , · gcm with x = . Ni with the magic proton shell P =
28 at a slightly higher density of 3 . · gcm . The drip elementis the Nickel isotope Ni with x = .
264 at a density of 3 , · gcm . The neutron gas at the dripconsists of 8 neutrons in a spherical cell with radius 46,39 fm. In contrast the the SkM* force theneutron drip densities of SkI3 ans SkI4 are both lower than predicted in [7]. With the estimates basedon the symmetry energy the drip occurs at a density of 4 , · gcm with x = .
396 using SKI3and 4 , · gcm with x = .
379 using SKI4.Skyrme force ρ pre − drip [ gcm ] pre-drip nuclei drip element ρ drip [ gcm ] E tot A [MeV] x = ZA µ n [MeV]SkI3 3 , ·
11 11634 Se Se , · -1,6701 0.288 2 , · − SkI4 3 , ·
11 9828 Ni Ni , · -1,538 0.264 1 , · − C. Skyrme force SIII
The Skyrme force SIII was proposed by Beiner et al. in 1975 [34]. It is one of the oldest Skyrmeforces still used today, which uses E s tdLS for the spin-orbit coupling. The authors have performeda detailed study of the influence of the force parameters on the binding energies, charge densities,radii and single-particle energies. They explored also magic nuclei. For the SIII force the pre-drip16lement is Kr at ρ = , · gcm with magic N=82 shell . The neutron drip starts at a densityof ρ = , · gcm with the Krypton isotope Kr and x = . , · gcm with x = . ρ pre − drip [ gcm ] pre-drip nuclei drip element ρ drip [ gcm ] E tot A [MeV] x = ZA µ n [MeV]SIII 3 , ·
11 11836 Kr Kr , · -1,5468 0.277 2 , · − D. Skyrme forces SLy4 and SLy6
The construction of the SLyx Skyrme forces was motivated by the most accurate description ofneutron-rich nuclei. The SLy6 force does not use the J -term. The correction due to the centerof mass motion is introduced in a self-consistent way using the full microscopic ansatz (details see inref. [35]). The neutron drip starts at a density of ρ = , · gcm with nucleus Kr and x = . Kr with a magic N=82 shell appears at a density of ρ = , · gcm . Hencethe neutron gas at the drip consists of 6 unbound neutrons in the spherical cell with a radius of 49 , h ˆ P cm i = ∑ k , k ′ ≥ ∑ m , m ′ ≥ p k , k ′ p m , m ′ h BCS | ˆ a † k ˆ a k ′ ˆ a † m ˆ a m ′ | BCS i (41)are taken into account, for details see ref. [21]. The pre-drip element is Sr with a magic N=82shell . The neutron gas at the drip point consists of 4 neutrons in a spherical shell of radius 50 ,
73 fm.The neutron drip starts at a density of ρ = , · gcm with drip element Sr and x = . K ∞ ,sum-rule enhancement factor κ or the asymmetry coefficient a sym are identical. There are only somedifferences in (a) the energy per nucleon (15,97 MeV for SLy4 and 15,90 MeV for SLy6), (b) thenuclear density n (0.160 fm − for SLy4 and 0,159 fm − for SLy6), (c) the effective mass m ∗ / m (0,695 for SLy4 and 0,690 for SLy6) and (d) the use of the spin-orbit coupling: Sly4 uses E ext LS ,whereas Sly6 uses E std LS . We found that SLy6 yields a slightly higher neutron drip density with Kr
17s drip-element. But actually the two drip elements of SLy4 and SLy6 are isobars. One more force ofthis family, SLy7, is using E ( J ) LS for the spin-orbit interaction which we do not include in our analysis.With the estimates based on the symmetry energy the drip occurs at a density of 4 , · gcm with x = .
388 using SLy4 and 4 , · gcm with the same value x = .
388 using SLy6.Skyrme force ρ pre − drip [ gcm ] pre-drip nuclei drip element ρ drip [ gcm ] E tot A [MeV] x = ZA µ n [MeV]SLy4 3 , ·
11 12038 Sr Sr , · -1,6490 0.306 2 , · − SLy6 3 , ·
11 11836 Kr Kr , · -1,521258 0.29 3 , · − E. Skyrme force SLy230a
The force SLy230a [35] is using the following parameters for symmetric nuclear matter:E/A [Mev] ρ fm K[MeV] E sym m ∗ / m -15.988 0.160 229.87 31.97 0.697This parametrization is quite common and matches very well to the binding energies and thecharge radii of doubly-magic nuclei. Pre-drip element at a density of ρ = , · gcm is Mo with a magic N=82 shell . In ref. [8] the authors found Molybdenum isotope to be the last pre-dripnucleus. Actually they found the sequence of pre-drip nuclei to be Mo -
Zr -
Sr -
Kr. Theneutron gas at the drip point with drip element Mo and x = .
333 is very diluted with 2 neutronsin a spherical cell of radius 53,784 fm. With the estimates based on the symmetry energy the dripoccurs at a density of 4 , · gcm with x = .
388 as SLy4 or SLy6.Skyrme force ρ pre − drip [ gcm ] pre-drip nuclei drip element ρ drip [ gcm ] E tot A [MeV] x = ZA µ n [MeV]SLy230a 2 , ·
11 12442 Mo Mo , · -1,6788 0.333 2 , · − The most important characteristics of the drip-point for the 7 selected Skyrme parametrizations aresummarized in the table below. One can see that the number of dripped neutrons in the cell varies18rom 2 to 12 and the neutron-gas has densities around n n ∼ − fm − , but the chemical potentialsof the dripped neutrons are almost of the same value: µ n ∼ ( , − , ) MeV. Although the protonnumbers of the drip-point elements vary from 22 to 42 the chemical potential of the electrons µ e hasalso nearly the same value for these Skyrme forces, around
25 MeV .Therefore the lower drip density is the most important point of our investigation. Obviously, fora lower drip point density the outer crust would be thinner. This may significantly change theobservable signatures of neutron stars as pointed out e.g. in [6].
Drip-point characteristics for 7 most popular Skyrme forces
Skyrme drip element n ND [ gcm ] n-gas h r i bar[fm] h r i W S [fm] n n [fm − ] µ n − m n c [MeV] µ e [MeV]SkI3 Se , · , · − , · − Ni , · , · − , · − ∗ Ti , · , · − , · − Kr , ·
12 4,745 50,10 3 , · − , · − Sr , · , · − , · − Kr , · , · − , · − Mo , · , · − , · − VI. DISCUSSION AND CONCLUSIONS
The presented calculations demonstrate that almost all obtained drip densities are lower than thecommonly accepted value of Baym, Bethe and Pethick, ref. [17]. Only the SKM ∗ -force matches quitewell this prediction. But on the other hand, the Wigner-Seitz cell is mainly filled by the relativisticelectron gas, only about of the volume of the WS-cell is occupied by the nucleus. Therefore,practically the whole pressure comes from the electrons. In our calculations the interaction of thenucleus with the surrounding electron gas, the so-called lattice energy, is taken into account in a newself-consistent way. At every baryon density we have calculated every single nucleus and obtainedprecise values of the neutron chemical potential and the neutron gas density in the WS cell. Theseobservables are calculated within a fully microscopic approach and directly linked to the parametersof the special Skyrme force. Within our approach we are able to describe exactly the nuclear changes19ith increasing baryon density. So, as a main result of our investigations, most of our drip densitiesare lower than the commonly accepted values predicted between 4 and 5 × gcm [36]. The firstcalculations of Baym, Pethik and Sutherland [7] predicted the neutron drip density 4 , · gcm .They found the electron chemical potential to be 26 . Kr as the drip element. In their calculation Kr remains thefavoured element even somewhat beyond that drip-point. Using a compressible liquid-drop modelBaym, Bethe and Pethik [17] confirmed the value 4 , · gcm as the density where the neutron-dripstarts. But even the drip densities based on the estimates using the symmetry energy of ref. [14] dosuggest lower densities of the neutron drip for most forces used (see appendix).However, at such high densities more accurate calculations of electron-nucleus interaction arerequired. Such calculations have been carried out in ref. [26] where the electron distributions wereobtained by solving the Poisson equation for electrostatic potential in combination with nuclearstructure calculations using a RMF model.Magic numbers play an important role in nuclear structures and stability. The pre-drip elements,most of them stabilized by magic neutron shells, are the last stable neutron-rich elements beforeneutrons drip out. Hence one can expect that magic neutron shells may stabilize the configurationseven at higher densities. Negele and Vautherin found Kr to be the last pre-drip element, stabilizedby the magic N=82 shell [8]. This magic shell appears also in the Skyrme parametrizations SKI3,SIII, SLy4, SLy6 and SLy230a. On the other hand, the Skyrme force SKM* belongs to the groupof parametrizations favouring the magic N=50 shell in the drip region. This shell survives untilthe density ρ = , · gcm . At higher densities the number of neutrons increases up to 60 forthe drip element Ti at ρ = , · gcm . According to ref. [37], there are mainly two types ofparametrizations, favouring either the N=50 shell or the
N=82 shell in the neutron-drip region. Itis also shown in ref. [37], that the nuclear symmetry energy has a great influence on the chemicalcomposition of the outer crust and therefore on the magic neutron shells in the drip-region. In Fig.2 one can observe that there is a local minimum at proton number Z =
50. But this minimum getsweaker with increasing density and is shifted to a proton number smaller than 50. Hence in ourcalculations the magic P=50 shell is quenched at the edge of the outer crust.Looking at the Skyrme forces used in our analysis, one can see that the values of E / A of the dripelements vary between − . − . E / A value is slightly higher ( − .
35 MeV). This can be explained by the smaller neutron numberof the magic shell appearing in the drip region. Accordingly, the drip element Ti obtained with20kM* has the smallest proton number ( Z =
22) of all Skyrme forces within our investigation. Buton the other hand, the atomic numbers of the drip elements obtained for other forces vary from 28(Ni) for SkI4 up to 42 (Mo) for SLy230a. The element
Mo appears also in the calculations ofNegele-Vautherin, ref. [8], but as a pre-drip element. Looking at the mass-energy density profiles ofneutron stars in ref. [5], one may conclude that the neutron drip-density is about ρ ∼ · gcm .This is a slightly lower than the value reported in ref. [17]. Using our method most drip pointdensities are lower than 4 · gcm . Applying the method described in [13] to all Skyrme forceslisted in [14] we have obtained the table presented in Appendix A which includes more than 200Skyrme forces. As one can see, drip densities vary from 2 , · gcm to 4 , · gcm . Obviously,a lower drip density corresponds to a thinner outer crust. These results may be useful for calculatingcooling rates and mechanical properties of the outer crust.Finally, we conclude with a few remarks:First, one can see that due to ρ ND the modern Skyrme parametrization SKM ∗ ( ρ ND = , · gcm )is indeed better parametrization than the older SKM force ( ρ ND = , · gcm ), supposing that theneutron drip density predicted by [17] is true.Second, the drip densities obtained with SkTK, SIII ∗ and the family ZRXx are too high andthird, the drip densities obtained with the family v070 - v110 are too low.Third, the Skyrme forces ZR3a, ZR3b and ZR3c with S < ppendix A: Neutron drip-point characteristics for different Skyrme forces Skyrme force E S δ x ρ ND [ gcm ] BSK1 -15.81 27.81 0.2524 0.3738 4 . · BSK2 -15.80 28.00 0.2507 0.3746 4 . · BSK2’ -15.79 28.00 0.2506 0.3747 4 . · BSK3 -15.81 27.93 0.2514 0.3743 4 . · BSK4 -15.77 28.00 0.2503 0.3748 4 . · BSK5 -15.80 28.70 0.2452 0.3774 4 . · BSK6 -15.75 28.00 0.25 0.375 4 . · BSK7 -15.76 28.00 0.2501 0.3749 4 . · BSK8 -15.83 28.00 0.2511 0.3744 4 . · BSK9 -15.92 30.00 0.2372 0.3814 4 . · BSK10 -15.91 30.00 0.2371 0.3814 4 . · BSK11 -15.86 30.00 0.2346 0.3818 4 . · BSK12 -15.86 30.00 0.2364 0.3818 4 . · BSK13 -15.86 30.00 0.2364 0.3818 4 . · BSK14 -15.85 30.00 0.2363 0.3818 4 . · BSK15 -16.04 30.00 0.2388 0.3806 4 . · BSK16 -16.05 30.00 0.2389 0.3805 4 . · BSK17 -16.06 30.00 0.2390 0.3805 4 . · BSK18 -16.06 30.00 0.2390 0.3805 4 . · BSK19 -16.08 30.00 0.2393 0.3803 4 . · BSK20 -16.08 30.00 0.2393 0.3803 4 . · BSK21 -16.05 30.00 0.2389 0.3805 4 . · E -16.13 27.66 0.2582 0.3709 4 . · Es -16.02 27.44 0.2599 0.37 4 . · f − -16.02 32.00 0.225 0.38750 4 . · f + -16.04 32.00 0.2252 0.3874 4 . · f -16.03 32.00 0.2251 0.3874 4 . · FPLyon -15.92 30.93 0.2307 0.3846 4 . · Gs -15.59 31.13 0.2307 0.3846 4 . · GS1 -16.03 28.86 0.2471 0.3764 4 . · E S δ x ρ ND [ gcm ] GS2 -16.01 25.96 0.2715 0.3642 4 . · GS3 -16.00 21.49 0.3208 0.3396 4 . · GS4 -15.96 12.83 0.4980 0.251 4 . · GS5 -15.91 18.70 0.3604 0.3198 4 . · GS6 -16.04 14.33 0.4558 0.2721 4 . · GSkI -16.02 32.03 0.2248 0.3876 4 . · GSkII -16.12 30.49 0.2364 0.3818 4 . · KDE -15.99 31.97 0.2248 0.3876 4 . · KDE0v -16.10 32.98 0.2199 0.39 4 . · KDE0v1 -16.23 34.58 0.2122 0.3939 4 . · LNS -15.32 33.43 0.2076 0.3962 3 . · MSk1 -15.83 30.00 0.2360 0.382 4 . · MSk2 -15.83 30.00 0.2360 0.382 4 . · MSk3 -15.79 28.00 0.2506 0.3747 4 . · MSk4 -15.79 28.00 0.2506 0.3747 4 . · MSk5 -15.79 28.00 0.2506 0.3747 4 . · MSk5 ∗ -15.78 28.00 0.2504 0.3748 4 . · MSk6 -15.79 28.00 0.2506 0.3747 4 . · MSk7 -15.80 27.95 0.2511 0.3744 4 . · MSk8 -15.80 27.93 0.2513 0.3743 4 . · MSk9 -15.80 28.00 0.002507 0.3746 4 . · MSkA -15.99 30.35 0.2356 0.3822 4 . · MSL0 -16.00 30.00 0.002383 0.03808 4 . · NRAPR -15.85 32.78 0.2180 0.391 4 . · PRC45 -15.82 51.01 0.1446 0.4277 4 . · RATP -16.05 29.26 0.2444 0.03778 4 . · Rs -15.59 30.82 0.2271 0.3864 4 . · Sefm068 -15.92 88.57 0.0862 0.4569 4 . · Sefm074 -15.81 33.40 0.2138 0.3931 4 . · Sefm081 -15.69 30.76 0.2288 0.3856 4 . · E S δ x ρ ND [ gcm ] Sefm09 -15.55 27.78 0.2489 0.3755 4 . · Sefm1 -15.40 24.81 0.2731 0.3634 3 . · SGI -15.89 28.33 0.002493 0.3753 4 . · SGII -15.60 26.83 0.2575 0.37125 4 . · SGOI -16.63 45.20 0.1696 0.4152 5 . · SGOII -16.70 93.98 0.0852 0.4574 5 . · SI -15.99 29.24 0.2437 0.3781 4 . · SII -15.99 34.16 0.2116 0.3942 4 . · SIII -15.85 28.16 0.2501 0.3749 4 . · SIII ∗ -16.07 31.97 0.2583 0.3708 . · SIV -15.96 31.22 0.2293 0.3853 4 . · Sk1’ -15.99 29.35 0.2429 0.3785 4 . · SK255 -16.33 37.40 0.1986 0.4007 4 . · SK272 -16.28 37.40 0.1980 0.401 4 . · SkA -15.99 32.91 0.19 0.3905 4 . · SkA25s20 -16.07 33.78 0.2148 0.3926 4 . · SkA35s15 -16.01 30.56 0.446 0.3827 4 . · SkA35s20 -16.08 33.57 0.2164 0.3918 4 . · SkA35s25 -16.14 36.98 0.1985 0.4007 4 . · SkA45s20 -16.08 33.39 0.2172 0.3914 4 . · SkB -15.99 23.88 0.2921 0.3539 4 . · SkI1 -15.95 37.53 0.1937 0.4031 4 . · SkI2 -15.78 33.37 0.2136 0.3922 4 . · SkI3 -15.98 34.83 0.2078 0.3961 4 . · SkI4 -15.95 29.50 0.2412 0.3794 4 . · SkI5 -15.85 36.64 0.1969 0.4015 4 . · SkI6 -15.89 29.90 0.2375 0.3812 4 . · SkM -15.77 30.75 0.2230 0.385 . · SkM ∗ -15.77 30.03 0.2349 0.38255 4 . · SkM1 -15.77 25.17 0.2753 0.3623 4 . · SkMP -15.56 29.89 0.2331 0.3834 4 . · E S δ x ρ ND [ gcm ] SkO -15.84 31.97 0.2289 0.3855 4 . · SkO’ -15.75 31.95 0.2219 0.3890 4 . · SkP -15.95 30.00 0.2376 0.3812 4 . · SkRA -15.78 31.32 0.2631 0.3868 . · SkS1 -15.86 28.75 0.2456 0.3772 4 . · SkS2 -15.89 29.23 0.2424 0.3788 4 . · SkS3 -15.88 28.84 0.2452 0.3774 4 . · SkS4 -15.88 28.35 0.2490 0.3755 4 . · SkSC1 -15.85 28.10 0.2506 0.3747 4 . · SkSC2 -15.90 24.74 0.2817 0.3591 4 . · SkSC3 -15.85 27.01 0.2597 0.3715 4 . · SkSC4 -15.87 28.80 0.2454 0.3773 4 . · SkSC4o -15.87 27.00 0.26 0.37 4 . · SkSC5 -15.85 30.99 0.2294 0.3853 4 . · SkSC6 -15.92 24.57 0.2837 0.2581 4 . · SkSC10 -15.96 22.83 0.3034 0.3483 4 . · SkSC11 -15.87 28.80 0.2454 0.3773 4 . · SkSC14 -15.92 30.00 0.2372 0.3814 4 . · SkSC15 -15.88 28.00 0.2518 0.3741 4 . · SkSP.1 -15.90 28.00 0.2521 0.3739 4 . · SkT -15.40 33.66 0.2073 0.3963 4 . · SkT1 -15.98 32.02 0.2244 0.3878 4 . · SkT2 -15.94 32.00 0.2240 0.388 4 . · SkT3 -15.95 31.50 0.2373 0.3863 4 . · SkT4 -15.96 35.24 0.2054 0.3973 4 . · SkT5 -16.00 37.00 0.68 0.4016 4 . · SkT6 -15.96 29.97 0.2379 0.3811 4 . · SkT7 -15.94 29.52 0.2409 0.3795 4 . · SkT8 -15.94 29.92 0.2380 0.381 4 . · SkT9 -15.88 29.76 0.2384 0.3808 4 . · E S δ x ρ ND [ gcm ] SkT1 ∗ -16.20 32.31 0.2253 0.3873 4 . · SkT3 ∗ -16.20 31.97 0.275 0.3862 4 . · SkT1a -15.98 32.02 0.2244 0.3878 4 . · SkT2a -15.94 32.00 0.2240 0.388 4 . · SkT3a -15.95 31.50 0.2273 0.38630 4 . · SkT4a -15.96 35.45 0.2042 0.3979 4 . · SkT5a -16.00 37.00 0.1968 0.4016 4 . · SkT6a -15.96 29.97 0.2379 0.3910 4 . · SkT7a -15.94 29.52 0.2409 0.3795 4 . · SkT8a -15.94 29.92 0.2381 0.381 4 . · SkT9a -15.88 29.76 0.2384 0.3808 4 . · SkTK -16.70 35.57 0.2122 0.3939 . · SKX -16.05 31.10 0.2313 0.3843 4 . · SKXce -15.86 30.15 0.2353 0.3823 4 . · SKXm -16.04 31.20 0.2305 0.3847 4 . · SKxs15 -15.76 31.88 0.2244 0.3878 4 . · SKxs20 -15.81 35.50 0.2022 0.3989 4 . · SKxs25 -15.87 39.60 0.1835 0.4082 4 . · SKz-1 -16.01 32.00 0.2249 0.3875 4 . · SKz0 -16.01 32.00 0.2249 0.3875 4 . · SKz1 -16.01 32.01 0.2248 0.3876 4 . · SKz2 -16.01 32.01 0.2248 0.3876 4 . · SKz3 -16.01 32.01 0.2248 0.3876 4 . · SKz4 -16.01 32.01 0.2248 0.3876 4 . · SLy0 -15.97 31.98 0.2245 0.3877 4 . · SLy1 -15.99 31.99 0.2247 0.3876 4 . · SLy2 -15.99 32.00 0.2246 0.3877 4 . · SLy230a -15.99 31.99 0.2247 0.3876 4 . · SLy230b -15.97 32.01 0.2243 0.3878 4 . · SLy3 -15.94 31.97 0.2242 0.3879 4 . · E S δ x ρ ND [ gcm ] SLy4 -15.97 32.00 0.2244 0.3878 4 . · SLy5 -15.99 32.01 0.2245 0.3877 4 . · SLy6 -15.92 31.96 0.2240 0.388 4 . · SLy7 -15.90 31.99 0.2235 0.3882 4 . · SLy8 -15.97 32.00 0.2244 0.3878 4 . · SLy9 -15.80 31.98 0.2223 0.3888 4 . · SLy10 -15.90 31.90 0.2241 0.3879 4 . · SQMC1 -14.00 29.68 0.2131 0.3934 . · SQMC2 -14.29 28.70 0.2239 0.388 . · SQMC3 -15.98 . · SQMC600 -15.74 34.38 0.2074 0.3963 4 . · SQMC650 -15.57 33.65 0.2094 0.3953 4 . · SQMC700 -15.49 33.47 0.2095 0.3952 4 . · SQMC750 -15.60 33.75 0.2092 0.3954 4 . · SSK -16.16 33.50 0.2175 0.3912 4 . · SV -16.05 32.82 0.2202 0.3899 4 . · SV-bas -15.91 30.00 0.2371 0.3814 4 . · SV-min -15.91 30.66 0.2324 0.3838 4 . · SVI -15.76 26.88 0.2595 0.3702 4 . · SVII -15.79 26.96 0.2592 0.3704 4 . · SV-K218 -15.90 30.00 0.2369 0.3815 4 . · SV-K226 -15.90 30.00 0.2369 0.3815 4 . · SV-K241 -15.91 30.00 0.2371 0.3814 4 . · SV-kap00 -15.90 30.00 0.2369 0.3815 4 . · SV-kap02 -15.90 30.00 0.2369 0.3815 4 . · SV-kap06 -15.91 30.00 0.2371 0.3814 4 . · SV-mas07 -15.89 30.00 0.2368 0.3816 4 . · SV-mas08 -15.90 30.00 0.2369 0.3815 4 . · SV-mas10 -15.91 30.00 0.2371 0.3814 4 . · SV-sym28 -15.47 28.47 0.2564 0.3718 4 . · E S δ x ρ ND [ gcm ] SV-sym32 -15.94 32.00 0.2240 0.388 4 . · SV-sym34 -15.97 34.00 0.2123 0.3938 4 . · SV-tls -15.89 30.00 0.2367 0.3816 4 . · T -15.93 28.35 0.2498 0.3751 4 . · T11 -16.01 32.00 0.2249 0.3975 4 . · T12 -16.00 32.00 0.2247 0.3876 4 . · T13 -16.00 32.00 0.2247 0.3876 4 . · T14 -15.99 32.00 0.2246 0.3877 4 . · T15 -16.01 32.00 0.2249 0.3875 4 . · T16 -16.01 32.00 0.2249 0.3875 4 . · T21 -16.03 32.00 0.2251 0.3874 4 . · T22 -16.02 32.00 0.225 0.3875 4 . · T23 -16.01 32.00 0.2249 0.3875 4 . · T24 -16.01 32.00 0.2249 0.3875 4 . · T25 -15.99 32.00 0.2246 0.3877 4 . · T26 -15.98 32.00 0.2249 0.3875 4 . · T31 -16.02 32.00 0.225 0.3875 4 . · T32 -16.03 32.00 0.2251 0.3874 4 . · T33 -16.02 32.00 0.225 0.3875 4 . · T34 -16.02 32.00 0.225 0.3875 4 . · T35 -16.00 32.00 0.2247 0.3876 4 . · T36 -15.99 32.00 0.2246 0.3877 4 . · T41 -16.06 32.00 0.2255 0.3872 4 . · T42 -16.05 32.00 0.2254 0.3873 4 . · T43 -16.04 32.00 0.2252 0.3874 4 . · T44 -16.02 32.00 0.225 0.3875 4 . · T45 -16.02 32.00 0.225 0.3875 4 . · T46 -16.00 32.00 0.2247 0.3876 4 . · T51 -16.05 32.00 0.2254 0.3873 4 . · E S δ x ρ ND [ gcm ] T52 -16.06 32.00 0.2255 0.3872 4 . · T53 -16.02 32.00 0.225 0.3875 4 . · T54 -16.03 32.00 0.2251 0.3874 4 . · T55 -16.03 32.00 0.2251 0.3874 4 . · T56 -16.01 32.00 0.2249 0.3875 4 . · T61 -16.07 32.00 0.2256 0.3872 4 . · T62 -16.07 32.00 0.2256 0.3872 4 . · T63 -16.06 32.00 0.2255 0.3872 4 . · T64 -16.03 32.00 0.2251 0.3874 4 . · T65 -16.04 32.00 0.2252 0.3874 4 . · T66 -16.02 32.00 0.225 0.3875 4 . · v070 -15.78 27.98 0.2506 0.3747 . · v075 -15.80 28.00 0.2507 0.3746 . · v080 -15.79 28.00 0.2506 0.3747 . · v090 -15.79 28.00 0.2506 0.3747 . · v100 -15.79 28.00 0.2506 0.3747 . · v105 -15.79 28.00 0.2506 0.3747 . · v110 -15.79 28.00 0.2506 0.3747 . · Z -15.97 26.82 0.2604 0.3698 4 . · ZR1a -16.99 9.84 0.65 0.75 . · ZR1b -16.99 18.50 0.385 0.3075 . · ZR1c -16.99 31.50 0.2407 0.3796 . · ZR2b -16.99 11.95 0.5562 0.2219 . · ZR2c -16.99 27.43 0.2725 0.3637 . · ZR3a -16.99 -138.96 -0.0631 0.53155 . · ZR3b -16.99 -100.46 -0.08847 0.5442 . · ZR3c -16.99 -42.71 -0.2239 0.61195 . · Zs -15.88 26.69 0.2629 0.3682 4 . · Zs ∗ -15.96 28.80 0.2466 0.3767 4 . · ppendix B: E tot / A over proton number for densities ( . − . ) · gcm calculated with SkM* Figure 1:
The binding energies per nucleon E/A of Fe isotopes calculated for baryon densities ( . − . ) · gcm with the SkM* force. As one can see, F e forms the ground state of the outer crust, atthese densities. ppendix C: E tot / A over proton number for densities ( . − . ) · gcm directly before the drip point(SkM*) Figure 2:
The binding energy per nucleon as a function of proton number in the drip-region of the outer crust,calculated for the SkM* force. The neutron number of every nucleus shown on the plot minimizes E/A of theisotopic chain. ppendix D: Shown are ground state elements over baryon density at the drip point (SkM*) Figure 3:
Neutron drip region calculated with SkM*. Shown are ground states elements over baryon densityat the drip point. The drip element T i consists of the nuclei T i and 2 free neutrons forming a gas. Withincreasing density the nucleus in the WS-cell gets more and more neutron-rich (decreasing x nuc , where onlythe bound neutrons are taken into account).
1] M. Camenzind,
Compact Objects in Astrophysics , Springer Verlag, Heidelberg (2007)[2] David Tsang, Jacelyn S. Read, Tanja Hinderer, Anthony L. Piro and Ruxandra Bondarescu, Phys. Rev.Lett.
2D Cooling of magnetized neutron stars , Astronomy andAstrophysics 486, 255-271 (2008)[4] W. Becker,
Neutron Stars and Pulsars , Springer, Berlin Heidelberg (2009)[5] N. K. Glendenning,
Compact Stars, Nuclear Physics, Particle Physics and general Relativity , Secondedition, Springer Press, Heidelberg (2000)[6] N. K. Glendenning,
Special and General Relativity With Applications to White Dwarfs, Neutron Starsand Black Holes , Springer Press, Heidelberg (2007)[7] G. Baym, C. Pethik and P. Sutherland,
The ground state of matter at high densities: Equation of stateand stellar models , The Astrophysical Journal 170, 299-317 (1971)[8] J.W. Negele and D. Vautherin,
Neutron star matter at sub-nuclear densities , Nucl.Phys. (A 207 298-320(1973)[9] E.P. Wigner, F. Seitz,
On the constitution of metallic sodium , Phys. Review 43, (1933)[10] E.P. Wigner, F. Seitz,
On the constitution of metallic sodium , Phys. Review 46, (1934)[11] C. J. Pethik and D. G. Ravenhall,
An Introduction to Matter at Sub-nuclear Densities , IN: J. Venturaand D. Pines,
Neutron Stars: Theory and Observation , Series C: Mathematical and Physical Sciences334, Springer Science and Business Media, B. V. (1991)[12] K. Bethge,
Kernphysik, Eine Einführung , Springer-Verlag, Berlin Heidelberg (1996)[13] P. Haensel, A.Y. Potekhin and D.G. Yakovlev,
Neutron Stars 1, Equation of State and Structure , Astro-physics and Science Library, Volume 326, Springer, New York (2007)[14] M. Dutra, O. Lourenço, J. S. Sá Martins, A. Delfino, J.R. Stone and P.D. Stevenson,
Skyrme interac-tion and nuclear matter constraints , Physical Review C 85, 035201 (2012)[15] Bhaskar Biswas, Rana Nandi, Prasanta Char, and Sukanta Bose, Phys. Rev. D 100, 044056 (20219)[16] Fabian Gittins, Niels Andersson, and Jonas P. Pereira, Phys. Rev. D 101,103025 (2020)[17] G. Baym, H.A. Bethe and C. Pethik,
The ground state of Matter at High Densities: Equation of Stateand Stellar Models , Astrophysical Journal 170, 299 (1971)
18] S. B. Ruester, M. Hempel and J. Schaffner-Bielich,
Outer crust of non-accreting cold neutron stars
Phys. Rev. C 73, 035804 (2006)[19] P. Ring and P. Schuck,
The nuclear many-body problem , 3rd printing, Springer-Verlag Heidelberg,New York (2004)[20] W. Greiner and J. Maruhn,
Kernmodelle , Verlag Harry Deutsch, Frankfurt am Main (1995)[21] M. Bender,
Exotische Atomkerne im Hartree-Fock Modell , Dissertation, Frankfurt am Main (1997)[22] W. Greiner,
Quantenmechanik: Spezielle Kapitel , Edition Harri Deutsch, Frankfurt am Main (1989)[23] C. Titin-Schnaider and P. Quentin,
Coulomb exchange contribution in the nuclear Hartree-Fock calcu-lations , Phys. Lett. 49B, 397-400 (1974)[24] P.-G.Reinhard and H. Flocard, nuclear forces and isotopes shifts , Nucl. Phys. A584, 467-488 (1995)[25] T.J. Bürvenich, I.N. Mishustin and W. Greiner,
Nuclei embedded in an electron gas , Physical ReviewC 76, 034310 (2007)[26] C. Ebel, T.J. Bürvenich and I.N. Mishustin,
Nuclear structure calculations for neutron-star crust , Eur.Phys. J. A 54, 27 (2018)[27] K. Langanke, J. A. Maruhn and S. E. Koonin (Eds.)
Computational nuclear physics 1, Nuclear struc-ture , Springer, Berlin Heidelberg (1991)[28] K. Rutz,
Struktur von Atomkernen im Relativistic-Mean-Field-Modell , Dissertation, Frankfurt amMain (1999)[29] F. Douchin and P. Haensel,
A unified equation of state of dense matter and neutron star structure , As-tronomy an Astrophysics 380, 151 (2001)[30] H. Krivine, J. Treiner and O. Bohigas,
Derivation of a fluid-dynamical Lagrangian and electric giantresonances , Nucl. Phys A 336 155 (1980)[31] J. Bartel, P. Quentin, M. Brack, C. Guet and H.-B. Kansson,
Towards a better parametrization ofSkyrme-like effective forces: A critical study of the SkM force , Nucl. Phys. A 386 79 (1982)[32] A.H. Wapstra and K. Bos,
Atomic Data and Nuclear Data Tables ,19 177 (1977)[33] J. Friedrich and P.-G. Reinhard,
Skyrme-force parametrization:Least-squares fit to nuclear ground-state properties , Phys. Rev. C 33 335 (1986)[34] M. Beiner, H. Flocard, Ngyen Van Giai and Quentin,
Ground-state properties and self-consistent cal-culations with Skyrme interaction , Nucl. Phys. A 238 29 (1975)[35] E. Chabanat, P. Bonche, P. Haensel, J. Meyer and F. Schaeffer,
Skyrme parametrization from subnu-clear to neutron star densities Part II: Nuclei far from stability , Nucl. Phys. A 635, 231 (1998)
36] F. Weber:
Private communications [37] X. Roco-Maza, J. Piekarewicz, T. Garcia-Galvez and M. Centelles,
Influence of the nuclear symmetryenergy on the structure and the composition of the outer crust , IN: C. Bertulani et al.,
Neutron Star Crust103, Nova Science Publishers, New York (2012)