Native Molecular Geometry Optimizer¶

vibe-qc ships a standalone molecular geometry optimizer that needs no ASE. It wraps analytic SCF nuclear gradients where they are validated and uses central finite differences of the full energy for wavefunction methods without an FD-tight analytic gradient. It is available both as a library function and through run_job. A conservative steepest-descent + Brent line-search backend is also available for systems where the Hessian approximation misbehaves.

Quick start¶

from vibeqc import Molecule, Atom
from vibeqc.molecular_optimize import optimize_molecule

mol = Molecule([
    Atom(8, [ 0.00,  0.00,  0.00]),
    Atom(1, [ 0.00,  1.43, -0.98]),
    Atom(1, [ 0.00, -1.43, -0.98]),
])

result = optimize_molecule(
    mol,
    basis_name="def2-svp",
    method="rks",
    functional="PBE",
)
print(result.system)           # optimised Molecule (bohr)
print(result.energy)           # final energy (Ha)
print(result.n_iter)           # number of BFGS steps
print(result.converged)        # True / False

Through `run_job`¶

Pass optimizer_backend="native" to bypass ASE:

from vibeqc import run_job

run_job(
    mol,
    basis="def2-svp",
    method="rks",
    functional="PBE",
    optimize=True,
    optimizer_backend="native",   # ← no ASE needed
    output="h2o_opt",
)

The default optimizer_backend="auto" prefers ASE when installed and falls back to the native path otherwise, so existing workflows are unchanged.

Backends¶

Three geometry-optimizer backends are available through optimizer_backend= (both via run_job and directly):

Backend	`optimizer_backend=`	How it works	Best for
ASE BFGS	`"ase"` (or `"auto"` when ASE installed)	ASE’s `BFGSLineSearch`, quasi-Newton with Wolfe-condition line search	Routine use; fastest convergence near minimum
Native L-BFGS-B	`"native"`	scipy’s L-BFGS-B, limited-memory BFGS with box constraints	No-ASE workflows; frozen atoms
Brent steepest-descent	`"brent"`	Steepest-descent direction + Brent 1-D line search per step	Flat/dispersion-bound PESs where Hessian extrapolation misbehaves

The Brent backend never takes an uphill step (each line search is a rigorous 1-D minimisation) and needs no Hessian approximation. Its cost is higher per geometry step (multiple energy evaluations per line search), but it is exceptionally robust on difficult surfaces. Use it when the quasi-Newton backends fail to converge.

# Via run_job:
run_job(mol, basis="def2-svp", method="rks", functional="PBE",
        optimize=True, optimizer_backend="brent", output="h2o_brent")

# Direct API:
from vibeqc.molecular_optimize import optimize_molecule_brent
result = optimize_molecule_brent(mol, basis_name="def2-svp", method="rks",
                                  functional="PBE")

The native optimizer family (`geom_opt=`)¶

Beyond the three backends above, run_job exposes the uniform vibeqc.geomopt framework through geom_opt=. Each optimizer walks the same provider surface (SCF energy + analytic gradient, with dispersion and solvation folded in when requested):

`geom_opt=`	Algorithm	Target
`"sd"`	Steepest descent + line search	minimum
`"cg"`	Conjugate gradient	minimum
`"bfgs"`	BFGS quasi-Newton	minimum
`"lbfgs"`	Limited-memory BFGS	minimum
`"trust"`	Trust-region Newton	minimum or TS
`"rfo"`	Rational function optimization	minimum
`"gdiis"`	Geometry DIIS	minimum
`"fire"`	Fast inertial relaxation (FIRE)	minimum
`"ef"`	Eigenvector following	minimum or TS
`"prfo"`	Partitioned RFO	transition state
`"dimer"`	Dimer method (gradient-only saddle search)	transition state

run_job(mol, basis="def2-svp", method="rks", functional="PBE",
        optimize=True, geom_opt="bfgs", fmax=0.01,
        output="h2o_native_bfgs")

# Transition-state search (ef / prfo / dimer / trust):
run_job(mol, basis="def2-svp", method="rks", functional="PBE",
        optimize=True, geom_opt="ef", geom_target="transition_state",
        output="ts_search")

geom_target="transition_state" selects saddle-point mode; geom_opt="dimer" requires it (a dimer search of a minimum is refused with a ValueError). Per-optimizer tuning knobs pass through geom_opt_options={...}, and the initial-Hessian policy through geom_hessian_init= / geom_hessian_update=.

Every optimizer streams a per-step progress table into the job .out while it runs, so long queue jobs can be monitored live:

  step            E (Ha)           dE      max|g|      |step|  conv
  evaluating initial energy and gradient ...
    -74.96517778       --      1.9577e-02       --     gmax=1.958e-02✗
    -74.96565587   -4.781e-04  6.9972e-03   2.891e-02  gmax=6.997e-03✗
    -74.96579842   -1.426e-04  4.8506e-03   1.633e-02  gmax=4.851e-03✗
    -74.96589991   -1.015e-04  1.3428e-03   2.670e-02  gmax=1.343e-03✗
    -74.96590107   -1.152e-06  3.2652e-04   1.673e-03  gmax=3.265e-04✓

(H2O/STO-3G RHF, geom_opt="bfgs".) The conv column shows each active convergence criterion with its current value and pass/fail state, plus rejected markers for trust-region step rejections; the optimization stops once all active criteria pass.

Coordinate systems (`geom_coords=`)¶

The native optimizer family walks Cartesian coordinates by default. geom_coords="dlc" selects delocalised internal coordinates (bonds / angles / torsions combined into a non-redundant set), which typically converge in fewer steps for covalently bonded molecules:

run_job(mol, basis="def2-svp", method="rks", functional="PBE",
        optimize=True, geom_opt="bfgs", geom_coords="dlc",
        output="h2o_dlc")

The DLC back-transform rebuilds the Wilson B-matrix pseudoinverse at the current geometry each Newton iteration, wraps torsion residuals onto the correct 2π branch, and seeds each step from the previous geometry, so large steps and near-planar torsions (±π branch cut) stay in the Newton convergence basin. The auto-generated primitive set does not yet guarantee completeness for every topology: when it spans fewer than the expected 3N−6 internal degrees of freedom, construction emits a RuntimeWarning and the optimization may not fully converge; use the default geom_coords="cartesian" for such systems.

Each step’s SCF warm-starts from the previous step’s converged density (mean-field methods, gas phase): the provider hands the prior result to the SCF as an initial_guess=READ restart, cutting the per-step iteration count substantially since optimizer steps are small. A warm-started SCF that fails to converge is retried once from the cold default guess; a still-nonconverged SCF raises a clear error instead of feeding an unconverged gradient to the optimizer. Opt out (e.g. to reproduce older runs’ iteration counts) with MolecularSCFProvider(..., warm_start=False) when driving run_geomopt directly.

Supported methods¶

Method	Gradient	Notes
`rhf`, `uhf`	Analytic	All-electron, closed- and open-shell
`rks`, `uks`	Analytic	All DFT functionals with analytic gradients
`selected_ci`, `dmrg`, `v2rdm`, `transcorrelated_ci`	Central FD	2-point finite difference on energy (`fd_step_bohr=0.005`)
`casci`	Central FD	2-point FD on energy
`casscf`	Analytic / Central FD	Validated analytic gradient for state-specific, closed-shell, default-`compute_wz` CASSCF (v0.15.0); SA-CASSCF, open-shell, and `compute_wz=True` fall back to full-energy FD (see the note below)
`caspt2`, `nevpt2`	Central FD	Relaxed full-energy central FD for CASSCF-referenced, gas-phase runs with `compute_corr_grad=True` on `caspt2_options` / `nevpt2_options`; otherwise the optimizer owns the same full-energy FD fallback (see the note below)

Dispersion corrections (D3-BJ) and implicit solvation (CPCM/COSMO) are folded into the energy and gradient automatically when passed. Implicit solvation composes with the mean-field methods only (rhf / uhf / rks / uks); requesting solvent=... with the CAS-family methods or rohf raises a clear ValueError (there is no CPCM composition for those methods, and silently optimizing the gas-phase surface instead would misreport the result).

On both backends every per-step energy evaluation receives the same solver options as the final single point, active_space, cas_reference, and the wavefunction option structs (selected_ci_options, dmrg_options, v2rdm_options, transcorrelated_options, casci_options, casscf_options, caspt2_options, nevpt2_options), so an SA-CASSCF optimization (casscf_options=CASSCFOptions(nroots=2)) walks the state-averaged surface its final energy is reported on, and a selected_ci optimization keeps its active-space truncation at every FD displacement instead of falling back to full-space CI.

As of v0.15.0 the state-specific CASSCF analytic gradient is a validated full-energy derivative (examples/regression/casscf_gradient_fd_reproducer.py passes to ~1e-7 Ha/bohr), and all three molecular backends (the L-BFGS-B primary, the Brent backend optimizer_backend="brent", and the uniform geom_opt framework) use it for state-specific, closed-shell, default- compute_wz CASSCF. They share one decision so they always walk the same surface. Outside that validated envelope they fall back to full-energy central FD: state-averaged CASSCF (casscf_options=CASSCFOptions(nroots=2), whose analytic gradient is only finiteness/translational-invariance checked), open-shell CASSCF (the kernel is the closed-shell RHF formalism), and the experimental compute_wz=True CP-MCSCF correction.

CASPT2/NEVPT2 optimizations use the runner-supplied relaxed full-energy central FD gradient when the run requests and can produce it: a CASSCF reference (casscf_options given), gas phase, and compute_corr_grad=True on caspt2_options / nevpt2_options. With the default compute_corr_grad=False the solver only computes the bare CASSCF reference gradient, so the optimizers fall back to their outer full-energy central FD rather than walk a surface inconsistent with the reported PT2 energy; the same fallback applies to CASCI-on-HF-referenced and solvated runs. All three backends share this decision too.

from vibeqc.solvers import CASPT2Options, CASSCFOptions

result = optimize_molecule(
    mol,
    basis_name="6-31g",
    method="caspt2",
    active_space=(2, 2),
    casscf_options=CASSCFOptions(),
    caspt2_options=CASPT2Options(compute_corr_grad=True),  # relaxed FD gradient
)

Trajectory collection¶

Set record_trajectory=True (the default) to collect per-step geometries and energies:

result = optimize_molecule(mol, basis_name="sto-3g", method="rhf")
for i, (frame, e) in enumerate(
    zip(result.trajectory_frames, result.trajectory_energies)
):
    print(f"  step {i}: E = {e:.8f} Ha")

When used through run_job(optimize=True, output_qvf=True), the trajectory is embedded in the QVF archive for vibe-view’s animation player, identical behaviour to the ASE backend.

Convergence control¶

Parameter	Default	Meaning
`conv_tol_grad`	`4.5e-4`	Gradient norm convergence (Ha/bohr)
`conv_tol_energy`	`1e-6`	Energy change tolerance (Ha)
`max_iter`	`100`	Maximum steps

For DFT jobs where the SCF may struggle at intermediate geometries (common with PBE + minimal basis sets), pass the appropriate rks_options / uks_options with increased max_iter:

from vibeqc import RKSOptions

rks_opts = RKSOptions()
rks_opts.max_iter = 80
rks_opts.use_diis = True

result = optimize_molecule(
    mol, basis_name="sto-3g", method="rks", functional="PBE",
    rks_options=rks_opts,
)

API reference¶

class vibeqc.molecular_optimize.MolecularOptimizeResult(system, energy, gradient, n_iter, converged, trajectory_frames=None, trajectory_energies=None)[source]¶

Bases: object

Container for molecular geometry optimization results.

Parameters:

system (Molecule)
energy (float)
gradient (np.ndarray)
n_iter (int)
converged (bool)
trajectory_frames (Optional[list[Molecule]])
trajectory_energies (Optional[list[float]])

vibeqc.molecular_optimize.optimize_molecule(molecule, basis_name, *, method='rhf', functional=None, rhf_options=None, uhf_options=None, rks_options=None, uks_options=None, cisd_options=None, selected_ci_options=None, dmrg_options=None, v2rdm_options=None, transcorrelated_options=None, casci_options=None, caspt2_options=None, nevpt2_options=None, casscf_options=None, active_space=None, cas_reference=None, max_iter=100, conv_tol_grad=0.00045, conv_tol_energy=1e-06, gradient_options=None, grid_options=None, dispersion_params=None, solvent=None, record_trajectory=True, progress=False, fd_step_bohr=0.005, freeze_indices=None)[source]¶

Relax molecular geometry using analytic gradients + L-BFGS-B.

Parameters:

molecule (vibeqc._vibeqc_core.Molecule) – Starting geometry (Cartesian coordinates in bohr).
basis_name (str) – Basis-set name (rebuilt at each geometry step).
method (str) – "rhf", "uhf", "rks", "uks", or a wavefunction method ("selected_ci", "dmrg", "v2rdm", "transcorrelated_ci", "casci", "casscf"). Wavefunction methods fall back to central finite differences on the energy.
functional (str | None) – XC functional string for "rks" / "uks" (e.g. "PBE").
uks_options (vibeqc._vibeqc_core.UKSOptions | None) – Per-method SCF options. If None, defaults are used.
/ (transcorrelated_options / casci_options / caspt2_options)
/
casscf_options (Any) – Wavefunction-solver options, forwarded to every per-step energy evaluation (the FD path) exactly as the final single point receives them – an SA-CASSCF optimization (casscf_options=CASSCFOptions(nroots=2)) walks the state-averaged surface it reports.
active_space (tuple[int, int] | None) – (n_active_orbitals, n_active_electrons) truncation for the wavefunction methods, applied at every per-step evaluation. Without it a selected_ci step would run full-space CI.
cas_reference (str | None) – Reference-orbital choice for the determinant solvers ("rhf" / "uhf" / "uno").
max_iter (int) – Maximum L-BFGS-B iterations.
conv_tol_grad (float) – Gradient convergence tolerance (Ha/bohr). Default 4.5e-4 corresponds to ~0.01 eV/Å – tight enough for routine use.
conv_tol_energy (float) – Energy convergence tolerance (Ha). Controls the scipy ftol parameter.
gradient_options (vibeqc._vibeqc_core.GradientOptions | None) – Options for the analytic gradient kernels (density fitting, COSX, etc.).
grid_options (vibeqc._vibeqc_core.GridOptions | None) – DFT integration grid options (RKS / UKS only).
dispersion_params (Any) – A D3BJParams instance – if provided, the D3-BJ energy and gradient are folded into the objective.
solvent (Any) – A SolventModel or preset string / dict for CPCM implicit solvation (v0.9.0).
record_trajectory (bool) – If True (default), collect per-step geometries and energies for downstream visualisation (QVF animation player).
progress (bool) – If True, print per-step energy and gradient norms to stdout.
fd_step_bohr (float) – Finite-difference step size for wavefunction-method gradients (bohr). Default 0.005 (≈ 0.0026 Å).
freeze_indices (Sequence[int] | None) – Atom indices to hold fixed during the relaxation. Implemented via per-coordinate L-BFGS-B (fixed, fixed) bounds, mirroring vibeqc.bipole_optimize.relax_atoms(). The SCF + gradient still see every atom; the optimizer simply cannot move the frozen ones, and the reported |grad| excludes them so the convergence metric reflects only the free degrees of freedom.
rhf_options (vibeqc._vibeqc_core.RHFOptions | None)
uhf_options (vibeqc._vibeqc_core.UHFOptions | None)
rks_options (vibeqc._vibeqc_core.RKSOptions | None)
uks_options
cisd_options (Any)
selected_ci_options (Any)
dmrg_options (Any)
v2rdm_options (Any)
transcorrelated_options (Any)
casci_options (Any)
caspt2_options (Any)
nevpt2_options (Any)

Return type:

MolecularOptimizeResult

vibeqc.molecular_optimize.optimize_molecule_brent(molecule, basis_name, *, method='rhf', functional=None, rhf_options=None, uhf_options=None, rks_options=None, uks_options=None, cisd_options=None, selected_ci_options=None, dmrg_options=None, v2rdm_options=None, transcorrelated_options=None, casci_options=None, caspt2_options=None, nevpt2_options=None, casscf_options=None, active_space=None, cas_reference=None, max_iter=100, conv_tol_grad=0.00045, gradient_options=None, grid_options=None, dispersion_params=None, solvent=None, record_trajectory=True, progress=False, fd_step_bohr=0.005, freeze_indices=None, line_search_step=0.05, line_search_tol=1e-05)[source]¶

Relax molecular geometry using steepest-descent + Brent line search.

At each geometry step the analytic (or finite-difference) gradient defines the steepest-descent direction. A 1-D line search using Brent’s method finds the optimal step length along that direction.

This is a conservative, gradient-driven optimiser that never takes uphill steps. Use optimizer_backend="brent" in run_job to select it from the top-level API.

Parameters:

molecule (vibeqc._vibeqc_core.Molecule)
basis_name (str)
method (str)
functional (str | None)
rhf_options (vibeqc._vibeqc_core.RHFOptions | None)
uhf_options (vibeqc._vibeqc_core.UHFOptions | None)
rks_options (vibeqc._vibeqc_core.RKSOptions | None)
uks_options (vibeqc._vibeqc_core.UKSOptions | None)
cisd_options (Any)
selected_ci_options (Any)
dmrg_options (Any)
v2rdm_options (Any)
transcorrelated_options (Any)
casci_options (Any)
caspt2_options (Any)
nevpt2_options (Any)
casscf_options (Any)
active_space (tuple[int, int] | None)
cas_reference (str | None)
max_iter (int)
conv_tol_grad (float)
gradient_options (vibeqc._vibeqc_core.GradientOptions | None)
grid_options (vibeqc._vibeqc_core.GridOptions | None)
dispersion_params (Any)
solvent (Any)
record_trajectory (bool)
progress (bool)
fd_step_bohr (float)
freeze_indices (Sequence[int] | None)
line_search_step (float)
line_search_tol (float)

Return type:

MolecularOptimizeResult

vibeqc.molecular_optimize.brent_minimize_1d(f, a, b, c, *, tol=1e-05, max_iter=100, progress=False)[source]¶

Brent’s 1-D minimisation without derivatives.

Finds a local minimum of the scalar function f within the bracketing triplet a < b < c where f(b) < f(a) and f(b) < f(c). The algorithm combines golden-section search with inverse parabolic interpolation.

This is the classic Brent (1973) algorithm as described in Numerical Recipes Sec. 10.2.

Returns (x_min, f_min, n_eval).

Parameters:

a (float)
b (float)
c (float)
tol (float)
max_iter (int)
progress (bool)

Return type:

tuple[float, float, int]

Native Molecular Geometry Optimizer¶

Quick start¶

Through run_job¶

Backends¶

The native optimizer family (geom_opt=)¶

Coordinate systems (geom_coords=)¶

Supported methods¶

Trajectory collection¶

Convergence control¶

API reference¶

Through `run_job`¶

The native optimizer family (`geom_opt=`)¶

Coordinate systems (`geom_coords=`)¶