Multi-k periodic SCF¶

vibe-qc treats multi-k periodic SCF as a native feature, not as a Python import wrapper around another quantum-chemistry code. PySCF and CRYSTAL are reference programs: run them out-of-process, parse their outputs, and compare energies, occupations, bands, and convergence traces. They are not vibe-qc backends.

Current Status¶

Native, usable today:

run_rhf_periodic_gamma_gdf(...) for Γ-only native Gaussian density fitting RHF, RKS, and hybrids across PeriodicSystem.dim = 1, 2, 3.
run_krhf_periodic_gdf(..., kmesh, ...) and run_krks_periodic_gdf(..., kmesh, ...) as the public k-point API for native GDF in 1D, 2D, and 3D. A Γ-only mesh (1, 1, 1) delegates to the Γ-GDF implementation; any other Monkhorst-Pack mesh now runs the full native multi-k GDF loop — k-dependent 3c/2c blocks, Bloch-phase assembly, multi-k DIIS. The multi-k Ewald-gauge fix brings a multi-k mesh to within 0.000 mHa of the Γ-point energy in the molecular limit. Also reachable through the high-level run_periodic_job(..., gdf_kmesh=(n, n, n)).
compute_2c_eri_lattice_blocks(...) and compute_3c_eri_lattice_blocks(...) expose translation-resolved periodic DF blocks in 1D, 2D, and 3D. Summing those blocks reproduces the Γ metric/tensor; the native multi-k GDF loop will apply Bloch phases to these blocks instead of rebuilding libint shell pairs. bloch_sum_2c_eri_blocks(...) and bloch_sum_3c_eri_blocks(...) pin that phase convention, and build_lpq_bloch_native(...) fits the phase-assembled tensor in the phase-assembled auxiliary metric.
run_rhf_periodic_scf(..., kmesh, opts) with opts.lattice_opts.coulomb_method = CoulombMethod.EWALD_3D.
run_rks_periodic_scf(..., kmesh, opts) with the same EWALD_3D route.
General lattice vectors: the code works with arbitrary 3D Bravais cells and uses the supplied lattice matrix to build reciprocal-space k-points and Bloch sums. Tests include skew/hexagonal cells so the implementation does not silently assume cubic axes.
Closed-shell finite-temperature Fermi-Dirac occupations for periodic multi-k Ewald RHF/RKS via opts.smearing_temperature (k_B T in Hartree). The native Γ-only GDF RHF/RKS driver and high-level run_periodic_job(smearing_temperature=...) support the same occupation model for Γ-point metal/debug work.

Not native yet:

Periodic UKS smearing is not implemented yet. A positive smearing_temperature on the UKS multi-k Ewald path raises NotImplementedError so metallic spin-polarized inputs do not run with silently integer occupations.
Production-scale all-electron MgO/pob-TZVP at dense k-meshes needs the native FFTDF/GDF path. The legacy Ewald RKS path is useful for debug and correctness work, but its periodic XC memory footprint is intentionally guarded.

Minimal RHF Example¶

import numpy as np
import vibeqc as vq

a = 10.0
system = vq.PeriodicSystem(
    3,
    np.eye(3) * a,
    [vq.Atom(1, [0.0, 0.0, 0.0]),
     vq.Atom(1, [0.0, 0.0, 1.4])],
)
basis = vq.BasisSet(system.unit_cell_molecule(), "sto-3g")
kmesh = vq.monkhorst_pack(system, [2, 2, 2])

opts = vq.PeriodicRHFOptions()
opts.lattice_opts.coulomb_method = vq.CoulombMethod.EWALD_3D
opts.lattice_opts.cutoff_bohr = 12.0
opts.lattice_opts.nuclear_cutoff_bohr = 15.0
opts.damping = 0.3
opts.use_diis = True

result = vq.run_rhf_periodic_scf(
    system, basis, kmesh, opts,
    omega=0.5,
    spacing_bohr=0.3,
)
print(result.energy, result.converged)

Minimal RKS Example¶

opts = vq.PeriodicKSOptions()
opts.functional = "pbe"
opts.lattice_opts.coulomb_method = vq.CoulombMethod.EWALD_3D
opts.lattice_opts.cutoff_bohr = 12.0
opts.lattice_opts.nuclear_cutoff_bohr = 15.0
opts.damping = 0.3
opts.use_diis = True

result = vq.run_rks_periodic_scf(
    system, basis, kmesh, opts,
    omega=0.5,
    spacing_bohr=0.3,
)
print(result.energy, result.e_xc, result.converged)

For MgO/pob-TZVP-class RKS, prefer this route only for small debug meshes until native FFTDF/GDF lands; the full periodic XC grid can be too memory-heavy for the legacy Ewald driver.

Smearing For Metals¶

Closed-shell periodic RHF/RKS support Fermi-Dirac smearing on the multi-k Ewald path and on the native Γ-only GDF path:

opts.smearing_temperature = vq.resolve_smearing_temperature(
    1000.0,
    unit="kelvin",
).temperature

The option stores k_B T in Hartree. You can set it directly, convert from Kelvin/eV/Ry, use a preset, or ask for a conservative automatic guess:

opts.smearing_temperature = vq.resolve_smearing_temperature("metal").temperature
opts.smearing_temperature = vq.resolve_smearing_temperature("0.1 eV").temperature
opts.smearing_temperature = vq.resolve_smearing_temperature(0.1, unit="eV").temperature
opts.smearing_temperature = vq.resolve_smearing_temperature(
    "auto",
    metallic=True,
).temperature

Numeric strings may carry their own unit ("1000 K", "0.1 eV", "0.01 Ha", "0.02 Ry"), which is often clearer in input files than a bare Hartree number. "auto" is intentionally cautious: if no metallic hint or prior band gap is supplied, it returns zero smearing rather than silently changing an insulating calculation. With smearing enabled, the drivers bisect the chemical potential so that the weighted k-mesh occupation satisfies the per-cell electron count:

sum_k w_k sum_i n_i(k) = N_electrons
n_i(k) = 2 / (1 + exp((eps_i(k) - mu) / T))

The result object reports:

result.energy: internal total energy.
result.free_energy: variational free energy E - T S.
result.entropy: dimensionless electronic entropy S/k_B.
result.fermi_level: chemical potential in Hartree.
result.occupations: per-k occupation arrays in [0, 2].

DIRECT_TRUNCATED currently rejects positive smearing_temperature; that backend would otherwise ignore the field and run a hard-Aufbau calculation.

For the high-level native-GDF job runner, pass the same value as a keyword:

result = vq.run_periodic_job(
    system,
    basis,
    method="RKS",
    functional="pbe",
    smearing_temperature="auto",
    smearing_metallic=True,
)

CRYSTAL And PySCF Parity¶

CRYSTAL and PySCF remain essential references, but the validation boundary is external:

Write a CRYSTAL or PySCF input next to the vibe-qc input.
Run the external program as a separate executable or script.
Parse its output into a comparison record.
Store method, basis, lattice, k-mesh, energy, occupation, and convergence metadata in the benchmark manifest.

This keeps licensing, dependency, and numerical-responsibility boundaries clean. It also makes it obvious whether a discrepancy is in vibe-qc or in an adapter/parser.

Native FFTDF/GDF Roadmap¶

The production target is a native periodic stack:

FFTDF first for speed and plane-wave-style Coulomb handling.
GDF for all-electron Gaussian density fitting and CRYSTAL/PySCF parity benchmarks.
Range-separated exchange screening for hybrids and HF in large cells.
Disk-backed Cholesky/DF tensors for correlated methods later.

The public run_krhf_periodic_gdf / run_krks_periodic_gdf names are the stable API. A Γ-only mesh runs through the native Γ-GDF bridge; any other Monkhorst-Pack mesh runs the native multi-k GDF loop. The remaining roadmap work is performance and dense-solid scaling (disk-backed tensors, IBZ-weighted J/K), not correctness of the multi-k path itself.