compatiability with pyscf 2.7.0 (#216)

* compatiability with pyscf 2.7.0

* relax test threshold

* relax hessian threshold

Author: wxj6000

gpu4pyscf/hessian/rhf.py

@@ -67,7 +67,7 @@ def hess_elec(hessobj, mo_energy=None, mo_coeff=None, mo_occ=None,
                                     max_memory, log)
     t1 = log.timer_debug1('hess elec', *t1)
     if h1mo is None:
-        h1mo = hessobj.make_h1(mo_coeff, mo_occ, hessobj.chkfile, atmlst, log)
+        h1mo = hessobj.make_h1(mo_coeff, mo_occ, None, atmlst, log)
         t1 = log.timer_debug1('making H1', *t1)
     if mo1 is None or mo_e1 is None:
         mo1, mo_e1 = hessobj.solve_mo1(mo_energy, mo_coeff, mo_occ, h1mo,
@@ -486,8 +486,7 @@ def gen_hop(hobj, mo_energy=None, mo_coeff=None, mo_occ=None, verbose=None):
                                  max_memory, log)
     de2 += hobj.hess_nuc()
 
-    # Compute H1 integrals and store in hobj.chkfile
-    hobj.make_h1(mo_coeff, mo_occ, hobj.chkfile, atmlst, log)
+    h1ao_cache = hobj.make_h1(mo_coeff, mo_occ, None, atmlst, log)
 
     aoslices = mol.aoslice_by_atom()
     s1a = -mol.intor('int1e_ipovlp', comp=3)
@@ -500,8 +499,7 @@ def h_op(x):
         s1ao = 0
         for ia in range(natm):
             shl0, shl1, p0, p1 = aoslices[ia]
-            h1ao_i = lib.chkfile.load(hobj.chkfile, 'scf_f1ao/%d' % ia)
-            h1ao += numpy.einsum('x,xij->ij', x[ia], h1ao_i)
+            h1ao += numpy.einsum('x,xij->ij', x[ia], h1ao_cache[ia])
             s1ao_i = numpy.zeros((3,nao,nao))
             s1ao_i[:,p0:p1] += s1a[:,p0:p1]
             s1ao_i[:,:,p0:p1] += s1a[:,p0:p1].transpose(0,2,1)
@@ -518,8 +516,7 @@ def h_op(x):
 
         for ja in range(natm):
             q0, q1 = aoslices[ja][2:]
-            h1ao = lib.chkfile.load(hobj.chkfile, 'scf_f1ao/%s'%ja)
-            hx[ja] += numpy.einsum('xpq,pq->x', h1ao, dm1) * 4
+            hx[ja] += numpy.einsum('xpq,pq->x', h1ao_cache[ja], dm1) * 4
             hx[ja] -= numpy.einsum('xpq,pq->x', s1a[:,q0:q1], dme1[q0:q1]) * 2
             hx[ja] -= numpy.einsum('xpq,qp->x', s1a[:,q0:q1], dme1[:,q0:q1]) * 2
         return hx.ravel()
@@ -677,7 +674,6 @@ def __init__(self, scf_method):
         self.stdout = scf_method.stdout
         self.mol = scf_method.mol
         self.base = scf_method
-        self.chkfile = None #scf_method.chkfile
         self.max_memory = self.mol.max_memory
         self.atmlst = range(self.mol.natm)
         self.de = numpy.zeros((0,0,3,3))  # (A,B,dR_A,dR_B)

gpu4pyscf/hessian/uhf.py

@@ -68,7 +68,7 @@ def hess_elec(hessobj, mo_energy=None, mo_coeff=None, mo_occ=None,
                                     max_memory, log)
     t1 = log.timer_debug1('hess elec', *t1)
     if h1mo is None:
-        h1mo = hessobj.make_h1(mo_coeff, mo_occ, hessobj.chkfile, atmlst, log)
+        h1mo = hessobj.make_h1(mo_coeff, mo_occ, None, atmlst, log)
         t1 = log.timer_debug1('making H1', *t1)
     if mo1 is None or mo_e1 is None:
         mo1, mo_e1 = hessobj.solve_mo1(mo_energy, mo_coeff, mo_occ, h1mo,
@@ -364,7 +364,7 @@ def _ao2mo(mat, mo, mocc):
         tmp = contract('xij,jo->xio', mat, mocc)
         return contract('xik,ip->xpk', tmp, mo)
     cupy.get_default_memory_pool().free_all_blocks()
-    
+
     avail_mem = get_avail_mem()
     blksize = int(avail_mem*0.4) // (8*3*nao*nao*4) // ALIGNED * ALIGNED
     blksize = min(8, blksize)
@@ -377,7 +377,7 @@ def _ao2mo(mat, mo, mocc):
     for p0,p1 in lib.prange(0,mol.natm,blksize):
         blk = sorted_idx[p0:p1]
         atom_groups.append(blk)
-        
+
     mo1sa = [None] * mol.natm
     mo1sb = [None] * mol.natm
     e1sa = [None] * mol.natm
@@ -404,7 +404,7 @@ def _ao2mo(mat, mo, mocc):
         tol = mf.conv_tol_cpscf
         mo1, e1 = ucphf.solve(fx, mo_energy, mo_occ, h1vo, s1vo,
                               max_cycle=max_cycle, level_shift=level_shift, tol=tol, verbose=verbose)
-        
+
         mo1a = mo1[0].reshape(-1,3,nao,nocca)
         mo1b = mo1[1].reshape(-1,3,nao,noccb)
         e1a = e1[0].reshape(-1,3,nocca,nocca)
@@ -478,8 +478,7 @@ def gen_hop(hobj, mo_energy=None, mo_coeff=None, mo_occ=None, verbose=None):
                                  max_memory, log)
     de2 += hobj.hess_nuc()
 
-    # Compute H1 integrals and store in hobj.chkfile
-    hobj.make_h1(mo_coeff, mo_occ, hobj.chkfile, atmlst, log)
+    h1ao_cache = hobj.make_h1(mo_coeff, mo_occ, None, atmlst, log)
 
     aoslices = mol.aoslice_by_atom()
     s1a = -mol.intor('int1e_ipovlp', comp=3)
@@ -492,8 +491,7 @@ def h_op(x):
         s1ao = 0
         for ia in range(natm):
             shl0, shl1, p0, p1 = aoslices[ia]
-            h1ao_i = lib.chkfile.load(hobj.chkfile, 'scf_f1ao/%d' % ia)
-            h1ao += numpy.einsum('x,xij->ij', x[ia], h1ao_i)
+            h1ao += numpy.einsum('x,xij->ij', x[ia], h1ao_cache[ia])
             s1ao_i = numpy.zeros((3,nao,nao))
             s1ao_i[:,p0:p1] += s1a[:,p0:p1]
             s1ao_i[:,:,p0:p1] += s1a[:,p0:p1].transpose(0,2,1)
@@ -510,8 +508,7 @@ def h_op(x):
 
         for ja in range(natm):
             q0, q1 = aoslices[ja][2:]
-            h1ao = lib.chkfile.load(hobj.chkfile, 'scf_f1ao/%s'%ja)
-            hx[ja] += numpy.einsum('xpq,pq->x', h1ao, dm1) * 4
+            hx[ja] += numpy.einsum('xpq,pq->x', h1ao_cache[ja], dm1) * 4
             hx[ja] -= numpy.einsum('xpq,pq->x', s1a[:,q0:q1], dme1[q0:q1]) * 2
             hx[ja] -= numpy.einsum('xpq,qp->x', s1a[:,q0:q1], dme1[:,q0:q1]) * 2
         return hx.ravel()

gpu4pyscf/lib/tests/test_to_gpu.py

@@ -64,12 +64,12 @@ def test_rhf(self):
     def test_rks(self):
         mf = rks.RKS(mol).to_gpu()
         e_tot = mf.to_gpu().kernel()
-        assert numpy.abs(e_tot - -74.73210527989748) < 1e-7
+        assert numpy.abs(e_tot - -74.73210527989748) < 1e-6
 
         mf = rks.RKS(mol).run()
         gobj = mf.nuc_grad_method().to_gpu()
         g = gobj.kernel()
-        assert numpy.abs(lib.fp(g) - -0.04340162663176693) < 1e-7
+        assert numpy.abs(lib.fp(g) - -0.04340162663176693) < 1e-6
 
         # RKS Hessian it not supported yet
         # mf = rks.RKS(mol).run()
@@ -114,7 +114,7 @@ def test_df_b3lyp(self):
     def test_df_RKS(self):
         mf = rks.RKS(mol, xc='wb97x').density_fit().to_gpu()
         e_tot = mf.to_gpu().kernel()
-        assert numpy.abs(e_tot - -75.30717654021076) < 1e-7
+        assert numpy.abs(e_tot - -75.30717654021076) < 1e-6
 
         mf = rks.RKS(mol, xc='wb97x').density_fit().run()
         gobj = mf.nuc_grad_method().to_gpu()
@@ -125,7 +125,7 @@ def test_df_RKS(self):
         mf.conv_tol_cpscf = 1e-7
         hobj = mf.Hessian().to_gpu()
         h = hobj.kernel()
-        assert numpy.abs(lib.fp(h) - 2.187025544697092) < 1e-4
+        assert numpy.abs(lib.fp(h) - 2.1858589608638384) < 1e-4
 
 if __name__ == "__main__":
     print("Full tests for to_gpu module")

gpu4pyscf/scf/uhf.py

@@ -20,6 +20,8 @@
 import cupy
 from pyscf.scf import uhf
 from pyscf import lib as pyscf_lib
+from pyscf import __config__
+
 from gpu4pyscf.scf.hf import eigh, damping, level_shift
 from gpu4pyscf.scf import hf
 from gpu4pyscf.lib import logger
@@ -151,7 +153,9 @@ def energy_elec(mf, dm=None, h1e=None, vhf=None):
 class UHF(hf.SCF):
     from gpu4pyscf.lib.utils import to_gpu, device
 
-    _keys = {'e_disp', 'screen_tol', 'conv_tol_cpscf', 'h1e', 's1e'}
+    init_guess_breaksym = getattr(__config__, 'scf_uhf_init_guess_breaksym', 1)
+
+    _keys = {'e_disp', 'screen_tol', 'conv_tol_cpscf', 'h1e', 's1e', "init_guess_breaksym"}
     def __init__(self, mol):
         hf.SCF.__init__(self, mol)
         self.nelec = None

gpu4pyscf/tests/test_dft.py

@@ -56,7 +56,7 @@ def test_b3lyp_with_d3bj():
     mf.conv_tol_cpscf = 1e-8
     mf.disp = 'd3bj'
     e_dft = mf.kernel()
-    assert np.abs(e_dft - -685.0326965348272) < 1e-8
+    assert np.abs(e_dft - -685.0326965348272) < 1e-7
 
     g = mf.nuc_grad_method().kernel()
     assert np.abs(cupy.linalg.norm(g) - 0.17498362161082373) < 1e-5
@@ -73,7 +73,7 @@ def test_b3lyp_d3bj():
     mf.conv_tol = 1e-10
     mf.conv_tol_cpscf = 1e-8
     e_dft = mf.kernel()
-    assert np.abs(e_dft - -685.0326965348272) < 1e-8
+    assert np.abs(e_dft - -685.0326965348272) < 1e-7
 
     g = mf.nuc_grad_method().kernel()
     assert np.abs(cupy.linalg.norm(g) - 0.17498362161082373) < 1e-5
@@ -91,7 +91,7 @@ def test_DFUKS():
     mf.conv_tol_cpscf = 1e-8
     mf.disp = 'd3bj'
     e_dft = mf.kernel()
-    assert np.abs(e_dft - -685.0326965349493) < 1e-8
+    assert np.abs(e_dft - -685.0326965349493) < 1e-7
 
     g = mf.nuc_grad_method().kernel()
     assert np.abs(cupy.linalg.norm(g) - 0.17498264516108836) < 1e-5
@@ -108,7 +108,7 @@ def test_RKS():
     mf.conv_tol = 1e-12
     mf.disp = 'd3bj'
     e_dft = mf.kernel()
-    assert np.abs(e_dft - -685.0325611822375) < 1e-8
+    assert np.abs(e_dft - -685.0325611822375) < 1e-7
 
     g = mf.nuc_grad_method().kernel()
     assert np.abs(cupy.linalg.norm(g) - 0.1750368231223345) < 1e-6
@@ -124,7 +124,7 @@ def test_DFRKS_with_SMD():
     mf.conv_tol_cpscf = 1e-8
     mf.disp = 'd3bj'
     e_dft = mf.kernel()
-    assert np.abs(e_dft - -685.0578838805443) < 1e-8
+    assert np.abs(e_dft - -685.0578838805443) < 1e-7
 
     g = mf.nuc_grad_method().kernel()
     assert np.abs(cupy.linalg.norm(g) - 0.16804945458657145) < 1e-5
@@ -143,7 +143,7 @@ def test_DFUKS_with_SMD():
     mf.conv_tol_cpscf = 1e-8
     mf.disp = 'd3bj'
     e_dft = mf.kernel()
-    assert np.abs(e_dft - -685.05788388063) < 1e-8
+    assert np.abs(e_dft - -685.05788388063) < 1e-7
 
     g = mf.nuc_grad_method().kernel()
     assert np.abs(cupy.linalg.norm(g) - 0.1680496465773684) < 1e-5

setup.py

@@ -138,7 +138,7 @@ def initialize_with_default_plat_name(self):
     ],
     cmdclass={'build_py': CMakeBuildPy},
     install_requires=[
-        'pyscf~=2.6.0',
+        'pyscf~=2.7.0',
         'pyscf-dispersion',
         f'cupy-cuda{CUDA_VERSION}',
         'geometric',