# Post-SCF methods: RI-MP2

[Hartree–Fock](../10-foundations/hartree-fock.md) is a mean-field method: each electron feels only the
*average* field of the others, so it misses the way electrons **instantaneously avoid** one another. The
energy of that missing effect is the **correlation energy**, and recovering it is the job of **post-SCF**
(post-Hartree–Fock) methods. This chapter covers the family qc-rs implements today: **RI-MP2** and its
spin-scaled variants.

## Theory: correlation and second-order perturbation theory

Recall the definition from [Part II](../10-foundations/hartree-fock.md):

$$
E_{\text{corr}} = E_{\text{exact}} - E_{\text{HF}} \;<\; 0 .
$$

It is a small fraction of the total energy but *chemically decisive* — bond energies, reaction barriers, and
dispersion all live there. The cheapest systematic way to estimate it is **second-order Møller–Plesset
perturbation theory (MP2)**: treat the difference between the true electron–electron repulsion and the HF
mean field as a perturbation, and take the second-order energy correction. For a closed shell it is

$$
E^{(2)} = \sum_{ij}^{\text{occ}} \sum_{ab}^{\text{virt}}
\frac{ (ia\,|\,jb)\,\big[\,2(ia\,|\,jb) - (ib\,|\,ja)\,\big] }
     { \varepsilon_i + \varepsilon_j - \varepsilon_a - \varepsilon_b } ,
$$

a sum over occupied ($i,j$) and virtual ($a,b$) orbitals of the two-electron integrals $(ia|jb)$ divided by
orbital-energy denominators. The two terms in the bracket are the **opposite-spin** and **same-spin**
contributions — qc-rs reports them separately (`e_os`, `e_ss`), which is what makes the spin-scaled variants
below possible.

### The RI ("density-fitting") approximation

The bottleneck is the four-index integral $(ia|jb)$. The **resolution-of-identity (RI)**, or density-fitting,
approximation fits the AO pair density into an **auxiliary basis** $\{\chi_P\}$ under the Coulomb metric,

$$
(\mu\nu\,|\,\lambda\sigma) \;\approx\; \sum_{PQ} (\mu\nu\,|\,P)\,\big(V^{-1}\big)_{PQ}\,(Q\,|\,\lambda\sigma),
\qquad
V_{PQ} = (P|Q),
$$

collapsing the four-index object into a product of three-index ones through a single **whitened factor**

$$
\boxed{\,B_{P,\mu\nu} = \sum_Q \big(V^{-1/2}\big)_{PQ}\,(Q\,|\,\mu\nu)\,}
\qquad\Longrightarrow\qquad
(\mu\nu\,|\,\lambda\sigma) \approx \sum_P B_{P,\mu\nu}\,B_{P,\lambda\sigma}.
$$

$V^{-1/2}$ comes from eigendecomposing the symmetric metric $V=U\Lambda U^\top$ and dropping near-linearly-
dependent modes below a cutoff (`ri.metric_cutoff`, default $10^{-10}$) — the only approximation RI makes,
controlled entirely by the auxiliary basis and that cutoff, identically for every qc-rs backend. Crucially,
$B$ depends only on integrals, never on the density, so it is built **once** and reused; MP2 needs it in the
molecular-orbital form $(ia|jb)\approx\sum_P B_{P,ia}B_{P,jb}$ (the same $B$ that RI-JK uses for the SCF,
half-transformed to occupied×virtual indices), which is why RI is a single unified machinery
([full derivation](../10-foundations/density-fitting-ri.md)) shared by the SCF's `eri="ri-*"` and MP2's
`ric=`. This is the *only* MP2 path in qc-rs — fast and accurate, the fitting error well below chemical
accuracy — which is why you supply a **correlation-fitting auxiliary basis** with `ric=` on
`qc.chk.new(...)`.

## Usage

MP2 is an **LCT** method (`lct`, for the correlation methods that sit on top of an SCF reference). Run the
SCF first, then add the `lct` step:

```python
import qc
water = "O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469"

hf = qc.chk.new(atom=water, ao="cc-pvdz", unit="angstrom",
                ric="cc-pvdz-ri/mp2fit").scf(ref="r").run()   # note the ric= aux basis
m  = hf.lct(method="mp2").run()

m.lct.energy    # -76.230747   total MP2 energy (E_HF + E_corr)
m.lct.e_corr    # -0.203953    the correlation energy
m.lct.e_os      # opposite-spin part
m.lct.e_ss      # same-spin part
```

The `ric=` auxiliary basis (here `cc-pvdz-ri/mp2fit`, the correlation-fitting set matched to `cc-pvdz`) is
required — RI-MP2 has no non-RI fallback. Choose the aux to match the orbital basis (`cc-pvtz` →
`cc-pvtz-ri/mp2fit`, def2 sets → their `def2-*/mp2fit`).

:::{tip} How good is RI-MP2?
The RI fitting error is tiny. For this water/cc-pVDZ case, qc's RI-MP2 correlation energy (−0.203953) agrees
with a conventional (non-RI) MP2 to **1.5 × 10⁻⁵ Ha** — far below chemical accuracy (~1.6 mHa).
:::

### Spin-component-scaled variants

Because MP2 slightly overcounts correlation, two empirical rescalings of the opposite-/same-spin parts often
improve accuracy at no extra cost — just change `method`:

```python
scs = hf.lct(method="scs-mp2").run().lct.energy    # -76.226839   SCS-MP2 (Grimme)
sos = hf.lct(method="sos-mp2").run().lct.energy    # -76.224885   SOS-MP2 (same-spin dropped, O(N⁴))
```

- **`scs-mp2`** scales opposite- and same-spin by different factors (⅚ and ⅓).
- **`sos-mp2`** keeps *only* the opposite-spin term (scaled 1.3), which enables an O(N⁴) Laplace algorithm.

### Open shells and frozen core

RI-MP2 works on **RHF, UHF, and ROHF** references (an open-shell radical uses its UHF/ROHF orbitals). The
**frozen-core** approximation (excluding chemically inert core orbitals from the correlation sum, the usual
default in the field) and other tuning are `iop` keys under `lct.*` — e.g. `iop={"lct.frozen_core": True}` —
documented in the [IOP reference](../40-reference/iop.md).

:::{note} What is implemented today
The **RI-MP2 family** (`mp2` / `scs-mp2` / `sos-mp2`) computes real energies. The other `lct` methods
(`cc2`) and the multireference methods (`casscf`, `caspt2`, `nevpt2`, `fci`, `dmrg`) and excited-state `td`
are recognized but **mock** (they set up the workflow without a real correlation energy yet). This chapter
covers the methods that produce real numbers.
:::

## Worked example: the correlation energy

```python
import qc
water = "O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469"

hf = qc.chk.new(atom=water, ao="cc-pvdz", unit="angstrom",
                ric="cc-pvdz-ri/mp2fit").scf(ref="r").run()
print(f"E(RHF)     = {hf.scf.energy:.6f}")            # -76.026794

for method in ("mp2", "scs-mp2", "sos-mp2"):
    e = hf.lct(method=method).run().lct.energy
    print(f"E({method:8}) = {e:.6f}")
# E(mp2     ) = -76.230747   (E_corr = -0.203953)
# E(scs-mp2 ) = -76.226839
# E(sos-mp2 ) = -76.224885
```

The RHF energy drops by ~0.20 Ha once correlation is added — small next to the total, but this is the part
that gets bond energies and reaction thermochemistry right.

:::{exercise}
:label: ex-postscf

1. You call `hf.lct(method="mp2").run()` but `hf` was built without `ric=`. What goes wrong, and why does MP2
   need it when the SCF did not?
2. From `m.lct.e_os` and `m.lct.e_ss`, how would you reconstruct the plain MP2 correlation energy? And what
   does SOS-MP2 do to `e_ss`?
:::

:::{solution} ex-postscf
:class: dropdown

1. RI-MP2 has **no non-RI path**, so it needs a correlation-fitting auxiliary basis; without `ric=` there is
   no aux to factor the four-index integrals through. The SCF did not need one because its J/K build used the
   4-center ERIs (or the JK-fit aux) — a *different* fitting set from the MP2 correlation aux.
2. Plain MP2: `e_corr = e_os + e_ss`. **SOS-MP2 drops `e_ss` entirely** (scales it by 0) and scales `e_os` by
   1.3, which is what lets it run in O(N⁴).
:::

RI-MP2 gave us a correlated energy at one geometry. To find equilibrium structures and forces, the
[next chapter](gradients-geomopt.md) turns to analytic gradients and geometry optimization.