2. Examples
This chapter introduces the primary use of MoHa calculation, and you can find those examples in the moha/data/examples directory. The accompanying user guide will present additional details on the various modules.
2.1. Hartree-Fock
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.
2.1.1. Plain Solver
A basic Hartree-Fock calculation with a plain solver can be performed by:
import moha
import numpy as np
mol,orbs = moha.IOSystem.from_file('../data/water.xyz','sto-3g.nwchem')
ham = moha.ChemicalHamiltonian.build(mol,orbs)
wfn = moha.HFWaveFunction(10,7,{'alpha':5,'beta':5})
hf_results = moha.PlainSCFSolver(ham,wfn).kernel()
2.1.2. DIIS Solver
The convergence of the plain Hartree-Fock solver is relatively slow. One of the most effective methods for speeding up such iterative sequences is the “Direct Inversion in the Iterative Subspace(DIIS)” approach.
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})
hf_results = DIISSCFSolver(ham,wfn).kernel()
2.2. Post-Hartree-Fock
Quantum chemistry is mature enough to have canonical approaches followed by the Hartree-Fock method. Their behavior, successes, and most importantly, their failures are widely known and accepted. With a couple of exceptions which we will not mention, they are:
Configuration Interaction(CI)
Coupled-Cluster(CC)
Perturbation Theory(PT)
2.2.1. Configuration Interaction
Configuration interaction approximates wave function by a linear expansion of N-electron basis functions made of a given one-electron basis. With CI wave function and CI basis set, the Schrödinger equation becomes a matrix-eigenvalue equation.
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.
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).
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()
2.2.2. Coupled-Cluster
Coupled cluster methods are among the most accurate electronic structure methods available today. Its wave function ansatz gives in exponential form \(e^{\hat{T}} \vert\Phi_0\rangle\)
Where the cluster operator \(\hat{T}\) is the sum of \(\{\hat{T}_{i}\}\) generates all possible determinants having the i-th excitations from the reference \(\vert\Phi_0\rangle\).
CCSD
Truncating the cluster operator \(\hat{T}\) at doubles leads to the coupled cluster singles and doubles method(CCSD). By keeping only \(\hat{T}_1\) and \(\hat{T}_2\), we reach the balance between accuracy and computational time.
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()
ccsolver = CCSDSolver(ham,wfn,hf_results)
cc_results = ccsolver.kernel()
CCSD(T)
Based on CCSD, the addition of perturbative triple correction gives the CCSD(T) method the “gold standard” in computational chemistry. It is now one of the most accurate methods applicable to reasonably large molecules.
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()
ccsolver = CCSD_TSolver(ham,wfn,hf_results)
cc_results = ccsolver.kernel()
2.2.3. Perturbation Theory
Perturbation theory is a collection of versatile methods used in many branches of science. It divides the system into a model part \(\hat{H}_0\) which is a known approximation to the real system \(\hat{H}\) and a perturbation part \(\lambda\hat{V}\), where \(\lambda\) is a parameter that is small enough to guarantee convergence.
Møller–Plesset perturbation theory is a particular case of perturbation theory, where we take the Fock operator \(\hat{F}\) as the model operator.
MP2
Depending on the choice of Møller–Plesset perturbation operator \(\hat{V}\), the first order correction energy \(E_{MP1}\) might be zero or none zero, but \(E_{MP0} + E_{MP1}\) is always equal to the Hartree–Fock energy \(E_{HF}\). Hence, the first meaningful correction in Møller–Plesset perturbation theory is second-order energy.
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()
ptsolver = MP2Solver(ham,wfn,hf_results)
pt_results = ptsolver.kernel()
MP3
After the second-order energy, the involvement of third-order correction still improves the Hartree-Fock method with cheap costs.
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()
ptsolver = MP3Solver(ham,wfn,hf_results)
pt_results = ptsolver.kernel()
2.3. Property
This section contains information about atomic and molecular properties implemented in MoHa.
2.3.1. Population Analyasis
Population analysis is the study of charge distribution in molecules. It requires the overlap integrals and the electron density and information about the number of basis functions centered on each atom.
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})
hf_solver = PlainSCFSolver(ham,wfn)
hf_results = hf_solver.kernel()
pa_m_results = PopulationAnalysisMulliken(mol,orbs,ham,wfn).kernel()
pa_l_results = PopulationAnalysisLowdin(mol,orbs,ham,wfn).kernel()
2.3.2. Multipole Moment
Multipole moments are the coefficients of a series expansion of a potential due to continuous or discrete sources. A multipole moment usually involves powers (or inverse powers) of the distance to the origin and some angular dependence.
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})
hf_solver = PlainSCFSolver(ham,wfn)
hf_results = hf_solver.kernel()
mm_results = MultipoleMoment(mol,orbs,wfn).kernel([[1,0,0],[0,1,0],[0,0,1]])
2.3.3. Excitation Energy
The simplest ab initio method to calculate excited electronic states is configuration interaction singles (CIS). The fundamental idea behind CIS is the representation of the excited-state wave functions as linear combinations of singly excited determinants relative to the Hartree-Fock reference wave function.
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})
hf_solver = PlainSCFSolver(ham,wfn)
hf_results = hf_solver.kernel()
ee_cis_results = ExcitationEnergyCIS(ham,wfn).kernel()