{'fraction_non_alphanumeric': 0.039177142857142855, 'fraction_numerical': 0.016091428571428572, 'mean_word_length': 4.22848948374761, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Quantum Monte Carlo study of the ground state of the two-dimensional Fermi fluid

We have used the variational and diffusion quantum Monte Carlo methods to calculate the energy, pair correlation function, static structure factor, and momentum density of the ground state of the two-dimensional homogeneous electron gas. We have used highly accurate Slater-Jastrow-backflow trial wave functions and twist averaging to reduce finite-size effects where applicable. We compare our results with others in the literature and construct a local-density-approximation exchange-correlation functional for 2D systems.


I. INTRODUCTION
The homogeneous electron gas (HEG) plays a key role in modern condensed-matter theory. It consists of a set of electrons moving in a uniform, inert, neutralizing background, and is the simplest fully-interacting quantum manybody model of condensed matter. The three-dimensional HEG models the behavior of the conduction electrons in metals and semiconductors, while the two-dimensional (2D) HEG models the behavior of electrons confined to surfaces and thin layers. Quantum Monte Carlo (QMC) methods 1,2 have long played an important role in establishing the ground-state properties of the HEG. Most effort has gone into calculating the ground-state energy of different phases as a function of density, in order to establish the zero-temperature phase diagram. [3][4][5] In this article, we report QMC calculations of some other properties of the 2D HEG of interest to condensed-matter physicists: the pair correlation function (PCF), static structure factor (SSF), and momentum density (MD). We also report energy data for highdensity HEG's. We have confined our attention to the fluid phase, which is the ground state at the densities typically encountered in experiments.
The PCF, especially the contact PCF g(0), is a key ingredient in generalized-gradient-approximation exchangecorrelation functionals for density functional theory (DFT) calculations. The PCF has been studied several times using QMC, 3,4,6,7 but the value of g(0) at low densities has proved controversial because electrons approach one another infrequently, and the QMC results disagree with values calculated using ladder theory. 8 Our QMC data, obtained using a different trial wave function from the earlier calculations, should help to clarify the situation. The SSF is related to the PCF by a Fourier transform. SSF data at small k are needed to establish the long-range behavior of the PCF.
The MD of the HEG is of considerable importance in Fermi liquid theory. To our knowledge, the only QMC MD data to have been published for the 2D HEG are those of Tanatar and Ceperley, 3 which used a relatively simple form of trial wave function. (A fit to QMC data generated by Conti 9 is reported in Ref. 10, but no details about the calculations are given.) Tanatar and Ceperley's low-density MD shows a very strange feature: the MD is lower at zero momentum than it is at the top of the Fermi edge. It is clearly important to provide new QMC MD data, in order to establish whether this is a genuine property of the HEG.
Finally, we report ground-state energy data for paramagnetic Fermi fluids, which we use to parameterize a localdensity-approximation exchange-correlation functional for use in DFT studies of 2D systems.
The rest of this article is arranged as follows. In Sec. II we describe the computational techniques used. In Sec. III we present the data we have generated. Finally, we draw our conclusions in Sec. IV. Densities are given in terms of the radius r s of the circle that contains one electron on average. We use Hartree atomic units ( = |e| = m e = 4πǫ 0 = 1) throughout this article. All our QMC calculations were performed using the casino code. 11

A. Trial wave functions
In the variational quantum Monte Carlo (VMC) method, expectation values are calculated with respect to an approximate trial wave function, the integrals being performed by a Monte Carlo technique. In diffusion quantum Monte Carlo 1,2 (DMC) the imaginary-time Schrödinger equation is used to evolve an ensemble of electronic configurations towards the ground state. The fermionic symmetry is maintained by the fixed-node approximation, 12 in which the   I: Energy, variance, and percentage of correlation energy retrieved using different methods for a 58-electron paramagnetic Fermi fluid of density parameter rs = 5 a.u. Twist averaging has not been used. "HF," "SJ-VMC," "SJB-VMC," "SJ-DMC," and "SJB-DMC" stand for Hartree-Fock theory, VMC with a Slater-Jastrow wave function, VMC with a Slater-Jastrow-backflow wave function, DMC with a Slater-Jastrow wave function, and DMC with a Slater-Jastrow-backflow wave function, respectively. The DMC energy data have been extrapolated to zero time step. The data marked with an asterisk were produced by Kwon et al., 19 while the data marked with a dagger were generated by Attaccalite et al. 21 at a time step of 0.1 a.u. (i.e. their data were not extrapolated to zero time step). The fraction of the correlation energy retrieved is computed on the assumption that our Slater-Jastrow-backflow DMC calculation retrieves 100% of the correlation energy.
nodal surface of the wave function is constrained to equal that of a trial wave function. The VMC algorithm generates electron configurations distributed according to the square of the trial wave function, while the DMC algorithm generates configurations distributed as the product of the trial wave function and its ground-state component. Our trial wave functions consisted of Slater determinants of plane-wave orbitals multiplied by a Jastrow correlation factor. The Jastrow factor contained polynomial and plane-wave expansions in electron-electron separation. 13 The orbitals in the Slater wave function were evaluated at quasiparticle coordinates related to the actual electron positions by backflow functions consisting of polynomial expansions in electron-electron separation. 14 The wave functions were optimized by variance minimization 15,16 and linear-least-squares energy minimization. 17 We simulated HEG's in finite, square cells subject to periodic boundary conditions. The many-body Bloch theorem 18 states that the wave function Ψ satisfies Ψ(r 1 , . . . , r i + R s , . . . , r N ) = exp(ik s · R s )Ψ(r 1 , . . . , r N ), where R s is a simulation-cell lattice vector and k s is the simulation-cell Bloch vector. In some of our calculations, and in previous QMC studies of the 2D HEG, 3,4,19 it has been assumed that k s = 0. However, in our calculations of the energy, PCF, and SSF we performed twist averaging, in which expectation values are averaged over k s in the first Brillouin zone of the simulation cell. 20 This procedure greatly reduces single-particle finite-size errors caused by shell-filling effects. The high quality of our trial wave functions is demonstrated in Table I, which shows QMC energies achieved using different levels of wave function for a 58-electron paramagnetic Fermi fluid of density parameter r s = 5 a.u. Backflow functions change the nodal surface of the trial wave function and can therefore improve the fixed-node DMC energy. In practice we find that backflow lowers the DMC energy substantially. Our VMC energies are significantly lower than those of Kwon et al., 19 as is our Slater-Jastrow-backflow DMC energy. On the other hand, the Slater-Jastrow-backflow DMC energy of Attaccalite et al. 21 is higher than that of Kwon et al. Our Slater-Jastrow DMC energy is slightly lower than that of Attaccalite et al., which in turn is lower than that of Kwon et al. Since the nodal surface is the same in the three calculations, these DMC energies really ought to agree. However, the trial wave function used by Kwon et al. is very much poorer than ours, as can be seen by comparing the VMC energies. Time-step and population-control biases in their DMC energies must be much greater, which may explain the discrepancy. The results of Attaccalite et al. have not been extrapolated to zero time step, hindering comparison. The VMC and DMC results of Rapisarda and Senatore 4 are in very close agreement with those of Kwon et al. 19 We have optimized a three-electron term in the Jastrow factor (together with the two-electron Jastrow terms and backflow functions) for a paramagnetic 58-electron HEG at r s = 5 a.u. The three-electron term lowered the nontwist-averaged VMC energy from −0.1493851 (6) to −0.1495111(5) a.u. per electron. The DMC energies at a time step of 0.1 a.u. without and with the three-body Jastrow factor are −0.149742(2) and −0.149740(2) a.u. per electron, respectively. As expected, the inclusion of the three-body term makes an insignificant difference to the DMC energy, because the DMC energy depends only on the nodal surface of the trial wave function, which is not directly affected by the Jastrow factor. We have therefore not used three-electron terms in our production calculations.

Evaluating the MD
Let Ψ(R) be the trial many-electron wave function, where R = (r 1 , . . . , r N ). Suppose the first N ↑ electrons are spin-up and the remainder are spin-down. The MD of spin-up electrons can be evaluated as where the angled brackets denote an average over the set of electron configurations generated in the VMC and DMC algorithms (which are distributed as |Ψ| 2 and |ΨΦ 0 |, respectively, where Φ 0 is the ground-state component of Ψ). We have restricted our attention to paramagnetic and fully ferromagnetic HEG's, so the total MD is equal to the spin-up MD. The integral in the expectation value of Eq. (2) is estimated by Monte Carlo sampling at each configuration R generated by the QMC algorithms, and the results are averaged. The use of a finite number of points in the evaluation of the integral at each R does not bias the QMC estimate of ρ(k). Suppose our finite simulation cell has area A and the simulation-cell Bloch vector is k s . We may write where the {G} are the simulation-cell reciprocal lattice points. Hence it is clear that ρ(k) is only nonzero if k = G+k s for some G. The MD is only defined for a discrete set of momenta at any given k s . One cannot twist average as such; instead, altering k s leads to the MD being defined at a different set of momenta. We simply report MD's obtained using k s = 0 (i.e., no twist was applied).

Evaluating the SSF
The SSF may be evaluated as is the Fourier transform of the density operator. S(k) is only nonzero at simulation-cell G vectors, even if the simulation-cell Bloch vector is nonzero. We can twist average when we calculate S(G).

Evaluating the PCF
The spherically averaged PCF is which can be evaluated by binning the electron-electron distances in the configurations generated by the QMC algorithms. Twist averaging introduces no complications.

Extrapolated estimation
IfÂ is an operator that does not commute with the Hamiltonian then the errors in the VMC and DMC estimates A VMC and A DMC of the expectation value ofÂ are linear in the error in the trial wave function; however, the error in the extrapolated estimate 2A DMC − A VMC is quadratic in the error in the trial wave function. 2 We have used extrapolated estimation in most of our calculations of expectation values. Examples of extrapolation are shown in Figs. 5 and 8, and the upper panel of Fig. 4. In each case the VMC, DMC, and extrapolated estimates are in good agreement, implying that the error resulting from the use of a DMC mixed estimate is small. Gori-Giorgi et al. 7 used reptation 22 QMC to accumulate the PCF and SSF, in which pure expectation values are obtained with respect to the fixed-node ground-state wave function, so that extrapolation is unnecessary.

C. Time-step and population-control biases
Finite-time-step errors in the twist-averaged DMC energy were removed by linear extrapolation to zero time step. An example is shown in Fig. 1; it can be seen that the time-step bias is in fact very small in any case. We checked that the other expectation values were converged with respect to the time step: see Figs. 4 and 5. We used a target population of 1600 configurations in all our DMC calculations, making population-control bias negligible.

D. Finite-size bias
Expectation values obtained in a finite N -electron cell subject to periodic boundary conditions differ from the corresponding infinite-system values because of "single-particle" shell-filling effects, as well as the neglect of longranged correlations and the compression of the exchange-correlation hole into the simulation cell. Single-particle finite-size effects can be removed from the energy, the SSF, and the PCF by twist averaging, as explained in Sec. II A. We have recently demonstrated that the finite-size error in the energy per particle in a 2D HEG falls off as N −5/4 , enabling us accurately to extrapolate QMC energies to infinite system size. 23 For the PCF, SSF, and MD we simply verified that the QMC data had converged with respect to system size (see Figs. 5 and 8, and the lower panel of Fig.  4).  Where N = ∞, the DMC energy has been extrapolated to infinite system size.

A. Energies
DMC energies of paramagnetic Fermi fluids at different densities and system sizes are shown in Table II. Our results for the energies of different phases of the 2D HEG at low density are reported elsewhere. 5 At r s = 5 and 10 a.u., Rapisarda and Senatore 4 obtained infinite-system energies of −0.1490(1) and −0.08512(2) a.u. per particle using DMC with a Slater-Jastrow wave function. Kwon et al. 19 obtained DMC energies of −0.2098(3), −0.1495(1), and −0.08536(2) a.u. per electron at r s = 1, 5, and 10, respectively, using a Slater-Jastrow-backflow wave function. Our DMC energies are somewhat lower than these data, as expected from the results shown in Table I.
Let the correlation energy per electron E c be the difference between the ground-state energy per electron and the Hartree-Fock energy. We fit the form proposed by Attaccalite et al. 21 to our correlation energies for paramagnetic HEG's: where A 0 = −0.1925, B 0 = √ 2(10 − 3π)/(3π), and D 0 = −A 0 H 0 . We fit to the infinite-system DMC energies shown in Table II and also to the DMC energies of low-density paramagnetic HEG's reported in Ref. 5 (at r s = 20, 25, 30, 35, and 40 a.u.). Our fitting parameters are shown in Table III (7) can be used as a local-density-approximation exchange-correlation functional in DFT calculations for 2D systems.
Unlike Attaccalite et al., 21 we fit Eq. (7) to paramagnetic data only; we do not attempt to calculate the spinpolarization-dependence of the energy of the HEG.

B. MD's
The MD's of paramagnetic HEG's are shown in Fig. 3, and a more detailed graph of the MD at r s = 5 a.u. is shown in Fig. 4. The upper panel of Fig. 4 demonstrates that the extrapolated estimate is accurate and that the DMC results are converged with respect to the time step. It is clear from the lower panel of Fig. 4 that, although backflow makes a significant improvement to the QMC energy estimates, it has very little effect on the MD. The inclusion of backflow results in a small transfer of weight to wavevectors above the Fermi edge, as expected, because a greater fraction of correlation energy is retrieved. It can also be seen that the MD's obtained at different system sizes are in agreement. We have therefore plotted data obtained at different system sizes together in Fig. 3. To our knowledge, the only previous QMC studies of the MD of the 2D HEG are those of Tanatar and Ceperley 3 and Conti. 9 At r s = 10 a.u., Tanatar and Ceperley found that the MD at small wave vectors is lower than the value near the Fermi edge. Our data do not show this unusual feature. Tanatar and Ceperley used a relatively inflexible Slater-Jastrow wave  (7) for the correlation energy of a paramagnetic HEG. Equation (7) was fitted to the infinite-system DMC energies given in Table II and the low-density energy data given in Ref. 5. It was found that fixing H0 = 0 did not affect the quality of the fit.  function, which may be the reason for the discrepancy. Giuliani and Vignale 10 quote a formula for the MD, which was obtained by fitting to QMC data generated by Conti. 9 We have fitted our MD's to a simplified version of the form suggested in Ref. 10: where x = r s k and the a i are fitting parameters. g(0) is the contact PCF, which we evaluated using Eq. (9). The fitted parameters are given in Table IV. Our values for the discontinuity at the Fermi edge are slightly smaller than those reported in Ref. 10.

C. SSF's
VMC and DMC SSF's of a 50-electron paramagnetic HEG at r s = 1 a.u. are shown in Fig. 5. It can be seen that the difference between the VMC and DMC data is in most cases smaller than the difference between the two sets of DMC data, implying that the errors due to extrapolated estimation are small. On the other hand, the difference between the data with k s = 0 and the twist-averaged data is significant. In particular, the former has some unusual features close to integer multiples of the Fermi wave vector, one of which is shown in the inset to Fig. 5. Elsewhere, twist averaging has only a small effect. At all densities we find that the twist-averaged SSF in a 50-or 58-electron cell is in agreement with the SSF in a 114-electron cell, as can be seen in Fig. 5. Since the statistical errors are less significant at the smaller system sizes, we have used N = 50 or 58 electrons in the high-density data reported below.       The SSF's of paramagnetic and ferromagnetic fluids at low density are shown in Fig. 6, while the SSF's of paramagnetic fluids at high density are shown in Fig. 7. It can be seen that a peak in the SSF at about k = 2.62k F appears at low density, perhaps due to incipient Wigner crystallization. Our SSF's are in good agreement with those of Gori-Giorgi et al. 7

D. PCF's
We compare PCF's obtained at different system sizes using different QMC methods for a paramagnetic HEG at r s = 5 a.u. in Fig. 8. The difference between the VMC and DMC PCF's, and the difference between extrapolated PCF's obtained with and without backflow is small, implying that the error in the extrapolated PCF is small. Twist averaging has a small effect on the PCF, but the twist-averaged PCF's at N = 58 and N = 114 electrons are very similar, implying that the finite-size error in the twist-averaged PCF at N = 58 is small. We have also verified that the PCF is converged with respect to the DMC time step.
Our PCF's are shown in Figs. 9 (high density) and 10 (low density), along with the results of Gori-Giorgi et al. 7 Our PCF's are in good agreement with those of Gori-Giorgi et al. (as expected, from the SSF results in Fig. 7).
The contact PCF g(0) is especially important in the construction of generalized-gradient-approximation exchangecorrelation functionals. 10 We give our g(0) values in Table V, and we plot r s g(0) against r s in Fig. 11, along with some other results in the literature. Because our PCF's have converged with respect to system size, and our VMC and DMC results agree with each other when backflow is used, we have simply averaged our VMC and DMC g(0) data at for paramagnetic HEG's at high density. kF = √ 2/rs is the Fermi wave vector. Slater-Jastrow-backflow wave functions and twist averaging were used. The system sizes are N = 50, 58, and 58 at rs = 1, 5, and 10 a.u., respectively. The curve marked "HF" shows the Hartree-Fock SSF. The dotted lines show the SSF obtained by taking the Fourier transform of the PCF data of Gori-Giorgi et al. 7 different system sizes in order to reduce the statistical noise in our estimate of the contact PCF. Our results are in reasonably good agreement with the fit to the earlier QMC data of Gori-Giorgi et al. 7 and also with the expression for g(0) obtained using ladder theory by Nagano et al. 24 Interestingly, our results clearly disagree with the more recent calculation of g(0) within ladder theory by Qian, 8 which involved fewer approximations than the work of Nagano et al. The close agreement between our results and those of Nagano et al. must therefore be regarded as a coincidence. The fact that our QMC calculations, using a different trial wave function, are consistent with the data of Gori-Giorgi et al. strongly suggests that the QMC results for g(0) are reliable, whereas ladder theory is of limited use at low densities. Our results are also in clear disagreement with the formula proposed by Polini et al., 25 which interpolates between the results of ladder theory at high-density (where it should be exact) and a partial-wave analysis at low density.
The fit to our g(0) data shown in Fig. 11 is 1 + Ar s + Br 2 s exp(−Er s ) r s ≥ 1 1 + ar s + br 2 s + cr 3 where A = −0.25724, B = 0.071116, and E = 0.98553 were obtained by fitting. Polini et al. 25 have shown that

IV. CONCLUSIONS
We have studied the ground-state properties of the fluid phases of the 2D HEG using QMC. We used highly accurate trial wave functions and dealt with finite-size effects by twist averaging. Twist averaging removes some strange features in the SSF, but our PCF's and SSF's are in good agreement with analytic fits to earlier QMC data, 7 confirming the accuracy of these formulas. Our MD's show some qualitative differences from earlier QMC results, 3 however; in particular, we do not observe an increase in the MD as the Fermi edge is approached at low density. Finally, we have reported DMC energy data for the high-density 2D HEG, which we used to construct a new exchange-correlation functional for 2D DFT calculations.

V. ACKNOWLEDGMENTS
We acknowledge financial support from Jesus College, Cambridge and the UK Engineering and Physical Sciences Research Council (EPSRC). Computing resources were provided by the Cambridge High Performance Computing   Table V), Gori-Giorgi et al., 7 Nagano et al., 24 Polini et al., 25 and Qian. 8