Update autoaux tests

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)

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]

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')