# Perturbation Calculation of the Uniform Electron Gas with a Transcorrelated
Hamiltonian

Hongjun Luo<EMAIL_ADDRESS>Max-Planck-Institute for Solid State Research,
Heisenbergstraße 1, 70569 Stuttgart, Germany  Ali Alavi<EMAIL_ADDRESS>Max-Planck-Institute for Solid State Research, Heisenbergstraße 1, 70569
Stuttgart, Germany  Department of Chemistry, University of Cambridge,
Lensfield Road, Cambridge, CB2 1EW, United Kingdom

###### Abstract

With a transcorrelated Hamiltonian, we perform a many body perturbation (MBPT)
calculation on the uniform electron gas in the high density regime. By using a
correlation factor optimised for a single determinant Jastrow ansatz, the
second order correlation energy is calculated as $\frac{1-\ln
2}{\pi^{2}}\ln(r_{s})-0.05075$. This already reproduces the exact logarithmic
term of the random phase approximation (RPA) result, while the constant term
is roughly $7\%$ larger than the RPA one. The close agreement with the RPA
method demonstrates that the transcorrelated method offers a viable and
potentially efficient method for treating metallic systems.

##  I Introduction

Electronic structure calculations usually suffer from numerical problems
originating from the Coulomb potential ${r_{12}^{-1}}$. On the one hand, the
short-range singularity leads to a non-smoothness in the many body wave
function, characterised by a cusp at the electron coalescence point $r_{12}=0$
Kato (1957), which causes a slow convergence with respect to the basis-set in
configuration expansions of the wave function . On the other hand, the slow
decay of the long-range tail of the Coulomb potential also causes a major
problem in calculations of metallic systems. Any straightforward perturbation
treatment leads to divergent values in the thermodynamic limit, and in order
to get a meaningful perturbation result one has to sum over some sub-sequences
of the perturbation expansion up to infinite order, a technique which is
usually referred to as the random phase approximation (RPA). Bohm and Pines
(1951); Pines and Bohm (1952); Bohm and Pines (1953); Gell-Mann and Brueckner
(1957); Ehrenreich and Cohen (1959)

One way to deal with the short-range singularity of the wave function is to
make use of the Jastrow ansatz Jastrow (1955)

$\Psi=e^{\tau}\Phi,$ (1)

to split the cusp into a pair correlation factor

$\tau(\mathbf{r}_{1},\cdots,\mathbf{r}_{N})=\frac{1}{2}\sum_{i,j}u(\mathbf{r}_{i},\mathbf{r}_{j}).$
(2)

With a factor $u$ fulfilling the asymptotic condition

$u(\mathbf{r}_{i},\mathbf{r}_{j})\sim\frac{1}{2}|\mathbf{r}_{i}-\mathbf{r}_{j}|,\quad,\mbox{when}\
|\mathbf{r}_{i}-\mathbf{r}_{j}|\rightarrow 0,$ (3)

the smoothness of the function $\Phi$ is one order higher than that of the
function $\Psi$ Fournais  _et al._ (2005). This product form is widely used in
quantum Monte Carlo methods Ceperley (1978); Umrigar, Nightingale, and Runge
(1993); Foulkes  _et al._ (2001), where the expectation value of the
Hamiltonian can be evaluated by random sampling. One can thus minimise this
expectation value of energy by adjusting the correlation factor $\tau$. Due to
the involved high dimensional integrals, it is difficult to use such kind of
variational treatment of the Jastrow wave function in conventional quantum
chemistry, where the function $\Phi$ needs to be further approximated by
configuration expansions. Alternatively, the transcorrelated (TC) method of
Boys and Handy Boys (1969); Boys and Handy (1969a, b, c, d); Handy (1969,
1972, 1973) offers a fairly efficient way to deal with the Jastrow ansatz.
With this method, the original Schrödinger equation is transformed into a non-
Hermitian eigenvalue problem

$\hat{H}_{\mbox{\scriptsize TC}}\Phi=E\Phi,\qquad\hat{H}_{\mbox{\scriptsize
TC}}=e^{-\tau}\hat{H}e^{\tau},$ (4)

via an exact similarity transformation, which removes the involved exponential
factor, and results in an exact effective Hamiltonian containing up to three
body terms:

$\displaystyle\hat{H}_{\mbox{\scriptsize TC}}$ $\displaystyle=$
$\displaystyle\hat{H}+[\hat{H},\tau]+\frac{1}{2}[[\hat{H},\tau],\tau]$ (5)
$\displaystyle=$
$\displaystyle\hat{H}-\sum_{i}\left(\frac{1}{2}\triangledown_{i}^{2}\tau+(\triangledown_{i}\tau)\cdot\triangledown_{i}+\frac{1}{2}(\triangledown_{i}\tau)^{2}\right).$

This method was initially studied by Boys and Handy for the single determinant
ansatz, and has only recently been combined with configuration expansion
methods of quantum chemistryTen-no (2000); Hino, Tanimura, and Ten-no (2001,
2002); Luo (2011); Luo and Alavi (2018); Cohen  _et al._ (2019); Guther  _et
al._ (2021); Liao  _et al._ (2021); Schraivogel  _et al._ (2021). These works
reveal a high level of efficiency of the transcorrelated Hamiltonian in
handling short range correlations, that the resolving of the cusp leads to
significant speed up of the basis convergence, and further, more speed up can
be achieved when more generalised factors are used (instead of the simple
$r_{12}$ type or $F(r_{12})$ type factors).

Regarding the study of long-range correlation with the transcorrelated method,
there exists only a few works dealing with the single determinant ansatz for
uniform electron gases (UEG) Armour (1980); N. Umezawa and S. Tsuneyuki
(2004); Luo (2012). By choosing a proper correlation factor $\tau$, the
leading order singularity of the Coulomb potential can be removed from the
effective Hamiltonian, both for the short-range and the long-range
singularities. We therefore expect that the divergence of the perturbation
energy for metallic systems can be cured to some extent. At least it should be
possible to get non-divergent results in the thermodynamic limit at low order
perturbation theory without any RPA-type treatment. It would then be
interesting to know whether this kind of non-RPA treatment leads to meaningful
results. In this work, we will test the performance of a MBPT treatment of the
transcorrelated Hamiltonian on UEG in the high density regime, which has a
well established RPA result Gell-Mann and Brueckner (1957); Onsager, Mittag,
and Stephen (1966)

$\frac{E_{c}}{N}=\frac{1-\ln 2}{\pi^{2}}\ln{r_{s}}-0.047+\cdots.$ (6)

In the next section, we describe the perturbation treatment of the TC
Hamiltonian, and in sections III and IV we present details of the calculations
of the first and the second order energies, where a Jastrow factor optimised
for the Hartree-Fock reference function is used. In section V, we also present
some calculation results with the widely used Coulomb-Yukawa factor, which
will be followed by some concluding remarks in section VI.

##  II Perurbation treatment of the transcorrelated Hamiltonian

Uniform electron gases are described by the Hamiltonian

$\hat{{\cal
H}}=-\frac{1}{2}\sum_{i}\triangledown_{i}^{2}+\sum_{i<j}\frac{1}{r_{ij}}.$ (7)

The system is characterized by the charge density $\rho$ (or equivalently the
Wigner–Seitz radius $r_{s}$, or the Fermi wavevector $k_{F}$). In order to get
the explicit dependence of the physical quantities on $r_{s}$, we re-scale the
spatial coordinates

${\bf r}_{i}=\alpha r_{s}{\bf
r}^{\prime}_{i},\quad\alpha=\left(\frac{4}{9\pi}\right)^{\frac{1}{3}},$ (8)

so that the system keeps a fixed density

$r_{s}^{\prime}=\frac{1}{\alpha},\
k_{F}^{\prime}=1,\rho^{\prime}=\frac{N}{\Omega}=\frac{1}{3\pi^{2}},$ (9)

but has a re-scaled Hamiltonian

$\displaystyle\hat{{\cal H}}$ $\displaystyle=$
$\displaystyle\frac{1}{\alpha^{2}r_{s}^{2}}\hat{H},$ (10)
$\displaystyle\hat{H}$ $\displaystyle=$
$\displaystyle-\frac{1}{2}\sum_{i}\triangledown_{i}^{2}+\alpha
r_{s}\sum_{i<j}\frac{1}{r_{ij}}.$ (11)

The formal volume $\Omega$ and the particle number $N$ will finally be taken
as $\infty$ in the thermodynamic limit. In this paper we will only use this
re-scaled form and thus will ignore the prime on the variables. We need only
to remember that the final energies should be divided by $\alpha^{2}r_{s}^{2}$
according to equation (10). Then we use Jastrow ansatz (1) for the ground
state wave function of $\hat{H}$, and the two-body correlation factor
$u(\mathbf{r}_{i},\mathbf{r}_{j})$ can be written as
$u(\mathbf{r}_{i}-\mathbf{r}_{j})$ due to the translational symmetry. The
transcorrelated Hamiltonian (5) can be derived straightforwardly in the second
quantisation formalism Luo (2012) in terms of plane wave orbitals

$\displaystyle\hat{H}_{\mbox{\scriptsize TC}}$ $\displaystyle=$ $\displaystyle
e^{-\tau}\hat{H}e^{\tau}=\hat{H}+\hat{K}+\hat{L},$ (12) $\displaystyle\hat{K}$
$\displaystyle=$
$\displaystyle\frac{1}{2\Omega}\sum_{\sigma\sigma^{\prime}}\sum_{kpq}\left(k^{2}\tilde{u}(\mathbf{k})-(\mathbf{p}-\mathbf{q})\cdot\mathbf{k}\
\tilde{u}(\mathbf{k})-{\left(\mathcal{F}(\triangledown
u)^{2}\right)}(\mathbf{k})\right)$ (13) $\displaystyle
a_{\mathbf{p}-\mathbf{k},\sigma}^{\dagger}a_{\mathbf{q}+\mathbf{k},\sigma^{\prime}}^{\dagger}a_{\mathbf{q},\sigma^{\prime}}a_{\mathbf{p},\sigma},$
$\displaystyle\hat{L}$ $\displaystyle=$
$\displaystyle-\frac{1}{2\Omega^{2}}\sum_{\sigma\sigma^{\prime}\sigma^{\prime\prime}}\sum_{kk^{\prime}pqs}\tilde{u}(\mathbf{k})\tilde{u}(\mathbf{k}^{\prime})\mathbf{k}^{\prime}\cdot\mathbf{k}\
a_{\mathbf{p}-\mathbf{k},\sigma}^{\dagger}a_{\mathbf{q}+\mathbf{k}^{\prime},\sigma^{\prime}}^{\dagger}a_{\mathbf{s}+\mathbf{k}-\mathbf{k}^{\prime},\sigma^{\prime\prime}}^{\dagger}a_{\mathbf{s},\sigma^{\prime\prime}}a_{\mathbf{q},\sigma^{\prime}}a_{\mathbf{p},\sigma},$
(14)

where $\tilde{u}(\mathbf{k})$ is the Fourier transformation of
$u(\mathbf{r})$. ${\left(\mathcal{F}(\triangledown u)^{2}\right)}(\mathbf{k})$
is the Fourier transformation of $(\triangledown u(\mathbf{r}))^{2}$ and can
be calculated via the convolution theorem

${\left(\mathcal{F}(\triangledown
u)^{2}\right)}(\mathbf{k})=-\frac{1}{(2\pi)^{3}}\int(\mathbf{k}-\mathbf{k}^{\prime})\cdot\mathbf{k}^{\prime}\tilde{u}(|\mathbf{k}-\mathbf{k}^{\prime}|)\tilde{u}(k^{\prime})d^{3}k^{\prime}.$
(15)

Here we have used the fact that in the thermodynamic limit, the function
$\tilde{u}(\mathbf{k})$ is actually a 1D function of $k$ due to rotational
symmetry.

It is of crucial importance to take a proper Jastrow factor $u$. It has to
fulfil asymptotic conditions at both of the limits Luo (2012)

$\displaystyle\tilde{u}(k)$ $\displaystyle\sim$
$\displaystyle-\frac{w(k)}{k^{2}}\qquad\mbox{when }k\sim\infty,$ (16)
$\displaystyle\tilde{u}(k)$ $\displaystyle\sim$
$\displaystyle-\frac{\sqrt{w(k)}}{\sqrt{\rho}k}\qquad\mbox{when }k\sim 0,$
(17)

where $w(k)=\frac{4\pi\alpha r_{s}}{k^{2}}$ is the Fourier transformation of
the Coulomb potential and $\rho$ is the constant density given in equation
(9). These conditions remove the leading order singular term of the Coulomb
potential and lead to a finite second order perturbation energy. However, they
do not guarantee the quality of the perturbation energy. Based on some recent
numerical studies of the TC method on small molecular systems Cohen  _et al._
(2019); Guther  _et al._ (2021), the correlation factor optimised for a single
determinant reference $\Phi$ (e.g., the Hartree-Fock wave function) also
performs well when $\Phi$ is further treated by configuration expansions, at
least for weakly correlated systems. For UEG in the high density limit, the
correlation factor optimised for the Hartree-Fock reference takes the
following form

$\tilde{u}(k)=\frac{k^{2}-\sqrt{k^{4}+4\rho w(k)k^{2}T^{2}_{2}(k)}}{2\rho
k^{2}T_{2}(k)}.$ (18)

$T_{2}(k)$ is a volume factor defined by

$T_{2}(k)\equiv\frac{2}{N}\sum_{\mathbf{p}}\Theta(1-p)\Theta(|\mathbf{p}-\mathbf{k}|-1),$
(19)

with $\Theta$ being the Heaviside step function. In the thermodynamic limit,
$T_{2}(k)$ can be calculated as

$T_{2}(k)=\begin{cases}1,&k>2,\\\ \frac{3}{4}k-\frac{1}{16}k^{3},&k\leq
2.\end{cases}$ (20)

The expression (18) of $\tilde{u}(k)$ is a solution to a quadratic equation

$\rho k^{2}T_{2}(k)\tilde{u}^{2}(k)-k^{2}\tilde{u}(k)-w(k)T_{2}(k)=0,$ (21)

which can can be derived both by a variational treatment of the Jastrow ansatz
Gaskell (1961, 1962); Talman (1974) and by transcorrelated method Armour
(1980); Luo (2012).

The Hamiltonian (5) is then split into two parts for a perturbation treatment

$\hat{H}_{\mbox{\scriptsize
TC}}=\hat{H}_{0}+\hat{W}_{1},\qquad\hat{W}_{1}=\hat{W}+\hat{K}+\hat{L},$ (22)

where for simplicity we take only the kinetic terms in $\hat{H}_{0}$

$\displaystyle\hat{H}_{0}$ $\displaystyle=$
$\displaystyle\sum_{p\sigma}\frac{1}{2}p^{2}\hat{n}_{\mathbf{p},\sigma},$ (23)
$\displaystyle\hat{W}$ $\displaystyle=$
$\displaystyle\frac{1}{2\Omega}\sum_{\sigma\sigma^{\prime}}\sum_{kpq}\ w(k)\
a_{\mathbf{p}-\mathbf{k},\sigma}^{\dagger}a_{\mathbf{q}+\mathbf{k},\sigma^{\prime}}^{\dagger}a_{\mathbf{q},\sigma^{\prime}}a_{\mathbf{p},\sigma}.$
(24)

With this partition of the TC Hamiltonian, the perturbation energies can be
written as

$\displaystyle E_{0}$ $\displaystyle=$
$\displaystyle\langle\Phi_{0}|\hat{H}_{0}|\Phi_{0}\rangle,$ (25)
$\displaystyle E_{1}$ $\displaystyle=$
$\displaystyle\langle\Phi_{0}|\hat{W}_{1}|\Phi_{0}\rangle,$ (26)
$\displaystyle E_{2}$ $\displaystyle=$
$\displaystyle-\langle\Phi_{0}|(\hat{W}_{1}-E_{1})\frac{1}{\hat{H}_{0}-E_{0}}(\hat{W}_{1}-E_{1})|\Phi_{0}\rangle,$
(27)

where $E_{0}$ is simply the total kinetic energy and $E_{1}$ contains the
exchange energy $\langle\Phi_{0}|\hat{W}|\Phi_{0}\rangle$ and the first order
correlation energy $\langle\Phi_{0}|\hat{K}+\hat{L}|\Phi_{0}\rangle$. The
second order energy contributes only to the correlation energy.

##  III The first order perturbation energy

The first order correlation energy is composed of the expectation values of
the two body operator $\hat{K}$ and the three body operator $\hat{L}$. For
calculations of these values, we need to find all possible contractions of the
creation and annihilation operators. From the two body operator we get a
direct and an exchange contribution

$\langle\Phi|\hat{K}|\Phi\rangle=-\frac{N^{2}}{2\Omega}{\left(\mathcal{F}(\triangledown
u)^{2}\right)}(\mathbf{0})+\frac{N}{2\Omega}\sum_{\mathbf{k}}{\left(\mathcal{F}(\triangledown
u)^{2}\right)}(\mathbf{k})\left(1-T_{2}(k)\right).$ (28)

Here only the ${\left(\mathcal{F}(\triangledown u)^{2}\right)}$ term of
$\hat{K}$ makes a non-vanishing contribution, while the other terms coming
from $[\hat{H},\tau]$ have no contribution to the expectation value. The
expectation value of the three body operator also has two possible
contractions

$<\Phi|\hat{L}|\Phi>=\frac{N^{2}}{2\Omega^{2}}\sum_{\mathbf{k}}\left(1-T_{2}(k)\right)k^{2}\tilde{u}^{2}(\mathbf{k})-\frac{N}{\Omega^{2}}\sum_{\mathbf{k}\mathbf{k}^{\prime}}\tilde{u}(\mathbf{k})\tilde{u}(\mathbf{k}^{\prime})\mathbf{k}\cdot\mathbf{k}^{\prime}O_{3}(\mathbf{k},\mathbf{k}^{\prime}).$
(29)

The volume factor $O_{3}$ is defined as

$O_{3}(\mathbf{k},\mathbf{k}^{\prime})=\frac{2}{N}\sum_{\mathbf{p}}\Theta(1-p)\Theta(1-|\mathbf{p}-\mathbf{k}|)\Theta(1-|\mathbf{p}-\mathbf{k}^{\prime}|),$
(30)

which can be calculated analytically in the thermodynamic limitLuo (2012).

In the thermodynamic limit the summation should be replaced by a integration

$\frac{1}{\Omega}\sum_{\mathbf{k}}\rightarrow\frac{1}{(2\pi)^{3}}\int d^{3}k.$
(31)

This replacement will be used in the following context without any
notification.

The first term in equation (29) is partly cancelled by the first term in
equation (28) and we have finally three terms in the expression of the first
order correlation energy

$\displaystyle E_{1}^{c}$ $\displaystyle=$ $\displaystyle
E_{11}^{c}+E_{12}^{c}+E_{13}^{c},$ (32) $\displaystyle E_{11}^{c}$
$\displaystyle=$ $\displaystyle-\frac{\rho N}{2(2\pi)^{3}}\int
T_{2}(k)\mathbf{k}^{2}\tilde{u}^{2}(\mathbf{k})d^{3}k,$ (33) $\displaystyle
E_{12}^{c}$ $\displaystyle=$ $\displaystyle\frac{N}{2(2\pi)^{3}}\int
d^{3}k\left({\cal F}((\triangledown
u)^{2})(\mathbf{k})\left(1-T_{2}(k)\right)\right),$ (34) $\displaystyle
E_{13}^{c}$ $\displaystyle=$ $\displaystyle-\frac{N}{(2\pi)^{6}}\int
d^{3}kd^{3}k^{\prime}\
\tilde{u}(\mathbf{k})\tilde{u}(\mathbf{k}^{\prime})\mathbf{k}\cdot\mathbf{k}^{\prime}O_{3}(\mathbf{k},\mathbf{k}^{\prime}).$
(35)

For a more transparent treatment, we reformulate the expression (18) of
$\tilde{u}(k)$ as

$\tilde{u}(k)=-\alpha r_{s}\frac{8\pi
T_{2}(k)}{k^{2}\left(k^{2}+\sqrt{k^{4}+\frac{16}{3\pi}\alpha
r_{s}T_{2}^{2}(k)}\right)}.$ (36)

The expression of the first order energy contains only quadratic terms of
$\tilde{u}(k)$ and the prefactor $\alpha r_{s}$ of $\tilde{u}$ in equation
(36) will be cancelled by the prefactor $\frac{1}{(\alpha r_{s})^{2}}$ in the
final energy expression according to equation (10). If $r_{s}$ is set to be
$0$ in the remaining integrands in the first order energy, it is found that
$E^{c}_{11}$ will diverge while the expression of $E^{c}_{12}$ and
$E^{c}_{13}$ still give finite results. The singular term in $E^{c}_{11}$
originates from the integration over small $k$. In order to get a analytic
result of the singular term we rearrange the expression of $E_{11}^{c}$ as

$\displaystyle\frac{E_{11}^{c}}{N}$ $\displaystyle=$
$\displaystyle-\frac{1}{6\pi^{2}(2\pi)^{3}}\int
T_{2}(k)k^{2}\tilde{u}^{2}(k)d^{3}k$ (38) $\displaystyle=$
$\displaystyle-\frac{3}{16}\int_{0}^{2}T_{2}(k)\left(\frac{k^{2}-\sqrt{k^{4}+\frac{16}{3\pi}\alpha
r_{s}T_{2}^{2}(k)}}{T_{2}(k)}\right)^{2}dk-\frac{4\pi}{6\pi^{2}(2\pi)^{3}}\int_{2}^{\infty}T_{2}(k)k^{4}\tilde{u}^{2}(k)dk$
$\displaystyle=$
$\displaystyle-\frac{3}{16}\int_{0}^{2}\frac{3k}{4}\left(\frac{k^{2}-\sqrt{k^{4}+\frac{16}{3\pi}\alpha
r_{s}(\frac{3}{4}k)^{2}}}{\frac{3}{4}k}\right)^{2}dk$
$\displaystyle-\frac{3}{16}\int_{0}^{2}\left[T_{2}(k)\left(\frac{\frac{16}{3\pi}\alpha
r_{s}T_{2}(k)}{k^{2}+\sqrt{k^{4}+\frac{16}{3\pi}\alpha
r_{s}T_{2}^{2}(k)}}\right)^{2}-\frac{3k}{4}\left(\frac{\frac{16}{3\pi}\alpha
r_{s}(\frac{3}{4}k)}{k^{2}+\sqrt{k^{4}+\frac{16}{3\pi}\alpha
r_{s}(\frac{3}{4}k)^{2}}}\right)^{2}\right]dk$
$\displaystyle-\frac{1}{12\pi^{4}}\int_{2}^{\infty}T_{2}(k)k^{4}\tilde{u}^{2}(k)dk.$

In the first line of expression (38) $T_{2}(k)$ is simplified with its leading
term $\frac{3}{4}k$ in the small k region, and this simplified integrand is
subtracted out from the original integrand in the second line. The simplified
integral in the first line can be calculated analytically as

$\displaystyle-\frac{3}{16}\int_{0}^{2}\frac{3k}{4}\left(\frac{k^{2}-\sqrt{k^{4}+\frac{16}{3\pi}\alpha
r_{s}(\frac{3}{4}k)^{2}}}{\frac{3}{4}k}\right)^{2}dk$ (39) $\displaystyle=$
$\displaystyle\alpha^{2}r_{s}^{2}\left[\frac{9}{32\pi^{2}}\ln(r_{s})+\frac{1}{16}\left(\frac{3}{\pi}\right)^{2}\left(\ln\left(2\sqrt{\frac{3\alpha}{\pi}}\right)-\ln
8\right)+\frac{1}{64}\left(\frac{3}{\pi}\right)^{2}\right]+o(r_{s}^{2}).$

In the second line, the singularities of the two terms cancel each other and
hence they contribute at most only to the constant term in the correlation
energy. For the calculation of the contribution to the constant, we can simply
take out the common $\alpha^{2}r_{s}^{2}$ factor and set $r_{s}=0$ in the
remaining integrand, which gives

$\displaystyle-\frac{3}{16}\int_{0}^{2}\left[T_{2}(k)\left(\frac{\frac{16}{3\pi}\alpha
r_{s}T_{2}(k)}{k^{2}+\sqrt{k^{4}+\frac{16}{3\pi}\alpha
r_{s}T_{2}^{2}(k)}}\right)^{2}-\frac{3k}{4}\left(\frac{\frac{16}{3\pi}\alpha
r_{s}(\frac{3}{4}k)}{k^{2}+\sqrt{k^{4}+\frac{16}{3\pi}\alpha
r_{s}(\frac{3}{4}k)^{2}}}\right)^{2}\right]dk$ (40) $\displaystyle=$
$\displaystyle-\alpha^{2}r_{s}^{2}\left(-\frac{8\cdot
3^{2}}{16^{2}\pi^{2}}+\frac{12}{16^{2}\pi^{2}}-\frac{1}{9\cdot
2^{5}\pi^{2}}\right)+o(r_{s}^{2}).$

In a similar way, the integral in the third line can be easily calculated as

$-\frac{1}{12\pi^{4}}\int_{2}^{\infty}T_{2}(k)k^{4}\tilde{u}^{2}(k)dk=-\frac{\alpha^{2}r_{s}^{2}}{18\pi^{2}}+o(r_{s}^{2}).$
(41)

$E_{12}^{c}$ and $E_{13}^{c}$ contribute also only at most to the constant
term of the correlation energy and can also be simplified with the above
method. After this simplification, it is still too difficult to be calculated
analytically due to the complicated integrand. There is, however, no problem
for a numerical calculation of such 1D integrals, which leads to

$\displaystyle\frac{E_{12}^{c}}{N}$ $\displaystyle=$ $\displaystyle
0.02352\alpha^{2}r_{s}^{2}+o(r_{s}^{2}),$ (42)
$\displaystyle\frac{E_{13}^{c}}{N}$ $\displaystyle=$
$\displaystyle-0.01310\alpha^{2}r_{s}^{2}+o(r_{s}^{2}).$ (43)

By adding all the results together, we get the final result for the first
order correlation energy

$\frac{E_{1}^{c}}{N}=\alpha^{2}r_{s}^{2}\left[\frac{9}{32\pi^{2}}\ln(r_{s})-0.05576\right]+o(r_{s}^{2}).$
(44)

The first order correlation energy is essentially the correlation energy
produced by the single determinant Jastrow ansatz. It is interesting to find
that the leading logarithmic term $\frac{9}{32\pi^{2}}\ln(r_{s})$ is the same
as that obtained by Talman based on linked cluster expansion Talman (1974).
This logarithmic term makes up roughly $90\%$ of the exact RPA logarithmic
term $\frac{1-\ln{2}}{\pi^{2}}\ln r_{s}$, which is also a typical portion of
correlated energy one may expect from variational quantum Monte Carlo
calculations based on single determinant Jastrow ansatz.

##  IV The second order perturbation energy

The transcorrelated effective potential $\hat{W}_{1}$ depends on $r_{s}$
nonlinearly. Besides the original Coulomb potential $w(k)$, which is
proportional to $r_{s}$, there are some terms which rely on $\tilde{u}(k)$
linearly and others quadratically. According to equation (36), in the high
density regime where $r_{s}\sim 0$, $\tilde{u}(k)$ scales roughly linearly
with $r_{s}$. To get an estimation of the contribution to the second order
energy $E_{2}$ from each quadratic term in the effective potential
$\hat{W}_{1}$, we can first set $r_{s}=0$ in the denominator of the expression
(36) of $\tilde{u}(k)$

$\tilde{u}(k)\approx-\frac{4\pi\alpha r_{s}T_{2}(k)}{k^{4}}.$ (45)

If this leads only to finite integrals in the calculation of $E_{2}$, we can
ignore this quadratic term for the current calculation, since it gives at most
a $O(r_{s}^{3})$ contribution to $E_{2}$ according to equation (27). Based on
this estimation we find there is only one quadratic term left, where the
approximation (45) leads to divergent integrals in the calculation of $E_{2}$.
This term looks like

$-\frac{\rho}{2\Omega}\sum_{\sigma\sigma^{\prime}}\sum_{kpq}k^{2}\tilde{u}^{2}(\mathbf{k})a_{\mathbf{p}-\mathbf{k},\sigma}^{\dagger}a_{\mathbf{q}+\mathbf{k},\sigma^{\prime}}^{\dagger}a_{\mathbf{q},\sigma^{\prime}}a_{\mathbf{p},\sigma},$
(46)

which is a contraction of the three body potential $\hat{L}$ in equation (14)
by taking $\mathbf{k}=\mathbf{k}^{\prime}$. Consequently for the calculation
of the two leading terms in the second order energy $E_{2}$, we need only to
take the following approximation of of the effective potential

$\displaystyle\hat{W}_{1}$ $\displaystyle\approx$
$\displaystyle\frac{1}{2\Omega}\sum_{\sigma\sigma^{\prime}}\sum_{kpq}\
G(\mathbf{p},\mathbf{q},\mathbf{k})\
a_{\mathbf{p}-\mathbf{k},\sigma}^{\dagger}a_{\mathbf{q}+\mathbf{k},\sigma^{\prime}}^{\dagger}a_{\mathbf{q},\sigma^{\prime}}a_{\mathbf{p},\sigma},$
(47) $\displaystyle G(\mathbf{p},\mathbf{q},\mathbf{k})$ $\displaystyle=$
$\displaystyle G_{0}(k)-(\mathbf{p}-\mathbf{q})\cdot\mathbf{k}\tilde{u}(k),$
(48) $\displaystyle G_{0}(k)$ $\displaystyle=$ $\displaystyle
w(k)+k^{2}\tilde{u}(k)-\rho k^{2}\tilde{u}^{2}(k).$ (49)

By taking this two body potential approximation of $\hat{W}_{1}$, the second
order energy consists two contributions

$\displaystyle E_{2}$ $\displaystyle=$
$\displaystyle-\langle\Phi_{0}|(\hat{W}_{1}-E_{1})\frac{1}{\hat{H}_{0}-E_{0}}(\hat{W}_{1}-E_{1})|\Phi_{0}\rangle$
(50) $\displaystyle\approx$ $\displaystyle E_{21}+E_{22},$

obtained by a direct contraction and an exchange contraction respectively.

The direct contraction term of $E_{2}$ can be written as

$\displaystyle E_{21}$ $\displaystyle=$
$\displaystyle-\frac{2}{4\Omega^{2}}\sum_{\sigma\sigma^{\prime}}\sum_{kpq}\frac{G(\mathbf{p},\mathbf{q},\mathbf{k})G(\mathbf{p}-\mathbf{k},\mathbf{q}+\mathbf{k},-\mathbf{k})}{\frac{1}{2}\left((\mathbf{p}-\mathbf{k})^{2}+(\mathbf{q}+\mathbf{k})^{2}-p^{2}-q^{2}\right)}\Theta(1-p)\Theta(1-q)\Theta(|\mathbf{p}-\mathbf{k}|-1)\Theta(|\mathbf{q}+\mathbf{k}|-1)$
(51) $\displaystyle=$
$\displaystyle-\frac{4}{\Omega^{2}}\sum_{kpq}\int_{0}^{\infty}dt\left\\{G(\mathbf{p},\mathbf{q},\mathbf{k})G(\mathbf{p}-\mathbf{k},\mathbf{q}+\mathbf{k},-\mathbf{k})\right.$
$\displaystyle\qquad\left.\Theta(1-p)\Theta(1-q)\Theta(|\mathbf{p}-\mathbf{k}|-1)\Theta(|\mathbf{q}+\mathbf{k}|-1)e^{-\left(2k^{2}-2(\mathbf{p}-\mathbf{q})\cdot\mathbf{k}\right)t}\right\\}$
$\displaystyle=$
$\displaystyle-\frac{4}{\Omega^{2}}\sum_{kpq}\int_{0}^{\infty}dt\left\\{\left(\left(G_{0}(k)-k^{2}\tilde{u}(k)\right)^{2}-\left(k^{2}-(\mathbf{p}-\mathbf{q})\cdot\mathbf{k}\right)^{2}\tilde{u}^{2}(k)\right)\right.$
$\displaystyle\qquad\left.\Theta(1-p)\Theta(1-q)\Theta(|\mathbf{p}-\mathbf{k}|-1)\Theta(|\mathbf{q}+\mathbf{k}|-1)e^{-\left(2k^{2}-2(\mathbf{p}-\mathbf{q})\cdot\mathbf{k}\right)t}\right\\}$
$\displaystyle=$
$\displaystyle-\frac{4\Omega}{(2\pi)^{3}}\int\int_{0}^{\infty}dtd^{3}k\left(\left(G_{0}(k)-k^{2}\tilde{u}(k)\right)^{2}Q^{2}(k,t)-\tilde{u}^{2}(k)A(k,t)\right).$

$Q(k,t)$ and $A(k,t)$ are integrals

$\displaystyle Q(k,t)$ $\displaystyle=$ $\displaystyle\frac{1}{(2\pi)^{3}}\int
d^{3}p\left\\{\exp(-\left(k^{2}-2\mathbf{p}\cdot\mathbf{k}\right)t)\Theta(1-p)\Theta(|\mathbf{p}-\mathbf{k}|-1)\right\\},$
(52) $\displaystyle A(k,t)$ $\displaystyle=$
$\displaystyle\frac{1}{(2\pi)^{6}}\int\int
d^{3}pd^{3}q\left\\{\left(k^{2}-(\mathbf{p}-\mathbf{q})\cdot\mathbf{k}\right)^{2}\right.$
(53)
$\displaystyle\left.\Theta(1-p)\Theta(1-q)\Theta(|\mathbf{p}-\mathbf{k}|-1)\Theta(|\mathbf{q}+\mathbf{k}|-1)e^{-\left(2k^{2}-2(\mathbf{p}-\mathbf{q})\cdot\mathbf{k}\right)t}\right\\}$
$\displaystyle=$ $\displaystyle\frac{1}{4}\frac{\partial^{2}}{\partial
t^{2}}\left(Q^{2}(k,t)\right).$

Detailed calculations dealing with $Q(k,t)$ and $A(k,t)$ are presented in
Appendix. The results are written as

$\displaystyle\int_{0}^{\infty}Q^{2}(k,t)dt$ $\displaystyle=$
$\displaystyle\frac{1-\ln 2}{3(2\pi)^{4}}k+X(k),$ (54)
$\displaystyle\int_{0}^{\infty}A(k,t)dt$ $\displaystyle=$
$\displaystyle\frac{1}{2}k^{2}\left(\frac{1}{6\pi^{2}}\right)^{2}T_{2}(k),$
(55)

where $X(k)$ is defined piecewisely for $k<2$

$\displaystyle X(k)$ $\displaystyle=$
$\displaystyle\frac{1}{2^{6}k(2\pi)^{4}}\frac{1}{5!}\left\\{-704k^{2}-48k^{4}-2048\ln
2\right.$ (56) $\displaystyle+\left[1024+960k-160k^{3}+12k^{5}\right]\ln(2+k)$
$\displaystyle\left.+\left[1024-960k+160k^{3}-12k^{5}\right]\ln(2-k)\right\\},\qquad\mbox{for
}k<2,$

and for $k\geq 2$

$\displaystyle X(k)$ $\displaystyle=$
$\displaystyle\frac{1}{2^{6}k(2\pi)^{4}}\frac{1}{5!}\left\\{128k^{3}+2816k+\left(2k\right)^{3}\left(160-8k^{2}\right)\ln(2k)\right.$
(57)
$\displaystyle+\left(2k-4\right)^{3}\left(\left(2k-4\right)^{2}+40k\right)\ln(2k-4)$
$\displaystyle\left.+\left(2k+4\right)^{3}\left(\left(2k+4\right)^{2}-40k\right)\ln(2k+4)\right\\}-\frac{1-\ln
2}{3(2\pi)^{4}}k,\qquad\mbox{for }k\geq 2.$

By using the expression (54) and (55) in equation (51) we have

$\displaystyle\frac{E_{21}^{c}}{N}$ $\displaystyle=$
$\displaystyle-\frac{4}{\rho^{\prime}(2\pi)^{3}}\int\int_{0}^{\infty}dtd^{3}k\left(\left(G_{0}(k)-k^{2}\tilde{u}(k)\right)^{2}Q^{2}(k,t)-\tilde{u}^{2}(k)A(k,t)\right)$
(58) $\displaystyle=$ $\displaystyle-\frac{3}{2\pi}\int
d^{3}k\left\\{\left(G_{0}(\mathbf{k})-k^{2}\tilde{u}(k)\right)^{2}\left(\frac{1-\ln
2}{3(2\pi)^{4}}k+X(k)\right)\right.$
$\displaystyle\left.-\frac{1}{2}k^{2}\left(\frac{1}{6\pi^{2}}\right)^{2}T_{2}(k)\tilde{u}^{2}(k)\right\\},$

where we find that the last term equals $-E_{11}^{c}/N.$ Therefore we have

$\frac{E_{21}^{c}}{N}=-\frac{1-\ln 2}{(2\pi)^{5}}\int
d^{3}k\left\\{\frac{k^{5}}{T_{2}^{2}(k)}\tilde{u}^{2}(k)\right\\}-\frac{3}{2\pi}\int
d^{3}k\left\\{\frac{k^{4}}{T_{2}^{2}(k)}\tilde{u}^{2}(k)X(k)\right\\}-\frac{E_{11}^{c}}{N},$
(59)

where we have used a simple relation

$G_{0}(k)-k^{2}\tilde{u}(k)=-\frac{k^{2}\tilde{u}(k)}{T_{2}(k)},$ (60)

which can be obtained from equation (49) and (21). As far as the singularity
at $k\sim 0$ is concerned, the first integrand is similar to the integrand of
$E_{11}^{c}$

$\frac{k^{5}}{T_{2}^{2}(k)}\tilde{u}^{2}(k)=\left(\frac{4}{3}\right)^{3}T_{2}(k)k^{2}\tilde{u}^{2}(k)+\tilde{u}^{2}(k)O(k^{4}).$
(61)

We can reformulate the expression of $E_{21}^{c}$as

$\displaystyle\frac{E_{21}^{c}}{N}$ $\displaystyle=$
$\displaystyle-\frac{1-\ln 2}{(2\pi)^{5}}\left(\frac{4}{3}\right)^{3}\int
d^{3}k\left\\{k^{2}T_{2}(k)\tilde{u}^{2}(k)\right\\}$
$\displaystyle-\frac{1-\ln 2}{(2\pi)^{5}}\int
d^{3}k\left\\{\frac{k^{5}}{T_{2}^{2}(k)}\tilde{u}^{2}(k)-\left(\frac{4}{3}\right)^{3}T_{2}(k)k^{2}\tilde{u}^{2}(k)\right\\}$
$\displaystyle-\frac{3}{2\pi}\int
d^{3}k\left\\{\frac{k^{4}}{T_{2}^{2}(k)}\tilde{u}^{2}(k)X(k)\right\\}-\frac{E_{11}^{c}}{N}$
$\displaystyle=$ $\displaystyle\frac{32}{9}(1-\ln
2)\frac{E_{11}^{c}}{N}-\frac{E_{11}^{c}}{N}$ $\displaystyle-\frac{1-\ln
2}{(2\pi)^{5}}\int
d^{3}k\left\\{\frac{k^{5}}{T_{2}^{2}(k)}\tilde{u}^{2}(k)-\left(\frac{4}{3}\right)^{3}T_{2}(k)k^{2}\tilde{u}^{2}(k)\right\\}$
$\displaystyle-\frac{3}{2\pi}\int
d^{3}k\left\\{\frac{k^{4}}{T_{2}^{2}(k)}\tilde{u}^{2}(k)X(k)\right\\},$

where the remaining integrals are not any more singular at $k\sim 0$ in the
high density limit. Then we can use the simple expression (45) and the
integrals can be numerically evaluated as

$\frac{E_{21}^{c}}{N}=\frac{32}{9}(1-\ln
2)\frac{E_{11}^{c}}{N}-\frac{E_{11}^{c}}{N}-0.00273\alpha^{2}r_{s}^{2},$ (62)

and the sum with the first order energy gives

$\frac{E_{1}^{c}+E_{21}^{c}}{N}=\alpha^{2}r_{s}^{2}\left[\frac{1-\ln
2}{\pi^{2}}\ln(r_{s})-0.0645\right].$ (63)

The exchange contraction term of $E_{2}$ can be written as

$\displaystyle E_{22}$ $\displaystyle=$ $\displaystyle
2\frac{1}{4\Omega^{2}}\sum_{\sigma}\sum_{kk^{\prime}p}\frac{G(\mathbf{p}-\mathbf{k},\mathbf{p}-\mathbf{k^{\prime}},\mathbf{k}^{\prime})G(\mathbf{p},\mathbf{p}-\mathbf{k}-\mathbf{k}^{\prime},\mathbf{k})}{\frac{1}{2}\left((\mathbf{p}-\mathbf{k})^{2}+(\mathbf{p}-\mathbf{k}^{\prime})^{2}-p^{2}-(\mathbf{p}-\mathbf{k}-\mathbf{k}^{\prime})^{2}\right)}$
(64)
$\displaystyle\qquad\Theta(1-p)\Theta(1-|\mathbf{p}-\mathbf{k}-\mathbf{k}^{\prime}|)\Theta(|\mathbf{p}-\mathbf{k}|-1)\Theta(|\mathbf{p}-\mathbf{k}^{\prime}|-1)$
$\displaystyle=$
$\displaystyle\frac{N}{2\Omega^{2}}\sum_{kk^{\prime}}\frac{\left(G_{0}(k)-(\mathbf{k}+\mathbf{k}^{\prime})\cdot\mathbf{k}\tilde{u}(k)\right)\left(G_{0}(k^{\prime})-(\mathbf{k}^{\prime}-\mathbf{k})\cdot\mathbf{k}^{\prime}\tilde{u}(k^{\prime})\right)}{-\mathbf{k}\cdot\mathbf{k}^{\prime}}T_{4}(\mathbf{k},\mathbf{k}^{\prime}),$

where

$T_{4}(\mathbf{k},\mathbf{k}^{\prime})\equiv\frac{2}{N}\sum_{p}\Theta(1-p)\Theta(1-|\mathbf{p}-\mathbf{k}-\mathbf{k}^{\prime}|)\Theta(|\mathbf{p}-\mathbf{k}|-1)\Theta(|\mathbf{p}-\mathbf{k}^{\prime}|-1),$
(65)

which can be calculated analytically. Luo (2012) Therefore in the
thermodynamic limit we have

$\displaystyle\frac{E_{22}}{N}$ $\displaystyle=$
$\displaystyle-\frac{1}{2(2\pi)^{6}}\int
d^{3}kd^{3}k^{\prime}\frac{\left(G_{0}(k)-(\mathbf{k}+\mathbf{k}^{\prime})\cdot\mathbf{k}\tilde{u}(k)\right)\left(G_{0}(k^{\prime})-(\mathbf{k}^{\prime}-\mathbf{k})\cdot\mathbf{k}^{\prime}\tilde{u}(k^{\prime})\right)}{\mathbf{k}\cdot\mathbf{k}^{\prime}}T_{4}(\mathbf{k},\mathbf{k}^{\prime})$
(66) $\displaystyle=$
$\displaystyle-\frac{1}{2(2\pi)^{6}}\int_{0}^{\infty}dk\int_{0}^{\infty}dk^{\prime}\int_{0}^{\pi}d\theta\left\\{8\pi^{2}k^{2}k^{\prime
2}\sin\theta T_{4}(k,k^{\prime},\theta)\right.$
$\displaystyle\left.\frac{\left(G_{0}(k)-(k^{2}+kk^{\prime}\cos\theta)\tilde{u}(k)\right)\left(G_{0}(k^{\prime})-(k^{\prime
2}-kk^{\prime}\cos\theta)\tilde{u}(k^{\prime})\right)}{kk^{\prime}\cos\theta}\right\\}$

which can be calculated numerically as

$\frac{E_{22}}{N}=0.01375\alpha^{2}r_{s}^{2}.$ (67)

The total correlation energy is finally estimated as

$\frac{E^{c}}{N}=\alpha^{2}r_{s}^{2}\left[\frac{1-\ln
2}{\pi^{2}}\ln(r_{s})-0.05075\right]+o(r_{s}^{2}).$ (68)

This result has precisely the same logarithmic term as that of the RPA result,
while the constant term is only roughly $7\%$ larger.

##  V Calculations with the Coulomb-Yukawa factor

The above calculations are based on the optimised Jastrow factor for the given
system, which is somewhat complicated for the analytic calculations. Can we
use some simpler factors instead, which still fulfil the asymptotic
conditions?

In quantum Monte Carlo, a Coulomb-Yukawa type function

$u(r)=-\frac{\alpha r_{s}}{f^{2}r}(1-e^{-fr}),\qquad
f^{2}=\sqrt{4\pi\rho\alpha r_{s}},$ (69)

is broadly used as the Jastrow factor. Here the expression is reformulated for
the re-scaled Hamiltonian (11). The Fourier transformation of $u(r)$ has the
following simple form

$\tilde{u}(k)=-\frac{4\pi\alpha
r_{s}}{f^{2}}\left(\frac{1}{k^{2}}-\frac{1}{k^{2}+f^{2}}\right)=-\frac{4\pi\alpha
r_{s}}{k^{2}(k^{2}+f^{2})},$ (70)

where the parameters are determined to fulfil the asymptotic conditions. This
factor is much simpler for the integral calculations and serves as a good
candidate for the test.

In order to get an idea about the performance of this factor in the second
order transcorrelated perturbation calculations, we focus only on the leading
terms of the perturbation energies. The leading term in the first order energy
is given by $E^{c}_{11}$ in equation (33)

$\displaystyle\frac{E_{11}^{c}}{N}$ $\displaystyle=$
$\displaystyle-\frac{\rho}{2(2\pi)^{3}}\int
T_{2}(k)\mathbf{k}^{2}\left(\frac{4\pi\alpha
r_{s}}{k^{2}(k^{2}+f^{2})}\right)^{2}d^{3}k$ (71) $\displaystyle\approx$
$\displaystyle-\frac{1}{6\pi^{2}(2\pi)^{3}}\int_{0}^{2}\frac{3}{4}k^{3}\left(\frac{4\pi\alpha
r_{s}}{k^{2}(k^{2}+f^{2})}\right)^{2}4\pi k^{2}dk$ $\displaystyle=$
$\displaystyle-\frac{\alpha^{2}r_{s}^{2}}{2\sqrt{3}\pi^{2}\alpha^{2}}\frac{1}{\sqrt{r_{s}}}+\cdots,$

which does not lead to a logarithmic term, but rather a more singular
$r^{-\frac{1}{2}}$ term.

According to equation (58), the $E_{11}^{c}$ term will be completely canceled
by the second order energy. The leading term of the total second order
perturbation energy can be evaluated by

$\displaystyle\frac{E_{11}^{c}+E_{21}^{c}}{N}$ $\displaystyle\approx$
$\displaystyle-\frac{3}{2\pi}\int
d^{3}k\left\\{\left(G_{0}(\mathbf{k})-k^{2}\tilde{u}(k)\right)^{2}\frac{1-\ln
2}{3(2\pi)^{4}}k\right\\}$ (72) $\displaystyle=$
$\displaystyle-\frac{3}{2\pi}\int d^{3}k\left\\{\left(\frac{4\pi\alpha
r_{s}}{k^{2}+f^{2}}+\frac{4\pi\alpha
r_{s}f^{2}}{(k^{2}+f^{2})^{2}}\right)^{2}\frac{1-\ln
2}{3(2\pi)^{4}}k\right\\}$ $\displaystyle\approx$
$\displaystyle\alpha^{2}r_{s}^{2}\frac{1-\ln 2}{2\pi}\ln r_{s},$

which reproduces only half of the exact logarithmic term.

This example shows that the perturbation results relies strongly on the
selection of the correlation factor. The current study suggests that a
correlation factor optimised for a single determinant reference function
(e.g., the Hartree-Fock reference function) serves as a good candidate. This
optimisation can be realized by solving either a variational equation, or a
transcorrelated equation for the factor. Luo (2012)

##  VI Conclusions and remarks

By using a proper Jastrow factor, the leading singularity of the Coulomb
potential can be removed from the transcorrelated Hamiltonian. This cures the
divergence problem of the MBPT method for metallic systems at low orders. Test
calculations on uniform electron gases shows that already at the second order
perturbation level, one can get meaningful results for the electron
correlations. We may expect that this method can be developed into an
efficient method for metals, curing problems both at the short-range and the
long-range.

The selection of the Jastrow factor is crucial for the calculations. The
necessary condition that the factor should fulfil the long- and short-range
asymptotic conditions does not guarantee performance. For example, the
Coulomb-Yukawa type factor Foulkes  _et al._ (2001), which is broadly used in
quantum Monte Carlo, does not lead to a meaningful perturbation result. We
also note that this method only cures the divergence at low orders; starting
at 4th order, there are still divergent integrals in the perturbation
energies.

##  Appendix A Calculation of $Q(k,t)$

The function $Q(k,t)$ can be calculated separately for $k\geq 2$ and $k<2$.

###  A.1 Case $k\geq 2$

For $k\geq 2$, the integration volume is a simple unit sphere

$\displaystyle Q(k,t)$ $\displaystyle=$
$\displaystyle\frac{e^{-k^{2}t}}{(2\pi)^{2}}\int_{0}^{1}\int_{0}^{\pi}\left\\{\exp(2pk\cos\theta\
t)p^{2}\sin\theta\right\\}d\theta dp$ (73) $\displaystyle=$
$\displaystyle\frac{e^{-k^{2}t}}{(2\pi)^{2}(2kt)^{3}}\left\\{2kt\left(e^{2kt}+e^{-2kt}\right)-\left(e^{2kt}-e^{-2kt}\right)\right\\}.$

In the calculation of $\int_{0}^{\infty}A(k,t)dt$ integrals we also need to
evaluate the leading orders of $Q$

$\displaystyle Q(k,t)$ $\displaystyle=$
$\displaystyle\frac{e^{-k^{2}t}}{(2\pi)^{2}(2kt)^{3}}\left\\{4kt\left(1+2k^{2}t^{2}+O(t^{4})\right)-\left(4kt+(2kt)^{3}/3+O(t^{5})\right)\right\\}$
$\displaystyle=$
$\displaystyle\frac{e^{-k^{2}t}}{6\pi^{2}}\left\\{1+O(t^{2})\right\\},$

and thus

$\displaystyle Q(k,0)$ $\displaystyle=$ $\displaystyle\frac{1}{6\pi^{2}},$
(74) $\displaystyle\frac{\partial}{\partial t}\left(Q^{2}(k,t)\right)_{t=0}$
$\displaystyle=$ $\displaystyle-2k^{2}\left(\frac{1}{6\pi^{2}}\right)^{2}.$
(75)

###  A.2 case $k<2$

For $k<2$, the integration volume is the unit sphere centered at $\mathbf{0}$
subtracts its intersection with the unit sphere centered at $\mathbf{k}$

$\displaystyle Q(k,t)$ $\displaystyle=$
$\displaystyle\frac{e^{-k^{2}t}}{(2\pi)^{2}(2kt)^{3}}\left\\{2kt\left(e^{2kt}+e^{-2kt}\right)-\left(e^{2kt}-e^{-2kt}\right)\right\\}$
(76)
$\displaystyle-\frac{1}{(2\pi)^{2}}\int_{0}^{\sqrt{1-k^{2}/4}}\int_{-\sqrt{1-r^{2}}+k/2}^{\sqrt{1-r^{2}}-k/2}dzdr\left\\{r\exp(2kzt)\right\\}$
$\displaystyle=$
$\displaystyle\frac{1}{(2\pi)^{2}(2kt)^{3}}\left\\{(2kt+1)\left(e^{-2kt-k^{2}t}-e^{-2kt+k^{2}t}\right)+2k^{2}t\right\\}.$

Again we need the leading terms of the Taylor expansion of $Q(k,t)$ with
respect to $t$, which can be calculated as

$Q(k,t)=\frac{1}{6\pi^{2}}\left\\{T_{2}(k)-k^{2}t+O(t^{2})\right\\}.$ (77)

Then we can get

$\displaystyle Q(k,0)$ $\displaystyle=$
$\displaystyle\frac{1}{6\pi^{2}}T_{2}(k),$ (78)
$\displaystyle\frac{\partial}{\partial t}\left(Q^{2}(k,t)\right)_{t=0}$
$\displaystyle=$
$\displaystyle-2k^{2}\left(\frac{1}{6\pi^{2}}\right)^{2}T_{2}(k).$ (79)

Here we find equation (75) and equation (79) are identical, since for $k>2$,
$T_{2}(k)=1.$

##  Appendix B Calculation of $\int_{0}^{\infty}Q^{2}(k,t)dt$ and
$\int_{0}^{\infty}A(k,t)dt$

###  B.1 case $k\geq 2$

According to equation (73)

$\displaystyle\int_{0}^{\infty}Q^{2}(k,t)dt$ $\displaystyle=$
$\displaystyle\frac{1}{2^{6}(2\pi)^{4}}\int_{0}^{\infty}dt\frac{e^{-2k^{2}t}}{(kt)^{6}}\left(2kt\left(e^{2kt}+e^{-2kt}\right)-\left(e^{2kt}-e^{-2kt}\right)\right)^{2}$
$\displaystyle=$
$\displaystyle\frac{1}{2^{6}k(2\pi)^{4}}\int_{0}^{\infty}dx\frac{e^{-2kx}}{x^{6}}\left(2x\left(e^{2x}+e^{-2x}\right)-\left(e^{2x}-e^{-2x}\right)\right)^{2}.$

This integral has the form $\int_{0}^{\infty}\frac{F(x)}{x^{6}}dx$. Since
$Q(k,t)$ is finite at $t=0,$ the leading term of $F(x)$ has the order $x^{6}$.
We can then use the simple relation

$\int_{0}^{\infty}\frac{F(x)}{x^{6}}dx=\int_{0}^{\infty}\frac{F^{(5)}(x)}{5!x}dx.$
(80)

The calculation with $F^{(5)}$ is lengthy but straightforward, and in the end
we are left with only the following type integrals

$\int_{0}^{\infty}xe^{-ax}dx=a^{-2},\quad\int_{0}^{\infty}e^{-ax}dx=a^{-1},\quad\int_{0}^{\infty}\frac{e^{-ax}-e^{-bx}}{x}dx=\ln\frac{b}{a}.$

The final result of the integral is

$\displaystyle\int_{0}^{\infty}Q^{2}(k,t)dt$ $\displaystyle=$
$\displaystyle\frac{1}{2^{6}k(2\pi)^{4}}\frac{1}{5!}\left\\{128k^{3}+2816k+\left(2k\right)^{3}\left(160-8k^{2}\right)\ln(2k)\right.$
$\displaystyle+\left(2k-4\right)^{3}\left(\left(2k-4\right)^{2}+40k\right)\ln(2k-4)$
$\displaystyle\left.+\left(2k+4\right)^{3}\left(\left(2k+4\right)^{2}-40k\right)\ln(2k+4)\right\\}$

The integral $\int_{0}^{\infty}A(k,t)dt$ can be simply calculated by partial
integration

$\displaystyle\int_{0}^{\infty}A(k,t)dt$ $\displaystyle=$
$\displaystyle\int_{0}^{\infty}dt\left\\{\frac{1}{4}\frac{\partial^{2}}{\partial
t^{2}}\left(Q^{2}(k,t)\right)\right\\}$ $\displaystyle=$
$\displaystyle-\frac{1}{4}\frac{\partial}{\partial
t}\left(Q^{2}(k,t)\right)_{t=0}$ $\displaystyle=$
$\displaystyle\frac{k^{2}}{2}\left(\frac{1}{6\pi^{2}}\right)^{2}.$

###  B.2 case $k<2$

In the same way as above, the integral $\int_{0}^{\infty}Q^{2}(k,t)dt$ can be
calculated as

$\displaystyle\int Q^{2}(k,t)dt$ $\displaystyle=$
$\displaystyle\frac{1}{(2\pi)^{4}}\int_{0}^{\infty}dt\frac{1}{(2kt)^{6}}\left\\{(2kt+1)\left(e^{-2kt-k^{2}t}-e^{-2kt+k^{2}t}\right)+2k^{2}t\right\\}^{2}$
$\displaystyle=$
$\displaystyle\frac{1}{2^{6}k(2\pi)^{4}}\frac{1}{5!}\left\\{1856k^{2}-48k^{4}-\left[2048+2560k^{2}\right]\ln
2\right.$ $\displaystyle+\left[1024+960k-160k^{3}+12k^{5}\right]\ln(2+k)$
$\displaystyle\left.+\left[1024-960k+160k^{3}-12k^{5}\right]\ln(2-k)\right\\}.$

For the calculation of the logarithmic term in $E_{2}^{c}$, we split out the
leading order term of the above integral and write it as

$\int_{0}^{\infty}Q^{2}(k,t)dt=\frac{1-\ln 2}{3(2\pi)^{4}}k+X(k),$

with $X(k)\sim O(k^{2})$ being the high order part

$\displaystyle X(k)$ $\displaystyle=$
$\displaystyle\frac{1}{2^{6}k(2\pi)^{4}}\frac{1}{5!}\left\\{-704k^{2}-48k^{4}-2048\ln
2\right.$ $\displaystyle+\left[1024+960k-160k^{3}+12k^{5}\right]\ln(2+k)$
$\displaystyle\left.+\left[1024-960k+160k^{3}-12k^{5}\right]\ln(2-k)\right\\}.$

The integral on $A$ can also be calculated simply by partial integration

$\displaystyle\int_{0}^{\infty}A(k,t)dt$ $\displaystyle=$
$\displaystyle\int_{0}^{\infty}dt\left\\{\frac{1}{4}\frac{\partial^{2}}{\partial
t^{2}}\left(Q^{2}(k,t)\right)\right\\}$ $\displaystyle=$
$\displaystyle-\frac{1}{4}\frac{\partial}{\partial
t}\left(Q^{2}(k,t)\right)_{t=0}$ $\displaystyle=$
$\displaystyle\frac{k^{2}}{2}\left(\frac{1}{6\pi^{2}}\right)^{2}T_{2}(k).$

## References

  * Kato (1957) T. Kato, “On the eigenfunctions of many-particle systems in quantum mechanics,” Commun. Pure Appl. Math. 10, 151–177 (1957).
  * Bohm and Pines (1951) D. Bohm and D. Pines, “A collective description of electron interactions. I. magnetic interactions,” Phys. Rev. 82, 625–634 (1951).
  * Pines and Bohm (1952) D. Pines and D. Bohm, “A collective description of electron interactions: II. collectivevsindividual particle aspects of the interactions,” Phys. Rev. 85, 338–353 (1952).
  * Bohm and Pines (1953) D. Bohm and D. Pines, “A collective description of electron interactions: III. coulomb interactions in a degenerate electron gas,” Phys. Rev. 92, 609–625 (1953).
  * Gell-Mann and Brueckner (1957) M. Gell-Mann and K. A. Brueckner, “Correlation energy of an electron gas at high density,” Phys. Rev. 106, 364–368 (1957).
  * Ehrenreich and Cohen (1959) H. Ehrenreich and M. H. Cohen, “Self-consistent field approach to the many-electron problem,” Phys. Rev. 115, 786–790 (1959).
  * Jastrow (1955) R. Jastrow, “Many-body problems with strong forces,” Phys. Rev. 98, 1479–1484 (1955).
  * Fournais  _et al._ (2005) S. Fournais, M. Hoffmann-Ostenhof, T. Hoffmann-Ostenhof, and T. Østergaard Sørensen, “Sharp regularity results for coulombic many-electron wave functions,” Commun. Math. Phys. 255, 183–227 (2005).
  * Ceperley (1978) D. Ceperley, “Ground state of the fermion one-component plasma: A Monte Carlo study in two and three dimensions,” Phys. Rev. 18, 3126–3138 (1978).
  * Umrigar, Nightingale, and Runge (1993) C. J. Umrigar, M. P. Nightingale, and K. J. Runge, “A diffusion monte-carlo algorithm with very small time-step errors,” J. Chem. Phys. 99, 2865–2890 (1993).
  * Foulkes  _et al._ (2001) W. M. C. Foulkes, L. Mitas, R. J. Needs, and G. Rajagopal, “Quantum Monte Carlo simulations of solids,” Rev. Mod. Phys. 73, 33–83 (2001).
  * Boys (1969) S. F. Boys, “Some bilinear convergence characteristics of the solutions of dissymmetric secular equations,” Proc. Roy. Soc. A 309, 195–208 (1969).
  * Boys and Handy (1969a) S. F. Boys and N. C. Handy, “A condition to remove the indeterminacy in interelectronic correlation functions,” Proc. Roy. Soc. A 309, 209–220 (1969a).
  * Boys and Handy (1969b) S. F. Boys and N. C. Handy, “The determination of energies and wavefunctions with full electronic correlation,” Proc. Roy. Soc. A 310, 43–61 (1969b).
  * Boys and Handy (1969c) S. F. Boys and N. C. Handy, “A calculation for the energies and wavefunctions for states of neon with full electronic correlation accuracy,” Proc. Roy. Soc. A 310, 63–78 (1969c).
  * Boys and Handy (1969d) S. F. Boys and N. C. Handy, “A first solution, for lih, of a molecular transcorrelated wave equation by means of restricted numerical integration,” Proc. Roy. Soc. A 311, 309–329 (1969d).
  * Handy (1969) N. C. Handy, “Energies and expectation values for be by transcorrelated method,” J. Chem. Phys. 51, 3205–3212 (1969).
  * Handy (1972) N. C. Handy, “The transcorrelated method for accurate correlation energies using gaussian-type functions,” Mol. Phys. 23, 1–27 (1972).
  * Handy (1973) N. C. Handy, “Towards and understanding of the form of correlated wavefunctions for atoms,” J. Chem. Phys. 58, 279–287 (1973).
  * Ten-no (2000) S. Ten-no, “A feasible transcorrelated method for treating electronic cusps using a fronzen Gaussina geminal,” Chem. Phys. Lett. 330, 169–174 (2000).
  * Hino, Tanimura, and Ten-no (2001) O. Hino, Y. Tanimura, and S. Ten-no, “Biorthogonal approach for explicitly correlated calculations using the transcorrelated hamiltonian,” J. Chem. Phys. 115, 7865–7871 (2001).
  * Hino, Tanimura, and Ten-no (2002) O. Hino, Y. Tanimura, and S. Ten-no, “Application of the transcorrelated hamiltonian to the linearized coupled cluster singles and doubles model,” Chem. Phys. Lett. 353, 317 (2002).
  * Luo (2011) H. Luo, “Complete optimisation of multi-configuration jastrow wave functions by variational transcorrelated method,” J. Chem. Phys. 135, 024109 (2011).
  * Luo and Alavi (2018) H. Luo and A. Alavi, “Combining the transcorrelated method with full configuration interaction quantum monte carlo: Application to the homogeneous electron gas,” J. Chem. Theory Comput. 14, 1403–1411 (2018).
  * Cohen  _et al._ (2019) A. J. Cohen, H. Luo, K. Guther, W. Dobrautz, D. P. Tew, and A. Alavi, “Similarity transformation of the electronic schrödinger equation via jastrow factorization,” J. Chem. Phys. 151, 061101 (2019).
  * Guther  _et al._ (2021) K. Guther, A. J. Cohen, H. Luo, and A. Alavi, “Binding curve of the beryllium dimer using similarity-transformed fciqmc: Spectroscopic accuracy with triple-zeta basis sets,” The Journal of Chemical Physics 155, 011102 (2021).
  * Liao  _et al._ (2021) K. Liao, T. Schraivogel, H. Luo, D. Kats, and A. Alavi, “Towards efficient and accurate ab initio solutions to periodic systems via transcorrelation and coupled cluster theory,” Phys. Rev. Research 3, 033072 (2021).
  * Schraivogel  _et al._ (2021) T. Schraivogel, A. J. Cohen, A. Alavi, and D. Kats, “Transcorrelated coupled cluster methods,” J. Chem. Phys. 155, 191101 (2021).
  * Armour (1980) E. A. G. Armour, “The calculation of the ground-state energy of the free-electon gas by the transcorrelated method,” J. Phys. C 13, 343–348 (1980).
  * N. Umezawa and S. Tsuneyuki (2004) N. Umezawa and S. Tsuneyuki, “Ground-state correlation energy for the homogeneous electron gas calculated by the transcorrelated method,” Phys. Rev. B 69, 165102 (2004).
  * Luo (2012) H. Luo, “Transcorrelated calculations of homogeneous electron gases,” J. Chem. Phys. 136, 224111 (2012).
  * Onsager, Mittag, and Stephen (1966) L. Onsager, L. Mittag, and M. J. Stephen, “Integrals in the theory of electron correlations,” Annalen der Physik 473, 71–77 (1966).
  * Gaskell (1961) T. Gaskell, “The collective treatment of a fermi gas: II,” Proc. Phys. Soc. 77, 1182–1192 (1961).
  * Gaskell (1962) T. Gaskell, “The collective treatment of many-body systems: III,” Proc. Phys. Soc. 80, 1091–1100 (1962).
  * Talman (1974) J. Talman, “Linked-cluster expansion for Jastrow-type wave functions and its application to the electron-gas problem,” Phys. Rev A 10, 1333–1344 (1974).