AICCM2026DEV-B: an independent finite-torus derivation

Status and scope

aiccm2026dev-b is an experimental, independently derived restricted and unrestricted Hartree–Fock and Kohn–Sham method with finite-torus MP2 and coupled-cluster extensions. It coexists with the earlier CCM implementation and is selected with

result = vibeqc.run_periodic_job(
    system,
    basis,
    method="RHF",
    jk_method="aiccm2026dev-b",
    aiccm_lattice_extension=(4, 4, 4),
    aiccm_backend="four_center",  # or "ri" / "rijcosx"
)

The implementation does not use vibeqc.periodic.ccm and does not copy the historical program. The earlier implementation, the historical program, and published CCM formulae are comparison targets. The present definition is the finite Born–von Karman (BvK) torus Hamiltonian restricted to translation-invariant determinants. It is evaluated in the finite translation group’s reciprocal irreducible representations. This is an exact change of basis, not an appeal to reciprocal-space convergence.

The implemented phase is deliberately narrow:

• RHF, RKS, UHF, and UKS (including hybrid functionals), integer occupations, neutral primitive cell;

• 1D, 2D, and 3D periodic cells, with an explicit dimensional Coulomb-gauge caveat below;

• direct corrected-gauge four-centre, pair-resolved RI, and RIJCOSX electron-repulsion backends;

• the existing periodic integral, GDF, Ewald, DIIS, basis, and lattice infrastructure;

• 3D canonical RHF RI-MP2, unrestricted RI-MP2, full-domain UCCSD(T), and restricted/unrestricted local-PNO routes on the exact real torus;

• gauge-invariant SCF one-particle properties; no charged-cell background, analytic gradients, Berry-phase polarization, or relaxed correlated density matrices yet.

The principal source studied was Peintinger and Bredow, J. Comput. Chem. 35, 839 (2014), DOI 10.1002/jcc.23550, together with Chapter 8 of Peintinger’s thesis and the earlier CCM papers cited there. The KS-ADFT CCM of Janetzko, Köster, and Salahub, J. Chem. Phys. 128, 024102 (2008), DOI 10.1063/1.2817582, was audited separately because it introduces two- and three-center multiplicative weights in deMon2k. Those sources motivate the questions; they do not define the implementation below.

1. Finite cyclic translation group

Let the primitive direct lattice be

[ A=(\mathbf a_1,\mathbf a_2,\mathbf a_3), \qquad \mathbf R_{\mathbf n}=A\mathbf n, ]

with only the first (d) directions active. The mathematical matrix (A) has lattice vectors as columns. The repository stores those vectors as rows of PeriodicSystem.lattice, so the implementation uses its transpose. This is immaterial for diagonal cells but essential for a skew lattice.

For positive integers (N_i), with (N_i=1) on inactive directions, define

[ \mathcal T_N = \mathbb Z_{N_1}\times\mathbb Z_{N_2} \times\mathbb Z_{N_3}, \qquad N_c=N_1N_2N_3. ]

The BvK identifications are

[ \mathbf R_{\mathbf n+N_i\mathbf e_i}\equiv\mathbf R_{\mathbf n}. ]

This quotient is the cyclic cluster. It is not a finite open cluster and it has no surface. Its superlattice vectors are (N_i\mathbf a_i). AO labels are ((\mu,\mathbf n)), with (\mu) in the primitive cell and (\mathbf n\in\mathcal T_N).

The primary user parameter is the real-space tuple ((N_1,N_2,N_3)), called the lattice extension. Alternatively, a shell radius (s_i) requests the odd extension

[ N_i=2s_i+1, ]

whose representatives run from (-s_i) through (+s_i). The Wigner–Seitz cell of the superlattice has a half-extent (N_i/2) in primitive lattice coordinates. This is analogous to k sampling because the reciprocal net below is its exact character group. It is not an additional real-space truncation: the selected Coulomb kernel remains BvK-periodic, and its Ewald or fitted long-range terms are not discarded outside a per-atom sphere.

The characters of this finite Abelian group are

[ \chi_{\mathbf m}(\mathbf n) =\exp\left(2\pi i\sum_i \frac{m_i n_i}{N_i}\right), \qquad m_i=0,\ldots,N_i-1, ]

corresponding to the Gamma-centred mesh

[ \mathbf k_{\mathbf m}=\sum_i \frac{m_i}{N_i}\mathbf b_i \pmod{\mathcal L^*}. ]

An even-mesh half shift is a different boundary condition and is therefore not permitted. Symmetry reduction is also not used: all (N_c) characters belong to the exact transform.

For any block-circulant AO matrix,

[ X_{\mu\mathbf0,\nu\mathbf R}=X_{\mu\nu}(\mathbf R), ]

the finite Fourier pair is

[ X_{\mu\nu}(\mathbf k) =\sum_{\mathbf R\in\mathcal T_N} e^{+i\mathbf k\cdot\mathbf R}X_{\mu\nu}(\mathbf R), ]

[ X_{\mu\nu}(\mathbf R) =\frac1{N_c}\sum_{\mathbf k} e^{-i\mathbf k\cdot\mathbf R}X_{\mu\nu}(\mathbf k). ]

Consequently, a Gamma-point calculation on the (N_c)-cell BvK supercell, constrained to commute with primitive translations, and a primitive-cell calculation on the above full mesh are unitarily identical. This equality holds for every finite (N); it is not merely a large-cluster limit.

2. Wigner–Seitz representatives and boundary weights

For centres (A) and (B) with intra-cell displacement (\boldsymbol\delta_{AB}=\boldsymbol\tau_B-\boldsymbol\tau_A), an element ([\mathbf n]\in\mathcal T_N) has infinitely many representatives

[ \boldsymbol\delta_{AB}+A\left(\mathbf n+ \operatorname{diag}(N_1,N_2,N_3)\mathbf z\right), \qquad \mathbf z\in\mathbb Z^d. ]

Let (M_{AB}([\mathbf n])) be the set that minimises this Cartesian distance from centre (A). If (m_{AB,[\mathbf n]}=|M_{AB}([\mathbf n])|), the weight of each tied representative is

[ w_{AB}(\mathbf r\mid[\mathbf n])=\frac1{m_{AB,[\mathbf n]}},\qquad \sum_{\mathbf r\in M_{AB}([\mathbf n])} w_{AB}(\mathbf r\mid[\mathbf n])=1. ]

For coincident basis offsets in a one-dimensional four-cell cluster, residue 2 has representatives (+2) and (-2), each with weight (1/2). In skew cells there may be higher-multiplicity edge or vertex ties. The code solves the closest-vector problem and stops only when a smallest-singular-value bound proves that no unexamined image can tie or improve the minimum. The implementation accepts the explicit centre-pair offset, so a non-Bravais basis is not reduced to the coincident-centre special case.

Where the weights enter

The weights select a representative of a translation equivalence class. They do not multiply an AO merely because its atom lies on a geometric boundary. For a periodised one-electron kernel,

[ h^N_{\mu\nu}([\mathbf R]) =\sum_{\mathbf r\in M([\mathbf R])} w(\mathbf r\mid[\mathbf R])h^N_{\mu\nu}(\mathbf r). ]

All terms on the right are identical under a superlattice translation, so the weights form a partition of unity. Similarly, after fixing the first AO at the origin, a periodised four-centre integral is

[ (\mu\mathbf0,\nu[\mathbf R]\mid \lambda[\mathbf S],\sigma[\mathbf T])N =\sum{\mathbf r,\mathbf s,\mathbf t} w_{\mathbf R}(\mathbf r)w_{\mathbf S}(\mathbf s)w_{\mathbf T}(\mathbf t) (\mu\mathbf0,\nu\mathbf r\mid \lambda\mathbf s,\sigma\mathbf t)_N. ]

Again, this is an average over equivalent representatives of the same periodised integral. It therefore leaves the value unchanged. Applying the same class rule to every symmetry-related term preserves

[ (12\mid34)=(21\mid34)=(12\mid43)=(34\mid12)^*. ]

These identities are necessary for the Coulomb and exchange matrices to be functional derivatives of one scalar RHF energy.

This differs fundamentally from multiplying an ordinary free-space ERI by products of pair-dependent boundary factors. Pair-product factors generally change under ((12\mid34)\leftrightarrow(34\mid12)) and therefore do not, without an additional derivation or symmetrisation, define the usual RHF energy functional.

3. Finite-torus Coulomb Hamiltonian

The Coulomb kernel is periodised on the BvK supercell. In three dimensions it is the zero-average solution

[ -\nabla^2G_N(\mathbf r) =4\pi\left(\sum_{\mathbf L_N}\delta(\mathbf r-\mathbf L_N) -\frac1{\Omega_N}\right), ]

or, equivalently,

[ G_N(\mathbf r)=\frac{4\pi}{\Omega_N} \sum_{\mathbf G_N\ne\mathbf0} \frac{e^{i\mathbf G_N\cdot\mathbf r}}{|\mathbf G_N|^2}. ]

The (\mathbf G=0) coefficient is fixed to zero. The constant background in the Poisson equation makes each separated charge component finite; for a neutral electron-plus-nucleus cell the arbitrary constant cancels in the total energy. The current implementation rejects charged primitive cells rather than silently choosing a jellium convention.

The electron–nucleus, nucleus–nucleus, and Hartree terms must use this same gauge. The implementation obtains all three from vibe-qc’s common periodic Ewald/GDF path. Exchange uses the same periodised interaction. Its finite mesh (q=0) singular channel is treated with the BvK-supercell Madelung (exxdiv="ewald") correction already implemented in the periodic GDF backend. Using the primitive-cell Madelung constant here would define a different finite Hamiltonian and is incorrect.

The finite-group construction itself applies unchanged in one, two, or three periodic dimensions. The currently executable Coulomb kernels differ by dimension. In 3D the four-center route uses the neutral Ewald J split and the Ewald exact-exchange correction above. In 1D and 2D it uses vibe-qc’s direct-truncated lower-dimensional lattice gauge, without a 3D spheropole or exchange-divergence term. The RI routes use their dimension-aware fitted Coulomb metric. Therefore the lower-dimensional four-center and RI routes are each defined and variational for their chosen kernel, but they are not yet a matched absolute-energy backend pair. A common 1D wire or 2D slab Green function must be established before backend differences can be called fitting errors.

4. RHF energy, Fock matrix, and double counting

Use a spin-summed density (D(\mathbf k)). With uniform weight (w_{\mathbf k}=1/N_c), the electron count per primitive cell is

[ N_e=\sum_{\mathbf k}w_{\mathbf k} \operatorname{Tr}[D(\mathbf k)S(\mathbf k)]. ]

The closed-shell projector condition is

[ D(\mathbf k)S(\mathbf k)D(\mathbf k)=2D(\mathbf k). ]

Define the periodised Coulomb and exchange contractions in the conventional way,

[ J_{12}[D]=\sum_{34}D_{34}(12\mid34)N, \qquad K{12}[D]=\sum_{34}D_{34}(13\mid24)_N. ]

Translation conservation makes these block diagonal in (\mathbf k), with the usual momentum-transfer sums implicit in the GDF contraction. The energy per primitive cell is

[ \mathcal E_N[D]=E_{NN}^N +\sum_{\mathbf k}w_{\mathbf k} \operatorname{Tr}[D(\mathbf k)h(\mathbf k)] +\frac12\sum_{\mathbf k}w_{\mathbf k} \operatorname{Tr}\left[D(\mathbf k) \left(J(\mathbf k)-\frac12K(\mathbf k)\right)\right]. ]

Its derivative is

[ F(\mathbf k)=h(\mathbf k)+J(\mathbf k)-\frac12K(\mathbf k). ]

The factors (1/2) and (1/2) are respectively the pair-counting factor and closed-shell exchange factor. No separate “central cell” energy is added, and no Wigner–Seitz factor is applied a second time. Each finite-group translation class occurs exactly once; boundary representative weights sum to one. This is the double-counting argument.

For RKS, let (E_{\mathrm{xc}}[\rho_D]) be the libxc functional and (\alpha) the functional’s exact-exchange fraction. The finite-torus energy is

[ \mathcal E_N^{\mathrm{RKS}}[D]=E_{NN}^N +\operatorname{Tr}_w[Dh] +\frac12\operatorname{Tr}_w[DJ] -\frac{\alpha}{4}\operatorname{Tr}w[DK] +E{\mathrm{xc}}[\rho_D], ]

where (\operatorname{Tr}_w) includes the uniform finite-character weights. Its derivative is

[ F^{\mathrm{RKS}}=h+J-\frac{\alpha}{2}K+V_{\mathrm{xc}}. ]

Thus pure semilocal DFT has (\alpha=0), whereas hybrids use the same periodised exchange tensor and finite-size convention as RHF.

5. Variational statement

For a fixed real-space extension, basis, periodised Coulomb gauge, and nuclear geometry, aiccm2026dev-b minimises the appropriate SCF functional over

[ \mathcal V_N={D:\ D^\dagger=D,\ DSD=2D,
\operatorname{Tr}(DS)=N_cN_e,\ [D,T_{\mathbf R}]=0}. ]

Equivalently, it minimises independently represented occupied subspaces at all finite-group (\mathbf k), subject to the common electron count. The generalised Roothaan equations are

[ F(\mathbf k)C(\mathbf k) =S(\mathbf k)C(\mathbf k)\varepsilon(\mathbf k). ]

DIIS and damping may alter the path to a stationary point but not the functional. A converged solution is stationary in (\mathcal V_N); as with ordinary RHF or approximate KS-DFT, SCF alone does not prove it is the global minimum. A fully unconstrained Gamma-supercell calculation has a larger determinant space and may break primitive translation symmetry. Such a broken-symmetry result is not required to equal the primitive-cell mesh result.

Unrestricted variational space

Let (P^\alpha) and (P^\beta) be spin densities with unit occupation and let (P=P^\alpha+P^\beta). The finite-torus UHF functional per primitive cell is

[ \mathcal E_N^{\mathrm{UHF}}[P^\alpha,P^\beta] =E_{NN}^N+\sum_{\sigma\in{\alpha,\beta}} \operatorname{Tr}_w[P^\sigma h] +\frac12\operatorname{Tr}w[PJ[P]] -\frac12\sum\sigma\operatorname{Tr}_w[P^\sigma K[P^\sigma]]. ]

Its independent functional derivatives are

[ F^\sigma=h+J[P]-K[P^\sigma]. ]

The minimized set is

[ \mathcal V_N^{\mathrm U}={(P^\alpha,P^\beta): (P^\sigma)^\dagger=P^\sigma,\quad P^\sigma S P^\sigma=P^\sigma,\quad [P^\sigma,T_{\mathbf R}]=0,\quad \operatorname{Tr}(P^\sigma S)=N_c n_\sigma}. ]

Here (n_\alpha=(N_e+2S)/2) and (n_\beta=(N_e-2S)/2) are integers fixed by the primitive-cell multiplicity. This is a genuine unrestricted variation: alpha and beta orbitals are not constrained to share a spatial subspace. UHF stationarity does not prove spin purity. The implementation therefore reports both (\langle S^2\rangle) and the ideal (S(S+1)).

For UKS with a hybrid fraction (a_x), replace the UHF exchange energy by

[ -\frac{a_x}{2}\sum_\sigma \operatorname{Tr}w[P^\sigma K[P^\sigma]] +E{\mathrm{xc}}[\rho_\alpha,\rho_\beta], ]

which gives

[ F^\sigma_{\mathrm{UKS}}=h+J[P]-a_xK[P^\sigma] +V_{\mathrm{xc}}^\sigma. ]

No new Wigner–Seitz weight appears in these equations. Both spin projectors use the same finite Hamiltonian and the same translation-class partition.

6. Algorithm

1. Validate neutral, integer-occupation, zero-temperature RHF/RKS/UHF/UKS input in one, two, or three periodic dimensions.

2. Form (\mathcal T_N), its exact Wigner–Seitz representative partition, and the unreduced Gamma-centred character mesh.

3. Build lattice-summed overlap and one-electron matrices with the common periodic gauge.

4. Build the periodised J and K matrices with one of the three backends defined below. Apply the BvK-supercell exchange-divergence correction only in the 3D Ewald gauge.

5. Solve the per-k generalised eigenproblems using the existing canonical linear-dependence handling and SCF accelerator.

6. Return energy per primitive cell and measure (\max_k\lVert D_kS_kD_k-2D_k\rVert_F), electron-count error, Wigner–Seitz partition error, and the imaginary residue of the inverse Bloch density.

The inverse Bloch transform stored in the implementation is also the bridge for future explicitly real-space post-HF work.

Integral backends

All backends act on the same finite-character density and differ only in their representation of the electron-repulsion tensor.

1. four_center uses the periodic BIPOLE engine: in 3D, direct four-centre short-range J/K plus the reciprocal long-range Hartree term and BvK exchange correction; in 1D/2D, its direct-truncated dimensional gauge.

2. ri uses pair-resolved RSGDF or MDF. With auxiliary functions (P,Q) and Coulomb metric (V_{PQ}=(P|Q)),

   [ (\mu\nu|\lambda\sigma){\mathrm{RI}} =\sum{PQ}(\mu\nu|P)(V^{-1})_{PQ}(Q|\lambda\sigma). ]

   The fit is applied consistently to J and K.

3. rijcosx uses the same pair-resolved RI J and evaluates exact exchange by periodic chain-of-spheres quadrature. It is therefore a controlled three-centre/seminumerical approximation, not a different cyclic Hamiltonian.

The q-only compcell fit is rejected: on tight cells it does not represent one consistent four-centre finite-torus tensor and is known to admit non-physical converged fixed points. The present native bridges also reject one-cell RI/RIJCOSX RKS and one-cell RIJCOSX RHF rather than silently substituting a different Gamma-only algorithm.

7. Canonical finite-torus RI-MP2

MP2 is built on the converged B-stream RI-RHF determinant, not on the legacy Gamma-supercell GDF driver. That distinction is numerical as well as formal: on the two-cell 8 x 12 x 12 bohr H2 control, the legacy Gamma-supercell HF reference differs from B RI-HF by 0.0046408868 Ha per cell. Reusing it would change the zeroth-order Hamiltonian.

Let (i,j) denote occupied bands, (a,b) virtual bands, and let (N_k=N_c). The pair-resolved canonical auxiliary factor is

[ \widetilde L^{\mathbf k_i\mathbf k_a}{P\mu\nu} =\left[U{\mathbf q}\lambda_{\mathbf q,+}^{-1/2} U_{\mathbf q}^{\dagger}T_{\mathbf q}\right]_{P\mu\nu}, \qquad \mathbf q=\mathbf k_a-\mathbf k_i. ]

Here only metric eigenvalues above the declared threshold enter. Keeping the factor in the original auxiliary-AO basis is important: the compact factor (\lambda^{-1/2}U^{\dagger}T) is invariant in an SCF contraction with its own adjoint, but independently diagonalised (\mathbf q) and (-\mathbf q) blocks have arbitrary eigenvector phases. Those phases need not cancel in an MP2 contraction between two blocks. The canonical matrix function above is unique on the retained subspace and obeys the required conjugation relation.

Transform one occupied and one virtual index,

[ L^{\mathbf k_i\mathbf k_a}{Pia} =\sum{\mu\nu}C_{\mu i}^{*}(\mathbf k_i) \widetilde L^{\mathbf k_i\mathbf k_a}{P\mu\nu} C{\nu a}(\mathbf k_a). ]

For every triple ((\mathbf k_i,\mathbf k_a,\mathbf k_j)), the fourth character is fixed exactly by

[ \mathbf k_b=\mathbf k_i-\mathbf k_a+\mathbf k_j+\mathbf G, ]

where (mathbf G) returns the point to the finite reciprocal group. The RI two-electron block in the normalization used here is

[ V_{ij}^{ab}(\mathbf k_i,\mathbf k_j,\mathbf k_a) =\frac1{N_k}\sum_P L^{\mathbf k_i\mathbf k_a}{Pia} L^{\mathbf k_j\mathbf k_b}{Pjb}. ]

With

[ \Delta_{ij}^{ab}=\varepsilon_{i\mathbf k_i} +\varepsilon_{j\mathbf k_j}-\varepsilon_{a\mathbf k_a} -\varepsilon_{b\mathbf k_b}, ]

the closed-shell spin components per primitive cell are

[ E_{\mathrm{OS}}^{(2)}=\frac1{N_k} \sum_{\mathbf k_i\mathbf k_j\mathbf k_a}\sum_{ijab} \frac{|V_{ij}^{ab}(\mathbf k_a)|^2}{\Delta_{ij}^{ab}}, ]

[ E_{\mathrm{SS}}^{(2)}=\frac1{N_k} \sum_{\mathbf k_i\mathbf k_j\mathbf k_a}\sum_{ijab} \frac{|V_{ij}^{ab}(\mathbf k_a)|^2 -V_{ij}^{ab}(\mathbf k_a)^{*}V_{ij}^{ba}(\mathbf k_b)} {\Delta_{ij}^{ab}}, ]

and (E_{\mathrm{MP2}}^{(2)}=E_{\mathrm{OS}}^{(2)}+ E_{\mathrm{SS}}^{(2)}). MP2 is not variational: it is the second-order Rayleigh–Schrodinger correction for the finite RI Hamiltonian around the stationary RHF determinant. The finite-group Fourier transform nevertheless makes it exactly equivalent to canonical RI-MP2 in the corresponding BvK supercell when the same Coulomb gauge and fit are used.

The implemented 3D H2/STO-3G two-point mesh agrees with an out-of-process PySCF KRHF/KMP2 calculation to (5.9\times10^{-10}) Ha in HF and (2.4\times10^{-10}) Ha in correlation energy per cell. The machine-readable record is docs/manuscripts/aiccm_comparison/data/h2_mp2_2026-06-21.json. MP2 currently fails closed in 1D and 2D until O1 supplies a common lower-dimensional Coulomb gauge.

Unrestricted MP2

For an unrestricted reference, occupied and virtual labels carry an explicit spin. With antisymmetrized same-spin integrals and ordinary opposite-spin integrals, the second-order correction is

[ E_{\mathrm{UMP2}}^{(2)}=E_{\alpha\alpha}^{(2)} +E_{\beta\beta}^{(2)}+E_{\alpha\beta}^{(2)}, ]

[ E_{\sigma\sigma}^{(2)}= \frac14\sum_{ij\in\sigma}\sum_{ab\in\sigma} \frac{|\langle i_\sigma j_\sigma\Vert a_\sigma b_\sigma\rangle|^2} {\varepsilon_i^\sigma+\varepsilon_j^\sigma -\varepsilon_a^\sigma-\varepsilon_b^\sigma}, ]

[ E_{\alpha\beta}^{(2)}= \sum_{i a\in\alpha}\sum_{j b\in\beta} \frac{|(i_\alpha a_\alpha|j_\beta b_\beta)|^2} {\varepsilon_i^\alpha+\varepsilon_j^\beta -\varepsilon_a^\alpha-\varepsilon_b^\beta}. ]

Every sum also contains the finite characters, with the same momentum conservation and one final division by (N_c) as the restricted expression. The implementation evaluates this equation after the exact inverse transform of the B RI Hamiltonian. Zero PNO threshold and canonical occupieds give the full UMP2 result. Independent unitary localization of the alpha and beta occupied projectors leaves both spin densities invariant; the coupled local MP2 equations are then required because each localized occupied Fock block is not diagonal.

8. Real-torus local correlation and coupled cluster

Exact inverse transform of the fitted Hamiltonian

Let the normalized finite-group Bloch AO be

[ |\mu\mathbf{k}\rangle =N_c^{-1/2}\sum_{\mathbf R\in\mathcal T_N} e^{+i\mathbf k\cdot\mathbf R}|\mu\mathbf R\rangle . ]

For any translation-invariant one-particle matrix, the full finite-torus AO matrix is therefore

[ X_{\mu\mathbf R,\nu\mathbf S} =\frac{1}{N_c}\sum_{\mathbf k} e^{+i\mathbf k\cdot\mathbf R}X_{\mu\nu}(\mathbf k) e^{-i\mathbf k\cdot\mathbf S}. ]

This formula is used for both overlap and Fock. It is a unitary change of representation, not a new Gamma-point SCF.

The RI factors need one more character. Let (\widetilde L^{\mathbf k_i\mathbf k_a}_{P,\mu\nu}) be the canonical auxiliary-AO factor of Section 7 and (\mathbf q=\mathbf k_a-\mathbf k_i). Define

[ B_{P\mathbf T,\mu\mathbf R,\nu\mathbf S} =\frac{1}{N_c^2}\sum_{\mathbf k_i\mathbf k_a} e^{+i\mathbf k_i\cdot\mathbf R} e^{-i\mathbf k_a\cdot\mathbf S} e^{+i\mathbf q\cdot\mathbf T} \widetilde L^{\mathbf k_i\mathbf k_a}_{P,\mu\nu}. ]

The corresponding real-torus MO coefficient is

[ C_{\mu\mathbf R,p\mathbf k} =N_c^{-1/2}e^{+i\mathbf k\cdot\mathbf R}C_{\mu p}(\mathbf k). ]

Contracting the previous two equations gives

[ B_{P\mathbf T,i\mathbf k_i,a\mathbf k_a} =N_c^{-1}e^{+i\mathbf q\cdot\mathbf T} L^{\mathbf k_i\mathbf k_a}_{Pia}. ]

Thus a product of factors with transfers (\mathbf q) and (-\mathbf q) satisfies

[ \sum_{P\mathbf T}B_{P\mathbf T,ia}B_{P\mathbf T,jb} =\frac{1}{N_c}\sum_P L_{Pia}L_{Pjb}, ]

while any nonconserving transfer vanishes by character orthogonality. This is exactly the normalization in the momentum-space MP2 integral. The real and character representations therefore describe one fitted finite Hamiltonian.

Local MP2 and CCSD(T)

The present local implementation deliberately builds the whole finite torus. It includes every occupied pair and every ordered occupied triple required by the spin-integrated equations, obtains a total correlation energy (E_{\mathrm{corr},N}), and reports

[ E_{\mathrm{corr}}^{\mathrm{cell}} =\frac{E_{\mathrm{corr},N}}{N_c}. ]

There is no additional boundary or pair multiplicity. A future translation- reduced implementation may retain one orbit under simultaneous translation of all occupied indices, but then its multiplicity must follow the orbit-stabilizer theorem. Selecting visually unique pairs and multiplying by an occurrence factor would recreate the historical ambiguity.

Occupied orbitals are localized with the Pipek–Mezey objective built from Mulliken populations. This is invariant to the choice of position origin and does not use the ill-defined ordinary position operator on a torus. Molecular Foster–Boys localization is rejected. The occupied rotation is unitary, so with complete pair domains, no PNO truncation, no pair screening, and complete occupied coupling, local MP2 and local CCSD reproduce their canonical finite-torus limits. With complete triple virtual spaces the same statement holds for the perturbative triples correction.

The current defaults allow PNO truncation but keep all occupied pairs and all occupied couplings. Distance-based pair screening and local auxiliary fitting are disabled because their existing molecular implementations use Euclidean distances across the cell boundary. They require minimum-image Wannier and auxiliary domains before they can be enabled without breaking translation symmetry.

Neither MP2 nor CCSD(T) is variational. CCSD is stationary only through its left-right coupled-cluster Lagrangian, and the perturbative triples energy has no upper-bound property. PNO/domain truncations add controlled numerical approximations, not a new Hamiltonian or a variational subspace theorem.

On the H2/STO-3G two-cell control, no-truncation real-torus local MP2 equals the character-space result within (7\times10^{-18}) hartree per cell. On the LiH/STO-3G two-cell control, no-truncation local and canonical finite-torus DF-CCSD correlation energies differ by (6.9\times10^{-12}) hartree total; their nonzero triples corrections differ by (1.4\times10^{-12}) hartree. These are finite-Hamiltonian exact-limit tests, not external infinite-crystal CCSD(T) benchmarks.

For UCCSD, use spin orbitals and

[ |\Psi_{\mathrm{CCSD}}\rangle=e^{T_1+T_2}|\Phi_0\rangle, \qquad \langle\Phi_\mu|e^{-T}He^T|\Phi_0\rangle=0 ]

for every single and double excitation (\mu). The energy is

[ E_{\mathrm{CCSD}}=\langle\Phi_0|e^{-T}He^T|\Phi_0\rangle. ]

The unrestricted PNO pilot expands each pair amplitude into the full finite-torus spin-orbital virtual space, evaluates the complete projected residual, and projects back. With zero PNO threshold this projection is a unitary coordinate change and recovers full-domain UCCSD. Perturbative triples are evaluated from the converged amplitudes. This establishes the equations and exact PNO limit, but the current implementation remains an explicitly cost-capped O(N^6) oracle: it does not yet skip translation-equivalent pairs or triples and therefore does not yet establish reduced scaling.

Post-SCF occupied localization

Let the complete finite character group contain (N_c) points. For each isolated occupied band, the canonical back transform is

[ |w_{n\mathbf R}^{0}\rangle =\frac{1}{\sqrt{N_c}}\sum_{\mathbf k} e^{-i\mathbf k\cdot\mathbf R}|\psi_{n\mathbf k}\rangle . ]

All (N_c n_{\mathrm{occ}}) functions are rotated together, (C^{\mathrm{loc}}=C^{0}U), with (U^\dagger U=I). Consequently

[ C^{\mathrm{loc}}C^{\mathrm{loc}\dagger} =C^0UU^\dagger C^{0\dagger}=C^0C^{0\dagger}. ]

This proves invariance of the finite-torus density and every RHF or RKS energy term. Localizing independent k blocks would not produce the required real-space family.

Exact Marzari–Vanderbilt localization requires cross-k matrix elements of the exponential position operator. Those integrals are not exposed by the current native API. The implemented Wannier route jointly diagonalizes Lowdin-projected AO-centre operators (\cos(2\pi r_\alpha/L_\alpha)) and (\sin(2\pi r_\alpha/L_\alpha)). If

[ z_{n\alpha}=\langle w_n|e^{2\pi i r_\alpha/L_\alpha}|w_n\rangle, \qquad \Omega_{n\alpha} =-\left(\frac{L_\alpha}{2\pi}\right)^2\log|z_{n\alpha}|^2, ]

the latter is the reported circular spread. It has the correct torus topology but is labelled a projected circular AO-centre approximation. Exact continuum spreads and Wannier90 parity remain open.

The IAO alternative constructs polarized intrinsic atomic orbitals from a minimal reference and maximizes squared atomic populations. Its present cross overlap is an explicit-supercell integral rather than a periodized torus integral. A four-cell H2 run gives translation covariance residuals of (1.17\times10^{-13}) for the Wannier path and (7.08\times10^{-6}) for IAO. Density invariance remains at machine precision. IAO translation maps are therefore diagnostic and are not used to skip correlated pairs.

Wraparound is flagged when the antipodal-shell probability exceeds 0.05 or when (\sqrt{\Omega_n}) exceeds one quarter of the shortest cyclic-cluster vector. A nonpositive indirect gap fails closed. Metallic-band disentanglement is not implemented.

PAOs, PNOs, and pair orbits

For an (S)-orthonormal occupied coefficient matrix (C_o), the coefficient projector out of the occupied space is

[ Q=I-C_oC_o^\dagger S,\qquad C_o^\dagger S Q=0. ]

If (E_D) selects AOs in a pair domain, let (\widetilde V_D=QE_D) and diagonalize its Gram matrix:

[ G_D=\widetilde V_D^\dagger S\widetilde V_D=XgX^\dagger, \qquad V_D=\widetilde V_DX_+g_+^{-1/2}. ]

Then (V_D^\dagger S V_D=I) and (C_o^\dagger S V_D=0). This removes both occupied leakage and redundant PAOs.

For pair amplitude matrix (T_{ij}), the model pair density is

[ D_{ij}=\frac{T_{ij}T_{ij}^\dagger+T_{ij}^\dagger T_{ij}} {1+\delta_{ij}}. ]

It is Hermitian positive semidefinite because (x^\dagger D_{ij}x=(|T_{ij}^\dagger x|^2+ |T_{ij}x|^2)/(1+\delta_{ij})\ge0). Its eigenvectors are PNOs. A zero threshold retains the complete PAO span, and positive thresholds give nested ranks. This proves rank nesting, not monotonic MP2 energy error: MP2 is nonvariational and convergence of the energy must be measured.

Translations act simultaneously on both indices of an unordered occupied pair. For measured localization permutations (p_g),

[ \mathcal O(i,j)={\operatorname{sort}(p_g(i),p_g(j)):g\in G_N}. ]

The per-cell coefficient is (|\mathcal O|/N_c), including smaller orbits with nontrivial stabilizers such as half-cell pairs on an even ring. The code verifies a disjoint partition of all unordered pairs, but it still evaluates all pairs. The orbit reduction factor is not a timing claim.

At complete PAO domains and zero PNO thresholds, the B MP2 wrapper also evaluates a rotation-invariant full-space DF-MP2 contraction. The raw iterative local-pair energy and the correction are both reported. Before the real local solver, a full-rank real time-reversal gauge is constructed and its occupied projector is checked against the complex localized projector. This avoids silently discarding imaginary coefficients. No correction is applied to truncated calculations.

Finite-cluster space-group symmetry

Write primitive lattice vectors as columns of (A), and let (N=\operatorname{diag}(N_1,N_2,N_3)). A crystallographic operation (g={W|\mathbf w}) preserves the finite torus exactly when

[ N^{-1}WN\in\mathbb Z^{3\times3}. ]

For primitive atom (a), retain both the atom permutation and integer image shift,

[ W\mathbf f_a+\mathbf w=\mathbf f_{p_g(a)}+\mathbf q_{g,a}, ]

so ((a,\mathbf r)) maps to

[ \left(p_g(a),W\mathbf r+\mathbf q_{g,a}\pmod{\mathbf N}\right). ]

The image shift is essential for screw and glide operations. Reciprocal characters transform as (\mathbf k’=W^{-T}\mathbf k), and a k-orbit has weight (|\mathcal O_k|/N_c). spglib supplies the crystallographic group; the B module builds its exact cluster-compatible subgroup.

Symmetry is diagnostic only. At Gamma, AO matrices can be checked with the Reynolds projector

[ \mathcal P_G[M]=\frac{1}{|G_N|}\sum_g U_g M U_g^T, ]

whose idempotency follows from group closure. At general k, nonsymmorphic sewing phases depend on each backend’s Bloch convention. Therefore symmetry_mode="integrals" fails closed until representative shell-pair and quartet scattering passes full-build Fock and energy parity. Primitive shell pairs and intrinsically eightfold-unique shell quartets are already partitioned into exact point-group orbits and reported as reduction counts; no libint call is skipped from those counts yet.

9. Limits and convergence

Exact finite-size identity

For the same finite-torus Hamiltonian and translation-invariant variational space,

[ E_N^{\Gamma\text{-supercell}}/N_c =E_N^{\text{Gamma-centred mesh}} ]

up to numerical error. This is the strongest internal reference and should hold at every mesh size.

Infinite-crystal limit

As all active (N_i\to\infty), the character sum approaches the Brillouin zone integral and the BvK Green’s function approaches the chosen infinite periodic Coulomb convention. Insulators with localised density matrices are expected to converge rapidly in short-range pieces; exact exchange and the Coulomb finite-size terms can converge algebraically. The sequence need not be monotone because each (N) defines a different finite Hamiltonian.

Molecular limit

For a one-cell cluster with increasing vacuum, image interactions and the BvK Madelung correction vanish, and the result approaches molecular RHF or RKS in the same basis. Basis linear dependence is handled per k by canonical orthogonalisation; dropping different ranks at different k is diagnosed and can spoil smooth convergence.

BSSE is not repaired by the cyclic construction. Counterpoise or basis extrapolation is a separate modelling choice.

10. Assumptions and explicit divergences from historical CCM

Audit of the 2008 deMon2k KS-ADFT CCM

Janetzko, Köster, and Salahub define the two-center occurrence weight

[ \omega_{MN}=\frac{1}{n_{MN}},\qquad T^{\mathrm{CCM}}=\sum_{\mu\nu}P^{\mathrm{CCM}}{\mu\nu} \omega{MN}\sum_{\nu’}T_{\mu\nu’}. ]

For a three-center quantity they do not merely multiply one selected integral. Their Eq. (20) averages both nested enumerations (M’\to N’) and (N’\to M’), divided by

[ n_{MNC}=n_{MN}(n_{MC}+n_{NC}),\qquad \omega_{MNC}=n_{MNC}^{-1}. ]

That explicit two-order average matters. It means the simple (M\leftrightarrow N) objection to the 2014 four-center factor cannot be transferred unchanged to the deMon2k three-center tensor. If both enumerations are complete, Eq. (20) is symmetric in the two AO centers by construction.

The two-center formula still requires the stronger Hermiticity identity

[ \omega_{MN}\sum_{\nu’}I_{\mu\nu’} =\omega_{NM}\sum_{\mu’}I_{\nu\mu’}. ]

The free-space identity (I_{\mu\nu}=I_{\nu\mu}), cited in the paper, does not by itself prove this equality for two different center-specific Wigner– Seitz enumerations. Likewise, variational auxiliary fitting requires the weighted three-center map and weighted auxiliary Coulomb metric to be adjoints under one finite inner product. The paper supplies a fitted-energy expression and its derivative, so the intended scalar functional is clearer than in the 2014 four-center derivation, but it does not prove this quotient-space adjoint condition for arbitrary non-Bravais or skew clusters.

There is also a separate finite-size issue. The paper explicitly states that its CCM and the corresponding periodic calculation are not fully equivalent because interactions outside the Wigner–Seitz regions are omitted. Its auxiliary-density lattice sum is truncated after one neighboring shell, and the XC quadrature is evaluated on a reference-cluster grid with translated auxiliary functions. Those can converge with cluster size, but they are not an exact finite-(N) BvK identity. Thus the deMon2k multiplicative weights remain a symmetry and variational audit target; the present derivation does not label them broken where Eq. (20) already enforces the relevant AO-pair symmetry.

1. Finite torus, not a weighted open cluster. The method is defined by a BvK Hamiltonian first. Wigner–Seitz geometry is a representation of its translation classes.

2. Boundary weights are a partition of unity. They average only tied representatives of the same periodised object. They are not attached to atoms, shells, or arbitrary AO pairs.

3. No historical product-weighted ERI. The published pair-product four-centre factors are not invariant under all ERI permutations in their displayed form. Without a separate scalar functional derivation, using them directly would make the variational and double-counting status unclear. aiccm2026dev-b therefore does not use them.

4. No Fock symmetrisation as a repair. Hermiticity follows from the Hamiltonian. A non-Hermitian raw Fock is treated as a bug.

5. Long range is part of the Hamiltonian. It is not an after-the-fact Madelung energy added to a short-range SCF.

6. Charged cells fail closed. No implicit uniform-background energy is reported.

7. Exchange is not independently truncated. It uses the finite-torus periodised kernel and BvK-size divergence correction.

8. The reciprocal implementation is exact, not a reference shortcut. It is the diagonal representation of the finite cyclic translation group.

9. Source agreement is secondary. Historical numerical agreement does not override idempotency, permutation symmetry, variational consistency, gauge consistency, or finite-supercell/mesh equivalence.

11. Validation ladder

The committed tests cover the selector, exact mesh, Wigner–Seitz ties, partition unity, inverse Fourier transform, metric idempotency, fail-closed domain, RHF/RKS runner dispatch, and native 1D/2D/3D SCF runs through all three integral backends. The benchmark harness adds:

1. 1D alternating hydrogen chain, increasing cyclic length;

2. direct comparison with both the historical union12 and separately developed symmetric aiccm2026dev-a weights in the in-repo CCM on identical geometry;

3. exact comparison with the native Gamma-centred k-mesh representation;

4. energy, convergence, idempotency, electron count, wall time, and implementation-to-implementation differences;

5. later 3D LiH, MgO, NaCl, and diamond ladders after each basis/gauge is pinned to an out-of-process reference input.

External CRYSTAL or PySCF comparisons remain out of process and are weak signals unless their cell, Coulomb, exchange-divergence, basis, and k-mesh conventions are matched explicitly.