Refactor API

This commit is contained in:
Max Luebke 2024-03-18 23:05:59 +01:00
parent 54136c8e0c
commit d92ab31765
4 changed files with 312 additions and 186 deletions

View File

@ -14,6 +14,7 @@
#include <map>
#include <string>
#include <sys/types.h>
#include <utility>
#include <vector>
enum { POET_SOL = 0, POET_EXCH, POET_KIN, POET_EQUIL, POET_SURF };
@ -33,42 +34,51 @@ public:
this->LoadDatabase(database.c_str());
this->RunFile(input_script.c_str());
this->parseInitValues();
this->parseInit();
}
std::vector<std::string> getInitNames() const {
std::vector<std::string> names;
IPhreeqcPOET(const std::string &database, const std::string &input_script,
const std::vector<std::string> &solutionInitVector,
std::uint32_t n_cells)
: IPhreeqc(), n_cells(n_cells), solutionInitVector(solutionInitVector) {
this->LoadDatabase(database.c_str());
this->RunFile(input_script.c_str());
for (auto &sol : this->initial_names) {
names.insert(names.end(), sol.begin(), sol.end());
if (n_cells > 1) {
std::string copy_string =
"COPY cell 1 2-" + std::to_string(n_cells) + "\n";
}
return names;
}
std::vector<std::vector<double>> getInitValues() const {
return this->initial_concentrations;
}
struct PhreeqcMat {
std::vector<std::string> names;
std::vector<int> ids;
std::vector<std::vector<double>> values;
};
std::vector<int> getSolutionIds() const { return this->solution_ids; }
PhreeqcMat getPhreeqcMat();
PhreeqcMat getPhreeqcMat(const std::vector<int> &ids);
std::map<int, std::string> raw_dumps() {
std::map<int, std::string> dumps;
this->SetDumpStringOn(true);
for (const auto &sol_id : this->solution_ids) {
for (const auto &[sol_id, unused] : this->raw_initials) {
std::string call_string = "DUMP\n -cells " + std::to_string(sol_id);
this->RunString(call_string.c_str());
dumps[sol_id] = this->GetDumpString();
}
this->SetDumpStringOn(false);
return dumps;
}
using essential_names = std::array<std::vector<std::string>, 5>;
using ModulesArray = std::array<std::uint32_t, 5>;
ModulesArray getModuleSizes() const {
ModulesArray module_sizes;
for (std::uint8_t i = 0; i < 5; i++) {
@ -79,13 +89,13 @@ public:
}
private:
using essential_names = std::array<std::vector<std::string>, 5>;
// required only for simulation
essential_names initial_names;
std::vector<std::vector<double>> initial_concentrations;
std::vector<int> solution_ids;
std::uint32_t n_cells;
std::vector<std::string> solutionInitVector;
void parseInitValues();
void parseInit();
essential_names dump_essential_names(std::size_t cell_number);
@ -110,13 +120,15 @@ private:
return Utilities::Rxn_find(this->PhreeqcPtr->Get_Rxn_surface_map(), n);
}
void resolveSolutionUseKW(
const std::vector<IPhreeqc::SolutionMapping> &unresolved,
std::map<int, std::pair<essential_names, std::vector<double>>>
&mapped_values);
using RawMap = std::map<int, std::pair<essential_names, std::vector<double>>>;
void valuesFromModule(const std::string &module_name, int cell_number,
essential_names &names, std::vector<double> &values);
essential_names union_raws(const RawMap &raws);
std::vector<std::vector<double>>
conc_from_essentials(const RawMap &raws, const essential_names &names);
// void valuesFromModule(const std::string &module_name, int cell_number,
// essential_names &names, std::vector<double> &values);
std::string subExchangeName(std::string name) {
for (const auto &species : this->PhreeqcPtr->Get_species_list()) {
@ -128,4 +140,6 @@ private:
}
return name;
};
RawMap raw_initials;
};

119
poet/src/GetSet.cpp Normal file
View File

@ -0,0 +1,119 @@
#include <IPhreeqcPOET.hpp>
IPhreeqcPOET::essential_names
IPhreeqcPOET::dump_essential_names(std::size_t cell_number) {
essential_names eNames;
// Solutions
if (this->Get_solution(cell_number) != NULL) {
std::vector<std::string> &eSolNames = eNames[POET_SOL];
this->Get_solution(cell_number)->dump_essential_names(eSolNames);
}
// Exchange
if (this->Get_exchange(cell_number) != NULL) {
std::vector<std::string> &eExchNames = eNames[POET_EXCH];
this->Get_exchange(cell_number)->dump_essential_names(eExchNames);
}
// Kinetics
if (this->Get_kinetic(cell_number) != NULL) {
std::vector<std::string> &eKinNames = eNames[POET_KIN];
this->Get_kinetic(cell_number)->dump_essential_names(eKinNames);
}
// PPassemblage
if (this->Get_equilibrium(cell_number) != NULL) {
std::vector<std::string> &eEquNames = eNames[POET_EQUIL];
this->Get_equilibrium(cell_number)->dump_essential_names(eEquNames);
}
// Surface
if (this->Get_surface(cell_number) != NULL) {
std::vector<std::string> &eSurfNames = eNames[POET_SURF];
this->Get_surface(cell_number)->dump_essential_names(eSurfNames);
}
return eNames;
}
std::vector<double>
IPhreeqcPOET::get_essential_values(std::size_t cell_number,
const std::vector<std::string> &order) {
std::vector<double> essentials;
// Solutions
if (this->Get_solution(cell_number) != NULL) {
std::vector<double> sol_values;
this->Get_solution(cell_number)->get_essential_values(sol_values, order);
essentials.insert(essentials.end(), sol_values.begin(), sol_values.end());
}
// Exchange
if (this->Get_exchange(cell_number) != NULL) {
std::vector<double> exch_values;
this->Get_exchange(cell_number)->get_essential_values(exch_values);
essentials.insert(essentials.end(), exch_values.begin(), exch_values.end());
}
// Kinetics
if (this->Get_kinetic(cell_number) != NULL) {
std::vector<double> kin_values;
this->Get_kinetic(cell_number)->get_essential_values(kin_values);
essentials.insert(essentials.end(), kin_values.begin(), kin_values.end());
}
// PPassemblage
if (this->Get_equilibrium(cell_number) != NULL) {
std::vector<double> equ_values;
this->Get_equilibrium(cell_number)->get_essential_values(equ_values);
essentials.insert(essentials.end(), equ_values.begin(), equ_values.end());
}
// Surface
if (this->Get_surface(cell_number) != NULL) {
std::vector<double> surf_values;
this->Get_surface(cell_number)->get_essential_values(surf_values);
essentials.insert(essentials.end(), surf_values.begin(), surf_values.end());
}
return essentials;
}
void IPhreeqcPOET::set_essential_values(std::size_t cell_number,
const std::vector<std::string> &order,
std::vector<double> &values) {
auto dump_it = values.begin();
// Solutions
if (this->Get_solution(cell_number) != NULL) {
this->Get_solution(cell_number)->set_essential_values(dump_it, order);
}
// Exchange
if (this->Get_exchange(cell_number) != NULL) {
this->Get_exchange(cell_number)->set_essential_values(dump_it);
}
// Kinetics
if (this->Get_kinetic(cell_number) != NULL) {
this->Get_kinetic(cell_number)->set_essential_values(dump_it);
}
// PPassemblage
if (this->Get_equilibrium(cell_number) != NULL) {
this->Get_equilibrium(cell_number)->set_essential_values(dump_it);
}
// Surface
if (this->Get_surface(cell_number) != NULL) {
this->Get_surface(cell_number)->set_essential_values(dump_it);
}
}

View File

@ -38,57 +38,58 @@ createConcVector(const std::vector<std::string> &conc_names,
return conc_vec;
}
void IPhreeqcPOET::valuesFromModule(const std::string &module_name,
int cell_number, essential_names &names,
std::vector<double> &values) {
std::size_t dest_module_i = 0;
std::vector<double> to_insert;
if (module_name == "exchange") { // 1
this->Get_exchange(cell_number)->dump_essential_names(names[POET_EXCH]);
this->Get_exchange(cell_number)->get_essential_values(to_insert);
dest_module_i = 1;
} else if (module_name == "kinetics") { // 2
this->Get_kinetic(cell_number)->dump_essential_names(names[POET_KIN]);
this->Get_kinetic(cell_number)->get_essential_values(to_insert);
dest_module_i = 2;
} else if (module_name == "surface") { // 4
this->Get_surface(cell_number)->dump_essential_names(names[POET_SURF]);
this->Get_surface(cell_number)->get_essential_values(to_insert);
dest_module_i = 4;
}
// void IPhreeqcPOET::valuesFromModule(const std::string &module_name,
// int cell_number, essential_names &names,
// std::vector<double> &values) {
// std::size_t dest_module_i = 0;
// std::vector<double> to_insert;
// if (module_name == "exchange") { // 1
// this->Get_exchange(cell_number)->dump_essential_names(names[POET_EXCH]);
// this->Get_exchange(cell_number)->get_essential_values(to_insert);
// dest_module_i = 1;
// } else if (module_name == "kinetics") { // 2
// this->Get_kinetic(cell_number)->dump_essential_names(names[POET_KIN]);
// this->Get_kinetic(cell_number)->get_essential_values(to_insert);
// dest_module_i = 2;
// } else if (module_name == "surface") { // 4
// this->Get_surface(cell_number)->dump_essential_names(names[POET_SURF]);
// this->Get_surface(cell_number)->get_essential_values(to_insert);
// dest_module_i = 4;
// }
std::size_t offset = 0;
for (std::size_t i = 0; i < dest_module_i; i++) {
offset += names[i].size();
}
// std::size_t offset = 0;
// for (std::size_t i = 0; i < dest_module_i; i++) {
// offset += names[i].size();
// }
values.insert(values.begin() + offset, to_insert.begin(), to_insert.end());
}
// values.insert(values.begin() + offset, to_insert.begin(), to_insert.end());
// }
void IPhreeqcPOET::resolveSolutionUseKW(
const std::vector<IPhreeqc::SolutionMapping> &unresolved,
std::map<int, std::pair<essential_names, std::vector<double>>>
&mapped_values) {
// void IPhreeqcPOET::resolveSolutionUseKW(
// const std::vector<IPhreeqc::SolutionMapping> &unresolved,
// std::map<int, std::pair<essential_names, std::vector<double>>>
// &mapped_values) {
for (const auto &input : unresolved) {
if (mapped_values.find(input.module_n) != mapped_values.end()) {
continue;
}
// for (const auto &input : unresolved) {
// if (mapped_values.find(input.module_n) != mapped_values.end()) {
// continue;
// }
essential_names new_conc_names;
new_conc_names[0] = mapped_values[input.sol_n].first[0];
// essential_names new_conc_names;
// new_conc_names[0] = mapped_values[input.sol_n].first[0];
const auto &curr_sol_vec = mapped_values[input.sol_n].second;
std::vector<double> new_conc_values(
curr_sol_vec.begin(), curr_sol_vec.begin() + new_conc_names[0].size());
// const auto &curr_sol_vec = mapped_values[input.sol_n].second;
// std::vector<double> new_conc_values(
// curr_sol_vec.begin(), curr_sol_vec.begin() +
// new_conc_names[0].size());
valuesFromModule(input.module_name, input.module_n, new_conc_names,
new_conc_values);
// valuesFromModule(input.module_name, input.module_n, new_conc_names,
// new_conc_values);
mapped_values[input.module_n] =
std::make_pair(new_conc_names, new_conc_values);
}
}
// mapped_values[input.module_n] =
// std::make_pair(new_conc_names, new_conc_values);
// }
// }
void IPhreeqcPOET::parseInitValues() {
@ -106,8 +107,8 @@ void IPhreeqcPOET::parseInitValues() {
curr_conc_names, this->get_essential_values(id, curr_conc_names[0]));
}
const auto unresolved_modules = this->getSolutionMapping();
resolveSolutionUseKW(unresolved_modules, init_values);
// const auto unresolved_modules = this->getSolutionMapping();
// resolveSolutionUseKW(unresolved_modules, init_values);
std::vector<int> ids_to_erase;
@ -182,122 +183,4 @@ void IPhreeqcPOET::parseInitValues() {
// solution_only, init_values[val.sol_n].second);
// }
// }
}
IPhreeqcPOET::essential_names
IPhreeqcPOET::dump_essential_names(std::size_t cell_number) {
essential_names eNames;
// Solutions
if (this->Get_solution(cell_number) != NULL) {
std::vector<std::string> &eSolNames = eNames[POET_SOL];
this->Get_solution(cell_number)->dump_essential_names(eSolNames);
}
// Exchange
if (this->Get_exchange(cell_number) != NULL) {
std::vector<std::string> &eExchNames = eNames[POET_EXCH];
this->Get_exchange(cell_number)->dump_essential_names(eExchNames);
}
// Kinetics
if (this->Get_kinetic(cell_number) != NULL) {
std::vector<std::string> &eKinNames = eNames[POET_KIN];
this->Get_kinetic(cell_number)->dump_essential_names(eKinNames);
}
// PPassemblage
if (this->Get_equilibrium(cell_number) != NULL) {
std::vector<std::string> &eEquNames = eNames[POET_EQUIL];
this->Get_equilibrium(cell_number)->dump_essential_names(eEquNames);
}
// Surface
if (this->Get_surface(cell_number) != NULL) {
std::vector<std::string> &eSurfNames = eNames[POET_SURF];
this->Get_surface(cell_number)->dump_essential_names(eSurfNames);
}
return eNames;
}
std::vector<double>
IPhreeqcPOET::get_essential_values(std::size_t cell_number,
const std::vector<std::string> &order) {
std::vector<double> essentials;
// Solutions
if (this->Get_solution(cell_number) != NULL) {
std::vector<double> sol_values;
this->Get_solution(cell_number)->get_essential_values(sol_values, order);
essentials.insert(essentials.end(), sol_values.begin(), sol_values.end());
}
// Exchange
if (this->Get_exchange(cell_number) != NULL) {
std::vector<double> exch_values;
this->Get_exchange(cell_number)->get_essential_values(exch_values);
essentials.insert(essentials.end(), exch_values.begin(), exch_values.end());
}
// Kinetics
if (this->Get_kinetic(cell_number) != NULL) {
std::vector<double> kin_values;
this->Get_kinetic(cell_number)->get_essential_values(kin_values);
essentials.insert(essentials.end(), kin_values.begin(), kin_values.end());
}
// PPassemblage
if (this->Get_equilibrium(cell_number) != NULL) {
std::vector<double> equ_values;
this->Get_equilibrium(cell_number)->get_essential_values(equ_values);
essentials.insert(essentials.end(), equ_values.begin(), equ_values.end());
}
// Surface
if (this->Get_surface(cell_number) != NULL) {
std::vector<double> surf_values;
this->Get_surface(cell_number)->get_essential_values(surf_values);
essentials.insert(essentials.end(), surf_values.begin(), surf_values.end());
}
return essentials;
}
void IPhreeqcPOET::set_essential_values(std::size_t cell_number,
const std::vector<std::string> &order,
std::vector<double> &values) {
auto dump_it = values.begin();
// Solutions
if (this->Get_solution(cell_number) != NULL) {
this->Get_solution(cell_number)->set_essential_values(dump_it, order);
}
// Exchange
if (this->Get_exchange(cell_number) != NULL) {
this->Get_exchange(cell_number)->set_essential_values(dump_it);
}
// Kinetics
if (this->Get_kinetic(cell_number) != NULL) {
this->Get_kinetic(cell_number)->set_essential_values(dump_it);
}
// PPassemblage
if (this->Get_equilibrium(cell_number) != NULL) {
this->Get_equilibrium(cell_number)->set_essential_values(dump_it);
}
// Surface
if (this->Get_surface(cell_number) != NULL) {
this->Get_surface(cell_number)->set_essential_values(dump_it);
}
}

110
poet/src/Init.cpp Normal file
View File

@ -0,0 +1,110 @@
#include <IPhreeqcPOET.hpp>
#include <algorithm>
#include <cstddef>
#include <tuple>
#include <vector>
static inline std::vector<std::string>
unionStringVectors(const std::vector<std::string> &a,
const std::vector<std::string> &b) {
std::vector<std::string> result;
std::set_union(a.begin(), a.end(), b.begin(), b.end(),
std::back_inserter(result));
return result;
}
static std::vector<double>
createConcVector(const std::vector<std::string> &conc_names,
const std::vector<double>::const_iterator &conc_it,
const std::vector<std::string> &new_names) {
std::vector<double> conc_vec;
for (const auto &i : new_names) {
auto it = std::find(conc_names.begin(), conc_names.end(), i);
auto conc_index = std::distance(conc_names.begin(), it);
if (it != conc_names.end()) {
conc_vec.push_back(*(conc_it + conc_index));
} else {
conc_vec.push_back(std::numeric_limits<double>::quiet_NaN());
}
}
return conc_vec;
}
void IPhreeqcPOET::parseInit() {
std::map<int, std::pair<essential_names, std::vector<double>>> init_values;
for (const auto &[id, val] : this->PhreeqcPtr->Get_Rxn_solution_map()) {
// A key less than zero indicates an internal solution
if (id < 0) {
continue;
}
const essential_names curr_conc_names = this->dump_essential_names(id);
auto pair = std::make_pair(
curr_conc_names, this->get_essential_values(id, curr_conc_names[0]));
raw_initials[id] = pair;
}
}
IPhreeqcPOET::essential_names IPhreeqcPOET::union_raws(const RawMap &raws) {
IPhreeqcPOET::essential_names names;
for (const auto &[sol_id, val] : raws) {
for (std::size_t i = 0; i < names.size(); i++) {
names[i] = unionStringVectors(names[i], val.first[i]);
}
for (auto &specie : names[POET_EXCH]) {
specie = this->subExchangeName(specie);
}
}
return names;
}
std::vector<std::vector<double>>
IPhreeqcPOET::conc_from_essentials(const RawMap &raws,
const essential_names &names) {
std::vector<std::vector<double>> values;
for (const auto &[key, val] : raws) {
std::vector<double> combined_conc_vec;
for (std::size_t i = 0, offset = 0; i < names.size(); i++) {
std::vector<double> union_vec =
createConcVector(val.first[i], val.second.begin() + offset, names[i]);
combined_conc_vec.insert(combined_conc_vec.end(), union_vec.begin(),
union_vec.end());
offset += val.first[i].size();
}
values.push_back(combined_conc_vec);
}
}
IPhreeqcPOET::PhreeqcMat IPhreeqcPOET::getPhreeqcMat() {
std::vector<std::vector<double>> values;
const IPhreeqcPOET::essential_names names =
this->union_raws(this->raw_initials);
const std::vector<std::vector<double>> conc_values =
this->conc_from_essentials(this->raw_initials, names);
std::vector<int> ids(this->raw_initials.size());
std::transform(this->raw_initials.begin(), this->raw_initials.end(),
ids.begin(), [](const auto &pair) { return pair.first; });
}