vibe-qc is an open-source quantum chemistry and solid-state code — C++17 backend, Python frontend via pybind11, ASE Calculator interface — being written in the open one session at a time. Day 1 shipped the molecular stack and the first Ewald infrastructure. Day 2 closed out v0.2.0 — end-to-end 3D Ewald dispatch, bulk benchmarks, multi-k DIIS, periodic Becke partition — and opened v0.2.5 with SYM3a orbit-reduced lattice-integral storage. Day 3 was narrower by design: a release-target decision and the first piece of the next milestone. 746 tests in, 754 out.
Setting the public-release bar at v0.4.0
The first question of the day was where to draw the public-release line. Three candidates.
v0.2.0 is technically coherent — “3D periodic bulk HF/DFT with quantitative Ewald Coulomb” is a complete story, and we are 95% of the way there (the remaining piece is tightening the CRYSTAL cross-check witness bounds on LiH, NaCl, MgO, and Si). But shipping there would hand a user a code that oscillates the moment they try a metal or a narrow-gap oxide. The DIIS accelerator that works beautifully on insulators stalls on systems with near-degenerate occupied and virtual orbitals near the Fermi level. That is not a missing feature — it is a known failure mode with known remedies, and shipping before those remedies exist would damage the project’s credibility faster than it builds it.
v0.3.0 adds visualisation polish: properties, cube files, better band plots. Useful, but it does not fix the convergence problem. A code that looks nicer while still oscillating on MgO is not a more useful code.
v0.4.0 is where vibe-qc becomes “I can do real solid-state chemistry with this.” The convergence-tooling track (C1: level shifting, Fermi-Dirac smearing, second-order SCF) handles the metals and tight-gap systems. Phase 14 adds effective-core potentials via libecpint, opening the pob-* ECP basis sets for Rb through Lu — metal oxides, perovskites, lanthanoids. Phase 15 adds periodic UHF/UKS for magnetic systems and spin-polarised transition-metal compounds. That combination covers the tutorial cases that currently bounce off the wall, and it matches the capability level where researchers in solid-state chemistry will actually trust the results enough to publish. That is the bar. The full roadmap is on vibe-qc.com.
Phase C1a: Saunders-Hillier level shift in periodic Ewald SCF
The first piece of the v0.4.0 convergence track shipped today. The Saunders-Hillier level shift (Saunders and Hillier, 1973; also Goedecker, 1996) modifies the Fock matrix before diagonalisation:
$$\tilde{\mathbf{F}} = \mathbf{F} + b\left(\mathbf{S} – \tfrac{1}{2}\mathbf{S}\mathbf{D}\mathbf{S}\right)$$
where $b$ is the shift parameter and $\mathbf{D}$ is the density matrix. The effect is to raise virtual MO eigenvalues by $b$ while leaving the occupied block untouched. This widens the effective HOMO-LUMO gap seen by the diagonaliser each iteration, suppressing the near-degenerate orbital swaps that cause oscillation in small-gap systems. The SCF fixed point is unchanged — at convergence, $\mathbf{F}\mathbf{D} = \mathbf{D}\mathbf{F}$ (in the $\mathbf{S}$-metric), so the shift term vanishes and the converged density is the same regardless of $b$.
That last point is the correctness contract that the test suite validates directly. A tight H₂ chain at $a = 8$ bohr (Γ point) and at $a = 10$ bohr with a $[2,2,2]$ k-mesh: $b = 0.3$ reproduces the $b = 0$ total energy to within $10^{-9}$ Ha. The reported mo_energies in the result object are the un-shifted physical orbital energies — the final self-consistency pass at convergence omits the shift, so the orbital spectrum is independent of the $b$ value chosen during iteration. Choosing $b$ too large just slows convergence; choosing it too small provides no benefit on the difficult cases. A value between 0.2 and 0.5 Ha covers most practical situations.
The implementation is wired into both the Γ-Ewald (run_rhf_periodic_gamma_ewald3d) and multi-k Ewald (run_rhf_periodic_multi_k_ewald3d) drivers via a new level_shift field on PeriodicRHFOptions, PeriodicSCFOptions, and PeriodicKSOptions. The default is 0.0 — the field is purely opt-in, and existing callers see bit-for-bit identical results. The SCF convergence section of the user guide covers when and how to use it. Eight new tests: field exposure on all three option structs, Γ-Ewald and multi-k inertness contracts, MO eigenvalues physical at convergence, default behaviour preserved. The molecular RHF/UHF/RKS/UKS level-shift wiring (C1a-2) needs a small parallel-code change in four C++ drivers and is deferred to a follow-up commit.
A segfault, triaged but not fixed
Honest engineering means mentioning this. A crash surfaced in the molecular gradient code: compute_gradient(Molecule, BasisSet, RHFResult) → overlap_gradient_contribution → null function pointer inside __kmp_invoke_microtask. Multiple OpenMP worker threads hitting pc = 0x0 is a strong signal of either a libint2 Engine dispatch table that was not initialized when freshly-copied engines raced into .compute() from worker threads, or an out-of-bounds access into the engine pool when nested or dynamic OpenMP ran more threads than the pool was sized for.
The existing 32 gradient tests all pass cleanly, so the crash triggers on a specific input pattern — most likely a basis set with shells exceeding the prototype’s max_l, or a thread-count change between pool construction and the parallel region. The crash binary’s offsets shifted after rebuild, making direct symbolisation unhelpful. This is deferred: waiting on a reproducer that pins the exact shell configuration. It will be fixed before v0.4.0 tags, because the gradient code is load-bearing for geometry optimisation and eventually periodic forces.
What is next
Two immediate priorities. First, C1b — Fermi-Dirac smearing with fractional occupations and an electronic-entropy contribution to the free energy. This is the single biggest credibility-unblocker for v0.4.0: it covers metals, oxide surfaces, and defect cells with mid-gap states that DIIS plus level shifting still cannot handle. Second, the segfault reproducer — once we have an input that reliably triggers the crash, the fix should be straightforward.
SYM3b (kernel-level compute reduction — the $|G|$-fold wall-clock speedup from space-group symmetry, not just the memory saving from SYM3a) remains open on the v0.2.5 track and will slot in once C1b is done. The CRYSTAL cross-check pass that formally closes v0.2.0 is also still outstanding.
The code is at vibe-qc.com. MPL 2.0.