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

feat: add H(0) and O(0) to solution list

Browse Source 
refactor: Remove knowledge on "fixed" solution entries from non-wrapper functions
This commit is contained in:
Max Lübke 2024-10-23 10:01:49 +02:00
parent c046e907b4
commit 26f6e11177
9 changed files with 88 additions and 74 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

9
poet/include/PhreeqcMatrix.hpp
View File

@ -166,19 +166,18 @@ public:
};
/**
* @brief Get all found solution names
* @brief Get all found solution names.
*
* @return std::vector<std::string> Vector containing all solution names.
*/
std::vector<std::string>
getSolutionNames(bool include_h_o_charge = false) const;
std::vector<std::string> getSolutionNames() const;
/**
* @brief Get solution total names of all found solutions (excluding H, O,
* Charge)
* Charge, H(0), O(0))
*
* @return std::vector<std::string> Names of all found solutions (excluding H,
* O, Charge)
* O, Charge, H(0), O(0))
*/
std::vector<std::string> getSolutionPrimaries() const;

15
poet/src/PhreeqcMatrix/Access.cpp
View File

@ -157,25 +157,18 @@ PhreeqcMatrix::STLExport PhreeqcMatrix::get(VectorExportType type,
return result;
}
std::vector<std::string>
PhreeqcMatrix::getSolutionNames(bool include_h_o_charge) const {
std::vector<std::string> PhreeqcMatrix::getSolutionNames() const {
std::vector<std::string> names;
if (include_h_o_charge) {
names.push_back("H");
names.push_back("O");
names.push_back("Charge");
}
const auto &first_element = _m_map.begin()->second;
for (std::size_t i = 3; i < _m_map.begin()->second.size(); i++) {
for (const auto &element : _m_map.begin()->second) {
// assuming the element vector always starts with the solution components
if (first_element[i].type != PhreeqcComponent::SOLUTION) {
if (element.type != PhreeqcComponent::SOLUTION) {
break;
}
names.push_back(first_element[i].name);
names.push_back(element.name);
}
return names;

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

@ -17,14 +17,12 @@
#include <utility>
#include <vector>
constexpr std::size_t SKIP_H_O_CB = 3;
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), SKIP_H_O_CB);
SolutionWrapper::names(solution, placeholder));
}
template <enum PhreeqcMatrix::PhreeqcComponent comp, class T>
@ -122,8 +120,7 @@ create_vector_from_phreeqc(Phreeqc *phreeqc, int id,
// Solution
SolutionWrapper sol_wrapper(
Utilities::Rxn_find(phreeqc->Get_Rxn_solution_map(), id),
{solution_names.begin() + SKIP_H_O_CB, solution_names.end()});
Utilities::Rxn_find(phreeqc->Get_Rxn_solution_map(), id), solution_names);
base_add_to_element_vector<PhreeqcMatrix::PhreeqcComponent::SOLUTION>(
sol_wrapper, solution_names, elements);

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

@ -1,11 +1,13 @@
#include "SolutionWrapper.hpp"
#include "NameDouble.h"
#include <set>
#include <vector>
SolutionWrapper::SolutionWrapper(
cxxSolution *soln, const std::vector<std::string> &solution_order_)
: solution(soln), solution_order(solution_order_) {
this->num_elements = solution_order.size() + NUM_ESSENTIALS;
cxxSolution *soln, const std::vector<std::string> &_solution_order)
: solution(soln), solution_order(_solution_order.begin() + NUM_ESSENTIALS,
_solution_order.end()) {
this->num_elements = _solution_order.size();
auto &totals = solution->Get_totals();
}
@ -14,6 +16,8 @@ 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) {
@ -29,6 +33,10 @@ void SolutionWrapper::get(std::span<LDBLE> &data) const {
void SolutionWrapper::set(const std::span<LDBLE> &data) {
std::size_t i = NUM_ESSENTIALS;
cxxNameDouble new_totals;
new_totals["H(0)"] = data[3];
new_totals["O(0)"] = data[4];
for (const auto &tot_name : solution_order) {
const double value = data[i++];
@ -46,18 +54,18 @@ SolutionWrapper::names(cxxSolution *solution,
std::vector<std::string> &solution_order) {
std::vector<std::string> names;
names.push_back("H");
names.push_back("O");
names.push_back("Charge");
names.insert(names.end(), ESSENTIALS.begin(), ESSENTIALS.end());
for (const auto &[tot_name, value] : solution->Get_totals()) {
if (tot_name == "H(0)" || tot_name == "O(0)") {
continue;
}
solution_order.push_back(tot_name);
std::set<std::string> names_set;
for (const auto &name : solution->Get_totals()) {
names_set.insert(name.first);
}
names.insert(names.end(), solution_order.begin(), solution_order.end());
for (const auto &to_erase : ESSENTIALS) {
// don't care if the element was not found
names_set.erase(to_erase);
}
names.insert(names.end(), names_set.begin(), names_set.end());
return names;
}

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

@ -3,7 +3,9 @@
#include "NameDouble.h"
#include "Solution.h"
#include "WrapperBase.hpp"
#include <array>
#include <cstddef>
#include <string>
#include <vector>
class SolutionWrapper : public WrapperBase {
@ -18,9 +20,14 @@ public:
static std::vector<std::string>
names(cxxSolution *solution, std::vector<std::string> &solution_order);
std::vector<std::string> getEssentials() const;
private:
cxxSolution *solution;
const std::vector<std::string> solution_order;
static constexpr std::size_t NUM_ESSENTIALS = 3;
static constexpr std::array<std::string, 5> ESSENTIALS = {"H", "O", "Charge",
"H(0)", "O(0)"};
static constexpr std::size_t NUM_ESSENTIALS = ESSENTIALS.size();
};

10
poet/test/testPhreeqcMatrix.cpp
View File

@ -20,8 +20,8 @@ POET_TEST(PhreeqcMatrixOneSolution) {
EXPECT_EQ(ids[0], 1);
PhreeqcMatrix::STLExport exported_init = pqc_mat.get();
// ID + H,O,Charge + 4 Solutions + 4 Equil incl. params
EXPECT_EQ(exported_init.names.size(), 12);
// ID + H,O,Charge,H(0),O(0) + 4 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) {
@ -58,7 +58,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 > 8 && i < 13) {
if (i > 10 && i < 15) {
EXPECT_TRUE(std::isnan(exported.values[i]));
continue;
}
@ -109,7 +109,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 > 8 && i < 13) {
if (i > 10 && i < 15) {
EXPECT_TRUE(std::isnan(exported.values[i]));
continue;
}
@ -134,7 +134,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 > 8 && i < 13) {
if (i > 10 && i < 15) {
EXPECT_TRUE(std::isnan(exported.values[i]));
continue;
}

54
poet/test/testPhreeqcMatrix.hpp.in
View File

@ -20,23 +20,26 @@ Dolomite 0.0 0
END)";
const std::vector<std::string> expected_names = {
"ID", "H", "O", "Charge", "Ca", "Cl",
"Mg", "Na", "Calcite_eq", "Calcite_si", "Dolomite_eq", "Dolomite_si"};
"ID", "H", "O", "Charge", "H(0)", "O(0)", "Ca",
"Cl", "Mg", "Na", "Calcite_eq", "Calcite_si", "Dolomite_eq", "Dolomite_si"};
const std::vector<double> expected_values = {1,
111.01243521533338,
111.01243522078481,
55.506218386370165,
-4.7941959748097226e-13,
0.1,
0.5,
0.1,
0.1,
-4.7940411972733642e-13,
1.2051199818688651e-25,
0,
0.10000000000000005,
0.49999844804496646,
0.10000000000000002,
0.10000000000000001,
1,
0,
0,
0};
const std::vector<double> expected_errors = {
0, 1e-3, 1e-3, 1e-15, 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
@ -45,27 +48,19 @@ 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",
"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", "H(0)",
"O(0)", "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,
@ -88,6 +83,8 @@ const std::vector<double> expected_errors = {
1e-5,
1e-5,
1e-5,
1e-5,
1e-5,
std::numeric_limits<double>::quiet_NaN(),
std::numeric_limits<double>::quiet_NaN(),
std::numeric_limits<double>::quiet_NaN(),
@ -96,10 +93,11 @@ const std::vector<double> expected_errors = {
0};
const std::vector<std::string> expected_names_erased = {
"ID", "H", "O", "Charge", "Ba", "Cl", "S(-2)",
"S(6)", "Sr", "Barite", "Barite_p1", "Celestite", "Celestite_p1"};
"ID", "H", "O", "Charge", "H(0)",
"O(0)", "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", "Ba", "Cl",
"S(-2)", "S(6)", "Sr", "Celestite_eq", "Celestite_si"};
"ID", "H", "O", "Charge", "H(0)", "O(0)", "Ba",
"Cl", "S(-2)", "S(6)", "Sr", "Celestite_eq", "Celestite_si"};
} // namespace barite_test

21
src/phreeqcpp/GFZ.cpp
View File

@ -1,14 +1,27 @@
#include "Phreeqc.h"
#include <set>
const std::set<std::string> to_ignore = {"H", "O", "Charge", "H(0)", "O(0)"};
std::vector<std::string> Phreeqc::find_all_valence_states(
const std::vector<std::string> &&solution_names, const std::size_t offset) {
std::vector<std::string> solution_with_valences(
solution_names.begin(), solution_names.begin() + offset);
const std::vector<std::string> &&solution_names) {
std::vector<std::string> solution_with_valences;
solution_with_valences.reserve(solution_names.size());
// to avoid duplicates store already evaluated master species
std::set<std::string> master_species_found;
for (std::size_t i = offset; i < solution_names.size(); i++) {
auto is_ignored = [](const std::string &in) {
return (to_ignore.find(in) != to_ignore.end());
};
for (std::size_t i = 0; i < solution_names.size(); i++) {
if (is_ignored(solution_names[i])) {
solution_with_valences.push_back(solution_names[i]);
continue;
}
const auto master_primary =
master_bsearch_primary(solution_names[i].c_str());

3
src/phreeqcpp/Phreeqc.h
View File

@ -1842,8 +1842,7 @@ protected:
public:
std::vector<std::string>
find_all_valence_states(const std::vector<std::string> &&solution_names,
const std::size_t offset);
find_all_valence_states(const std::vector<std::string> &&solution_names);
static const class const_iso iso_defaults[];
static const int count_iso_defaults;