Second-order Møller–Plesset perturbation theory (MP2)

Modules: mp, pbc.mp

The MP2 and coupled-cluster functionalities of PySCF are similar. See also Coupled-cluster theory.

Introduction

Second-order Møller–Plesset perturbation theory (MP2) [1] is a post-Hartree–Fock method. MP2 calculations can be performed with or without density fitting, depending on the initial SCF calculation.

A simple example (see examples/mp/00-simple_mp2.py) of running an MP2 calculation is

#!/usr/bin/env python
#
# Author: Qiming Sun <osirpt.sun@gmail.com>
#

'''
A simple example to run MP2 calculation.
'''

import pyscf

mol = pyscf.M(
    atom = 'H 0 0 0; F 0 0 1.1',
    basis = 'ccpvdz')

mf = mol.RHF().run()

mf.MP2().run()

which outputs

converged SCF energy = -99.9873974403487
E(MP2) = -100.198764900659  E_corr = -0.211367460310054

namely, the Hartree–Fock energy, the MP2 energy and their difference, the MP2 correlation energy.

Note

The last line in the code example above could have been replaced by

pyscf.mp.MP2(mf).kernel()

or

pyscf.mp.MP2(mf).run()

for the same result.

Spin symmetry

The MP2 module in PySCF supports a number of reference wavefunctions with broken spin symmetry. In particular, MP2 can be performed with a spin-restricted, spin-unrestricted, and general (spin-mixed) Hartree-Fock solution, leading to the RMP2, UMP2, and GMP2 methods.

The module-level mp.MP2(mf) constructor can infer the correct method based on the level of symmetry-breaking in the mean-field argument. For more explicit control or inspection, the respective classes and functions can be found in mp2.py (restricted), ump2.py (unrestricted), and gmp2.py (general).

For example, a spin-unrestricted calculation on triplet oxygen can be performed as follows:

from pyscf import gto, scf, mp
mol = gto.M(
    atom = 'O 0 0 0; O 0 0 1.2',  # in Angstrom
    basis = 'ccpvdz',
    spin = 2
)
mf = scf.HF(mol).run() # this is UHF
mymp = mp.MP2(mf).run() # this is UMP2
print('UMP2 total energy = ', mymp.e_tot)

Properties

A number of properties are available at the MP2 level.

Unrelaxed 1- and 2-electron reduced density matrices can be calculated. They are returned in the MO basis:

dm1 = mymp.make_rdm1()
dm2 = mymp.make_rdm2()

Analytical nuclear gradients can be calculated [2][3][4]

mygrad = mymp.nuc_grad_method().run()

Frozen orbitals

By default, MP2 calculations in PySCF correlate all electrons in all available orbitals. To freeze the lowest-energy core orbitals, use the frozen keyword argument:

mymp = mp.MP2(mf, frozen=2).run()

To freeze occupied and/or unoccupied orbitals with finer control, a list of 0-based orbital indices can be provided as the frozen keyword argument:

# freeze 2 core orbitals
mymp = mp.MP2(mf, frozen=[0,1]).run()
# freeze 2 core orbitals and 3 unoccupied orbitals
mymp = mp.MP2(mf, frozen=[0,1,16,17,18]).run()

Job control

Avoid t2 storage

If the t2 amplitudes are not required after the MP2 calculation, they don’t need to be saved in memory:

mymp = mp.MP2(mf)
# by default, with_t2=True
mymp.kernel(with_t2=False)

References

1

C. Møller and M. S. Plesset. Note on an approximation treatment for many-electron systems. Phys. Rev., 46:618, Oct 1934. doi:10.1103/PhysRev.46.618.

2

J. A. Pople, R. Krishnan, H. B. Schlegel, and J. S. Binkley. Derivative studies in hartree-fock and møller-plesset theories. Int. J. Quantum Chem., 16:225, 1979. doi:10.1002/qua.560160825.

3

N. C. Handy, R. D. Amos, J. F. Gaw, J. E. Rice, and E. D. Simandiras. The elimination of singularities in derivative calculations. Chem. Phys. Lett., pages 151, 1985. doi:10.1016/0009-2614(85)87031-7.

4

Y. Yamaguchi and H. F. Schaefer III. Analytic Derivative Methods in Molecular Electronic Structure Theory: A New Dimension to Quantum Chemistry and its Applications to Spectroscopy, chapter Fundamentals and Theory, pages 325. John Wiley & Sons, Ltd., 2011. doi:10.1002/9780470749593.hrs006.