Geometric analysis

Geometric analysis#

Not every molecular property comes from the electron density. This last properties chapter covers qc.prop’s purely geometric descriptors — quantities computable from the nuclear positions (and, for a couple, the van der Waals radii) alone, without touching the wavefunction at all. They connect qc-rs to rotational spectroscopy, packing/solvation-accessible-surface reasoning, and radial structural analysis.

Theory: rotational constants from the moment of inertia#

A rigid molecule’s rotational energy levels — what a microwave/rotational spectrum measures — are set by its moment-of-inertia tensor, built purely from nuclear masses and positions relative to the centre of mass:

\[ I_{\alpha\beta} = \sum_A m_A\big(|\mathbf r_A|^2\delta_{\alpha\beta} - r_{A,\alpha}r_{A,\beta}\big). \]

Diagonalizing \(\mathbf I\) gives the three principal moments of inertia \(I_a\le I_b\le I_c\) (in \(\text{amu}\cdot\text{Å}^2\)), and each converts to a rotational constant

\[ X = \frac{h}{8\pi^2 c\,I_X}, \qquad X \in \{A,B,C\}, \]

conventionally reported in GHz or cm⁻¹. A linear molecule has \(I_a=0\) (no moment of inertia about its own axis) and only two independent rotational constants; an asymmetric top (the general case, e.g. water) has all three distinct; a symmetric top has two equal. These constants are exactly what a rotational or microwave spectrum measures directly — computing them from a qc-rs geometry is the bridge from “a computed structure” to “a predicted rotational spectrum,” and (with the analogous mass-weighted machinery from Hessians & thermochemistry) the rotational partition function that feeds thermochemistry.

Usage#

import qc
water = "O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469"
m = qc.chk.new(atom=water, ao="sto-3g", unit="angstrom").scf(ref="r").run()

rc = qc.prop.geom.rotational_constants(m)
rc["constants_ghz"]           # [822.15, 437.53, 285.56]  A, B, C (GHz)
rc["moments_amu_angstrom2"]   # [0.615, 1.155, 1.770]      I_a, I_b, I_c

Water’s three distinct rotational constants confirm it as an asymmetric top — consistent with its bent, low-symmetry shape. This needs only the converged geometry, not the wavefunction, so it works identically on an optimized or an as-typed structure (though only the former corresponds to a real equilibrium spectrum).

Other geometric descriptors#

Three further leaves round out qc.prop.geom, all built from nuclear positions plus van der Waals radii rather than the wavefunction:

surface_area(mychk) — the molecular (solvent-accessible-style) surface area, from the union of scaled van der Waals spheres, with a per-atom breakdown. Useful alongside dispersion/solvation reasoning, where buried vs. exposed surface matters.
free_volume(mychk) — partitions a bounding box into occupied (inside any atom’s vdW sphere) and free volume, reporting the free fraction — a packing-density descriptor relevant to crystals and cavities.
rdf(mychk) — the radial distribution of atoms around a chosen centre, binned by distance: a structural fingerprint useful for clusters and disordered/extended systems where “distance to the nearest neighbour of each type” is the natural first question.

sa = qc.prop.geom.surface_area(m)
sa["area_angstrom2"], sa["per_atom_angstrom2"]     # total + per-atom surface area

fv = qc.prop.geom.free_volume(m)
fv["free_fraction"]                                 # fraction of the bounding box that is empty

Worked example#

import qc
water = "O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469"
m = qc.chk.new(atom=water, ao="sto-3g", unit="angstrom").scf(ref="r").run()

rc = qc.prop.geom.rotational_constants(m)
print("A, B, C (GHz):", [round(x, 2) for x in rc["constants_ghz"]])
print("I_a, I_b, I_c (amu·Å²):", [round(x, 4) for x in rc["moments_amu_angstrom2"]])
# A, B, C (GHz): [822.15, 437.53, 285.56]
# I_a, I_b, I_c (amu·Å²): [0.6147, 1.1551, 1.7698]

Exercise 19

Linear CO₂ and bent SO₂ are both triatomic. Without computing anything, predict how many distinct rotational constants each will have, and why.
rotational_constants needs only the geometry — no SCF, no basis set choice beyond what you already used to optimize the structure. Why does it not need the wavefunction at all?
You want to know how “crowded” the interior of a molecular crystal’s unit cell is. Which qc.prop.geom leaf gives you a single descriptive number for that?

Solution to Exercise 19

Linear CO₂: \(I_a=0\) about its own axis, so it has only 2 distinct rotational constants (\(B=C\) from the two equal perpendicular moments). Bent SO₂: a genuine asymmetric top with all three moments distinct, so it has 3 distinct constants — exactly the water pattern in this chapter.
The moment-of-inertia tensor is built purely from nuclear masses and positions — no electron density or orbitals enter the formula at all. The geometry you feed it may have come from a full SCF optimization, but the rotational-constant calculation itself is classical mechanics on point masses.
qc.prop.geom.free_volume(mychk) — its free_fraction is exactly a single number summarizing how much of a bounding region is empty vs. occupied by atomic volume.

That completes the molecular-properties suite tour — from charges and bonds, through topology and weak interactions, to reactivity, orbital localization, electrostatics, and now geometry. Together with the SCF, correlation, gradients, and environment chapters earlier in Part III, you have the full day-to-day qc-rs toolkit.