diff --git a/include/gauxc/enums.hpp b/include/gauxc/enums.hpp
index 76d4500c..ed09068d 100644
--- a/include/gauxc/enums.hpp
+++ b/include/gauxc/enums.hpp
@@ -35,7 +35,17 @@ enum class AtomicGridSizeDefault {
   UltraFineGrid,  ///< Ultrafine grid (appropriate accuracy)
   SuperFineGrid,  ///< Superfine grid (most accurate)
   GM3,            ///< Treutler-Ahlrichs GM3
-  GM5             ///< Treutlet-Ahlrichs GM5
+  GM5,            ///< Treutler-Ahlrichs GM5
+  PySCF0,         ///< PySCF default level 0
+  PySCF1,         ///< PySCF default level 1
+  PySCF2,         ///< PySCF default level 2 (angular points ~ fine grid)
+  PySCF3,         ///< PySCF default level 3
+  PySCF4,         ///< PySCF default level 4 (radial points ~ fine grid, angular points ~ ultrafine grid)
+  PySCF5,         ///< PySCF default level 5
+  PySCF6,         ///< PySCF default level 6 (radial points ~ ultrafine grid, angular points ~ superfine grid)
+  PySCF7,         ///< PySCF default level 7
+  PySCF8,         ///< PySCF default level 8
+  PySCF9          ///< PySCF default level 9 (radial points ~ superfine grid)
 };
 
 /**
diff --git a/src/molgrid_defaults.cxx b/src/molgrid_defaults.cxx
index 61b66830..a40858e0 100644
--- a/src/molgrid_defaults.cxx
+++ b/src/molgrid_defaults.cxx
@@ -15,6 +15,40 @@
 
 namespace GauXC {
 
+namespace detail {
+
+inline RadialSize get_pyscf_radial_size( AtomicNumber Z, int level ) {
+  if (level < 0 || level > 8) 
+    GAUXC_GENERIC_EXCEPTION("Invalid PySCF grid level: " + std::to_string(level) + ". Valid levels are 0-8.");
+
+  if ( Z.get() <= 2 ) {
+    if (level == 0) return RadialSize( 10 );
+    return RadialSize( 20 + 10*level );  // ..., 30, 40, 50, ...
+  } else if ( Z.get() <= 10 ) {
+    if (level == 0) return RadialSize( 15 );
+    if (level == 1) return RadialSize( 40 );
+    return RadialSize( 30 + 15 * level );  // ..., 60, 75, 90, ...
+  } else if ( Z.get() <= 18 ) {
+    if (level == 0) return RadialSize( 20 );
+    return RadialSize( 35 + 15 * level );  // ..., 50, 65, 80, ...
+  } else if ( Z.get() <= 36 ) {
+    if (level == 0) return RadialSize( 30 );
+    return RadialSize( 45 + 15 * level );  // ..., 60, 75, 90, ...
+  } else if ( Z.get() <= 54 ) {
+    if (level == 0) return RadialSize( 35 );
+    return RadialSize( 50 + 15 * level );  // ..., 65, 80, 95, ...
+  } else if ( Z.get() <= 86 ) {
+    if (level == 0) return RadialSize( 40 );
+    return RadialSize( 55 + 15 * level );  // ..., 70, 85, 100, ...
+  } else if ( Z.get() <= 118 ) {
+    if (level == 0) return RadialSize( 50 );
+    return RadialSize( 60 + 15 * level );  // ..., 75, 90, 105, ...
+  } else {
+    GAUXC_GENERIC_EXCEPTION("Z > 118 Not Supported for PySCF Grid Defaults");
+  }
+}
+}
+
 RadialScale default_mk_radial_scaling_factor( AtomicNumber _Z ) {
   auto Z = _Z.get();
   switch(Z) {
@@ -117,6 +151,36 @@ std::tuple<RadialSize,AngularSize>
     case AtomicGridSizeDefault::GM5:
       return std::make_tuple( RadialSize(50), AngularSize(302) );
 
+    case AtomicGridSizeDefault::PySCF0:
+      return std::make_tuple( detail::get_pyscf_radial_size(Z, 0), Z.get() <= 2 ? AngularSize(50) : (Z.get() <= 10 ? AngularSize(86) : AngularSize(110)));
+
+    case AtomicGridSizeDefault::PySCF1:
+      return std::make_tuple( detail::get_pyscf_radial_size(Z, 1), Z.get() <= 2 ? AngularSize(110) : AngularSize(194));
+
+    case AtomicGridSizeDefault::PySCF2:
+      return std::make_tuple( detail::get_pyscf_radial_size(Z, 2), Z.get() <= 2 ? AngularSize(194) : AngularSize(302));
+
+    case AtomicGridSizeDefault::PySCF3:
+      return std::make_tuple( detail::get_pyscf_radial_size(Z, 3), Z.get() <= 10 ? AngularSize(302) : AngularSize(434));
+
+    case AtomicGridSizeDefault::PySCF4:
+      return std::make_tuple( detail::get_pyscf_radial_size(Z, 4), Z.get() <= 2 ? AngularSize(434) : AngularSize(590));
+
+    case AtomicGridSizeDefault::PySCF5:
+      return std::make_tuple( detail::get_pyscf_radial_size(Z, 5), Z.get() <= 2 ? AngularSize(590) : AngularSize(770));
+
+    case AtomicGridSizeDefault::PySCF6:
+      return std::make_tuple( detail::get_pyscf_radial_size(Z, 6), Z.get() <= 2 ? AngularSize(770) : AngularSize(974));
+
+    case AtomicGridSizeDefault::PySCF7:
+      return std::make_tuple( detail::get_pyscf_radial_size(Z, 7), Z.get() <= 2 ? AngularSize(974) : AngularSize(1202));
+
+    case AtomicGridSizeDefault::PySCF8:
+      return std::make_tuple( detail::get_pyscf_radial_size(Z, 8), AngularSize(1202));
+
+    case AtomicGridSizeDefault::PySCF9:
+      return std::make_tuple( RadialSize(200), AngularSize(1454) );
+
     case AtomicGridSizeDefault::FineGrid:
       return std::make_tuple( RadialSize(75), AngularSize(302) );
 
diff --git a/tests/molgrid_test.cxx b/tests/molgrid_test.cxx
index 1de1be98..e51379ca 100644
--- a/tests/molgrid_test.cxx
+++ b/tests/molgrid_test.cxx
@@ -219,6 +219,258 @@ TEST_CASE("MolGrid Defaults", "[molgrid]") {
     }
   }
 
+  SECTION("PySCF0 Grid Size") {
+    int64_t refgr[] = { -1,
+      10, 10 ,
+      15, 15, 15, 15, 15, 15, 15, 15 ,
+      20, 20, 20, 20, 20, 20, 20, 20,
+      30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 30, 35, 35,
+      35, 35, 35, 35, 35, 35, 35, 35, 35, 35, 35, 35, 35, 35, 35, 35,
+      40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40, 40,
+      50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50, 50
+    };
+    int64_t refga[] = { -1,
+      50, 50 ,
+      86, 86, 86, 86, 86, 86, 86, 86 ,
+      110, 110, 110, 110, 110, 110, 110, 110,
+      110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110,
+      110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110,
+      110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110,
+      110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110
+    };
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF0 );
+      REQUIRE( gr.get() == refgr[Z] );
+      REQUIRE( ga.get() == refga[Z] );
+    }
+  }
+
+  SECTION("PySCF1 Grid Size") {
+    int64_t refgr[] = { -1,
+      30, 30 ,
+      40, 40, 40, 40, 40, 40, 40, 40 ,
+      50, 50, 50, 50, 50, 50, 50, 50,
+      60, 60, 60, 60, 60, 60, 60, 60, 60, 60, 60, 60, 60, 60, 60, 60, 60, 60, 65, 65,
+      65, 65, 65, 65, 65, 65, 65, 65, 65, 65, 65, 65, 65, 65, 65, 65,
+      70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70,
+      75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75
+    };
+    int64_t refga[] = { -1,
+      110, 110 ,
+      194, 194, 194, 194, 194, 194, 194, 194 ,
+      194, 194, 194, 194, 194, 194, 194, 194,
+      194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194,
+      194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194,
+      194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194,
+      194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194, 194
+    };
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF1 );
+      REQUIRE( gr.get() == refgr[Z] );
+      REQUIRE( ga.get() == refga[Z] );
+    }
+  }
+
+  SECTION("PySCF2 Grid Size") {
+    int64_t refgr[] = { -1,
+      40, 40 ,
+      60, 60, 60, 60, 60, 60, 60, 60 ,
+      65, 65, 65, 65, 65, 65, 65, 65,
+      75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 75, 80, 80,
+      80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80,
+      85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85, 85,
+      90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90
+    };
+    int64_t refga[] = { -1,
+      194, 194 ,
+      302, 302, 302, 302, 302, 302, 302, 302 ,
+      302, 302, 302, 302, 302, 302, 302, 302,
+      302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302,
+      302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302,
+      302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302,
+      302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302, 302
+    };
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF2 );
+      REQUIRE( gr.get() == refgr[Z] );
+      REQUIRE( ga.get() == refga[Z] );
+    }
+  }
+
+  SECTION("PySCF3 Grid Size") {
+    int64_t refgr[] = { -1,
+      50, 50 ,
+      75, 75, 75, 75, 75, 75, 75, 75 ,
+      80, 80, 80, 80, 80, 80, 80, 80,
+      90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 90, 95, 95,
+      95, 95, 95, 95, 95, 95, 95, 95, 95, 95, 95, 95, 95, 95, 95, 95,
+      100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100,
+      105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105
+    };
+    int64_t refga[] = { -1,
+      302, 302 ,
+      302, 302, 302, 302, 302, 302, 302, 302 ,
+      434, 434, 434, 434, 434, 434, 434, 434,
+      434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434,
+      434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434,
+      434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434,
+      434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434, 434
+    };
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF3 );
+      REQUIRE( gr.get() == refgr[Z] );
+      REQUIRE( ga.get() == refga[Z] );
+    }
+  }
+
+  SECTION("PySCF4 Grid Size") {
+    int64_t refgr[] = { -1,
+      60, 60 ,
+      90, 90, 90, 90, 90, 90, 90, 90 ,
+      95, 95, 95, 95, 95, 95, 95, 95,
+      105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 105, 110, 110,
+      110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110, 110,
+      115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115, 115,
+      120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120
+    };
+    int64_t refga[] = { -1,
+      434, 434 ,
+      590, 590, 590, 590, 590, 590, 590, 590 ,
+      590, 590, 590, 590, 590, 590, 590, 590,
+      590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590,
+      590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590,
+      590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590,
+      590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590, 590
+    };
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF4 );
+      REQUIRE( gr.get() == refgr[Z] );
+      REQUIRE( ga.get() == refga[Z] );
+    }
+  }
+
+  SECTION("PySCF5 Grid Size") {
+    int64_t refgr[] = { -1,
+      70, 70 ,
+      105, 105, 105, 105, 105, 105, 105, 105 ,
+      110, 110, 110, 110, 110, 110, 110, 110,
+      120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 120, 125, 125,
+      125, 125, 125, 125, 125, 125, 125, 125, 125, 125, 125, 125, 125, 125, 125, 125,
+      130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130, 130,
+      135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135
+    };
+    int64_t refga[] = { -1,
+      590, 590 ,
+      770, 770, 770, 770, 770, 770, 770, 770 ,
+      770, 770, 770, 770, 770, 770, 770, 770,
+      770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770,
+      770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770,
+      770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770,
+      770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770, 770
+    };
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF5 );
+      REQUIRE( gr.get() == refgr[Z] );
+      REQUIRE( ga.get() == refga[Z] );
+    }
+  }
+
+  SECTION("PySCF6 Grid Size") {
+    int64_t refgr[] = { -1,
+      80, 80 ,
+      120, 120, 120, 120, 120, 120, 120, 120 ,
+      125, 125, 125, 125, 125, 125, 125, 125,
+      135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 135, 140, 140,
+      140, 140, 140, 140, 140, 140, 140, 140, 140, 140, 140, 140, 140, 140, 140, 140,
+      145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145, 145,
+      150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150
+    };
+    int64_t refga[] = { -1,
+      770, 770 ,
+      974, 974, 974, 974, 974, 974, 974, 974 ,
+      974, 974, 974, 974, 974, 974, 974, 974,
+      974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974,
+      974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974,
+      974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974,
+      974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974, 974
+    };
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF6 );
+      REQUIRE( gr.get() == refgr[Z] );
+      REQUIRE( ga.get() == refga[Z] );
+    }
+  }
+
+  SECTION("PySCF7 Grid Size") {
+    int64_t refgr[] = { -1,
+      90, 90 ,
+      135, 135, 135, 135, 135, 135, 135, 135 ,
+      140, 140, 140, 140, 140, 140, 140, 140,
+      150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 150, 155, 155,
+      155, 155, 155, 155, 155, 155, 155, 155, 155, 155, 155, 155, 155, 155, 155, 155,
+      160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160, 160,
+      165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165
+    };
+    int64_t refga[] = { -1,
+      974, 974 ,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202 ,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202
+    };
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF7 );
+      REQUIRE( gr.get() == refgr[Z] );
+      REQUIRE( ga.get() == refga[Z] );
+    }
+  }
+
+  SECTION("PySCF8 Grid Size") {
+    int64_t refgr[] = { -1,
+      100, 100 ,
+      150, 150, 150, 150, 150, 150, 150, 150 ,
+      155, 155, 155, 155, 155, 155, 155, 155,
+      165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 165, 170, 170,
+      170, 170, 170, 170, 170, 170, 170, 170, 170, 170, 170, 170, 170, 170, 170, 170,
+      175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175, 175,
+      180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180, 180
+    };
+    int64_t refga[] = { -1,
+      1202, 1202 ,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202 ,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202,
+      1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202, 1202
+    };
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF8 );
+      REQUIRE( gr.get() == refgr[Z] );
+      REQUIRE( ga.get() == refga[Z] );
+    }
+  }
+
+  SECTION("PySCF9 Grid Size") {
+    for(auto Z = 1; Z < 119; ++Z) {
+      auto [gr, ga] = default_grid_size( AtomicNumber(Z), RadialQuad::MuraKnowles,
+        AtomicGridSizeDefault::PySCF9 );
+      REQUIRE( gr.get() == 200 );
+      REQUIRE( ga.get() == 1454 );
+    }
+  }
+
 }