Autogenerated protonic basis

Author: susilehtola

src/contrib/hneoci.cpp

@@ -106,7 +106,34 @@ int main_guarded(int argc, char **argv) {
   BasisSetLibrary baslib;
   baslib.load_basis(settings.get_string("Basis"));
   BasisSetLibrary pbaslib;
-  pbaslib.load_basis(settings.get_string("ProtonBasis"));
+  if(settings.get_string("ProtonBasis").size())
+    pbaslib.load_basis(settings.get_string("ProtonBasis"));
+  else {
+    // Automatic basis. Get element
+    ElementBasisSet orbbas(baslib.get_element("H"));
+    orbbas.decontract();
+
+    // New basis
+    ElementBasisSet pbas("H", 0);
+    // Loop over angular momentum
+    for(int am=0;am<=orbbas.get_max_am();am++) {
+      arma::vec exps;
+      arma::mat coeffs;
+      orbbas.get_primitives(exps,coeffs,am);
+
+      // Scale exponents
+      exps *= sqrt(proton_mass);
+
+      for(auto expn: exps) {
+        FunctionShell fn(am);
+        fn.add_exponent(1.0,expn);
+        pbas.add_function(fn);
+      }
+    }
+    pbaslib.add_element(pbas);
+
+    printf("Autogenerated protonic basis by scaling exponents by %e\n",sqrt(proton_mass));
+  }
 
   std::vector<atom_t> atoms(1);
   atoms[0].el="H";