max/iphreeqc
1
0
Fork 0
mirror of https://git.gfz-potsdam.de/naaice/iphreeqc.git synced 2025-12-15 16:18:22 +01:00
Code Issues Packages Projects Releases Wiki Activity

Merge branch 'ml/make-H0-O0-optional' into 'poet'

Browse Source 
feat: add option to include H(0) and O(0) in PhreeqcMatrix constructor and related functions

See merge request naaice/iphreeqc!26
This commit is contained in:
Max Lübke 2024-12-18 17:10:49 +01:00
commit c6e14a0ce9
8 changed files with 56 additions and 45 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

10
poet/include/PhreeqcMatrix.hpp
View File

@ -44,10 +44,12 @@ public:
* Construct a new Phreeqc Matrix object by reading the database and input
* script already present as a string.
*
* @param database
* @param input_script
* @param database Phreeqc database as a string.
* @param input_script Phreeqc input script as a string.
* @param with_h0_o0 Whether to include H(0) and O(0) in the output or not.
*/
PhreeqcMatrix(const std::string &database, const std::string &input_script);
PhreeqcMatrix(const std::string &database, const std::string &input_script,
bool with_h0_o0 = false);
/**
* @brief Construct a new Phreeqc Matrix object
@ -299,4 +301,6 @@ private:
std::shared_ptr<PhreeqcKnobs> _m_knobs;
std::string _m_database;
bool _m_with_h0_o0;
};

4
poet/src/PhreeqcMatrix/Ctor.cpp
View File

@ -9,8 +9,8 @@
#include <string>
PhreeqcMatrix::PhreeqcMatrix(const std::string &database,
const std::string &input_script)
: _m_database(database) {
const std::string &input_script, bool with_h0_o0)
: _m_database(database), _m_with_h0_o0(with_h0_o0) {
this->_m_pqc = std::make_shared<IPhreeqc>();
this->_m_pqc->LoadDatabaseString(database.c_str());

6
poet/src/PhreeqcMatrix/Init.cpp
View File

@ -18,12 +18,14 @@
#include <utility>
#include <vector>
bool include_h0_o0 = false;
static std::vector<std::string> dump_solution_names(cxxSolution *solution,
Phreeqc *phreeqc) {
std::vector<std::string> placeholder;
return phreeqc->find_all_valence_states(
SolutionWrapper::names(solution, placeholder));
SolutionWrapper::names(solution, include_h0_o0, placeholder));
}
template <enum PhreeqcMatrix::PhreeqcComponent comp, class T>
@ -205,6 +207,8 @@ void PhreeqcMatrix::initialize() {
this->_m_pqc->RunString("RUN_CELLS\n-cells 1\nEND");
}
include_h0_o0 = this->_m_with_h0_o0;
std::vector<std::string> solutions = find_all_solutions(phreeqc);
for (auto &[id, solution] : phreeqc->Get_Rxn_solution_map()) {

17
poet/src/Wrapper/SolutionWrapper.cpp
View File

@ -16,8 +16,6 @@ void SolutionWrapper::get(std::span<LDBLE> &data) const {
data[0] = solution->Get_total_h();
data[1] = solution->Get_total_o();
data[2] = solution->Get_cb();
data[3] = solution->Get_total("H(0)");
data[4] = solution->Get_total("O(0)");
std::size_t i = NUM_ESSENTIALS;
for (const auto &tot_name : solution_order) {
@ -38,9 +36,6 @@ void SolutionWrapper::set(const std::span<LDBLE> &data) {
const double &total_o = data[1];
const double &cb = data[2];
new_totals["H(0)"] = data[3];
new_totals["O(0)"] = data[4];
for (const auto &tot_name : solution_order) {
const double value = data[i++];
@ -54,19 +49,25 @@ void SolutionWrapper::set(const std::span<LDBLE> &data) {
}
std::vector<std::string>
SolutionWrapper::names(cxxSolution *solution,
SolutionWrapper::names(cxxSolution *solution, bool include_h0_o0,
std::vector<std::string> &solution_order) {
std::vector<std::string> names;
names.insert(names.end(), ESSENTIALS.begin(), ESSENTIALS.end());
if (include_h0_o0) {
names.push_back("H(0)");
names.push_back("O(0)");
}
std::set<std::string> names_set;
for (const auto &name : solution->Get_totals()) {
names_set.insert(name.first);
}
for (const auto &to_erase : ESSENTIALS) {
// don't care if the element was not found
// remove H(0) and O(0) from the set as they are already in the vector (if)
for (const auto &to_erase : {"H(0)", "O(0)"}) {
names_set.erase(to_erase);
}

7
poet/src/Wrapper/SolutionWrapper.hpp
View File

@ -17,7 +17,8 @@ public:
void set(const std::span<LDBLE> &data);
static std::vector<std::string>
names(cxxSolution *solution, std::vector<std::string> &solution_order);
names(cxxSolution *solution, bool include_h0_o0,
std::vector<std::string> &solution_order);
std::vector<std::string> getEssentials() const;
@ -25,8 +26,8 @@ private:
cxxSolution *solution;
const std::vector<std::string> solution_order;
static constexpr std::array<const char *, 5> ESSENTIALS = {"H", "O", "Charge",
"H(0)", "O(0)"};
static constexpr std::array<const char *, 3> ESSENTIALS = {"H", "O",
"Charge"};
static constexpr std::size_t NUM_ESSENTIALS = ESSENTIALS.size();
};

43
poet/test/testInput.hpp.in
View File

@ -24,17 +24,14 @@ RUN_CELLS
END)";
const std::vector<std::string> expected_names = {
"ID", "H", "O", "Charge", "H(0)", "O(0)",
"C(-4)", "C(4)", "Ca", "Cl", "Mg", "Na",
"Calcite_eq", "Calcite_si", "Dolomite_eq", "Dolomite_si"};
"ID", "H", "O", "Charge", "C(-4)", "C(4)", "Ca",
"Cl", "Mg", "Na", "Calcite_eq", "Calcite_si", "Dolomite_eq", "Dolomite_si"};
const std::vector<double> expected_values = {1,
111.01243522078478,
55.506323405120348,
-4.7919630342304069e-13,
0,
2.0180693312370871e-15,
0,
3.500625006800175e-05,
0.12131244646561848,
0.49999844804496646,
@ -46,8 +43,7 @@ const std::vector<double> expected_values = {1,
0};
const std::vector<double> expected_errors = {
0, 1e-3, 1e-3, 1e-15, 1e-15, 1e-15, 1e-5, 1e-5,
1e-5, 1e-5, 1e-5, 1e-5, 1e-5, 0, 1e-5, 0};
0, 1e-3, 1e-3, 1e-15, 1e-5, 1e-5, 1e-5, 1e-5, 1e-5, 1e-5, 1e-5, 0, 1e-5, 0};
const std::string phreeqc_database = R"database(@POET_PHREEQCDAT_DB@)database";
} // namespace base_test
@ -56,19 +52,27 @@ namespace barite_test {
const std::string script = R"barite(@POET_BARITE_PQI@)barite";
const std::string database = R"barite(@POET_BARITE_DB@)barite";
const std::vector<std::string> expected_names = {
"ID", "H", "O", "Charge", "H(0)",
"O(0)", "Ba", "Cl", "S(-2)", "S(6)",
"Sr", "Barite", "Barite_p1", "Celestite", "Celestite_p1",
"Celestite_eq", "Celestite_si"};
const std::vector<std::string> expected_names = {"ID",
"H",
"O",
"Charge",
"Ba",
"Cl",
"S(-2)",
"S(6)",
"Sr",
"Barite",
"Barite_p1",
"Celestite",
"Celestite_p1",
"Celestite_eq",
"Celestite_si"};
const std::vector<double> expected_values_line_one = {
1,
111.01243359383071,
55.508698688362124,
-1.2153078399577636e-09,
0,
8.6371063688066983e-15,
1.0000001848805677e-12,
1.0000000116187218e-12,
0,
@ -86,8 +90,6 @@ const std::vector<double> expected_errors = {
1e-3,
1e-3,
1e-15,
1e-15,
1e-15,
1e-5,
1e-5,
1e-5,
@ -101,13 +103,12 @@ const std::vector<double> expected_errors = {
0};
const std::vector<std::string> expected_names_erased = {
"ID", "H", "O", "Charge", "H(0)",
"O(0)", "Ba", "Cl", "S(-2)", "S(6)",
"Sr", "Barite", "Barite_p1", "Celestite", "Celestite_p1"};
"ID", "H", "O", "Charge", "Ba", "Cl", "S(-2)",
"S(6)", "Sr", "Barite", "Barite_p1", "Celestite", "Celestite_p1"};
const std::vector<std::string> expected_names_subset = {
"ID", "H", "O", "Charge", "H(0)", "O(0)", "Ba",
"Cl", "S(-2)", "S(6)", "Sr", "Celestite_eq", "Celestite_si"};
"ID", "H", "O", "Charge", "Ba", "Cl",
"S(-2)", "S(6)", "Sr", "Celestite_eq", "Celestite_si"};
} // namespace barite_test
namespace test_engine {

10
poet/test/testPhreeqcMatrix.cpp
View File

@ -23,8 +23,8 @@ POET_TEST(PhreeqcMatrixOneSolution) {
EXPECT_FALSE(pqc_mat.checkIfExists(2));
PhreeqcMatrix::STLExport exported_init = pqc_mat.get();
// ID + H,O,Charge,H(0),O(0) + 6 Solutions + 4 Equil incl. params
EXPECT_EQ(exported_init.names.size(), 16);
// ID + H,O,Charge + 6 Solutions + 4 Equil incl. params
EXPECT_EQ(exported_init.names.size(), 14);
EXPECT_EQ(exported_init.names, base_test::expected_names);
for (std::size_t i = 0; i < exported_init.values.size(); ++i) {
@ -73,7 +73,7 @@ POET_TEST(PhreeqcMatrixMultiSolution) {
EXPECT_EQ(exported.names, barite_test::expected_names);
for (std::size_t i = 0; i < exported.names.size(); i++) {
if (i > 10 && i < 15) {
if (i > 8 && i < 13) {
EXPECT_TRUE(std::isnan(exported.values[i]));
continue;
}
@ -124,7 +124,7 @@ POET_TEST(PhreeqcMatrixCtor) {
EXPECT_EQ(exported.names, barite_test::expected_names);
for (std::size_t i = 0; i < exported.names.size(); i++) {
if (i > 10 && i < 15) {
if (i > 8 && i < 13) {
EXPECT_TRUE(std::isnan(exported.values[i]));
continue;
}
@ -149,7 +149,7 @@ POET_TEST(PhreeqcMatrixOperator) {
EXPECT_EQ(exported.names, barite_test::expected_names);
for (std::size_t i = 0; i < exported.names.size(); i++) {
if (i > 10 && i < 15) {
if (i > 8 && i < 13) {
EXPECT_TRUE(std::isnan(exported.values[i]));
continue;
}

4
poet/test/testPhreeqcRunner.cpp
View File

@ -48,10 +48,10 @@ POET_TEST(PhreeqcRunnerSimulation) {
const bool is_first_half = cell_index < half_cells;
if (is_first_half) {
EXPECT_EQ(simulationInOut[cell_index][0], expected_value_first_half);
EXPECT_TRUE(std::isnan(simulationInOut[cell_index][11]));
EXPECT_TRUE(std::isnan(simulationInOut[cell_index][9]));
} else {
EXPECT_EQ(simulationInOut[cell_index][0], expected_value_second_half);
EXPECT_FALSE(std::isnan(simulationInOut[cell_index][11]));
EXPECT_FALSE(std::isnan(simulationInOut[cell_index][9]));
}
EXPECT_NEAR(simulationInOut[cell_index][1], 111, 1);