# Symmetry & group theory

Molecular symmetry is not a bookkeeping convenience — it is a mathematical guarantee that lets you
**block-diagonalize** the SCF eigenvalue problem before you ever build a Fock matrix, cutting a size-$n$
diagonalization into several much smaller ones and, as a side effect, giving every molecular orbital a
label (`A1`, `B2`, …) instead of a bare index. This chapter derives the machinery — the group-theoretic
projector, the symmetry-adapted linear combination (SALC) basis, and its interaction with overlap
orthogonalization — that [the SCF chapter's `symmetry=` option](../20-guide/scf.md) invokes, grounded in
`.design/24qc.scf-salc.md`.

## Point groups and the symmetry operation on AOs

A rigid molecule's **point group** $G$ is the finite set of proper/improper rotations that map the nuclear
framework onto itself. Every element $g\in G$ acts on 3-D space by an orthogonal matrix $R(g)$, and this
same $g$ permutes atomic centers: $A\mapsto gA$. To use $G$ in the SCF, you need to know how $g$ acts on a
contracted Gaussian AO, not merely on a point in space.

An AO is a product of a radial Gaussian and an angular part centered on atom $A$:

$$
\chi_{A,l,k,m}(\mathbf r) = N_{l,k}\, p_{lm}(\mathbf r - \mathbf R_A)\, e^{-\alpha_k|\mathbf r-\mathbf
R_A|^2},
$$

where $l$ is the angular momentum, $k$ the contraction/primitive index, and $m$ the angular function index
($m=-l,\dots,l$ for spherical harmonics, or a Cartesian monomial index). Because $R(g)$ is a rigid rotation
of space, it maps the radial Gaussian to itself (rotation-invariant) and mixes only the angular functions
of the same $l$ on the image center $gA$:

$$
\hat U(g)\,\chi_{A,l,k,m} = \sum_{m'} T^{(l)}_{m'm}(g)\,\chi_{gA,l,k,m'}.
$$

$T^{(l)}(g)$ is the **Wigner-style representation matrix** of $g$ restricted to angular momentum $l$ — for
$l=1$ (Cartesian $x,y,z$) it is exactly the $3\times3$ point-group matrix $R(g)$; for higher $l$ it is
obtained as the $l$-th symmetric-power representation of $R(g)$ acting on the monomials $x^ay^bz^c$
($a+b+c=l$), or, for spherical harmonics, via the Cartesian-to-spherical transform $C_l$ (the same
transform that converts Cartesian to spherical AOs elsewhere in the code):

$$
T^{\text{sph},l}(g) = C_l\,T^{\text{cart},l}(g)\,C_l^{+},
$$

with $C_l^+$ the Moore–Penrose pseudo-inverse. The full **AO action matrix** $U_g$ assembles these
shell-local blocks at their AO offsets:

$$
(U_g)_{\,o(gA,l,k)+m',\ o(A,l,k)+m} = T^{(l)}_{m'm}(g).
$$

$U_g$ is a permutation-times-rotation matrix — sparse at the shell level, never built densely in
production code — with two defining properties that every implementation must verify numerically:
$T^{(l)}(g)T^{(l)}(h)\approx T^{(l)}(gh)$ (representation homomorphism) and $T^{(l)}(E)\approx I$
(identity). A scalar one-electron operator such as the overlap or the Fock matrix is symmetric under every
group element by construction of the molecular Hamiltonian, so $U_gSU_g^{\mathsf T}=S$ and
$U_gFU_g^{\mathsf T}=F$ for the *exact* Fock matrix — any deviation observed on a computed $F$ is numerical
noise from screening/direct integral evaluation, not physics, which motivates the Fock symmetrization step
below.

## The symmetry projector and character theory

qc-rs targets **Abelian point groups** in its production symmetry path
($C_1,C_s,C_i,C_2,C_{2v},C_{2h},D_2,D_{2h}$, plus axis-specific variants) — every irrep is
one-dimensional and every character is $\pm1$, which makes the projector construction and its numerical
verification simple and robust. (`qc-mol`'s tables in fact cover 65 point groups including non-Abelian
ones such as $C_{3v}$ or $T_d$, for labeling/detection purposes; the production **Fock-blocking** SCF path
falls back to $C_1$ whenever the requested group is non-Abelian or the atom-mapping/projector checks
fail — see the failure-policy table below.)

For an Abelian group, the projector onto irrep $\Gamma$ is the group average of $U_g$ weighted by the
character $\chi_\Gamma(g)$:

$$
P_\Gamma = \frac{1}{|G|}\sum_{g\in G}\chi_\Gamma(g)\,U_g,
$$

where $|G|$ is the number of symmetry operations. (The general, non-Abelian formula includes an irrep
dimension factor $d_\Gamma$, $P_\Gamma=\tfrac{d_\Gamma}{|G|}\sum_g\chi_\Gamma(g)^*U_g$, which reduces to
the above when $d_\Gamma=1$.) A correctly constructed set of projectors satisfies three group-theoretic
identities that qc-rs's test suite checks directly:

$$
P_\Gamma^2 \approx P_\Gamma \quad\text{(idempotent)}, \qquad
P_\Gamma P_\Lambda \approx 0\ (\Gamma\ne\Lambda) \quad\text{(orthogonal)}, \qquad
\sum_\Gamma P_\Gamma \approx I \quad\text{(complete)}.
$$

Completeness can fail — geometry not exactly on the requested symmetry axis (tolerance is
$10^{-4}$–$10^{-3}$ bohr for atom mapping), or a basis that isn't itself symmetry-complete — in which case
qc-rs's default policy (`symmetry="auto"`) silently falls back to $C_1$; an explicit `symmetry=True` treats
the same failure as an error instead.

## Building the SALC basis

$P_\Gamma$ is a projector onto a subspace of the full AO space, not yet a basis for it. Diagonalizing it,

$$
P_\Gamma = V_\Gamma\, n_\Gamma\, V_\Gamma^{\mathsf T},
$$

and keeping the columns whose eigenvalue exceeds a cutoff $\tau_{\text{salc}}$ (default $10^{-8}$) gives
the **symmetry-adapted linear combination (SALC)** basis for that irrep:

$$
Q_\Gamma = V_\Gamma[:,\, i \mid n_i > \tau_{\text{salc}}].
$$

Stacking every irrep's SALC columns side by side gives the full transform from the AO basis to the
block-labeled SAO basis, $Q = [\,Q_{\Gamma_1}\ Q_{\Gamma_2}\ \cdots\,]$, and any AO-basis operator becomes
block-diagonal in it: $A^{\text{SAO}} = Q^{\mathsf T}A^{\text{AO}}Q$. (An alternative construction — rank-
revealing QR/SVD on projected trial vectors $P_\Gamma e_i$ — trades a slightly more involved implementation
for a more direct rank decision; qc-rs's v1 uses the projector-eigendecomposition approach above for its
simplicity.)

:::{prf:algorithm} SALC-blocked SCF eigenproblem
:label: alg-salc-block

**Input:** AO overlap $S$, AO Fock $F$, point group $G$ with irreps $\{\Gamma\}$, cutoff $\tau_{\text{salc}}$.
**Output:** per-irrep orbital coefficients, energies, and an MO→irrep label for every orbital.

1. For each irrep $\Gamma$, build the projector $P_\Gamma=\tfrac1{|G|}\sum_g\chi_\Gamma(g)U_g$ and diagonalize
   it; keep eigenvectors with eigenvalue $>\tau_{\text{salc}}$ as $Q_\Gamma$.
2. Verify $P_\Gamma^2\approx P_\Gamma$, $P_\Gamma P_\Lambda\approx0$, $\sum_\Gamma P_\Gamma\approx I$; on
   failure, fall back to $C_1$ (`symmetry="auto"`) or raise (`symmetry=True`).
3. For each irrep block: form $S_\Gamma=Q_\Gamma^{\mathsf T}SQ_\Gamma$, diagonalize
   $S_\Gamma=U_\Gamma s_\Gamma U_\Gamma^{\mathsf T}$, and drop eigenvalues below the orthogonalization
   cutoff (default $10^{-8}$) to get the non-redundant, orthonormal transform
   $X_\Gamma = s_{\Gamma,\text{nr}}^{-1/2}U_{\Gamma,\text{nr}}^{\mathsf T}Q_\Gamma^{\mathsf T}$.
4. Form the block Fock $F_\Gamma' = X_\Gamma F X_\Gamma^{\mathsf T}$ and diagonalize it independently:
   $F_\Gamma'C_\Gamma' = C_\Gamma'\varepsilon_\Gamma$.
5. Back-transform each block's MOs to the AO basis, $C_\Gamma^{\text{AO}} = C_\Gamma'^{\mathsf T}X_\Gamma$,
   and label every orbital in that block with irrep $\Gamma$.
6. Merge all blocks' orbitals for the Aufbau/occupation step, ordering by orbital energy across blocks.
:::

Step 3 is exactly the non-redundant overlap orthogonalization from [the linear-algebra
foundations](variational-lcao-basis.md), applied *inside* each symmetry block instead of once over the
full AO space — symmetry blocking and non-redundancy handling compose rather than compete, because each
irrep subspace is overlap-closed ($U_gSU_g^{\mathsf T}=S$ guarantees $S_\Gamma$ doesn't mix with
$S_\Lambda$). In the trivial $C_1$ case, $Q=I$ and this reduces exactly to ordinary non-redundant
orthogonalization — `symmetry=False` and a symmetry-enabled run on a molecule with no symmetry beyond
identity must and do agree bit-for-bit.

```{mermaid}
flowchart LR
    AO["AO basis"] -->|"project by irrep<br/>P_Gamma"| SALC["SALC blocks Q_Gamma"]
    SALC -->|"S_Gamma = Q^T S Q<br/>block-diagonalize"| ORTH["orthonormal blocks X_Gamma"]
    ORTH -->|"F'_Gamma = X F X^T"| DIAG["diagonalize each block<br/>independently"]
    DIAG -->|"back-transform + label"| MO["MOs with irrep labels"]
```

Verified example — water in $C_{2v}$, STO-3G (7 AOs total):

```python
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").set_sym("C2v")
m.symmetry_block_dims()
# {'A1': 4, 'A2': 0, 'B1': 1, 'B2': 2}     -- sums to 7, matches the AO count

m = m.scf(ref="r", symmetry=True).run()
m.scf.energy, m.scf.converged, m.scf.state_irrep
# (-74.9629466565387, True, 'A1')
m.current_mo_irreps[:5]
# ['A1', 'A1', 'B2', 'A1', 'B1']    -- the 5 occupied MOs, in canonical Mulliken labels
```

The block dimensions ($4+0+1+2=7$) sum exactly to the AO count — the projector-completeness identity
$\sum_\Gamma P_\Gamma\approx I$ made concrete — and the totally-symmetric $A_1$ irrep dominates because the
oxygen $2s$/$2p_z$ and both hydrogen $1s$ combinations are symmetric under both the $C_2$ rotation and both
mirror planes.

## Fock symmetrization

Even when the input density is exactly symmetry-adapted, a computed AO Fock matrix can carry small
off-block numerical noise from screened/direct integral evaluation — the exact operator commutes with
every $U_g$, but floating-point evaluation does not do so to machine precision automatically. qc-rs's
default `symmetrize_fock=true` removes this by group-averaging the *scalar* (totally-symmetric) AO
operators after each Fock build:

$$
\mathcal S[A] = \frac{1}{|G|}\sum_{g\in G}U_g A U_g^{\mathsf T}.
$$

For an operator that itself transforms as a non-trivial irrep $\Lambda$ (a dipole component, an external
field perturbation — relevant to linear-response theory, not the ground-state SCF loop), the analogous
character-weighted average is

$$
\mathcal S_\Lambda[A] = \frac{1}{|G|}\sum_{g\in G}\chi_\Lambda(g)\,U_g A U_g^{\mathsf T}.
$$

The ground-state SCF Fock matrix is always totally symmetric ($\Lambda=A_1/A_g$), so `symmetrize_fock`
only ever needs $\mathcal S[\cdot]$; the $\Lambda\ne A_1$ form is future machinery for symmetry-selected
response properties.

## MO irrep assignment, degeneracy, and restart

Because {prf:ref}`alg-salc-block` diagonalizes each irrep block *independently*, every resulting MO
carries an unambiguous irrep label straight from which block produced it — no post-hoc analysis is needed,
unlike the labeling schemes some other codes bolt on after a full-space diagonalization. This label
survives into `current_mo.irreps` in the checkpoint and is what `mychk.current_mo_irreps` and
`mychk.scf.state_irrep` (the many-electron state's own overall irrep, the product of the singly/partially
occupied orbitals' irreps for a closed-shell determinant this is just $A_1$/$A_g$) expose.

The label is not indelible, though — it is metadata *about how the orbitals were produced*, and several
things invalidate it:

| Situation | `irreps` after |
|---|---|
| Fresh SALC-block diagonalization | `Some([...])`, one label per MO |
| A rotation confined within one irrep block | labels preserved |
| Localization mixing MOs across irreps (Boys, Pipek–Mezey, IBO) | `None` — a localized orbital has no single irrep |
| `guess("read", irreps="ignore")` | `None`, dropped on import |
| `guess("read", irreps="preserve")` on a mismatched point group/AO layout | an error, not a silent drop |
| `guess("read", irreps="auto")` | preserved if point group + SALC layout hash match, else dropped |
| `spin_break="mix"`/`"afm"` (intentionally symmetry-broken UHF) | dropped — the whole point is to leave the irrep-pure solution |

Two more subtleties matter once you go beyond simple closed-shell molecules. First, **near-degenerate
orbitals of different irreps** are not actually mixed by the block-diagonal solver — even if two orbitals
from different blocks land at nearly the same orbital energy (common for, e.g., the near-degenerate $\pi$
orbitals of a linear or highly symmetric fragment reduced to a lower-symmetry Abelian subgroup), each stays
confined to its own block and its own label; only the *merge-and-order-by-energy* step (
{prf:ref}`alg-salc-block` step 6) interleaves them for Aufbau occupation. Second, **non-Abelian groups have
genuinely multi-dimensional irreps** ($E$, $T$, …) whose partner functions are not distinguished by a label
alone; qc-rs's `IrrepLabel` metadata already carries an optional partner index for this case, but the
production Fock-blocking path does not yet exploit multi-dimensional projectors — non-Abelian requests
fall back to the $C_1$ path today (see the table below), and a molecule whose true symmetry is, say, $C_{3v}$
or $T_d$ is run with symmetry disabled or reduced to an Abelian subgroup.

## Failure policy

Symmetry detection and adaptation can fail for purely numerical reasons (geometry not exactly on-axis,
floating-point atom-mapping mismatch) as well as structural ones (non-Abelian group requested). qc-rs's
policy distinguishes `symmetry="auto"` (silent, safe fallback) from an explicit `symmetry=True` (fail
loud, because the user asked for symmetry specifically and a silent fallback could hide a geometry bug):

| Failure | `symmetry="auto"` | `symmetry=True` |
|---|---|---|
| Point-group table has no usable character data | fall back to $C_1$ | error |
| Atom mapping under a symop fails | fall back to $C_1$ | error |
| Projector idempotency check fails | fall back to $C_1$ | error |
| SALC rank sum $\ne n_{\text{ao}}$ | fall back to $C_1$ | error |
| Overlap block has a negative eigenvalue | error | error |
| Non-redundancy cutoff empties every block | error | error |

(A negative overlap eigenvalue or a completely emptied block are numerical/basis-linear-dependence
failures unrelated to symmetry per se — they error under both policies, matching the same failures
[the linear-algebra foundations chapter](variational-lcao-basis.md) describes for the unsymmetrized case.)

:::{exercise}
:label: ex-symmetry-theory

1. Water's $C_{2v}$ STO-3G run above gives block dimensions $\{A_1{:}4,\ A_2{:}0,\ B_1{:}1,\ B_2{:}2\}$.
   Water has 7 AOs (O: $1s,2s,2p_x,2p_y,2p_z$; each H: $1s$). Using the fact that $A_2$ transforms as
   $R_z$ (a rotation about the $C_2$ axis, antisymmetric under both mirror planes) and that none of the
   valence AOs listed transform that way at this geometry, explain in one sentence why the $A_2$ block is
   empty.
2. Why does `symmetrize_fock` only ever need the totally-symmetric average $\mathcal S[\cdot]$ for a
   ground-state SCF loop, never the irrep-weighted $\mathcal S_\Lambda[\cdot]$?
3. A user requests `symmetry="C3v"` for ammonia and it silently falls back to $C_1$ under `symmetry="auto"`,
   giving the right energy but no irrep labels. What is the qc-rs-specific (not general chemistry) reason
   this happens, and what would `symmetry=True` have done instead?
:::

:::{solution} ex-symmetry-theory
:class: dropdown

1. None of water's AOs (s and p functions centered on the atoms, in the standard $yz$-molecular-plane
   orientation) transform like a rotation about the $C_2$ axis — s functions are always totally symmetric,
   and each p function transforms like $x$, $y$, or $z$ (i.e. like $B_1$, $B_2$, or $A_1$ in the standard
   $C_{2v}$ character table), none of which is $A_2$. Since the projector $P_{A_2}$ only has something to
   project *onto* if some AO combination already carries $A_2$ character, and none does at this
   angular-momentum level, $P_{A_2}$ has rank 0 on this basis.
2. The ground-state electronic Hamiltonian (and hence the exact SCF Fock operator) is invariant under
   every point-group operation by construction — it has no preferred direction beyond the molecule's own
   symmetry — so it is always the totally symmetric irrep, $A_1$/$A_g$. The irrep-weighted average
   $\mathcal S_\Lambda[\cdot]$ is only needed for operators that are *not* totally symmetric, such as an
   individual dipole component or an external perturbing field, which arise in linear-response/property
   calculations, not in the SCF energy loop itself.
3. qc-rs's production SALC/Fock-blocking machinery targets **Abelian** point groups only ($C_1, C_s, C_i,
   C_2, C_{2v}, C_{2h}, D_2, D_{2h}$); $C_{3v}$ has a genuinely two-dimensional irrep ($E$) that the current
   one-dimensional-character projector construction cannot block correctly, so the failure-policy table's
   "point-group table has no usable character data" / structural-mismatch row applies and `"auto"` falls
   back to $C_1$ — the energy is still exact (a $C_1$ SCF is always correct, just unblocked and unlabeled),
   only the speedup and the labels are lost. `symmetry=True` would have raised an error instead of silently
   dropping to $C_1$, surfacing the non-Abelian-group limitation immediately rather than after the fact.
:::

The projector-and-block machinery here is also the foundation two later chapters build on directly:
[Density-fitting/RI](density-fitting-ri.md) and the post-HF chapters treat symmetry as an orthogonal axis
(RI never assumes symmetry), while a future CASSCF active-space selection and TD-DFT symmetry-forbidden
transition analysis (noted as future work in `.design/24qc.scf-salc.md`) will reuse exactly this
irrep-labeled MO basis.
