From 60c1eae2e093e08246579bd65cd7b13c173080bc Mon Sep 17 00:00:00 2001
From: "Josef M. Gallmetzer" <64498081+galjos@users.noreply.github.com>
Date: Tue, 7 Jul 2026 14:13:42 +0200
Subject: [PATCH 1/3] Add in-memory thermochemistry from PySCF results
Add pyscf_thermo(mean_field, ...) which reads the geometry (Bohr ->
Angstrom) and energy (mf.e_tot) from a converged PySCF mean-field object
and runs the RRHO thermochemistry, taking vibrational frequencies directly
or deriving them from a Hessian via pyscf.hessian.thermo. Optional pyscf
extra; pyscf is imported lazily so the core install is unaffected.
Closes #102
---
ThermoScreening/thermo/api.py | 98 ++++++++++++++++++++++++++++++++++
docs/api.rst | 1 +
docs/usage.rst | 14 +++++
pyproject.toml | 4 ++
tests/calculator/test_pyscf.py | 60 +++++++++++++++++++++
5 files changed, 177 insertions(+)
create mode 100644 tests/calculator/test_pyscf.py
diff --git a/ThermoScreening/thermo/api.py b/ThermoScreening/thermo/api.py
index aefdfc6..76d2785 100644
--- a/ThermoScreening/thermo/api.py
+++ b/ThermoScreening/thermo/api.py
@@ -568,6 +568,104 @@ def cclib_thermo(
)
+def _pyscf_frequencies_from_hessian(mol, hessian):
+ """Vibrational frequencies (cm^-1) from a PySCF Hessian, via pyscf."""
+ try:
+ from pyscf.hessian import thermo as pyscf_thermo_module
+ except ImportError as exc:
+ raise TSValueError(
+ "pyscf is required to compute frequencies from a Hessian; install it "
+ "with 'pip install thermoscreening[pyscf]'."
+ ) from exc
+ freq_info = pyscf_thermo_module.harmonic_analysis(mol, hessian) # pragma: no cover
+ return np.real(freq_info["freq_wavenumber"]) # pragma: no cover
+
+
+def pyscf_thermo(
+ mean_field,
+ frequencies=None,
+ hessian=None,
+ energy=None,
+ temperature=298.15,
+ pressure=101325,
+ charge=0.0,
+ spin=None,
+ quasi_rrho=False,
+):
+ """
+ Run the thermochemistry from an in-memory PySCF calculation.
+
+ PySCF is a Python library, so its results live in memory rather than a file.
+ Pass a converged mean-field object; its geometry (Bohr -> Angstrom) and
+ energy (``mean_field.e_tot``) are used, and the vibrational frequencies are
+ taken from ``frequencies`` or computed from ``hessian`` via
+ ``pyscf.hessian.thermo``.
+
+ Parameters
+ ----------
+ mean_field : object
+ A converged PySCF mean-field object (has ``.mol`` and ``.e_tot``).
+ frequencies : array_like, optional
+ Vibrational frequencies in cm^-1. If omitted, ``hessian`` is used.
+ hessian : array_like, optional
+ A PySCF Hessian (from ``mf.Hessian().kernel()``); frequencies are derived
+ from it. Requires the ``pyscf`` extra. Ignored if ``frequencies`` is given.
+ energy : float, optional
+ Electronic energy in Hartree. Defaults to ``mean_field.e_tot``.
+ temperature : float
+ Temperature in K. Default 298.15.
+ pressure : float
+ Pressure in Pa. Default 101325.
+ charge : float
+ System charge. Default 0.0.
+ spin : float, optional
+ Spin quantum number S. Defaults to the minimum-spin electron-count guess.
+ quasi_rrho : bool
+ If True, use Grimme's quasi-RRHO vibrational entropy. Default False.
+
+ Returns
+ -------
+ Thermo
+ The thermo calculation object.
+
+ Raises
+ ------
+ TSValueError
+ If neither ``frequencies`` nor ``hessian`` is given (or pyscf is missing
+ when a ``hessian`` is used).
+ """
+ from ase import Atoms
+ from ase.units import Bohr
+
+ mol = mean_field.mol
+ symbols = [mol.atom_pure_symbol(index) for index in range(mol.natm)]
+ positions = np.asarray(mol.atom_coords()) * Bohr # Bohr -> Angstrom
+ atoms = Atoms(symbols=symbols, positions=positions)
+
+ if energy is None:
+ energy = float(mean_field.e_tot)
+
+ if frequencies is None:
+ if hessian is None:
+ raise TSValueError(
+ "Provide frequencies=... or hessian=... to pyscf_thermo."
+ )
+ frequencies = _pyscf_frequencies_from_hessian(mol, hessian)
+ frequencies = np.asarray(frequencies, dtype=float)
+
+ return run_thermo(
+ vibrational_frequencies=frequencies,
+ atoms=atoms,
+ energy=energy,
+ temperature=temperature,
+ pressure=pressure,
+ charge=charge,
+ spin=spin,
+ engine="dftb+",
+ quasi_rrho=quasi_rrho,
+ )
+
+
def execute(input_file: str) -> Thermo:
"""
Execute the thermo calculation. Returns thermo calculation object.
diff --git a/docs/api.rst b/docs/api.rst
index d19c187..e266f58 100644
--- a/docs/api.rst
+++ b/docs/api.rst
@@ -9,6 +9,7 @@ High-level pipeline
.. autofunction:: ThermoScreening.thermo.api.xtb_cli_thermo
.. autofunction:: ThermoScreening.thermo.api.orca_thermo
.. autofunction:: ThermoScreening.thermo.api.cclib_thermo
+.. autofunction:: ThermoScreening.thermo.api.pyscf_thermo
.. autofunction:: ThermoScreening.thermo.api.run_thermo
.. autofunction:: ThermoScreening.thermo.api.execute
diff --git a/docs/usage.rst b/docs/usage.rst
index dfb034c..74b24b5 100644
--- a/docs/usage.rst
+++ b/docs/usage.rst
@@ -135,6 +135,20 @@ which auto-detects the program via `cclib `_:
``read_cclib`` returns the parsed ``(atoms, frequencies, energy)`` if you want
them directly; ``energy=`` overrides the parsed energy.
+**PySCF** results live in memory, so pass the mean-field object directly (with
+the ``pyscf`` extra, ``pip install thermoscreening[pyscf]``):
+
+.. code-block:: python
+
+ from ThermoScreening.thermo.api import pyscf_thermo
+
+ mf = mol.RKS(xc="b3lyp").run()
+ hessian = mf.Hessian().kernel()
+ thermo = pyscf_thermo(mf, hessian=hessian) # or frequencies=
+
+The geometry and energy (``mf.e_tot``) come from the mean field; frequencies are
+taken from ``frequencies=`` or derived from ``hessian=`` via ``pyscf.hessian``.
+
Temperature scans
-----------------
diff --git a/pyproject.toml b/pyproject.toml
index da94e9e..bfd1cf1 100644
--- a/pyproject.toml
+++ b/pyproject.toml
@@ -49,6 +49,10 @@ Issues = "https://github.com/MolarVerse/ThermoScreening/issues"
qm = [
"cclib >= 1.8, < 2.0",
]
+# In-memory thermochemistry from PySCF results
+pyscf = [
+ "pyscf",
+]
test = [
"pytest",
"coverage",
diff --git a/tests/calculator/test_pyscf.py b/tests/calculator/test_pyscf.py
new file mode 100644
index 0000000..27a330d
--- /dev/null
+++ b/tests/calculator/test_pyscf.py
@@ -0,0 +1,60 @@
+import math
+import sys
+
+import numpy as np
+import pytest
+
+import ThermoScreening.thermo.api as api
+from ThermoScreening.thermo.api import pyscf_thermo, _pyscf_frequencies_from_hessian
+from ThermoScreening.exceptions import TSValueError
+
+
+class _FakeMol:
+ natm = 3
+
+ def atom_coords(self): # Bohr
+ return np.array([[0.0, 0.0, 0.221722],
+ [0.0, 1.430901, -0.886569],
+ [0.0, -1.430901, -0.886569]])
+
+ def atom_pure_symbol(self, index):
+ return ["O", "H", "H"][index]
+
+
+class _FakeMeanField:
+ def __init__(self, e_tot=-76.4):
+ self.mol = _FakeMol()
+ self.e_tot = e_tot
+
+
+def test_pyscf_thermo_from_frequencies():
+ thermo = pyscf_thermo(_FakeMeanField(), frequencies=[1600.0, 3700.0, 3800.0])
+ assert thermo.electronic_energy() == pytest.approx(-76.4) # from mf.e_tot
+ assert math.isfinite(thermo.total_EeGtot())
+
+
+def test_pyscf_thermo_energy_override():
+ thermo = pyscf_thermo(_FakeMeanField(), frequencies=[1600.0, 3700.0, 3800.0], energy=-77.0)
+ assert thermo.electronic_energy() == pytest.approx(-77.0)
+
+
+def test_pyscf_thermo_requires_frequencies_or_hessian():
+ with pytest.raises(TSValueError, match="Provide frequencies"):
+ pyscf_thermo(_FakeMeanField())
+
+
+def test_pyscf_thermo_hessian_path(monkeypatch):
+ monkeypatch.setattr(
+ api, "_pyscf_frequencies_from_hessian",
+ lambda mol, hessian: np.array([1600.0, 3700.0, 3800.0]),
+ )
+ thermo = pyscf_thermo(_FakeMeanField(), hessian=object())
+ assert math.isfinite(thermo.total_EeGtot())
+
+
+def test_pyscf_frequencies_from_hessian_missing_dependency(monkeypatch):
+ # force the pyscf import inside the helper to fail, regardless of whether
+ # pyscf happens to be installed in the test environment (None -> ImportError)
+ monkeypatch.setitem(sys.modules, "pyscf.hessian.thermo", None)
+ with pytest.raises(TSValueError, match="pyscf is required"):
+ _pyscf_frequencies_from_hessian(object(), object())
From f0c5497e7fdfb3e59ede7d49ab8553144408b1b6 Mon Sep 17 00:00:00 2001
From: "Josef M. Gallmetzer" <64498081+galjos@users.noreply.github.com>
Date: Tue, 7 Jul 2026 14:16:05 +0200
Subject: [PATCH 2/3] Keep PySCF imaginary modes as negative reals; real
integration test
harmonic_analysis defaults to imaginary_freq=True (complex wavenumbers),
which np.real would flatten to ~0 and silently turn a saddle point into a
minimum. Pass imaginary_freq=False so imaginary modes are negative reals,
matching the other engines. Add a skippable real-pyscf test validating the
harmonic_analysis path returns the 3N-6 real vibrational modes for water.
---
ThermoScreening/thermo/api.py | 9 +++++++--
tests/calculator/test_pyscf.py | 19 +++++++++++++++++++
2 files changed, 26 insertions(+), 2 deletions(-)
diff --git a/ThermoScreening/thermo/api.py b/ThermoScreening/thermo/api.py
index 76d2785..79518f9 100644
--- a/ThermoScreening/thermo/api.py
+++ b/ThermoScreening/thermo/api.py
@@ -577,8 +577,13 @@ def _pyscf_frequencies_from_hessian(mol, hessian):
"pyscf is required to compute frequencies from a Hessian; install it "
"with 'pip install thermoscreening[pyscf]'."
) from exc
- freq_info = pyscf_thermo_module.harmonic_analysis(mol, hessian) # pragma: no cover
- return np.real(freq_info["freq_wavenumber"]) # pragma: no cover
+ # imaginary_freq=False stores imaginary modes as negative real wavenumbers
+ # (the convention the other engines use), so a saddle point isn't silently
+ # flattened into a minimum by dropping a complex frequency.
+ freq_info = pyscf_thermo_module.harmonic_analysis( # pragma: no cover
+ mol, hessian, imaginary_freq=False
+ )
+ return np.asarray(freq_info["freq_wavenumber"], dtype=float) # pragma: no cover
def pyscf_thermo(
diff --git a/tests/calculator/test_pyscf.py b/tests/calculator/test_pyscf.py
index 27a330d..aad7398 100644
--- a/tests/calculator/test_pyscf.py
+++ b/tests/calculator/test_pyscf.py
@@ -1,9 +1,12 @@
+import importlib.util
import math
import sys
import numpy as np
import pytest
+_HAS_PYSCF = importlib.util.find_spec("pyscf") is not None
+
import ThermoScreening.thermo.api as api
from ThermoScreening.thermo.api import pyscf_thermo, _pyscf_frequencies_from_hessian
from ThermoScreening.exceptions import TSValueError
@@ -58,3 +61,19 @@ def test_pyscf_frequencies_from_hessian_missing_dependency(monkeypatch):
monkeypatch.setitem(sys.modules, "pyscf.hessian.thermo", None)
with pytest.raises(TSValueError, match="pyscf is required"):
_pyscf_frequencies_from_hessian(object(), object())
+
+
+@pytest.mark.skipif(not _HAS_PYSCF, reason="pyscf is not installed")
+def test_pyscf_frequencies_from_hessian_real():
+ # real pyscf: harmonic_analysis on a water HF/STO-3G Hessian returns the
+ # 3N-6 vibrational modes as a real (cm^-1) array (imaginary_freq=False)
+ from pyscf import gto, scf
+
+ mol = gto.M(atom="O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469",
+ basis="sto-3g", verbose=0)
+ mf = scf.RHF(mol).run()
+ hessian = mf.Hessian().kernel()
+
+ frequencies = _pyscf_frequencies_from_hessian(mol, hessian)
+ assert frequencies.dtype == float # real (negative for any imaginary mode)
+ assert len(frequencies) == 3 # 3N - 6 for a bent triatomic
From bc26455e37411de5b0e8436ef957d3412189e19c Mon Sep 17 00:00:00 2001
From: "Josef M. Gallmetzer" <64498081+galjos@users.noreply.github.com>
Date: Tue, 7 Jul 2026 14:23:08 +0200
Subject: [PATCH 3/3] Make the pyscf missing-dependency test order-robust
Intercept the import itself instead of a sys.modules sentinel, which a
prior import of pyscf.hessian.thermo (e.g. the real integration test) would
defeat via the fromlist hasattr shortcut.
---
tests/calculator/test_pyscf.py | 17 +++++++++++++----
1 file changed, 13 insertions(+), 4 deletions(-)
diff --git a/tests/calculator/test_pyscf.py b/tests/calculator/test_pyscf.py
index aad7398..674c4d2 100644
--- a/tests/calculator/test_pyscf.py
+++ b/tests/calculator/test_pyscf.py
@@ -1,6 +1,5 @@
import importlib.util
import math
-import sys
import numpy as np
import pytest
@@ -56,9 +55,19 @@ def test_pyscf_thermo_hessian_path(monkeypatch):
def test_pyscf_frequencies_from_hessian_missing_dependency(monkeypatch):
- # force the pyscf import inside the helper to fail, regardless of whether
- # pyscf happens to be installed in the test environment (None -> ImportError)
- monkeypatch.setitem(sys.modules, "pyscf.hessian.thermo", None)
+ # force the pyscf import inside the helper to fail, robustly -- independent
+ # of whether pyscf is installed or already imported (which would defeat a
+ # sys.modules sentinel) -- by intercepting the import itself
+ import builtins
+
+ real_import = builtins.__import__
+
+ def fake_import(name, *args, **kwargs):
+ if name.startswith("pyscf"):
+ raise ImportError("no pyscf")
+ return real_import(name, *args, **kwargs)
+
+ monkeypatch.setattr(builtins, "__import__", fake_import)
with pytest.raises(TSValueError, match="pyscf is required"):
_pyscf_frequencies_from_hessian(object(), object())