From c7ca579ac43fe506bc96518338e08f84da8a1ca0 Mon Sep 17 00:00:00 2001 From: Qiming Sun Date: Wed, 28 Aug 2024 13:13:10 -0700 Subject: [PATCH] Update autoaux tests --- pyscf/df/addons.py | 11 ++++++++--- pyscf/df/autoaux.py | 18 +++++++++++++++++- pyscf/df/test/test_addons.py | 29 ++++++++++++++++++++++------- 3 files changed, 47 insertions(+), 11 deletions(-) diff --git a/pyscf/df/addons.py b/pyscf/df/addons.py index 3a22528009..91884d3c45 100644 --- a/pyscf/df/addons.py +++ b/pyscf/df/addons.py @@ -31,6 +31,11 @@ ETB_BETA = getattr(__config__, 'df_addons_aug_dfbasis', 2.0) FIRST_ETB_ELEMENT = getattr(__config__, 'df_addons_aug_start_at', 36) # 'Rb' +# TODO: Switch to other default scheme for auxiliary basis generation. +# The auxiliary basis set generated by version 2.6 (and earlier) lacks compact +# functions. It may cause higher errors in ERI integrals. +USE_VERSION_26_AUXBASIS = True + # Obtained from http://www.psicode.org/psi4manual/master/basissets_byfamily.html DEFAULT_AUXBASIS = { # AO basis JK-fit MP2-fit @@ -94,9 +99,9 @@ def _aug_etb_element(nuc_charge, basis, beta): # 1: H - Be, 2: B - Ca, 3: Sc - La, 4: Ce - max_shells = 4 - conf.count(0) - if 0: - # This is the method that PySCF-2.6 (and older versions) generates - # auxiliary basis. It estimates the exponents ranges by geometric average. + if USE_VERSION_26_AUXBASIS: + # This is the method that version 2.6 (and earlier) generates auxiliary + # basis. It estimates the exponents ranges by geometric average. # This method is not recommended because it tends to generate diffused # functions. Important compact functions might be improperly excluded. l_max = min(l_max, max_shells) diff --git a/pyscf/df/autoaux.py b/pyscf/df/autoaux.py index af5ef50416..4c3268b95c 100644 --- a/pyscf/df/autoaux.py +++ b/pyscf/df/autoaux.py @@ -53,9 +53,18 @@ def _primitive_emin_emax(basis): r_exp = np.einsum('pi,pq,qi->i', cs, r_ints, cs) k = 2**(2*l+1) * factorial(l+1)**2 / factorial(2*l+2) - # Eq (9) in the paper. r_exp**2 in the denominator may be a bug in bse + # Eq (9) in the paper #e_eff = 2 * k**2 / (np.pi * r_exp) + # r_exp**2 is applied in bse. For primitive functions, this leads to + # e_eff = exponent of the basis e_eff = 2 * k**2 / (np.pi * r_exp**2) + # For primitive functions, e_eff may be slightly different to the + # exponent due to the rounding errors in gaussian_int function. + # When a particular shell has only one primitive function, one auxiliary + # function should be generated. This error can introduce an additional + # auxiliary function. + # Slightly reduce e_eff to remove the extra auxiliary functions. + e_eff -= 1e-8 e_eff_by_l[l] = max(e_eff.max(), e_eff_by_l[l]) return np.array(emax_by_l), np.array(emin_by_l), np.array(e_eff_by_l) @@ -131,6 +140,9 @@ def expand(symb): raise RuntimeError(f'Failed to generate even-tempered auxbasis for {symb}') uniq_atoms = {a[0] for a in mol._atom} + if bse.basis_set_exchange is None: + return {symb: expand(symb) for symb in uniq_atoms} + if isinstance(mol.basis, str): try: elements = [gto.charge(symb) for symb in uniq_atoms] @@ -160,6 +172,10 @@ def autoabs(mol): http://doi.org/10.1063/1.2752807 ''' from pyscf.gto.basis import bse + if bse is None: + print('Package basis-set-exchange not available') + raise ImportError + uniq_atoms = {a[0] for a in mol._atom} if isinstance(mol.basis, str): elements = [gto.charge(symb) for symb in uniq_atoms] diff --git a/pyscf/df/test/test_addons.py b/pyscf/df/test/test_addons.py index 32046f1b4c..6ba601961e 100644 --- a/pyscf/df/test/test_addons.py +++ b/pyscf/df/test/test_addons.py @@ -35,12 +35,19 @@ def test_aug_etb(self): basis = 'cc-pvdz', ) auxbasis = df.addons.aug_etb(mol) - self.assertEqual(len(auxbasis['O']), 42) - self.assertEqual(len(auxbasis['H']), 13) + if df.addons.USE_VERSION_26_AUXBASIS: + self.assertEqual(len(auxbasis['O']), 36) + self.assertEqual(len(auxbasis['H']), 12) + else: + self.assertEqual(len(auxbasis['O']), 42) + self.assertEqual(len(auxbasis['H']), 13) auxbasis = df.addons.autoaux(mol) self.assertEqual(len(auxbasis['O']), 37) - self.assertEqual(len(auxbasis['H']), 13) + # If basis-set-exchange pacakge is installed, this test will fail. + # bse-0.9 produces 13 shells due to round-off errors. + # The correct number of generated basis functions should be 12. + self.assertEqual(len(auxbasis['H']), 12) mol = gto.M( verbose = 0, @@ -50,8 +57,12 @@ def test_aug_etb(self): basis = ('cc-pvdz', [[4, (1., 1.)]]) ) auxbasis = df.addons.aug_etb(mol) - self.assertEqual(len(auxbasis['O']), 59) - self.assertEqual(len(auxbasis['H']), 16) + if df.addons.USE_VERSION_26_AUXBASIS: + self.assertEqual(len(auxbasis['O']), 36) + self.assertEqual(len(auxbasis['H']), 12) + else: + self.assertEqual(len(auxbasis['O']), 59) + self.assertEqual(len(auxbasis['H']), 16) auxbasis = df.addons.autoaux(mol) self.assertEqual(len(auxbasis['O']), 47) @@ -80,8 +91,12 @@ def test_make_auxbasis(self): 'O': ('631g', [[0, 0, (1., 1.)]])} ) auxbasis = df.addons.make_auxbasis(mol) - self.assertEqual(len(auxbasis['O']), 35) - self.assertEqual(len(auxbasis['H']), 3) + if df.addons.USE_VERSION_26_AUXBASIS: + self.assertEqual(len(auxbasis['O']), 32) + self.assertEqual(len(auxbasis['H']), 3) + else: + self.assertEqual(len(auxbasis['O']), 35) + self.assertEqual(len(auxbasis['H']), 3) def test_default_auxbasis(self): mol = gto.M(atom='He 0 0 0; O 0 0 1', basis='ccpvdz')