refactor: clean up and try CI build?

pyxdh/DerivOnce/deriv_once_df.py

@@ -1,20 +1,18 @@
+# basic utilities
 import numpy as np
-import opt_einsum
+from opt_einsum import contract as einsum
+from scipy.linalg import solve_triangular
+# python utilities
 from abc import ABC
 from functools import partial
-import os
-from scipy.linalg import solve_triangular
-
+# pyscf utilities
 from pyscf import gto
 from pyscf.df.grad.rhf import _int3c_wrapper as int3c_wrapper
-
+# pyxdh utilities
 from pyxdh.DerivOnce import DerivOnceSCF, DerivOnceMP2
-
-MAXMEM = float(os.getenv("MAXMEM", 2))
+from pyxdh.Utilities import cached_property
+# simplification
 st = partial(solve_triangular, lower=True)
-np.einsum = partial(np.einsum, optimize=["greedy", 1024 ** 3 * MAXMEM / 8])
-einsum = opt_einsum.contract
-np.set_printoptions(8, linewidth=1000, suppress=True)
 
 
 class DerivOnceDFSCF(DerivOnceSCF, ABC):
@@ -23,10 +21,12 @@ def __init__(self, config):
         super(DerivOnceDFSCF, self).__init__(config)
         self.aux_jk = self.scf_eng.with_df.auxmol  # type: gto.Mole
 
-    def _get_eri0_ao(self):
+    @cached_property
+    def eri0_ao(self):
         raise AssertionError("eri0 should not be called in density fitting module!")
 
-    def _get_eri1_ao(self):
+    @cached_property
+    def eri1_ao(self):
         raise AssertionError("eri1 should not be called in density fitting module!")
 
     @staticmethod
@@ -50,41 +50,30 @@ def __init__(self, config):
         self._Y_ao_ri = NotImplemented  # type: np.ndarray
         self._Y_ia_ri = NotImplemented  # type: np.ndarray
 
-    @property
+    @cached_property
     def int2c2e_ri(self):
-        if self._int2c2e_ri is NotImplemented:
-            return self.aux_ri.intor("int2c2e")
-        return self._int2c2e_ri
+        return self.aux_ri.intor("int2c2e")
 
-    @property
+    @cached_property
     def int3c2e_ri(self):
-        if self._int3c2e_ri is NotImplemented:
-            return int3c_wrapper(self.mol, self.aux_ri, "int3c2e", "s1")()
-        return self._int3c2e_ri
+        return int3c_wrapper(self.mol, self.aux_ri, "int3c2e", "s1")()
 
-    @property
+    @cached_property
     def L_ri(self):
-        if self._L_ri is NotImplemented:
-            return np.linalg.cholesky(self.int2c2e_ri)
-        return self._L_ri
+        return np.linalg.cholesky(self.int2c2e_ri)
 
-    @property
+    @cached_property
     def L_inv_ri(self):
-        if self._L_inv_ri is NotImplemented:
-            return np.linalg.inv(self.L_ri)
-        return self._L_inv_ri
+        return np.linalg.inv(self.L_ri)
 
-    @property
+    @cached_property
     def Y_ao_ri(self):
-        if self._Y_ao_ri is NotImplemented:
-            return einsum("μνQ, PQ -> μνP", self.int3c2e_ri, self.L_inv_ri)
-        return self._Y_ao_ri
+        return einsum("μνQ, PQ -> μνP", self.int3c2e_ri, self.L_inv_ri)
 
-    @property
+    @cached_property
     def Y_ia_ri(self):
-        if self._Y_ia_ri is NotImplemented:
-            return einsum("μνP, μi, νa -> iaP", self.Y_ao_ri, self.Co, self.Cv)
-        return self._Y_ia_ri
+        return einsum("μνP, μi, νa -> iaP", self.Y_ao_ri, self.Co, self.Cv)
 
-    def _get_t_iajb(self):
+    @cached_property
+    def t_iajb(self):
         return einsum("iaP, jbP -> iajb", self.Y_ia_ri, self.Y_ia_ri) / self.D_iajb

pyxdh/DerivTwice/__init__.py

@@ -5,8 +5,7 @@
     "DipDerivSCF", "DipDerivNCDFT", "DipDerivMP2", "DipDerivXDH",
     "DerivTwiceUSCF", "DerivTwiceUMP2",
     "HessUSCF", "HessUMP2",
-    "PolarUSCF",
-    "PolarUMP2",
+    "PolarUSCF", "PolarUMP2",
 ]
 
 from pyxdh.DerivTwice.deriv_twice_r import DerivTwiceSCF, DerivTwiceNCDFT, DerivTwiceMP2, DerivTwiceXDH
@@ -15,6 +14,5 @@
 from pyxdh.DerivTwice.dipderiv_r import DipDerivSCF, DipDerivNCDFT, DipDerivMP2, DipDerivXDH
 
 from pyxdh.DerivTwice.deriv_twice_u import DerivTwiceUSCF, DerivTwiceUMP2
-from pyxdh.DerivTwice.hess_uscf import HessUSCF, HessUMP2
-from pyxdh.DerivTwice.polar_uscf import PolarUSCF
-from pyxdh.DerivTwice.polar_ump2 import PolarUMP2
+from pyxdh.DerivTwice.hess_u import HessUSCF, HessUMP2
+from pyxdh.DerivTwice.polar_u import PolarUSCF, PolarUMP2

pyxdh/DerivTwice/hess_uscf.py → pyxdh/DerivTwice/hess_u.py

@@ -6,7 +6,7 @@
 from pyxdh.Utilities import timing, cached_property
 # pytest
 from pyscf import gto, scf
-from pyxdh.DerivOnce import GradUSCF, GradUMP2
+from pyxdh.DerivOnce import GradUSCF
 from pyxdh.Utilities import FormchkInterface
 from pkg_resources import resource_filename
 

pyxdh/DerivTwice/polar_u.py

@@ -0,0 +1,48 @@
+# basic utilities
+import numpy as np
+from opt_einsum import contract as einsum
+# pyxdh utilities
+from pyxdh.DerivTwice import DerivTwiceUSCF, PolarSCF, DerivTwiceUMP2
+# pytest
+from pyscf import gto, scf
+from pyxdh.DerivOnce import DipoleUSCF, DipoleUMP2
+from pyxdh.Utilities import FormchkInterface
+from pkg_resources import resource_filename
+
+
+class PolarUSCF(DerivTwiceUSCF, PolarSCF):
+    pass
+
+
+class PolarUMP2(DerivTwiceUMP2, PolarUSCF):
+
+    def _get_E_2_MP2_Contrib(self):
+        A, B = self.A, self.B
+        so, sv = self.so, self.sv
+
+        E_2_MP2_Contrib = (
+            + einsum("xBpq, xApq -> AB", B.pdA_D_r_oovv, A.B_1)
+            + einsum("Aai, Bai -> AB", A.U_1[0, :, sv[0], so[0]], self.RHS_B[0])
+            + einsum("Aai, Bai -> AB", A.U_1[1, :, sv[1], so[1]], self.RHS_B[1])
+            + einsum("xpq, xABpq -> AB", self.D_r, self.pdB_F_A_mo)
+        )
+        return E_2_MP2_Contrib
+
+
+class TestPolarU:
+
+    mol = gto.Mole(atom="C 0. 0. 0.; H 1. 0. 0.; H 0. 2. 0.; H 0. 0. 1.5", basis="6-31G", spin=1, verbose=0).build()
+
+    def test_u_uhf_polar(self):
+        scf_eng = scf.UHF(self.mol); scf_eng.conv_tol = 1e-12; scf_eng.conv_tol_grad = 1e-10; scf_eng.max_cycle = 256; scf_eng.run()
+        diph = DipoleUSCF({"scf_eng": scf_eng})
+        polh = PolarUSCF({"deriv_A": diph})
+        formchk = FormchkInterface(resource_filename("pyxdh", "Validation/gaussian/CH3-HF-freq.fchk"))
+        assert np.allclose(- polh.E_2, formchk.polarizability(), atol=1e-6, rtol=1e-4)
+
+    def test_u_mp2_polar(self):
+        scf_eng = scf.UHF(self.mol); scf_eng.conv_tol = 1e-12; scf_eng.conv_tol_grad = 1e-10; scf_eng.max_cycle = 256; scf_eng.run()
+        diph = DipoleUMP2({"scf_eng": scf_eng, "cphf_tol": 1e-12})
+        polh = PolarUMP2({"deriv_A": diph})
+        formchk = FormchkInterface(resource_filename("pyxdh", "Validation/gaussian/CH3-MP2-freq.fchk"))
+        assert np.allclose(- polh.E_2, formchk.polarizability(), atol=1e-6, rtol=1e-4)