# Towards efficient and accurate _ab initio_ solutions to periodic systems via
transcorrelation and coupled cluster theory

Ke Liao<EMAIL_ADDRESS>Max Planck Institute for Solid State Research,
Heisenbergstrasse 1, 70569 Stuttgart, Germany  Thomas Schraivogel  Max Planck
Institute for Solid State Research, Heisenbergstrasse 1, 70569 Stuttgart,
Germany  Hongjun Luo  Max Planck Institute for Solid State Research,
Heisenbergstrasse 1, 70569 Stuttgart, Germany  Daniel Kats<EMAIL_ADDRESS>Max Planck Institute for Solid State Research, Heisenbergstrasse 1, 70569
Stuttgart, Germany  Ali Alavi<EMAIL_ADDRESS>Max Planck Institute for
Solid State Research, Heisenbergstrasse 1, 70569 Stuttgart, Germany  Yusuf
Hamied Department of Chemistry, University of Cambridge, Lensfield Road,
Cambridge CB2 1EW, United Kingdom

###### Abstract

We propose a streamlined combination scheme of the transcorrelation (TC) and
coupled cluster (CC) theory, which not only increases the convergence rate
with respect to the basis set, but also extends the applicability of the
lowest order CC approximations to strongly correlated regimes in the three
dimensional uniform electron gas (3D UEG). With the correct physical insights
built into the correlator used in TC, highly accurate ground state energies
with errors $\leq 0.001$ a.u./electron relative to the state-of-the-art
quantum Monte Carlo results can be obtained across a wide range of densities.
The greatly improved efficiency and accuracy of our methods hold great promise
for strongly correlated solids where many other methods fail.

##  I Introduction

The coupled cluster (CC) methodologies Čížek (1966); Čížek and Paldus (1971,
1980) at the level of singles and doubles (CCSD) and perturbative triples
(CCSD(T)) Raghavachari  _et al._ (1989) have become the _de facto_ standard of
single-reference _ab initio_ quantum chemistry, and can be applied to systems
consisting of hundreds of electrons Riplinger  _et al._ (2016); Schmitz and
Hättig (2016); Schwilk  _et al._ (2017); Nagy  _et al._ (2018). In the past
few years, these methods have also shown promise in applications to the solid
state Booth  _et al._ (2013); Liao and Grüneis (2016); McClain  _et al._
(2017); Gruber  _et al._ (2018); Liao  _et al._ (2019); Gao  _et al._ (2020),
although significant challenges remain before they can be routinely applied,
as for example density functional theories are. On the one hand, because of
quite steep computational scaling (${\cal O}(N^{6})$ and ${\cal O}(N^{7})$ for
CCSD and CCSD(T) respectively), it is desirable to keep the methods at the
lowest possible CC level, namely CCSD, whilst maintaining accuracy. The more
accurate CCSD(T), as a perturbative correction to CCSD, additionally fails for
metals Shepherd and Grüneis (2013). It is also desirable that the CC methods
can be extended to more strongly correlated systems, where the single
reference nature of these approximations breaks down. There have been various
attempts to develop modified CCSD methods with a higher accuracy for weakly
Meyer (1971); Nooijen and Le Roy (2006); Neese  _et al._ (2009); Huntington
and Nooijen (2010) and strongly Paldus  _et al._ (1984); Piecuch and Paldus
(1991); Bartlett and Musiał (2006); Robinson and Knowles (2011); Small  _et
al._ (2014) correlated systems. The distinguishable cluster singles and
doubles (DCSD) Kats and Manby (2013a); Kats (2014) is one such method, which
has shown promise in improving CCSD in weakly and strongly correlated
molecular systems Kats  _et al._ (2015); Rishi  _et al._ (2017); Tsatsoulis
_et al._ (2017).

In a separate development, there has been renewed interest in so-called
transcorrelated (TC) methods Boys and Linnett (1969); Boys  _et al._ (1969);
Boys and Handy (1969); Hino  _et al._ (2001, 2002); Umezawa and Tsuneyuki
(2004); Sakuma and Tsuneyuki (2006); Ochi and Tsuneyuki (2015); Ochi  _et al._
(2016, 2017); Luo and Alavi (2018); Jeszenszki  _et al._ (2018); Dobrautz  _et
al._ (2019); Cohen  _et al._ (2019); Jeszenszki  _et al._ (2020), based on
Jastrow factorisation of the electronic wavefunction, which result in
effective similarity transformed (ST) Hamiltonians Luo and Alavi (2018); Cohen
_et al._ (2019). Although TC methods were originally proposed as a way to
accelerate basis set convergence in electronic wavefunctions, it has become
apparent that such similarity transformations can also be extremely helpful in
the context of strongly correlated systems. For example, in the repulsive 2D
Fermi-Hubbard model, it was found that, with a suitable Gutzwiller correlator,
extremely compact forms of ground state right eigenvectors of the ST Hubbard
Hamiltonian could be obtained Dobrautz  _et al._ (2019), dominated by the
Hartree-Fock (HF) determinant. Since single-reference coupled cluster methods
work best when the wavefunction is dominated by the HF determinant, and
furthermore, since the coupled cluster method can itself be cast in terms of a
similarity transformation of the Hamiltonian, it is natural to ask if the two
concepts - coupled cluster and transcorrelation - can be usefully combined
into a single framework, whereby the compactification generated by the TC
method is exploited by the CC method, to extend its range of applicability
into more strongly correlated systems. The purpose of this paper is to report
such an investigation, applied to the uniform electron gas (UEG), over a broad
range of densities.

The 3D UEG model assumes that the background is evenly and positively charged,
and that the electrons interact with each other via the Coulomb interaction.
As simple as it is, the UEG possesses an intricate phase diagram Wigner
(1934); Ceperley and Alder (1980), which can only be accurately described by
theories that perform consistently well over a broad range of densities.
Historically, the UEG model has also played an important role in the
development of many useful approximations. For example, several successful
local and gradient-corrected density functionals Kohn and Sham (1965); Perdew
and Zunger (1981); Perdew  _et al._ (1996) are based on the UEG; the random
phase approximation (RPA) Bohm and Pines (1953); Gell-Mann and Brueckner
(1957) was developed in a pursuit of understanding metals using the UEG as a
model. In recent years, the UEG has attracted studies from various highly
accurate _ab initio_ methods and spurred the development of several new
methods Umezawa and Tsuneyuki (2004); Shepherd  _et al._ (2012a, b); Shepherd
and Grüneis (2013); Shepherd  _et al._ (2014); Neufeld and Thom (2017);
Ruggeri  _et al._ (2018); Luo and Alavi (2018); Irmler  _et al._ (2019).

When applying CC to the UEG, we work in a plane wave basis; momentum
conservation then excludes all single excitations from the CC ansatz, greatly
simplifying the resulting amplitude equations. As a result, the TC Hamiltonian
can be treated with relative ease, allowing us to investigate whether the CC
method can be beneficially applied to the TC Hamiltonian. We will investigate
the CCD and DCD approximations, in the context of the TC Hamiltonian, and show
that with a suitable form of the correlator, the basis set convergence can be
greatly accelerated (as expected), but in addition highly accurate energies
can be obtained across a broad range of densities $0.5\leq r_{s}\leq 50$,
covering both the weakly and strongly correlated regimes. This gives us
confidence that the method, once suitably generalised to real systems (which
will need to include the singles contribution), will allow a highly accurate
yet efficient methodology for the solid state.

In the rest part of this paper, we review the UEG model, CC and TC theories in
Sec. II.1, II.2 and II.3, respectively; in Sec. II.4 we discuss the important
approximations made to the TC Hamiltonian; we demonstrate our scheme for
choosing the optimal parameters in the correlator in Sec. III; we showcase and
discuss our TC-CCD/DCD results in comparison with benchmark data in Sec. IV;
and finally we conclude the paper in Sec. V with some outlooks for future
directions.

##  II Theory

###  II.1 Three Dimensional Uniform Electron Gas

The 3D UEG is the simplest model for realistic periodic solids, of which the
Hamiltonian in real space reads

$\hat{H}=-\sum_{i}\frac{1}{2}\nabla^{2}_{i}+\sum_{i\neq j}\frac{1}{|{\bf
r_{i}-r_{j}}|}+{\rm const.},$ (1)

where the const. includes the interactions between electrons and the
homogenous positive background charge, and the interactions between the
electrons and their own periodic images, which is termed as the Madelung
constant and will disappear as the size of the simulation cell goes to
infinity. Atomic units are used to simplify the equations. When plane wave
basis functions and a simple cubic simulation cell of volume $\Omega=L^{3}$
are used, we can reformulate the Hamiltonian in a second-quantized form,

$\hat{H}=\sum_{p}\sum_{\sigma}\frac{1}{2}{\bf
k}_{p}^{2}a^{\dagger}_{p,\sigma}a_{p,\sigma}+\frac{1}{2}\sum_{pqrs}\sum_{\sigma\sigma^{\prime}}V^{rs}_{pq}a^{\dagger}_{{p},\sigma}a^{\dagger}_{{q},\sigma^{\prime}}a_{{s},\sigma^{\prime}}a_{{r},\sigma},$
(2)

where for simplicity we use $p,q,r,s\dots$ indices as a compact form for the
general momentum (plane wave basis function) indices ${\bf k}_{p},{\bf
k}_{q},{\bf k}_{r},{\bf k}_{s}\dots$ and hereon we use the two terms plane
wave basis function and orbital equivalently. We stress that due to momentum
conservation, i.e. ${\bf k}\equiv{\bf k}_{r}-{\bf k}_{p}={\bf k}_{q}-{\bf
k}_{s}$, there are only three free indices among $pqrs$, and the interactions
with the homogenous positive background charge are cancelled by the divergent
Coulomb potential at $k=0$, which is defined as $V^{rs}_{pq}\equiv V({\bf
k})=\frac{4\pi}{\Omega{\bf k}^{2}}$. We also ignore the Madelung contribution
in the Hamiltonian which can be added posteriorly to the ground state energy.
The electron density of the system can be described by the Wigner-Seitz radius
$r_{s}=\left(\frac{3}{4\pi N}\right)^{1/3}L$, where $N$ is the number of
electrons.

###  II.2 Coupled/Distinguishable Cluster Doubles

In the CC ansatz, we let the many-electron ground state wavefunction to be

$\Phi=e^{\hat{T}}\Phi_{0}$ (3)

where $\Phi_{0}$ is the HF wavefunction and $\hat{T}$ is a cluster operator.
In the case of the UEG, we work in a plane wave basis and $\Phi_{0}$ is given
by the Fermi sphere. We will investigate the CCD approximation and its
distinguishable variant (DCD) Kats and Manby (2013b, a), which is based on a
modification of the CCD amplitude equations by neglecting inter-cluster
exchange diagrams and ensuring the particle-hole symmetry and exactness for
two electrons. Alternatively, DCD can be derived from screened Coulomb
considerations Kats (2016). We start with the canonical CCD and later
highlight the differences between CCD and DCD.

In CCD, the full cluster operator is approximated by the doubles excitations
only,

$\hat{T}\approx\hat{T}_{2}=\frac{1}{2}\sum_{ijab}T_{ab}^{ij}\sum_{\sigma\sigma^{\prime}}a^{\dagger}_{{a},\sigma}a^{\dagger}_{{b},\sigma^{\prime}}a_{{j},\sigma^{\prime}}a_{{i},\sigma},$
(4)

where $T^{ij}_{ab}$ are the doubles amplitudes. Following convention, we use
$i,j,k\dots$ and $a,b,c\dots$ to represent occupied and unoccupied orbitals in
the reference determinant, respectively. Again, the momentum conservation
ensures that only 3 indices of the amplitude tensor are free, saving a great
deal in storing them in the computer memory.

The $T_{2}$ amplitudes are obtained by solving the projective doubles
amplitude equations,

$\left<\Phi_{ij}^{ab}\left|e^{-\hat{T}_{2}}\hat{H}e^{\hat{T}_{2}}\right|\Phi_{0}\right>=0,$
(5)

where $\Phi_{ij}^{ab}$ are doubly substituted determinants. To be specific, a
functional form of the residual, which unifies CCD and DCD for closed shell
systems, can be written as

$\displaystyle R^{ij}_{ab}=$ $\displaystyle
V_{ab}^{ij}+V_{ab}^{cd}T_{cd}^{ij}+I_{kl}^{ij}T_{ab}^{kl}+X_{al}^{cj}T_{cb}^{il}+\tilde{T}_{ac}^{ik}V_{kl}^{cd}\tilde{T}_{db}^{lj}$
(6)
$\displaystyle+\mathcal{\hat{P}}(ia;jb)\left[x_{a}^{c}T_{cb}^{ij}-x_{k}^{i}T_{ab}^{kj}+\chi_{al}^{ci}(T_{bc}^{lj}-T_{cb}^{lj})\right.$
$\displaystyle\left.-V_{ka}^{ic}T_{cb}^{kj}-V_{kb}^{ic}T_{ac}^{kj}+\tilde{T}_{ac}^{ik}V_{kb}^{cj}\right],$

where we define the permutation operator
$\mathcal{\hat{P}}(ia;jb)T_{ab}^{ij}\equiv T_{ab}^{ij}+T_{ba}^{ji}$ and the
following intermediates,

$\displaystyle I_{kl}^{ij}$
$\displaystyle=\begin{cases}V_{kl}^{ij}+V_{kl}^{cd}T_{cd}^{ij},&\quad\text{CCD}\\\
V_{kl}^{ij},&\quad\text{DCD}\end{cases}$ (7) $\displaystyle X_{al}^{cj}$
$\displaystyle=\begin{cases}V_{kl}^{cd}T_{ad}^{kj},&\quad\text{CCD}\\\
0,&\quad\text{DCD}\end{cases}$ (8) $\displaystyle x_{a}^{c}$
$\displaystyle=\begin{cases}f_{a}^{c}-\tilde{T}_{ad}^{kl}V_{lk}^{dc},&\quad\text{CCD}\\\
f_{a}^{c}-\frac{1}{2}\tilde{T}_{ad}^{kl}V_{lk}^{dc},&\quad\text{DCD}\end{cases}$
(9) $\displaystyle x_{k}^{i}$
$\displaystyle=\begin{cases}f_{k}^{i}+\tilde{T}_{cd}^{il}V_{lk}^{dc},&\quad\text{CCD}\\\
f_{k}^{i}+\frac{1}{2}\tilde{T}_{cd}^{il}V_{lk}^{dc},&\quad\text{DCD}\end{cases}$
(10) $\displaystyle\chi_{al}^{ci}$
$\displaystyle=\begin{cases}V_{kl}^{cd}T_{da}^{ki},&\quad\text{CCD}\\\
0,&\quad\text{DCD}\end{cases}$ (11) $\displaystyle\tilde{T}_{ab}^{ij}$
$\displaystyle=2T_{ab}^{ij}-T_{ba}^{ij}.$ (12)

We note that in this case the Fock matrix $f^{p}_{q}$ is diagonal, with the
diagonal elements being the orbital energies $\epsilon_{p}$. A straightforward
way to update the $T_{2}$ amplitudes at iteration $n+1$ will be

$\displaystyle\Delta_{ab}^{ij}$
$\displaystyle=\frac{R_{ab}^{ij}}{\epsilon_{i}+\epsilon_{j}-\epsilon_{a}-\epsilon_{b}},$
(13) $\displaystyle T_{ab}^{ij}(n+1)$
$\displaystyle=T_{ab}^{ij}(n)+\Delta_{ab}^{ij}.$ (14)

Of course, more advanced iterative schemes can be used, e.g. DIIS Pulay (1980,
1982), to accelerate convergence rate.

Using the converged $T_{2}$ amplitudes, the correlation energy is expressed as

$E_{\rm c}=\tilde{T}_{ab}^{ij}V_{ij}^{ab},$ (15)

and the total energy is

$E=E_{\rm HF}+E_{\rm c},$ (16)

where $E_{\rm HF}=\left<\Phi_{0}\right|\hat{H}\left|\Phi_{0}\right>$.

###  II.3 Transcorrelation

In the transcorrelation framework the many-electron wavefunction is written as

$\Psi=e^{\hat{\tau}}\Phi$ (17)

where $\hat{\tau}=\frac{1}{2}\sum_{i\neq j}u({\bf r_{i},r_{j}})$ is a
correlator consisting of pair correlations $u({\bf r_{i},r_{j}})$, whose form
will be discussed later. $\Phi$ should satisfy the similarity-transformed
eigenvalue equation

$\hat{H}_{\rm tc}\Phi=E\Phi,\quad\hat{H}_{\rm
tc}=e^{-\hat{\tau}}\hat{H}e^{\hat{\tau}}.$ (18)

It is worth pointing out that at this stage, no approximations have been made,
and the spectra $E$ of $\hat{H}_{\rm tc}$ are the same as of the original
Hamiltonian. As shown in Ref. Luo and Alavi (2018), the second quantized form
of the $\hat{H}_{\rm tc}$ in a plane wave basis is

$\displaystyle\hat{H}_{\rm tc}$
$\displaystyle=\hat{H}+\frac{1}{2}\sum_{\sigma\sigma^{\prime}}\sum_{pqrs}{\omega}_{pq}^{rs}a^{\dagger}_{{p},\sigma}a^{\dagger}_{{q},\sigma^{\prime}}a_{{s},\sigma^{\prime}}a_{{r},\sigma}$
(19)
$\displaystyle+\frac{1}{2}\sum_{\sigma\sigma^{\prime}\sigma^{\prime\prime}}\sum_{pqorst}{\omega}^{rst}_{pqo}a^{\dagger}_{{p},\sigma}a^{\dagger}_{{q},\sigma^{\prime}}a^{\dagger}_{{o},\sigma^{\prime\prime}}a_{{t},\sigma^{\prime\prime}}a_{{s},\sigma^{\prime}}a_{{r},\sigma},$

where momentum conservation requires ${\bf k}\equiv{\bf k}_{r}-{\bf k}_{p}$,
${\bf k^{\prime}}\equiv{\bf k}_{q}-{\bf k}_{s}$ and ${\bf k}_{o}={\bf
k}_{t}+{\bf k-k^{\prime}}$, and we define

$\displaystyle{\omega}^{rs}_{pq}$
$\displaystyle=\frac{1}{\Omega}\left[k^{2}\tilde{u}({\bf k})-({\bf k}_{r}-{\bf
k}_{s})\cdot{\bf k}\tilde{u}({\bf k})\right]$ (20)
$\displaystyle\quad+\frac{1}{\Omega}\sum_{\bf k^{\prime}}({\bf
k-k^{\prime}})\cdot{\bf k^{\prime}}\tilde{u}({\bf k-k^{\prime}})\tilde{u}({\bf
k^{\prime}}),$ $\displaystyle{\omega}^{rst}_{pqo}$
$\displaystyle=-\frac{1}{\Omega^{2}}\tilde{u}({\bf k})\tilde{u}({\bf
k^{\prime}}){\bf k}\cdot{\bf k^{\prime}}.$ (21)

The TC Hamiltonian has additional 2-body and 3-body interactions. Due to one
of the additional 2-body interactions, the TC Hamiltonian is non-hermitian.
This fact can pose some difficulties for variational methods, but not so for
projection methods such as full configuration interaction Monte Carlo (FCIQMC)
Booth  _et al._ (2009, 2013) and CC.

###  II.4 Approximations to the 3-body operator

The additional 3-body operator when treated without approximations will
increase the computational scaling of CCD or DCD from $N^{6}$ to $N^{7}$. To
seek a good balance between the computational cost and the accuracy, we
include only up to effective 2-body operators arising from normal-ordering the
3-body operator. In this approximation, only the normal-ordered 3-body
operator is excluded. We can justify this approximation by analogy to the HF
approximation, which constructs a mean-field solution by including only the
single and double contractions from the Coulomb operator. In cases where the
mean-field approximation is reasonably good, the contribution of the missing
normal-ordered Coulomb operator is small, compared to the single and double
contractions. In contrast to the HF approximation, the parameters in the
correlator in general allow a tuning of the strength of the missing normal-
ordered 3-body operator, which we will discuss in the next section.

In general, we can write our approximated Hamiltonian as

$\displaystyle\hat{H}_{\rm tc}$ $\displaystyle=\tilde{E}_{\rm
HF}+\sum_{\sigma}\sum_{p}\tilde{\epsilon}_{p}\\{a^{\dagger}_{{p},\sigma}a_{{p},\sigma}\\}$
(22)
$\displaystyle+\frac{1}{2}\sum_{\sigma\sigma^{\prime}}\sum_{pqrs}\tilde{V}^{rs}_{pq}\\{a^{\dagger}_{{p},\sigma}a^{\dagger}_{{q},\sigma^{\prime}}a_{{s},\sigma^{\prime}}a_{{r},\sigma}\\}$
$\displaystyle+E_{\rm
T}+\sum_{\sigma}\sum_{p}\tilde{\omega}_{p}\\{a^{\dagger}_{{p},\sigma}a_{{p},\sigma}\\}$
$\displaystyle+\frac{1}{2}\sum_{\sigma\sigma^{\prime}}\sum_{pqrs}\tilde{w}^{rs}_{pq}\\{a^{\dagger}_{{p},\sigma}a^{\dagger}_{{q},\sigma^{\prime}}a_{{s},\sigma^{\prime}}a_{{r},\sigma}\\},$

where $E_{T}$ refers to the triply-contracted 3-body operator contribution,
$\tilde{\omega}_{p}$ is the doubly-contracted 3-body integral and
$\tilde{\omega}^{rs}_{pq}$ is the singly-contracted 3-body integral. The curly
brackets indicate that the operators are normal-ordered with respect to the HF
vacuum (Fermi sphere). We emphasize that in Eq. (22) $\tilde{E}_{\rm HF}$ and
$\tilde{\epsilon}_{p}$ are calculated now with the modified 2-body integrals
$\tilde{V}_{pq}^{rs}=w^{rs}_{pq}+V^{rs}_{pq}$. For clarity, we outline the
procedures for our TC-CCD/DCD framework.

  1. 1.

Evaluating $\omega^{rs}_{pq}$ and $V^{rs}_{pq}$ and combining them into
$\tilde{V}_{pq}^{rs}\leftarrow w^{rs}_{pq}+V^{rs}_{pq}$;

  2. 2.

Calculating $\tilde{\epsilon}_{p}=\frac{{\bf
k}_{p}^{2}}{2}+\sum_{i}(2\tilde{V}^{pi}_{pi}-\tilde{V}^{pi}_{ip})$;

  3. 3.

Calculating $\tilde{E}_{\rm
HF}=2\sum_{i}^{N/2}\tilde{\epsilon}_{i}-\sum_{ij}(2\tilde{V}^{ij}_{ij}-\tilde{V}^{ij}_{ji})$
and $E_{\rm T}$;

  4. 4.

Evaluating $\tilde{\omega}_{p}$, and defining
$\epsilon_{p}\leftarrow\tilde{\epsilon}_{p}+\tilde{\omega}_{p}$;

  5. 5.

Evaluating the singly-contracted 3-body integral $\tilde{w}^{rs}_{pq}$ and
redefining $V^{rs}_{pq}\leftarrow\tilde{V}_{pq}^{rs}+\tilde{w}^{rs}_{pq}$;

  6. 6.

Solving the usual CCD/DCD amplitude equations using $\epsilon_{p}$ and
$V^{rs}_{pq}$ for $T_{2}$ and obtaining $E_{\rm c}$;

  7. 7.

The total energy is $E=\tilde{E}_{\rm HF}+E_{\rm T}+E_{\rm c}$.

For details on the mathematical expressions for the contractions of the 3-body
operator, we refer to the Supplementary Information.

##  III Choice and optimization of the correlator

Past experience with the TC method has shown that the form of the correlator
$\hat{\tau}$ is of extreme importance in the TC method, otherwise the benefit
of transcorrelation is lost – $\Phi$ can be simpler than $\Psi$ only if the
correlator captures the correct physics of the pair correlations in the
system. An inappropriate correlator can actually lead to a harder problem than
the original Schrödinger equation. In our previous study of the exact TC
method in the UEG, we proposed a form of the correlator (shown below) which
was found to work successfully in accelerating convergence to the basis set
limit, without changing the correlation that could be captured with the basis
set by a FCI level $\Phi$ function. In the present study, since we will be
approximating the $\Phi$ with the CC ansatz, we additionally require the
correlator $\hat{\tau}$ to capture some of the correlation inside the Hilbert
space.

Here we propose a physically motivated correlator that mimics the behavior of
the correlation hole between two unlike-spin electrons as $r_{s}$ varies in 3D
UEG. The correlation hole can be examined by the pair-correlation function
$g(r_{12})$ in real space, as studied in Ref. Gori-Giorgi  _et al._ (2000),
which shows that the correlation hole between two unlike-spin electrons grows
deeper and wider as the Wigner-Seitz radius $r_{s}$ increases or as the
electron density decreases. Fig. 1 provides a sketch of the Jastrow factor
with our proposed correlator $u(r_{12})$ as the exponent, which captures the
desired behavior.

Figure 1: A sketch of the Jastrow factor with the proposed correlator as the
exponent.

We point out that the functional form of this correlator, which reads in real
and reciprocal space respectively as

$\displaystyle u(r)$ $\displaystyle=-\frac{r}{\pi}\left({\rm
si}(k_{c}r)+\frac{{\rm cos}(k_{c}r)}{k_{c}r}+\frac{{\rm
sin}(k_{c}r)}{(k_{c}r)^{2}}\right),$ (23) $\displaystyle\tilde{u}({\bf k})$
$\displaystyle=\begin{cases}-\frac{4\pi}{k^{4}},\quad|{\bf k}|>k_{c},\\\
0,\quad|{\bf k}|\leq k_{c},\end{cases}$ (24)

where ${\rm si}(x)=-\int_{x}^{\infty}\frac{\sin(t)}{t}\mathrm{d}t$, was first
reported in Ref. Luo and Alavi (2018) to satisfy the cusp condition between
two electrons with opposite spins at short inter-electron distance and its
influence is reduced to nonexistence as the complete basis set (CBS) limit is
reached. This was done by choosing $k_{c}$ to be the same as the plane wave
cutoff momentum, $k_{\rm F}$, which defines how many plane waves are included
as basis functions. In contrast, to mimic the behavior of the pair-correlation
function, as a first attempt in the present study we choose the parameter in
this correlator such that the first nonzero root of Eq. (23) is fixed to be at
$r_{s}$, irrespective of the basis set. This is achieved by setting

$k_{c}=\frac{R_{1}}{r_{s}},$ (25)

where $R_{1}\approx 2.322502989$. This choice can be rationalized by the
physical picture that at lower densities, electrons prefer to stay further
away from each other. Furthermore, this correlator, regardless of the choice
of $k_{c}$, retains the cusp condition for two electrons with unlike-spins at
$r=0$ Luo and Alavi (2018) and should increase the convergence rate of the
computed energies with respect to the employed basis set towards the CBS
limit.

To further justify the choice of this correlator, we show that for UEG with 14
electrons at $r_{s}=5$, where traditional CCD exhibits a large error, the most
compact expansion of the wavefunction in Slater determinant space is reached
at this value of $k_{c}$. In Fig. 2, we show the weights of the HF determinant
extracted from TC-FCIQMC simulations using different $k_{c}$ values, without
making approximations to the 3-body operator. We note that due to the discrete
momentum mesh as a result of using a finite simulation cell, the possible
choices are $k_{c}=\frac{2\pi\sqrt{n}}{L}$, $n\in\mathbb{N}$
111$n=n_{1}^{2}+n_{2}^{3}+n_{3}^{2}$, where $n_{1},n_{2},n_{3}$ are the
components of a ${\bf k}$ vector in units of $\frac{2\pi}{L}$. So the possible
values of n are 0, 1, 2, 3, 4, 5, 6, 8, …, where $L$ is the length of the
cubic cell. In this case, $k_{c}=\frac{R_{1}}{r_{s}}$ is equivalent to
$k_{c}=\frac{2\sqrt{2}\pi}{L}$, and for this choice of $k_{c}$ the exact
ground state wavefunction of the transcorrelated Hamiltonian has the highest
weight on the HF determinant.

Figure 2:  The weights of the HF determinant (Fermi sphere) as a function of
$k_{c}$ extracted from their corresponding TC-FCIQMC simulations (red dashed
line), and that of the normal FCIQMC simulation (horizontal red solid line).
The system consists of 14 electrons with $r_{s}=5$ and a basis set including
57 plane waves. $5\times 10^{8}$ walkers are used in all simulations and the
initiator threshold is set to 3. No approximations are made to the 3-body
interactions.

However, we find that this intuitive choice of $k_{c}$ is not always the
optimal, especially at extremely low density regimes. It is reasonable to
expect that the optimal $k_{c}$ for those systems should deviate slightly from
$\frac{R_{1}}{r_{s}}$. So we scan a range of $k_{c}$ values around it to
locate the one that minimizes the norm of the closed-shell amplitudes for
double excitations of electrons with opposite spins,
$\|T_{2}^{\downarrow\uparrow}\|$, in the TC-CCD/DCD calculations with a small
basis set, see Fig. 3 222Since the $\|T_{2}^{\downarrow\uparrow}\|$ in TC-CCD
and TC-DCD show the same trend as a function of $k_{c}$, we show only the
latter in the figure..

Ideally, two separate correlators should be used for electrons with parallel
and anti-parallel spins, and their parameters should be optimized
simultaneously using the norm of the full amplitudes in a similar manner. For
the present study, we argue that the correlations between two parallel-spin
electrons are dominated by the exchange effects, which are already captured by
the anti-symmetry in the Slater determinants. Therefore, we focus on capturing
the correct physics between electrons with opposite spins in the correlator,
i.e. the changing depth and width of the correlation hole as a function of
$r_{s}$ Gori-Giorgi  _et al._ (2000), and minimizing the corresponding
amplitudes in the CC ansätze. Indeed, we found in practice the minima in
$\|T_{2}^{\downarrow\uparrow}\|$ as a function of $k_{c}$ are more pronounced,
and thus easier to spot than those in the norm of the full amplitudes,
$\|T_{2}\|$. We stress that this compact form of wavefunction at the optimal
choice of $k_{c}$ should greatly benefit approximate methods like CCD and DCD,
whose accuracy relies on the assumption that the true ground state
wavefunction is compact around the reference determinant, which is normally
chosen to be the HF determinant.

(a) (b)

Figure 3:  The norm of the amplitudes for excitations of unlike-spin electron
pairs, $\|T_{2}^{\downarrow\uparrow}\|$, as a function of $k_{c}$, calculated
by TC-DCD method. The 14-electron and 54-electron systems use a basis set
including 57 and 257 plane waves, respectively. All possible contractions from
the 3-body interactions are included, excluding the normal-ordered 3-body
interactions. The solid horizontal color lines represent the
$\|T_{2}^{\downarrow\uparrow}\|$ in the canonical DCD calculations with the
same settings as their TC counterparts.

By including the most important contractions, the error of neglecting the rest
of the 3-body interactions are well under control in that it scales
approximately as $\tilde{u}^{2}(k)k^{2}\sim\frac{1}{k^{6}}$. We note in
passing that if we choose correlators that do not truncate at small $k$, such
as the Yukawa-Coulomb correlator in Ref. Grüneis  _et al._ (2013) or the
Gaskell correlator in Ref. Gaskell (1961), the iterative solution of the
amplitude equations becomes too unstable to converge at low densities. We
attribute this instability in these cases to the large missing normal-ordered
3-body interactions, similar to the instability in a HF self-consistent
solution when the missing normal-ordered Coulomb interactions are large.

The FCIQMC calculations are carried out using the NECI program Guther  _et
al._ (2020). The CCD and DCD along with the TC integrals are implemented in a
Python program using the automatic tensor contraction engine CTF Solomonik
_et al._ (2014) and the NumPy package Harris  _et al._ (2020).

##  IV Results on the 3D UEG

$r_{s}$ |  $k_{c}$ ($\frac{2\pi}{L}$) | CCD | DCD | TC-CCD | TC-DCD | TC-FCIQMC | BF-DMC
---|---|---|---|---|---|---|---
0.5 | 1 | 3.41278 | 3.41252 | 3.41258 | 3.41244 | 3.41241(1) | 3.41370(2)
1.0 | 1 | 0.56975 | 0.56909 | 0.56891 | 0.56859 | 0.56861(1) | 0.56958(1)
2.0 | $\sqrt{2}$ | -0.00623 | -0.00748 | -0.00707 | -0.00800 | -0.00868(2) | -0.007949(7)
5.0 | $\sqrt{2}$ | -0.07618 | -0.07788 | -0.07816 | -0.07929 | -0.08002(2) | -0.079706(3)
10.0 | $\sqrt{2}$ | -0.05137 | -0.05289 | -0.05420 | -0.05509 | N/A | -0.055160(2)
20.0 | 2 | -0.02924 | -0.03035 | -0.03136 | -0.03201 | N/A | -0.0324370(8)
50.0 | $\sqrt{6}$ | -0.01261 | -0.01323 | -0.01350 | -0.01384 | N/A | -0.0146251(3)
Table 1: Total energy (a.u./electron), including the Madelung constant, of the 14-electron 3D UEG using different methods. $r_{s}$ |  $k_{c}$ ($\frac{2\pi}{L}$) | CCD | DCD | TC-CCD | TC-DCD | TC-FCIQMC | BF-DMC
---|---|---|---|---|---|---|---
0.5 | $\sqrt{2}$ | 3.22079 | 3.22052 | 3.22077 | 3.22071 | 3.22042(2) | 3.22112(4)
1.0 | $\sqrt{2}$ | 0.53069 | 0.53001 | 0.52982 | 0.52968 | 0.52973(3) | 0.52989(4)
2.0 | 2 | -0.01162 | -0.01286 | -0.01324 | -0.01379 | N/A | -0.01311(2)
5.0 | $\sqrt{5}$ | -0.07492 | -0.07655 | -0.07750 | -0.07837 | N/A | -0.079036(3)
10.0 | $2\sqrt{2}$ | -0.05016 | -0.05157 | -0.05230 | -0.05322 | N/A | -0.054443(2)
20.0 | 3 | -0.02846 | -0.02925 | -0.03055 | -0.03113 | N/A | -0.032047(2)
50.0 | 4 | -0.01223 | -0.01267 | -0.01263 | -0.01281 | N/A | N/A
Table 2: Total energy (a.u./electron), including the Madelung constant, of the
54-electron 3D UEG using different methods.

(a) (b)

Figure 4: Total energy per electron relative to the extrapolated value
retrieved as a function of $1/M$, where $M$ is the number of plane wave basis
functions for 3D UEG with 14 electrons. Linear extrapolations using the left
most two points are used in all cases.

(a) (b)

Figure 5: Errors in energy per electron relative to benchmark data on the 3D
UEG using TC-CCD, TC-DCD, CCD and DCD methods. (a) Results for 14 electrons.
For $r_{s}=0.5-5$, TC-FCIQMC data Luo and Alavi (2018) are used and for the
rest BF-DMC data Shepherd  _et al._ (2014) are used for benchmark. (b) Results
for 54 electrons. For $r_{s}=0.5-1$, TC-FCIQMC data Luo and Alavi (2018) are
used and for the rest BF-DMC data López Ríos  _et al._ (2006) are used for
benchmark. The grey shaded areas stand for the $\pm$0.001 a.u./electron
accuracy relative to the reference data. The colorful shaded areas reflect the
uncertainties in the TC-CCD and TC-DCD energies due to slightly different
choices of the $k_{c}$ values.

We first examine the basis set convergence behavior of TC-CCD/DCD compared to
the canonical ones. In Fig. 4 we present the total energy per electron
relative to the extrapolated value for each method, retrieved as a function of
the inverse of the employed number of plane waves, $1/M$. As mentioned before,
our correlator satisfies the cusp condition at the coalescence point of two
electrons with opposite spins. So the accelerated convergence behavior in the
TC methods compared to the canonical ones is not surprising. The acceleration
is the most obvious at high density regimes, since at low densities the
required number of basis functions to reach convergence in both the TC and
non-TC methods is relatively small. These observations are consistent with
those of the TC-FCIQMC reported in Ref. Luo and Alavi (2018).

The optimal $k_{c}$ values, (TC-)CCD/DCD energies at CBS and the benchmark
data are listed in Table 1 and 2 for the 14- and 54-electron 3D UEG,
respectively. In Fig. 5 we present the errors of total energies per electron
calculated by TC-CCD, TC-DCD, CCD and DCD relative to the most accurate FCIQMC
Luo and Alavi (2018); Shepherd  _et al._ (2012b); Neufeld and Thom (2017) and
backflow DMC (BF-DMC) López Ríos  _et al._ (2006); Shepherd  _et al._ (2014)
results on the 14- and 54-electron 3D UEG. The finite basis set errors in our
methods have been carefully eliminated by extrapolation to the CBS limit.

In general the accuracies of the TC methods are greatly improved compared to
their canonical counterparts, especially in regions ($r_{s}=5-50$) where the
latter exhibit the largest errors. More importantly, the improved accuracies
are retained when going from the 14- to the 54-electron system. We highlight
that the TC-DCD achieves an accuracy of $\leq 0.001$ a.u./electron across a
wide range of densities, i.e. $r_{s}=0.5-50$ for the 14-electron and
$r_{s}=0.5-20$ for the 54-electron 3D UEG, with an exception at $r_{s}=10$ for
the latter where it drops slightly out of the 0.001 a.u./electron accuracy. We
argue that with the next possible smaller value of $k_{c}=\sqrt{6}$, which
yields a marginally lower $\|T_{2}^{\downarrow\uparrow}\|$, instead of the
current choice of $k_{c}=2\sqrt{2}$, the 0.001 a.u./electron accuracy at
$r_{s}=10$ can be regained. The discrete grid of the $k$-mesh makes it hard to
pick the optimal $k_{c}$ in Fig. 3(b). However, as the system gets larger and
the $k$-mesh gets finer, the $\|T_{2}^{\downarrow\uparrow}\|$ as a function of
$k_{c}$ will also be smoother, and the choice of the optimal $k_{c}$ will
become more definite. We use colorful shaded areas in Fig. 5 to reflect the
uncertainties due to the possible choices of $k_{c}$ which yield similar
$\|T_{2}^{\downarrow\uparrow}\|$ values in Fig. 3.

At high densities, i.e. $r_{s}=0.5-2$, the canonical DCD is already very
accurate, and the main benefit from TC there is in accelerating the basis set
convergence. Overall, DCD exhibits smaller errors than CCD, which agrees with
earlier comparative studies between DCSD and CCSD Kats and Manby (2013a); Kats
(2014); Kats  _et al._ (2015).

##  V Conclusions

We demonstrated that the correlator Eq. (23), used with transcorrelated
coupled cluster theory, drastically improves the accuracy of approximate
methods, i.e. CCD and DCD, for 3D UEG across a wide range of densities. The
basis set convergence rate is also improved thanks to the fact that the
correlator satisfies the cusp condition at the coalescence point of two
unlike-spin electrons. We have explored the mechanism behind the improved
accuracy of the TC-CCD and TC-DCD methods, which is related to a
compactification of the many-electron wavefunction in Slater determinant space
when the dominant pair correlations between electrons with unlike spins are
directly included in the correlator. The optimization of the parameter in the
correlator is seamlessly incorporated within the TC-CCD/DCD framework, without
requiring an external algorithm. We notice that a range-separation scheme of
CCD can also achieve similar accuracy in 3D UEG, but without improving the
basis set convergence rate Shepherd  _et al._ (2014). Comparatively speaking,
our methods are systematically improvable, in that a more flexible form of the
correlator can be designed by a combination of a series of functions Cohen
_et al._ (2019) or by a general function approximator, e.g. an artificial
neural network, to include further correlation effects such as nuclear-
electron correlations and correlations between parallel pairs of electrons.
Other systematic ways of optimizing the correlator in combination with VMC
Cohen  _et al._ (2019) can also be explored. When going to real periodic
solids, TC-CCSD and TC-DCSD will be needed; extra efforts are also required to
compute the additional integrals besides the Coulomb integrals, where the most
computationally demanding part is the singly-contracted 3-body integrals which
scales like $\mathcal{O}(N_{o}^{2}N_{v}^{4})$, where $N_{o}$ and $N_{v}$ are
the number of occupied and unoccupied orbitals, respectively. Fortunately, the
computation of the extra integrals scales no worse than the CCSD or DCSD
algorithm and it can be compensated by the accelerated convergence rate
towards CBS limit in the TC framework. These perspectives will be important in
extending the encouraging performance of the current TC-CCD and TC-DCD methods
from the UEG to real periodic solids with moderate to strong correlation.

Competing Interests The authors declare that they have no competing financial
interests.

Acknowledgments We thank A. Grüneis and N. Masios for helpful discussions and
P. López Ríos for providing the BF-DMC data for the 14 electron UEG.

## References

  * Čížek (1966) J. Čížek, J. Chem. Phys. 45, 4256 (1966).
  * Čížek and Paldus (1971) J. Čížek and J. Paldus, International Journal of Quantum Chemistry 5, 359 (1971).
  * Čížek and Paldus (1980) J. Čížek and J. Paldus, Phys. Scr. 21, 251 (1980).
  * Raghavachari  _et al._ (1989) K. Raghavachari, G. W. Trucks, J. A. Pople, and M. Head-Gordon, Chem. Phys. Lett. 157, 479 (1989).
  * Riplinger  _et al._ (2016) C. Riplinger, P. Pinski, U. Becker, E. F. Valeev, and F. Neese, J. Chem. Phys. 144, 024109 (2016).
  * Schmitz and Hättig (2016) G. Schmitz and C. Hättig, J. Chem. Phys. 145, 234107 (2016).
  * Schwilk  _et al._ (2017) M. Schwilk, Q. Ma, C. Köppl, and H.-J. Werner, J. Chem. Theory Comput. 13, 3650 (2017).
  * Nagy  _et al._ (2018) P. R. Nagy, G. Samu, and M. Kállay, J. Chem. Theory Comput. 14, 4193 (2018).
  * Booth  _et al._ (2013) G. H. Booth, A. Grüneis, G. Kresse, and A. Alavi, Nature 493, 365 (2013).
  * Liao and Grüneis (2016) K. Liao and A. Grüneis, Journal of Chemical Physics 145, 0 (2016).
  * McClain  _et al._ (2017) J. McClain, Q. Sun, G. K.-L. Chan, and T. C. Berkelbach, J. Chem. Theory Comput. 13, 1209 (2017).
  * Gruber  _et al._ (2018) T. Gruber, K. Liao, T. Tsatsoulis, F. Hummel, and A. Grüneis, Phys. Rev. X 8, 021043 (2018).
  * Liao  _et al._ (2019) K. Liao, X.-Z. Li, A. Alavi, and A. Grüneis, npj Comput Mater 5, 110 (2019).
  * Gao  _et al._ (2020) Y. Gao, Q. Sun, J. M. Yu, M. Motta, J. McClain, A. F. White, A. J. Minnich, and G. K.-L. Chan, Phys. Rev. B 101, 165138 (2020).
  * Shepherd and Grüneis (2013) J. J. Shepherd and A. Grüneis, Physical Review Letters 110 (2013), 10.1103/PhysRevLett.110.226401, arXiv:1310.6059 .
  * Meyer (1971) W. Meyer, Int. J. Quantum Chem.: Quantum Chem. Symp. 5, 341 (1971).
  * Nooijen and Le Roy (2006) M. Nooijen and R. J. Le Roy, J. Mol. Struc. 768, 25 (2006).
  * Neese  _et al._ (2009) F. Neese, F. Wennmohs, and A. Hansen, J. Chem. Phys. 130, 114108 (2009).
  * Huntington and Nooijen (2010) L. M. J. Huntington and M. Nooijen, J. Chem. Phys. 133, 184109 (2010).
  * Paldus  _et al._ (1984) J. Paldus, J. Čížek, and M. Takahashi, Phys. Rev. A 30, 2193 (1984).
  * Piecuch and Paldus (1991) P. Piecuch and J. Paldus, Int. J. Quantum Chem. 40, 9 (1991).
  * Bartlett and Musiał (2006) R. J. Bartlett and M. Musiał, J. Chem. Phys. 125, 204105 (2006).
  * Robinson and Knowles (2011) J. B. Robinson and P. J. Knowles, J. Chem. Phys. 135, 044113 (2011).
  * Small  _et al._ (2014) D. W. Small, K. V. Lawler, and M. Head-Gordon, J. Chem. Theory Comput. 10, 2027 (2014).
  * Kats and Manby (2013a) D. Kats and F. R. Manby, J. Chem. Phys. 139, 021102 (2013a).
  * Kats (2014) D. Kats, The Journal of Chemical Physics 141, 061101 (2014).
  * Kats  _et al._ (2015) D. Kats, D. Kreplin, H.-J. Werner, and F. R. Manby, J. Chem. Phys. 142, 064111 (2015).
  * Rishi  _et al._ (2017) V. Rishi, A. Perera, M. Nooijen, and R. J. Bartlett, J. Chem. Phys. 146, 144104 (2017).
  * Tsatsoulis  _et al._ (2017) T. Tsatsoulis, F. Hummel, D. Usvyat, M. Schütz, G. H. Booth, S. S. Binnie, M. J. Gillan, D. Alfé, A. Michaelides, and A. Grüneis, J. Chem. Phys. 146, 204108 (2017).
  * Boys and Linnett (1969) S. F. Boys and J. W. Linnett, Proceedings of the Royal Society of London. A. Mathematical and Physical Sciences 309, 195 (1969).
  * Boys  _et al._ (1969) S. F. Boys, N. C. Handy, and J. W. Linnett, Proceedings of the Royal Society of London. A. Mathematical and Physical Sciences 310, 43 (1969).
  * Boys and Handy (1969) S. F. Boys and N. C. Handy, Proceedings of the Royal Society A: Mathematical, Physical and Engineering Sciences 309, 209 (1969).
  * Hino  _et al._ (2001) O. Hino, Y. Tanimura, and S. Ten-no, J. Chem. Phys. 115, 7865 (2001).
  * Hino  _et al._ (2002) O. Hino, Y. Tanimura, and S. Ten-no, Chemical Physics Letters 353, 317 (2002).
  * Umezawa and Tsuneyuki (2004) N. Umezawa and S. Tsuneyuki, Phys. Rev. B 69, 165102 (2004).
  * Sakuma and Tsuneyuki (2006) R. Sakuma and S. Tsuneyuki, J. Phys. Soc. Jpn. 75, 103705 (2006).
  * Ochi and Tsuneyuki (2015) M. Ochi and S. Tsuneyuki, Chemical Physics Letters 621, 177 (2015).
  * Ochi  _et al._ (2016) M. Ochi, Y. Yamamoto, R. Arita, and S. Tsuneyuki, J. Chem. Phys. 144, 104109 (2016).
  * Ochi  _et al._ (2017) M. Ochi, R. Arita, and S. Tsuneyuki, Phys. Rev. Lett. 118, 026402 (2017).
  * Luo and Alavi (2018) H. Luo and A. Alavi,  (2018), 10.1021/acs.jctc.7b01257.
  * Jeszenszki  _et al._ (2018) P. Jeszenszki, H. Luo, A. Alavi, and J. Brand, Phys. Rev. A 98, 053627 (2018).
  * Dobrautz  _et al._ (2019) W. Dobrautz, H. Luo, and A. Alavi, Physical Review B 99 (2019), 10.1103/PhysRevB.99.075119.
  * Cohen  _et al._ (2019) A. J. Cohen, H. Luo, K. Guther, W. Dobrautz, D. P. Tew, and A. Alavi, J. Chem. Phys. 151, 061101 (2019).
  * Jeszenszki  _et al._ (2020) P. Jeszenszki, U. Ebling, H. Luo, A. Alavi, and J. Brand, Phys. Rev. Research 2, 043270 (2020).
  * Wigner (1934) E. Wigner, Phys. Rev. 46, 1002 (1934).
  * Ceperley and Alder (1980) D. M. Ceperley and B. J. Alder, Phys. Rev. Lett. 45, 566 (1980).
  * Kohn and Sham (1965) W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965).
  * Perdew and Zunger (1981) J. P. Perdew and A. Zunger, Phys. Rev. B 23, 5048 (1981).
  * Perdew  _et al._ (1996) J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996).
  * Bohm and Pines (1953) D. Bohm and D. Pines, Phys. Rev. 92, 609 (1953).
  * Gell-Mann and Brueckner (1957) M. Gell-Mann and K. A. Brueckner, Physical Review 106, 364 (1957).
  * Shepherd  _et al._ (2012a) J. J. Shepherd, G. Booth, A. Grüneis, and A. Alavi, Physical Review B - Condensed Matter and Materials Physics 85, 1 (2012a), arXiv:1109.2635 .
  * Shepherd  _et al._ (2012b) J. J. Shepherd, G. H. Booth, and A. Alavi, J. Chem. Phys. 136, 244101 (2012b).
  * Shepherd  _et al._ (2014) J. J. Shepherd, T. M. Henderson, and G. E. Scuseria, Phys. Rev. Lett. 112, 133002 (2014).
  * Neufeld and Thom (2017) V. A. Neufeld and A. J. W. Thom, J. Chem. Phys. 147, 194105 (2017).
  * Ruggeri  _et al._ (2018) M. Ruggeri, P. L. Ríos, and A. Alavi, Phys. Rev. B 98, 161105 (2018).
  * Irmler  _et al._ (2019) A. Irmler, F. Hummel, and A. Grüneis, arXiv:1903.05458 [cond-mat, physics:physics] (2019), arXiv:1903.05458 [cond-mat, physics:physics] .
  * Kats and Manby (2013b) D. Kats and F. R. Manby, J. Chem. Phys. 138, 144101 (2013b).
  * Kats (2016) D. Kats, J. Chem. Phys. 144, 044102 (2016).
  * Pulay (1980) P. Pulay, Chemical Physics Letters 73, 393 (1980).
  * Pulay (1982) P. Pulay, Journal of Computational Chemistry 3, 556 (1982).
  * Booth  _et al._ (2009) G. H. Booth, A. J. Thom, and A. Alavi, Journal of Chemical Physics 131, 054106 (2009).
  * Gori-Giorgi  _et al._ (2000) P. Gori-Giorgi, F. Sacchetti, and G. B. Bachelet, Phys. Rev. B 61, 7353 (2000), arXiv:cond-mat/9909448 .
  * Note (1) $n=n_{1}^{2}+n_{2}^{3}+n_{3}^{2}$, where $n_{1},n_{2},n_{3}$ are the components of a ${\bf k}$ vector in units of $\frac{2\pi}{L}$. So the possible values of n are 0, 1, 2, 3, 4, 5, 6, 8, …
  * Note (2) Since the $\|T_{2}^{\downarrow\uparrow}\|$ in TC-CCD and TC-DCD show the same trend as a function of $k_{c}$, we show only the latter in the figure.
  * Grüneis  _et al._ (2013) A. Grüneis, J. J. Shepherd, A. Alavi, D. P. Tew, and G. H. Booth, The Journal of Chemical Physics 139, 084112 (2013).
  * Gaskell (1961) T. Gaskell, Proc. Phys. Soc. 77, 1182 (1961).
  * Guther  _et al._ (2020) K. Guther, R. J. Anderson, N. S. Blunt, N. A. Bogdanov, D. Cleland, N. Dattani, W. Dobrautz, K. Ghanem, P. Jeszenszki, N. Liebermann, G. L. Manni, A. Y. Lozovoi, H. Luo, D. Ma, F. Merz, C. Overy, M. Rampp, P. K. Samanta, L. R. Schwarz, J. J. Shepherd, S. D. Smart, E. Vitale, O. Weser, G. H. Booth, and A. Alavi, J. Chem. Phys. 153, 034107 (2020).
  * Solomonik  _et al._ (2014) E. Solomonik, D. Matthews, J. R. Hammond, J. F. Stanton, and J. Demmel, Journal of Parallel and Distributed Computing Domain-Specific Languages and High-Level Frameworks for High-Performance Computing, 74, 3176 (2014).
  * Harris  _et al._ (2020) C. R. Harris, K. J. Millman, S. J. van der Walt, R. Gommers, P. Virtanen, D. Cournapeau, E. Wieser, J. Taylor, S. Berg, N. J. Smith, R. Kern, M. Picus, S. Hoyer, M. H. van Kerkwijk, M. Brett, A. Haldane, J. F. del Río, M. Wiebe, P. Peterson, P. Gérard-Marchant, K. Sheppard, T. Reddy, W. Weckesser, H. Abbasi, C. Gohlke, and T. E. Oliphant, Nature 585, 357 (2020).
  * López Ríos  _et al._ (2006) P. López Ríos, A. Ma, N. D. Drummond, M. D. Towler, and R. J. Needs, Physical Review E 74 (2006), 10.1103/PhysRevE.74.066701.