Reference system approach within the white-dwarfs theory

Reference system approach of non-relativistic electron fluid theory was adapted for calculation of characteristics of the electron-nuclear model at the densities typical of degenerate dwarfs. Two- and three correlation functions of degenerate relativistic electron gas have been calculated in the momentum-frequency representation in the local field approximation. Main contributions of the Coulomb interactions to the energy and equation of state of the model at $T=0$ K have been calculated in the adiabatic approximation.


Introduction
A hundred years have passed since the discovery of degenerate dwarfs [1]. The theory of an internal structure of cold dwarfs was developed by Chandrasekhar in 1940s, and it was based on the equation of state of ideal relativistic electron gas at T = 0 K [2,3]. Generalization of this theory followed in the next decades, when in the works of many authors there were investigated effects of such important factors as axial rotation [4,5], Coulomb interactions [6], incomplete degeneration of an electron subsystem [7,8], effects of magnetic fields [9,10], effects of general relativity [5,11], processes of neutronization [12], etc. Interpretation of the whole diversity of properties of the dwarfs obtained from the observations of space [13,14], requires constructing a general theory that also takes into account the effects of the above mentioned factors, among which there are the competing ones.
The effect of interactions that play an important role in determining the structure of dwarfs at different masses and luminosities, and especially for the case of massive cold dwarfs, is the least studied. Basing on Wigner-Seitz, Thomas-Fermi approximations and non-relativistic random phase approximation corresponding to the approximate accounting of two-particle electron correlations, Salpeter [6] showed that Coulomb interactions lead to a small decrease of pressure of degenerate relativistic electron gas at T = 0 K, which is still considered to be the basis of Chandasekhar's theory [10].
Due to the high density (∼ 10 5 g/cm 3 ), the matter in the cores of degenerate dwarfs has a metal electron structure with completely collectivized electrons, and the Fermi momentum is of the order m 0 c. Thus, the non-relativistic approach is not applicable, and it complicates the calculations, because the correlation functions of the reference system (of an interacting relativistic electron gas) are less investigated than the functions of the analogical non-relativistic model. On the other hand, the electron-nuclear model at "dwarfs" densities is slightly non-ideal, which allows one to use a random phase approximation in order to account for the interactions. A consistent approach that corresponds to the modern metals many-electron theory is proposed in [15]. A disadvantage of this work is the use of approximate expressions for twoand three particle correlation functions of ideal relativistic degenerate electron gas. In the present work we use an accurate static two-electron correlation function and a three-particle correlation function in the long-wavelength approximation, which improves the reliability of the results for the energy of the ground state of the model and equation of state.

General relations
We consider a more realistic, compared to [15], spatially homogeneous electrically neutral model, which consists of N e electrons, N 1 nuclei of charge z 1 and N 2 nuclei of charge z 2 in the volume V in thermodynamic limit N e , V → ∞, N e /V = const at the temperatures much lower than the temperature of the electron subsystem degeneration. Generalization for a larger number of nuclei species is obvious.
For the Hamiltonian of the model we use a secondary quantization representation for the electrons and a coordinate one for the nuclear subsystem:Ĥ 0 = k,s E k a + k,s a k,s (2.2) is the Hamiltonian of free electrons (E k = [(m 0 c 2 ) 2 + 2 k 2 c 2 ] 1/2 − m 0 c 2 ), is the operator of electron-electron interactions, is the operator of electron interactions with i-th nuclear subsystem, is the sum of direct nuclear interactions. Here, is the structure factor of i-th nuclear subsystem, a + k,s , a k,s are the creation and annihilation operators of electrons in quantum states with the given vector k and the spin variable s, R i l is the radius-vector of l-th nucleus with the charge z i . For the calculation of a partition function of an electron subsystem in the fixed nuclei field in the grand canonical ensemble we summarize the reference system approach developed in [16,17] for a description of non-relativistic model of electron liquid. In the formula (2.7), N e is the operator of the number of electrons, µ is the 43001-2 variable of chemical potential. As in [16], let us move to the "frequency" representation of electron operators . . . . In the "frequency" representation andĤ µ =Ĥ 0 − µN e . The operator T β is a generalized ordering operator, moreover, ν = 2πnβ −1 ; n = 0, ±1, ±2, . . . , and the calculation of the average of the product of operators a k,s (ν * ) is performed according to the rule [18,19] where G k,s (ν * ) = (iν * − E k + µ) −1 exp (iδν * ) is the spectral image of one-electron Green's function of an ideal system (δ → +0), which is the reference system for a description of the interacting electron gas model. We use model with the Hamiltonian H 0 +V ee as a reference system to calculate Z(µ): The symbol Â 0 denotes a statistical averaging over the states, and Â e for the states of the reference system; is the grand potential of an ideal system of electrons, and Ω e (µ) is the grand potential of the interacting electron gas. Expanding the operatorŜ en (ν) in the formula (2.12) in a power series, averaging by states for the reference system and presenting the result in an exponential form, we have obtained contributions to the grand potential of electron-nuclear interactions ×μ n (q 1 , q 2 , . . . , q n |0)δ q 1 +...+q n ,0 , (2.14) is a cumulant part of the average of the product of n operatorsρ q,0 . This is a static limit of n-particles correlation function of the reference system µ n (y 1 , . . . , y n ) =μ n (q 1 , . . . , q n |ν 1 , . . . , ν n ) = β −1 T β {ρ q 1 ,ν 1ρ q 2 ,ν 2 . . .ρ q n ,ν n } c e . (2.16) The functions (2.15), (2.16) are a generalization of "many-tails" of non-relativistic theory of the electronion model of metals [16,17] and have a well-defined physical meaning. In particular, the functions (2.16) are a spectral representation of n-particles correlation functions, which are given in the coordinate space. For example, the binary distribution function of the interacting electron gas F 2 (r) is associated with the functionμ 2 (q, −q|ν, −ν) of the expression The dynamic correlation functions of n-particles (2.16) are determined by the polarization operators M n (y 1 , . . . , y n ) of n-particles [17] These ratios generalize the well-known random phase approximation, in which M 2 (y, −y) ⇒μ 0 2 (y, −y), M n (y 1 , . . . , y n ) ⇒μ 0 n (y 1 , . . . , y n ), wherẽ is the spectral image of correlation functions of an ideal electron gas model, M n (y 1 , . . . , y n ) =μ 0 n (y 1 , . . . , y n ) for the n 3, so the problem of calculating the functions (2.16) is reduced to the calculation of correlation functions of the ideal relativistic gas model (2.19) and the local field correction function G(y) for the relativistic interacting electron gas.

43001-4
According to the definition (2.19), the functionsμ 0 n (y 1 , . . . , y n ) are a convolution of one-particle Green's functions, the same as in a non-relativistic case: Factorizing the products of Green's functions and using the ratio we obtain representation functionsμ 0 n (y 1 , . . . , y n ):

Two-particle correlation function
Rewriting the sum over vector k via integral using a spherical coordinate system and integrating by the angular variables, we obtain a representation: Here, there appear dimensionless variables In the static case, As in the non-relativistic case, the function J 2 (q * , 0) has a weak logarithmic peculiarity of type (x − 1 2 q * ) ln |x − 1 2 q * | in the vicinity of q * = 2x (|q| = 2k F ). In general, the correlation functions of a relativistic model are similar to the corresponding functions of a non-relativistic model, which is clearly visible from the asymptotic of the functionμ 0 2 (q, −q|0): in the relativistic case,μ . A peculiarity of the functionμ 0 2 (q, −q|0) is its strong dependence on the relativistic parameter, as shown in figure 1.
In the dynamic case, theμ 0 2 (y, −y) can be represented as the approximatioñ

43001-6
Reference system approach that only slightly deviates near the maximum from the numerical result, calculated by formula (3.4). In the non-relativistic limit and the long-wavelength approximation C(q * |x) → 2, v → ν(2ε F q/k F ) −1 , and the function (3.9) coincides with the known expression for a non-relativistic functionμ 0 2 (y, −y) [15]. The expression (3.9) is here only for completeness. We used the result of numerical calculations for the dynamic functionμ 0 2 (y, −y) to obtain the correlation energy of the reference system (the model of interacting electron gas).

Three-particle correlation function
The functionμ 0 3 (y, −y, 0), which is a partial case of three-particle dynamic functions when q 2 = −q 1 , ν 2 = −ν 1 , at a full degeneration has an exact analytical image: (3.10) In the static limit,μ In the formula (3.11), a "non-relativistic" scale was used for the wave vector (q ≡ |q|/k F ). In the long-wavelength limit,μ (3.12) Dependence of the dimensionless factor J 3 (q, −q, 0|x) on the wave vector and the relativistic parameter is illustrated in figure 2. As in non-relativistic case, the function (3.11) has a logarithmic peculiarity at q = 2k F . The formulae (3.6), (3.9) and (3.11), and figures 1, 2 reveal the general property of correlation functionsμ 0 n (y 1 , . . . , y n ) -a steep decrease for large wave vectors (|q i | > 2k F ), providing a convergence of integrals in the series (2.14).
In general, for a rough estimation of convergence of series (2.14) we consider a chemically homogeneous model (z 1 = z 2 = z, N 1 + N 2 = N n ), constraining the integration for each independent vector q i of the area |q i | < 2k F , neglecting the screening interactions, replacing the product of structural factors S q 1 S q 2 . . . S q n with N n , we replace the functionsμ 0 , which approximately corresponds to the long-wavelength asymptotic. For the magnitude of n-th member of series (2.14) it can be estimated as where α 0 = e 2 / c is the fine structure parameter. Hence, the series (2.14) is an expansion for a dimensionless parameter zα 0 , which varies from 0.014 (helium dwarf) to 0.19 (iron dwarf). For the typical dwarfs, mainly consisting of nitrogen and oxygen, zα 0 ≈ 0.1. That expansion parameter is a small value, which makes it possible to restrict the consideration to two-and three-electron correlations (we note, that correlation energy of the reference system is of the order of α 2 0 ). Moreover, for the three-electron functionμ 0 3 (q 1 , q 2 , −q 1 − q 2 |0), there can be used an approximate analytical expression, because the main 43001-7 contributions provide the two-electron correlations, and the contribution of three-electron correlations play a role of a correction. The calculation of correlation functionsμ 0 3 (y 1 , y 2 , −y 1 − y 2 ) andμ 0 4 (y 1 , −y 1 , y 2 , −y 2 ) in the static or dynamic cases in the non-relativistic theory is based on Feynman identity [22], which allows one to integrate over angular variables of vector k at a fixed configuration of vectors q 1 and q 2 in terms of γ 3 (y i , y j ). Unfortunately, this identity cannot be used in the exact calculation due to the complex dependence of relativistic electron energy on the wave vector.
In order to make an approximate calculation using the identity transformation, we present γ 3 (q i , q j ) as which is asymptotically correct both at small and at large q i . According to the Feynman identity, Rewriting the sum over vector k via integral, we use dimensionless variablek = |k|/k F ,q i = |q i |/k F , and the spherical coordinate system, the Oz axis of which coincides with the vector ρ i j , we perform integration over the angular variables, reducing γ 3 (q i , q j ) to one-dimensional integral: Here, the following notations are introduced: , q R is the radius of the circle, circumscribing the triangle constructed on the vectors q 1 , q 2 , −q 1 − q 2 ; t i j is a cosine of the angle between the vectors q i , q j . In the formulae (3.17), (3.18), the variables k and q i , q j are dimensionless (in unit k F ). The substitution C(q i |k) → 2(1 + x 2 k 2 ) 1/2 corresponds to long-wavelength approximation, which allows one to rewriteμ 0 3 (q 1 , q 2 , q 3 |0) in a compact form: In this formula, ; Dependence of functionμ 0 3 (q 1 , q 2 , −q 1 − q 2 |0)(3N e ) −1 (m 0 c 2 x 2 ) 2 on wave vectors for different t are illustrated in figures 3, 4.

The local field correction function
It is well known from the non-relativistic electron fluid theory [16] that the local field correction function (LFCF) in the weak non-ideal model is a universal function of the variable y = (q, ν). It does not depend on any parameters and corresponds to the approximation (3.21) We have calculated the LFCF, using the correlation functions of the relativistic ideal electron model. Summarizing over the frequencies ν 1 and ν * (which appears inμ 0 4 (y, −y, y 1 , −y 1 ) and using the formula (3.2), we obtain a representation: The behavior of relativistic LFCF at small and large values of vector q is the same as that of non-relativistic case [23]: . . by q ≫ k F and any ν. For the numerical calculation at the absolute zero temperature, we have used a cylindrical coordinate system for the vectors k 1 , k 2 [k j = (z j , ρ j , ϕ j ), k 2 j = z 2 j + ρ 2 j ], in which (k j , q) = qz j , (k 1 − k 2 ) 2 = ρ 2 1 + ρ 2 2 + (z 1 − z 2 ) 2 − 2ρ 1 ρ 2 cos (ϕ 1 − ϕ 2 ), (k 1 , k 2 ) = z 1 z 2 + ρ 1 ρ 2 cos (ϕ 1 − ϕ 2 ). Integrating over the angular variables ϕ 1 , ϕ 2 , we reduce G id (y) to the following 4-dimensional integral: Here, the following notations are used: For a comparison G id (q, ν) with LFCF of non-relativistic theory, there were used the variables q = |q|k −1 F , ν = x 2ν ,ν = ν( 2 k 2 F /m) −1 . In figure 5 there is shown a function G id (q, ν), calculated by the formula (3.24) atν = 0.01, which is very close to the static limit. It is obvious that the asymptotic behavior G id (q, ν) for small and large values of the |q| almost does not depend on the relativistic parameter. The deviation of a relativistic correction from non-relativistic one is significant near its maximum, which monotonously decreases with an increase of relativistic parameter. Figure 6 illustrates the behavior of G id (q, ν) at a very high value of frequency (ν = 0.5). In general, the behavior of the correction corresponds to the one of non-relativistic theory, and some deviations are caused by frequency renormalization at x 1.

The energy of ground state
As we know from the theory of non-relativistic electron fluid, the transition from thermodynamic description to the quantum mechanics is performed by the following procedure: instead of the chemical potential of interacting system µ there should be used a chemical potential of the ideal system µ 0 , and in the Ω(µ) there should be considered only the so-called "normal" components, because the chemical potential shift µ − µ 0 is compensated with contributions of the so-called "abnormal" components, which arise in the correlation functionsμ 0 n (y 1 , . . . , y n ) at n 4. In the Green's function approach, this fact was first discovered in [24] and in the reference system approach -in [17] within the perturbation theory.
Taking into account a weak non-ideal matter of dwarfs, we will consider only two-and three-particle electron correlations. In this approach, the ground state energy of the model is given as follows: . . S (i n ) −q n δ q 1 +...+q n ,0μn (q 1 . . . q n |0). Here, E e is the ground state energy of the reference system (an interacting relativistic electron gas), which can be calculated using the expression whereμ λ 2 (y, −y) is the two-particle dynamic correlation function of the auxiliary model of electrons with the potential of interactions λV q , and E 0 is the energy of the ideal electron system Extracting from the second term of formula (4.2) the contribution of ideal correlation 2V q 0 k,s V q n k+q/2,s n k−q/2,s (4.4) and considering, that G id (y) does not depend on the "coupling constant" λ, we have represented E e in a traditional form

43001-12
where is the so-called correlation energy (the contribution of non-ideal correlations). In units m 0 c 2 we obtain where is the contribution of an ideal system per one electron, −3α 0 x(4π) −1 is the contribution of interactions of the Hartree-Fock approximation, α 2 0 ε c (x) is the correlation contribution. According to our calculations, the ε c (x) can be approximated with the following expression: At a 1 = a 2 = d 0 = 0, the expression matches the approximation for the correlation energy [25], which is calculated using the Monte-Carlo method [26] ε MC c (x). As we can see in figure 7, in the range x 1, the expression (4.9) is close to ε MC c (x), and the deviation ε In order to calculate the contributions of electron-nuclear interactions in the products of structure factors in the formula (4.1), we have selected one-particle and two-particle sums by the coordinates of nuclei, and the three-nuclear effective interactions were neglected. In this approach, pol (x), 10 · ε pol (x) and ε c (x) on the relativistic parameter.

43001-13
where E pol is the polarization energy of electron fluid by nuclei, which does not depend on the structure of the nuclear subsystem, (4.11) The configuration energy is determined by the structure of a nuclear subsystem, and expressed through an effective two-particle potential of interactions of the nuclei, which is formed by two-and three-electron correlations: (4.12)

The approximation of two-electron correlations
Let us rewrite the component E (2) pol , calculated in the local field approximation, in the form where the dimensionless function ε (2) pol (x) is of the same order as ε c (x), and z n = (z n 1 N 1 + z n 2 N 2 )(N 1 + N 2 ) −1 . The function ε (2) pol (x) can be approximated by the following expression: (4.14) As it is shown in figure 7, ε (2) pol (x) has linear asymptotics at x ≫ 1, as well as ε c (x). However, the polarization energy E (2) pol exceeds the correlation energy of the reference system by about z α −1/2 0 ≈ 10 z times, and for z ∼ 10 it is comparable with E HF .
The configuration energy was calculated in the coordinate representation and, introducing the effective two-nuclear interactions, yielded

43001-14
In the formula (4.15) the sum of the vector q includes a component with q = 0. Therefore, To simplify the calculation of the lattice sum, we adopt a simple model of nuclei distribution: 17) where N j is the number of all knots on the j-th coordination sphere, and N (i) j is the number of the knots occupied by nuclei with charge z i (i = 1, 2). In this model, where R j is the radius of the j-th coordination sphere. The effective two-particle potential is screened, and at small and medium distances between nuclei it is close to the expression 19) and the screening radius is of the order 0.1a B (where a B = 2 /m 0 e 2 is the Bohr radius). At large distances, V 2 (R) oscillates, but with a small amplitude, The configuration energy for a simple cubic lattice of nuclei is calculated numerically and can be represented as with a dimensionless factor approximated by the expression where all the coefficients a 1 , . . . , a 6 are the functions of z , that is The cofficients a i j of the formula (4.24) are listed in table 1. In figure 8 there is shown a dependence of the configuration energy on the relativistic parameter.

The effect of three-electron correlations
According to the formula (4.11), the contribution of three-particle correlations is represented as The result of numerical calculations is illustrated in figure 7, which shows that for x > 1 the ratio ε (3) pol (x) ≈ 0.1ε (2) pol (x) is satisfied. Therefore, at sufficiently large values of nuclei charges, E (3) pol is not less than the correlation energy of the reference system: at z 6, the contribution E pol is close to the correlation energy, at z 12, it exceeds the correlation energy by 5 times, and at z = 26 -more than by 20 times. The result of numerical calculation ε (3) pol (x) is approximated by the expression Similar to the formula (4.15), we have introduced a correction to the effective two-nuclear potential by the three-electron correlations (4.27) where the sum over the vector q includes the component q = 0. As shown in figure 9, V 2 (R) is a weak attracting potential of the type of quantum package screening potential, which is close to the expression The functions A(x) and γ(x) are determined as follows: The contribution to the configuration energy of the model by the three-particle correlations in the model (4.17) takes the form, similar to the formula (4.18): This contribution is calculated for a simple cubic lattice of nuclei and is represented by L (x| z ).  The dependence of a dimensionless factor ε (3) L (x| z ) is illustrated in figure 10. At sufficiently large nuclei charges z and x 2, the function ε L (x| z ), but it has a positive sign. It is approximated by the expression The coefficients of the formula (4.32) are shown in table 2.

43001-18
x Here, is the contribution of the ideal degenerate relativistic spatially homogeneous electron gas; is the contribution of Coulomb interactions in the two-electron correlations approximation; is the contribution of the three-particle electron correlations. In the region x 1, all contributions [with the exception of ε (3) L (x| z )] to the model energy caused by interactions are negative monotonously decreasing functions of the relativistic parameter. In the twoelectron correlations approximation, the equation of state (5.2) numerically is very close to the result of Salpeter [6].

Conclusions
Reference system approach was adapted for the description of a degenerate relativistic electron subsystem. There were investigated the features of two-and three-particle correlation functions in a wide domain of a relativistic parameter, as well as there was obtained an exact expression for the static twoparticle correlation function. For the first time, the expressions for the three-particle correlation function were approximated and the local field correction of interacting relativistic electron gas was studied, which is the basis for the calculation of the energy and structural characteristics of degenerated dwarfs. The energy of the ground state of the electron-nuclear model, as well as the equation of state of the model, have been calculated in a wide range of the relativistic parameter at absolute zero temperature. As it was  Figure 11. The ratio of pressure P(x) to the pressure of the ideal relativistic electron gas P 0 (x) as function of the relativistic parameter and nuclear charge (curve 1z 1 = z 2 = 2, 2z 1 = z 2 = 6, 3z 1 = 7; z 2 = 8; 4z 1 = z 2 = 12).
shown in our calculation, the contributions of Coulomb interactions to the energy of the ground state and pressure, caused by two-electron correlations are important and increase with an increase of nuclear charge. The contributions caused by three-electron correlations are much smaller but they exceed the contribution of correlation energy of the electron fluid, especially at large values of the nuclear charge.