Add chemistry initialization to InitialList class

bench/dolo/het/dolo_200.R

@@ -28,8 +28,23 @@ diffusion_setup <- list(
     alpha_y = 1e-6
 )
 
+chemistry_setup <- list(
+    dht_species = c(
+        "H" = 3,
+        "O" = 3,
+        "Charge" = 3,
+        "C(4)" = 6,
+        "Ca" = 6,
+        "Cl" = 3,
+        "Mg" = 5,
+        "Calcite" = 4,
+        "Dolomite" = 4
+    )
+)
+
 # Define a setup list for simulation configuration
 setup <- list(
     Grid = grid_setup, # Parameters related to the grid structure
-    Diffusion = diffusion_setup # Parameters related to the diffusion process
+    Diffusion = diffusion_setup, # Parameters related to the diffusion process
+    Chemistry = chemistry_setup # Parameters related to the chemistry process
 )

src/Init/ChemistryInit.cpp

@@ -1,6 +1,15 @@
 #include "InitialList.hpp"
 
 namespace poet {
+
+void InitialList::initChemistry(const Rcpp::List &chem) {
+  this->dht_defined = chem.containsElementNamed("dht_species");
+
+  if (this->dht_defined) {
+    this->dht_species = Rcpp::as<NamedVector<uint32_t>>(chem["dht_species"]);
+  }
+}
+
 InitialList::ChemistryInit InitialList::getChemistryInit() const {
   ChemistryInit chem_init;
 
@@ -12,6 +21,9 @@ InitialList::ChemistryInit InitialList::getChemistryInit() const {
 
   chem_init.pqc_sol_order = pqc_sol_order;
 
+  chem_init.dht_defined = dht_defined;
+  chem_init.dht_species = dht_species;
+
   return chem_init;
 }
 } // namespace poet

src/Init/InitialList.cpp

@@ -1,4 +1,5 @@
 #include "InitialList.hpp"
+#include "DataStructures/NamedVector.hpp"
 #include <Rcpp/internal/wrap.h>
 #include <Rcpp/proxy/ProtectedProxy.h>
 #include <Rcpp/vector/instantiation.h>
@@ -10,6 +11,7 @@ namespace poet {
 void InitialList::initializeFromList(const Rcpp::List &setup) {
   initGrid(setup[grid_key]);
   initDiffusion(setup[diffusion_key]);
+  initChemistry(setup[chemistry_key]);
 }
 
 void InitialList::importList(const Rcpp::List &setup) {
@@ -51,6 +53,11 @@ void InitialList::importList(const Rcpp::List &setup) {
       setup[static_cast<int>(ExportList::CHEM_PQC_IDS)]);
   this->pqc_sol_order = Rcpp::as<std::vector<std::string>>(
       setup[static_cast<int>(ExportList::CHEM_PQC_SOL_ORDER)]);
+
+  this->dht_defined =
+      Rcpp::as<bool>(setup[static_cast<int>(ExportList::CHEM_DHT_DEFINED)]);
+  this->dht_species = Rcpp::as<NamedVector<uint32_t>>(
+      setup[static_cast<int>(ExportList::CHEM_DHT_SPECIES)]);
 }
 
 Rcpp::List InitialList::exportList() {
@@ -78,6 +85,10 @@ Rcpp::List InitialList::exportList() {
   out[static_cast<int>(ExportList::CHEM_PQC_IDS)] = Rcpp::wrap(this->pqc_ids);
   out[static_cast<int>(ExportList::CHEM_PQC_SOL_ORDER)] =
       Rcpp::wrap(this->pqc_sol_order);
+  out[static_cast<int>(ExportList::CHEM_DHT_DEFINED)] =
+      Rcpp::wrap(this->dht_defined);
+  out[static_cast<int>(ExportList::CHEM_DHT_SPECIES)] =
+      Rcpp::wrap(this->dht_species);
 
   return out;
 }

src/Init/InitialList.hpp

@@ -1,10 +1,12 @@
 #pragma once
 
+#include "Chemistry/ChemistryDefs.hpp"
+#include "DataStructures/NamedVector.hpp"
+#include <Rcpp/vector/instantiation.h>
 #include <tug/Boundary.hpp>
 
 #include <RInside.h>
 #include <Rcpp.h>
-#include <Rcpp/vector/instantiation.h>
 #include <array>
 #include <cstddef>
 #include <cstdint>
@@ -12,7 +14,6 @@
 #include <IPhreeqcPOET.hpp>
 #include <map>
 #include <string>
-#include <utility>
 #include <vector>
 
 #include <DataStructures/Field.hpp>
@@ -50,6 +51,8 @@ private:
     CHEM_PQC_SCRIPTS,
     CHEM_PQC_IDS,
     CHEM_PQC_SOL_ORDER,
+    CHEM_DHT_DEFINED,
+    CHEM_DHT_SPECIES,
     ENUM_SIZE
   };
 
@@ -150,7 +153,7 @@ private:
   // Chemistry Members
   static constexpr const char *chemistry_key = "Chemistry";
 
-  void initChemistry();
+  void initChemistry(const Rcpp::List &chem_input);
 
   std::string database;
   std::vector<std::string> pqc_scripts;
@@ -158,6 +161,9 @@ private:
 
   std::vector<std::string> pqc_sol_order;
 
+  bool dht_defined;
+  NamedVector<std::uint32_t> dht_species;
+
 public:
   struct ChemistryInit {
     uint32_t total_grid_cells;
@@ -166,6 +172,9 @@ public:
     std::vector<std::string> pqc_scripts;
     std::vector<int> pqc_ids;
     std::vector<std::string> pqc_sol_order;
+
+    bool dht_defined;
+    NamedVector<std::uint32_t> dht_species;
   };
 
   ChemistryInit getChemistryInit() const;