more nactive options

ExaChem · Jan 10, 2025 · 860166f · 860166f
1 parent dc30374
commit 860166f
Show file tree

Hide file tree

Showing 24 changed files with 87 additions and 56 deletions.
diff --git a/ci/reference_output/co.cc-pvdz.gfccsd.json b/ci/reference_output/co.cc-pvdz.gfccsd.json
@@ -540,7 +540,6 @@
       "writet": "true",
       "writet_iter": 5,
       "ccsd_maxiter": 50,
-      "nactive": 0,
       "debug": false,
       "profile_ccsd": false,
       "balance_tiles": "true",

diff --git a/ci/reference_output/h2o_ducc.cc-pvdz.ducc.json b/ci/reference_output/h2o_ducc.cc-pvdz.ducc.json
@@ -272,7 +272,7 @@
       "balance_tiles": "true"
     },
     "DUCC": {
-      "nactive": 2
+      "nactive_va": 2
     }
   }
 }
diff --git a/docs/schema/input_schema.json b/docs/schema/input_schema.json
@@ -296,7 +296,16 @@
         "ccsd_maxiter": {
           "type": "number"
         },
-        "nactive": {
+        "nactive_oa": {
+          "type": "number"
+        },
+        "nactive_ob": {
+          "type": "number"
+        },
+        "nactive_va": {
+          "type": "number"
+        },
+        "nactive_vb": {
           "type": "number"
         },
         "ducc_lvl": {

diff --git a/docs/user_guide/coupledcluster.rst b/docs/user_guide/coupledcluster.rst
@@ -121,11 +121,17 @@ The double unitary CC formalism (DUCC) is described in the following paper.
 .. code-block:: json
 
  "CC": {
-   "nactive"  : 0,
+   "nactive_oa"  : 0,
+   "nactive_ob"  : 0,
+   "nactive_va"  : 0,
+   "nactive_vb"  : 0,
    "ducc_lvl" : 2
  }
 
-:nactive: ``[default=0]`` An integer that specifies the number of active virtual orbitals included in the Hamiltonian H. All occupied orbitals are included.
+:nactive_oa: ``[default=0]`` An integer that specifies the number of active occupied alpha orbitals included in the Hamiltonian H.
+:nactive_ob: ``[default=0]`` An integer that specifies the number of active occupied beta orbitals included in the Hamiltonian H.
+:nactive_va: ``[default=0]`` An integer that specifies the number of active virtual alpha orbitals included in the Hamiltonian H.
+:nactive_vb: ``[default=0]`` An integer that specifies the number of active virtual beta orbitals included in the Hamiltonian H.
 
 :ducc_lvl: ``[default=2]`` An integer that specifies the level of DUCC theory.
 

diff --git a/exachem/cc/ducc/ducc-t_ccsd.cpp b/exachem/cc/ducc/ducc-t_ccsd.cpp
@@ -1623,7 +1623,10 @@ bool ducc_tensors_exist(std::vector<Tensor<T>>& rwtensors, std::vector<std::stri
 
 void reset_ducc_runcontext(ExecutionContext& ec, ChemEnv& chem_env) {
   CCSDOptions& ccsd_options = chem_env.ioptions.ccsd_options;
-  chem_env.run_context["ducc"]["nactive"]  = ccsd_options.nactive;
+  chem_env.run_context["ducc"]["nactive_oa"]  = ccsd_options.nactive_oa;
+  chem_env.run_context["ducc"]["nactive_ob"]  = ccsd_options.nactive_ob;
+  chem_env.run_context["ducc"]["nactive_va"]  = ccsd_options.nactive_va;
+  chem_env.run_context["ducc"]["nactive_vb"]  = ccsd_options.nactive_vb;
   chem_env.run_context["ducc"]["ducc_lvl"] = ccsd_options.ducc_lvl;
   chem_env.run_context["ducc"]["level0"]   = false;
   chem_env.run_context["ducc"]["level1"]   = false;
@@ -1659,7 +1662,7 @@ void DUCC_T_CCSD_Driver(ChemEnv& chem_env, ExecutionContext& ec, const TiledInde
 
   SystemData& sys_data = chem_env.sys_data;
   const size_t nelectrons_alpha = sys_data.nelectrons_alpha;
-  const size_t nactv = chem_env.ioptions.ccsd_options.nactive;
+  const size_t nactv = chem_env.ioptions.ccsd_options.nactive_va;
   const int ducc_lvl       = chem_env.ioptions.ccsd_options.ducc_lvl;
 
   const TiledIndexSpace& O = MO("occ");
@@ -1778,7 +1781,7 @@ void DUCC_T_CCSD_Driver(ChemEnv& chem_env, ExecutionContext& ec, const TiledInde
   bool      dtensors_exist = ducc_tensors_exist(ducc_tensors, dt_files, drestart);
 
   if(drestart && dtensors_exist) {
-    const size_t orig_nactv = chem_env.run_context["ducc"]["nactive"];
+    const size_t orig_nactv = chem_env.run_context["ducc"]["nactive_va"];
     if(orig_nactv != nactv) { dtensors_exist = false; }
     const int orig_ducc_lvl = chem_env.run_context["ducc"]["ducc_lvl"];
     if(orig_ducc_lvl > ducc_lvl) { dtensors_exist = false; }
@@ -2146,7 +2149,7 @@ void DUCC_T_CCSD_Driver(ChemEnv& chem_env, ExecutionContext& ec, const TiledInde
   // TODO: FCIDUMP
   // if(sys_data.nbf_orig <= 20) {
   //   int              nfrzc = 0;
-  //   std::vector<int> symvec(chem_env.ioptions.ccsd_options.nactive);
+  //   std::vector<int> symvec(chem_env.ioptions.ccsd_options.nactive_va);
   //   std::fill(symvec.begin(), symvec.end(), 1);
 
   //   double            enuc = sys_data.results["output"]["SCF"]["nucl_rep_energy"];

diff --git a/exachem/cholesky/cholesky_2e.cpp b/exachem/cholesky/cholesky_2e.cpp
@@ -125,17 +125,23 @@ std::tuple<TiledIndexSpace, TAMM_SIZE> setup_mo_red(ExecutionContext& ec, ChemEn
 }
 
 std::tuple<TiledIndexSpace, TAMM_SIZE> setupMOIS(ExecutionContext& ec, ChemEnv& chem_env,
-                                                 bool triples, int nactv) {
+                                                 bool triples) {
   SystemData& sys_data    = chem_env.sys_data;
   TAMM_SIZE   n_occ_alpha = sys_data.n_occ_alpha;
   TAMM_SIZE   n_occ_beta  = sys_data.n_occ_beta;
 
   const int rank         = ec.pg().rank().value();
   auto      ccsd_options = chem_env.ioptions.ccsd_options;
 
+  const int nactv = chem_env.ioptions.ccsd_options.nactive_va;
+
   Tile tce_tile      = ccsd_options.tilesize;
   bool balance_tiles = ccsd_options.balance_tiles;
 
+  // TODO: Implement check for UHF
+  if(nactv > sys_data.n_vir_alpha && sys_data.is_restricted)
+    tamm_terminate("[DUCC ERROR]: nactive_va > n_vir_alpha");
+
   const std::string jkey    = "tilesize";
   const bool        user_ts = chem_env.jinput["CC"].contains(jkey) ? true : false;
 
@@ -281,14 +287,7 @@ void update_sysdata(ExecutionContext& ec, ChemEnv& chem_env, TiledIndexSpace& MO
     }
     sys_data.update();
     if(!is_mso) std::tie(MO, total_orbitals) = setup_mo_red(ec, chem_env);
-    else {
-      const int nactive = chem_env.ioptions.ccsd_options.nactive;
-      if(nactive > 0) {
-        // Only DUCC uses this right now
-        std::tie(MO, total_orbitals) = cholesky_2e::setupMOIS(ec, chem_env, false, nactive);
-      }
-      else std::tie(MO, total_orbitals) = cholesky_2e::setupMOIS(ec, chem_env);
-    }
+    else { std::tie(MO, total_orbitals) = cholesky_2e::setupMOIS(ec, chem_env); }
   }
 }
 

diff --git a/exachem/cholesky/cholesky_2e.hpp b/exachem/cholesky/cholesky_2e.hpp
@@ -30,7 +30,7 @@ std::tuple<TiledIndexSpace, TAMM_SIZE> setup_mo_red(ExecutionContext& ec, ChemEn
                                                     bool triples = false);
 
 std::tuple<TiledIndexSpace, TAMM_SIZE> setupMOIS(ExecutionContext& ec, ChemEnv& chem_env,
-                                                 bool triples = false, int nactv = 0);
+                                                 bool triples = false);
 
 void update_sysdata(ExecutionContext& ec, ChemEnv& chem_env, TiledIndexSpace& MO,
                     bool is_mso = true);

diff --git a/exachem/cholesky/cholesky_2e_driver.cpp b/exachem/cholesky/cholesky_2e_driver.cpp
@@ -14,7 +14,7 @@ using namespace exachem::scf;
 void exachem::cholesky_2e::cholesky_2e_driver(ExecutionContext& ec, ChemEnv& chem_env) {
   using T = double;
 
-  if(!chem_env.scf_context.scf_converged) scf::scf_driver(ec, chem_env);
+  if(!chem_env.scf_context.skip_scf) scf::scf_driver(ec, chem_env);
 
   Tensor<T> C_AO      = chem_env.scf_context.C_AO;
   Tensor<T> C_beta_AO = chem_env.scf_context.C_beta_AO;
@@ -49,16 +49,8 @@ void exachem::cholesky_2e::cholesky_2e_driver(ExecutionContext& ec, ChemEnv& che
   // TODO: MO here is MSO, rename MO to MSO
   TiledIndexSpace MO;
   TAMM_SIZE       total_orbitals;
-  if(ccsd_options.nactive > 0) {
-    // Only DUCC uses this right now
-    std::tie(MO, total_orbitals) =
-      cholesky_2e::setupMOIS(ec, chem_env, false, ccsd_options.nactive);
-    // TODO: Implement check for UHF
-    if(ccsd_options.nactive > sys_data.n_vir_alpha && sys_data.is_restricted)
-      tamm_terminate("[DUCC ERROR]: nactive > n_vir_alpha");
-  }
-  else std::tie(MO, total_orbitals) = cholesky_2e::setupMOIS(ec, chem_env);
-  chem_env.is_context.MSO = MO;
+  std::tie(MO, total_orbitals) = cholesky_2e::setupMOIS(ec, chem_env);
+  chem_env.is_context.MSO      = MO;
 
   if(ccsd_options.skip_ccsd) {
     Tensor<T>::deallocate(C_AO, F_AO);

diff --git a/exachem/common/chemenv.cpp b/exachem/common/chemenv.cpp
@@ -176,7 +176,10 @@ void ChemEnv::write_json_data(const std::string cmodule) {
     results["input"][cmodule]["writet"]        = str_bool(ccsd.writet);
     results["input"][cmodule]["writet_iter"]   = ccsd.writet_iter;
     results["input"][cmodule]["ccsd_maxiter"]  = ccsd.ccsd_maxiter;
-    results["input"][cmodule]["nactive"]       = ccsd.nactive;
+    results["input"][cmodule]["nactive_oa"]    = ccsd.nactive_oa;
+    results["input"][cmodule]["nactive_ob"]    = ccsd.nactive_ob;
+    results["input"][cmodule]["nactive_va"]    = ccsd.nactive_va;
+    results["input"][cmodule]["nactive_vb"]    = ccsd.nactive_vb;
     results["input"][cmodule]["debug"]         = ccsd.debug;
     results["input"][cmodule]["profile_ccsd"]  = ccsd.profile_ccsd;
     results["input"][cmodule]["balance_tiles"] = str_bool(ccsd.balance_tiles);
@@ -212,7 +215,10 @@ void ChemEnv::write_json_data(const std::string cmodule) {
 
   if(cmodule == "DUCC") {
     // DUCC options
-    results["input"]["DUCC"]["nactive"] = ccsd.nactive;
+    results["input"]["DUCC"]["nactive_oa"] = ccsd.nactive_oa;
+    results["input"]["DUCC"]["nactive_ob"] = ccsd.nactive_ob;
+    results["input"]["DUCC"]["nactive_va"] = ccsd.nactive_va;
+    results["input"]["DUCC"]["nactive_vb"] = ccsd.nactive_vb;
   }
 
   if(cmodule == "EOMCCSD") {

diff --git a/exachem/common/context/scf_context.hpp b/exachem/common/context/scf_context.hpp
@@ -23,7 +23,8 @@ class SCFContext {
   bool                     no_scf;
   bool                     do_df;
 
-  bool scf_converged{false};
+  // bool scf_converged{false};
+  bool skip_scf{false};
 
   void update(double hf_energy, std::vector<size_t> shell_tile_map, tamm::Tensor<TensorType> C_AO,
               tamm::Tensor<TensorType> F_AO, tamm::Tensor<TensorType> C_beta_AO,
@@ -36,4 +37,4 @@ class SCFContext {
     this->F_beta_AO      = F_beta_AO;
     this->no_scf         = no_scf;
   }
-};
+};
diff --git a/exachem/common/options/input_options.cpp b/exachem/common/options/input_options.cpp
@@ -126,7 +126,7 @@ void CCSDOptions::print() {
   std::cout << " ndiis                = " << ndiis << std::endl;
   std::cout << " threshold            = " << threshold << std::endl;
   std::cout << " tilesize             = " << tilesize << std::endl;
-  if(nactive > 0) std::cout << " nactive              = " << nactive << std::endl;
+  // if(nactive > 0) std::cout << " nactive              = " << nactive << std::endl;
   if(pcore > 0) std::cout << " pcore                = " << pcore << std::endl;
   std::cout << " ccsd_maxiter         = " << ccsd_maxiter << std::endl;
   txt_utils::print_bool(" freeze_atomic       ", freeze_atomic);
@@ -226,7 +226,10 @@ void CCSDOptions::initialize() {
   tilesize       = 40;
   ndiis          = 5;
   lshift         = 0;
-  nactive        = 0;
+  nactive_oa     = 0;
+  nactive_ob     = 0;
+  nactive_va     = 0;
+  nactive_vb     = 0;
   ducc_lvl       = 2;
   ccsd_maxiter   = 50;
   freeze_core    = 0;

diff --git a/exachem/common/options/input_options.hpp b/exachem/common/options/input_options.hpp
@@ -116,7 +116,7 @@ class CCSDOptions: public CommonOptions {
   std::pair<bool, double> tamplitudes{false, 0.05};
   std::vector<int>        cc_rdm{};
 
-  int  nactive;
+  int  nactive_oa, nactive_ob, nactive_va, nactive_vb;
   int  ducc_lvl;
   int  ccsd_maxiter;
   bool freeze_atomic{false};

diff --git a/exachem/common/options/parse_ccsd_options.cpp b/exachem/common/options/parse_ccsd_options.cpp
@@ -25,7 +25,8 @@ void ParseCCSDOptions::parse_check(json& jinput) {
       "comments", "threshold",    "tilesize",     "ducc_lvl",
       "lshift",   "ndiis",        "ccsd_maxiter",   "freeze",       "PRINT",
       "readt",    "writet",       "writet_iter",  "debug",
-      "nactive",  "profile_ccsd", "balance_tiles",  "ext_data_path"};
+      "nactive_oa", "nactive_ob", "nactive_va", "nactive_vb", 
+      "profile_ccsd", "balance_tiles",  "ext_data_path"};
   // clang-format on
   for(auto& el: jinput["CC"].items()) {
     if(std::find(valid_cc.begin(), valid_cc.end(), el.key()) == valid_cc.end())
@@ -39,7 +40,10 @@ void ParseCCSDOptions::parse(ChemEnv& chem_env) {
   update_common_options(chem_env);
 
   parse_option<int>(cc_options.ndiis, jcc, "ndiis");
-  parse_option<int>(cc_options.nactive, jcc, "nactive");
+  parse_option<int>(cc_options.nactive_oa, jcc, "nactive_oa");
+  parse_option<int>(cc_options.nactive_ob, jcc, "nactive_ob");
+  parse_option<int>(cc_options.nactive_va, jcc, "nactive_va");
+  parse_option<int>(cc_options.nactive_vb, jcc, "nactive_vb");
   parse_option<int>(cc_options.ducc_lvl, jcc, "ducc_lvl");
   parse_option<int>(cc_options.ccsd_maxiter, jcc, "ccsd_maxiter");
   parse_option<double>(cc_options.lshift, jcc, "lshift");

diff --git a/exachem/scf/scf_hartree_fock.cpp b/exachem/scf/scf_hartree_fock.cpp
@@ -1237,7 +1237,7 @@ void exachem::scf::SCFHartreeFock::scf_hf(ExecutionContext& exc, ChemEnv& chem_e
   chem_env.is_context.AO_tis   = scf_vars.tAOt;
   chem_env.is_context.AO_ortho = tAO_ortho;
 
-  chem_env.scf_context.scf_converged = true;
+  // chem_env.scf_context.scf_converged = true;
 
   chem_env.scf_context.update(ehf, scf_vars.shell_tile_map, C_alpha_tamm, Fa_global, C_beta_tamm,
                               Fb_global, chem_env.ioptions.scf_options.noscf);

diff --git a/exachem/task/ec_task.cpp b/exachem/task/ec_task.cpp
@@ -91,7 +91,8 @@ void execute_task(ExecutionContext& ec, ChemEnv& chem_env, std::string ec_arg2)
   chem_env.shells           = chem_env.ec_basis.shells;
   chem_env.sys_data.has_ecp = chem_env.ec_basis.has_ecp;
 
-  if(chem_env.atoms.size() <= 30) exachem::task::geometry_analysis(ec, chem_env);
+  if(chem_env.atoms.size() > 1 && chem_env.atoms.size() <= 30)
+    exachem::task::geometry_analysis(ec, chem_env);
 
   print_geometry(ec, chem_env);
 

diff --git a/exachem/task/numerical_gradients.cpp b/exachem/task/numerical_gradients.cpp
@@ -1,7 +1,7 @@
 /*
  * ExaChem: Open Source Exascale Computational Chemistry Software.
  *
- * Copyright 2023-2024 Pacific Northwest National Laboratory, Battelle Memorial Institute.
+ * Copyright 2023-2025 Pacific Northwest National Laboratory, Battelle Memorial Institute.
  *
  * See LICENSE.txt for details
  */
@@ -25,8 +25,8 @@ double get_task_energy(ExecutionContext& ec, ChemEnv& chem_env) {
 }
 
 double compute_energy(ExecutionContext& ec, ChemEnv& chem_env, std::string ec_arg2) {
+  // std::cout << "computing energy" << std::endl;
   execute_task(ec, chem_env, ec_arg2);
-  chem_env.scf_context.scf_converged = false;
   return get_task_energy(ec, chem_env);
 }
 
@@ -145,4 +145,4 @@ Matrix compute_gradients(ExecutionContext& ec, ChemEnv& chem_env, const std::vec
   return gradients;
 }
 
-} // namespace exachem::task
+} // namespace exachem::task
diff --git a/exachem/task/numerical_gradients.hpp b/exachem/task/numerical_gradients.hpp
@@ -1,7 +1,7 @@
 /*
  * ExaChem: Open Source Exascale Computational Chemistry Software.
  *
- * Copyright 2023-2024 Pacific Northwest National Laboratory, Battelle Memorial Institute.
+ * Copyright 2023-2025 Pacific Northwest National Laboratory, Battelle Memorial Institute.
  *
  * See LICENSE.txt for details
  */
@@ -18,4 +18,4 @@ Matrix compute_gradients(ExecutionContext& ec, ChemEnv& chem_env, const std::vec
 
 double compute_energy(ExecutionContext& ec, ChemEnv& chem_env, std::string ec_arg2);
 
-} // namespace exachem::task
+} // namespace exachem::task
diff --git a/inputs/aminoprop.json b/inputs/aminoprop.json
@@ -62,7 +62,6 @@
     "lshift": 0.0,
     "ndiis": 5,
     "ccsd_maxiter": 100,
-    "nactive": 0,
     "freeze":{
       "atomic":false,
       "core": 0,

diff --git a/inputs/caffeine.json b/inputs/caffeine.json
@@ -45,7 +45,7 @@
   },
   "CC": {
     "threshold": 1e-6,
-    "nactive": 36,
+    "nactive_va": 36,
     "CCSD(T)": {
       "ccsdt_tilesize": 28
     },
@@ -61,4 +61,4 @@
     "ccsd": false,
     "ccsd_t": false
   }  
-}
+}
diff --git a/inputs/ci/h2o_ducc.json b/inputs/ci/h2o_ducc.json
@@ -38,8 +38,13 @@
       "ccsdt_tilesize": 28
     },
     "debug": false,
-    "nactive": 2,
-    "ducc_lvl": 2
+    "nactive_va": 2,
+    "ducc_lvl": 2,
+    "freeze": {
+      "atomic": false,
+      "core": 0,
+      "virtual": 0
+    }    
   },
   "TASK": {
     "ducc": true  

diff --git a/inputs/ci/h2o_eom.json b/inputs/ci/h2o_eom.json
@@ -37,12 +37,12 @@
       "ccsdt_tilesize": 28
     },
     "debug": false,
-    "nactive": 0,
+    "nactive_va": 0,
     "EOMCCSD": {
       "eom_nroots": 4
     }
   },
   "TASK": {
     "eom_ccsd": true
   }  
-}
+}
diff --git a/inputs/example.json b/inputs/example.json
@@ -133,7 +133,10 @@
         "lshift": 0,
         "ndiis": 5,
         "ccsd_maxiter": 50,
-        "nactive": 0,
+        "nactive_oa": 0,
+        "nactive_ob": 0,
+        "nactive_va": 0,
+        "nactive_vb": 0,
         "ducc_lvl": 2,
         "freeze": {
           "atomic": false,