Fix periodic PySCF calculations (dftd4#254)

* fix periodic pyscf calculations; ported fix from dftd3/simple-dftd3#74 * add test for pbc --------- Co-authored-by: Alex Izvorski <[email protected]>
grimme-lab · Aug 7, 2024 · 97b23b4 · 97b23b4
1 parent 4aa9e79
commit 97b23b4
Show file tree

Hide file tree

Showing 2 changed files with 38 additions and 1 deletion.
diff --git a/python/dftd4/pyscf.py b/python/dftd4/pyscf.py
@@ -100,10 +100,18 @@ def kernel(self) -> Tuple[float, np.ndarray]:
         """
         mol = self.mol
 
+        lattice = None
+        periodic = None
+        if hasattr(mol, 'lattice_vectors'):
+            lattice = mol.lattice_vectors()
+            periodic = np.array([True, True, True], dtype=bool)
+
         disp = DispersionModel(
             np.asarray([gto.charge(sym) for sym in mol.elements]),
             mol.atom_coords(),
             mol.charge,
+            lattice=lattice,
+            periodic=periodic,
         )
 
         param = DampingParam(

diff --git a/python/dftd4/test_pyscf.py b/python/dftd4/test_pyscf.py
@@ -23,7 +23,7 @@
 try:
     import dftd4.pyscf as disp
     import pyscf
-    from pyscf import gto, scf
+    from pyscf import gto, scf, pbc
 except ModuleNotFoundError:
     pyscf = None
 
@@ -151,3 +151,32 @@ def test_gradient_hf() -> None:
     )
     grad = disp.energy(scf.RHF(mol)).run().nuc_grad_method()
     assert grad.kernel() == approx(ref, abs=1.0e-7)
+
+
+@pytest.mark.skipif(pyscf is None, reason="requires pyscf")
+def test_pbc():
+    mol = gto.M(
+        atom="""
+             O  -1.6256  -0.0413   0.3705
+             H  -0.7061  -0.0938   0.0934
+             H  -2.0618  -0.7328  -0.1359
+             """,
+        basis="def2-tzvp",
+    )
+
+    pmol = pbc.gto.M(
+        atom="""
+             O  -1.6256  -0.0413   0.3705
+             H  -0.7061  -0.0938   0.0934
+             H  -2.0618  -0.7328  -0.1359
+             """,
+        basis="def2-tzvp",
+        a=[[3, 0, 0], [0, 3, 0], [0, 0, 3]],
+    )
+
+    xc = 'pbe'
+
+    e_mol_disp = disp.DFTD4Dispersion(mol, xc=xc).kernel()[0]
+    e_pbc_disp = disp.DFTD4Dispersion(pmol, xc=xc).kernel()[0]
+
+    assert e_mol_disp != e_pbc_disp