diff --git a/gpu4pyscf/_patch_pyscf.py b/gpu4pyscf/_patch_pyscf.py
index bb4206553..068a2ba20 100644
--- a/gpu4pyscf/_patch_pyscf.py
+++ b/gpu4pyscf/_patch_pyscf.py
@@ -16,7 +16,9 @@
 import cupy
 import pyscf
 
-if int(pyscf.__version__.split('.')[1]) <= 10:
+pyscf_version = int(pyscf.__version__.split('.')[1])
+
+if pyscf_version <= 10:
     def _fftdf_to_gpu(self):
         from gpu4pyscf.pbc.df.fft import FFTDF
         return FFTDF(self.cell, self.kpts)
@@ -35,26 +37,6 @@ def _gdf_to_gpu(self):
     from pyscf.pbc.df.df import GDF
     GDF.to_gpu = _gdf_to_gpu
 
-    from pyscf.solvent.grad import pcm as pcm_grad
-    if hasattr(pcm_grad, 'WithSolventGrad'):
-        def _pcm_grad_to_gpu(self):
-            from pyscf.tdscf.rhf import TDBase
-            from pyscf.solvent.pcm import PCM
-            assert isinstance(self.base.with_solvent, PCM)
-            if isinstance(self, TDBase):
-                raise NotImplementedError('.to_gpu() for PCM-TDDFT')
-            return self.base.to_gpu().PCM().Gradients()
-        pcm_grad.WithSolventGrad.to_gpu = _pcm_grad_to_gpu
-
-    from pyscf.solvent.hessian import pcm as pcm_hess
-    if hasattr(pcm_hess, 'WithSolventHess'):
-        def _pcm_hessian_to_gpu(self):
-            from pyscf.tdscf.rhf import TDBase
-            if isinstance(self, TDBase):
-                raise NotImplementedError('.to_gpu() for PCM-TDDFT')
-            return self.base.to_gpu().PCM().Hessian()
-        pcm_hess.WithSolventHess.to_gpu = _pcm_hessian_to_gpu
-
     # patch PySCF Cell class, updating lattice parameters is not avail in pyscf 2.10
     from pyscf import lib
     from pyscf.lib import logger
@@ -182,3 +164,93 @@ def find_boundary(a):
         return np.asarray(Ls, order='C')
     # Patch the get_lattice_Ls for pyscf-2.9 or older
     Cell.get_lattice_Ls = get_lattice_Ls
+
+if pyscf_version <= 11:
+    # In pyscf-2.11, the auxbasis_response attribute is not registered in the
+    # df.Hessian._keys. Consequently, this key is excluded by the conversion in
+    # utils.to_cpu()
+    from pyscf.df.hessian import rhf, rks, uhf, uks
+    rhf.Hessian._keys = uhf.Hessian._keys = \
+            rks.Hessian._keys = uks.Hessian._keys = {'auxbasis_response',}
+
+    from pyscf.lib import misc
+    misc._ATTRIBUTES_IN_NPARRAY = {
+        'kpt', 'kpts', '_kpts', 'kpts_band', 'mesh', 'frozen'}
+    def to_gpu(method, out=None):
+        '''Convert a method to its corresponding GPU variant, and recursively
+        converts all attributes of a method to cupy objects or gpu4pyscf objects.
+        '''
+        # If a GPU class inherits a CPU code, the "to_gpu" method may be resolved
+        # and available in the GPU class. Skip the conversion in this case.
+        if method.__module__.startswith('gpu4pyscf'):
+            return method
+
+        if out is None:
+            if isinstance(method, (misc.SinglePointScanner, misc.GradScanner)):
+                method = method.undo_scanner()
+
+            from importlib import import_module
+            mod = import_module(method.__module__.replace('pyscf', 'gpu4pyscf'))
+            try:
+                cls = getattr(mod, method.__class__.__name__)
+            except AttributeError:
+                if hasattr(cls, 'from_cpu'):
+                    # the customized to_gpu function can be accessed at module
+                    # levelin gpu4pyscf.
+                    return cls.from_cpu(method)
+                raise
+
+            # Allow gpu4pyscf to customize the to_gpu method for PySCF classes.
+            if hasattr(mod, 'from_cpu'):
+                return mod.from_cpu(method)
+
+            # A temporary GPU instance. This ensures to initialize private
+            # attributes that are only available for GPU code.
+            cls = getattr(mod, method.__class__.__name__)
+            out = method.view(cls)
+
+        elif hasattr(out, 'from_cpu'):
+            out.__dict__.update(out.__class__.from_cpu(method).__dict__)
+            return out
+
+        # Convert only the keys that are defined in the corresponding GPU class
+        cls_keys = [getattr(cls, '_keys', ()) for cls in out.__class__.__mro__[:-1]]
+        out_keys = set(out.__dict__).union(*cls_keys)
+        # Only overwrite the attributes of the same name.
+        keys = out_keys.intersection(method.__dict__)
+
+        for key in keys:
+            val = getattr(method, key)
+            if isinstance(val, np.ndarray):
+                if key not in misc._ATTRIBUTES_IN_NPARRAY:
+                    val = cupy.asarray(val)
+            elif hasattr(val, 'to_gpu'):
+                val = val.to_gpu()
+            setattr(out, key, val)
+        if hasattr(out, 'reset'):
+            try:
+                out.reset()
+            except NotImplementedError:
+                pass
+        return out
+    misc.to_gpu = to_gpu
+
+    from pyscf.solvent.grad import pcm as pcm_grad
+    if hasattr(pcm_grad, 'WithSolventGrad'):
+        def _pcm_grad_to_gpu(self):
+            from pyscf.tdscf.rhf import TDBase
+            from pyscf.solvent.pcm import PCM
+            assert isinstance(self.base.with_solvent, PCM)
+            if isinstance(self, TDBase):
+                raise NotImplementedError('.to_gpu() for PCM-TDDFT')
+            return misc.to_gpu(self, self.base.to_gpu().Gradients())
+        pcm_grad.WithSolventGrad.to_gpu = _pcm_grad_to_gpu
+
+    from pyscf.solvent.hessian import pcm as pcm_hess
+    if hasattr(pcm_hess, 'WithSolventHess'):
+        def _pcm_hessian_to_gpu(self):
+            from pyscf.tdscf.rhf import TDBase
+            if isinstance(self, TDBase):
+                raise NotImplementedError('.to_gpu() for PCM-TDDFT')
+            return misc.to_gpu(self, self.base.to_gpu().Hessian())
+        pcm_hess.WithSolventHess.to_gpu = _pcm_hessian_to_gpu
diff --git a/gpu4pyscf/cc/tests/test_ccsd.py b/gpu4pyscf/cc/tests/test_ccsd.py
index 91ca21d0e..a6a2e570d 100644
--- a/gpu4pyscf/cc/tests/test_ccsd.py
+++ b/gpu4pyscf/cc/tests/test_ccsd.py
@@ -2,7 +2,6 @@
 import numpy as np
 import pyscf
 import pytest
-from packaging import version
 from gpu4pyscf.cc import ccsd_incore
 
 def setUpModule():
@@ -44,7 +43,6 @@ def test_ccsd_incore_update_amps(self):
         self.assertAlmostEqual(abs(r1 - t1).max(), 0, 9)
         self.assertAlmostEqual(abs(r2 - t2).max(), 0, 9)
 
-    @pytest.mark.skipif(version.parse(pyscf.__version__) <= version.parse('2.4.0'), reason='requires pyscf 2.5 or higher')
     def test_ccsd_incore_kernel(self):
         ref = mf.CCSD().run()
         mcc = ccsd_incore.CCSD(mf.to_gpu()).run()
diff --git a/gpu4pyscf/df/df.py b/gpu4pyscf/df/df.py
index 46aa44a90..9946fc8fc 100644
--- a/gpu4pyscf/df/df.py
+++ b/gpu4pyscf/df/df.py
@@ -40,6 +40,8 @@
 
 class DF(lib.StreamObject):
 
+    use_gpu_memory = True
+
     _keys = {'intopt', 'nao', 'naux', 'cd_low', 'mol', 'auxmol', 'use_gpu_memory'}
 
     def __init__(self, mol, auxbasis=None):
@@ -47,7 +49,6 @@ def __init__(self, mol, auxbasis=None):
         self.stdout = mol.stdout
         self.verbose = mol.verbose
         self.max_memory = mol.max_memory
-        self.use_gpu_memory = True
         self._auxbasis = auxbasis
 
         self.auxmol = None
diff --git a/gpu4pyscf/df/hessian/rhf.py b/gpu4pyscf/df/hessian/rhf.py
index 417b63b9c..a7e080215 100644
--- a/gpu4pyscf/df/hessian/rhf.py
+++ b/gpu4pyscf/df/hessian/rhf.py
@@ -645,7 +645,7 @@ class Hessian(rhf_hess.Hessian):
 
     _keys = {'auxbasis_response',}
 
-    auxbasis_response = 1
+    auxbasis_response = 2
     partial_hess_elec = partial_hess_elec
     make_h1 = make_h1
     get_jk_mo = _get_jk_mo
diff --git a/gpu4pyscf/df/hessian/rks.py b/gpu4pyscf/df/hessian/rks.py
index 17e881fe2..a41aa6824 100644
--- a/gpu4pyscf/df/hessian/rks.py
+++ b/gpu4pyscf/df/hessian/rks.py
@@ -107,7 +107,7 @@ class Hessian(rks_hess.Hessian):
 
     _keys = {'auxbasis_response',}
 
-    auxbasis_response = 1
+    auxbasis_response = 2
     partial_hess_elec = partial_hess_elec
     make_h1 = make_h1
     get_jk_mo = df_rhf_hess._get_jk_mo
diff --git a/gpu4pyscf/df/hessian/tests/test_df_rks_hessian.py b/gpu4pyscf/df/hessian/tests/test_df_rks_hessian.py
index f737e92ab..61f9ce2cb 100644
--- a/gpu4pyscf/df/hessian/tests/test_df_rks_hessian.py
+++ b/gpu4pyscf/df/hessian/tests/test_df_rks_hessian.py
@@ -66,6 +66,7 @@ def test_df_lda(self):
 
         mf = mf.to_gpu()
         hessobj = mf.Hessian()
+        hessobj.auxbasis_response = 1
         hess_gpu = hessobj.kernel()
         assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
 
@@ -81,6 +82,7 @@ def test_df_gga(self):
 
         mf = mf.to_gpu()
         hessobj = mf.Hessian()
+        hessobj.auxbasis_response = 1
         hessobj.base.cphf_grids = hessobj.base.grids
         hess_gpu = hessobj.kernel()
         assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
@@ -97,6 +99,7 @@ def test_df_mgga(self):
 
         mf = mf.to_gpu()
         hessobj = mf.Hessian()
+        hessobj.auxbasis_response = 1
         hessobj.base.cphf_grids = hessobj.base.grids
         hess_gpu = hessobj.kernel()
         assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
@@ -104,4 +107,3 @@ def test_df_mgga(self):
 if __name__ == "__main__":
     print("Full Tests for DF RKS Hessian")
     unittest.main()
-    
\ No newline at end of file
diff --git a/gpu4pyscf/df/hessian/tests/test_df_uks_hessian.py b/gpu4pyscf/df/hessian/tests/test_df_uks_hessian.py
index 31804715d..8ba5180c1 100644
--- a/gpu4pyscf/df/hessian/tests/test_df_uks_hessian.py
+++ b/gpu4pyscf/df/hessian/tests/test_df_uks_hessian.py
@@ -67,6 +67,7 @@ def test_df_lda(self):
 
         mf = mf.to_gpu()
         hessobj = mf.Hessian()
+        hessobj.auxbasis_response = 1
         hess_gpu = hessobj.kernel()
         assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
 
@@ -82,6 +83,7 @@ def test_df_gga(self):
 
         mf = mf.to_gpu()
         hessobj = mf.Hessian()
+        hessobj.auxbasis_response = 1
         hessobj.base.cphf_grids = hessobj.base.grids
         hess_gpu = hessobj.kernel()
         assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
@@ -98,6 +100,7 @@ def test_df_mgga(self):
 
         mf = mf.to_gpu()
         hessobj = mf.Hessian()
+        hessobj.auxbasis_response = 1
         hessobj.base.cphf_grids = hessobj.base.grids
         hess_gpu = hessobj.kernel()
         assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
@@ -105,4 +108,3 @@ def test_df_mgga(self):
 if __name__ == "__main__":
     print("Full Tests for DF UKS Hessian")
     unittest.main()
-    
diff --git a/gpu4pyscf/df/hessian/uhf.py b/gpu4pyscf/df/hessian/uhf.py
index 872fcf572..ab56ecd24 100644
--- a/gpu4pyscf/df/hessian/uhf.py
+++ b/gpu4pyscf/df/hessian/uhf.py
@@ -710,7 +710,7 @@ class Hessian(uhf_hess.Hessian):
 
     _keys = {'auxbasis_response',}
 
-    auxbasis_response = 1
+    auxbasis_response = 2
     partial_hess_elec = partial_hess_elec
     make_h1 = make_h1
     get_jk_mo = _get_jk_mo
diff --git a/gpu4pyscf/df/hessian/uks.py b/gpu4pyscf/df/hessian/uks.py
index 5b20c1711..54c990c81 100644
--- a/gpu4pyscf/df/hessian/uks.py
+++ b/gpu4pyscf/df/hessian/uks.py
@@ -135,7 +135,7 @@ class Hessian(uks_hess.Hessian):
 
     _keys = {'auxbasis_response',}
 
-    auxbasis_response = 1
+    auxbasis_response = 2
     partial_hess_elec = partial_hess_elec
     make_h1 = make_h1
     get_jk_mo = df_uhf_hess._get_jk_mo
diff --git a/gpu4pyscf/df/tests/test_df_tdrhf_grad.py b/gpu4pyscf/df/tests/test_df_tdrhf_grad.py
index 919b6f9b1..a3f274145 100644
--- a/gpu4pyscf/df/tests/test_df_tdrhf_grad.py
+++ b/gpu4pyscf/df/tests/test_df_tdrhf_grad.py
@@ -19,7 +19,6 @@
 from pyscf import scf, dft, tdscf
 import gpu4pyscf
 from gpu4pyscf import scf as gpu_scf
-from packaging import version
 
 atom = """
 O       0.0000000000     0.0000000000     0.0000000000
@@ -27,8 +26,6 @@
 H       0.0000000000     0.7570000000     0.5870000000
 """
 
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse("2.5.0")
-
 bas0 = "def2svpd"
 
 
diff --git a/gpu4pyscf/df/tests/test_df_tdrks_grad.py b/gpu4pyscf/df/tests/test_df_tdrks_grad.py
index 5ebaa015a..da414bb27 100644
--- a/gpu4pyscf/df/tests/test_df_tdrks_grad.py
+++ b/gpu4pyscf/df/tests/test_df_tdrks_grad.py
@@ -19,7 +19,6 @@
 from pyscf import scf, dft, tdscf
 import gpu4pyscf
 from gpu4pyscf import scf as gpu_scf
-from packaging import version
 
 atom = """
 O       0.0000000000     0.0000000000     0.0000000000
@@ -27,8 +26,6 @@
 H       0.0000000000     0.7570000000     0.5870000000
 """
 
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse("2.5.0")
-
 bas0 = "def2svpd"
 
 
diff --git a/gpu4pyscf/df/tests/test_df_tduhf_grad.py b/gpu4pyscf/df/tests/test_df_tduhf_grad.py
index 08af7b5c3..f2a10297f 100644
--- a/gpu4pyscf/df/tests/test_df_tduhf_grad.py
+++ b/gpu4pyscf/df/tests/test_df_tduhf_grad.py
@@ -19,7 +19,6 @@
 from pyscf import scf, dft, tdscf
 import gpu4pyscf
 from gpu4pyscf import scf as gpu_scf
-from packaging import version
 
 atom = """
 O       0.0000000000     0.0000000000     0.0000000000
@@ -27,8 +26,6 @@
 H       0.0000000000     0.7570000000     0.5870000000
 """
 
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse("2.5.0")
-
 bas0 = "cc-pvdz"
 
 
diff --git a/gpu4pyscf/df/tests/test_df_tduks_grad.py b/gpu4pyscf/df/tests/test_df_tduks_grad.py
index 96d2825f6..96510bbeb 100644
--- a/gpu4pyscf/df/tests/test_df_tduks_grad.py
+++ b/gpu4pyscf/df/tests/test_df_tduks_grad.py
@@ -19,7 +19,6 @@
 from pyscf import scf, dft, tdscf
 import gpu4pyscf
 from gpu4pyscf import scf as gpu_scf
-from packaging import version
 
 atom = """
 O       0.0000000000     0.0000000000     0.0000000000
@@ -27,8 +26,6 @@
 H       0.0000000000     0.7570000000     0.5870000000
 """
 
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse("2.5.0")
-
 bas0 = "cc-pvdz"
 
 
diff --git a/gpu4pyscf/dft/gks.py b/gpu4pyscf/dft/gks.py
index 8d83e3f92..48daa8cb0 100644
--- a/gpu4pyscf/dft/gks.py
+++ b/gpu4pyscf/dft/gks.py
@@ -23,7 +23,7 @@ class GKS(gks.GKS, GHF):
     def __init__(self, mol, xc='LDA,VWN'):
         raise NotImplementedError
 
-    reset = rks.RKS.reset
+    reset = rks.KohnShamDFT.reset
     energy_elec = rks.RKS.energy_elec
     get_veff = NotImplemented
     nuc_grad_method = NotImplemented
diff --git a/gpu4pyscf/dft/rks.py b/gpu4pyscf/dft/rks.py
index a6c62fd75..5826e1540 100644
--- a/gpu4pyscf/dft/rks.py
+++ b/gpu4pyscf/dft/rks.py
@@ -245,7 +245,21 @@ def dump_flags(self, verbose=None):
         log.info('small_rho_cutoff = %g', self.small_rho_cutoff)
         return self
 
-    reset = rks.KohnShamDFT.reset
+    def reset(self, mol=None):
+        hf.SCF.reset(self, mol)
+        self.grids.reset(mol)
+        self.nlcgrids.reset(mol)
+        self._numint.reset()
+        # The cphf_grids attribute is not available in the PySCF CPU version.
+        # In PySCF's to_gpu() function, this attribute is not initialized.
+        if hasattr(self, 'cphf_grids'):
+            self.cphf_grids.reset(self.mol)
+        else:
+            cphf_grids = self.cphf_grids = gen_grid.Grids(self.mol)
+            cphf_grids.prune = gen_grid.sg1_prune
+            cphf_grids.atom_grid = (50,194)
+        return self
+
     do_nlc = rks.KohnShamDFT.do_nlc
 
 hf.KohnShamDFT = KohnShamDFT
@@ -264,19 +278,6 @@ def dump_flags(self, verbose=None):
         hf.RHF.dump_flags(self, verbose)
         return KohnShamDFT.dump_flags(self, verbose)
 
-    def reset(self, mol=None):
-        hf.SCF.reset(self, mol)
-        self.grids.reset(mol)
-        self.nlcgrids.reset(mol)
-        self._numint.reset()
-        # The cphf_grids attribute is not available in the PySCF CPU version.
-        # In PySCF's to_gpu() function, this attribute is not properly
-        # initialized. mol of the KS object must be used for initialization.
-        if mol is None:
-            mol = self.mol
-        self.cphf_grids.reset(mol)
-        return self
-
     def nuc_grad_method(self):
         from gpu4pyscf.grad import rks as rks_grad
         return rks_grad.Gradients(self)
diff --git a/gpu4pyscf/dft/uks.py b/gpu4pyscf/dft/uks.py
index eb412ccc8..4eeef1e2a 100644
--- a/gpu4pyscf/dft/uks.py
+++ b/gpu4pyscf/dft/uks.py
@@ -119,17 +119,8 @@ def __init__(self, mol, xc='LDA,VWN'):
     energy_elec = energy_elec
     energy_tot = hf.RHF.energy_tot
     init_guess_by_vsap = uks_cpu.UKS.init_guess_by_vsap
-
     to_hf = NotImplemented
 
-    def reset(self, mol=None):
-        hf.SCF.reset(self, mol)
-        self.grids.reset(mol)
-        self.nlcgrids.reset(mol)
-        self.cphf_grids.reset(mol)
-        self._numint.reset()
-        return self
-
     def nuc_grad_method(self):
         from gpu4pyscf.grad import uks as uks_grad
         return uks_grad.Gradients(self)
diff --git a/gpu4pyscf/grad/tests/test_grid_response.py b/gpu4pyscf/grad/tests/test_grid_response.py
index 3918aa32c..2bf9d4eb1 100644
--- a/gpu4pyscf/grad/tests/test_grid_response.py
+++ b/gpu4pyscf/grad/tests/test_grid_response.py
@@ -19,19 +19,16 @@
 import pytest
 from pyscf.dft import rks as cpu_rks
 from gpu4pyscf.dft import rks as gpu_rks
-from packaging import version
 
-atom = '''
-O       0.0000000000    -0.0000000000     0.1174000000
-H      -0.7570000000    -0.0000000000    -0.4696000000
-H       0.7570000000     0.0000000000    -0.4696000000
-'''
-
-bas0='def2-tzvpp'
-grids_level = 5
-nlcgrids_level = 3
 def setUpModule():
     global mol_sph, mol_cart
+    atom = '''
+    O       0.0000000000    -0.0000000000     0.1174000000
+    H      -0.7570000000    -0.0000000000    -0.4696000000
+    H       0.7570000000     0.0000000000    -0.4696000000
+    '''
+
+    bas0='def2-tzvpp'
     mol_sph = pyscf.M(atom=atom, basis=bas0, max_memory=32000,
                       output='/dev/null', verbose=1)
 
diff --git a/gpu4pyscf/grad/tests/test_rhf_grad.py b/gpu4pyscf/grad/tests/test_rhf_grad.py
index 8a5773587..db9a8b2ee 100644
--- a/gpu4pyscf/grad/tests/test_rhf_grad.py
+++ b/gpu4pyscf/grad/tests/test_rhf_grad.py
@@ -23,20 +23,16 @@
 from pyscf.grad import rhf as rhf_grad_cpu
 from gpu4pyscf.grad import rhf as rhf_grad_gpu
 from gpu4pyscf.lib.multi_gpu import num_devices
-from packaging import version
-
-atom = '''
-O       0.0000000000    -0.0000000000     0.1174000000
-H      -0.7570000000    -0.0000000000    -0.4696000000
-H       0.7570000000     0.0000000000    -0.4696000000
-'''
-
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse('2.5.0')
-
-bas0='cc-pvtz'
 
 def setUpModule():
     global mol_sph, mol_cart
+    atom = '''
+    O       0.0000000000    -0.0000000000     0.1174000000
+    H      -0.7570000000    -0.0000000000    -0.4696000000
+    H       0.7570000000     0.0000000000    -0.4696000000
+    '''
+    bas0='cc-pvtz'
+
     mol_sph = pyscf.M(atom=atom, basis=bas0, max_memory=32000,
                       output='/dev/null', verbose=1)
 
@@ -70,11 +66,9 @@ def test_grad_rhf(self):
     def test_grad_cart(self):
         _check_grad(mol_cart, tol=1e-6)
 
-    @pytest.mark.skipif(pyscf_25, reason='requires pyscf 2.6 or higher')
     def test_grad_d3bj(self):
         _check_grad(mol_sph, tol=1e-6, disp='d3bj')
 
-    @pytest.mark.skipif(pyscf_25, reason='requires pyscf 2.6 or higher')
     def test_grad_d4(self):
         _check_grad(mol_sph, tol=1e-6, disp='d4')
 
diff --git a/gpu4pyscf/grad/tests/test_rks_grad.py b/gpu4pyscf/grad/tests/test_rks_grad.py
index 3368785aa..75af61a6a 100644
--- a/gpu4pyscf/grad/tests/test_rks_grad.py
+++ b/gpu4pyscf/grad/tests/test_rks_grad.py
@@ -18,7 +18,6 @@
 import pytest
 from pyscf.dft import rks as cpu_rks
 from gpu4pyscf.dft import rks as gpu_rks
-from packaging import version
 
 atom = '''
 O       0.0000000000    -0.0000000000     0.1174000000
diff --git a/gpu4pyscf/grad/tests/test_tddft_opt.py b/gpu4pyscf/grad/tests/test_tddft_opt.py
index 533c458fc..ffe78e7cd 100644
--- a/gpu4pyscf/grad/tests/test_tddft_opt.py
+++ b/gpu4pyscf/grad/tests/test_tddft_opt.py
@@ -20,7 +20,6 @@
 from pyscf.geomopt.geometric_solver import optimize
 import gpu4pyscf
 from gpu4pyscf import scf as gpu_scf
-from packaging import version
 from gpu4pyscf.lib.multi_gpu import num_devices
 
 atom = """
diff --git a/gpu4pyscf/grad/tests/test_tdrhf_grad.py b/gpu4pyscf/grad/tests/test_tdrhf_grad.py
index 52430d8b1..208661c81 100644
--- a/gpu4pyscf/grad/tests/test_tdrhf_grad.py
+++ b/gpu4pyscf/grad/tests/test_tdrhf_grad.py
@@ -19,7 +19,6 @@
 from pyscf import scf, dft, tdscf
 import gpu4pyscf
 from gpu4pyscf import scf as gpu_scf
-from packaging import version
 from gpu4pyscf.lib.multi_gpu import num_devices
 
 atom = """
@@ -28,8 +27,6 @@
 H       0.0000000000     0.7570000000     0.5870000000
 """
 
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse("2.5.0")
-
 bas0 = "cc-pvdz"
 
 def diagonalize(a, b, nroots=5):
diff --git a/gpu4pyscf/grad/tests/test_tdrks_grad.py b/gpu4pyscf/grad/tests/test_tdrks_grad.py
index 6c285b87a..cea339380 100644
--- a/gpu4pyscf/grad/tests/test_tdrks_grad.py
+++ b/gpu4pyscf/grad/tests/test_tdrks_grad.py
@@ -19,7 +19,6 @@
 from pyscf import scf, dft, tdscf
 import gpu4pyscf
 from gpu4pyscf import scf as gpu_scf
-from packaging import version
 from gpu4pyscf.lib.multi_gpu import num_devices
 
 atom = """
@@ -28,8 +27,6 @@
 H       0.0000000000     0.7570000000     0.5870000000
 """
 
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse("2.5.0")
-
 bas0 = "cc-pvdz"
 
 def diagonalize(a, b, nroots=5):
diff --git a/gpu4pyscf/grad/tests/test_tduhf_grad.py b/gpu4pyscf/grad/tests/test_tduhf_grad.py
index 68e67a08b..0c4202815 100644
--- a/gpu4pyscf/grad/tests/test_tduhf_grad.py
+++ b/gpu4pyscf/grad/tests/test_tduhf_grad.py
@@ -19,7 +19,6 @@
 from pyscf import scf, dft, tdscf
 import gpu4pyscf
 from gpu4pyscf import scf as gpu_scf
-from packaging import version
 from gpu4pyscf.lib.multi_gpu import num_devices
 
 atom = """
@@ -28,8 +27,6 @@
 H       0.0000000000     0.7570000000     0.5870000000
 """
 
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse("2.5.0")
-
 bas0 = "cc-pvdz"
 
 def diagonalize(a, b, nroots=5):
diff --git a/gpu4pyscf/grad/tests/test_tduks_grad.py b/gpu4pyscf/grad/tests/test_tduks_grad.py
index 6b40a5bbb..6bf67499c 100644
--- a/gpu4pyscf/grad/tests/test_tduks_grad.py
+++ b/gpu4pyscf/grad/tests/test_tduks_grad.py
@@ -19,7 +19,6 @@
 from pyscf import scf, dft, tdscf
 import gpu4pyscf
 from gpu4pyscf import scf as gpu_scf
-from packaging import version
 from gpu4pyscf.lib.multi_gpu import num_devices
 
 atom = """
@@ -28,8 +27,6 @@
 H       0.0000000000     0.7570000000     0.5870000000
 """
 
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse("2.5.0")
-
 bas0 = "cc-pvdz"
 
 def diagonalize(a, b, nroots=5):
diff --git a/gpu4pyscf/grad/tests/test_uhf_grad.py b/gpu4pyscf/grad/tests/test_uhf_grad.py
index 9e072f492..b17fa3acb 100644
--- a/gpu4pyscf/grad/tests/test_uhf_grad.py
+++ b/gpu4pyscf/grad/tests/test_uhf_grad.py
@@ -18,7 +18,6 @@
 import pytest
 from pyscf import lib, gto
 from gpu4pyscf import scf
-from packaging import version
 
 atom = '''
 O       0.0000000000    -0.0000000000     0.1174000000
@@ -26,8 +25,6 @@
 H       0.7570000000     0.0000000000    -0.4696000000
 '''
 
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse('2.5.0')
-
 bas0='cc-pvtz'
 
 def setUpModule():
@@ -67,12 +64,10 @@ def test_grad_cart(self):
         print('---- testing UHF Cart -------')
         _check_grad(mol_cart, tol=1e-10)
 
-    @pytest.mark.skipif(pyscf_25, reason='requires pyscf 2.6 or higher')
     def test_grad_d3bj(self):
         print('---- testing UHF with D3(BJ) ----')
         _check_grad(mol_sph, tol=1e-6, disp='d3bj')
 
-    @pytest.mark.skipif(pyscf_25, reason='requires pyscf 2.6 or higher')
     def test_grad_d4(self):
         print('------- UHF with D4 -----')
         _check_grad(mol_sph, tol=1e-6, disp='d4')
diff --git a/gpu4pyscf/grad/tests/test_uks_grad.py b/gpu4pyscf/grad/tests/test_uks_grad.py
index 17f192f62..e1e601f4a 100644
--- a/gpu4pyscf/grad/tests/test_uks_grad.py
+++ b/gpu4pyscf/grad/tests/test_uks_grad.py
@@ -17,7 +17,6 @@
 import unittest
 import pytest
 from gpu4pyscf.dft import uks
-from packaging import version
 
 atom = '''
 O       0.0000000000    -0.0000000000     0.1174000000
diff --git a/gpu4pyscf/gto/tests/test_int1e_grids_ip.py b/gpu4pyscf/gto/tests/test_int1e_grids_ip.py
index 4e46382b2..bb10aa71a 100644
--- a/gpu4pyscf/gto/tests/test_int1e_grids_ip.py
+++ b/gpu4pyscf/gto/tests/test_int1e_grids_ip.py
@@ -65,6 +65,7 @@ class KnownValues(unittest.TestCase):
     def test_int1e_grids_ip_full_tensor_cart(self):
         mol = mol_cart
         fakemol = gto.fakemol_for_charges(grid_points)
+        fakemol.cart = True
 
         int3c2e_ip1 = mol._add_suffix('int3c2e_ip1')
         cintopt = gto.moleintor.make_cintopt(mol._atm, mol._bas, mol._env, int3c2e_ip1)
@@ -183,6 +184,7 @@ def test_int1e_grids_ip_contracted_cart(self):
 
         mol = mol_cart
         fakemol = gto.fakemol_for_charges(grid_points)
+        fakemol.cart = True
 
         int3c2e_ip1 = mol._add_suffix('int3c2e_ip1')
         cintopt = gto.moleintor.make_cintopt(mol._atm, mol._bas, mol._env, int3c2e_ip1)
diff --git a/gpu4pyscf/gto/tests/test_int1e_grids_ipip.py b/gpu4pyscf/gto/tests/test_int1e_grids_ipip.py
index 388d721dd..936db6e1a 100644
--- a/gpu4pyscf/gto/tests/test_int1e_grids_ipip.py
+++ b/gpu4pyscf/gto/tests/test_int1e_grids_ipip.py
@@ -68,6 +68,7 @@ def test_int1e_grids_ipip1_charge_contracted_cart(self):
 
         mol = mol_cart
         fakemol = gto.fakemol_for_charges(grid_points)
+        fakemol.cart = True
 
         int3c2e_ipip1 = mol._add_suffix('int3c2e_ipip1')
         cintopt = gto.moleintor.make_cintopt(mol._atm, mol._bas, mol._env, int3c2e_ipip1)
@@ -170,6 +171,7 @@ def test_int1e_grids_ipvip1_charge_contracted_cart(self):
 
         mol = mol_cart
         fakemol = gto.fakemol_for_charges(grid_points)
+        fakemol.cart = True
 
         int3c2e_ipvip1 = mol._add_suffix('int3c2e_ipvip1')
         cintopt = gto.moleintor.make_cintopt(mol._atm, mol._bas, mol._env, int3c2e_ipvip1)
@@ -272,6 +274,7 @@ def test_int1e_grids_ip1ip2_charge_contracted_cart(self):
 
         mol = mol_cart
         fakemol = gto.fakemol_for_charges(grid_points)
+        fakemol.cart = True
 
         int3c2e_ip1ip2 = mol._add_suffix('int3c2e_ip1ip2')
         cintopt = gto.moleintor.make_cintopt(mol._atm, mol._bas, mol._env, int3c2e_ip1ip2)
@@ -374,6 +377,7 @@ def test_int1e_grids_ipip2_charge_contracted_cart(self):
 
         mol = mol_cart
         fakemol = gto.fakemol_for_charges(grid_points)
+        fakemol.cart = True
 
         # Note: we cannot compute ipip2 (dCdC) directly due to numerical problems,
         #       pyscf treat a point charge as a sharp Gaussian, and we cannot take 2nd derivative of it.
diff --git a/gpu4pyscf/hessian/rhf.py b/gpu4pyscf/hessian/rhf.py
index b0a56722e..2ac74c1f3 100644
--- a/gpu4pyscf/hessian/rhf.py
+++ b/gpu4pyscf/hessian/rhf.py
@@ -953,6 +953,12 @@ def dump_flags(self, verbose=None):
                  self.max_memory, lib.current_memory()[0])
         return self
 
+    def reset(self, mol=None):
+        if mol is not None:
+            self.mol = mol
+        self.base.reset(mol)
+        return self
+
 class Hessian(HessianBase):
     '''Non-relativistic restricted Hartree-Fock hessian'''
 
diff --git a/gpu4pyscf/lib/tests/test_to_gpu.py b/gpu4pyscf/lib/tests/test_to_gpu.py
index dd11b9fce..5d1e9554d 100644
--- a/gpu4pyscf/lib/tests/test_to_gpu.py
+++ b/gpu4pyscf/lib/tests/test_to_gpu.py
@@ -19,24 +19,17 @@
 import pytest
 from pyscf import scf, lib
 from pyscf.dft import rks
-from packaging import version
-
-atom = '''
-O       0.0000000000    -0.0000000000     0.1174000000
-H      -0.7570000000    -0.0000000000    -0.4696000000
-H       0.7570000000     0.0000000000    -0.4696000000
-'''
-
-bas='sto3g'
-grids_level = 1
-pyscf_24 = version.parse(pyscf.__version__) <= version.parse('2.4.0')
 
 def setUpModule():
     global mol
-    mol = pyscf.M(atom=atom, basis=bas, max_memory=32000)
-    mol.output = '/dev/null'
-    mol.build()
-    mol.verbose = 1
+    atom = '''
+    O       0.0000000000    -0.0000000000     0.1174000000
+    H      -0.7570000000    -0.0000000000    -0.4696000000
+    H       0.7570000000     0.0000000000    -0.4696000000
+    '''
+    bas='sto3g'
+    mol = pyscf.M(atom=atom, basis=bas, max_memory=32000, output = '/dev/null',
+                  verbose=6)
 
 def tearDownModule():
     global mol
@@ -44,7 +37,6 @@ def tearDownModule():
     del mol
 
 class KnownValues(unittest.TestCase):
-    @pytest.mark.skipif(pyscf_24, reason='requires pyscf 2.5 or higher')
     def test_rhf(self):
         mf = scf.RHF(mol).to_gpu()
         e_tot = mf.to_gpu().kernel()
@@ -59,7 +51,6 @@ def test_rhf(self):
         # h = mf.Hessian().to_gpu()
         # h.kernel()
 
-    @pytest.mark.skipif(pyscf_24, reason='requires pyscf 2.5 or higher')
     def test_rks(self):
         mf = rks.RKS(mol).to_gpu()
         e_tot = mf.to_gpu().kernel()
@@ -75,7 +66,6 @@ def test_rks(self):
         # h = mf.Hessian().to_gpu()
         # h.kernel()
 
-    @pytest.mark.skipif(pyscf_24, reason='requires pyscf 2.5 or higher')
     def test_df_RHF(self):
         mf = scf.RHF(mol).density_fit().to_gpu()
         e_tot = mf.to_gpu().kernel()
@@ -89,10 +79,10 @@ def test_df_RHF(self):
         mf = scf.RHF(mol).density_fit().run()
         mf.conv_tol_cpscf = 1e-7
         hobj = mf.Hessian().to_gpu()
+        hobj.auxbasis_response = 1
         h = hobj.kernel()
         assert numpy.abs(lib.fp(h) - 2.198079352288524) < 1e-4
 
-    @pytest.mark.skipif(pyscf_24, reason='requires pyscf 2.5 or higher')
     def test_df_b3lyp(self):
         mf = rks.RKS(mol, xc='b3lyp').density_fit().to_gpu()
         e_tot = mf.to_gpu().kernel()
@@ -106,10 +96,10 @@ def test_df_b3lyp(self):
         mf = rks.RKS(mol, xc='b3lyp').density_fit().run()
         mf.conv_tol_cpscf = 1e-7
         hobj = mf.Hessian().to_gpu()
+        hobj.auxbasis_response = 1
         h = hobj.kernel()
         assert numpy.abs(lib.fp(h) - 2.1527804103141848) < 1e-4
 
-    @pytest.mark.skipif(pyscf_24, reason='requires pyscf 2.5 or higher')
     def test_df_RKS(self):
         mf = rks.RKS(mol, xc='wb97x').density_fit().to_gpu()
         e_tot = mf.to_gpu().kernel()
@@ -123,10 +113,10 @@ def test_df_RKS(self):
         mf = rks.RKS(mol, xc='wb97x').density_fit().run()
         mf.conv_tol_cpscf = 1e-7
         hobj = mf.Hessian().to_gpu()
+        hobj.auxbasis_response = 1
         h = hobj.kernel()
         assert numpy.abs(lib.fp(h) - 2.1858589608638384) < 1e-4
 
 if __name__ == "__main__":
     print("Full tests for to_gpu module")
     unittest.main()
-    
diff --git a/gpu4pyscf/mp/dfmp2.py b/gpu4pyscf/mp/dfmp2.py
index 23ca83a3f..b337a099e 100644
--- a/gpu4pyscf/mp/dfmp2.py
+++ b/gpu4pyscf/mp/dfmp2.py
@@ -163,3 +163,5 @@ def update_amps(self, t2, eris):
 
     def init_amps(self, mo_energy=None, mo_coeff=None, eris=None, with_t2=WITH_T2):
         return kernel(self, mo_energy, mo_coeff, eris, with_t2)
+
+DFRMP2 = DFMP2
diff --git a/gpu4pyscf/mp/tests/test_mp2.py b/gpu4pyscf/mp/tests/test_mp2.py
index b51278166..1c944945e 100644
--- a/gpu4pyscf/mp/tests/test_mp2.py
+++ b/gpu4pyscf/mp/tests/test_mp2.py
@@ -47,9 +47,6 @@ def tearDownModule():
     mol.stdout.close()
     del mol, mf
 
-import pyscf
-from packaging import version
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse('2.5.0')
 class KnownValues(unittest.TestCase):
     def test_mp2(self):
         nocc = mol.nelectron//2
@@ -139,7 +136,6 @@ def test_to_cpu(self):
         e_cpu = pt.kernel()[0]
         assert abs(e_cpu - e_gpu) < 1e-6
 
-    @pytest.mark.skipif(pyscf_25, reason='requires pyscf 2.6 or higher')
     def test_to_gpu(self):
         pt = mp_cpu.mp2.MP2(mf)
         e_cpu = pt.kernel()[0]
diff --git a/gpu4pyscf/pbc/dft/tests/test_multigrid.py b/gpu4pyscf/pbc/dft/tests/test_multigrid.py
index a283047c7..299945963 100644
--- a/gpu4pyscf/pbc/dft/tests/test_multigrid.py
+++ b/gpu4pyscf/pbc/dft/tests/test_multigrid.py
@@ -81,13 +81,11 @@ class KnownValues(unittest.TestCase):
     def test_get_pp(self):
         ref = MultiGridNumInt_cpu(cell_orth).get_pp()
         out = multigrid.MultiGridNumInt(cell_orth).get_pp().get()
-        self.assertEqual(out.shape, ref.shape)
         self.assertAlmostEqual(abs(ref-out).max(), 0, 9)
 
     def test_get_nuc(self):
         ref = MultiGridNumInt_cpu(cell_orth).get_nuc()
         out = multigrid.MultiGridNumInt(cell_orth).get_nuc().get()
-        self.assertEqual(out.shape, ref.shape)
         self.assertAlmostEqual(abs(ref-out).max(), 0, 9)
 
     def test_eval_nucG(self):
@@ -132,7 +130,10 @@ def test_get_rho(self):
         np.random.seed(2)
         dm = np.random.random((nao,nao)) - .5
         dm = dm.dot(dm.T)
-        ref = MultiGridNumInt_cpu(cell_orth).get_rho(dm)
+        if hasattr(multigrid_cpu, 'MultiGridNumInt'):
+            ref = multigrid_cpu.MultiGridNumInt(cell_orth).get_rho(cell_orth, dm, None)
+        else:
+            ref = multigrid_cpu.MultiGridFFTDF(cell_orth).get_rho(dm)
         out = multigrid.MultiGridNumInt(cell_orth).get_rho(dm).get()
         self.assertAlmostEqual(abs(ref-out).max(), 0, 9)
 
diff --git a/gpu4pyscf/solvent/tests/test_pcm.py b/gpu4pyscf/solvent/tests/test_pcm.py
index 3c0a08ce8..7aee174ae 100644
--- a/gpu4pyscf/solvent/tests/test_pcm.py
+++ b/gpu4pyscf/solvent/tests/test_pcm.py
@@ -20,15 +20,6 @@
 from pyscf import gto
 from gpu4pyscf import scf, dft
 from gpu4pyscf.solvent import pcm
-from packaging import version
-try:
-    # Some PCM methods are registered when importing the CPU version.
-    # However, pyscf-2.7 does note automatically import this module.
-    from pyscf.solvent import pcm as pcm_on_cpu
-except ImportError:
-    pass
-
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse('2.5.0')
 
 def setUpModule():
     global mol, epsilon, lebedev_order
@@ -118,51 +109,28 @@ def test_dfuks(self):
         e_tot = _energy_with_solvent(dft.UKS(mol, xc='b3lyp').density_fit(), 'IEF-PCM')
         print(f"Energy error in DFUKS with IEF-PCM: {numpy.abs(e_tot - -71.67135250643567)}")
         assert numpy.abs(e_tot - -71.67135250643567) < 1e-5
-        
-    def test_to_cpu(self):
-        mf = dft.RKS(mol, xc='b3lyp')
-        e_gpu = mf.kernel()
-        mf = mf.to_cpu()
-        e_cpu = mf.kernel()
-        assert abs(e_cpu - e_gpu) < 1e-8
 
-        mf = dft.RKS(mol, xc='b3lyp').density_fit()
-        e_gpu = mf.kernel()
-        chg = mf.analyze()[0][1]
-        mf_cpu = mf.to_cpu()
-        e_cpu = mf_cpu.kernel()
-        assert abs(e_cpu - e_gpu) < 1e-8
-        chg_ref = mf_cpu.analyze()[0][1]
-        assert abs(chg - chg_ref).max() < 1e-5
-
-    @pytest.mark.skipif(pyscf_25, reason='requires pyscf 2.6 or higher')
     def test_to_gpu(self):
-        import pyscf
-        mf = pyscf.dft.RKS(mol, xc='b3lyp').PCM()
-        e_cpu = mf.kernel()
-        mf = mf.to_gpu()
-        e_gpu = mf.kernel()
-        assert abs(e_cpu - e_gpu) < 1e-8
-
-        mf = pyscf.dft.RKS(mol, xc='b3lyp').density_fit().PCM()
+        mf = mol.RKS(xc='b3lyp').PCM()
         e_cpu = mf.kernel()
         mf = mf.to_gpu()
         e_gpu = mf.kernel()
         assert abs(e_cpu - e_gpu) < 1e-8
-
-    @pytest.mark.skipif(pyscf_25, reason='requires pyscf 2.6 or higher')
-    def test_to_cpu_1(self):
-        mf = dft.RKS(mol, xc='b3lyp').PCM()
-        e_gpu = mf.kernel()
         mf = mf.to_cpu()
         e_cpu = mf.kernel()
         assert abs(e_cpu - e_gpu) < 1e-8
 
-        mf = dft.RKS(mol, xc='b3lyp').density_fit().PCM()
-        e_gpu = mf.kernel()
-        mf = mf.to_cpu()
+        mf = mol.RKS(xc='b3lyp').density_fit().PCM()
         e_cpu = mf.kernel()
+        mf = mf.to_gpu()
+        e_gpu = mf.kernel()
         assert abs(e_cpu - e_gpu) < 1e-8
+        chg = mf.analyze()[0][1]
+        mf_cpu = mf.to_cpu()
+        e_cpu = mf_cpu.kernel()
+        assert abs(e_cpu - e_gpu) < 1e-8
+        chg_ref = mf_cpu.analyze()[0][1]
+        assert abs(chg - chg_ref).max() < 1e-5
 
     def test_df_and_pcm(self):
         mol = gto.M(atom='H 0 0 0; H 0 0 1')
diff --git a/gpu4pyscf/solvent/tests/test_pcm_grad.py b/gpu4pyscf/solvent/tests/test_pcm_grad.py
index c17e05f3c..d728e9379 100644
--- a/gpu4pyscf/solvent/tests/test_pcm_grad.py
+++ b/gpu4pyscf/solvent/tests/test_pcm_grad.py
@@ -21,9 +21,6 @@
 from gpu4pyscf import scf
 from gpu4pyscf.solvent import pcm
 from gpu4pyscf.solvent.grad import pcm as pcm_grad
-from packaging import version
-
-pyscf_25 = version.parse(pyscf.__version__) <= version.parse('2.5.0')
 
 def setUpModule():
     global mol, epsilon, lebedev_order
@@ -191,7 +188,6 @@ def test_uhf_grad_IEFPCM(self):
         print(f"Gradient error in UHF with IEFPCM: {numpy.linalg.norm(g0 - grad)}")
         assert numpy.linalg.norm(g0 - grad) < 1e-6
 
-    @pytest.mark.skipif(pyscf_25, reason='requires pyscf 2.6 or higher')
     def test_to_cpu(self):
         mf = scf.RHF(mol).PCM()
         mf.verbose = 0
@@ -215,7 +211,6 @@ def test_to_cpu(self):
         grad_cpu = gradobj.kernel()
         assert numpy.linalg.norm(grad_gpu - grad_cpu) < 1e-8
 
-    @pytest.mark.skipif(pyscf_25, reason='requires pyscf 2.6 or higher')
     def test_to_gpu(self):
         mf = pyscf.scf.RHF(mol).PCM()
         mf.verbose = 0
diff --git a/gpu4pyscf/solvent/tests/test_pcm_hessian.py b/gpu4pyscf/solvent/tests/test_pcm_hessian.py
index e962c344e..ba78bfdc0 100644
--- a/gpu4pyscf/solvent/tests/test_pcm_hessian.py
+++ b/gpu4pyscf/solvent/tests/test_pcm_hessian.py
@@ -20,13 +20,10 @@
 from pyscf import gto
 from gpu4pyscf.solvent import pcm
 from gpu4pyscf import scf, dft
-from packaging import version
 from gpu4pyscf.solvent.hessian.pcm import analytical_grad_vmat, analytical_hess_nuc, analytical_hess_solver, analytical_hess_qv
 from gpu4pyscf.lib.cupy_helper import contract
 from gpu4pyscf.lib.multi_gpu import num_devices
 
-pyscf_211 = version.parse(pyscf.__version__) <= version.parse('2.11.0')
-
 def setUpModule():
     global mol, epsilon, lebedev_order, eps, xc, tol
     mol = gto.Mole()
@@ -325,7 +322,6 @@ def test_hess_solver_ssvpe(self):
 
         assert abs(ref_grad_vmat - test_grad_vmat).max() < 1e-9
 
-    @pytest.mark.skipif(pyscf_211, reason='requires pyscf 2.12 or higher')
     def test_to_gpu_to_cpu(self):
         mol = gto.Mole()
         mol.atom = '''
@@ -346,10 +342,10 @@ def test_to_gpu_to_cpu(self):
         hess_gpu = hessobj.kernel()
         hessobj = hessobj.to_cpu()
         hess_cpu = hessobj.kernel()
-        assert np.linalg.norm(hess_cpu - hess_gpu) < 1e-5
+        assert np.linalg.norm(hess_cpu - hess_gpu) < 2e-6
         hessobj = hessobj.to_gpu()
         hess_gpu = hessobj.kernel()
-        assert np.linalg.norm(hess_cpu - hess_gpu) < 1e-5
+        assert np.linalg.norm(hess_cpu - hess_gpu) < 2e-6
 
         mf = dft.RKS(mol, xc='b3lyp').density_fit().PCM()
         mf.conv_tol = 1e-12
@@ -357,15 +353,14 @@ def test_to_gpu_to_cpu(self):
         mf.grids.atom_grid = (50,194)
         mf.kernel()
         hessobj = mf.Hessian()
-        # The auxbasis_response attribute was not handled in pyscf-2.11
-        hessobj.auxbasis_response = 2
+        hessobj.auxbasis_response = 1
         hess_gpu = hessobj.kernel()
         hessobj = hessobj.to_cpu()
         hess_cpu = hessobj.kernel()
-        assert np.linalg.norm(hess_cpu - hess_gpu) < 1e-5
+        assert np.linalg.norm(hess_cpu - hess_gpu) < 2e-6
         hessobj = hessobj.to_gpu()
         hess_gpu = hessobj.kernel()
-        assert np.linalg.norm(hess_cpu - hess_gpu) < 1e-5
+        assert np.linalg.norm(hess_cpu - hess_gpu) < 2e-6
 
 if __name__ == "__main__":
     print("Full Tests for Hessian of PCMs")
diff --git a/gpu4pyscf/solvent/tests/test_smd_grad.py b/gpu4pyscf/solvent/tests/test_smd_grad.py
index ad8024967..ac1209324 100644
--- a/gpu4pyscf/solvent/tests/test_smd_grad.py
+++ b/gpu4pyscf/solvent/tests/test_smd_grad.py
@@ -246,46 +246,30 @@ def test_Br(self):
         _check_grad(atom, solvent='water')
         _check_grad(atom, solvent='toluene')
 
-    @pytest.mark.skipif(pyscf_211, reason='requires pyscf 2.12 or higher')
     def test_to_gpu(self):
-        import pyscf
-        mf = pyscf.dft.RKS(mol, xc='b3lyp').SMD()
+        mf = mol.RKS(xc='b3lyp').SMD()
         mf.conv_tol = 1e-12
         mf.kernel()
         gradobj = mf.nuc_grad_method()
         g_cpu = gradobj.kernel()
         gradobj = gradobj.to_gpu()
         g_gpu = gradobj.kernel()
-        assert numpy.linalg.norm(g_cpu - g_gpu) < 1e-5
-
-        mf = pyscf.dft.RKS(mol, xc='b3lyp').density_fit().SMD()
-        mf.conv_tol = 1e-12
-        mf.kernel()
-        gradobj = mf.nuc_grad_method()
-        g_cpu = gradobj.kernel()
-        gradobj = gradobj.to_gpu()
-        g_gpu = gradobj.kernel()
-        assert numpy.linalg.norm(g_cpu - g_gpu) < 1e-5
-
-    @pytest.mark.skipif(pyscf_211, reason='requires pyscf 2.12 or higher')
-    def test_to_cpu(self):
-        mf = dft.RKS(mol, xc='b3lyp').SMD()
-        mf.conv_tol = 1e-12
-        mf.kernel()
-        gradobj = mf.nuc_grad_method()
-        g_gpu = gradobj.kernel()
+        assert numpy.linalg.norm(g_cpu - g_gpu) < 1e-6
         gradobj = gradobj.to_cpu()
         g_cpu = gradobj.kernel()
-        assert numpy.linalg.norm(g_cpu - g_gpu) < 1e-5
+        assert numpy.linalg.norm(g_cpu - g_gpu) < 1e-6
 
-        mf = dft.RKS(mol, xc='b3lyp').density_fit().SMD()
+        mf = mol.RKS(xc='b3lyp').density_fit().SMD()
         mf.conv_tol = 1e-12
         mf.kernel()
         gradobj = mf.nuc_grad_method()
+        g_cpu = gradobj.kernel()
+        gradobj = gradobj.to_gpu()
         g_gpu = gradobj.kernel()
+        assert numpy.linalg.norm(g_cpu - g_gpu) < 1e-6
         gradobj = gradobj.to_cpu()
         g_cpu = gradobj.kernel()
-        assert numpy.linalg.norm(g_cpu - g_gpu) < 1e-5
+        assert numpy.linalg.norm(g_cpu - g_gpu) < 1e-6
 
 if __name__ == "__main__":
     print("Full Tests for Gradient of SMD")
diff --git a/gpu4pyscf/solvent/tests/test_smd_hessian.py b/gpu4pyscf/solvent/tests/test_smd_hessian.py
index 566b9786c..03605f80c 100644
--- a/gpu4pyscf/solvent/tests/test_smd_hessian.py
+++ b/gpu4pyscf/solvent/tests/test_smd_hessian.py
@@ -22,9 +22,6 @@
 from gpu4pyscf.solvent.hessian import smd as smd_hess
 from gpu4pyscf.solvent.grad import smd as smd_grad
 from gpu4pyscf.solvent import smd
-from packaging import version
-
-pyscf_211 = version.parse(pyscf.__version__) <= version.parse('2.11.0')
 
 def setUpModule():
     global mol
@@ -226,7 +223,6 @@ def test_Br(self):
         _check_hess(atom, solvent='water')
         _check_hess(atom, solvent='toluene')
 
-    @pytest.mark.skipif(pyscf_211, reason='requires pyscf 2.12 or higher')
     def test_to_gpu_to_cpu(self):
         mf = dft.RKS(mol, xc='b3lyp').SMD()
         mf.conv_tol = 1e-12
@@ -236,24 +232,24 @@ def test_to_gpu_to_cpu(self):
         hess_gpu = hessobj.kernel()
         hessobj = hessobj.to_cpu()
         hess_cpu = hessobj.kernel()
-        assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
+        assert numpy.linalg.norm(hess_cpu - hess_gpu) < 2e-6
         hessobj = hessobj.to_gpu()
         hess_gpu = hessobj.kernel()
-        assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
+        assert numpy.linalg.norm(hess_cpu - hess_gpu) < 2e-6
 
         mf = dft.RKS(mol, xc='b3lyp').density_fit().SMD()
         mf.conv_tol = 1e-12
         mf.conv_tol_cpscf = 1e-7
         mf.kernel()
         hessobj = mf.Hessian()
-        hessobj.auxbasis_response = 2
+        hessobj.auxbasis_response = 1
         hess_gpu = hessobj.kernel()
         hessobj = hessobj.to_cpu()
         hess_cpu = hessobj.kernel()
-        assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
+        assert numpy.linalg.norm(hess_cpu - hess_gpu) < 2e-6
         hessobj = hessobj.to_gpu()
         hess_gpu = hessobj.kernel()
-        assert numpy.linalg.norm(hess_cpu - hess_gpu) < 1e-5
+        assert numpy.linalg.norm(hess_cpu - hess_gpu) < 2e-6
 
     def test_cds(self):
         from gpu4pyscf.solvent.hessian.smd import get_cds as get_cds_hess