# Single-particle excitations in the uniform electron gas by diagrammatic
Monte Carlo

Kristjan Haule  Department of Physics and Astronomy, Rutgers University,
Piscataway, New Jersey 08854, USA<EMAIL_ADDRESS>Kun Chen
Department of Physics and Astronomy, Rutgers University, Piscataway, New
Jersey 08854, USA  Center for Computational Quantum Physics, Flatiron
Institute, 162 5th Avenue, New York, NY 10010. The Flatiron Institute is a
division of the Simons Foundation

###### Abstract

We calculate the single-particle excitation spectrum and the Landau Fermi
liquid parameters for the archetypal model of solids, the three-dimensional
uniform electron gas, with the numerically exact variational diagrammatic
Monte Carlo method. In the metallic range of density, we establish benchmark
values for the wave-function renormalization factor $Z$, the effective mass
$m^{*}/m$, and the Landau parameters $F_{0}^{s}$ and $F_{0}^{a}$ with
unprecedented accuracy, and we resolve the long-standing puzzle of non-
monotonic dependence of mass on density. We also exclude the possibility that
experimentally measured large reduction of bandwidth in Na metal can originate
from the charge and spin fluctuations contained in the model of the uniform
electron gas.

## Introduction

The uniform electron gas (UEG) is the most fundamental model for understanding
the electronic properties of metallic materials. The ground-state properties
of the model have been very precisely calculated by quantum Monte Carlo
methods [1], and this allowed one to build approximate density functionals [2,
3], which are at the heart of the ab-initio approaches in material science and
modern theory-driven materials design. The knowledge of the low energy
excitations of the same model remain challenging to evaluate accurately [4, 5,
6, 7, 8, 9, 10], even though such calculations are important for building more
sophisticated density functionals [11, 12, 13], and these excitations are
directly measured in experiments on simple metals, such as alkaline materials.
Some aspects of the excitation spectra, such as the quasiparticle
renormalization amplitude, were recently determined by extention of the
variational Monte Carlo method in Ref. [14], which turn out to be in
remarkable good agreement with our current results.

In the metallic regime, the low-energy properties of the electron liquid are
dominated by the long-lived quasiparticles near the Fermi surface, and their
dynamics is described by a handful of the Fermi liquid parameters. These
parameters completely characterize the low energy excitation spectra of the
metallic state. Unfortunately, they are very challenging to calculate by a
first principle approach, therefore they are usually treated as
phenomenological parameters requiring input from experiments.

Here we develop an extension of the recently introduced variational
diagrammatic Monte Carlo (VDMC) method [15], which fills this void, and allows
us to determine the single-particle excitations of UEG with unprecedented
accuracy. In this letter, we calculate the single-particle excitation spectra,
and in particular, we give controlled values of the wave-function
renormalization factor $Z$, the quasiparticle effective mass ratio $m^{*}/m$
and also the Landau Fermi liquid parameters $F_{0}^{a}$ and $F_{0}^{s}$. Our
computed values are free of systematic error, and their uncertainty is mainly
controlled by the statistical error, and hence our established value can be
used as a precise benchmark for new method development. Moreover, these
precise Fermi liquid parameters are also useful for building more
sophisticated density functionals. Finally, the method we develop here can be
used to solve more sophisticated models, and can also be used in the ab-initio
framework on models of realistic materials, a development which is currently
underway [16].

## Results

### The Feynman expansion algorithm

The VDMC method [15] is a flavor of diagrammatic Monte Carlo method (DMC) [17,
18, 19, 20, 21, 22, 23, 24], which samples high-order Feynman diagrams with a
Monte Carlo importance sampling. The novelty of VDMC is two-fold: i) it
optimizes the starting point of the perturbative expansion in such a way that
the expansion converges very rapidly with the increasing perturbation order.
ii) it efficiently combines an exponentially large number of Feynman diagrams,
which mostly cancel among themselves due to alternating fermionic sign so that
the groups of diagrams can be efficiently sampled with the Monte Carlo
importance sampling hence avoiding the explosion of statistical error with
perturbative order.

In Ref. 15 we computed the spin and the charge response functions of the UEG
model with VDMC by evaluating the Feynman diagrams for the polarization
function. A similar type of Feynman expansion in terms of non-interacting
single-particle Green’s function, and statically screened Coulomb interaction
does not converge rapidly enough to establish a reliable infinite order
result, hence we here develop an alternative approach.

Figure 1:  Feynman diagrams for the self-energy in terms of the three leg
vertex $\Gamma_{3}$, which is expanded in bare series in terms of
$G_{\mathbf{k}}^{0}$ and partially screened interaction
$v_{\mathbf{q}}=\frac{8\pi}{q^{2}+\lambda}$ and counter-terms
$(\frac{\lambda}{8\pi})^{N}(\frac{8\pi}{q^{2}+\lambda})^{N+1}$.

In this work, we show that extremely rapid convergence with perturbation order
can be achieved by using a Hedin-type equation, in which we first compute the
numerically exact screened interaction $W_{\mathbf{q}}$ (previously developed
in Ref. 15), and we then expand only the three-point vertex function
$\Gamma_{3}$ in powers of the bare electron propagator $G^{0}_{\mathbf{k}}$,
and statically screened interaction $v_{q}(\lambda)$, with proper counter
terms defined in the Method section. Here the screened Coulomb interaction
$v_{q}(\lambda)$ has a Yukawa form, characterized by the inverse screening
length $\lambda$. This screening parameter has to be determined by the
principle of minimal sensitivity in order to achieve rapid convergence of the
perturbative series, so that the extrapolation to infinite order is possible.
Fig. 1 shows the sketch of the corresponding Feynman diagrams up to the third
order. Below we apply the algorithm to the UEG model, although the method is
completely general and could as well be carried out for realistic material in
the ab-initio framework.

### The single particle excitations

Figure 2:  The wave-function renormalization factor $Z$ versus screening
parameter $\lambda$ for various perturbation orders $N=1...5$ and for
$r_{s}=1,2,3$ and $4$. The insets show the convergence of $Z$ with
perturbation order $N$ when its value is taken at the extremal $\lambda$. The
numbers next to each point show the value of $\lambda$ used for each
calculated point. Panel d) compares current VDMC results with prior Monte
Carlo results from Ref. [14] and G0W0 from Ref. [25].

We first present the single-particle excitation spectral results. Fig. 2a-c
show how the wave-function renormalization factor $Z$ depends on the screening
parameter $\lambda$ in our theory. To determine the optimized parameter
$\lambda$, we scan $Z(\lambda)$ for each $r_{s}$, and determine it with the
principle of minimal sensitivity. For efficiency, we here sample the self-
energy only at the Fermi wave vector $k_{F}$ and at the two lowest Matsubara
frequencies, which is sufficient to determine $Z$. We notice that for the
first two orders, no counter term in the parameter $\lambda$ occurs, therefore
the curve $Z(\lambda)$ displayed in Fig. 2 does not have extremum, while all
higher-order terms have a well-defined maximum, which broadens and develops
into a broad plateau with increasing order. The insets of Figs. 2a-c show
optimized $Z$ versus perturbation order, where the first two orders are
evaluated at the optimal $\lambda$ of the third order, and for later orders,
we take the value in the maximum. We also display the value of $\lambda$ used
at each order. From Fig. 2 it is apparent that beyond order three the rate of
convergence to limiting value of $Z$ is extremely fast, and therefore we can
confidently determine the first three digits of $Z$. The values and the
estimated error-bar from the extrapolation and statistical errors are shown in
Table 1.

$r_{s}$ | $Z$ | $m^{*}/m$ | $F_{0}^{a}$ | $F_{0}^{s}$
---|---|---|---|---
1 | 0.8725(2) | 0.955(1) | -0.171(1) | -0.209(5)
2 | 0.7984(2) | 0.943(3) | -0.271(2) | -0.39(1)
3 | 0.7219(2) | 0.965(3) | -0.329(3) | -0.56(1)
4 | 0.6571(2) | 0.996(3) | -0.368(4) | -0.83(2)
Table 1:  Landau liquid parameters: The wave-function renormalization factor
$Z$, effective mass $m^{*}/m$, and the Landau parameters $F_{0}^{a}$,
$F_{0}^{s}$ for various values of the density parameter $r_{s}$ with the
estimated error.

In Fig. 2d we compare our computed $Z(r_{s})$ with the previous best available
estimates, obtained by various flavors of Monte Carlo (MC) methods, which are
reproduced from Ref. [14]. Note that all these published MC methods rely on
fixed node approximation and the thermodynamic limit extrapolation, hence they
have an inherent systematic error, nevertheless they turn out to be in very
good agreement with current VDMC results. Our current work based on VDMC has
only statistical error, and a small error in extrapolating in perturbation
order, and is thus far more precise than previous best results. We notice that
previous MC results are broadly consistent with our results, with SJ-VMC
method predicting slightly too large and BF-VMC and BF-RMC slightly too small
value. It is also well known that G0W0 predicts quite accurate $Z$ values,
however, we can now confidently claim that in the range of metallic densities,
G0W0 consistently underestimates $Z$.

Figure 3:  Electron effective mass: The upper panel shows our calculated
effective mass versus perturbation order for $r_{s}=1-4$. The lower panel
compares the $r_{s}$ dependence of the effective mass of this work (VDMC) with
the prior analytic and numeric work from Ref. 26.

Once the extremal value of $\lambda$ is determined, we compute the entire
momentum and frequency dependence of the self-energy, which allows us to
determine also the momentum derivative of the self-energy, and hence the
effective mass of the electron through the relation

$\displaystyle\frac{m}{m^{*}}=Z\left(1+\frac{m}{k_{F}}\frac{d\Sigma(k_{F},\omega=0)}{dk}\right)$
(1)

The convergence of the effective mass ratio $m^{*}/m$ with perturbation order
is shown in Fig. 3a, and its dependence on $r_{s}$ is displayed in Fig. 3b.

The dependence of the effective mass $m^{*}/m$ on $r_{s}$ has been
controversial for many decades. Some theories predict that the ratio is
monotonically decreasing with increasing $r_{s}$ [7, 27], while others predict
the existence of a turning point $r_{s}^{*}$ [26, 28, 29, 30, 31] at which the
trend is reversed. Our controlled results confirm the correctness of the later
theories. Furthermore, we compare our controlled VDMC results with previous
best estimates, which are based on the theory of many-body local field factors
[26]. This theory includes vertex corrections associated with charge and spin
fluctuations, extracted from available Monte Carlo data. We notice that G0W0
overestimates the effective mass in the entire range of metallic densities.
The density fluctuations beyond RPA are included in theory with $G_{+}$ local
field corrections, which reduce the mass substantially and bring it very close
to our VDMC results at small $r_{s}$. However, beyond $r_{s}>3$ our VDMC
results are closer to the theory which contains both the charge and the spin
fluctuations ($G_{+}\&G_{-}$), hence we can infer that at moderate
correlations strength, the spin fluctuations start to play an important role,
and charge fluctuations are no longer sufficient in determining the mass of
the electron gas.

### The Landau liquid parameters

With precisely calculated effective mass, as well as the spin and charge
susceptibility determined in our previous work [15], we can calculate Landau
parameters $F_{0}^{a}$ and $F_{0}^{s}$, which are obtained from
$\frac{\chi_{s}}{\chi_{s}^{0}}=\frac{m^{*}}{m}\frac{1}{1+F_{0}^{a}}$ and
$\frac{P_{q=0}}{P_{q=0}^{0}}=\frac{m^{*}}{m}\frac{1}{1+F_{0}^{s}}$ . Here
$\chi_{s}$ and $P_{q}$ are the spin susceptibility and charge polarization,
while $\chi_{s}^{0}$ and $P_{q}^{0}$ are their non-interacting analogues. In
table 1 we list our calculated Landau parameters $F^{a}_{0}$ and $F^{s}_{0}$,
together with the estimation of their error, which mostly comes from error in
determining spin and charge susceptibility in Ref. 15. While the Landau
parameters, which determine the interaction between quasiparticle, have been
estimated by various approximate numerical methods before [7], to our
knowledge their numerically controlled value has not be obtained before.

### The spectral function and the bandwidth

Figure 4:  The spectral function and $\Sigma_{k=0}(\omega)$ at $r_{s}=4$ and
$k=0$ as relevant for bandwidth of Na metal.

The present VDMC algorithm also allows us to compute a numerically controlled
value for the dynamic self-energy on the imaginary axis. Analytic continuation
is needed to obtain the self-energy on the real frequency axis. We use the
maximum entropy method to compute the quasiparticle energy at the $k=0$ point,
which determines the bandwidth of the electron dispersion, i.e., the energy
difference between the Fermi level and the lowest possible quasiparticle
energy. In Fig. 4 we display the self-energy, as well as the spectral function
at momentum $k=0$ and finite frequency. We notice that the imaginary part of
the self-energy starts to grow rapidly when the energy of the single-particle
excitations exceeds the plasma frequency $\omega_{p}\approx 1.881E_{F}$.
Consequently, there appears a strong pole at $\omega\approx-2.4E_{F}$ due to
such plasma excitations, and makes quasiparticle approximation invalid at a
frequency below $\omega<-E_{F}$, as the real part of the self-energy is no
longer a linear function of frequency. However, around $E_{F}$ the real-part
of $\Sigma$ is still quite close to a linear function, and only minor
deviations are noticed. Consequently, the renormalization of the dispersion
can not substantially deviate from our earlier estimation of $m^{*}/m$, which
is valid at the Fermi level. Our numerical estimation based on the
analytically continued self-energy is that the spectral function at $r_{s}=4$
and $k=0$ has a peak around $-0.96E_{F}$, which deviates from the non-
interacting value for only 4%, hence the bandwidth reduction due to
interactions at $r_{s}=4$ is only of the order of 4%. This value is much
smaller than the experimental estimation of the bandwidth reduction in Na
metal, in which the measured ARPES bandwidth appears to be renormalized for
about 18-25% [32, 33]. However, our estimated bandwidth is definitely not
substantially larger as compared to the non-interacting bandwidth, in contrast
to several other many-body calculations [34, 35], and is neither substantially
smaller as in early GW calculations [36] or GW with paramagnon vertex
corrections [8]. Based on our very precise estimation of the single-particle
self-energy, we can confidently exclude a possibility of such a dramatic
reduction of the bandwidth in the model of electron gas due to correlation
effects at the density corresponding to Na metal. This large reduction of the
effective mass in ARPES thus requires an alternative explanation, which was
assigned to the interaction in the final states [37, 35] in ARPES, surface
effects [38], and possibly the lattice effects, i.e, deviation of Na metal
from the continuous model of the uniform electron gas.

In summary, we established the low energy excitation spectrum of the uniform
electron gas at metallic density using recently developed VDMC. Controlled
values of $Z$, $m^{*}/m$, $F_{0}^{s}$, and $F_{0}^{a}$ are given, which agree
with the state of the art calculations in the field, but here we provide much
more precise values than previously known.

## Methods

The Hamiltonian of UEG problem is

$\displaystyle\hat{H}=\sum_{{\mathbf{k}}\sigma}\left({{\mathbf{k}}^{2}}-\mu\right)\hat{\psi}^{\dagger}_{{\mathbf{k}}\sigma}\hat{\psi}_{{\mathbf{k}}\sigma}+$
$\displaystyle\frac{1}{2V}\sum_{\begin{subarray}{c}{\mathbf{q}}\neq 0\\\
{\mathbf{k}}{\mathbf{k}}^{\prime}\sigma\sigma^{\prime}\end{subarray}}\frac{8\pi}{q^{2}}\hat{\psi}^{\dagger}_{{\mathbf{k}}+{\mathbf{q}}\sigma}\hat{\psi}^{\dagger}_{{\mathbf{k}}^{\prime}-{\mathbf{q}}\sigma^{\prime}}\hat{\psi}_{{\mathbf{k}}^{\prime}\sigma^{\prime}}\hat{\psi}_{{\mathbf{k}}\sigma},$
(2)

where $\hat{\psi}$/$\hat{\psi}^{\dagger}$ are the annihilation/creation
operator of an electron, $\mu$ is the chemical potential controlling the
density of the electrons in the system, and the long-range Coulomb repulsion
is $8\pi/q^{2}$, as we measure the energy in units of Rydbergs, and the wave
number $k,q$ in units of inverse Bohr radius.

The expansion in terms of the bare interaction is divergent, therefore we
first transform the original problem into an equivalent but a more appropriate
problem for power expansion, which describes the emergent degrees of freedom
at the lowest order, and the corrections are perturbatively included with very
rapid convergence. Motivated by the well-known fact that the long-range
Coulomb interaction is screened in the solid and that the effective potential
of emerging quasiparticles differs from the bare potential, we introduce the
screening parameter $\lambda_{\mathbf{q}}$ and an electron potential
$v_{\mathbf{k}}$ into the quadratic part of the emergent Lagrangian $L_{0}$ of
the form

$\displaystyle
L_{0}=\sum_{{\mathbf{k}}\sigma}\psi^{\dagger}_{{\mathbf{k}}\sigma}\left(\frac{\partial}{\partial\tau}-\mu+{\mathbf{k}}^{2}+v_{\mathbf{k}}(\xi=1)\right)\psi_{{\mathbf{k}}\sigma}$
$\displaystyle+\sum_{{\mathbf{q}}\neq
0}\phi_{-{\mathbf{q}}}\frac{q^{2}+\lambda_{\mathbf{q}}}{8\pi}\phi_{\mathbf{q}},.$
(3)

We then add the following interacting part to the Lagrangian

$\displaystyle\Delta
L=-\sum_{{\mathbf{k}}\sigma}\psi^{\dagger}_{{\mathbf{k}}\sigma}v_{\mathbf{k}}(\xi)\psi_{{\mathbf{k}}\sigma}-\xi\sum_{{\mathbf{q}}\neq
0}\phi_{-{\mathbf{q}}}\frac{\lambda_{\mathbf{q}}}{8\pi}\phi_{\mathbf{q}}$
$\displaystyle+\sqrt{\xi}\frac{i}{\sqrt{2V}}\sum_{{\mathbf{q}}\neq
0}\left(\phi_{{\mathbf{q}}}\rho_{-{\mathbf{q}}}+\rho_{{\mathbf{q}}}\phi_{-{\mathbf{q}}}\right).$
(4)

so that, when the number $\xi$ is set to unity, $L(\xi)=L_{0}(\xi)+\Delta
L(\xi)$ is Lagrangian of UEG. Indeed integrating out the bosonic fields
$\phi_{\mathbf{q}}$ from Lagrangian $L$, we get the Lagrangian corresponding
to the original Hamiltonian Eq. (2). Here $\rho_{\mathbf{q}}$ is the density
fluctuation of the problem
$\rho_{\mathbf{q}}=\sum_{{\mathbf{k}}\sigma}\psi^{\dagger}_{{\mathbf{k}}\sigma}\psi_{{\mathbf{k}}+q\sigma}$.
Note that the first two terms in $\Delta L$ are the counterterms [39] which
exactly cancel the two terms we added to $L_{0}$ above. We use the number
$\xi$ to track the order of the Feynman diagrams so that order $N$
contribution sums up all diagrams carrying the factor $\xi^{N}$. We set
$\xi=1$ once we enumerate all the diagrams of a certain order.

The emergent screening length $\lambda_{\mathbf{q}}$ and effective potential
$v_{\mathbf{k}}$ are not a-priory known and need to be properly optimized to
achieve an optimal speed of convergence. We note in passing that determining
those parameters self-consistently, i.e., $\lambda_{\mathbf{q}}$ from the
self-consistent polarization, and $v_{\mathbf{k}}$ from the single-particle
self-energy, is not the most optimal choice for the speed of convergence.
Determining them by the principal of minimal sensitivity is a much better
choice, as pointed out by Kleinert and Feynman [40, 41, 42, 43]. They showed
that when an effective parameter of a theory is optimized with this principle,
the perturbative expansion converges very fast, and can succeed even when the
interaction is strong, and regular perturbation theory fails.

To make algorithm sufficiently simple to implement, we take
$\lambda_{\mathbf{q}}$ to be ${\mathbf{q}}$ independent constant ($\lambda$),
which is already sufficient for rapid convergence of the series. We emphasize
that for any choice of these parameters we are guaranteed to converge to the
same answer, provided that the series converges. Furthermore, we found that
the convergence of the expansion is best when the Fermi surface of both the
dressed $G_{\mathbf{k}}$ and the bare $G^{0}_{\mathbf{k}}$ Green’s function at
each order is fixed with the Luttinger’s theorem so that the density and the
Fermi surface volume is not changed with the increasing perturbation order. We
therefore, expand $v_{\mathbf{k}}$ in power series
$v_{\mathbf{k}}=(\Sigma^{x}_{\mathbf{k}}(\lambda)-\Sigma^{x}_{k_{F}}(\lambda))+\xi\,s_{1}+\xi^{2}\,s_{2}\cdots$,
and we determine $s_{N}$ so that all contributions at the order $\xi^{N}$ do
not alter the physical volume of the Fermi surface. Similarly to optimizing
$\lambda_{\mathbf{q}}$, we found that it is sufficient to take $s_{N}$
constants independent of the momentum. Since the exchange
($\Sigma^{x}_{\mathbf{k}}$) is static and is typically large, we accommodate
it at the zeroth-order into the effective potential, so that at the first
order we recover the GW type self-energy with $G_{\mathbf{k}}$ at the screened
Hartree-Fock (screened by screening length $\lambda$) and exact
$W_{\mathbf{q}}$.

As mentioned before, the algorithm depicted in Fig. 1 needs a numerically
exact (converged) $W_{\mathbf{q}}$, which is first computed with the algorithm
of Ref. 15. It was shown in Ref. 15 that the most rapidly converging scheme
for charge and spin-susceptibilities is the so-called vertex correction
scheme, in which we add an infinite sum of ladder diagrams on both sides of a
polarization Feynman diagram. To do that, we first precompute the three-point
ladder vertex and then attach it to both sides of a polarization Feynman
diagram while the diagrams are sampled, and at the same time, we eliminate all
ladder-type diagrams from the sampling, to avoid double-counting of diagrams.
Next, we use Hedin’s type equation depicted in Fig. 1 in which one fermion
propagator is dressed and requires self-consistent $G$. It easy to see that it
is sufficient to use bold $G$ of the lower order $N-1$ when evaluating self-
energy at order $N$, to avoid the expensive self-consistent calculation.
Finally, we use the finite temperature imaginary-time formalism, and we set
the temperature to $T=0.04\,E_{F}$, which is sufficiently below the Fermi
liquid scale, so that is essentially equivalent to zero temperature.

## References

  * [1] Ceperley, D. M. & Alder, B. J.  Ground state of the electron gas by a stochastic method.  _Phys. Rev. Lett._ 45, 566–569, DOI: 10.1103/PhysRevLett.45.566 (1980).
  * [2] Jones, R. O. & Gunnarsson, O.  The density functional formalism, its applications and prospects.  _Rev. Mod. Phys._ 61, 689–746, DOI: 10.1103/RevModPhys.61.689 (1989).
  * [3] Perdew, J. P., Burke, K. & Ernzerhof, M.  Generalized gradient approximation made simple.  _Phys. Rev. Lett._ 77, 3865–3868, DOI: 10.1103/PhysRevLett.77.3865 (1996).
  * [4] Simion, G. E. & Giuliani, G. F.  Many-body local fields theory of quasiparticle properties in a three-dimensional electron liquid.  _Physical Review B_ 77, 035131 (2008).
  * [5] Takada, Y.  Inclusion of vertex corrections in the self-consistent calculation of quasiparticles in metals.  _Phys. Rev. Lett._ 87, 226402, DOI: 10.1103/PhysRevLett.87.226402 (2001).
  * [6] Yasuhara, H., Yoshinaga, S. & Higuchi, M.  Why is the bandwidth of sodium observed to be narrower in photoemission experiments?  _Phys. Rev. Lett._ 83, 3250–3253, DOI: 10.1103/PhysRevLett.83.3250 (1999).
  * [7] Yasuhara, H. & Ousaka, Y.  Effective mass, landau interaction function and self-energy of an electron liquid.  _International Journal of Modern Physics B_ 06, 3089–3145, DOI: 10.1142/S0217979292002401 (1992).  https://doi.org/10.1142/S0217979292002401.
  * [8] Zhu, X. & Overhauser, A. W.  Plasmon-pole and paramagnon-pole model of an electron liquid.  _Phys. Rev. B_ 33, 925–936, DOI: 10.1103/PhysRevB.33.925 (1986).
  * [9] Richardson, C. F. & Ashcroft, N. W.  Dynamical local-field factors and effective interactions in the three-dimensional electron liquid.  _Phys. Rev. B_ 50, 8170–8181, DOI: 10.1103/PhysRevB.50.8170 (1994).
  * [10] Kukkonen, C. A. & Overhauser, A. W.  Electron-electron interaction in simple metals.  _Phys. Rev. B_ 20, 550–557, DOI: 10.1103/PhysRevB.20.550 (1979).
  * [11] Perdew, J. P., Kurth, S., Zupan, A. c. v. & Blaha, P.  Accurate density functional with correct formal properties: A step beyond the generalized gradient approximation.  _Phys. Rev. Lett._ 82, 2544–2547, DOI: 10.1103/PhysRevLett.82.2544 (1999).
  * [12] Tran, F. & Blaha, P.  Accurate band gaps of semiconductors and insulators with a semilocal exchange-correlation potential.  _Phys. Rev. Lett._ 102, 226401, DOI: 10.1103/PhysRevLett.102.226401 (2009).
  * [13] Sun, J., Ruzsinszky, A. & Perdew, J. P.  Strongly constrained and appropriately normed semilocal density functional.  _Phys. Rev. Lett._ 115, 036402, DOI: 10.1103/PhysRevLett.115.036402 (2015).
  * [14] Holzmann, M. _et al._ Momentum distribution of the homogeneous electron gas.  _Phys. Rev. Lett._ 107, 110402, DOI: 10.1103/PhysRevLett.107.110402 (2011).
  * [15] Chen, K. & Haule, K.  A combined variational and diagrammatic quantum Monte Carlo approach to the many-electron problem.  _Nature Communications_ 10, 3725, DOI: 10.1038/s41467-019-11708-6 (2019).
  * [16] Haule, K. & Mandal, S.  All electron gw with linearized augmented plane waves for metals and semiconductors.  _arXiv preprint arXiv:2008.07727_ (2020).
  * [17] Prokof’ev, N. V. & Svistunov, B. V.  Polaron problem by diagrammatic quantum monte carlo.  _Physical review letters_ 81, 2514 (1998).
  * [18] Prokof’ev, N. & Svistunov, B.  Fermi-polaron problem: Diagrammatic monte carlo method for divergent sign-alternating series.  _Physical Review B_ 77, 020408 (2008).
  * [19] Van Houcke, K. _et al._ Feynman diagrams versus fermi-gas feynman emulator.  _Nature Physics_ 8, 366 (2012).
  * [20] Houcke, K. V., Kozik, E., Prokof’ev, N. & Svistunov, B.  Diagrammatic monte carlo.  _Physics Procedia_ 6, 95 – 105, DOI: https://doi.org/10.1016/j.phpro.2010.09.034 (2010).  Computer Simulations Studies in Condensed Matter Physics XXI.
  * [21] Kozik, E. _et al._ Diagrammatic monte carlo for correlated fermions.  _EPL (Europhysics Letters)_ 90, 10004 (2010).
  * [22] Deng, Y., Kozik, E., Prokof’ev, N. V. & Svistunov, B. V.  Emergent bcs regime of the two-dimensional fermionic hubbard model: Ground-state phase diagram.  _EPL (Europhysics Letters)_ 110, 57001 (2015).
  * [23] Rossi, R.  Determinant diagrammatic monte carlo algorithm in the thermodynamic limit.  _Phys. Rev. Lett._ 119, 045701, DOI: 10.1103/PhysRevLett.119.045701 (2017).
  * [24] Rossi, R., Ohgoe, T., Van Houcke, K. & Werner, F.  Resummation of diagrammatic series with zero convergence radius for strongly correlated fermions.  _arXiv preprint arXiv:1802.07717_ (2018).
  * [25] Hedin, L.  New method for calculating the one-particle green’s function with application to the electron-gas problem.  _Phys. Rev._ 139, A796–A823, DOI: 10.1103/PhysRev.139.A796 (1965).
  * [26] Simion, G. E. & Giuliani, G. F.  Many-body local fields theory of quasiparticle properties in a three-dimensional electron liquid.  _Phys. Rev. B_ 77, 035131, DOI: 10.1103/PhysRevB.77.035131 (2008).
  * [27] Lam, J.  Thermodynamic properties of the electron gas at metallic densities.  _Phys. Rev. B_ 5, 1254–1260, DOI: 10.1103/PhysRevB.5.1254 (1972).
  * [28] Hubbard, J.  The description of collective motions in terms of many-body perturbation theory.  _Proc. R. Soc. Lond. A_ 240, 539–560, DOI: http://doi.org/10.1098/rspa.1957.0106 (1957).
  * [29] Rice, T.  The effects of electron-electron interaction on the properties of metals.  _Annals of Physics_ 31, 100 – 129, DOI: https://doi.org/10.1016/0003-4916(65)90234-4 (1965).
  * [30] Hedin, L. & Lundqvist, S.  Effects of electron-electron and electron-phonon interactions on the one-electron states of solids.  vol. 23 of _Solid State Physics_ , 1 – 181, DOI: https://doi.org/10.1016/S0081-1947(08)60615-3 (Academic Press, 1970).
  * [31] MacDonald, A. H., Dharma-wardana, M. W. C. & Geldart, D. J. W.  Density functional approximation for the quasiparticle properties of simple metals. i. theory and electron gas calculations.  _Journal of Physics F: Metal Physics_ 10, 1719–1736, DOI: 10.1088/0305-4608/10/8/010 (1980).
  * [32] Jensen, E. & Plummer, E. W.  Experimental band structure of na.  _Phys. Rev. Lett._ 55, 1912–1915, DOI: 10.1103/PhysRevLett.55.1912 (1985).
  * [33] Lyo, I.-W. & Plummer, E. W.  Quasiparticle band structure of na and simple metals.  _Phys. Rev. Lett._ 60, 1558–1561, DOI: 10.1103/PhysRevLett.60.1558 (1988).
  * [34] Maezono, R., Towler, M. D., Lee, Y. & Needs, R. J.  Quantum monte carlo study of sodium.  _Phys. Rev. B_ 68, 165103, DOI: 10.1103/PhysRevB.68.165103 (2003).
  * [35] Takada, Y.  Inclusion of vertex corrections in the self-consistent calculation of quasiparticles in metals.  _Phys. Rev. Lett._ 87, 226402, DOI: 10.1103/PhysRevLett.87.226402 (2001).
  * [36] Northrup, J. E., Hybertsen, M. S. & Louie, S. G.  Theory of quasiparticle energies in alkali metals.  _Phys. Rev. Lett._ 59, 819–822, DOI: 10.1103/PhysRevLett.59.819 (1987).
  * [37] Shung, K. W. K. & Mahan, G. D.  Calculated photoemission spectra of na.  _Phys. Rev. Lett._ 57, 1076–1079, DOI: 10.1103/PhysRevLett.57.1076 (1986).
  * [38] Shung, K. W. K., Sernelius, B. E. & Mahan, G. D.  Self-energy corrections in photoemission of na.  _Phys. Rev. B_ 36, 4499–4502, DOI: 10.1103/PhysRevB.36.4499 (1987).
  * [39] Wu, W., Ferrero, M., Georges, A. & Kozik, E.  Controlling feynman diagrammatic expansions: Physical nature of the pseudogap in the two-dimensional hubbard model.  _Phys. Rev. B_ 96, 041105, DOI: 10.1103/PhysRevB.96.041105 (2017).
  * [40] Stevenson, P. M.  Gaussian effective potential: Quantum mechanics.  _Physical Review D_ 30, 1712 (1984).
  * [41] Stevenson, P. M.  Gaussian effective potential. ii. $\lambda$ $\varphi$ 4 field theory.  _Physical Review D_ 32, 1389 (1985).
  * [42] Stevenson, P. M. & Tarrach, R.  The return of $\lambda$$\varphi$4\.  _Physics Letters B_ 176, 436–440 (1986).
  * [43] Kleinert, H.  _Path integrals in quantum mechanics statistics and polymer physics_ (World Scientific, 1995).

## Acknowledgements

This work was supported by the Simons Foundation through Simons fellowship and
Simons Collaboration on the Many Electron Problem. KH acknowledges supported
of NSF DMR-1709229. The Flatiron Institute is a division of the Simons
Foundation.

## Author contributions statement

K.C. and K.H. both developed two independent codes to crosschecked the
results, K.H. wrote the first draft of the manuscript and K.C. and K.H.
contributed to the manuscript and approved the text.

## Additional information

Accession codes: The source code has been made available for download under
gnu license at: https://github.com/haulek/VDMC

Competing interests: The authors declare no competing interests.