diff --git a/poet/include/PhreeqcEngine.hpp b/poet/include/PhreeqcEngine.hpp index 3b131aa6..fd455d90 100644 --- a/poet/include/PhreeqcEngine.hpp +++ b/poet/include/PhreeqcEngine.hpp @@ -17,18 +17,27 @@ #include #include #include -#include -#include #include -#include #include -#include #include -#include #include +/** + * @brief Class for running Phreeqc wrappped in POET + * + * Direct interface to Phreeqc, without utilizing *_MODIFY keywords/scripts to + * set new values. Use with already initialized Phreeqc config. + * + */ class PhreeqcEngine : public IPhreeqc { public: + /** + * @brief Construct a new Phreeqc Engine object + * + * Construct a new Phreeqc Engine object by previously initialized POETConfig. + * + * @param config Holds the configuration for the Phreeqc engine. + */ PhreeqcEngine(const POETConfig &config) : IPhreeqc() { this->LoadDatabaseString(config.database.c_str()); this->RunString(config.input_script.c_str()); @@ -36,27 +45,19 @@ public: this->init_wrappers(config.cell); } + /** + * @brief Siimulate a cell for a given time step + * + * @param cell_values Vector containing the input values for the cell + * (**including the ID**). Output values are written back in place to this + * vector! + * @param time_step Time step to simulate in seconds + */ void runCell(std::vector &cell_values, double time_step); - void run(int start_cell, int end_cell) { - const std::string runs_string = "RUN_CELLS\n -cells " + - std::to_string(start_cell) + "-" + - std::to_string(end_cell) + "\nEND\n"; - this->RunString(runs_string.c_str()); - } - - void run(int start_cell, int end_cell, double time_step) { - std::stringstream time_ss; - time_ss << std::fixed << std::setprecision(20) << time_step; - - const std::string runs_string = - "RUN_CELLS\n -cells " + std::to_string(start_cell) + "-" + - std::to_string(end_cell) + "\n -time_step " + time_ss.str() + "\nEND\n"; - this->RunString(runs_string.c_str()); - } - private: void init_wrappers(const POETInitCell &cell); + void run(double time_step); cxxSolution *Get_solution(std::size_t n) { return Utilities::Rxn_find(this->PhreeqcPtr->Get_Rxn_solution_map(), n); diff --git a/poet/include/PhreeqcInit.hpp b/poet/include/PhreeqcInit.hpp index b81b2d56..100992d4 100644 --- a/poet/include/PhreeqcInit.hpp +++ b/poet/include/PhreeqcInit.hpp @@ -19,34 +19,88 @@ enum { POET_SOL = 0, POET_EXCH, POET_KIN, POET_EQUIL, POET_SURF }; +/** + * @brief Class for initializing Phreeqc and create input for POET + * + * This class is used to initialize Phreeqc by reading a database and an input + * Phreeqc script. It then extracts the 'plain' initial values for each cell and + * module, which are then used for further processing in POET. + */ class PhreeqcInit : public IPhreeqc { public: + /** + * @brief Construct a new Phreeqc Init object + * + * @param database String containing the content of the database file + * @param input_script String containing the content of the input script + */ PhreeqcInit(const std::string &database, const std::string &input_script); + /** + * Hold column names, indexes and values for each cell defined in the input + * script + */ struct PhreeqcMat { std::vector names; std::vector ids; std::vector> values; }; + /** + * @brief Get the PhreeqcMatrix + * + * @return PhreeqcMat Value matrix of the parsed Phreeqc input script + */ PhreeqcMat getPhreeqcMat() const { return pqc_mat; }; + /** + * @brief Get the *_RAW input scripts for each cell + * + * @return std::map Contains the raw input scripts for each + * cell + */ std::map raw_dumps(); using essential_names = std::array, 5>; using ModulesArray = std::array; + /** + * @brief Get the actual reactant sizes for given cells + * + * @param cell_ids Vector of cell ids to return reactant sizes for + * @return ModulesArray Combined sizes of all cells for each module + */ ModulesArray getModuleSizes(const std::vector &cell_ids); + /** + * @brief Get solution primaries found by PhreeqcInit + * + * @return std::vector Vector containing only the solution + * primary names of the solution defined by the input script + */ std::vector getSolutionPrimaries() { return std::vector(this->surface_primaries.begin(), this->surface_primaries.end()); } + /** + * @brief Get the solution names for a given cell + * + * @param cell_id ID of the cell to get the solution names for + * @return std::vector Whole vector of solution names for the + * cell + */ std::vector getSolutionNames(int cell_id) { return this->raw_initials[cell_id].first[POET_SOL]; } + /** + * @brief Get the exchange names for a given cell + * + * @param cell_id ID of the cell to get the exchange names for + * @return std::vector Whole vector of exchange names for the + * cell. Empty if no exchange is defined. + */ std::vector getExchanger(int id) { if (this->exchanger.contains(id)) { return this->exchanger[id]; @@ -55,6 +109,13 @@ public: return {}; } + /** + * @brief Get the kinetics names for a given cell + * + * @param cell_id ID of the cell to get the kinetics names for + * @return std::vector Whole vector of kinetics names for the + * cell. Empty if no kinetics are defined. + */ std::vector getKineticsNames(int id) { if (this->kinetics.contains(id)) { return this->kinetics[id]; @@ -63,6 +124,13 @@ public: return {}; } + /** + * @brief Get the equilibrium names for a given cell + * + * @param cell_id ID of the cell to get the equilibrium names for + * @return std::vector Whole vector of equilibrium names for the + * cell. Empty if no equilibrium is defined. + */ std::vector getEquilibriumNames(int id) { if (this->equilibrium.contains(id)) { return this->equilibrium[id]; @@ -71,6 +139,13 @@ public: return {}; } + /** + * @brief Get the surface component names for a given cell + * + * @param cell_id ID of the cell to get the surface component names for + * @return std::vector Whole vector of surface component names + * for the cell. Empty if no surface is defined. + */ std::vector getSurfaceCompNames(int id) { if (this->surface_comps.contains(id)) { return this->surface_comps[id]; @@ -79,6 +154,13 @@ public: return {}; } + /** + * @brief Get the surface charge names for a given cell + * + * @param cell_id ID of the cell to get the surface charge names for + * @return std::vector Whole vector of surface charge names for + * the cell. Empty if no surface is defined. + */ std::vector getSurfaceChargeNames(int id) { if (this->surface_charge.contains(id)) { return this->surface_charge[id]; diff --git a/poet/src/Engine.cpp b/poet/src/Engine.cpp index 2f327bbf..538cf435 100644 --- a/poet/src/Engine.cpp +++ b/poet/src/Engine.cpp @@ -1,6 +1,8 @@ #include "PhreeqcEngine.hpp" #include +#include #include +#include void PhreeqcEngine::init_wrappers(const POETInitCell &cell) { @@ -215,6 +217,15 @@ void PhreeqcEngine::runCell(std::vector &cell_values, std::span cell_data{cell_values.begin() + 1, cell_values.end()}; this->set_essential_values(cell_data); - this->run(1, 1, time_step); + this->run(time_step); this->get_essential_values(cell_data); +} + +void PhreeqcEngine::run(double time_step) { + std::stringstream time_ss; + time_ss << std::fixed << std::setprecision(20) << time_step; + + const std::string runs_string = + "RUN_CELLS\n -cells 1\n -time_step " + time_ss.str() + "\nEND\n"; + this->RunString(runs_string.c_str()); } \ No newline at end of file