A comprehensive list of theoretical approximations that are used in computational chemistry

I am having trouble recognizing all the approximations that are used in computational chemistry. I would like to start an answer list (similar as to the list of resources for learning chemistry) that addresses this question.

I will try to answer this question myself. I am hoping that anyone who is knowledgeable in the topic will correct and/or contribute. I am planning to start from general (e.g. Born Oppenheimer, LCAO) to specific (e.g. pseudopotential, functional). I am also planning to include why this approximation is necessary.

• Asking for a comprehensive list makes this question impossible to answer between reasonable limits, as it effectively asks for a list that includes all the tricks used from integration, basis set approximations, through diagonalization to whatever that is included every single QChem software that meant for production and not just simple demonstrations. – Greg Sep 25 '16 at 17:34
• @Greg I will try to avoid the math. I know that it is important and that it is sometimes difficult to separate the physics/chemistry from math. Otherwise it is too large, as you say. For example, I will not delve into how Frank Boys figured out that STO could be represented with GTO. Nor will I explain the different types of diagonalization (conjugate gradient, Davidson, etc). I will try not go to basis sets either unless there is a physical/chemical reasoning. Additionally, the reason why some basis sets are better or worse can be illustrated by first showing what has been sheared away first. – CoffeeIsLife Sep 26 '16 at 2:26
• Overlooked in answers here is the reality of classical simulation, encompassing huge fields of monte carlo and molecular dynamics, which seeks to define mathematical fits to highly accurate quantum phenomena, and thereby runs a lot faster than solving the Schrodinger equation. There is no fundamental theory that proves this can work in all cases, but it has doubtless been useful: github.com/khavernathy/mcmd – khaverim Oct 17 '17 at 18:46
• @khaverim I suspect that this bias is primarily due to the site's users; most of us use or develop only the electronic structure theory side of things, so we have far fewer experts in MM/MD/MC, which means many related questions go completely unanswered. – pentavalentcarbon Oct 17 '17 at 20:40
• @khaverim if you have anything to add, you can always post it as a new answer. Even if you don't think you have a full answer, you could start a community wiki answer that could help this post become a more comprehensive reference for other users. – Tyberius Oct 17 '17 at 23:27

The goal of computational chemistry is to obtain the properties of a system. This is done by solving Dirac's Equations.

Treating particles as point particles with mass

In most computational software, particles are treated as points with some mass. Neutrons and protons may be lumped into a nucleus. However, this is not true in all cases as pointed out in this question: Does computational chemistry include neutrons in ground state calculations? .

Electrons, protons, and neutrons are not simply point particles with some mass and charge. The theory and forces can get quite complex. I haven't been able to find journals that address the effect of this approximation. But I will try to post something if I come across it.

Neglecting Relativistic Effects

Wien2k has a nice summary of relativistic effects that need to be considered [8]. The relativistic effects that need to be included are:

1. Mass-velocity correction
2. The Darwin term
3. Spin-orbit coupling
4. Indirect relativistic effect

However, relativistic effects usually become important for elements further down the fifth row in the periodic table. This is because relativistic effects are dependent on nuclear charge. The velocity of the electron as well as spin orbit coupling increases as the nuclear charge increases.

This does not mean that relativistic effects do not affect "light" atoms. A good example that shows this is the sodium doublet which is a result of spin-orbit coupling. Nevertheless, solving Schrodinger's equation is enough in many cases.

Most computational chemistry software have the option of scalar relativistic calculations. The scalar relativistic technique was developed by Koelling and Harmon. These are calculations that include the Darwin term and the mass velocity correction. According to Koelling and Harmon's paper:

We present a technique for the reduction of the Dirac equation which initially omits the spin-orbit interaction (thus keeping spin as a good quantum number), but retains all other relativistic kinematic effects such as mass-velocity, Darwin,and higher order terms.[9]

Some offer fully relativistic calculations, which include all (or most) relativistic effects. But these are rarer to find, and are only available for certain cases.

So why are relativistic effects neglected? Well, relativistic effects are small for light molecules and relativistic calculations are expensive

This is so because relativistic calculations need self-consistent solutions of about twice as many states as non-relativistic ones. [10]

Thus, scalar relativistic calculations offer a nice middle ground between efficiency and accuracy.

The Born-Oppenheimer Approximation

I am not sure what motivated Born and Oppenheimer to use this approximation. The seminal paper is in German.[1] It appears that the motivation was to simplify the Hamiltonian. According to the introductory course at MIT,[2]

it allows one to compute the electronic structure of a molecule without saying anything about the quantum mechanics of the nuclei

And according to Wikipedia,[3] it reduces the amount of computations:

For example, the benzene molecule consists of 12 nuclei and 42 electrons. The time independent Schrödinger equation, which must be solved to obtain the energy and wavefunction of this molecule, is a partial differential eigenvalue equation in 162 variables—the spatial coordinates of the electrons and the nuclei.

As many approximations, it does not hold true in all scenarios. Cases where the Born-Oppenheimer approximation fail are:[4]

ion-molecule, charge transfer, and other reactions involving obvious electronic curve crossings

Qualitatively, the Born-Oppenheimer approximation says that the nuclei are so slow moving that we can assume them to be fixed when describing the behavior of electrons. Mathematically(?), the Born-Oppenheimer approximation allows to treat the electrons and protons independently. This does not imply that the nuclei and electrons are independent of each other. In other words, it does not mean that the nuclei are not influenced by the motion of electrons. The nuclei still feel the motion of the electrons. In addition, the Born-Oppenheimer approximation does not say that the nuclei does not move. It only means that when describing the motion of electrons, we assume that the nuclei are fixed.

No Analytical Solution to Dirac/Schrodinger Equation

Unfortunately there is no analytic solution to the Dirac equation for any atom that has more than one electron even after the Born-Oppenheimer approximation (a list of quantum-mechanical problems that have an analytical solution is available on Wikipedia[5]). Many texts state that the reason as to why the Schrödinger equation is not exactly solvable for more than one electron is due to the Coulomb repulsion between electrons.[6]
However, this is not entirely true. A counterargument is Hooke's atom. The Hamiltonian for Hooke's atom has an Coulomb electron-electron repulsion term. However, it has an exact solution for more than one electron under certain circumstatnces.[7]

The true reason as to why the Schrödinger equation is not solvable for multi-electron atoms is due to the fact that the motion of electrons cannot be decoupled from each other. In other words, the Hamiltonian is not separable for a multi electron system. If we were to get rid of the electron-electron Coulomb repulsion, the motion of the electrons can be decoupled. This may be the reason as to why the electron-electron Coulomb repulsion (a.k.a. electron correlation) is used as the reason why the Schrödinger equation is not exactly solvable.

From non-interacting to the real thing

Since the Dirac equation cannot be solved analytically, we must make models that are solvable and add approximations to it. These approximations are further refined to get more accurate results.

The most simple model (and the foundation for computational chemistry) is the system of non-interacting electrons. As the name suggests, the electrons do not interact with other electrons. This allows us to write the Hamiltonian for all electrons as the sum of one-electron Hamiltonians.
The one electron Hamiltonian consists of a kinetic energy term and a potential energy term.

$$\mathcal{H}=\sum^N_i \left(\frac{\hbar^2}{2m_i\nabla^2 }+V_i\right)$$ The solutions for the non-interacting system of electrons are analytical. And they give a starting point for other calculations. Hartree-Fock, DFT, and solid state physics use this simplified model.

Note: it appears that there are several terms that say something similar to this. There is the independent electron approximation and the central field approximation. I'll dig into literature to see what are the differences between these three terms.

Variations of the non-interacting system of electrons

There are several variations of non-interacting electron systems. In all these, the electrons do not interact with each other. What sets them apart from each other is the potential energy term of the Hamiltonian $$V_i$$. This potential energy term is often called the effective potential.

Free electron model

This model is used as the starting point in solid state physics. The free electron model describes the behavior of valence electrons in a metal or semiconductor. Here, the potential $$V_i$$ is equal to 0 for all electrons. It does not feel the effect of other electrons or nuclei. The wavefunctions of this model are plane waves.

Nearly free electron model

The nearly free electron includes a weak periodic potential. In other words $$V_i(r)=V_i(r+k)$$ It is used to describe periodic systems such as ideal crystals. The wavefunctions of the nearly free electron model are Bloch waves. Bloch waves are plane waves that are multiplied by a periodic function.

Hartree method

In the Hartree method, the potential considers repulsion between two electrons and attraction to the nucleus. In this model,

Hartree assumed that each electron moves in the averaged potential of the electrostatic interactions with surrounding electrons [11]

Hartree replaced the electron-electron interaction with an effective potential that only depended on the coordinates of the $$i$$th electron. This effective potential describes an electron interacting with an electron cloud.

From one-electron orbitals to the multielectron orbital

As mentioned above, the total Hamiltonian can be approximated as the sum of one-electron Hamiltonians for non-interacting electrons. However, the wavefunction has to satisfy antisymmetry.

Exact exchange

Electron correlation

Some define electron correlation as everything that the Hartree-Fock method leaves out.

Periodic systems and pseudopotentials

References

1. M. Born and R. Oppenheimer, Ann. Phys. 1927, 389, 457–484.
doi: 10.1002/andp.19273892002
2. Born Oppenheimer Approximation. Open Courseware MIT:Introductory Quantum Mechanics. Fall 2005. Section 12 Lecture. (https://ocw.mit.edu/courses/chemistry/5-73-introductory-quantum-mechanics-i-fall-2005/lecture-notes/sec12.pdf) (https://ocw.mit.edu/courses/chemistry/5-73-introductory-quantum-mechanics-i-fall-2005/)
3. Born-Oppenheimer approximation. Wikipedia (https://en.wikipedia.org/wiki/Born%E2%80%93Oppenheimer_approximation)
4. L. J. Butler, Annu. Rev. Phys. Chem. 1998, 49, 125-71.
PMID: 15012427 doi: 10.1146/annurev.physchem.49.1.125
5. List of quantum-mechanical systems with analytical solutions. Wikipedia
6. LibreTexts: 9.1: The Schrödinger Equation For Multi-Electron Atoms
7. Hooke's atom. Wikipedia
8. Summary of relativistic effects. Wien2k (http://www.wien2k.at/reg_user/textbooks/WIEN2k_lecture-notes_2013/Relativity-NCM.pdf)
9. A technique for relativistic spin-polarised calculations. Journal of Physics C: Solid State Physics, Volume 10, Number 16 (http://iopscience.iop.org/article/10.1088/0022-3719/10/16/019/meta)
10. The Scalar Relativistic Approximation. Takeda, T. Z Physik B (1978) 32: 43. doi:10.1007/BF01322185 (http://link.springer.com/article/10.1007/BF01322185)
11. Density Functional Theory in Quantum Chemistry. Tsuneda, T. 2014. ISBN: 978-4-431-54824-9. Page 36.
• @Rodriguez I will look into that. Also, this is an ongoing answer. It is still not complete! Nevertheless I would appreciate any feedback. – CoffeeIsLife Sep 25 '16 at 16:52
• 1. Schroedinger equation is approximation, better go with Dirac to have relativity (and spin). 2. Born-Oppenheimer is two-step process. First is separation of total wavefunction into product of nuclear and electronic part, second is setting the nuclear kinetic energy to zero. (see en.wikipedia.org/wiki/Born%E2%80%93Oppenheimer_approximation) – ssavec Sep 25 '16 at 18:55
• Please use manual markup for links, i.e. [link](http://...). It looks better, is easier readable, and will not lead into nirvana if interpreted wrong by the page. Just for funs, you can check what I mean by clicking the links in revision 8. – Martin - マーチン Sep 26 '16 at 8:34
• 1. As @ssavec already mentioned, the first usual approximation is neglecting relativistic effects. 2. Big mistake in the very first sentence: the goal is to obtain properties of a molecular system. Then, depending on the formalism used (MM, WFT, DFT, ...) you are looking for a mathematical entity that can somehow spits out the values of the properties. And this entity is not always the wave function. 3. Don't get me wrong, but if you don't quite understand the purpose of arguably the central approximation in present day QC (BO one), it might be a bit to early for you to answer the question. – Wildcat Sep 26 '16 at 10:55
• @QuantumAMERICCINO There is some stuff I would explain differently. I think I will write something about the Variational Principle and (Post)-Hartree-Fock as a new answer. Maybe even some DFT. Then we can see how to connect it best to your content. – Feodoran May 5 '17 at 8:44

Variational Principle

The Variational Principle states that any approximate solution to the true wave function $\Psi$ will be higher in its total energy. The total energy is calculated as the expectation value of the Hamiltonian $E=\langle\Psi|\hat H|\Psi\rangle$. Thus we can set up a trial wave function $\tilde\Psi$ and vary its parameters until its corresponding energy $\tilde E=\langle\tilde\Psi|\hat H|\tilde\Psi\rangle$ reaches a minimum. This solution provides an upper limit for the true (exact) energy ($\tilde E \le E$) as well as the best possible approximation to the wave function.

When comparing variational methods we can directly judge which result is better by comparing their total energy. Examples for variational methods are Hartree-Fock and Configuration Interaction. Non-variational methods are Coupled Cluster, Density Functional Theory and Pertubation Theory.

Note that the total energy by itself has not much physical meaning and depends on different numerical parameters. Of practical interest are always difference between total energies, for example binding energies.

Hartree-Fock

The main issue of solving the Schrödinger is the Coulomb interaction between electrons. There is no general analytic solution for more than 1 electron. However, for the non-interacting system, which ignore electron-electron interaction completely, one can employ separation of variables for the electronic coordinates. This allows us to express the total electronic wave function as a product of one-electron wave functions (orbitals). This is called the Hartree-Product and when additionally accounting for the anti-symmetry as required by the Pauli principle, it leads to the Slater-Determinant.

In the Hartree-Fock method the Slater-Determinant, which is only exact for the non-interacting system, is combined with the Hamilonian including the electron-electron interation term. Further applying the Variational Principle leads to the Hartree-Fock equations, which may be solved numerically. The physical interpretation of this approximation is, that the electrons feel each other only in an averaged way, hence the name mean-field approximation.

For the numerical solution one needs to diagonalize the Fock matrix resulting in its eigenvalues and eigenvectors. Since the Fock matrix depends on its own eigenvectors, this is solved iteratively: After choosing an initial guess the result is used to update the Fock matrix, which is then diagonalized again. After each iteration an improved result is obtained. This is known as self-consistent field method.

The SCF algorithm scales with $M^3$, where $M$ is the number of basis sets. However, the preceding calculation of all required one- and two-electron integrals is usually more time consuming than the actual HF calculation.

Electron Correlation & Post-Hartree-Fock

The correlation energy is defined as the difference between the Hartree-Fock energy and the exact energy. Electron Correlation is the error made in the Hartree-Fock approximation, since the Ansatz made for the electronic wave function does not reflect its true form. The numerically exact solution can be obtained with the Full Configuration Interaction method.

Configuration Interaction

The electronic wave function, like any other wave function in quantum mechanics, can be expanded in an arbitrary basis set. In Configuration Interaction it is recognized that the electron configurations, which can be created based on the molecular orbitals obtained from a Hartree-Fock calculation, can be used as such a basis.

Note that this basis of electron configurations is a different one, than the one-electron basis set which is used to to represent the molecular orbitals.

A Configuration Interaction calculation is thus preceded by a Hartree-Fock calculation. The Slater-Determinant used in the Hartree-Fock calculation is the first, and usually most important, configuration. Further excited determinants (configurations) are generated by promoting electrons from occupied orbitals to virtual (in the HF configuration unoccupied) ones. Applying again the Variational Principle one can calculate that linear combination of such configurations that minimized the electronic energy.

Computationally one needs to set up the matrix representation of the Hamiltonian (in the basis of the chosen basis of the configurations) and diagonalize it. In contrast to the Fock matrix it does not depend on its own solution, so no iterative procedure is required here. However, this is commonly speed up by iterative approaches which are much less time consuming than an exact diagoanlization.

In total $\binom{M}{N}$ configurations can be generated, where $M$ is the number of orbitals and $N$ is the number of electrons. Since this number increases exponentially, only the smallest molecules can be calculated even on the fastest super computers (of course this also depends on various other numerical parameters). Therefore different truncation schemes exist to only include some of the configurations.

• Full CI includes all possible configurations and is the numerically exact limit.
• Hartree-Fock from the perpective of CI it is known the one-determinant approximation.
• CIS includes only configrations where 1 electron is excited with respect to the HF configuration (Single excitations). Due to Brillouin's theorem, no improvement of the HF ground state energy is made, but rough approximations for excited states can be obtained.
• CID includes only configrations where 2 electrons are excited with respect to the HF configuration (Double excitations). First improvements to the ground state energy.
• CISD since Single excitation couple with Double and those in turn with the HF ground state, this is an improvement over CID. Furthermore, there are much less Single excitations than Double excitations, which means computationally CISD has no considerable additional efford over CID.
• CISDT... The higher the excitation degree the smaller the correction to the electronic energy, but the higher the computational demand.
• CASCI stands for Complete Active Space CI. The active space is a chosen set of orbitals within which all possible configurations are considered. Orbitals with energy below the active space orbitals are (doubly) occupied for all configurations, orbitals above in energy are always left empty.

Multi-configurational self-consistent field

The MCSCF approach combines a CI calculation with the HF method in the sense that the coefficients of the chosen configurations (CI optimization) and the molecular orbital coefficient (SCF optimization) are optimized simultaneously. A common example is the CASSCF method.

Coupled Cluster

Coupled Cluster is a reparametrization of the electronic wave function. In the CI approach the wave function can be written as an excitation operator $\hat T$ acting directly on the HF determinant and creating all possible excited configurations. The CC wave function is obtained by having the exponential operator $\exp(\hat T)$ acting on the HF determinant.

The effect is that for example CCSD will not only include the Single and Double excitations, but certain Triple and Quadruple excitations as well, making the method size consitent.

Unfortunately the arising Coupled Cluster equations cannot be solved by employing the Variational Principle, hence the resulting total energy may be lower than the exact FCI energy. This also means that a CC calculation including all possible configurations is not equivalent to FCI. Although CCSD seems to be an exception here. CC is computationally more expensive than a CI calculation of the same excitation operator, but results are more accurate.

CCSD(T) is considered the gold standard of quantum chemistry and includes pertubative estimates of the Triple excitations. It's overall scaling is in the order of $N^7$.

Multi Reference Methods

Instead of creating excitated configurations only based on the HF configuration, one can also have multiple reference configurations. Commonly one starts with a MCSCF wave function and its configurationsi as reference. In a second step one can do for example do a MRCI or MRCC calculation.

This approach is required for strongly correlated systems, where other methods (DFT, CCSD(T)) are failing.

References

• A. Szabo and N. S. Ostlund. Modern Quantum Chemistry: Introduction to Advanced Electronic Structure Theory (Dover Books on Chemistry). Dover Publications, 1996.
• T. Helgaker, P. Jorgensen, and J. Olsen. Molecular electronic-structure theory. 1. Auflage. New York: Wiley, 2000.

Speeding up two-electron integral evaluation through approximate methods

Remarks on Notation (using Mulliken Notation)

• $\mu,\nu,\kappa,\lambda:$ atomic orbital (AO) basis functions
• $P,Q,R,S:$ auxiliary basis functions
• $(\mu\nu|\kappa\lambda) = \int\int\phi_\mu^*(r_1)\phi_\nu(r_1)\frac{1}{r_{12}} \phi_\kappa^*(r_2) \phi_\lambda(r_2)d\tau_1d\tau_2 = \left<\mu\kappa|\nu\lambda \right>$
• $(\mu\nu|1|P) = (\mu\nu P)$
• $(\mu\nu|g_{12}|P) = (\mu\nu|P)$ where $g_{12}$ is the two-electron interaction, almost always $= \frac{1}{r_{12}}$
• $S_{PQ} = (PQ)$
• $V_{PQ} = (P|Q)$

Resolution of the Identity approximations

Computation of 2-electron 4-centre Integrals $(\mu\nu|\lambda\sigma)$ can be a significant bottleneck in electronic structure calculations.

$\rightarrow$ Idea of RI is to avoid such integrals

Proposed solution: Reexpansion of pair products $|\mu\nu)$ with

$$|\mu\nu) \approx |\widetilde{\mu\nu}) = \sum_P C_{\mu\nu}^P |P)$$

Where we choose the auxiliary basis set $\{P\}$ to be $\sum_P |P)(P| \approx 1$ (Hence the name Resolution of the Identity)

If we approximate the unity operator $\hat{1}$ in the Coulomb metric $\hat{1} \approx \sum_{PQ} |\chi_P)(P|Q)^{-1} (\chi_Q|$

The integral $(\mu\nu|\kappa\lambda)$ can then be approximated as $(\mu\nu|\kappa\lambda) \approx \sum_{PQ}(\mu\nu|P)(P|Q)^{-1}(Q|\kappa\lambda)$.

The same formulation is obtained if the the integral error $\Delta (\mu\nu - \widetilde{\mu\nu}|\kappa\lambda - \widetilde{\kappa\lambda})$ is minimized with respect to the coefficients $C_{\mu\nu}^P$.

This leads to a reduction of computation cost since we need only to compute three-index integrals. These are much faster to compute than four-index integrals. Additionally, a lot fewer three-index index integrals have to be computed than four-index integrals would be.

The RI approximation is used to accelerate the calculation of the Coulomb operator $\hat{J}$ in the RI-J approximation Here the electron density is expanded in the auxiliary set ($D_{\kappa\lambda}$ is the density matrix):

$$J_{\mu\nu} = (\mu\nu|\rho) \approx \sum_{PQ} \sum_{\kappa\lambda} (\mu\nu|P)(P|Q)^{-1}(Q|\kappa\lambda)D_{\kappa\lambda}$$

If you implement this in your program in such a way that you calculate all the quantities strictly from right to left, the scaling is brought down from $N^4$ to $N^3$.

It can also be used to approximate the Exchange operator $\hat{K}$. Here no reduction of the scaling is possible, but the prefactor is reduced.

RI approximations are essential also for MP2/CC2-R12/F12 methods. There you need the robust expression $(\mu\nu|\kappa\lambda) \approx (\widetilde{\mu\nu|\kappa\lambda}) = (\widetilde{\mu\nu}|\kappa\lambda) + (\mu\nu|\widetilde{\kappa\lambda}) - (\widetilde{\mu\nu}|\widetilde{\kappa\lambda})$ for example for the operators $r_{12}$ and $[\hat{T}_{12},r_{12}]$.

The auxiliary basis sets can be used freely, but they have to be fitted for their purpose (meaning there are individual auxiliary basis sets for every operator, often denoted $jbas$ (for pair products $(ij|$), $cbas$ (for pair products $(ia|$), $i,j$ being occupied, $a$ being virtual orbitals. I know of no aux basis sets for doubly occupied pair products $(ab|$, so there are certain integrals where RI/DF can't be used afaik.

Cholesky Decomposition

The two-electron integrals $(\mu\nu|\kappa\lambda)$ can be expressed in a positive definite Hermitian $V_{\mu\nu,\kappa\lambda}$.

Since this matrix is positive definite it can be decomposed via a Cholesky decomposition $\mathbf{V} = \mathbf{LL}^{\dagger}$. For a detailed explanation read the original publication by Beebe and Linderberg.

In this approximation the basis set $\{P\}$ is conceived via the orbital basis, so no generation of a basis set is needed.

While CD is computationally more demanding than RI, it is numerically robust and can be used for more precise results. Since the only constricting value is the threshhold at which the expansion is stopped, the precision of the CD in theory can be chosen arbitrarily and systematically.

References

• P. Merlot, T. Kjaergaard, T. Helgaker, R. Lindh, F. Aquilante, S. Reine, T. B. Pedersen, J. Comp. Chem. 2013, 34, 1486 - 1496

• O. Vahtras, J. Almloef. Chem. Phys. Letters 1993, 213, 514-518

• F. Weigend, M. Kattanek, R. Ahlrichs, J. Chem. Phys., 2009, 130, 164106-1

• N. H. F. Beebe, J. Linderberg, Int. J. Q. Chem., 1977, 12, 683-705