Coupled cluster
"RI-CC" here means the two-electron integrals in the amplitude equations come from density fitting (three-centre B tensors from an RI auxiliary basis), not from exact four-centre integrals. Everything below needs an orbital basis and an RI auxiliary basis. Formal cost is O(N⁶) for CCSD and O(N⁷) for (T).
CCSD and CCSD(T)
What it is. RI-CCSD and the perturbative triples correction (T), in two implementations:
- Spin-adapted, closed shell (
ccsd_closed_shell,ccsd_t_closed_shell): amplitudes over spatial orbitals, the algorithm PySCF'scc.CCSDuses. This is what both the CLI and Python run for an RHF reference. Its VVVV block is 16× smaller than the spin-orbital one; measured about 8–10× faster than the spin-orbital CCSD at cc-pVDZ, and the (T) step 9.6–42× faster. - Spin-orbital (
ccsd,ccsd_t): built from the same RHF spatial orbitals and kept as the cross-check of the spin-adapted solvers, so it is closed shell only as well. Its (T) streams one occupied triple at a time, so memory is O(no·nv³)-class rather than the dense six-index tensor.
No open-shell CCD, CCSD or CCSD(T) exists in ferric: every CC solver reads
restricted orbitals and refuses an unrestricted reference. The only open-shell
coupled-cluster-type method is LinLCCD(hh): the CLI linlccd kind is
closed shell, and its open-shell (UHF) version is library-only
(ferric_cc::linlccd_u::u_linlccd).
Run it.
- CCSD:
method.kind = "ccsd"(examples/water-ccsd.toml, water/cc-pVDZ); Pythonferric.run_ccsd(mol, bs, aux). - CCSD(T):
method.kind = "ccsd(t)"(examples/water-ccsd-t.toml); Pythonferric.run_ccsd_t(mol, bs, aux). - CCD:
method.kind = "ccd"(examples/water-ccd.toml); Pythonferric.run_ccd(mol, bs, aux). - All three kinds are closed shell only and refuse
multiplicity > 1.
Aux basis. The RI error is not negligible at CC accuracy: on water /
cc-pVDZ with cc-pvdz-ri, RI-CCSD differs from exact-integral CCSD by about
1.3e-4 Ha. Choose the aux basis for the accuracy you need.
Accuracy. ccsd is Proven.
| Quantity | System / basis | Reference | Agreement | Pinned by |
|---|---|---|---|---|
| CCSD correlation energy | H2 / STO-3G | exact-integral numpy | −0.02052453 Ha | ferric-cc/src/ccsd.rs::test_ccsd_h2_sto3g |
| Closed-shell (T) | H2O / cc-pVDZ | PySCF ccsd_t() | ~1e-6 Ha (test asserts 1e-4) | ferric-cc/src/ccsd_t_closed_shell.rs::closed_shell_t_h2o_ccpvdz_matches_pyscf |
| Streaming vs dense spin-orbital (T) | H2O / cc-pVDZ | ferric's former dense path | 5e-16 Ha | ferric-cc/src/ccsd_t.rs::streaming_matches_dense_h2o_ccpvdz |
Limits. Closed-shell entry points only in the CLI and Python. The spin-adapted (T) rejects spin-orbital amplitudes with a typed error rather than mixing conventions.
Particle-particle ladder: dense or RI
The particle–particle ladder Σcd Wabcd τijcd
of closed-shell CCSD can be evaluated two ways, chosen by [cc] ladder
(Python run_ccsd(ladder=) / run_ccsd_t(ladder=)):
denseholds the chemist (ab|cd) block in core and builds Wabcd every iteration: nv⁴ memory, a few times over.riholds neither the (ab|cd) block nor Wabcd. It rebuilds the (ac|bd) entries from the RI 3-index tensors every iteration in ablk·nv³ panels: at most one nv⁴-sized panel at full width (dense holds four nv⁴ buffers at its peak), and narrower panels when the budget requires. It is the same algebra over the same tensors, so it is exact: on water/cc-pVDZ |ΔE| and max|Δt| are ≲ 3e-17 (round-off), and the alkane energies below agree to all 10 printed digits.auto(the default) runsdensewhenever its memory plan fits the budget andriotherwise, so it is never refused whereriwould run. When it falls back toriit says so on stderr. The route taken is printed (ladder = dense (auto)), logged, and returned (CcResult::ladder_route; PythonCcResult.ladder_route).
Reproducibility. Under auto, a closed-shell CCSD job whose dense plan
exceeds the memory budget runs on the RI ladder instead of being refused.
Because auto decides against the budget, a run without an explicit
[memory] budget_gb (auto-detected as 0.8 × the RAM available at start-up)
can take different routes on the same input from one run to the next. The
energies agree to round-off either way; the time does not. Set ladder
explicitly, or fix the budget, for reproducible timings.
Measured 2026-10-26, alkane/cc-pVDZ (cc-pvdz-ri), OPENBLAS_NUM_THREADS=1,
RAYON_NUM_THREADS=4, routes interleaved one per process, load average 3–6
from other jobs (benchmarks/harness/examples/ri_ccsd_memory.rs). RI ran at
full block width (ablk = nv, the default-budget
choice). Per-iteration time is (time for k iterations − setup time) / k,
repeated 4× per route for water, 3× for C3H8 and 2× for C4H10/C5H12;
ratios are of the medians. Peak RSS is one converged run per route.
| System | nv | Per iteration, dense | Per iteration, RI | RI / dense | Peak RSS, dense | Peak RSS, RI |
|---|---|---|---|---|---|---|
| water | 19 | 7.8 ms | 8.3 ms | 1.06 | — | — |
| C3H8 | 69 | 1.27 s | 1.64 s | 1.30 | 1.50 GB | 0.70 GB |
| C4H10 | 89 | 3.2–3.7 s | 5.2–5.3 s | 1.50 | 3.79 GB | 1.67 GB |
| C5H12 | 109 | 8.8–10.8 s | 14.5–14.6 s | 1.48 | 7.90 GB | 3.40 GB |
RI is 1.3–1.5× slower per iteration from nv ≈ 70 up and uses
2.1–2.3× less memory, which is why auto takes it only when dense does
not fit. The RI route still holds the no·nv³
(ia|bc) block in core, so it lowers the nv⁴ ceiling but does
not remove every large tensor.
LinLCCD and ωB97X-L-V
LinLCCD(hh) is linearized coupled-cluster doubles with the hole–hole ladder
kept to all orders, closed shell only. The ladder keeps the correlation energy
finite as the HOMO–LUMO gap closes, where MP2 diverges.
method.kind = "linlccd" (examples/water-linlccd.toml). Proven (narrow):
no other quantum chemistry code implements LinLCCD(hh), so its energy is
checked against an independent numpy solve of the same equations on PySCF
density-fitted integrals (H2O and NH3 through the closed-shell path, and UHF
OH through the library-only open-shell path; agreement ≤1.2e-12 Ha). With
the ladder off it reduces exactly to RI-MP2, and with exact integrals its
driver terms reproduce canonical MP2.
[mp2] linlccd_variant (Python run_linlccd(variant=...)) selects the
method's ladder terms: hh (default), drivers-only (no ladder, equal to
RI-MP2) or full (hole–hole plus particle–particle, with CCD-like VVVV
memory). Every variant is computed exactly by default. All three variants
match an independent numpy solve on the same density-fitted integrals
(H2O and NH3, 6-31G and cc-pVDZ, ≤ 4.4e-13 Ha); drivers-only also matches
PySCF DFMP2.
Local LinLCCD. With [local] scheme = "amplitude-threshold" and eps
(examples/water-linlccd-local.toml; Python
run_linlccd(..., local="amplitude-threshold", eps=1e-4)) the same variant
is solved in the Boys-localized basis with pair amplitudes at or below eps
dropped (see Exact and local correlation).
eps has no default and is printed and logged with the kept fraction;
eps = 0 reproduces the exact LinLCCD of the same variant and the numpy
reference (≤ 4.4e-13 Ha, every variant), which is how the local path is
anchored. [local] reference = true also runs the exact
LinLCCD and prints the local error. Closed shell and energy only.
ωB97X-L-V is a double-hybrid functional that uses short-range LinLCCD(hh)
instead of MP2 for its correlation term. It converges its own ωB97X-L
Kohn–Sham reference (a non-converged reference is an error), then adds the
LinLCCD(hh) correction on those orbitals. method.kind = "wb97x-l-v"
(examples/water-wb97xlv.toml). [dft] lambda and omega override the
published 0.6 and 0.1 Bohr⁻¹; omitting them gives the published values.
λ enters the amplitude equations as well as the energy (the paper's eqn 22), so
the correlation term is quadratic in λ at leading order. Validated against PySCF
with the published parameters plus a numpy LinLCCD(hh) on water and OH, and
against the paper's Be₂ bond energy (2.3 kcal/mol; ferric 2.299).
MP2-based double hybrids
B2PLYP and DSD-PBEP86: a KS reference with weighted exchange and
correlation components, plus scaled (SCS-)RI-MP2 correlation.
method.kind = "b2plyp" / "dsd-pbep86" (examples/water-b2plyp.toml);
Python ferric.run_double_hybrid(mol, bs, aux, kind="b2plyp"). Spike: no
comparison to a reference code yet.
Implementation
The equations each solver implements, with file and line, are in
crates/ferric-cc/theory.md.
All contractions go through einsum!, a macro that maps tensor
contractions onto BLAS3 GEMMs. The permutation copies that feed those GEMMs are
parallelized, because for a strided permutation the copy can dominate the
contraction it feeds: measured at 47% at nv = 40 and 70% at
nv = 80. The copies are bit-identical regardless of thread count,
since a permutation writes each output element once; a test pins this.
Memory
The amplitude tensors dominate and grow as \( n_o^2 n_v^2 \), or \( (2n_o)^2 (2n_v)^2 \) in the spin-orbital drivers. Memory budgets are enforced: an oversized job is refused with a breakdown naming the dominant term instead of being OOM-killed partway through.
Cite
CCSD: Scuseria, Janssen & Schaefer 1988; spin-adapted closed-shell CCSD equations: Hirata et al. 2004. (T): Raghavachari et al. 1989; closed-shell (T) algorithm: Rendell, Lee & Komornicki 1991. LinLCCD(hh): Carter-Fenk 2025. ωB97X-L-V: Ransford & Carter-Fenk 2026. B2PLYP: Grimme 2006. DSD-PBEP86: Kozuch & Martin 2011. Review: Bartlett & Musiał 2007. Full entries in References.