When the tutorial is the product: how vibe-qc treats documentation as a first-class deliverable

Split view: theory section of a vibe-qc tutorial on the left, running Python code on the right

vibe-qc is an open-source quantum chemistry and solid-state code: a C++17 numerical backend, a Python frontend via pybind11, and an ASE Calculator interface that connects it to the wider atomistic simulation ecosystem. The code is being written in the open, one session at a time, with the cyclic cluster model — a method for treating solid-state defects at post-HF accuracy — as its eventual destination. But the thing that distinguishes vibe-qc from other projects at this stage is not its roadmap length. It is that the tutorials are treated as first-class deliverables, written concurrently with the code they document, and held to the same quality gate as the implementation itself.

This post is about why that choice was made and what it looks like in practice.

What a tutorial is, in this project

Every tutorial in vibe-qc follows the same structure: a working Python script, representative output, an embedded figure where the calculation produces something visual, a theory section (roughly half a page, three to six equations), a references block, and a pointer to the next tutorial in the sequence. None of these sections are optional. A tutorial without a theory section is a recipe. A tutorial without a working script is an essay. The combination is what makes it an educational artefact rather than documentation.

The theory sections are typeset with explicit equations. Since peintinger.com does not currently run a MathJax or KaTeX plugin, the equations in the docs themselves are rendered by MathJax inside the Sphinx build at vibe-qc.com. For example, the molecular Hartree-Fock tutorial works through the Roothaan equations in matrix form, $\mathbf{F}\mathbf{C} = \mathbf{S}\mathbf{C}\boldsymbol{\varepsilon}$, where $\mathbf{F}$ is the Fock matrix, $\mathbf{C}$ the MO coefficient matrix, $\mathbf{S}$ the overlap matrix, and $\boldsymbol{\varepsilon}$ the diagonal matrix of orbital energies:

Split view: theory section of a vibe-qc tutorial on the left, running Python code on the right
The dual character of a vibe-qc tutorial: theory on the left (Bloch theorem, dispersion relation, DOS formula, foundation references), running code on the right. Both sections are required; neither is optional.

The citation discipline matters as much as the equations. Foundation-era papers — Roothaan 1951, Bloch 1928, Becke 1988, Pulay 1980 — are cited with full journal-volume-page triples drawn directly from the original. Post-2010 results that are harder to verify from memory carry a [Ref: verify] marker in the draft and are only promoted to full citation after the volume and page have been confirmed against the journal. This is not the standard practice for informal software documentation. It is standard practice for a journal article, and the tutorials are being written to that standard deliberately.

The parity matrix as a dual roadmap

The roadmap page carries two distinct structures. The first is an engineering milestone list tracking Ewald phases, symmetry sub-phases, and versioned capabilities. The second is a tutorial parity matrix that pins vibe-qc against ORCA’s molecular tutorial set and CRYSTAL’s solid-state tutorial set — the two reference programs that cover the two halves of what vibe-qc is trying to do.

Table showing vibe-qc tutorial parity against ORCA and CRYSTAL reference programs
The tutorial parity matrix from vibe-qc.com/roadmap.html. Green means a tutorial covering this capability is shipped and cross-checked against the reference program. Yellow means partial. Red means the capability is on the roadmap but the tutorial does not yet exist.

The value of the parity matrix is that it answers two different questions from a single table. For a user evaluating whether vibe-qc can support their work, it gives a one-page answer: “can I reproduce the ORCA tutorial on geometry optimisation with this code?” For a developer deciding what to build next, it gives a prioritised queue: the red cells in the CRYSTAL column are where the periodic capability is missing and the tutorials cannot yet be written. The engineering roadmap and the documentation roadmap are the same document.

The framing is deliberately not competitive. The parity matrix does not claim vibe-qc exceeds ORCA or CRYSTAL at anything. It says, plainly, where vibe-qc is today relative to what those programs can teach, and what it would take to close the gap. That honesty is part of the project’s credibility.

Fourteen tutorials, all shipped today

As of the current tutorial index, fourteen tutorials are live. They run from molecular Hartree-Fock — the entry point, which teaches the Roothaan equations, basis set input, and how to read an SCF trace — through open-shell methods, geometry optimisation, vibrational frequencies, thermodynamics at finite temperature, orbital visualisation, basis set convergence, and dispersion corrections, and then into the periodic side: periodic HF on a 1D hydrogen chain, periodic KS-DFT, Madelung constants via Ewald summation, and band structure and density of states.

The band structure tutorial is worth looking at as the canonical example of what the format can carry. It delivers a full Γ-to-X k-path calculation on the hydrogen chain, plots bands and DOS in a single combined figure, explains the Bloch theorem, $\psi_{n\mathbf{k}}(\mathbf{r}) = e^{i\mathbf{k}\cdot\mathbf{r}}u_{n\mathbf{k}}(\mathbf{r})$, and the origin of band dispersion in terms of tight-binding bonding and antibonding combinations, cites Bloch 1928 and Monkhorst-Pack 1976, and ends with a pointer forward to the Peierls distortion. The code block is short enough to read in under a minute. The theory section is long enough to understand what the code is computing and why the result looks the way it does.

That rhythm — code, output, figure, theory, references, next — is the same across all fourteen tutorials. A student who works through them in order gets a coherent course in molecular and periodic quantum chemistry. A researcher who arrives at a specific tutorial looking for a recipe gets that too, because the code blocks are self-contained. The two use cases are not in tension; they are served by the same document structure.

The no-code-required backlog

One discipline the roadmap makes explicit is the distinction between tutorials that require new code and tutorials that require only documentation work on existing capability. Several entries in the tutorial backlog fall into the second category: the capability is already in the code, the API is stable, and the only work is writing the tutorial itself — the explanation, the example script, the figure, the theory section, the references.

This is worth naming because it is easy to defer. The code works; the tutorial can wait. What the project has committed to instead is treating those deferred tutorials as a first-class backlog item, tracked in the parity matrix, with the same visibility as unimplemented features. A shipped capability without a tutorial is a capability that most users will never find.

The user guide plays a complementary role here. It covers the reference-level detail — how to specify k-point meshes, how to configure Ewald parameters, how to read the memory budget report — that would clutter a tutorial but that a user needs once they move past the introductory examples. The split between tutorial and user guide is deliberate: tutorials teach phenomena and procedures, the user guide documents options and behaviour.

Plots as the lesson, not decoration

Every tutorial that produces a figure ships that figure embedded in the page. The figures are not screenshots pasted in after the fact. They are generated by scripts in the repository, regenerable by any user who runs the tutorial code, and updated whenever the underlying calculation changes.

This matters because the figures are doing pedagogical work that the code blocks cannot. A code block teaches syntax and API. A band structure plot teaches what band dispersion looks like, where the gap sits, how the DOS integrates the band structure into a density, and what it means for a state at the zone boundary to be antibonding. The band structure tutorial includes both panels — bands and DOS in a single combined figure — because showing them side by side is how a reader learns to interpret both simultaneously.

The same logic applies to the dispersion tutorial, which shows a binding curve with and without D3-BJ correction, making visible the underbinding failure of uncorrected GGAs on weakly-bound systems. The figure is the argument. The code block shows how to reproduce it.

The doc build as a gate

The documentation is built with sphinx-build -W --keep-going: warnings are errors, and the build does not finish if any cross-reference is broken or any code block fails to parse. This is the same discipline applied to the test suite — 746 passing tests as of the most recent commit — extended to the documentation layer. A tutorial that links to a non-existent API page, or a code block that uses a function name that was renamed, fails the build and blocks the merge.

This is not the default posture for scientific software projects. The default is to build docs as a separate step, treat failures as non-blocking, and accumulate rot over time. The cost of running docs under the same gate as code is low — a few minutes added to CI. The benefit is that the tutorials stay synchronized with the implementation without anyone having to remember to update them.

Where this fits in the larger project

vibe-qc is aiming at the cyclic cluster model: a method that treats a finite cluster of atoms embedded in a periodic Madelung field, enabling post-HF correlation calculations on solid-state defects at a cost that scales with cluster size rather than unit cell size. The cyclic cluster model is v2.0 on the roadmap. Every tutorial from the molecular HF entry point through the periodic band structure is load-bearing preparation for the moment when a researcher can open a CIF file, define a vacancy cluster, run CCSD(T) on it, and understand what the result means.

The tutorials exist not to fill a documentation requirement but because the project’s audience includes people who need to understand what they are computing, not just how to invoke the code. That is a different audience from “users who already know periodic HF/DFT and just need the API reference.” The user guide and API reference serve that second audience. The tutorials serve everyone, and they are the harder document to write well.

The code lives at vibe-qc.com. MPL 2.0.

vibe-qc day 2: closing out Ewald and laying the symmetry foundation

vibe-qc day 2: closing out Ewald and laying the symmetry foundation

vibe-qc is an open-source quantum chemistry and solid-state code with a C++17 numerical backend, a Python frontend via pybind11, and an ASE Calculator interface for geometry optimization and structure I/O. It is being written in the open, one coding session at a time, with Claude Opus as the implementation engine and me steering at the architecture level. Day 1 shipped the molecular HF/DFT/MP2 stack (v0.1.0, April 18) and then pushed into periodic territory: one-electron lattice integrals, Gamma-only and multi-k RHF, multi-k KS-DFT, and the first Ewald infrastructure through Madelung-cancellation helpers. Day 2 closed out the v0.2.0 Ewald milestone and planted the first stake for v0.2.5 space-group symmetry. Four deliverables, 26 new tests (720 to 746 total), Sphinx docs building clean.

Phase 12e-c-4: end-to-end EWALD_3D dispatch

The previous session had the Ewald machinery in place but not yet wired to the user-facing entry points. Today that changed. Both run_rhf_periodic_scf(...) and run_rhf_periodic_gamma_scf(...) now branch on options.lattice_opts.coulomb_method: pass DIRECT_TRUNCATED and you get the existing C++ driver unchanged; pass EWALD_3D and you get the new Python driver, with nuclear_repulsion_per_cell automatically rerouted through Ewald summation as well. That last part matters: an Ewald electronic energy paired with a direct-truncation nuclear repulsion would be self-consistent in neither the long-range cancellation nor the ω-dependence, and the resulting total energy would be meaningless. Electronic and nuclear sides now always agree on which summation scheme is in use.

The benchmark suite that shipped alongside covers H₂ crystal, LiH rocksalt, MgO rocksalt, and Ne FCC — 11 tests in total, checking convergence, $\omega$-invariance (energy independence of the Ewald splitting parameter), an equation-of-state scan, and the equivalence of a [1,1,1] multi-k mesh with a Gamma-only calculation. All pass.

One finding worth documenting: atom positions relative to the FFT Poisson grid origin are not numerically irrelevant. Atoms sitting at box corners, rather than centered in the cell, inflate the $\omega$-invariance residual by roughly two orders of magnitude. The physics is correct either way, but the numerical conditioning is not. This is now flagged in the docs and the benchmark suite uses centered geometries as the reference.

Phase 12e-c-4c-iii-c: multi-k Pulay DIIS

The multi-k Ewald driver shipped without DIIS. For loose cells and [1,1,1] meshes that was acceptable — pure density damping converged, slowly. For tight cells with denser k-meshes it was not. An H₂ chain at a = 10 bohr with a [2,2,2] mesh simply plateaued under damping and never converged.

The fix is Pulay DIIS extended to the Brillouin zone. The commutator error vector at each k-point is:

$$\mathbf{e}(\mathbf{k}) = \mathbf{F}(\mathbf{k})\mathbf{D}(\mathbf{k})\mathbf{S}(\mathbf{k}) – \mathbf{S}(\mathbf{k})\mathbf{D}(\mathbf{k})\mathbf{F}(\mathbf{k})$$

and the Pulay B matrix couples iterations i and j with a k-weighted inner product:

$$B_{ij} = \sum_{\mathbf{k}} w_{\mathbf{k}}\, \mathrm{Re}\,\mathrm{tr}\!\left[\mathbf{e}_i(\mathbf{k})^\dagger\, \mathbf{e}_j(\mathbf{k})\right]$$

The k-weights are the standard Monkhorst-Pack weights that already appear in the energy and density assembly, so no new parameters appear. The DIIS coefficients solve the same linear system as in the molecular case; the only change is that the scalar error measure is now a sum over the zone rather than a single matrix norm.

The numbers are clean. H₂ chain, [1,1,1]: 13 iterations under damping, 7 with DIIS. H₂ chain, [2,2,2]: previously non-converging under damping, now 4 iterations with DIIS. This is the same qualitative behavior the molecular DIIS delivered on H₂O back on day 1 (52 iterations to 9), just extended to k-space.

Phase 12f: periodic Becke partition for tight-cell DFT

The molecular Becke fuzzy-cell partition (J. Chem. Phys. 88, 2547, 1988) assigns each grid point a weight that partitions unity across atomic cells using a smooth step function of interatomic distances. In a periodic system with a tight unit cell, image atoms from neighboring cells fall within the smoothing radius of the step function, and ignoring them breaks the partition: the weights no longer sum to the cell volume.

The fix is to extend the partition denominator over home atoms plus image atoms within a user-specified radius. The implementation is wired into run_rks_periodic via two new options: PeriodicKSOptions.use_periodic_becke (boolean, default False) and becke_image_radius_bohr (float). Default behavior is unchanged from v0.1.0; the periodic extension is an explicit opt-in for users who know they are working in a tight cell.

The empirical witness is a 5-bohr cubic H₂ cell. Without periodic Becke, the molecular partition produces a total grid weight of 18,526 — wildly wrong for a cell with volume 125 bohr³. With periodic Becke, the total grid weight comes out as 126, matching V_cell = 125 bohr³ to within the grid quadrature error. The discrepancy between the two modes (18,526 vs 126) is not a subtle numerical issue; it is a categorical correctness failure of the molecular partition in tight periodic geometry. Any DFT calculation on a tight crystal without this correction is integrating a density that does not integrate to the correct electron count per cell.

Phase SYM3a: orbit-reduced storage for lattice integrals

The v0.2.5 symmetry milestone starts here. The one-electron lattice integrals $S(\mathbf{g})$, $T(\mathbf{g})$, $V(\mathbf{g})$ are indexed by a cell offset vector $\mathbf{g} = (g_1, g_2, g_3)$. For a system with space-group symmetry, many (g, atom-pair) combinations are related by point-group operations and carry the same numerical content up to a known rotation. SYM3a exploits this to compress storage.

The machinery builds on SYM2c, which already identifies atom-pair and cell-index orbits for arbitrary (including non-origin-fixed) structures like NaCl. SYM3a takes those orbits and stores one representative sub-block of shape $(n_{\mathrm{AO},a},\, n_{\mathrm{AO},b})$ per orbit, rather than a full $(n_{\mathrm{bf}} \times n_{\mathrm{bf}})$ block at every cell triple. A round-trip compress/reconstruct test on simple cubic Pm-3m He at STO-3G, with a 10 bohr cutoff, partitions 33 cell triples into 5 orbits and achieves 6.6× memory reduction. Reconstruction is exact to 10⁻¹⁶.

What SYM3a does not yet do is reduce compute: the integral kernels still evaluate every shell-pair times cell-triple combination before handing off to the compressor. That is SYM3b, the next item on the list. SYM3b is where the |G|-fold reduction hits wall-clock time rather than just peak memory, because the kernel skips shell-pair/cell-triple combinations that are equivalent under the point group instead of evaluating and then discarding them. The storage handle SYM3a provides is what SYM3b will plug into.

Where v0.2.0 stands

The v0.2.0 contract is “ω-invariant 3D periodic SCF that reduces to molecular HF in the loose-box limit and integrates the unit-cell volume correctly in tight cells.” As of today that contract is end-to-end testable with two option flags:

options.lattice_opts.coulomb_method = EWALD_3D
options.ks_opts.use_periodic_becke = True

Internal-consistency benchmarks are in place across four crystal systems. The remaining gating item for the v0.2.0 tag is a cross-check pass against published CRYSTAL reference energies on LiH, NaCl, MgO, and Si at published geometries. That is a validation exercise against an external reference rather than new implementation work, and it is next on the list.

What is next

Immediately: SYM3b, the kernel-level compute reduction that makes the |G|-fold saving from space-group symmetry show up in wall-clock time. The substrate — Wigner D-matrices (SYM1), AO permutation matrices and orbit identification (SYM2a/b), atom-pair-resolved orbits for non-origin-fixed structures (SYM2c), and now the compressed storage layer (SYM3a) — is all in place. SYM3b is the payoff.

Further out, the headline feature of the entire project remains the cyclic cluster model. The roadmap runs from v1.0 (feature-complete molecular and periodic HF/DFT/MP2) through v2.0 (HF-CCM in 3D) and on through a sequence of independently publishable steps: MP2-CCM, local-MP2-CCM, CCSD-CCM, and eventually projection-based embedding where a CCM-correlated central cluster sits inside a DFT-treated periodic environment. The capstone is multi-scale defect chemistry at chemical accuracy on systems with hundreds of atoms in the cluster. Every piece shipped in the first two days of this project is load-bearing for that goal.

The code is at vibe-qc.com. Licensed MPL 2.0.