Symmetry & group theory

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 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.)

Algorithm 2 (SALC-blocked SCF eigenproblem)

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.

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\).
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).
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}\).
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\).
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\).
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, 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.

        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):

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 Algorithm 2 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 ( Algorithm 2 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 describes for the unsymmetrized case.)

Exercise 5

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.
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]\)?
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 to Exercise 5

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.
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.
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 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.