# Weak interactions (NCI, IGM)

QTAIM and ELF described the *strong* bonds within a molecule. Chemistry also runs on **weak, non-covalent
interactions** — hydrogen bonds, van der Waals contacts, steric clashes — that hold dimers, host–guest
complexes, and folded biomolecules together. They are too weak to show up as ordinary bonds, but they leave a
clear fingerprint in the density. This chapter maps them with **NCI** and **IGM**.

## Theory: interactions hide in the low-density, low-gradient regions

A non-covalent contact appears where two fragments' densities gently overlap — a region of **low density**
*and* **low reduced density gradient**. The **non-covalent interaction (NCI)** analysis exposes exactly these
regions using two fields:

- the **reduced density gradient** $s = \dfrac{|\nabla\rho|}{2(3\pi^2)^{1/3}\rho^{4/3}}$, which **spikes down
  toward zero** wherever an interaction flattens the density, and
- **$\text{sign}(\lambda_2)\,\rho$**, the density signed by the sign of the second Hessian eigenvalue
  $\lambda_2$, which tells you the interaction **type**: $\lambda_2 < 0$ = **attractive** (H-bond),
  $\lambda_2 \approx 0$ = **van der Waals**, $\lambda_2 > 0$ = **repulsive** (steric clash).

A 2-D plot of $s$ against $\text{sign}(\lambda_2)\rho$ shows a spike per interaction, coloured by type. The
**independent gradient model (IGM)** isolates the interaction signal more directly, by comparing the *true*
gradient to a **hypothetical non-interacting** one built from the same fragments:

$$
g^{\text{IGM}} = \Big\|\Big(\sum_A|\partial_x\rho_A|,\ \sum_A|\partial_y\rho_A|,\ \sum_A|\partial_z\rho_A|\Big)\Big\|,
\qquad
\delta g = g^{\text{IGM}} - |\nabla\rho| .
$$

Because $|\nabla\rho|=\big|\sum_A\nabla\rho_A\big|$ while $g^{\text{IGM}}$ sums each fragment's gradient
*component-by-component before* taking the magnitude, the triangle inequality guarantees
$g^{\text{IGM}}\ge|\nabla\rho|$ always — with equality wherever the fragments' gradients do not oppose each
other. $\delta g$ is therefore **exactly the cancellation** between fragments' gradients — nonzero only where
two atomic densities push against one another, which is precisely an interaction region. **IGMH** applies the
same idea with Hirshfeld-partitioned real (not promolecular free-atom) densities, giving a sharper,
SCF-density-based signal.

## Usage

These are `qc.prop.mesh` fields — computed on a grid around the molecule. A **hydrogen-bonded water dimer**
is the canonical example:

```python
import qc
dimer = ("O 0 0 0; H 0.76 0 0.59; H -0.76 0 0.59; "
         "O 0 0 2.98; H 0 0.76 3.57; H 0 -0.76 3.57")
m = qc.chk.new(atom=dimer, ao="cc-pvdz", unit="angstrom").scf(ref="r").run()   # E = -152.060006

nci = qc.prop.mesh.nci_data(m)     # {'sign_lambda2_rho': ..., 'rdg': ...}  the two NCI axes
igm = qc.prop.mesh.igm(m)          # {'origin','spacing','shape','rho','delta_g','sign_lambda2_rho','rdg'}
igmh = qc.prop.mesh.igmh(m)        # the Hirshfeld-partitioned IGM
```

`nci_data` returns the two point-wise fields (`sign_lambda2_rho`, `rdg`) that make the NCI scatter plot;
`igm`/`igmh` add the gridded density, the interaction descriptor `delta_g`, and the grid geometry
(`origin`/`spacing`/`shape`) for a 3-D isosurface.

### Seeing it

The whole point is the picture. From the [visualization chapter](../visualization.md):

```python
m.plot_nci()             # the 2-D NCI diagram: s vs sign(λ₂)ρ, with the interaction spikes
m.view3d("nci")          # the 3-D NCI isosurface, coloured by interaction type
```

For the water dimer, `plot_nci()` shows a spike at negative $\text{sign}(\lambda_2)\rho$ — the signature of
the **O–H···O hydrogen bond** — plus a broad van der Waals feature near zero.

## Intrinsic bond strength (IBSI)

To put a *number* on a contact's strength, `qc.prop.bond.ibsi` computes the **intrinsic bond-strength index**
(built on the IGM $\delta g$), which correlates with interaction energy across bond types — covalent and
non-covalent alike:

```python
qc.prop.bond.ibsi(m)     # per-atom-pair intrinsic bond-strength indices
```

## Worked example

```python
import qc
dimer = ("O 0 0 0; H 0.76 0 0.59; H -0.76 0 0.59; "
         "O 0 0 2.98; H 0 0.76 3.57; H 0 -0.76 3.57")
m = qc.chk.new(atom=dimer, ao="cc-pvdz", unit="angstrom").scf(ref="r").run()

print("dimer energy:", round(m.scf.energy, 6))     # -152.060006
nci = m.prop.mesh.nci_data()
print("NCI fields   :", list(nci.keys()))          # ['sign_lambda2_rho', 'rdg']
m.plot_nci()                                        # visualize the H-bond spike
```

:::{exercise}
:label: ex-nci

1. In an NCI plot, an interaction shows a spike at $\text{sign}(\lambda_2)\rho \approx -0.03$ a.u. Attractive
   or repulsive? What kind of contact is it likely to be?
2. Why do NCI/IGM need a *grid* (making them slower than a Mulliken charge), while a bond order does not?
3. You have two conformers of a molecule and want to know which has a stronger internal hydrogen bond. Name
   two `qc.prop` tools from this chapter you could use, and what each would tell you.
:::

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

1. **Attractive** ($\lambda_2 < 0$). A spike at low negative density is the signature of a **hydrogen bond**
   (a stronger attractive contact than a pure van der Waals interaction, which sits nearer zero).
2. NCI/IGM are defined **point-by-point in real space** — they evaluate $\rho$, $\nabla\rho$, and the Hessian
   on a grid to find the low-density/low-gradient regions. A Mayer/Wiberg bond order is an **algebraic**
   contraction of the density matrix, needing no grid.
3. `m.plot_nci()` (compare the depth/position of the H-bond spike between conformers) and
   `qc.prop.bond.ibsi(m)` (a numerical intrinsic bond-strength index for the H-bonded atom pair) — the
   conformer with the deeper spike / larger IBSI has the stronger hydrogen bond.
:::

Next, [aromaticity](aromaticity.md) quantifies a different collective property — the special stability of
conjugated rings.