Periodic HF on a 1D H chain

The simplest non-trivial periodic system: a one-dimensional chain of H₂ dimers.

import numpy as np
from vibeqc import (
    Atom, BasisSet, PeriodicSystem, PeriodicSCFOptions,
    monkhorst_pack, run_rhf_periodic,
)

# Unit cell: H2 oriented along the chain, 4.0 bohr between unit cells.
unit_cell = [Atom(1, [0, 0, 0]), Atom(1, [0, 0, 1.4])]
lattice = np.diag([4.0, 30.0, 30.0])   # 1D axis on x; y, z are vacuum

sysp = PeriodicSystem(
    dim=1,
    lattice=lattice,
    unit_cell=unit_cell,
)

basis = BasisSet(sysp.unit_cell_molecule(), "pob-tzvp")
kmesh = monkhorst_pack(sysp, [6, 1, 1])

opts = PeriodicSCFOptions()
opts.conv_tol_energy = 1e-10
opts.lattice_opts.cutoff_bohr = 12.0         # μν overlap cutoff
opts.lattice_opts.nuclear_cutoff_bohr = 40.0  # longer-range nuclear sum

result = run_rhf_periodic(sysp, basis, kmesh, opts)
print(f"E_RHF per cell = {result.energy:.10f}  ({result.n_iter} iters)")

What the options mean

  • PeriodicSystem.dim — number of periodic axes. dim=1 makes only lattice.col(0) a lattice vector; the other two columns are interpreted as vacuum with enough separation (30 bohr here) that images don’t overlap. dim=2 treats cols 0 and 1 as in-plane lattice vectors. dim=3 is full bulk.

  • lattice_opts.cutoff_bohr — the real-space distance beyond which (μν) pair integrals are dropped from the lattice sum. Bump it up until your energy stops changing.

  • lattice_opts.nuclear_cutoff_bohr — separate, usually larger, cutoff for the nuclear-attraction term (1/r tail decays more slowly than AO overlap).

Choosing a k-mesh

For the SCF to represent an insulator correctly, the mesh needs to sample enough of the Brillouin zone that the occupied manifold stays separated from the virtual. Start with [N, 1, 1] for 1D, [N, N, 1] for 2D, [N, N, N] for 3D, and increase N until the energy per cell stops changing. For small H-chain tests a 4×1×1 mesh is usually plenty; for covalent 3D crystals in TZ basis, 4×4×4 is a reasonable starting point.

What’s validated today

  • Molecular limit: a big-box (50 bohr) unit cell produces the same energy as molecular RHF — machine precision across dim ∈ {1, 2, 3} and multiple k-mesh sizes.

  • Bloch-sum machinery: kinetic-only folding equivalence between 1-cell × K k-points and K-cell × Γ, machine precision.

  • Non-trivial 1D SCF convergence at tight cell spacings.

3D bulk calculations with this direct-truncated Coulomb sum are cutoff-dependent (the infinite Coulomb series is conditionally convergent in 3D); quantitative benchmarking against CRYSTAL waits for the Ewald splitting in Phase 12e. See the roadmap.