In a recent paper, Lucco Castello et al. (arXiv:2107.03537) provided an accurate parameterization of classical one-component plasma bridge functions that was embedded in a novel dielectric scheme for strongly coupled electron liquids. Here, this approach is rigorously formulated, its set of equations is formally derived, and its numerical algorithm is scrutinized. A systematic comparison with available and new path integral Monte Carlo simulations reveals a rather unprecedented agreement especially in terms of the interaction energy and the long wavelength limit of the static local field correction.

The quantum one-component plasma or jellium (also referred to as uniform electron gas or homogeneous electron gas) is universally acclaimed as one of the most important idealized systems in condensed matter physics,1,2 statistical mechanics,3–6 and quantum chemistry.7,8 The jellium can serve as a model for many metals under the simplifying assumption that the charge density of the positive ionic background is uniformly smeared out.2–6 Hence, initial investigations focused on its ground state at metallic densities.9 The ground-state studies have led to remarkable insights such as the BCS superconductivity theory,10 Landau’s Fermi-liquid theory,11 the Bohm–Pines quasi-particle picture of collective excitations,12 Wigner’s crystallization paradigm,13 and even the Kohn–Sham formulation of density functional theory.14 

Increasing interest in warm dense matter,15,16 an exotic state of high-temperature highly compressed matter that is encountered in dense astrophysical objects (giant planet interiors, brown or white dwarfs, and neutron star crusts) and ultra-fast laser heating or shock compression of metals,17,18 has provided the impetus for the worldwide intense research activity targeted at the high-density finite temperature jellium.19 Only recently, an accurate description of the uniform electron gas has been achieved in this regime on the basis of a combination of novel path integral Monte Carlo (PIMC) methods.16,20–23

Owing to the fact that low-density finite temperature inhomogeneous electron systems are inaccessible even to the contemporary state-of-the-art experiments, much less attention has been paid to the strongly coupled uniform electron liquid. In fact, there have been rather few computational and theoretical studies, despite the numerous intriguing physical phenomena that have been speculated to manifest themselves at strong coupling. These include the possibility of a charge-density wave instability,24–26 the emergence of a spin-density wave,26–28 the possibility of a continuous paramagnetic to ferromagnetic transition,29–31 and the emergence of an incipient excitonic mode.32,33 It should be pointed out that this constitutes a particularly challenging regime for theoretical approaches due to the complex interplay between quantum effects (exchange degeneracy and diffraction), thermal excitations, and strong Coulomb correlations.

The thermodynamic properties and static structure of strongly coupled electron liquids were recently probed by extensive PIMC simulations at finite temperatures.34 A dielectric formalism scheme that handles quantum mechanical effects at the random phase approximation level, assumes a frequency independent local field correction, and treats strong Coulomb correlations within the classical hypernetted-chain (HNC) approximation35 was revealed to provide consistently accurate predictions for the interaction energies, although there was ample room for improvement regarding structural properties.34 From the integral equation theory of classical liquids,36–38 it is known that as the coupling increases, the bridge diagrams, which are neglected in the HNC approximation, have an increasing importance for two-particle correlations. Hence, the HNC-based scheme is expected to gradually become inaccurate as crystallization is approached.

In a recent letter,39 we extracted the bridge functions of the classical one-component plasma at multiple thermodynamic state points, spanning the whole dense liquid region, from specially designed molecular dynamics simulations. With this input, we constructed a very accurate closed-form bridge function parameterization that covers the entire non-trivial range. This analytic description was incorporated into a novel dielectric scheme that naturally extends the HNC dielectric scheme by including the exact Coulomb bridge function. A comparison with available and new PIMC simulations of the strongly coupled electron liquid revealed that the novel scheme leads to significant improvements over the HNC-based scheme. It should be noted that this marked the first incorporation of bridge functions into a dielectric scheme, although bridge functions had been earlier embedded in variants of the so-called classical mapping method.40–42 

In the present work, the integral equation theory based dielectric scheme is formally introduced and its emerging self-consistent set of equations is derived. The basic assumptions and main drawbacks of this dielectric scheme are discussed together with possible further refinements. Mathematical techniques are presented that decrease the computational cost and facilitate convergence. Extensive PIMC simulations are analyzed that extend our current picture of the finite temperature strongly coupled electron liquid. A systematic comparison is carried out with other dielectric schemes and “exact” PIMC simulations in terms of interaction energies, static structure factors, static local field corrections, and static density responses.

The uniform electron fluid (UEF) is a homogeneous quantum mechanical system that consists of electrons that are immersed in a rigid charge neutralizing background.43 In other words, the UEF constitutes the quantum mechanical analog of the classical one-component plasma (OCP). Since the electronic interactions involve not only Coulomb interactions but also exchange effects, the distinction between spin-up and spin-down electrons makes the UEF essentially a two-component system. As a consequence, the thermodynamic state points of the UEF are specified by three dimensionless parameters:16,19,43 the quantum coupling or Brueckner parameterrs = d/aB, where d is the Wigner–Seitz radius d = (4πn/3)−1/3 and aB = 2/(mee2) is the first Bohr radius; the degeneracy parameterθ = kBT/EF, where EF=[(6π2n)2/3/2](2/me) is the Fermi energy with respect to the Fermi wave-vector of spin-up electrons kF=(6π2n)1/3; and the spin polarization parameterξ = (nn)/n that is constrained within 0 ≤ ξ ≤ 1 under the spin-up choice nn. In the above description, is the reduced Planck constant, kB is the Boltzmann constant, e is the electron charge, me is the electron mass, n = n + n is the electron density, and T is the temperature. In contrast to the UEF, the thermodynamic state points of the classical OCP are specified by one dimensionless parameter: the classical coupling parameter Γ = e2/(dkBT)19,43–45 for which Γ = 2λ2rs/θ with λ3=(kFd)3=4/(9π).

In what follows, we shall restrict ourselves to the spin-unpolarized (paramagnetic) case, n = n or ξ = 0. The ideal (Lindhard) density response of free electrons, i.e., in the absence of Coulomb interactions, has the form19,43,46

χ0(k,ω)=2d3q(2π)3f0qf0q+kω+ϵ(q)ϵ(q+k)+ıη,
(1)

with ϵ(q) = 2q2/(2me) being the electron kinetic energy, η → 0+ owing to the adiabatic switching of the perturbation at t → − that ensures causality, and f0(q) being the Fermi–Dirac distribution function that is given by

f0q=1exp2q22mekBTμ̄+1,
(2)

with μ̄=βμ being the reduced chemical potential (μ is the chemical potential and 1/β = kBT) that is determined by the normalization condition [d3q/(2π)3]f0q=n/2.

Finally, we shall mainly focus on the high-degeneracy moderate-density range that lies beyond the warm dense matter regime15,16,19 but simultaneously prior to the Wigner crystallization.13,47 Roughly demarcating the upper warm dense matter boundary with rs ∼ 10,16,19 the above implies that we are interested in rs ∼ 10–100 and θ ∼ 1. In this range, the UEF behaves as an electron liquid whose properties are determined by the interplay between strong Coulomb correlations and quantum effects.

The dielectric formalism is based on the following: (i) the quantum fluctuation dissipation theorem that connects the dynamic structure factor (DSF) S(k, ω) with the imaginary part of the exact density–density response function χ(k, ω); (ii) the general relation that expresses χ(k, ω) in terms of the ideal (Lindhard) density response χ0(k, ω), which introduces the unknown dynamic local field correction (LFC) G(k, ω); and (iii) a functional relation connecting G(k, ω) with the static structure factor (SSF) S(k) that is obtained by perturbative quantum/classical BBGKY approaches or non-perturbative integral equation approaches.19,43 The first two building blocks are exact and common to all the dielectric schemes, whereas the third building block is approximate and differs.

The quantum fluctuation dissipation theorem has the general form46,48

S(k,ω)=π11eβωI{χ(k,ω)}.

Integration over the frequency domain together with the zero frequency moment rule S(k) = ∫S(k, ω) leads to S(k)=[/(2πn)]coth(βω)/2I{χ(k,ω)}dω.19,43 The extension of the χ(k, ω) domain with the aid of analytic continuation leads to the complex valued χ̃(k,z) and allows for integral computation with contour integration techniques.49 Owing to the infinitely many poles of the integrand, the integral is converted to an infinite sum, i.e.,

S(k)=1nβl=χ̃(k,ıωl),
(3)

where ωl = 2πl/(βℏ) are the so-called bosonic Matsubara frequencies.19,49

In addition, within the linear response theory, the polarization potential approach reveals that χ(k, ω) can always be expressed in terms of χ0(k, ω) and the unknown LFC as50–52 

χ(k,ω)=χ0(k,ω)1U(k)1G(k,ω)χ0(k,ω),
(4)

with U(k) = 4πe2/k2 being the regularized Fourier transform of the Coulomb pair interaction energy.

Furthermore, all dielectric formalism schemes approximate the LFC as a SSF functional, leading to closed self-consistent approaches.19,43 Most rigorous approaches treat quantum effects on the random phase approximation (RPA) level and correlation effects in a classical manner, such as the STLS scheme,53,54 VS scheme,55,56 CA scheme,57,58 and HNC scheme.34,35 These lead to a frequency independent LFC and closures of the form G(k) ≡ F[S(k)]. The same applies for the recently developed semi-empirical ESA scheme that utilizes established asymptotic limits and incorporates exact MC simulation results for the warm dense matter regime in order to construct an accurate frequency independent LFC.59,60 The qSTLS scheme constitutes a notable exception since it captures beyond-RPA quantum effects by truncating the quantum BBGKY hierarchy within the Wigner representation at its first member with the introduction of the standard STLS closure condition.61,62 This yields a dynamic LFC and a G(k, ω) ≡ F[S(k), ω] closure.

Combining the above expressions, regardless of the approximation scheme, the dielectric formalism ultimately leads to a functional equation of the form

S(k)=l=χ̃0(k,ıωl)/(nβ)1U(k)1F[S,ıωl]χ̃0(k,ıωl),
(5)

which can be solved in an iterative manner.

For one-component pair-interacting classical systems, the integral equation theory (IET) of liquids consists of two formally exact equations, the Ornstein–Zernike (OZ) integral equation and the non-linear closure equation,36–38 

h(r)=c(r)+nc(r)h(|rr|)d3r,
(6)
g(r)=expβu(r)+h(r)c(r)+b(r),
(7)

where g(r) is the radial distribution function or pair correlation function, h(r) = g(r) − 1 is the total correlation function, c(r) is the direct correlation function, and b(r) is the bridge function. It is noted that capital notations H(k), C(k), B(k) are reserved for the respective Fourier transforms. Fourier transformed radial distribution functions do not emerge in what follows and, thus, should not be confused with the frequency independent LFC G(k).

The bridge function is formally defined by a virial-type expansion whose unknown coefficients are given by multi-dimensional integrals with the integrands being products of Mayer functions f(r) =  exp[ −βu(r)] − 1. The operation of topological reduction allows us to recast the f-bond expansion as an h-bond expansion,36,37 a convenient resummation for long range interactions. However, slow convergence implies that rigorous computation of the exact b[h] functional is impossible. Hence, most theoretical attempts have focused on the formulation of phenomenological b[h] closures, typically of the form b(hc).63 The hypernetted-chain (HNC) approximation, which assumes that b[h] ≡ 0, serves as a prominent example. For the classical OCP, the HNC yields near-exact structural and thermodynamic results for Γ ≲ 1,64 but its accuracy strongly degrades as crystallization is approached.65–67 

A computationally costly but far more accurate alternative concerns the indirect bridge function extraction from computer simulations.68–70 Extraction efforts for the classical OCP date back to the 1980s,71–73 but their deficiencies have been documented.39,74,75 A very recent study managed to obtain very accurate OCP bridge functions for 17 state points Γ = 10–170 and construct an accurate analytic representation valid along the short range and intermediate range.39,75 In reduced units x = r/d, this parameterization reads as

b(x,Γ)=1f(x)bS(x,Γ)+f(x)bI(x,Γ),bS(x,Γ)=s0(Γ)+i=25si(Γ)xi,
(8)
bI(x,Γ)=l0(Γ)Γ5/6expl1(Γ)(x1.44)0.3x2×cosl2(Γ)(x1.44)+l3(Γ)exp3.5(x1.44),f(x)=0.51+erf5.0x1.5,

where bS(x, Γ) is the short range monotonic bridge function, bI(x, Γ) is the intermediate range oscillatory decaying bridge function, f(x) is a sigmoid switching function, si(Γ)=j=03sijΓ(lnΓ)j,li(Γ)=j=04lijΓ1/6(lnΓ)j are monotonic functions of Γ, and the sij,lij coefficients are tabulated.39,75 The closed system of Eqs. (6)(8) leads to a near-exact description of strong correlations in the classical OCP.

The derivation of the functional closure of the IET-based scheme is based on the technique originally developed by Tanaka for the HNC-based scheme.35 The difference is that the IET non-linear closure equation for an arbitrary known b[h] ≡ b(r) function is used instead of the HNC non-linear closure for b(r) ≡ 0. In this section, the basic assumptions and mathematical steps will be outlined in sufficient detail for this article to remain self-contained.

The starting point is the classical (θ) fluctuation dissipation theorem, which reads as37,48

S(k,ω)=1πnβI{χ(k,ω)}ω.

The integration over the frequency domain, together with the zero frequency moment rule S(k) = ∫S(k, ω) and the χ(k,0)=π1[I{χ(k,ω)}/ω]dω Kramers–Kronig causality relation, leads to76,77

S(k)=1nβχ(k,0).

Substituting for the general form of the density–density response function χ(k, ω), approximating the LFC with a frequency independent value G(k, ω) ≡ G(k), employing the Maxwellian result for the static limit of the ideal classical density–density response χ0(k, 0) = −, using the connecting relation S(k) = 1 + nH(k), and invoking the Fourier transformed OZ equation H(k) = C(k)/[1 − nC(k)], one obtains77–79 

βU(k)G(k)=C(k)+βU(k).

Having established the above linear GF[C] functional, the objective is to exploit the two exact IET equations to derive an integral GF[S] functional. Application of the gradient operator to the non-linear closure equation and the OZ equation [see Eqs. (6) and (7)] yields

rh(r)=h(r)+1βru(r)+rb(r)+nd3rrc(rr)h(r),
rh(r)=rc(r)+nd3rrc(rr)h(r).

Equating the right-hand sides of the above equations, two of the convolution terms cancel out each other, resulting in

rc(r)+βru(r)=rβu(r)b(r)h(r)+rb(r)+nh(r)d3rrc(rr)h(r),

whose Fourier transform (after using the multiplication, convolution, and differentiation properties) reads as

kC(k)+βU(k)=kB(k)+d3q(2π)3qH(kq)×B(q)βU(q)1+nH(q)+nH(q)C(q)+βU(q).

We proceed with substituting the linear GF[C] functional relation on both sides, operating with (k·) on both sides, and solving for G(k), which leads to

G(k)=B(k)βU(k)+d3q(2π)3kqk2U(q)U(k)H(kq)×B(q)βU(q)1+nH(q)+nH(q)G(q).

Using again the connecting relation S(k) = 1 + nH(k) in order to dispose of the Fourier transformed total correlation function and substituting for U(k) = 4πe2/k2 ultimately result in the sought-for IET functional

GIET(k)=B(k)βU(k)1nd3q(2π)3kqq2[S(kq)1]×1B(q)βU(q)G(q)1S(q)1.
(9)

As expected, the IET functional collapses to the HNC functional for b(r) ≡ 0, which explicitly reads as34,35

GHNC(k)=1nd3q(2π)3kqq2[S(kq)1]×1G(q)1S(q)1.
(10)

It is rather straightforward to express the IET functional in terms of the HNC functional,

GIET(k)=1nd3q(2π)3kqq2B(q)βU(q)[S(kq)1]+B(k)βU(k)+GHNC(k).

Considering the short range of the OCP bridge function and the long range of Coulomb interactions, this relation suggests that the IET scheme results in small corrections to the LFC of the HNC scheme at the long and intermediate wavelength regions.

Considering the isotropy of homogeneous electron liquids, the IET closure functional GIET(k) is formally expressed as a triple integral of the type ∫d3q(k · q)I(|kq|)J(|q|); see Eq. (9). Such a triple integral type is also encountered in the course of the Wertheim–Thiele derivation of the exact solution of the Percus–Yevick approximation for hard spheres,80,81 the Laplace transform based derivation of the exact solution of the soft mean spherical approximation for classical plasmas,82,83 and the mathematical treatment of the asymptotic memory kernel in mode coupling theories of classical supercooled liquids.84–86 This triple integral can be directly converted to a double integral with the introduction of azimuthally expanded two-center bipolar coordinates, which is equivalent to the sequential transformations q = p + (1/2)k, spherical coordinates for p assuming kẑ without loss of generality, and u = |p + (1/2)k|, w = |p − (1/2)k| with a surface element identity p2 sin θdθdp = (uw/k)dudw. Ultimately, the IET closure functional becomes

G(k)=B(k)βU(k)1n18π2k0B(u)βU(u)+1G(u)1×S(u)1udu|uk|u+ku2w2+k2u2wS(w)1dw.

Utilizing dimensionless variables for all the wave-vectors, i.e., yu/kF, zw/kF, and xk/kF, and rearranging the integrands, the double integral expression for the IET closure functional reads as

G(x)=B(x)βU(x)+38x0B(y)βU(y)+1G(y)1S(y)1×dyy|yx|y+xz2y2x2zS(z)1dz.

To our knowledge, bipolar coordinates have never been utilized in earlier applications of the CA scheme57,58 and the HNC scheme,34,35 whose LFC closure functional is also formally expressed as a triple integral of the type ∫d3q(k · q)I(|kq|)J(|q|), although this would lead to a drastic reduction in the computational cost without invoking extra approximations. This comes in stark contrast to the so-called modified versions of the CA scheme and the HNC scheme,35,87 where the triple integral is converted to a single integral after replacing S(kq) by an ad hoc Yukawa screening function of a characteristic screening wave-number that is determined by an interaction energy constraint.

Independent of the dielectric scheme (STLS, CA, qSTLS, HNC, and IET), the Matsubara summation of Eq. (5) is slowly converging, especially for small values of the degeneracy parameter θ. The implementation of mathematical tricks that speeds up the convergence rate is especially important for computationally heavy schemes, such as the HNC, IET, and qSTLS.

A significant speed-up of the rate of convergence can be achieved by splitting the Hartree–Fock SSF, i.e., the SSF in the absence of Coulomb pair interactions U(k) = 0, since its respective Matsubara summation can be calculated exactly.49,54,62,88 Introducing the auxiliary complex function Φ(k,z)=(2EF)/(3n)χ̃0(k,z) and normalized wave-vectors xk/kF, Eq. (5) becomes

S(x)=SHF(x)l=6πλrsθ1x21G(x)Φ2(x,l)1+4πλrs1x21G(x)Φ(x,l),

where, courtesy of logarithm properties and the product n=+a2+π2n2/b2+π2n2=sinh2(a)/sinh2(b), the Hartree–Fock SSF SHF(x)=(3/2)θl=Φ(x,l) is given by the following integral:49,54,88

SHF(x)=13θ4x0ydyln1+expμ̄(yx)2θ1+expμ̄(y+x)2θexpy2θμ̄+1.

Convergence can be further accelerated by splitting the square of the high frequency–short wavelength asymptotic form (l, q) of the auxiliary complex function Φ(x,l)=(4/3)x2/[x4+(2πlθ)2]+O(x4,l4) since its respective Matsubara summation can also be calculated exactly.49,54,62 This leads to

S(x)=SHF(x)S(x)6πλrsθ1x21G(x)×l=Φ2(x,l)1+4πλrs1x21G(x)Φ(x,l)Φ2(x,l),

where, by differentiating both sides of the known formula n=+(x4+a2n2)1=[π/(ax2)]cothπx2/a with respect to x, the residual SSF correction that is defined by S(x)=(6/π)λrsθ(1/x2)1G(x)l=Φ2(x,l) is given by49,54

S(x)=43πλrsθ1G(x)x2csch2x22θ+2θx2cothx22θ,

where coth(·) is the hyperbolic cotangent and csch(·) is the hyperbolic cosecant.

On the practical side, with an accuracy goal of 0.001% in the interaction energies, the above tricks speed up the convergence rate by more than two orders of magnitude. To be more specific, for rs = 100 and for θ = 0.5 or 4, this accuracy goal is well achieved within l = 512 when these two splitting procedures are implemented. On the other hand, for rs = 100, the interaction energy accuracy at l = 51 200 is merely 0.1% for θ = 4 and 0.7% for θ = 0.5 when no splitting is implemented.

From the general SSF decomposition in terms of SHF(x) and S(x), it is rather straightforward to prove that the frequency independent LFC assumption G(k, ω) ≡ G(k) implies the (3π/8)x4×[1 − S(x)] = λrs[1 − G(x)] asymptotic limit.49 When combined with Kimball’s expression89g(r̃=0)/r̃=(3π/8)limxx4[1S(x)] with r̃=rkF and the cusp relation89,90g(r̃=0)/r̃=λrsg(0), this yields the general asymptotic condition G(x) = 1 − g(0),49,91 where g(0) is the contact or on-top value of the radial distribution function. This asymptotic condition is often considered as a self-consistency condition since it solely originates from the first two building blocks of the dielectric formalism and needs to be independently satisfied by the third building block (closure functional). This has been shown to be valid in the case of the STLS scheme49 and can also be confirmed for other dielectric schemes (VS, CA, and HNC) and for the IET scheme, given the rapid B(x)/βU(x) decay to zero. Note that in the strongly coupled electron liquid regime, g(0) ≃ 0 is expected,59,92 leading to G(x) ≃ 1.

Nevertheless, the implicit form of the IET LFC can strongly inhibit its proper numerical convergence at short wavelengths. In particular, the converged LFC solution has been consistently observed to exhibit a short wavelength dependence, which abruptly drops from near-unity to near-zero. The unphysical asymptotic behavior always takes place close to the wavenumber cutoff considered in the integrations, regardless of its actual value. It translates to interaction energy errors of the order of 0.005%, which exceed our 0.001% accuracy goal.

In order to achieve proper convergence of the IET local field correction in the short wavelength limit, it is highly beneficial to split the STLS LFC from the IET LFC. The STLS LFC emerges by setting B(k) ≡ 0 and G(k) ≡ 0 in the right-hand side of the IET LFC. The standard single integral form of the STLS LFC is recovered from the double integral form of the IET LFC by employing the change in variables (y, z) → (t, s) of the form y=x2+s22xst, z = s and carrying out the t− integration. This leads to

G(x)=G1(x)+G2(x),G1(x)=340s2S(s)11+x2s22xslnx+sxsds,G2(x)=B(x)βU(x)38x0B(y)βU(y)+G(y)1S(y)1×dyy|yx|y+xz2y2x2zS(z)1dz.

In the above equation, G1(x) denotes the STLS LFC with a numerical asymptotic behavior G1(x) ≃ 1 and G2(x) denotes the residual IET LFC with a numerical asymptotic behavior G2(x) = 0. Thus, the IET LFC short wavelength limit now correctly converges to a value close to unity.

The closed normalized set of equations for the IET-based dielectric scheme comprises of the following: the normalization condition of the Fermi–Dirac energy distribution function [see Eq. (11), Sec. II A], the ideal Lindhard density response expressed through the auxiliary complex function Φ(x, l) and evaluated at the imaginary Matsubara frequencies ωl including the static limit [see Eqs. (12) and (13), Sec. II A], the Fourier transform of the classical OCP bridge function [see Eq. (14), Sec. II C], the infinite Matsubara summation expression for the SSF after separating the Hartree–Fock SSF and also the asymptotic component [see Eq. (15), Sec. III B], and the IET LFC double integral expression after utilizing bipolar coordinates and splitting the STLS LFC single integral expression [see Eq. (16), Secs. III A and III C],

0zdzexpzμ̄+1=23θ3/2,
(11)
Φ(x,l)=12x0yexpy2θμ̄+1lnx2+2xy2+2πlθ2x22xy2+2πlθ2dy,
(12)
Φ(x,0)=1θx0yexpy2θμ̄expy2θμ̄+12×y2x24ln2y+x2yx+xydy,
(13)
B(q)βU(q)=9π8θrsq0yby,2λ2rsθsinqλydy,
(14)
S(x)=SHF(x)S(x)6πλrsθ1x21G(x)×l=Φ2(x,l)1+4πλrs1x21G(x)Φ(x,l)Φ2(x,l),
(15)
G(x)=340s2S(s)11+x2s22xslnx+sxsds+B(x)βU(x)38x0B(y)βU(y)+G(y)1S(y)1×dyy|yx|y+xz2y2x2zS(z)1dz.
(16)

Concerning the origin of Eq. (14), within the classical x = r/d or q = kd normalization, it is rather straightforward to show that B(q)/βU(q)=(q/Γ)0xb(x,Γ)sin(qx)dx for the ratio of spatial Fourier transforms, where b(x, Γ) is directly adopted from Eq. (8). Naturally, the above expression needs to be translated to the quantum x = rkF or q = k/kF normalization, which is formally equivalent to the substitution qq/λ. Finally, it is apparent that the utilization of classical OCP bridge functions necessitates a mapping of the quantum states (rs, θ) to classical states (Γ) via Γ = 2λ2(rs/θ). The B(q)/βU(q) contribution is illustrated in Fig. 1.

FIG. 1.

The non-trivial −B(q)/βU(q) contribution in the normalized wavenumber range q = k/kF ≤ 4.5 for different values of the classical coupling parameter Γ = 2λ2(rs/θ). Results for Γ = 10–220. It should be noted that the classical OCP bridge function extrapolates smoothly in the range Γ = 170–220. For the computation, Eq. (14) was combined with Eq. (8). The numerical integration was performed with the doubly adaptive Clenshaw–Curtis quadrature method.

FIG. 1.

The non-trivial −B(q)/βU(q) contribution in the normalized wavenumber range q = k/kF ≤ 4.5 for different values of the classical coupling parameter Γ = 2λ2(rs/θ). Results for Γ = 10–220. It should be noted that the classical OCP bridge function extrapolates smoothly in the range Γ = 170–220. For the computation, Eq. (14) was combined with Eq. (8). The numerical integration was performed with the doubly adaptive Clenshaw–Curtis quadrature method.

Close modal

The accuracy goal was set to 0.001% with respect to the interaction energies. All improper integrals were numerically evaluated with the doubly adaptive Clenshaw–Curtis quadrature method, as implemented in the GNU Scientific Library, with a 0.1kF grid resolution and with a 40kF upper cutoff. The only exception was the complete Fermi–Dirac integral I1/2(μ̄), which is implemented as a special function in the GNU Scientific Library. It should be noted that non-adaptive quadrature rules were confirmed to require much denser grid spacing <0.005kF in order to satisfy the accuracy goal. The infinite Matsubara summation was truncated at |l| = 512. Convergence studies were carried out, for representative quantum coupling parameters and degeneracy parameters, in order to ensure that the chosen resolution and cutoffs do not affect the thermodynamic and structural results.

The iteration cycle proceeds as follows: (i) The reduced chemical potential μ̄ is calculated from Eq. (11) using a bisection root-finding algorithm. (ii) The Lindhard density responses (0 ≤ l ≤ 512) are computed from Eqs. (12) and (13). (iii) The Fourier transform of the classical OCP bridge function is computed from Eq. (14) and stored. (iv) With the RPA’s LFC as an initial guess, the SSF is evaluated from Eq. (15). (v) These LFC and SSF are substituted in the right-hand side of Eq. (16) for an initial evaluation of the IET LFC. (vi) The last two steps are repeated until the absolute relative difference between two successive IET LFC evaluations is smaller than 10−5 for all the grid points. Starting from the RPA solution, convergence is typically reached within 200 iterations. Especially for state points with rs > 100 and θ < 1.0, Broyles’s technique of mixing iterates65,93 was necessary to speed up and sometimes even to achieve convergence.

To compute accurate benchmark data for our new IET scheme, we have carried out direct PIMC simulations94–96 of the UEF without any nodal restrictions.97 As a consequence, our simulations are afflicted with the notorious fermion sign problem,98 which, in general, leads to an exponential increase in the compute time by increasing the system-size N or decreasing the temperature T; the reader is addressed to Refs. 99 and 100 for topical and accessible review articles. In practice, however, the sign problem is not severe as quantum exchange effects are effectively reduced by the strong Coulomb repulsion in the strongly coupled electron liquid regime.

The basic idea behind the PIMC method is to evaluate the (canonical, i.e., system-size N, volume V, and inverse temperature β = 1/kBT are fixed) partition function in coordinate space, which gives

Zβ,N,V=1N!N!σSNσSNsgnσ,σ×dRReβĤπ̂σπ̂σR.
(17)

Specifically, R=(r1,,rN)T contains the coordinates of all N = N + N electrons, and due to the antisymmetry of the fermionic density matrix under the exchange of particle coordinates, we have to explicitly take the sums over all elements σi of the respective permutation group SNi, i ∈ {↑, ↓}, where the sign function sgn(σ, σ) gives positive (negative) unity for an even (odd) number of pair permutations. Furthermore, the operators π̂σi realize the particular permutations for a corresponding element σi.

The problem with Eq. (17) is that the matrix elements of the density operator eβĤ cannot be readily evaluated since the kinetic (K̂) and potential (V̂) contributions to the Hamiltonian Ĥ=K̂+V̂ do not commute, eβĤeβK̂eβV̂. In order to overcome this obstacle, we utilize the exact semi-group property of ρ̂,

eβĤ=α=0P1eϵĤ,
(18)

where the definition ϵ = β/P has been employed. Applying Eq. (18) to Eq. (17) and inserting P − 1 unity operators of the form 1̂=dRα|RαRα| lead to the intermediate result

Zβ,N,V=1N!N!σSNσSNsgnσ,σ×dR0dRP1R0|eϵĤ|R0×R1|eϵĤ|R1RP1|eϵĤ|π̂σπ̂σR0,
(19)

which is still exact. Evidently, Eq. (19) requires the evaluation of P density matrices, but at P times the temperature. For a sufficiently large P, each of these factors can be straightforwardly evaluated using a suitable high-temperature approximation, such as the primitive factorization eϵĤeϵK̂eϵV̂. In fact, the factorization error of the latter decays as ∼1/P2,101 and the convergence in the limit of P is ensured by the celebrated Trotter formula,102 

limPeϵK̂eϵV̂P=eβ(K̂+V̂).
(20)

In practice, we find P = 200 sufficient to reduce the factorization error substantially below the noise level of the Monte Carlo simulation. For completeness, we note that Eq. (20) only holds for the case of potentials V̂ that are bounded from below, which is, indeed, the case for the UEF. Attractive potentials such as the Coulomb interaction between positive and negative charges require a modified procedure such as the pair approximation, which is discussed in detail, for instance, in Ref. 103. In addition, we note that higher-order factorizations of the thermal density matrix that converge as ∼1/P4104,105 and even ∼1/P6−8106 have been presented, although we do not find them necessary for the present conditions.

The final result for the PIMC partition function can then be written in the abbreviated form as

Z=dXW(X),
(21)

where the integration over the X=(R0,,RP1)T meta-variable also contains the summation over all the possible permutations. Furthermore, the weight function W(X) can be readily evaluated for each individual configuration X and contains contributions from the potential energy and from the free thermal density matrix; see, for instance, Ref. 94 for an accessible review article. Evidently, Eq. (21) requires the evaluation of a d = 3PN-dimensional integral, which easily leads to d ∼ 104 in the present study. While this renders the application of standard quadrature methods impractical due to the well-known curse of dimensionality, the Metropolis Monte Carlo method107 does not suffer from this drawback.

A graphical illustration of the PIMC approach is shown in Fig. 2. Specifically, we show a configuration X of N = 4 electrons with P = 6 high-temperature factors. It is evident that each particle is represented by a closed (i.e., β-periodic) path of particle coordinates in the imaginary time τ ∈ [0, β] (strictly speaking, it is τ/i[0, β], but it is conventional to drop the pre-factor for simplicity). This, in turn, corresponds to the famous classical isomorphism,108 where the complicated quantum many-body system of interest is mapped onto an effective classical system of interacting ring-polymers.

FIG. 2.

Illustration of the PIMC method. A configuration of N = 4 particles is shown in the τx plane for P = 6 high-temperature factors. Note the single pair exchange of the two particles in the center, leading to a negative sign for this particular configuration, W(X) < 0.

FIG. 2.

Illustration of the PIMC method. A configuration of N = 4 particles is shown in the τx plane for P = 6 high-temperature factors. Note the single pair exchange of the two particles in the center, leading to a negative sign for this particular configuration, W(X) < 0.

Close modal

An additional difficulty arises due to the indistinguishable nature of electrons of the same spin-species, which requires us to sample all the possible permutations of particle coordinates. Within the PIMC picture, the latter manifest as the so-called permutation cycles,109 which are trajectories comprising more than a single particle. An example for such a permutation cycle can be identified at the center of Fig. 2. As a consequence of the depicted pair exchange, we have to move twice through the imaginary time to return to the point of origin. This, in turn, results in a negative sign of the corresponding configuration weight, W(X) < 0, and thereby contributes to the aforementioned fermion sign problem.

In practice, we construct a Markov chain of random configurations {Xi} using a canonical adaption110 of the worm algorithm by Boninsegni et al.111,112

Let us next consider the PIMC estimation of the interaction energy per particle VN/N, which is shown in Fig. 3 for rs = 125 and θ = 2. In particular, PIMC simulations are, by default, only possible for a finite number of electrons N. The raw PIMC data for three different system sizes are shown as the red circles in Fig. 3 and exhibit a significant dependence on N. In practice, however, we are interested in the thermodynamic limit, i.e., in the limit of an infinite number of particles with the density being constant,

v=limNVNNrs,θ.
(22)

To eliminate the difference between the PIMC data for VN/N and Eq. (22), we use the finite-size correction by Chiesa et al.,113 which has subsequently been adapted to finite temperatures in Ref. 114,

Δv(N)=ωp4Ncothβωp2,
(23)

where the plasma frequency is given by ωp=3/rs3 in Hartree atomic units. To be more specific, Eq. (23) is based on the insight that the system-size dependence of VN/N is mainly the consequence of the approximation of a continuous integral [see Eq. (25)] by the sum over reciprocal lattice vectors due to the momentum quantization in a finite simulation cell. To the first order, this discretization error can be approximated by utilizing the exact long wave-length limit of the static structure factor115 

limk0S(k)=k22ωpcothβωp2.
(24)

A detailed derivation is beyond the scope of the present work and was presented by Drummond et al.116 Adding the finite-size correction given in Eq. (23) to the PIMC data leads to the green crosses in Fig. 3. Evidently, the bulk of the finite-size errors have been removed, and the corrected data points fall into an interval of 0.01% (horizontal light gray lines) around their common average value (horizontal blue line). All the PIMC results for the interaction energy that are shown in this work have been obtained by following this procedure.

FIG. 3.

The finite-size extrapolation of the PIMC data. The UEF interaction energy per particle is shown at rs = 125 and θ = 2 as a function of the inverse electron number. The red circles depict the raw PIMC results for different numbers of electrons N, and the green crosses are obtained by adding onto the former the finite-size correction from Eq. (23). The horizontal blue dotted line shows the average value of the corrected data points, and the light gray dotted lines indicate a deviation of 0.01%.

FIG. 3.

The finite-size extrapolation of the PIMC data. The UEF interaction energy per particle is shown at rs = 125 and θ = 2 as a function of the inverse electron number. The red circles depict the raw PIMC results for different numbers of electrons N, and the green crosses are obtained by adding onto the former the finite-size correction from Eq. (23). The horizontal blue dotted line shows the average value of the corrected data points, and the light gray dotted lines indicate a deviation of 0.01%.

Close modal

For completeness, it should be mentioned that finite-size effects are substantially more pronounced in the warm dense matter regime and that the simple first-order correction from Eq. (23) breaks down.19,117 A recent investigation of finite-size effects of a uniform electron gas at extreme densities and temperatures was presented by Dornheim and Vorberger.118 

Here, we compare the paramagnetic electron liquid interaction energies and static properties as computed from different dielectric formalism schemes with their “exact” counterparts as extracted from our PIMC simulations. It has been demonstrated that beyond the warm dense matter regime and especially for rs > 30, the STLS and the VS schemes yield increasingly inaccurate results.34 Only the IET, HNC, and qSTLS schemes will be numerically solved herein since the inclusion of the STLS and VS schemes has been judged to be meaningless. The HNC and IET comparison will lead to direct conclusions regarding the impact of the classical OCP bridge function inclusion, while the qSTLS and IET comparison will lead to indirect conclusions regarding the significance of the beyond-RPA quantum effects. Our HNC algorithm can be essentially obtained from the IET algorithm described in Sec. III D by setting B(q)/[βU(q)] ≡ 0 and has been systematically validated against the published results of Tanaka.35 On the other hand, our qSTLS algorithm has marked differences from the IET algorithm described in Sec. III D owing to the dynamic nature of the LFC and has been benchmarked against the published results of Schweng and Böhm.62 

The interaction energy per particle is generally obtained by U=(1/2)[d3k/(2π)3]U(k)S(k)1. Substituting for the Coulomb potential energy, introducing the (rs, θ) variables, and employing normalized wavenumbers lead to the standard UEF expression19,43

ũ(rs,θ)=1πλrs0S(x)1dx,
(25)

where ũ(rs,θ) denotes the interaction energy per particle normalized by the Hartree energy Eh = e2/aB.

The normalized interaction energies extracted from the PIMC simulations and computed with the qSTLS, HNC, and IET schemes are listed in Table I. The absolute relative deviations between the theoretical and “exact” interaction energies are also tabulated therein. (i) The qSTLS interaction energies are relatively inaccurate with 5.0%–8.7% relative deviations from PIMC results. The errors systematically decrease as θ increases and increase as rs increases. (ii) The HNC interaction energies are very accurate with 0.47%–1.37% relative deviations from PIMC results. The errors systematically decrease as θ increases (diminishing quantum effects) and increase as rs increases (stronger beyond-HNC classical pair correlations65–67). (iii) The IET interaction energies are revealed to be the most accurate with merely 0.05%–0.68% relative deviations from the PIMC results. No systematic tendencies have been observed in the errors with respect to either θ or rs. (iv) The IET scheme provides the most accurate interaction energy predictions for all 20 thermodynamic state points investigated, with appreciable improvements over the HNC interaction energy predictions. (v) Despite their very high accuracy, the IET (as well as the HNC) interaction energies cannot reproduce the exact θ-dependence of ũ(rs,θ) within the highly degenerate θ ≲ 1 range. To be more specific, the PIMC results reveal a monotonic |ũ| decrease as θ increases, while the IET (and HNC) results exhibit a monotonic |ũ| increase as θ increases within θ ≲ 1 and the correct monotonic |ũ| decrease only as θ increases within θ ≳ 1. For this reason, we did not construct an exchange–correlation free energy parameterization via the adiabatic connection formula.

TABLE I.

The interaction energy ũ(rs,θ) per particle (expressed in Hartree units) of the paramagnetic electron liquid: a comparison of the finite-size corrected “exact” PIMC results with the predictions of the qSTLS, HNC, and IET dielectric formalism schemes. The PIMC data for the first six state points are adapted from Ref. 34, while the PIMC data for the remaining 14 state points are new. The absolute relative deviations between the dielectric scheme results and the PIMC results are also reported.

rsθΓũPIMCũqSTLSeqSTLS (%)ũHNCeHNC (%)ũIETeIET (%)
100     0.50     108.6     −0.008 255 00     −0.007 623 82 7.646 −0.008 158 66 1.167 −0.008 221 81 0.402 
100     0.75     72.40     −0.008 245 70     −0.007 647 55 7.254 −0.008 164 90 0.980 −0.008 225 44 0.246 
100     1.00     54.30     −0.008 234 90     −0.007 667 12 6.895 −0.008 166 18 0.834 −0.008 225 59 0.113 
100     2.00     27.15     −0.008 176 50     −0.007 688 92 5.963 −0.008 129 05 0.580 −0.008 190 66 0.173 
100     4.00     13.58     −0.008 006 23     −0.007 603 08 5.035 −0.007 968 33 0.473 −0.008 031 43 0.315 
50     0.50     54.30     −0.016 007 00     −0.014 998 07 6.303 −0.015 898 41 0.678 −0.016 035 10 0.176 
60     0.50     65.16     −0.013 453 10     −0.012 560 41 6.636 −0.013 348 04 0.781 −0.013 460 14 0.052 
70     0.50     76.02     −0.011 611 75     −0.010 807 32 6.928 −0.011 509 38 0.882 −0.011 603 90 0.068 
80     0.50     86.88     −0.010 219 37     −0.009 485 45 7.182 −0.010 120 12 0.971 −0.010 201 49 0.175 
90     0.50     97.74     −0.009 128 62     −0.008 452 87 7.403 −0.009 032 93 1.048 −0.009 104 15 0.268 
110     0.50     119.5     −0.007 526 42     −0.006 943 43 7.746 −0.007 440 12 1.147 −0.007 496 75 0.394 
125     0.50     135.8     −0.006 654 21     −0.006 124 35 7.963 −0.006 573 77 1.209 −0.006 622 68 0.474 
125     0.75     90.50     −0.006 650 53     −0.006 143 35 7.626 −0.006 578 38 1.085 −0.006 625 56 0.442 
125     1.00     67.88     −0.006 643 36     −0.006 159 66 7.281 −0.006 579 99 0.954 −0.006 626 47 0.254 
125     1.50     45.25     −0.006 625 35     −0.006 179 40 6.731 −0.006 574 32 0.770 −0.006 621 12 0.064 
125     2.00     33.94     −0.006 602 98     −0.006 184 76 6.334 −0.006 559 00 0.666 −0.006 607 12 0.063 
150     0.50     162.9     −0.005 581 77     −0.005 119 15 8.288 −0.005 508 21 1.318 −0.005 547 97 0.606 
150     1.00     81.45     −0.005 571 34     −0.005 148 88 7.583 −0.005 513 37 1.040 −0.005 551 32 0.359 
200     0.50     217.2     −0.004 222 44     −0.003 855 65 8.687 −0.004 164 45 1.373 −0.004 193 73 0.680 
200     1.00     108.6     −0.004 217 10     −0.003 878 21 8.036 −0.004 168 13 1.161 −0.004 195 59 0.510 
rsθΓũPIMCũqSTLSeqSTLS (%)ũHNCeHNC (%)ũIETeIET (%)
100     0.50     108.6     −0.008 255 00     −0.007 623 82 7.646 −0.008 158 66 1.167 −0.008 221 81 0.402 
100     0.75     72.40     −0.008 245 70     −0.007 647 55 7.254 −0.008 164 90 0.980 −0.008 225 44 0.246 
100     1.00     54.30     −0.008 234 90     −0.007 667 12 6.895 −0.008 166 18 0.834 −0.008 225 59 0.113 
100     2.00     27.15     −0.008 176 50     −0.007 688 92 5.963 −0.008 129 05 0.580 −0.008 190 66 0.173 
100     4.00     13.58     −0.008 006 23     −0.007 603 08 5.035 −0.007 968 33 0.473 −0.008 031 43 0.315 
50     0.50     54.30     −0.016 007 00     −0.014 998 07 6.303 −0.015 898 41 0.678 −0.016 035 10 0.176 
60     0.50     65.16     −0.013 453 10     −0.012 560 41 6.636 −0.013 348 04 0.781 −0.013 460 14 0.052 
70     0.50     76.02     −0.011 611 75     −0.010 807 32 6.928 −0.011 509 38 0.882 −0.011 603 90 0.068 
80     0.50     86.88     −0.010 219 37     −0.009 485 45 7.182 −0.010 120 12 0.971 −0.010 201 49 0.175 
90     0.50     97.74     −0.009 128 62     −0.008 452 87 7.403 −0.009 032 93 1.048 −0.009 104 15 0.268 
110     0.50     119.5     −0.007 526 42     −0.006 943 43 7.746 −0.007 440 12 1.147 −0.007 496 75 0.394 
125     0.50     135.8     −0.006 654 21     −0.006 124 35 7.963 −0.006 573 77 1.209 −0.006 622 68 0.474 
125     0.75     90.50     −0.006 650 53     −0.006 143 35 7.626 −0.006 578 38 1.085 −0.006 625 56 0.442 
125     1.00     67.88     −0.006 643 36     −0.006 159 66 7.281 −0.006 579 99 0.954 −0.006 626 47 0.254 
125     1.50     45.25     −0.006 625 35     −0.006 179 40 6.731 −0.006 574 32 0.770 −0.006 621 12 0.064 
125     2.00     33.94     −0.006 602 98     −0.006 184 76 6.334 −0.006 559 00 0.666 −0.006 607 12 0.063 
150     0.50     162.9     −0.005 581 77     −0.005 119 15 8.288 −0.005 508 21 1.318 −0.005 547 97 0.606 
150     1.00     81.45     −0.005 571 34     −0.005 148 88 7.583 −0.005 513 37 1.040 −0.005 551 32 0.359 
200     0.50     217.2     −0.004 222 44     −0.003 855 65 8.687 −0.004 164 45 1.373 −0.004 193 73 0.680 
200     1.00     108.6     −0.004 217 10     −0.003 878 21 8.036 −0.004 168 13 1.161 −0.004 195 59 0.510 

It is important to point out that the employed classical OCP bridge function parameterization is strictly valid for 10 ≤ Γ ≤ 170.39,75 This upper threshold is surpassed by the state point (rs, θ) = (200, 0.50), which corresponds to Γ = 217.2. Since all the si(Γ) and li(Γ) coefficients involved in Eq. (8) are monotonic functions of Γ, it can be expected from the well-known continuity between the stable and metastable liquid states that the analytic classical OCP bridge function expression can be extrapolated without large errors toward the supercooled liquid regime Γ > 171.8. The same applies for extrapolations in the weakly interacting regime Γ < 10, but these are less significant due to the diminished bridge function impact.

Furthermore, we compare with the interaction energies as computed from three accurate parameterizations of the exchange–correlation free energy through the expression

ũ(rs,θ)=rsf̃xc(rs,θ)rs+2f̃xc(rs,θ),
(26)

which is acquired after the differentiation of the thermodynamic formula f̃xc(rs,Θ)=rs20rsrsũint(rs,Θ)drs.19 In particular, we consider the parameterization by Groth et al. (GDSMFB) that is based on simulation results obtained by various novel PIMC methods within 0.1 ≤ rs ≤ 20 and 0.5 ≤ θ ≤ 8,119 the parameterization by Karasiev et al. (KSDT) that is based on simulation results obtained by the restricted PIMC method within 1.0 ≤ rs ≤ 40 and 0.062 5 ≤ θ ≤ 8,120 and the corrected version of the parameterization by Karasiev et al. (corrKSDT) that is also based on the above restricted PIMC data.121 Despite some documented deficiencies of the restricted PIMC data input114 (concerning the uncontrolled fixed node approximation122 and the unsatisfactory treatment of finite-size effects117) and despite a procedural mistake in the original fitting procedure (concerning the utilization of an analytic ground-state fit instead of the actual ground-state MC data123), the KSDT interaction energies exhibit 0.39%–1.63% relative deviations from PIMC results, whereas the GDSMFB interaction energies exhibit 1.08%–4.94% relative deviations and the corrKSDT interaction energies exhibit 1.58%–6.90% relative deviations from PIMC results. Therefore, extrapolated GDSMFB, KSDT, and corrKSDT interaction energies are more accurate than the qSTLS results but less accurate than the IET results. More specifically, the KSDT interaction energies are even more accurate than the HNC results at some investigated states but never more accurate than the IET results. Hence, even though the GDSMFB and corrKSDT exchange–correlation free energy parameterization is much more accurate than the KSDT parameterization within the warm dense matter regime,119 it can be concluded that the latter extrapolates better within the strongly coupled electron regime, at least as far as interaction energies are concerned. For completeness, we point out that the GDSMFB, KSDT, and corrKSDT discrepancies in warm dense matter ranges do not impact density functional theory calculations.124 

Finally, we also compare with the interaction energies as computed from the direct parameterization of the interaction energy by Ichimaru et al. (IIT), which is based on the STLS interaction energies after their correction, to comply with variable coupling QMC simulations for the ground-state limit (θ → 0) and with variable coupling MC simulations for the classical limit (θ), via the implementation of an ad hoc θ− interpolation function.43,125,126 Although the interpolation function accuracy is unclear for intermediate degeneracies,19 the incorporation of exact strong coupling ground-state results127 and classical results128,129 hints that the IIT parameterization might be very accurate. This is verified by the comparison that reveals that the IIT interaction energies exhibit 0.06%–0.63% relative deviations from PIMC results. Remarkably, IIT interaction energies are much more accurate than the qSTLS, more accurate than the HNC, and even as accurate as the IET results. In particular, the IET interaction energies are more accurate for 11 and the IIT interaction energies more accurate for nine of the studied states, while the mean absolute relative errors are 0.29% for the IET and 0.35% for the IIT.

Characteristic static structure factors extracted from our PIMC simulations and computed with the three dielectric formalism schemes are illustrated in Figs. 4 and 5. For all the 20 investigated state points, the magnitudes and positions of the first SSF peak resulting from the PIMC simulations and from the qSTLS, HNC, and IET schemes are listed in Table II. The absolute relative deviations between the theoretical and the “exact” SSF peak values are also tabulated therein.

FIG. 4.

Dependence of the paramagnetic electron liquid static structure factor on the quantum coupling parameter rs. The results from the IET scheme (red solid line), the HNC scheme (blue dashed line), the qSTLS scheme (green dotted-dashed line), and PIMC simulations (black crosses) for θ = 0.50 and varying rs = 50, 100, 150, 200 in the normalized wavenumber range k ≤ 4.5kF. The superiority of the new IET scheme within the long wavelength range, in the maximum vicinity, and within the short wavelength range is obvious.

FIG. 4.

Dependence of the paramagnetic electron liquid static structure factor on the quantum coupling parameter rs. The results from the IET scheme (red solid line), the HNC scheme (blue dashed line), the qSTLS scheme (green dotted-dashed line), and PIMC simulations (black crosses) for θ = 0.50 and varying rs = 50, 100, 150, 200 in the normalized wavenumber range k ≤ 4.5kF. The superiority of the new IET scheme within the long wavelength range, in the maximum vicinity, and within the short wavelength range is obvious.

Close modal
FIG. 5.

Dependence of the paramagnetic electron liquid static structure factor on the quantum degeneracy parameter θ. The results from the IET scheme (red solid line), the HNC scheme (blue dashed line), the qSTLS scheme (green dotted-dashed line), and PIMC simulations (black crosses) for rs = 125 and varying θ = 0.5, 1.0, 1.5, 2.0 in the normalized wavenumber range k ≤ 4.5kF. The superiority of the new IET scheme within the long wavelength range, in the maximum vicinity, and within the short wavelength range is obvious.

FIG. 5.

Dependence of the paramagnetic electron liquid static structure factor on the quantum degeneracy parameter θ. The results from the IET scheme (red solid line), the HNC scheme (blue dashed line), the qSTLS scheme (green dotted-dashed line), and PIMC simulations (black crosses) for rs = 125 and varying θ = 0.5, 1.0, 1.5, 2.0 in the normalized wavenumber range k ≤ 4.5kF. The superiority of the new IET scheme within the long wavelength range, in the maximum vicinity, and within the short wavelength range is obvious.

Close modal
TABLE II.

The peak magnitude Smax and peak position argqSmax of the static structure factor S(k/kF; rs, θ) of the paramagnetic electron liquid: a comparison of the “exact” PIMC results with the predictions of the qSTLS, HNC, and IET dielectric formalism schemes. The PIMC data for the first six state points are adopted from Ref. 34, while the PIMC data for the remaining 14 state points are new. The absolute relative deviations between the dielectric scheme results and the PIMC results are also reported.

rsθSPIMCmaxSqSTLSmaxeqSTLS (%)SHNCmaxeHNC (%)SIETmaxeIET (%)argqSPIMCmaxargqSqSTLSmaxeqSTLS (%)argqSHNCmaxeHNC (%)argqSIETmaxeIET (%)
100 0.50 1.494 1.224 18.074 1.202 19.546 1.306 12.598 2.248 1.770 21.247 2.220 1.225 2.200 2.115 
100 0.75 1.449 1.191 17.779 1.203 16.932 1.301 10.159 2.248 1.770 21.247 2.220 1.225 2.210 1.670 
100 1.00 1.406 1.166 17.127 1.199 14.724 1.290 8.285 2.191 1.790 18.288 2.230 1.797 2.220 1.341 
100 2.00 1.285 1.102 14.204 1.164 9.366 1.227 4.461 2.248 1.870 16.798 2.270 0.999 2.250 0.110 
100 4.00 1.166 1.051 9.912 1.105 5.263 1.138 2.381 2.303 2.040 11.422 2.340 1.605 2.320 0.736 
50 0.50 1.252 1.107 11.569 1.093 12.677 1.136 9.211 2.248 1.910 15.018 2.290 1.889 2.250 0.110 
60 0.50 1.327 1.131 14.767 1.115 15.956 1.171 11.802 2.231 1.860 16.635 2.260 1.294 2.230 0.051 
70 0.50 1.351 1.156 14.450 1.137 15.822 1.205 10.861 2.188 1.830 16.355 2.250 2.843 2.220 1.471 
80 0.50 1.391 1.180 15.174 1.159 16.674 1.238 10.990 2.144 1.800 16.029 2.230 4.030 2.210 3.097 
90 0.50 1.438 1.202 16.418 1.181 17.905 1.272 11.551 2.231 1.790 19.772 2.220 0.499 2.210 0.947 
110 0.50 1.481 1.248 15.766 1.223 17.456 1.339 9.591 2.144 1.750 18.362 2.210 3.097 2.200 2.631 
125 0.50 1.538 1.288 16.310 1.254 18.502 1.390 9.661 2.188 1.720 21.383 2.200 0.557 2.200 0.557 
125 0.75 1.557 1.243 20.160 1.256 19.369 1.384 11.132 2.231 1.730 22.461 2.210 0.947 2.200 1.395 
125 1.00 1.501 1.209 19.455 1.251 16.671 1.369 8.765 2.188 1.750 20.011 2.220 1.471 2.210 1.014 
125 1.50 1.436 1.163 19.000 1.231 14.258 1.331 7.338 2.231 1.790 19.772 2.230 0.051 2.220 0.499 
125 2.00 1.365 1.133 17.039 1.209 11.464 1.292 5.360 2.274 1.820 19.952 2.250 1.040 2.230 1.919 
150 0.50 1.695 1.355 20.022 1.305 23.018 1.474 13.022 2.144 1.700 20.694 2.200 2.631 2.200 2.631 
150 1.00 1.586 1.255 20.916 1.300 18.038 1.447 8.768 2.231 1.720 22.909 2.210 0.947 2.210 0.947 
200 0.50 1.817 1.436 20.962 1.403 22.807 1.638 9.833 2.231 1.700 23.806 2.190 1.844 2.200 1.395 
200 1.00 1.738 1.337 23.112 1.395 19.757 1.600 7.943 2.231 1.700 23.806 2.200 1.395 2.200 1.395 
rsθSPIMCmaxSqSTLSmaxeqSTLS (%)SHNCmaxeHNC (%)SIETmaxeIET (%)argqSPIMCmaxargqSqSTLSmaxeqSTLS (%)argqSHNCmaxeHNC (%)argqSIETmaxeIET (%)
100 0.50 1.494 1.224 18.074 1.202 19.546 1.306 12.598 2.248 1.770 21.247 2.220 1.225 2.200 2.115 
100 0.75 1.449 1.191 17.779 1.203 16.932 1.301 10.159 2.248 1.770 21.247 2.220 1.225 2.210 1.670 
100 1.00 1.406 1.166 17.127 1.199 14.724 1.290 8.285 2.191 1.790 18.288 2.230 1.797 2.220 1.341 
100 2.00 1.285 1.102 14.204 1.164 9.366 1.227 4.461 2.248 1.870 16.798 2.270 0.999 2.250 0.110 
100 4.00 1.166 1.051 9.912 1.105 5.263 1.138 2.381 2.303 2.040 11.422 2.340 1.605 2.320 0.736 
50 0.50 1.252 1.107 11.569 1.093 12.677 1.136 9.211 2.248 1.910 15.018 2.290 1.889 2.250 0.110 
60 0.50 1.327 1.131 14.767 1.115 15.956 1.171 11.802 2.231 1.860 16.635 2.260 1.294 2.230 0.051 
70 0.50 1.351 1.156 14.450 1.137 15.822 1.205 10.861 2.188 1.830 16.355 2.250 2.843 2.220 1.471 
80 0.50 1.391 1.180 15.174 1.159 16.674 1.238 10.990 2.144 1.800 16.029 2.230 4.030 2.210 3.097 
90 0.50 1.438 1.202 16.418 1.181 17.905 1.272 11.551 2.231 1.790 19.772 2.220 0.499 2.210 0.947 
110 0.50 1.481 1.248 15.766 1.223 17.456 1.339 9.591 2.144 1.750 18.362 2.210 3.097 2.200 2.631 
125 0.50 1.538 1.288 16.310 1.254 18.502 1.390 9.661 2.188 1.720 21.383 2.200 0.557 2.200 0.557 
125 0.75 1.557 1.243 20.160 1.256 19.369 1.384 11.132 2.231 1.730 22.461 2.210 0.947 2.200 1.395 
125 1.00 1.501 1.209 19.455 1.251 16.671 1.369 8.765 2.188 1.750 20.011 2.220 1.471 2.210 1.014 
125 1.50 1.436 1.163 19.000 1.231 14.258 1.331 7.338 2.231 1.790 19.772 2.230 0.051 2.220 0.499 
125 2.00 1.365 1.133 17.039 1.209 11.464 1.292 5.360 2.274 1.820 19.952 2.250 1.040 2.230 1.919 
150 0.50 1.695 1.355 20.022 1.305 23.018 1.474 13.022 2.144 1.700 20.694 2.200 2.631 2.200 2.631 
150 1.00 1.586 1.255 20.916 1.300 18.038 1.447 8.768 2.231 1.720 22.909 2.210 0.947 2.210 0.947 
200 0.50 1.817 1.436 20.962 1.403 22.807 1.638 9.833 2.231 1.700 23.806 2.190 1.844 2.200 1.395 
200 1.00 1.738 1.337 23.112 1.395 19.757 1.600 7.943 2.231 1.700 23.806 2.200 1.395 2.200 1.395 

For the state points of interest, the liquid character of the UEF becomes apparent from the relatively large magnitude of the first SSF peak (which is generally around ∼1.5 and even reaches 1.82), the well-resolved first SSF trough around 3kF (especially for rs ≳ 100), and the well-resolved second SSF peak above 4kF (only for rs ≳ 100). (i) The IET scheme always generates the most accurate SSF across the entire interval: within the long wavelength range of k ≲ 2kF, in the Lorentzian shaped region that surrounds the first maximum k ∼ 2 − 2.5kF, and within the short wavelength range of k ≳ 2.5kF. (ii) The qSTLS scheme has the worst performance. It strongly underestimates the position of the first SSF peak by around 20%. The same behavior was earlier observed for the STLS and the VS scheme.34 These results can be anticipated from classical OCP liquids, where it has been established that stronger Coulomb correlations (neglected in STLS, VS, and qSTLS) not only increase the SSF peak magnitude but also slightly displace it toward longer wavenumbers. (iii) On the other hand, the IET and HNC schemes provide very accurate predictions for the SSF peak positions, 0.05%–3.10% and 0.05%–4.03% relative deviations from the PIMC results, respectively, with the IET scheme having the slight edge. Thus, it can be concluded that the first SSF peak position is mainly controlled by strong correlations. (iv) Regardless of the state, the IET scheme greatly improves the HNC prediction for the SSF peak magnitude, with the relative deviations from PIMC results being 2.38%–13.02% and 5.26%–23.02%, respectively. (v) For all states, the IET SSF is characterized by a very accurate long wavelength behavior, especially for kkF.

Despite the large improvements over the HNC scheme, the IET scheme seems to generate SSFs that are not accurate enough to justify the very accurate predictions of interaction energies. Detailed inspection of Figs. 4 and 5 reveals that the very accurate interaction energies are the result of favorable error cancellation in Eq. (25) since the IET SSF tends to be slightly too large within kFk ≤ 2kF, becomes clearly too small within 2.0kFk ≤ 2.5kF, and tends to be somewhat too large within 2.5kFk ≤ 4.0kF. The same reasoning applies for the HNC scheme. Note that a similar favorable error cancellation is responsible for the success of STLS generated interaction energies within the warm dense matter regime.19,126

Let us briefly discuss the radial distribution functions (RDFs) g(r) extracted from PIMC simulations and computed with the three dielectric formalism schemes. Since the RDFs encode the same correlation information as SSFs, it is unsurprising that the IET scheme greatly improves the predictions of the HNC scheme for the first maximum and first non-zero minimum magnitude, that the IET and HNC schemes provide very accurate similar predictions for the first maximum and first non-zero minimum positions, and that the qSTLS scheme has the worst performance. Finally, it is important to discuss the un-physical negative RDF region near the origin, a common pathology of all dielectric formalism schemes stemming from the approximate treatment of quantum effects.19 The IET and HNC schemes have negative regions of similar extent and signal power, as expected from the RPA-level treatment of quantum effects. On the other hand, the qSTLS scheme has a negative region of less extent and power, as expected from its beyond-RPA quantum features.

Characteristic static local field corrections and static density responses extracted from our PIMC simulations and computed with the three dielectric formalism schemes are illustrated in Figs. 6 and 7.

FIG. 6.

Dependence of the paramagnetic electron liquid static local field correction (main figure) and static density response (inset) on the quantum coupling parameter rs. The results from the IET scheme (red solid line), the HNC scheme (blue dashed line), the qSTLS scheme (green dotted-dashed line), and our PIMC simulations (black crosses) for θ = 1.00 and varying rs = 50, 100, 150, 200 in the normalized wavenumber range k ≤ 4.5kF. For both quantities, the superiority of the new IET scheme within the long wavelength range and in the vicinity of the maximum is obvious.

FIG. 6.

Dependence of the paramagnetic electron liquid static local field correction (main figure) and static density response (inset) on the quantum coupling parameter rs. The results from the IET scheme (red solid line), the HNC scheme (blue dashed line), the qSTLS scheme (green dotted-dashed line), and our PIMC simulations (black crosses) for θ = 1.00 and varying rs = 50, 100, 150, 200 in the normalized wavenumber range k ≤ 4.5kF. For both quantities, the superiority of the new IET scheme within the long wavelength range and in the vicinity of the maximum is obvious.

Close modal
FIG. 7.

Dependence of the paramagnetic electron liquid static local field correction (main figure) and static density response (inset) on the quantum degeneracy parameter θ. The results from the IET scheme (red solid line), the HNC scheme (blue dashed line), the qSTLS scheme (green dotted-dashed line), and our PIMC simulations (black crosses) for rs = 100 and varying θ = 0.5, 1, 2, 4 in the normalized wavenumber range k ≤ 4.5kF. For both quantities, the superiority of the new IET scheme within the long wavelength range and in the vicinity of the maximum is obvious.

FIG. 7.

Dependence of the paramagnetic electron liquid static local field correction (main figure) and static density response (inset) on the quantum degeneracy parameter θ. The results from the IET scheme (red solid line), the HNC scheme (blue dashed line), the qSTLS scheme (green dotted-dashed line), and our PIMC simulations (black crosses) for rs = 100 and varying θ = 0.5, 1, 2, 4 in the normalized wavenumber range k ≤ 4.5kF. For both quantities, the superiority of the new IET scheme within the long wavelength range and in the vicinity of the maximum is obvious.

Close modal

Let us first discuss the static local field correction: (i) The PIMC LFC exhibits a pronounced first local maximum of magnitude ∼1.10–1.25 located at k ∼ 2.35–2.60kF, which is followed by a shallow local minimum located at k ∼ 3.05–3.45kF and is ultimately succeeded by a rather sharp short wavelength increase. The asymptotic behavior does not contradict the short wavelength limit G(k) = 1 − g(0)49,91 that was discussed earlier since this is valid for frequency independent LFCs (as confirmed by the HNC and the IET asymptotics), but not for the exact static LFC. In the ground state θ → 0, the asymptotic behavior of the exact static LFC is described by the Holas expansion G(k, 0) = B + Ck2, which predicts a parabolic divergence.130 Such a parabolic asymptote has been empirically observed to persist also when θ ∼ 1,34,131 but not in the classical limit θ ≫ 1, where one gets G(k, 0) = 1.43 The above observation is confirmed by the present PIMC LFC since as θ increases, the asymptotic limit gradually switches from parabolically diverging to nearly unity in the depicted wavenumber range; see Figs. 7(a)7(d). (ii) Regardless of the state point of interest, the IET scheme generates the most accurate static LFC and the qSTLS scheme generates the least accurate. The IET scheme improves the HNC prediction for the LFC first peak magnitude, with the relative deviations from PIMC results being 1.4%–10.1% and 2.6%–12.0%, respectively. However, the HNC prediction for the LFC first peak position is always more accurate. (iii) Perhaps, the most remarkable feature of the IET scheme’s LFC concerns its very accurate long wavelength behavior that extends up to nearly k ≲ 2kF. This does not necessarily guarantee that the IET scheme satisfies to a high degree the compressibility sum rule (CSR), an exact relation connecting the long wavelength static LFC to the second density derivative of the exchange–correlation free energy, i.e., in (n, T) thermodynamic variables and cgs units,19,43,51,52

limk0G(k)k2=14πe22n2nfex(n,T),
(27)

or in (rs, θ) thermodynamic variables and Hartree units,

limx0G(x)x2=π12λrs4Θ22Θ2+rs22rs2+4Θrs2Θrs2ΘΘ2rsrsf̃xc(rs,θ).
(28)

This is simply because each dielectric scheme satisfies its individual CSR that involves its individual exchange–correlation free energy and not the exact exchange–correlation free energy. However, when combined with the very accurate IET interaction energies, the very accurate IET LFC long wavelength limit makes it very likely that the IET scheme satisfies the CSR to a large degree. We postpone such an investigation to a future work owing to its high computational cost since as Eq. (28) suggests, it requires an accurate parameterization of the IET exchange–correlation free energy. Nevertheless, this IET feature can be exploited to provide an alternative route for the accurate calculation of the paramagnetic electron liquid’s isothermal compressibility without involving neither PIMC simulations nor thermodynamic integration.

Let us now discuss the static density response function defined by χ(k) ≡ χ(k, ω = 0). (i) The qualitative behavior of the “exact” χ(k) of the paramagnetic electron liquid has been discussed earlier34 and is confirmed by our PIMC results. It is evident that the magnitude of the pronounced χ(k) extremum increases as rs decreases (see Fig. 6) and as θ decreases (see Fig. 7) and that the width of the χ(k) extremum becomes more sharp as rs increases (see Fig. 6) and especially as θ decreases (see Fig. 7). (ii) Since χ(k) is obtained by setting ω = 0 to Eq. (4) and only involves the common among schemes χ0(k) and the varying among schemes G(k), a comparison with respect to χ(k) should reflect the comparison with respect to G(k). In fact, the IET scheme generates the most accurate χ(k) and the qSTLS scheme, the least accurate χ(k). The IET scheme strongly improves the HNC prediction for the extremum magnitude, but the HNC prediction for the extremum position is somewhat more accurate. (iii) Finally, the IET scheme’s static density response is again nearly exact for long wavelengths up to nearly k ≲ 2kF.

The static approximation utilizes the exact static LFC PIMC results to close the set of equations of the dielectric formalism.33,132 In other words, this approach replaces the exact dynamic LFC G(k, ω) with its exact static limit G(k) = G(k, ω = 0) that is accessible from PIMC simulations. The approach has been demonstrated to yield basically exact results for the dynamic structure factor and other related dynamic quantities in their entire non-trivial range, provided that rs ≲ 4, and to even reproduce the most prominent characteristics of the same quantities even at stronger coupling.132 The approach has also been revealed to lead to accurate results for the SSF only for k ≲ 2.5kF but to overestimate short-range correlations, a shortcoming that stems from the divergence of the exact static LFC short wavelength limit.59 The SSF-related accuracy of the static approximation is expected to deteriorate further at strong coupling, where the dynamic nature of the LFC has been documented to be more important even in the classical limit.52,78,79 In order to draw more quantitative conclusions, the exact static LFC PIMC results were employed as the closure of the exact two building blocks of the dielectric formalism; see Eqs. (3) and (4). It is evident that the static approximation has a low computational cost similar to that of the RPA.

Systematic comparisons revealed the following: (i) Regardless of the state point, the static approximation begins to overestimate the SSF prior to the main peak. For some state points, the overestimation takes place already at k ≲ 1.5kF, i.e., right after the long wavelength range. (ii) As expected, the static approximation is much less accurate in the strongly coupled regime than in the warm dense matter regime. However, as θ increases, the accuracy level becomes significantly higher and the deviations become restricted at progressively shorter wavelengths, as revealed in Fig. 8. (iii) Since the static approximation either reproduces or overestimates the SSF in the entire wavenumber range, there is no possibility of favorable error cancellation during the computation of the interaction energy. As a result, the absolute relative deviations in the interaction energies always exceed 10% and can even reach 40%.

FIG. 8.

The paramagnetic electron liquid’s static structure factor for rs = 100 and θ = 0.50 (a), θ = 0.75 (b), θ = 1.00 (c), and θ = 2.00 (d). The results from the static approximation (red solid line and red circles) and PIMC simulations (black crosses) in the normalized wavenumber range k ≤ 4.5kF. The static approximation consistently overestimates the SSF, especially at small values of the degeneracy parameter, which leads to highly inaccurate interaction energies.

FIG. 8.

The paramagnetic electron liquid’s static structure factor for rs = 100 and θ = 0.50 (a), θ = 0.75 (b), θ = 1.00 (c), and θ = 2.00 (d). The results from the static approximation (red solid line and red circles) and PIMC simulations (black crosses) in the normalized wavenumber range k ≤ 4.5kF. The static approximation consistently overestimates the SSF, especially at small values of the degeneracy parameter, which leads to highly inaccurate interaction energies.

Close modal

Moreover, given the success of the IIT interaction energy parameterization in the strongly coupled regime (see Sec. V A), it is tempting to check whether an extrapolation of the effective static approximation (ESA)59,60 scheme would also prove to be successful. The ESA utilizes an analytical representation for the frequency independent LFC,60 which is based on a neural net representation of PIMC static LFC results,131 but also incorporates the static LFC exact long wavelength limit as expressed by the CSR (employing the GDSMFB exchange–correlation free energy parametrization119) and the frequency independent LFC exact short wavelength limit as expressed by the asymptotic self-consistency condition [employing a g(0) parameterization based on extrapolated restricted PIMC data59]. Owing to the latter feature, the ESA scheme does not suffer from the aforementioned drawback of the static approximation. In fact, the ESA was recently demonstrated to be a fast and reliable tool for the computation of numerous thermodynamic, static, and dynamic quantities within the warm dense matter relevant range of 0.7 ≤ rs ≤ 20 and 0 ≤ θ ≤ 4.59,60 Taking into account the satisfactory performance of the extrapolation of the GDSMFB parameterization toward low densities and the very small (albeit unphysical negative) values of the extrapolated g(0) parameterization toward low densities, the ESA performance for the paramagnetic electron liquid mainly depends on whether the analytical representation of the intermediate wavelength PIMC static LFC results can be successfully extrapolated toward high coupling parameters. In order to draw conclusions, the analytical ESA LFC was used as the closure of the exact two building blocks of the dielectric formalism; see Eqs. (3) and (4). It is evident that the ESA scheme also has a low computational cost similar to that of the RPA.

Systematic comparisons revealed the following: (i) As expected, the ESA LFC always extrapolates very well toward the long and the short wavelength ranges. (ii) At least for 0.75 ≤ θ ≤ 2.00 and regardless of rs, the ESA well describes the position and the magnitude of the LFC peak. On the other hand, the ESA strongly overestimates the magnitude of the LFC peak for θ = 0.5 especially at strong coupling and underestimates the magnitude of the LFC peak for θ = 4.0. (iii) Regardless of the state point of interest, the main drawback of the ESA LFC lies at its transition to the short wavelength regime. This transition is always too abrupt and is not accompanied by small oscillations, as one would expect from a frequency independent LFC at strong coupling (see the IET LFC and the HNC LFC). These characteristics are not present in the warm dense matter regime or in dielectric schemes that are only equipped to describe weak correlations (see the qSTLS LFC). This drawback can be traced back to the implementation of a sigmoid activation function during the construction of ESA.60 (iv) Within 1 ≤ θ ≤ 2 and regardless of rs, the ESA accurately estimates the position and magnitude of the SSF peak. On the other hand, the ESA strongly overestimates the magnitude of the SSF peak for θ = 0.5, 0.75 especially at strong coupling and underestimates the magnitude of the SSF peak for θ = 4.0. (v) Regardless of the state point, the ESA SSF directly transitions to its unity short wavelength limit after the Lorentzian peak and does not feature a shallow minimum followed by a second maximum. (vi) For most of the state points, the ESA SSF overestimates the PIMC SSF not only prior to the peak but also after. As a consequence, it leads to very inaccurate interaction energies. To sum up, the ESA does not extrapolate well to strong coupling, which can mostly be attributed to the abrupt transition from the LFC peak to the short wavelength range. Thus, it can be concluded that the future utilization of an ESA-like scheme at low densities would require a more involved activation function. Some characteristic examples are illustrated in Fig. 9.

FIG. 9.

The paramagnetic electron liquid’s static local field correction (main) and static structure factor (inset) for rs = 125 and θ = 0.75 (a), θ = 1.00 (b), and θ = 1.50 (c). The results from the IET scheme (red solid line), the ESA scheme (cyan dashed line), and PIMC simulations (black crosses) in the normalized wavenumber range k ≤ 4.5kF. The ESA shortcoming that concerns the abrupt LFC (and thus SSF) transition from the peak vicinity to the short wavelength limit is apparent.

FIG. 9.

The paramagnetic electron liquid’s static local field correction (main) and static structure factor (inset) for rs = 125 and θ = 0.75 (a), θ = 1.00 (b), and θ = 1.50 (c). The results from the IET scheme (red solid line), the ESA scheme (cyan dashed line), and PIMC simulations (black crosses) in the normalized wavenumber range k ≤ 4.5kF. The ESA shortcoming that concerns the abrupt LFC (and thus SSF) transition from the peak vicinity to the short wavelength limit is apparent.

Close modal

In this work, we investigated the thermodynamic and the structural properties of the paramagnetic electron liquid with different dielectric theories and with ab initio path integral Monte Carlo simulations. First and foremost, we formulated a novel IET scheme that combines the dielectric formalism of many fermion systems with the integral equation theory of classical liquids. Essentially, the IET scheme incorporates a near-exact recent parameterization of the classical OCP bridge function into the existing HNC scheme. In addition, we carried out extensive PIMC simulations of the paramagnetic electron liquid for 16 state points. Combined with the recent availability of PIMC simulations for 20 state points, these new results substantially extend our current picture of the finite temperature UEF in the strongly coupled regime.

After an extensive comparison with “exact” PIMC results, various dielectric formalism schemes, and numerous extrapolated high-density parameterizations, we have demonstrated the following: (i) The IET scheme yields more accurate results for the interaction energy, static structure factor, static local field correction, and static density response of the paramagnetic electron liquid than all other known schemes (HNC, qSTLS, VS, and STLS). (ii) The IET interaction energies exhibit a remarkable agreement with the PIMC interaction energies having an accuracy within 0.68% and an average accuracy of 0.29%. This has been attributed to a favorable error cancellation in the course of the static structure factor integration. (iii) The IET long wavelength static local field correction, which is connected with the isothermal compressibility through the eponymous sum rule, is nearly indistinguishable from its “exact” PIMC counterpart up to nearly k ≲ 2kF. (iv) The utilization of the IIT interaction energy parameterization leads to very accurate results for the paramagnetic electron liquid (similar accuracy level to the IET), in contrast to extrapolations of the GDSMFB, KSDT, and corrKSDT exchange–correlation free energy parameterizations. This has been attributed to the incorporation of exact ground-state and classical limit results and to the introduction of an accurate ad hoc θ− interpolation function. (v) The static approximation is much less accurate in the strongly coupled regime than in the warm dense matter regime. This confirms that the dynamic nature of the local field correction is an essential ingredient for the paramagnetic electron liquid. On the other hand, the analytic effective static approximation, which is very reliable in the warm dense matter regime, cannot be extrapolated to the strongly coupled regime primarily due to its direct transition to the short wavelength limit after the peak of the frequency independent local field correction.

At this point, it is important to bring forth the two main drawbacks of the IET scheme. First, the classical OCP bridge function is introduced as an analytic b(r, Γ) parameterization and not as a b[h] functional. As a consequence, the bridge function only reacts to the classical Coulomb interactions and not to quantum mechanical interactions (exchange degeneracy and diffraction effects). This shortcoming manifests itself in the lack of an appropriate ground-state limit, given the Γ = 2λ2(rs/θ) mapping of quantum to classical states that stems from the classical coupling parameter definition Γ = e2/(dkBT). A remedy can be found by recalling the ground-state quantum coupling parameter definition Γq = e2/(dEF) that leads to Γq = 2λ2rs. This suggests an effective coupling parameter definition Γeff=e2/{d[(kBT)2+EF2]1/2}16 that leads to the new Γ=2λ2rs/1+θ2 mapping. The utilization of the effective mapping in lieu of the classical mapping in the IET scheme led to nearly identical results for the state points investigated here. Nevertheless, this should be further tested for low-density state points near the ground-state (θ → 0). It should be mentioned that different empirical mappings could also be utilized based on enforcing consistency with the compressibility sum rule, based on enforcing consistency with some PIMC results (e.g., the static structure factor peak magnitude) or based on adopting empirical relations for a quantum temperature from the classical mapping method.133–135 Such options would probably lead to an improved IET scheme that is even more accurate, but they are not desirable because they add a strong phenomenological element to an otherwise rigorous theoretical approach. Second, quantum effects are treated on the random phase approximation level. A more sophisticated quantum treatment can be achieved by essentially combining the qSTLS and IET schemes. Given the studied static structure factors of these schemes, this seems to be a very promising strategy that could improve the peak magnitude prediction of the IET scheme. Unfortunately, it also adds an inconsistent element to the approach because the bridge function incorporation technique is based on the classical fluctuation dissipation theorem for a frequency independent local field correction, while the qSTLS scheme necessarily leads to a dynamic local field correction.

Future work will focus on the exploration of the aforementioned possibilities to remedy the drawbacks of the current version of the IET scheme. Moreover, the characterization of the IET scheme’s self-consistency concerning the compressibility sum rule and the application of the IET scheme for different values of the spin polarization parameter will also be pursued in the future.

Finally, benefitting from the recent availability of extensive path integral Monte Carlo simulation results for the finite temperature θ ∼ 1 paramagnetic electron fluid (non-ideal gas and liquid), the following practical guidelines can be formulated concerning the optimal dielectric formalism scheme. Within the warm dense matter regime of rs ≲ 20, the semi-empirical ESA constitutes the most accurate scheme. At the moderate coupling 20 ≲ rs ≲ 50 regime, the first-principle HNC scheme is the most accurate. At the strong coupling regime of rs ≳ 50, the novel first-principle IET scheme is the most accurate. Naturally, the exact quantum coupling parameter boundaries between these three regimes depend on the specific degenerate parameter value.

The present work was partly funded by the Swedish National Space Agency under Grant No. 143/16. The present work was also partly funded by the Center for Advanced Systems Understanding (CASUS), which is financed by Germany’s Federal Ministry of Education and Research (BMBF) and by the Saxon Ministry for Science, Culture and Tourism (SMWK) with tax funds on the basis of the budget approved by the Saxon State Parliament. All the PIMC simulations were carried out at the Norddeutscher Verbund für Hoch-und Höchstleistungsrechnen (HLRN) under Grant No. shp00026 and on a Bull Cluster at the Center for Information Services and High Performance Computing (ZIH) at the Technische Universität Dresden. The IET, HNC, and qSTLS schemes were numerically solved on resources provided by the Swedish National Infrastructure for Computing (SNIC) at the NSC (Linköping University) that is partially funded by the Swedish Research Council under Grant Agreement No. 2018-05973.

The authors have no conflicts to disclose.

The data that support the findings of this study are available from the corresponding author upon reasonable request and within the article.

1.
C.
Kittel
,
Quantum Theory of Solids
(
John Wiley & Sons
,
New York
,
1987
).
2.
D.
Pines
,
Elementary Excitations in Solids
(
Perseus Books Publishing
,
Boston, MA
,
1999
).
3.
E. K. U.
Gross
,
E.
Runge
, and
O.
Heinonen
,
Many-Particle Theory
(
Adam Hilger
,
Bristol
,
1991
).
4.
G. D.
Mahan
,
Many-Particle Physics
(
Kluwer Academic; Plenum Publishers
,
New York
,
2000
).
5.
D.
Kremp
,
M.
Schlanges
, and
W.-D.
Kraeft
,
Quantum Statistics of Nonideal Plasmas
(
Springer-Verlag
,
Berlin
,
2005
).
6.
M.
Bonitz
,
Quantum Kinetic Theory
(
Springer International Publishing
,
Switzerland
,
2016
).
7.
J. P.
Lowe
and
K. A.
Peterson
,
Quantum Chemistry
(
Elsevier Academic Press
,
San Diego
,
2006
).
8.
R. G.
Parr
and
W.
Yang
,
Density-Functional Theory of Atoms and Molecules
(
Clarendon Press
,
Oxford
,
1989
).
9.
D.
Pines
and
P.
Nozieres
,
The Theory of Quantum Liquids
(
Taylor & Francis Group, LLC
,
New York
,
1989
), Vol. I.
10.
J.
Bardeen
,
L. N.
Cooper
, and
J. R.
Schrieffer
,
Phys. Rev.
108
,
1175
(
1957
).
11.
A. A.
Abrikosov
and
I. M.
Khalatnikov
,
Rep. Prog. Phys.
22
,
329
(
1959
).
12.
D.
Bohm
and
D.
Pines
,
Phys. Rev.
85
,
338
(
1952
).
13.
E.
Wigner
,
Phys. Rev.
46
,
1002
(
1934
).
14.
W.
Kohn
and
L. J.
Sham
,
Phys. Rev.
140
,
A1133
(
1965
).
15.
F.
Graziani
,
M. P.
Desjarlais
,
R.
Redmer
, and
S. B.
Trickey
,
Frontiers and Challenges in Warm Dense Matter
(
Springer International
,
Switzerland
,
2014
).
16.
M.
Bonitz
,
T.
Dornheim
,
Zh. A.
Moldabekov
,
S.
Zhang
,
P.
Hamann
,
H.
Kählert
,
A.
Filinov
,
K.
Ramakrishna
, and
J.
Vorberger
,
Phys. Plasmas
27
,
042710
(
2020
).
17.
B. A.
Remington
,
R. P.
Drake
, and
D. D.
Ryutov
,
Rev. Mod. Phys.
78
,
755
(
2006
).
18.
V. E.
Fortov
,
Phys. Usp.
52
,
615
(
2009
).
19.
T.
Dornheim
,
S.
Groth
, and
M.
Bonitz
,
Phys. Rep.
744
,
1
(
2018
).
20.
T.
Schoof
,
M.
Bonitz
,
A.
Filinov
,
D.
Hochstuhl
, and
J. W.
Dufty
,
Contrib. Plasma Phys.
51
,
687
(
2011
).
21.
T.
Dornheim
,
S.
Groth
,
A.
Filinov
, and
M.
Bonitz
,
New J. Phys.
17
,
073017
(
2015
).
22.
F. D.
Malone
,
N. S.
Blunt
,
E. W.
Brown
,
D. K. K.
Lee
,
J. S.
Spencer
,
W. M. C.
Foulkes
, and
J. J.
Shepherd
,
Phys. Rev. Lett.
117
,
115701
(
2016
).
23.
T.
Dornheim
,
S.
Groth
,
F. D.
Malone
,
T.
Schoof
,
T.
Sjostrom
,
W. M. C.
Foulkes
, and
M.
Bonitz
,
Phys. Plasmas
24
,
056303
(
2017
).
24.
K. S.
Singwi
, “
Correlations and the possibility of a charge-density-wave instability in quantum electron liquids
,” in
Recent Progress in Many-Body Theories
, edited by
H.
Kümmel
and
M. L.
Ristig
(
Springer
,
Berlin
,
1984
), pp.
219
228
.
25.
H.
Iyetomi
,
K.
Utsumi
, and
S.
Ichimaru
,
Phys. Rev. B
24
,
3226
(
1981
).
26.
A. W.
Overhauser
,
Phys. Rev.
167
,
691
(
1968
).
27.
S.
Ichimaru
and
S.
Tanaka
,
Phys. Rev. B
36
,
6182
(
1987
).
28.
A. W.
Overhauser
,
Phys. Rev.
128
,
1437
(
1962
).
29.
W. M. C.
Foulkes
,
L.
Mitas
,
R. J.
Needs
, and
G.
Rajagopal
,
Rev. Mod. Phys.
73
,
33
(
2001
).
30.
G.
Ortiz
,
M.
Harris
, and
P.
Ballone
,
Phys. Rev. Lett.
82
,
5317
(
1999
).
31.
M.
Holzmann
and
S.
Moroni
,
Phys. Rev. Lett.
124
,
206404
(
2020
).
32.
Y.
Takada
,
Phys. Rev. B
94
,
245106
(
2016
).
33.
T.
Dornheim
,
S.
Groth
,
J.
Vorberger
, and
M.
Bonitz
,
Phys. Rev. Lett.
121
,
255001
(
2018
).
34.
T.
Dornheim
,
T.
Sjostrom
,
S.
Tanaka
, and
J.
Vorberger
,
Phys. Rev. B
101
,
045129
(
2020
).
35.
S.
Tanaka
,
J. Chem. Phys.
145
,
214104
(
2016
).
36.
H. L.
Frisch
and
J. L.
Lebowitz
,
The Equilibrium Theory of Classical Fluids
(
Benjamin
,
New York
,
1964
).
37.
J. P.
Hansen
and
I. R.
McDonald
,
Theory of Simple Liquids
(
Academic Press
,
London
,
2006
).
38.
A.
Santos
,
A Concise Course on the Theory of Classical Liquids
(
Springer
,
Heidelberg
,
2016
).
39.
F.
Lucco Castello
,
P.
Tolias
, and
T.
Dornheim
, “
Classical bridge functions in classical and quantum plasma liquids
,” arXiv:2107.03537.
40.
F.
Perrot
and
M. W. C.
Dharma-wardana
,
Phys. Rev. Lett.
87
,
206404
(
2001
).
41.
M. W. C.
Dharma-wardana
,
Phys. Rev. B
82
,
195303
(
2010
).
42.
Y.
Liu
and
J.
Wu
,
J. Chem. Phys.
140
,
084103
(
2014
).
43.
S.
Ichimaru
,
H.
Iyetomi
, and
S.
Tanaka
,
Phys. Rep.
149
,
91
(
1987
).
44.
M.
Baus
and
J.-P.
Hansen
,
Phys. Rep.
59
,
1
(
1980
).
45.
T.
Ott
,
H.
Thomsen
,
J.-W.
Abraham
,
T.
Dornheim
, and
M.
Bonitz
,
Eur. Phys. J. D
72
,
84
(
2018
).
46.
G. F.
Giuliani
and
G.
Vignale
,
Quantum Theory of The Electron Liquid
(
Cambridge University Press
,
New York
,
2005
).
47.
N. D.
Drummond
,
Z.
Radnai
,
J. R.
Trail
,
M. D.
Towler
, and
R. J.
Needs
,
Phys. Rev. B
69
,
085116
(
2004
).
48.
R.
Balescu
,
Equilibrium and Nonequilibrium Statistical Mechanics
(
John Wiley & Sons
,
New York
,
1975
).
49.
S.
Tanaka
and
S.
Ichimaru
,
J. Phys. Soc. Jpn.
55
,
2278
(
1986
).
50.
A. A.
Kugler
,
J. Stat. Phys.
12
,
35
(
1975
).
51.
S.
Ichimaru
,
Rev. Mod. Phys.
54
,
1017
(
1982
).
52.
S.
Ichimaru
,
Statistical Plasma Physics
(
CRC Press
,
Boca Raton
,
2018
), Vols. I & II.
53.
K. S.
Singwi
,
M. P.
Tosi
,
R. H.
Land
, and
A.
Sjölander
,
Phys. Rev.
176
,
589
(
1968
).
54.
S.
Tanaka
,
S.
Mitake
, and
S.
Ichimaru
,
Phys. Rev. A
32
,
1896
(
1985
).
55.
P.
Vashishta
and
K. S.
Singwi
,
Phys. Rev. B
6
,
875
(
1972
).
56.
T.
Sjostrom
and
J.
Dufty
,
Phys. Rev. B
88
,
115123
(
2013
).
57.
S.
Ichimaru
,
Phys. Rev. A
2
,
494
(
1970
).
58.
S.
Tanaka
and
S.
Ichimaru
,
Phys. Rev. B
39
,
1036
(
1989
).
59.
T.
Dornheim
,
A.
Cangi
,
K.
Ramakrishna
,
M.
Böhme
,
S.
Tanaka
, and
J.
Vorberger
,
Phys. Rev. Lett.
125
,
235001
(
2020
).
60.
T.
Dornheim
,
Z. A.
Moldabekov
, and
P.
Tolias
,
Phys. Rev. B
103
,
165102
(
2021
).
61.
A.
Holas
and
S.
Rahman
,
Phys. Rev. B
35
,
2720
(
1987
).
62.
H. K.
Schweng
and
H. M.
Böhm
,
Phys. Rev. B
48
,
2037
(
1993
).
63.
J. M.
Bomont
,
Adv. Chem. Phys.
139
,
1
(
2008
).
64.
J.-M.
Caillol
and
D.
Gilles
,
J. Phys. A: Math. Theor.
43
,
105501
(
2010
).
65.
K.-C.
Ng
,
J. Chem. Phys.
61
,
2680
(
1974
).
66.
P.
Tolias
and
F.
Lucco Castello
,
Phys. Plasmas
26
,
043703
(
2019
).
67.
F.
Lucco Castello
and
P.
Tolias
,
Contrib. Plasma Phys.
61
,
e202000105
(
2021
).
68.
M.
Llano-Restrepo
and
W. G.
Chapman
,
J. Chem. Phys.
97
,
2046
(
1992
).
69.
J.
Kolafa
,
S.
Labík
, and
A.
Malijevský
,
Mol. Phys.
100
,
2629
(
2002
).
70.
F.
Lucco Castello
,
P.
Tolias
, and
J. C.
Dyre
,
J. Chem. Phys.
154
,
034501
(
2021
).
71.
P. D.
Poll
,
N. W.
Ashcroft
, and
H. E.
DeWitt
,
Phys. Rev. A
37
,
1672
(
1988
).
72.
H.
Iyetomi
,
S.
Ogata
, and
S.
Ichimaru
,
Phys. Rev. A
46
,
1051
(
1992
).
73.
S.
Ogata
,
Phys. Rev. E
53
,
1094
(
1996
).
74.
Y.
Rosenfeld
,
Phys. Rev. A
46
,
1059
(
1992
).
75.
F.
Lucco Castello
and
P.
Tolias
, “
Bridge functions of classical one-component plasmas
,” arXiv:2108.09574.
76.
A. A.
Kugler
,
J. Stat. Phys.
8
,
107
(
1973
).
77.
K.
Tago
,
K.
Utsumi
, and
S.
Ichimaru
,
Prog. Theor. Phys.
65
,
54
(
1981
).
78.
M. S.
Murillo
,
Phys. Plasmas
7
,
33
(
2000
).
79.
P.
Tolias
and
F.
Lucco Castello
,
Phys. Plasmas
28
,
034502
(
2021
).
80.
M. S.
Wertheim
,
J. Math. Phys.
5
,
643
(
1964
).
81.
E.
Thiele
,
J. Chem. Phys.
39
,
474
(
1963
).
82.
R. G.
Palmer
and
J. D.
Weeks
,
J. Chem. Phys.
58
,
4171
(
1973
).
83.
P.
Tolias
,
S.
Ratynskaia
, and
U.
de Angelis
,
Phys. Rev. E
90
,
053101
(
2014
).
84.
U.
Bengtzelius
,
W.
Götze
, and
A.
Sjölander
,
J. Phys. C: Solid State Phys.
17
,
5915
(
1984
).
85.
E.
Flenner
and
G.
Szamel
,
Phys. Rev. E
72
,
031508
(
2005
).
86.
F.
Lucco Castello
and
P.
Tolias
,
Molecules
26
,
669
(
2021
).
87.
X.-Z.
Yan
and
S.
Ichimaru
,
J. Phys. Soc. Jpn.
56
,
3853
(
1987
).
88.
P.
Arora
,
K.
Kumar
, and
R. K.
Moudgil
,
Eur. Phys. J. B
90
,
76
(
2017
).
89.
J. C.
Kimball
,
Phys. Rev. A
7
,
1648
(
1973
).
90.
G.
Vignale
,
Phys. Rev. B
38
,
6445
(
1988
).
91.
W.
Stolzmann
and
M.
Rösler
,
Contrib. Plasma Phys.
41
,
203
(
2001
).
92.
H.
Yasuhara
,
Solid State Commun.
11
,
1481
(
1972
).
93.
A. A.
Broyles
,
J. Chem. Phys.
33
,
456
(
1960
).
94.
D. M.
Ceperley
,
Rev. Mod. Phys.
67
,
279
(
1995
).
95.
M. F.
Herman
,
E. J.
Bruskin
, and
B. J.
Berne
,
J. Chem. Phys.
76
,
5150
(
1982
).
96.
M.
Takahashi
and
M.
Imada
,
J. Phys. Soc. Jpn.
53
,
963
974
(
1984
).
97.
D. M.
Ceperley
,
J. Stat. Phys.
63
,
1237
1267
(
1991
).
98.
M.
Troyer
and
U.-J.
Wiese
,
Phys. Rev. Lett.
94
,
170201
(
2005
).
99.
T.
Dornheim
,
Phys. Rev. E
100
,
023307
(
2019
).
100.
T.
Dornheim
,
J. Phys. A: Math. Theor.
54
,
335001
(
2021
).
101.
H.
Kleinert
,
Path Integrals in Quantum Mechanics, Statistics, Polymer Physics, and Financial Markets
(
World Scientific
,
Singapore
,
2009
).
102.
H.
De Raedt
and
B.
De Raedt
,
Phys. Rev. A
28
,
3575
(
1983
).
103.
B.
Militzer
,
Comput. Phys. Commun.
204
,
88
96
(
2016
).
104.
K.
Sakkos
,
J.
Casulleras
, and
J.
Boronat
,
J. Chem. Phys.
130
,
204109
(
2009
).
105.
T.
Dornheim
,
S.
Groth
, and
M.
Bonitz
,
Contrib. Plasma Phys.
59
,
e201800157
(
2019
).
106.
R. E.
Zillich
,
J. M.
Mayrhofer
, and
S. A.
Chin
,
J. Chem. Phys.
132
,
044103
(
2010
).
107.
N.
Metropolis
,
A. W.
Rosenbluth
,
M. N.
Rosenbluth
,
A. H.
Teller
, and
E.
Teller
,
J. Chem. Phys.
21
,
1087
(
1953
).
108.
D.
Chandler
and
P. G.
Wolynes
,
J. Chem. Phys.
74
,
4078
(
1981
).
109.
T.
Dornheim
,
S.
Groth
,
A. V.
Filinov
, and
M.
Bonitz
,
J. Chem. Phys.
151
,
014108
(
2019
).
110.
F.
Mezzacapo
and
M.
Boninsegni
,
Phys. Rev. A
75
,
033201
(
2007
).
111.
M.
Boninsegni
,
N. V.
Prokofév
, and
B. V.
Svistunov
,
Phys. Rev. E
74
,
036701
(
2006
).
112.
M.
Boninsegni
,
N.
Prokofév
, and
B.
Svistunov
,
Phys. Rev. Lett.
96
,
070601
(
2006
).
113.
S.
Chiesa
,
D. M.
Ceperley
,
R. M.
Martin
, and
M.
Holzmann
,
Phys. Rev. Lett.
97
,
076404
(
2006
).
114.
E. W.
Brown
,
B. K.
Clark
,
J. L.
DuBois
, and
D. M.
Ceperley
,
Phys. Rev. Lett.
110
,
146405
(
2013
).
115.
A. A.
Kugler
,
Phys. Rev. A
1
,
1688
(
1970
).
116.
N. D.
Drummond
,
R. J.
Needs
,
A.
Sorouri
, and
W. M. C.
Foulkes
,
Phys. Rev. B
78
,
125106
(
2008
).
117.
T.
Dornheim
,
S.
Groth
,
T.
Sjostrom
,
F. D.
Malone
,
W. M. C.
Foulkes
, and
M.
Bonitz
,
Phys. Rev. Lett.
117
,
156403
(
2016
).
118.
T.
Dornheim
and
J.
Vorberger
,
J. Chem. Phys.
154
,
144103
(
2021
).
119.
S.
Groth
,
T.
Dornheim
,
T.
Sjostrom
,
F. D.
Malone
,
W. M. C.
Foulkes
, and
M.
Bonitz
,
Phys. Rev. Lett.
119
,
135001
(
2017
).
120.
V. V.
Karasiev
,
T.
Sjostrom
,
J. W.
Dufty
, and
S. B.
Trickey
,
Phys. Rev. Lett.
112
,
076403
(
2014
).
121.
V. V.
Karasiev
,
J. W.
Dufty
, and
S. B.
Trickey
,
Phys. Rev. Lett.
120
,
076401
(
2018
).
122.
T.
Schoof
,
S.
Groth
,
J.
Vorberger
, and
M.
Bonitz
,
Phys. Rev. Lett.
115
,
130402
(
2015
).
123.
G. G.
Spink
,
R. J.
Needs
, and
N. D.
Drummond
,
Phys. Rev. B
88
,
085121
(
2013
).
124.
V. V.
Karasiev
,
S. B.
Trickey
, and
J. W.
Dufty
,
Phys. Rev. B
99
,
195134
(
2019
).
125.
S.
Tanaka
,
Contrib. Plasma Phys.
57
,
126
(
2017
).
126.
S.
Groth
,
T.
Dornheim
, and
M.
Bonitz
,
Contrib. Plasma Phys.
57
,
137
(
2017
).
127.
D. M.
Ceperley
and
B. J.
Alder
,
Phys. Rev. Lett.
45
,
566
(
1980
).
128.
W. L.
Slattery
,
G. D.
Doolen
, and
H. E.
DeWitt
,
Phys. Rev. A
21
,
2087
(
1980
).
129.
W. L.
Slattery
,
G. D.
Doolen
, and
H. E.
DeWitt
,
Phys. Rev. A
26
,
2255
(
1982
).
130.
A.
Holas
, “
Exact asymptotic expression for the static dielectric function of a uniform electron liquid at large wave vector
,” in
Strongly Coupled Plasma Physics
, edited by
F. J.
Rogers
and
H. E.
DeWitt
(
Plenum
,
New York
,
1987
), pp.
463
482
.
131.
T.
Dornheim
,
J.
Vorberger
,
S.
Groth
,
N.
Hoffmann
,
Zh. A.
Moldabekov
, and
M.
Bonitz
,
J. Chem. Phys.
151
,
194104
(
2019
).
132.
P.
Hamann
,
T.
Dornheim
,
J.
Vorberger
,
Z. A.
Moldabekov
, and
M.
Bonitz
,
Phys. Rev. B
102
,
125150
(
2020
).
133.
M. W. C.
Dharma-wardana
and
F.
Perrot
,
Phys. Rev. Lett.
84
,
959
(
2000
).
134.
S.
Dutta
and
J.
Dufty
,
Europhys. Lett.
102
,
67005
(
2013
).
135.
Y.
Liu
and
J.
Wu
,
J. Chem. Phys.
141
,
064115
(
2014
).