92 lines
10 KiB
TeX
92 lines
10 KiB
TeX
%****************************************************************
|
|
\section{History}
|
|
%****************************************************************
|
|
It is well known that highly-accurate wave functions require the fulfilment (or near-fulfilment) of the electron-electron cusp conditions \cite{Kato51, Kato57, Pack66, Morgan93, Myers91, Tew08, QuasiExact09, ExSpherium10, QR12, Kurokawa13, Kurokawa14, Gruneis17}.
|
|
For correlated wave functions expanded in terms of products of one-electron Gaussian basis functions, the energy converges as $\order{L^{-3}}$, where $L$ is the maximum angular momentum of the basis set \cite{Kutzelnigg85}.
|
|
This slow convergence can be tracked down to the inadequacy of these products to properly model the Coulomb correlation hole \cite{Hattig12, Kong12}.
|
|
|
|
In the late 20's, Hylleraas solved this issue for the helium atom by introducing explicitly the interelectronic distance $\ree = \abs{\br_1 - \br_2}$ as an additional two-electron basis function \cite{Hylleraas28, Hylleraas29}.
|
|
As Kutzelnigg later showed, this leads to a prominent improvement of the energy convergence from $\order{L^{-3}}$ to $\order{L^{-7}}$ \cite{Kutzelnigg85}.
|
|
|
|
Around the same time, Slater, while studying the Rydberg series of helium \cite{Slater28}, suggested a new correlation factor known nowadays as a Slater geminal:
|
|
\begin{equation}
|
|
S_{12} = \exp(-\lambda_{12}\, \ree ).
|
|
\end{equation}
|
|
Unfortunately, the increase in mathematical complexity brought by $\ree$ or $S_{12}$ has been found to be rapidly computationally overwhelming.\footnote{Note that, although Slater was the first to propose such correlation factor, he suggested to set $\lambda_{12} = -1/2$ in order to ensure that the wave function fulfils Kato's electron-electron cusp condition. However, Hartree and Ingman \cite{Hartree33} found this correlation factor physically unreasonable due to its behaviour at large $\ree$, and suggested that a correlation factor of the form $1 - c \exp(-\lambda_{12}\, \ree)$ (with $\lambda_{12} >0$) would be more appropriate. We refer the interested reader to the review of Hattig et al.~\cite{Hattig12} for a detailed historical overview.}
|
|
|
|
In 1960, Boys \cite{Boys60} and Singer \cite{Singer60} independently proposed to use the Gaussian geminal (GG) correlation factor
|
|
\begin{equation}
|
|
\label{eq:gg}
|
|
G_{12} = \exp(-\lambda_{12}\, \ree^2 ),
|
|
\end{equation}
|
|
as \cite{Boys60}
|
|
|
|
\begin{quote}
|
|
``\textit{there are explicit formulas for all of the necessary many-dimensional integrals}".
|
|
\end{quote}
|
|
|
|
Interestingly, in the same article, a visionary Boys argued that, even if GGs do not fulfil the electron-electron cusp conditions, they could be used to fit $S_{12}$.
|
|
|
|
During the following years, variational calculations involving GGs flourished, giving birth to various methods, such as the exponentially-correlated Gaussian method \cite{RychBook, Bubin2005, Bubin2009, Szalewicz2010}.
|
|
However, this method was restricted to fairly small systems as it requires the optimisation of a large number of non-linear parameters.
|
|
In the 70's, the first MP2 calculations including GGs appeared thanks to the work by Pan and King \cite{Pan70, Pan72}, Adamowicz and Sadlej \cite{Adamowicz77a, Adamowicz77b, Adamowicz78}, and later Szalewicz \textit{et al.} \cite{Szalewicz82, Szalewicz83}.
|
|
Even if these methods represented a substantial step forward in terms of computational burden, they still require the optimization of a large number of non-linear parameters.
|
|
|
|
In 1985, Kutzelnigg derived a first form of the MP2-R12 equations using $\ree$ as a correlation factor \cite{Kutzelnigg85}.
|
|
Kutzelnigg's idea, which was more formally described together with Klopper in 1987 \cite{Klopper87}, dredged up an old problem: in addition to two-electron integrals (traditional ones and new ones), three-electron and four-electron integrals were required.
|
|
At that time, the only way to evaluate them would have been via an expensive one- or two-dimensional Gauss-Legendre quadrature \cite{Preiskorn85, Clementi89}.
|
|
Additionally, citing Kutzelnigg and Klopper \cite{Kutzelnigg91},
|
|
\begin{quote}
|
|
\textit{``even if fast procedures for the evaluation of these integrals were available, one would have to face the problem of the large number of these integrals; while that of two-electron integrals is $\sim N^{4}$, there are $\sim N^{6}$ three-electron and $\sim N^{8}$ four-electron integrals.
|
|
The storing and manipulating of these integrals could be handled only for extremely small basis sets.''}
|
|
\end{quote}
|
|
|
|
Undoubtedly, in the late 80's, the two-electron integrals technology was still in development \cite{MD78, PH78, King76, Dupuis76, Rys83, OS1, HGP, Gill94b}.
|
|
Nowadays, though still challenging, these integrals could be computed much more effectively via judicious recursive schemes, designing the quadrature only to the fundamental integrals \cite{3ERI1}.
|
|
Another important remark is that the actual number of \textit{significant} (i.e.~greater than a given threshold) three- and four-electron integrals in a large system, is, at worst, $\order{N^{3}}$ or $\order{N^{4}}$.
|
|
These kinds of scaling are achievable, for example, by exploiting robust density fitting \cite{Womack14} or upper bound-based screening methods, as discussed below \cite{3ERI2}.
|
|
|
|
Nevertheless, the success of the R12 method was triggered by the decision to avoid the computation of these three- and four-electron integrals through the use of the \textit{resolution of the identity} (RI) approximation \cite{Kutzelnigg91, Hattig12, Werner07}.
|
|
In this way, three- and four-electron integrals are approximated as linear combinations of products of two-electron integrals.
|
|
Several key developments and improvements of the original MP2-R12 approach have been proposed in the last decade \cite{MOLPRO2011, Hattig12, Werner07, Kong12}.
|
|
Of course, the accuracy of the RI approximation relies entirely on the assumption that the auxiliary basis set is sufficiently large, i.e.~$N_\text{RI} \gg N$, where $N$ and $N_\text{RI}$ are the number of basis functions in the primary and auxiliary basis sets, respectively.
|
|
|
|
The use of RI as method of choice does not seem definitive to us.
|
|
In fact, eschewing the RI approximation would offer at least two advantages:
|
|
i) smaller one-electron basis as the larger auxiliary basis set would not be required anymore;
|
|
ii) the three- and four-electron integrals would be computed exactly.
|
|
Moreover, one could avoid the commutator rearrangements involved in the computation of integrals over the kinetic energy operator \cite{Rohse93}.
|
|
|
|
In 1996, Persson and Taylor killed two birds with one stone.
|
|
Using a pre-optimized GG expansion fitting $\ree$
|
|
\begin{equation}
|
|
\label{eq:r12fit}
|
|
\ree \approx \sum_{k} a_{k} \qty[1-\exp(-\lambda_{k} \ree^{2}) ],
|
|
\end{equation}
|
|
they avoided the non-linear optimisation, and eschewed the RI approximation thanks to the analytical integrability of three- and four-electron integrals over GGs \cite{Persson96}.
|
|
They were able to show that a six- or ten-term fit introduces a $0.5$ mE$_\text{h}$ or $20$ $\mu$E$_\text{h}$ error, respectively \cite{Persson96}.
|
|
Unfortunately, further improvements were unsuccessful due to the failure of $\ree$ in modelling the correct behaviour of the wave function for intermediate and large $\ree$ \cite{Tew05, Tew06}.
|
|
In fact, Ten-no showed that $S_{12}$ is near-optimal at describing the correlation hole, and that a 10-term GG fit of $S_{12}$ yields very similar results.
|
|
This suggests that, albeit not catching the cusp per se, the Coulomb correlation hole can be accurately represented by GGs \cite{Tenno04a, Tenno07, May04, May05}.
|
|
|
|
Methods for evaluating many-electron integrals involving GGs have already been developed.
|
|
As mentioned previously, Persson and Taylor \cite{Persson97} derived recurrence relations based on Hermite Gaussians, analogously to the work of McMurchie and Davidson for two-electron integrals \cite{MD78}.
|
|
These recurrence relations were implemented by Dahle \cite{DahleThesis, Dahle2007, Dahle2008}.
|
|
Saito and Suzuki \cite{Saito01} also proposed an approach based on the work by Obara and Saika \cite{OS1, OS2}.
|
|
More recently, a general formulation using Rys polynomials \cite{King76, Dupuis76, Rys83} was published by Komornicki and King \cite{Komornicki11}.
|
|
Even if limited to the three-centre case, it is worth mentioning that May has also developed recurrence relations for two types of three-electron integrals \cite{MayThesis}.
|
|
These recurrence relations were implemented by Womack using automatically-generated code \cite{WomackThesis}.
|
|
Recently, we have developed recurrence relations for three- and four-electron integrals for generic correlation factors \cite{3ERI1, 4ERI1}.
|
|
All these recursive schemes have variable computational cost depending on the degree of contraction of the integral class to be computed.
|
|
Unsurprisingly, these algorithms are, to a large extent, complementary, and none of the algorithms has proven optimal under all circumstances.
|
|
|
|
A major limitation of all these approaches is that they do not include any integral screening.
|
|
Indeed, a remarkable consequence of the short-range nature of the Slater and Gaussian correlation factors is that, even if formally scaling as $\order{N^{6}}$ and $\order{N^{8}}$, there are only $\order{N^{2}}$ \textit{significant} (i.e.~greater than a given threshold) three- and four-electron integrals in a large system \cite{3ERI1, 3ERI2}.
|
|
Therefore, it is paramount to devise rigorous upper bounds to avoid computing the large number of negligible integrals.
|
|
|
|
This chapter is organised as follows.
|
|
First, we discuss Gaussian basis functions, many-electron integrals and the structure of the three- and four-electron operators considered here.
|
|
In the next three sections, we introduce the main ingredients for the efficient computation of three- and four-electron integrals involving GGs: i) fundamental integrals (FIs), ii) upper bounds (UBs), and iii) recurrence relations (RRs).
|
|
Finally, we give an overall view of our algorithm which is an extension of the late-contraction path of PRISM, also known as the Head-Gordon-Pople (HGP) algorithm (see Refs.~\cite{HGP, Gill94b} and references therein).
|
|
Note that the HGP algorithm corresponds to the Obara-Saika scheme where one has introduced additional RRs and provides a precise description of how these RRs can be used in tandem to good effects.
|