n-electron valence state perturbation theory

A possible approach is to define the perturber wavefunctions into Hilbert spaces ${\displaystyle S_{l}^{k}}$ defined by those determinants with given k and l labels. The determinants characterizing these spaces can be written as a partition comprising the same inactive (core + virtual) part ${\displaystyle \Phi _{l}^{-k}}$ and all possible valence (active) parts ${\displaystyle \Psi _{I}^{k}}$

${\displaystyle S_{l}^{k}\ {\stackrel {\mathrm {def} }{=}}\ \{\Phi _{l}^{-k}\Psi _{I}^{k}\}}$

The full dimensionality of these spaces can be exploited to obtain the definition of the perturbers, by diagonalizing the Hamiltonian inside them

${\displaystyle {\hat {\mathcal {P}}}_{S_{l}^{k}}{\hat {\mathcal {H}}}{\hat {\mathcal {P}}}_{S_{l}^{k}}\left|\Phi _{l}^{-k}\Psi _{\mu }^{v+k}\right\rangle =E_{l,\mu }\left|\Phi _{l}^{-k}\Psi _{\mu }^{v+k}\right\rangle }$

This procedure is impractical given its high computational cost: for each ${\displaystyle S_{l}^{k}}$ space, a diagonalization of the true Hamiltonian must be performed. Computationally, is preferable to improve the theoretical development making use of the modified Dyall's Hamiltonian ${\displaystyle {\hat {\mathcal {H}}}^{D}}$. This Hamiltonian behaves like the true Hamiltonian inside the CAS space, having the same eigenvalues and eigenvectors of the true Hamiltonian projected onto the CAS space. Also, given the decomposition for the wavefunction defined before, the action of the Dyall's Hamiltonian can be partitioned into

${\displaystyle {\hat {\mathcal {H}}}^{D}\left|\Phi _{l}^{-k}\Psi _{\mu }^{v+k}\right\rangle =E_{l,\mu }^{k}\left|\Phi _{l}^{-k}\Psi _{\mu }^{v+k}\right\rangle }$

stripping out the constant contribution of the inactive part and leaving a subsystem to be solved for the valence part

${\displaystyle {\hat {\mathcal {H}}}_{v}^{D}\left|\Psi _{\mu }^{v+k}\right\rangle =E_{\mu }^{k}\left|\Psi _{\mu }^{v+k}\right\rangle }$

The total energy ${\displaystyle E_{l,\mu }^{k}}$ is the sum of ${\displaystyle E_{\mu }^{k}}$ and the energies of the orbitals involved in the definition of the inactive part ${\displaystyle \Phi _{l}^{-k}}$. This introduces the possibility to perform a single diagonalization of the valence Dyall's Hamiltonian on the CASCI zero-order wavefunction and evaluate the perturber energies using the property depicted above.

A different choice in the development of the NEVPT approach is to choose a single function for each space ${\displaystyle S_{l}^{k}}$, leading to the Strongly Contracted (SC) scheme. A set of perturbative operators are used to produce a single function for each space, defined as the projection inside each space ${\displaystyle {\hat {\mathcal {P}}}_{S_{l}^{k}}}$ of the application of the Hamiltonian to the contracted zero order wavefunction. In other words,

${\displaystyle \Psi _{l}^{k}={\hat {\mathcal {P}}}_{S_{l}^{k}}{\hat {\mathcal {H}}}\Psi _{m}^{(0)}}$

where ${\displaystyle {\hat {\mathcal {P}}}_{S_{l}^{k}}}$ is the projector onto the subspace. This can be equivalently written as the application of a specific part of the Hamiltonian to the zero-order wavefunction

${\displaystyle \Psi _{l}^{k}=V_{l}^{k}\Psi _{m}^{(0)}}$

For each space, appropriate operators can be devised. We will not present their definition, as it could result overkilling. Suffice to say that the resulting perturbers are not normalized, and their norm

${\displaystyle N_{l}^{k}=\left\langle \Psi _{l}^{k}\left.\right|\Psi _{l}^{k}\right\rangle =\left\langle \Psi _{m}^{(0)}\left|\left(V_{l}^{k}\right)^{+}V_{l}^{k}\right|\Psi _{m}^{(0)}\right\rangle }$

plays an important role in the Strongly Contracted development. To evaluate these norms, the spinless density matrix of rank not higher than three between the ${\displaystyle \Psi _{m}^{(0)}}$ functions are needed.

An important property of the ${\displaystyle \Psi _{l}^{k}}$ is that any other function of the space ${\displaystyle S_{l}^{k}}$ which is orthogonal to ${\displaystyle \Psi _{l}^{k}}$ do not interact with the zero-order wavefunction through the true Hamiltonian. It is possible to use the ${\displaystyle \Psi _{l}^{k}}$ functions as a basis set for the expansion of the first-order correction to the wavefunction, and also for the expression of the zero-order Hamiltonian by means of a spectral decomposition

${\displaystyle {\hat {\mathcal {H}}}_{0}=\sum _{lk}\left|\Psi _{l}^{k}{}^{\prime }\right\rangle E_{l}^{k}\left\langle \Psi _{l}^{k}{}^{\prime }\right\rangle +\sum _{m}\left|\Psi _{m}^{(0)}\right\rangle E_{m}^{(0)}\left\langle \Psi _{m}^{(0)}\right|}$

where ${\displaystyle \left|\Psi _{l}^{k}{}^{\prime }\right\rangle }$ are the normalized ${\displaystyle \left|\Psi _{l}^{k}\right\rangle }$.

The expression for the first-order correction to the wavefunction is therefore

${\displaystyle \Psi _{m}^{(1)}=\sum _{kl}\left|\Psi _{l}^{k}{}^{\prime }\right\rangle {\frac {\left\langle \Psi _{l}^{k}{}^{\prime }\left|{\hat {\mathcal {H}}}\right|\Psi _{m}^{(0)}\right\rangle }{E_{m}^{(0)}-E_{l}^{k}}}=\sum _{kl}\left|\Psi _{l}^{k}{}^{\prime }\right\rangle {\frac {\sqrt {N_{l}^{k}}}{E_{m}^{(0)}-E_{l}^{k}}}}$

and for the energy is

${\displaystyle E_{m}^{(2)}=\sum _{kl}{\frac {\left|\left\langle \Psi _{l}^{k}{}^{\prime }\left|{\hat {\mathcal {H}}}\right|\Psi _{m}^{(0)}\right\rangle \right|^{2}}{E_{m}^{(0)}-E_{l}^{k}}}=\sum _{kl}{\frac {N_{l}^{k}}{E_{m}^{(0)}-E_{l}^{k}}}}$

This result still misses a definition of the perturber energies ${\displaystyle E_{l}^{k}}$, which can be defined in a computationally advantageous approach by means of the Dyall's Hamiltonian

${\displaystyle E_{l}^{k}={\frac {1}{N_{l}^{k}}}\left\langle \Psi _{l}^{k}\left|{\hat {\mathcal {H}}}^{D}\right|\Psi _{l}^{k}\right\rangle }$

leading to

${\displaystyle N_{l}^{k}E_{l}^{k}=\left\langle \Psi _{m}^{(0)}\left|\left(V_{l}^{k}\right)^{+}{\hat {\mathcal {H}}}^{D}V_{l}^{k}\right|\Psi _{m}^{(0)}\right\rangle =\left\langle \Psi _{m}^{(0)}\left|\left(V_{l}^{k}\right)^{+}V_{l}^{k}{\hat {\mathcal {H}}}^{D}\right|\Psi _{m}^{(0)}\right\rangle +\left\langle \Psi _{m}^{(0)}\left|\left(V_{l}^{k}\right)^{+}\left[{\hat {\mathcal {H}}}^{D},V_{l}^{k}\right]\right|\Psi _{m}^{(0)}\right\rangle }$

Developing the first term and extracting the inactive part of the Dyall's Hamiltonian it can be obtained

${\displaystyle E_{l}^{k}=E_{m}^{(0)}+\Delta \epsilon _{l}+{\frac {1}{N_{l}^{k}}}\left\langle \Psi _{m}^{(0)}\left|\left(V_{l}^{k}\right)^{+}\left[{\hat {\mathcal {H}}}_{v},V_{l}^{k}\right]\right|\Psi _{m}^{(0)}\right\rangle }$

with ${\displaystyle \Delta \epsilon _{l}}$ equal to the sum of the orbital energies of the newly occupied virtual orbitals minus the orbital energies of the unoccupied core orbitals.

The term that still needs to be evaluated is the bracket involving the commutator. This can be obtained developing each ${\displaystyle V}$ operator and substituting. To obtain the final result it is necessary to evaluate Koopmans matrices and density matrices involving only active indexes. An interesting case is represented by the contribution for the ${\displaystyle V_{ijrs}^{(0)}}$ case, which is trivial and can be demonstrated identical to the Møller–Plesset second-order contribution

${\displaystyle E_{m}^{(2)}\left(S_{rsij}^{0}\right)=-{\frac {N_{rsij}^{0}}{\epsilon _{r}+\epsilon _{s}-\epsilon _{i}-\epsilon _{j}}}}$

NEVPT2 can therefore be seen as a generalized form of MP2 to multireference wavefunctions.

An alternative approach, named Partially Contracted (PC) is to define the perturber wavefunctions in a subspace ${\displaystyle {\overline {S}}_{l}^{k}}$ of ${\displaystyle S_{l}^{k}}$ with dimensionality higher than one (like in case of the Strongly Contracted approach). To define this subspace, a set of functions ${\displaystyle \Phi }$ is generated by means of the ${\displaystyle V_{l}^{k}}$ operators, after decontraction of their formulation. For example, in the case of the ${\displaystyle V_{rsi}^{-1}}$ operator

${\displaystyle V_{rsi}^{-1}=\gamma _{rs}\sum _{a}\left(\left\langle rs\left.\right|ia\right\rangle E_{ri}E_{sa}+\left\langle sr\left.\right|ia\right\rangle E_{si}E_{ra}\right)\quad r\leq s}$

The Partially Contracted approach makes use of functions ${\displaystyle \Phi _{risa}=E_{ri}E_{sa}\Psi _{m}^{(0)}}$ and ${\displaystyle \Phi _{risa}=E_{si}E_{ra}\Psi _{m}^{(0)}}$. These functions must be orthonormalized and purged of linear dependencies which may arise. The resulting set spans the ${\displaystyle {\overline {S}}_{rsi}^{-1}}$ space.

Once all the ${\displaystyle {\overline {S}}_{l}^{k}}$ spaces have been defined, we can obtain as usual a set of perturbers from the diagonalization of the Hamiltonian (true or Dyall) inside this space

${\displaystyle {\hat {\mathcal {P}}}_{{\overline {S}}_{l}^{k}}{\hat {\mathcal {H}}}{\hat {\mathcal {P}}}_{{\overline {S}}_{l}^{k}}\left|\Psi _{l\mu }^{k}\right\rangle =E_{l,\mu }^{k}\left|\Psi _{l\mu }^{k}\right\rangle }$

As usual, the evaluation of the Partially Contracted perturbative correction by means of the Dyall Hamiltonian involves simply manageable entities for nowadays computers.

Although the Strongly Contracted approach makes use of a perturbative space with very low flexibility, in general it provides values in very good agreement with those obtained by the more decontracted space defined for the Partially Contracted approach. This can be probably explained by the fact that the Strongly Contracted perturbers are a good average of the totally decontracted perturbative space.

The Partially Contracted evaluation has a very little overhead in computational cost with respect to the Strongly Contracted one, therefore they are normally evaluated together.