Cumulant Structure Factors of Jellium

For the ground state of the homogeneous electron gas (jellium), it is shown how the cumulant decomposition of the 2-matrix leads to the cumulant decomposition of the structure factors S a,p ( q ) for the antiparallel (a) and parallel (p) spin pairs and how it simultaneously allows one to derive the momentum distribution n ( k ), which is a one-body quantity [Phys. Rev. A 86, 012508 (2012)]. The small- q and large- q behavior of S a,p ( q ), and their normalizations are derived and compared with the results of P. Gori-Giorgi et al . [Physica A 280, 199 (2000) and Phys. Rev. B 61, 7353 (2000)].


Introduction
This paper deals with the ground state of an extended spin-unpolarized (paramagnetic) homogeneous electron gas (HEG), which is one of the most widely studied systems with correlated electrons [1]. The gas is assumed to be uniform in space with the electron density v r e r = [2], the momentum distribution ( ) n k , and the static structure factors (SFs) ( ) a,p S q with "a" standing for electron pairs with antiparallel spins and "p" for pairs with parallel spins. For recent papers on the HEG, also known as jellium, see for example [3]- [15]. Recent GWA and QMC calculations for ( ) n k are [3] [4]. (For comparison with experiment see [5]. The value of F 0. 6 z ≈ for r s = 3.9 agrees with the results, obtained in Refs. [6] [7].) Overhauser considered the pair densities PDs ( ) a,p g r , which follow from the SFs by Fourier transform, and parameterized them in terms of certain 2-body wave functions (known as geminals) [16]. More general geminals should allow one to calculate not only the 2-body quantities ( ) a,p S q , but also the 1-body quantity ( ) n k . Searching in this direction, one immediately encounters the hierarchy of the reduced density matrices (RDMs) and their cumulant decomposition. (The concept of RDM is of great importance in quantum chemistry, see [17]- [22].) Here it is assumed that the 2-body RDM (2-matrix) is available from perturbation theory or otherwise (hierarchy truncation etc.) and, with the spin structure as a decisive point, it is elucidated how the contraction has to be performed for an extended system in such a way that the 1-body density ( ) n k is obtained. In view of the thermodynamic limit, this procedure calls for the cumulant decomposition of the 2-matrix into its Hartree-Fock part HF γ and its nonreducible remainder, the cumulant 2-matrix c γ , which proves to be the source of both ( ) n k and ( ) a,p S q [10]. Diagonalization of c γ yields cumulant geminals and the corresponding weights (i.e. the spectral resolution). A side product of this geminal analysis is the following. Usually the partitioning p 0 p S S S = + ∆ involves an unperturbed (ideal Fermi gas) term S 0 and the interaction-induced term p S ∆ arising from the pairwise Coulomb repulsion [15]. In contrast to that approach, an alternative partitioning is presented here as a consequence of the cumulant decomposition of the 2-matrix, namely p HF p S S C = − , with the Hartree-Fock term HF S [constructed from the correlated ( ) n k ] and a nonreducible remainder, i.e. the cumulant parallel-spin SF p C . The advantage of such a partitioning stems from the fact that the cumulant SFs a,p C are given within perturbation theory by linked (and therefore size-extensive) diagrams with two open particle-hole lines and arbitrary numbers of closed loops. As usual, the divergent terms are eliminated by a partial summation known as "the random phase approximation" (RPA) with the physical meaning of long-range correlations (screening), which becomes exact at the high-density limit of 0 s r → [23] [24]. RPA is a general method to treat the quantum-mechanical many-body problem, successfully applied not only to the HEG. It is now implemented in electronic structure codes for both solids and molecules [25]. A recent discussion concerns its failure when applied to the dissociation of 2 H : it is exact at the dissociation limit, but very inaccurate at intermediate internuclear distances, see Figure 1 in [26]. This RPA has been applied to the HEG quantities e, v , n(k), and ( ) a,p S q . There have been recent efforts to obtain n(k) and S(q) beyond RPA [7] [15]. In those publications, analytical constraints and fitting procedures to quantum-Monte-Carlo data (as "amends" for experimental data) have been combined for the best possible results. Other attempts to go beyond RPA and quantum-Monte-Carlo calculations are the operatorproduct expansion (OPE) [27] and the machine learning methods [28]. In the following text, the RPA is combined with the cumulant analysis described above. Altogether, the available data on the static SFs ( ) a,p S q of the HEG are compiled and analyzed from a new viewpoint, namely the cumulant decomposition, continuing earlier attempts in that direction [7] [10] [29]- [32]. From a formal point of view, the present study is closely related to the work of P. Gori-Giorgi et al. [15].
The paper is organized as follows. In Section 2, the 2-matrix 2 γ and its cumulant decomposition is intro- g r via the Fourier transform is elucidated in Appendix A, where the cusp and curvature theorems are summarized as well. In Section 3, the antiparallel-spin SFs are derived within RPA [the necessary technical details (particlehole propagator, etc.) being given in Appendix B]. The normalization of ( ) a S q yields the short-range correlation parameter ( ) a 0 g . In Section 4, the cumulant decomposition of the parallel-spin SFs with a complicated exchange integral is presented. The fourth moment of ( ) p S q yields another short-range correlation parameter, namely ( ) p 0 g′′ . Table 1 and Table 2 summarize the small-and large-q behavior and the normalization of the SFs. Section 5 concludes the paper with a summary and an outlook.

Basic Equations
]. The GGSB analysis [15] is based upon the decomposition ( ) In addition to the normalization integrals of the last column, it holds that vanishing at the high-density limit 0 s r → , which corresponds to "no Coulomb repulsion". The linear-cubic combination 3 3 4 16 q q − that appears several times in the following derivations may be regarded as the ideal-Fermi-gas term.
( ) 0 S q has a simple geometrical meaning: it is the average volume accessible to the occupied momentum when excited by a transfer momentum of amplitude q [33]. The interaction energy in the lowest order [i.e. with ( ) ( ) It is interesting to note that the plasmon sum rule [34] requires ( ) S q to possess the small-q asymptotics of where q is the momentum transfer, η is the energy or frequency.
In this paper, an alternative is presented, namely Note that what is called , reducible to the 1-body quantity ( ) n k , which follows from the exact 1-matrix, and the non-reducible "cumulant" component ( ) In this cumulant decomposition, the Coulomb repulsion 2 s r q enters the SF ( ) S q at two places: 1) the "Fock term" ( ) F q via the correlated momentum distribution ( ) n k and 2) the cumulant SF ( ) C q via linked diagrams, whereas in the GGSB decomposition the Coulomb repulsion is hidden only in ]. For plots of ( ) F q see Figure 1 in paper I of [35]. The Fock term ( ) F q brings in the momentum distribution ( ) n k with its nonanalytical peculiarities, caused by the RPA long-range correlation. The complexity of the HEG manifests itself in nonanalyticities or singularities of the type ln s s r r or ln q q , see e.g. ( ) This gives rise to the k-or q-dependences in relation to the dependence on r s . For example, according to the perturbation theory, ( ) n k should behave like 2 s r for small r s . This is true for k values away from the Fermi edge (i.e. 1 k  and 1 k  ).
However, at the Fermi edge ( ) n k involves terms like , see [29]. The cumulant decomposition seems to be more complicated than its conventional counterpart but, on the other hand, it may be considered a more natural partitioning. This is so because, within the many-body perturbation theory, the cumulant SF ( ) C q is given by linked diagrams, whereas ( ) S q ∆ contains also unlinked diagrams.
As already mentioned above, the spin-summed SF is given by ( ) The short hands The hermiticity of 2 γ is obvious and the permutational antisymmetry means that 2 2 γ γ → − when 1 and 2 or 1' and 2' are interchanged.
holds for the normalization. The corresponding contraction sum rule This matrix obviously has the spin structure respectively. Here and in the following, the arguments ( ) | , | ′ ′ r r r r of these matrices are suppressed for the sake of convenience. Their diagonals (for which 1 1 ′ = r r and 2 2 ′ = r r ) are denoted by see also Appendix A. A consequence of Equation (2.14) and Equations (A.6)-(A.8) is that the short-range correlation determines the asymptotic large-q behavior of the SFs and their normalizations.
Where do the cumulant 2-matrices χ ± come from? As mentioned above, they are given by linked diagrams.
with two open particle-hole lines (one running from 1 ′ r to 1 r and another one from 2 ′ r to 2 r ) is associated with an exchange diagram Note the slightly different definitions of the spin-summed quantities χ , h , and g on one hand, and C and S on the other.
which can be written as p

Antiparallel Spins: The Coulomb Hole
We move through the chain d d d a a h C C S χ → → → → , where "a" means antiparallel spins and "d" means "direct" diagrams (of the type depicted in Figure 2(a)). In its Cartesian space representation, dr χ is given by at the high-density limit [38]- [40]. Actually, it is exactly this high-density behavior of ( ) h r and ( ) g r for 0 r → that results from the ladder theory, which constitutes the best method to treat the short-range correlation [41]- [46] (for the recent developments see [27]). The correctness of Equations  How does in agreement with an ansatz of Iwamoto [47]. The coefficients 4,5 d are due to multipair and quasiparticle-qua- 0.02875 d = × = [48]. Note that the region c q q  shrinks with decreasing s r and finally vanishes for 0 s r → . As ( ) a S q contains (in the paramagnetic gas) half of the plasmon term, the other half has to come from ( ) p S q , which also has to compensate the ( s r -independent) ideal-Fermi-gas term.
In the transitional region between the small and large values of q there is a peculiar point. Namely, upon approaching from below the transition momenta of 2 = q , the topology changes from two overlapping to two non-overlapping Fermi spheres. This topology change causes a jump discontinuity in the ellipsis representing the terms beyond RPA. In the Cartesian space these jump discontinuities cause the Friedel oscillations in ( ) a g r . They make interpolations between 0 q → and q → ∞ with Padé approximants [49] or the robust interpolation scheme [50] Are there perhaps corrections due to terms beyond RPA? Indeed, the electron-electron coalescence cusp theorem

Parallel Spins: The Fermi Hole
For the parallel-spin, SF it is appropriate to write it as p F ∆ is known to some extent (at least sufficient for a qualitative discussion). As it can be seen from Figure 1 of [35], it is The quantity c is referred to as the Löwdin parameter or the index of non-idempotency. [P.
-O. Löwdin was first to address the meaning of the trace of the squared 1-matrix.] It is a correlation parameter that measures the correlation-induced non-idempotency of the momentum distribution ( ) n k and vanishes for 0 s r → . For the dependence of c on s r see [7]. For the kinetic part of the correlation energy The small-and large-q asymptotics of ( ) respectively. The small q behavior has been demonstrated with the decomposition ( ) ( ) ( ) Together with ( ) x p h C C χ → → → , where "p" means parallel spins and "x" stands for exchange. According to Figure 3(a), the exchange counterpart to dr χ is given by  (4.6) and (after a tedious derivation) its Fourier transform turns out to be equal 2  3  3  pl  xr  1  2  1 2  xr  2  2  2  2  2  2 2  2  2 2  2  1  2   1,2  1,2  1  2  1  2  1   2  2  1  2   d d  d  3  ,  π  ,  4π  4π 3 1, 1, Unfortunately this complicated integral is not known so far as an explicit function of q and s r . If it is approximately substituted by ( ) x1 C q , then the integration over u yields the well-known energy denominator for the exchange term: ( ) x I q is referred to as the Gutlé function [33]; for its definition and properties see Equations (B.9), (B.10), and , v q η is necessary to remove the long-range divergences of the bare interaction ( ) v q , this is not the case for x2 v as it follows from This is the famous result of Onsager, Mittag, and Stephen [54] [55]. The total interaction energy equals with 0 3 10 t = . Next the small-q behavior of ( ) xr C q is studied. Since the integral that enters Equation (4.7) is not known in a closed form, understanding the small-q behavior of p C is not as straightforward as that of a where the first ellipsis stands for higher-order terms [of the order ( ) 4 O q or ( ) 5 O q ] of Equation (3.3) and the second ellipsis means the beyond-RPA terms What contributes to the second term as a function of q and s r remains to be studied. If we assume that its power expansion begins with ( ) 4 O q then, in analogy with Equation (3.7), one obtains and thus (again for in agreement with the plasmon sum rule and [15]. Note that ( ) | , | χ ′ ′ r r r r (see Figure 3(a) and also Figure 3(a) in [10]). Upon carrying out the contraction procedure 2 2 ′ = r r , followed by integration over 2 r and the Fourier transform for 1 1 ′ − r r , the source quantity for the momentum distribution ( ) n k emerges, see [10]. On the other hand, upon taking the diagonal elements 1 1 ′ = r r and 2 2 ′ = r r and carrying out the Fourier transform for 12  perhaps with an additional correlation parameter 2 c . In comparison with Equation (3.10), the 6 1 q tail is strongly enhanced (i.e. four times larger). Note that the difference of the prefactors is just The asymptotics of the SFs are summarized in Table 1 and  , is in agreement with the slope of the ideal-Fermi-gas term at 0 q → .

Summary and Outlook
The focus of this paper is on the mathematics of the HEG, which on one hand is only a marginal point in the broad realm of correlated systems lacking rigorous solutions, but on the other constitutes an archetype of an extended many-body Fermi system. The small-q and large-q behaviors of the cumulant structure factors ( ) a,p C q and the structure factors ( ) a,p S q as well as their normalizations are summarized in Table 1 and Table 2, respectively. In these tables, "a" stands for the electron pairs with antiparallel spins and "p" stands for the pairs with parallel spins. This analysis is based upon rigorous constraints such as the perfect screening sum rule (or the charge neutrality condition), the plasmon sum rule with its inflexion-point trajectory, the on-top (the zero electron-electron distance) theorems for the pair densities (cusp for "a" and curvature for "p") and the "old" RPA (with the Feynman diagrams like those depicted in Figure 2 . A formal analysis is presented with the spin structure (2.10) as a decisive point, leaving the task to parameterize ( ) a,p C q analytically as functions of q and s r in a manner similar to that carried out in the GGSB paper [15]. The agreement of Equation (3.5) with the results of Macke [23] and Gell-Mann/Brueckner [24] and of Equation (4.10) with the calculations of Onsager et al. [54] also confirms the expressions (3.1) for "dr" and (4.5) for "xr", where "r" means RPA and "d" and "x" mean direct and exchange Feynman diagrams, respectively, see Figure 2(a) and Figure 3(a). A complete and self-consistent RPA description needs both ( ) dr C q , which is available from Equation S S S = + , leaving twice the half plasmon term. One may ask how the Fock term F ∆ , respectively the combination x C F − ∆ influences the PD ( ) p g r ? It may be that in the present paper a deeply lying confrontation emerges, namely between 1) perturbation theory, which is linear in the Feynman diagrams, and 2) the cumulant decomposition, for which higher-order reduced density matrices in terms of products of lower-order ones are characteristic.
Going back to the starting point, which rests upon the belief that the geminals and the corresponding weights that diagonalize the cumulant 2-matrix, are the most "natural" 2-body quantities to describe an extended manybody Fermi system, establishes the direction of the future studies. If these quantities are known, then they determine the structure factors, the interaction energy v and-according to [10]-also the momentum distribution ( ) n k in terms of closed-form expressions. The next task would be to derive an effective approximate 2-body scheme (a modified 2-body Schrödinger equation or the Bethe-Salpeter equation?) from the hierarchy of contracted Schrödinger equations such that it yields the cumulant geminals and their weights.