MP2-F12: Explicitly Correlated MP2

Code author: Erica Mitchell

Section author: Erica Mitchell

Module: Keywords, PSI Variables, MP2-F12

Introduction

Convergence of the one-particle basis set to the electronic cusp is slow with traditional methods. This slow convergence is inherent in all traditional electronic structure methods and has led to the construction of what are called explicitly correlated methods. Explicitly correlated methods include terms that have the interelectronic distance within them significantly increasing the speed at which convergence to the electronic cusp is achieved. These explicitly correlated terms contain the correlation factor, $f_{12}$, which describes the behavior of the electronic cusp. The correlation factor is the defining part of explicitly correlated methods.

The inclusion of $f_{12}$ into the equations makes the prefactor for the computation much larger than traditional methods. When running MP2-F12, the 3C(FIX) Ansatz is utilized along with the complementary auxiliary basis set (CABS) approach resulting in an $\mathcal{O}(N^5)$ scaling separate from the prefactor. To lower the prefactor, density-fitting may be utilized but scaling will remain the same.

An example input file is:

molecule h2o {
0 1
O
H 1 1.0
H 1 1.0 2 104.5
symmetry c1
}

set {
  basis         cc-pvtz-f12
  cabs_basis    cc-pvtz-f12-optri
  df_basis_f12  aug-cc-pvtz-ri
  freeze_core   True
  mp2_type      df
  f12_beta      1.0
  cabs_singles  True
}

energy('mp2-f12')

The energy('mp2-f12') call to energy() executes the predefined MP-F12 procedure, first calling the SCF module with a default RHF reference. The two-electron integrals may be solved either with or without DF, as specified with MP2_TYPE. The general outline of an MP2-F12 computation is as follows:

1. Form the CABS AO to MO transformation matrix

2. Compute the one- and two-electron integrals

3. Construct the F12 intermediates, $V$, $X$, $C$, and $B$

4. Solve the F12 energy equation for each occupied orbital pair

The results will then look something like:

===> DF-MP2-F12/3C(FIX) Energies <===

 Total DF-MP2-F12/3C(FIX) Energy:      -76.362440627264
    RHF Reference Energy:              -76.059038661557
    MP2 Correlation Energy:             -0.276433879145
    F12/3C(FIX) Correlation Energy:     -0.026029803210
    CABS Singles Correction:            -0.000938283352

The theory and common keywords used in MP2-F12 are presented below.

Theory

In MP2-F12, the second-order MBPT energy may be determined by minimizing the Hylleraas functional,

(1)$$E^{(2)}=⟨{\psi }^{(1)}|{\hat{H}}^{(0)}-E^{(0)}|{\psi }^{(1)}⟩+2⟨{\psi }^{(1)}|{\hat{H}}^{(1)}|{\psi }^{(0)}⟩==\underset{|{\psi }^{(1)}⟩}{min}$$

The resulting linear first-order equations can then be broken down into the conventional residual and the explicitly correlated residual,

(2)$$\begin{array}{r}{R}_{ab}^{ij}=G_{ab}^{ij}+({ϵ}_a+{ϵ}_b-{ϵ}_i-{ϵ}_j)T_{ab}^{ij}+T_{mn}^{ij}{C}_{ab}^{mn}\\ R_{kl}^{ij}=V_{kl}^{ij}+C_{ab}^{kl}{T}_{ab}^{ij}+[B_{kl,mn}-({ϵ}_i+{ϵ}_j)X_{kl,mn}]T_{mn}^{ij}\end{array}$$

where the last term of the $R_{ab}^{ij}$ couples the two residuals. Within the residuals, the correlation factor, $f_{12}$, and the projector, ${\hat{Q}}_{12}$, are present in what are known as the F12 intermediates,

(3)$$\begin{array}{rl}{V}_{kl}^{ij}& =⟨ij|r_{12}^{-1}{\hat{Q}}_{12}{f}_{12}|kl⟩\\ X_{kl,mn}& =⟨kl|f_{12}{\hat{Q}}_{12}{f}_{12}|mn⟩\\ C_{ab}^{kl}& =⟨kl|f_{12}{\hat{Q}}_{12}({\hat{F}}_1+{\hat{F}}_2)|ab⟩\\ B_{kl,mn}& =⟨kl|f_{12}{\hat{Q}}_{12}({\hat{F}}_1+{\hat{F}}_2){\hat{Q}}_{12}{f}_{12}|mn⟩\end{array}$$

The correlation factor and projector are what hold the majority of the approximations made in F12 theory. Firstly, the correlation factor is what describes the electronic cusp but, similar to how a nuclear cusp is simulated by contracted Gaussian-type orbitals, the $f_{12}$ is composed of multiple Gaussian-type geminals to imitate a Slater-type geminal,

(4)\[f_{12} = \frac{1 - e^{-\beta r_{12}}{\beta} \approx \sum^6_{n=1} c_n e^{-\alpha_n r_{12}^2}\]

where the coefficients and exponents are taken from Tew and Klopper [Tew:2005:074101] and can be scaled using the parameter $\beta$. Secondly, the projector is what holds the chosen Ansatz. In other words, the projector determines the choice of orthogonalization between the explicitly correlated functions and to the conventional doubly substituted determinants. The strongest orthogonality is observed with

(5)$${\hat{Q}}_{12}=(1-{\hat{O}}_1)(1-{\hat{O}}_2)-{\hat{V}}_1{\hat{V}}_2$$

where $\hat{O}$ represents $|i⟩⟨i|$ and $\hat{V}$ represents $|a⟩⟨a|$. This form of the projector is known as Ansatz 3.

Within modern F12 theory, the ansatz inherently holds the resolution-of-identity (RI). The RI can be taken as functions from the orbital basis set (OBS) and a chosen auxiliary basis set (ABS). It was observed by Valeev [Valeev:2004:190] that the choice of ABS can affect the error due to the use of a finite RI. He proposed that the ABS include the OBS explicitly so that the exact orthogonal complement projector ${\hat{Q}}_{12}=1-{\hat{P}}_{12}$ is properly approximated. This lead to what is known as the complementary auxiliary basis set or CABS. The expanded Ansatz 3 projector with the RI and CABS is then

(6)$${\hat{Q}}_{12}=1-|a^{\prime }j⟩⟨a^{\prime }j|-|i{b}^{\prime }⟩⟨i{b}^{\prime }|-|rs⟩⟨rs|$$

When applying the projector within the CABS approach to the F12 intermediates, the $C$ and $B$ intermediates contain terms which require commutation of the correlation factor, $f_{12}$, and the Fock operator, $\hat{F}$. To factorize these, there are three different approaches, A, B, and C. In general, Ansatz 3 with approximation C, or Ansatz 3C, is computationally cheaper than ans{"a}tze 3A and 3B. In Ansatz 3C, the unit operator in the projector is treated with the commutator approach while the rest utilizes direct RI expansions. Within Ansatz 3C, the F12 intermediates then take the form

(7)$$\begin{array}{rl}{V}_{kl}^{ij}& =(\overline{gf}{)}_{kl}^{ij}-{\bar{g}}_{a^{\prime }n}^{ij}{\bar{f}}_{kl}^{a^{\prime }n}-{\bar{g}}_{m{b}^{\prime }}^{ij}{\bar{f}}_{kl}^{m{b}^{\prime }}-{\bar{g}}_{rs}^{ij}{\bar{f}}_{kl}^{rs}\\ X_{kl,mn}& ={\left(f^2\right)}_{mn}^{kl}-{\bar{f}}_{i{b}^{\prime }}^{kl}{\bar{f}}_{mn}^{i{b}^{\prime }}-{\bar{f}}_{a^{\prime }j}^{kl}{\bar{f}}_{mn}^{a^{\prime }j}-{\bar{f}}_{pq}^{kl}{\bar{f}}_{mn}^{pq}\\ C_{ab}^{kl}& ={\bar{f}}_{a{b}^{\prime }}^{kl}{F}_b^{b^{\prime }}+{\bar{f}}_{a^{\prime }b}^{kl}{F}_a^{a^{\prime }}\\ B_{kl,mn}& =⟨kl|[{\hat{f}}_{12},[{\hat{T}}_{12},{\hat{f}}_{12}]]|mn⟩+{\hat{S}}_{mn,kl}[(\overline{f^2}{)}_{mn}^{p^{\prime }l}(F+K{)}_{p^{\prime }}^k+(\overline{f^2}{)}_{mn}^{k{q}^{\prime }}(F+K{)}_{q^{\prime }}^l]\\ & -{\bar{f}}_{mn}^{p^{\prime }{q}^{\prime }}{K}_{p^{\prime }}^{r^{\prime }}{\bar{f}}_{r^{\prime }{q}^{\prime }}^{kl}-{\bar{f}}_{mn}^{p^{\prime }j}{F}_{p^{\prime }}^{r^{\prime }}{\bar{f}}_{r^{\prime }j}^{kl}+{\bar{f}}_{mn}^{i{b}^{\prime }}{F}_i^j{\bar{f}}_{j{b}^{\prime }}^{kl}-{\bar{f}}_{mn}^{pb}{F}_p^r{\bar{f}}_{rb}^{kl}\\ & -2{\hat{S}}_{mn,kl}[{\bar{f}}_{mn}^{i{b}^{\prime }}{F}_i^{p^{\prime }}{\bar{f}}_{p^{\prime }{b}^{\prime }}^{kl}+{\bar{f}}_{mn}^{a^{\prime }b}{F}_{a^{\prime }}^q{\bar{f}}_{qb}^{kl}]\end{array}$$

where the following integral types are revealed,

(8)$$\begin{array}{rl}& ⟨p^{\prime }{q}^{\prime }|r_{12}^{-1}|r^{\prime }{s}^{\prime }⟩=g_{r^{\prime }{s}^{\prime }}^{p^{\prime }{q}^{\prime }}\\ & ⟨p^{\prime }{q}^{\prime }|{\hat{f}}_{12}|r^{\prime }{s}^{\prime }⟩=f_{r^{\prime }{s}^{\prime }}^{p^{\prime }{q}^{\prime }}\\ & ⟨p^{\prime }{q}^{\prime }|{\hat{f}}_{12}^2|r^{\prime }{s}^{\prime }⟩=(f^2{)}_{r^{\prime }{s}^{\prime }}^{p^{\prime }{q}^{\prime }}\\ & ⟨p^{\prime }{q}^{\prime }|{\hat{f}}_{12}{r}_{12}^{-1}|r^{\prime }{s}^{\prime }⟩=(gf{)}_{r^{\prime }{s}^{\prime }}^{p^{\prime }{q}^{\prime }}\\ & ⟨p^{\prime }{q}^{\prime }|[{\hat{f}}_{12},[{\hat{T}}_{12},{\hat{f}}_{12}]]|r^{\prime }{s}^{\prime }⟩\end{array}$$

Before the energy expression can be revealed, the choice of explicitly correlated amplitudes must be taken into account. Although the residuals can both be set to zero to solve for the minimized amplitudes, it has been determined that the best approach is by utilizing the electronic cusp conditions so that the s- and p- wave coalesence conditions are met [Tenno:2004:117]. This results in only the conventional residual, $R_{ab}^{ij}$, being minimized and the explicitly correlated amplitudes being

(9)$$T_{kl}^{ij}=\frac{3}{8}{\delta }_{ik}\delta jl+\frac{1}{8}{\delta }_{jk}\delta il$$

which is denoted as the 3C(FIX) Ansatz. The 3C(FIX) Ansatz has the benefits that it is free from the geminal basis set superposition error, orbital invariant, and diagonal. Since only diagonal terms are needed, the computational complexity reduces from $\mathcal{O}(N^6)$ with the minimized amplitudes to $\mathcal{O}(N^5)$ with fixed amplitudes.

Finally, the F12/3C(FIX) energy correction can be expressed as

(10)$$\mathrm{\Delta }{E}_{\mathrm{F}12/3\mathrm{C}(\mathrm{FIX})}=T_{kl}^{ij}(2\text{Tilde}{V}_{kl}^{ij}+\text{Tilde}{B}_{kl,mn}^{ij}{T}_{mn}^{ij})$$

where $\text{Tilde}{V}_{kl}^{ij}$ and $\text{Tilde}{B}_{kl,mn}^{ij}$ are

(11)$$\begin{array}{rl}\text{Tilde}{V}_{kl}^{ij}& =V_{kl}^{ij}-\frac{C_{ab}^{kl}{G}_{ab}^{ij}}{{ϵ}_a+{ϵ}_b-{ϵ}_i-{ϵ}_j}\\ \text{Tilde}{B}_{kl,mn}^{ij}& =B_{kl,mn}-({ϵ}_i+{ϵ}_j)X_{kl,mn}-\frac{C_{ab}^{mn}{C}_{ab}^{kl}}{{ϵ}_a+{ϵ}_b-{ϵ}_i-{ϵ}_j}\end{array}$$

With the slow convergence of the conventional doubles energy handled, the primary source of error then comes from the basis set incompleteness error from the HF reference. When the CABS is utilized, the relaxation of the HF orbitals in the presence of the CABS needs to be accounted for. This is done using the what’s termed the CABS singles correction. The CABS singles correction is reminiscent of the singles correction in MP2 without a RHF or UHF reference,

(12)$$\mathrm{\Delta }{E}_{\mathrm{CABS}\mathrm{Singles}}=\frac{|F_i^{\alpha }{|}^2}{{ϵ}_{\alpha }-{ϵ}_i}$$

where $\alpha$ runs over all virtual orbitals, OBS virtuals and CABS.

Therefore, the most accurate MP2-F12/3C(FIX) energy is given by

(13)$$E_{\mathrm{MP}2-\mathrm{F}12/3\mathrm{C}(\mathrm{FIX})}=E_{\mathrm{MP}2}+\mathrm{\Delta }{E}_{\mathrm{F}12}+\mathrm{\Delta }{E}_{\mathrm{CABS}\mathrm{Singles}}$$

Recommendations

MP2-F12 is a complicated theory but the implemented module should be simple to use. All keywords are appropriately documented in the Appendix F12 and care should be taken to read through them before changing defaults. Some basic recommendations for running the module are included below:

• MP2-F12 runs much faster with density-fitted integrals. Unless high accuracy is needed, DF requires much less computational resources and still provides quality results. If DF is chosen, the HF and MP2 modules are also set to run DF. Choose a DF basis set with at least triple-zeta quality as PSI4‘s DF-SCF module performs poorly with lower quality DF basis sets. This module utilizes robust density-fitting [Dunlap:2000:2113] so energies may differ from other programs. The default MP2_TYPE is DF.

• The length-scale parameter F12_BETA has been fiddled with in many past studies and they have determined that the exponents 0.9, 1.0, and 1.0 $a_0^{-1}$ are recommended for cc-pVXZ-F12 where X = D, T, and Q, respectively, while those for aug-cc-pVXZ where X = D, T, and Q are 1.1, 1.2, and 1.4 $a_0^{-1}$. The default value is set to 1.0 $a_0^{-1}$ to allow flexibility of use.

• MP2-F12 likes threads. The module was developed with threading in mind and the formation of the integrals and the tensor contractions rely on the threads for efficiency.

• MP2-F12 is memory greedy. The core algorithm has been developed to compute the integrals on-the-fly but the size of some integrals cannot be avoided. Care should be taken to choose whether the core or disk algorithm is chosen, F12_SUBTYPE .

• Like most wavefunction methods, freezing the core is good for both efficiency and correctness purposes.

• MP2-F12 is not symmetry aware.

• MP2-F12 cannot use Psi4’s Schwarz screening due to the presence of mixed basis sets. This is a work in progress but siginificantly increases the computational prefactor.

• At the moment, the MP2-F12 code is only compatible with RHF references.