3. Configuration Interaction

Configuration interaction approximates wave function by a linear expansion of N–electron basis functions with a given one–electron basis:

\[\vert \Psi_{CI}\rangle = \sum_{\lbrace i \rbrace}c_i \vert \Phi_i\rangle\]

where \(\Phi_i\) is the i–th configuration and \(c_i\) is the corresponding coefficient.

To optimize the configuration interaction wave function and minimize the energy with the normalized constrain \(\langle \Psi \vert \Psi \rangle = 1\). Variational principle is applied to varying the coefficients:

\[\mathop{min}\limits_{\{c_i\}}(\langle\Psi_{CI}\vert H \vert \Psi_{CI}\rangle-E\langle\Psi_{CI}\vert\Psi_{CI}\rangle)\]

The optimaztion using variational principle results in the Schrödinger equation becomes a matrix–eigenvalue equation:

\[\mathbf{H}\mathbf{C} = E \mathbf{C}\]

3.1. Algorithm

  • Step 1: Build the Integral

The pre–computed quantities for the SCF algorithm include:

  • The nuclear repulsion energy \(E_{nuc}\).

  • One–electron integrals overlap \(\mathbf{S}\), kinetic energy \(\mathbf{T}\), and nuclear attraction \(\mathbf{V}\).

  • Two–electron integral, the electron repulsion \((\mu\nu\vert\lambda\sigma)\).

Please check out the Hamiltonian section for further information to learn more about the pre–computed quantities.

  • Step 2: Optimization of Molecular Orbital Coefficients

The Hartree-Fock method is an uncorrelated mean-field theory that offers a qualitative description of chemical systems. Although Hartree-Fock theory is only qualitatively correct, it forms the basis for more accurate models and becomes the cornerstone of ab initio quantum chemistry.

Please check out the Hartree-Fock section for further information.

  • Step 3: Transformation of Atomic Orbital to Molecular Orbital

With the optimized LCAO_MO coefficients, we can transform the operators from the atomic orbital basis to the molecular orbital basis.

For the one electron operators:

\[h_{ij} = \sum_{\mu\nu}C_{\mu}^jC_{\nu}^i \langle\phi_{\mu}\vert\hat{h}\vert\phi_{\nu}\rangle = \sum_{\mu\nu} C_{\mu}^j C_{\nu}^i F_{\mu\nu}\]

For the two electron operators:

\[\langle pq \vert rs\rangle = \sum_{\mu} C_{\mu}^p \left[ \sum_{\nu} C_{\nu}^q \left[ \sum_{\lambda} C_{\lambda}^r \left[ \sum_{\sigma} C_{\sigma}^s \langle\mu\nu\vert\lambda\sigma\rangle\right]\right]\right]\]

  • Step 4: Construction of N–electron Basis Set

Configuration interaction method need an N–electron basis set \(\{\Phi_i\}\), here denotes explicitly a set of orthonormal Slater determinants made of a set of molecular spin orbitals \(\{\chi_i\}\).

\[\langle \Phi_i \vert \Phi_j \rangle = \delta_{ij}\]

For an N electron and K spatial orbitals system, the maximum N–electron basis set contains all Slater determinants generated by distributing N electrons among 2K spin orbitals. The number of Slater determinants gives:

\[N_{det} = \sum_{n=0}^{min(N,2K-N)}\binom{N}{n}\binom{2K-N}{n}\]

Without constraints, the number of Slater determinants is enormous. Therefore, to control the size of the N–electron basis set, a reasonable approach starts from a reference configuration. Usually, it is the ground state \(\vert\Phi_0\rangle\), which is formed by the N lowest energy spin orbitals, which results in Hartree–Fock calculation.

Then we append the other configuration by excitation level relative to the reference state. Thus, we build an N–electron basis set \(\{ \vert\Phi_0\rangle, \vert S\rangle,\vert D\rangle, \vert T \rangle \dots\}\) with the desired size. Where \(\vert S\rangle\) represents the single excitations configurations, \(\vert D\rangle\) represents the double excitations configurations, and so on.

  • Step 5: Calculation of Hamiltonian Matrix Elements

With N–electron basis set \(\{\Phi_i\}\), the Hamiltonian operator \(\hat{H}\) can be represent by a matrix,

\[\begin{split}\left[ \begin{split} \langle\Phi_0\vert\hat{H}\vert\Phi_0\rangle&\ &\langle\Phi_0\vert\hat{H}\vert\mathbf{h}\rangle\\ \langle\mathbf{h}\vert\hat{H}\vert\Phi_0\rangle&\ &\langle\mathbf{h}\vert\hat{H}\vert\mathbf{h}\rangle\\ \end{split} \right] \left[ \begin{split} 1\\ \mathbf{C}\\ \end{split} \right] = \left[ \begin{split} 1\\ \mathbf{C}\\ \end{split} \right] E\end{split}\]

where \(\vert\mathbf{h}\rangle\) indicates the set of excitation configurations.

Using the Slater–Condon rules, we can express the matrix elements \(H_{ij}\) in terms of one– and two– electron molecular integrals.

\[H_{ij}= \langle \Phi_i\vert\hat{H}\vert \Phi_j\rangle\]

  • First we evaluate the number of sphin orbitals substution between two determinants.

According to the Slater–Condon rules, all the determinants that differ by more than two spin orbitals results in zero, hence and we can only focus on diterminants that differ by zero– one– and two sphin orbitals.

  • Then evaluate the matrix elements by the number of sphin-orbitals substitution.

For identical determinants:

\[\langle\Phi_i \vert\hat{H}\vert \Phi_i\rangle = \sum_m^N \langle m \vert\hat{h}\vert m \rangle + \sum_{m>n}^N \langle mn \vert\vert mn \rangle\]

For determinants that differ by one spin orbitals:

\[\begin{split}\vert \Phi_i \rangle &= \vert \dots mn \dots \rangle\\ \vert \Phi_j \rangle &= \vert \dots pn \dots \rangle\\ \langle\Phi_i \vert\hat{H}\vert \Phi_j\rangle &= \langle m \vert\hat{h}\vert p \rangle + \sum_{n}^N \langle mn \vert\vert pn \rangle\end{split}\]

For determinants that differ by two spin orbitals:

\[\begin{split}\vert \Phi_i \rangle &= \vert \dots mn \dots \rangle\\ \vert \Phi_j \rangle &= \vert \dots pq \dots \rangle\\ \langle\Phi_i \vert\hat{H}\vert \Phi_j\rangle &= \langle mn \vert\vert pq \rangle\end{split}\]

  • Finally compute the phase factor.

Since the Slater determinant is antisymmetrized, when we evaluating the overlap between two determinants \(\langle\Phi_i\vert\Phi_j\rangle\) a phase factor.

When evaluting the Hamiltonian matrix element \(H_{ij}\), the phase is calculated as \(-1^{N_{perm}}\), where \(N_{perm}\) is the number permutations necessary to bring the spin-orbitals on which the holes are made to the positions of the praticles.

  • Step 6: Solution of the Matrix eigenvalue problem for the desired state

The optimaztion using variational principle results in the Schrödinger equation becomes a matrix–eigenvalue equation:

\[\mathbf{H}\mathbf{C} = E \mathbf{C}\]

Thus, the last step is to find the eigenvalues and eigenvectors of matrix \(\mathbf{H}\) with numerical method.

3.2. Examples

  • Full CI

The full configuration interaction(FCI) method assumes that all electrons are correlated among all orbitals in a given system. Hence it provides numerically exact solutions (within the infinitely flexible complete basis set) to the electronic time-independent, non-relativistic Schrödinger equation.

/data/examples/ci/fci.py
from moha import *

mol,orbs = IOSystem.from_file('../data/water.xyz','sto-3g.nwchem')
ham = ChemicalHamiltonian.build(mol,orbs)
wfn = HFWaveFunction(10,7,{'alpha':5,'beta':5})

scfsolver = PlainSCFSolver(ham,wfn)
hf_results = scfsolver.kernel()

cisolver = FullCISolver(ham,wfn,hf_results)
ci_results = cisolver.kernel()

  • CISD

FCI method is exact in a given atomic orbital basis but prohibitively expensive. Therefore, to balance accuracy and computational time, we truncate the CI space by excitation level relative to the reference state. The most common truncation of the CI space expansion is CI singles and doubles (CISD).

/data/examples/ci/cisd.py
from moha import *

mol,orbs = IOSystem.from_file('../data/water.xyz','sto-3g.nwchem')
ham = ChemicalHamiltonian.build(mol,orbs)
wfn = HFWaveFunction(10,7,{'alpha':5,'beta':5})

scfsolver = PlainSCFSolver(ham,wfn)
hf_results = scfsolver.kernel()

cisolver = CISDSolver(ham,wfn,hf_results)
ci_results = cisolver.kernel()