MDL/golem: added bunch of getMatrix* getters and poet/test/testGetters.cpp

poet/CMakeLists.txt

@@ -61,5 +61,8 @@ if (BUILD_TESTING AND STANDALONE_BUILD)
     gtest_discover_tests(poet_test)
 endif()
 
-add_executable(golemrunner test/testGolemRunner.cpp)
+add_executable(testGolemRunner test/testGolemRunner.cpp)
-target_link_libraries(golemrunner IPhreeqcPOET)
+target_link_libraries(testGolemRunner IPhreeqcPOET)
+
+add_executable(testGetters     test/testGetters.cpp)
+target_link_libraries(testGetters IPhreeqcPOET)

poet/include/PhreeqcMatrix.hpp

@@ -313,6 +313,48 @@ public:
    */
   bool withRedox() const { return _m_with_redox; }
 
+  // MDL
+  /**
+   * @brief Returns all column names of the Matrix pertaining to KINETICS
+   *
+   * This function returns a string vector.
+   *
+   * @return std::vector<std::string> Whole vector of names. Empty if no KINETICS
+   * is defined
+   */
+   std::vector<std::string> getMatrixKinetics() const;
+   
+  /**
+   * @brief Returns all column names of the Matrix pertaining to EQUILIBRIUM
+   *
+   * This function returns a string vector.
+   *
+   * @return std::vector<std::string> Whole vector of names. Empty if no EQUILIBRIUM
+   * is defined
+   */
+   std::vector<std::string> getMatrixEquilibrium() const;
+
+      
+   /*
+   
+   @brief Returns all column names of the Matrix pertaining to
+   quantities that must be transported
+   
+   @return std::vector<std::string> vector of names
+
+   */
+   std::vector<std::string> getMatrixMustTransport() const;
+
+   /*
+
+   @brief Returns all column names of the Matrix pertaining to
+   quantities that must NOT be transported but have to be included in
+   the output
+   
+   @return std::vector<std::string> vector of names
+   */
+   std::vector<std::string> getMatrixOutOnly() const;
+
 private:
   std::map<int, std::vector<element>> _m_map;
   std::map<int, std::vector<base_names>> _m_internal_names;
@@ -330,4 +372,4 @@ private:
 
   bool _m_with_h0_o0;
   bool _m_with_redox;
-};
+};

poet/src/PhreeqcMatrix/Access.cpp

@@ -240,3 +240,99 @@ double PhreeqcMatrix::operator()(int cell_id, const std::string &name) const {
 
   return it->value;
 }
+
+
+// MDL
+std::vector<std::string> PhreeqcMatrix::getMatrixKinetics() const {
+    std::vector<std::string> names;
+    
+    auto n = this->getIds().size();
+    for (auto i = 0; i<n; ++i) {
+	auto pqc_kinnames = this->getKineticsNames(i);
+	for (auto nam : pqc_kinnames ) {
+	    for (auto mat_name : this->get().names){
+		if (mat_name.starts_with(nam)) {
+		    names.push_back(mat_name);
+		}
+	    }
+	}
+    }
+    std::sort(names.begin(), names.end());
+    std::vector<std::string>::iterator it;
+    it = std::unique(names.begin(), names.end());
+    names.resize(std::distance(names.begin(),it) );
+    return names;
+}
+
+
+std::vector<std::string> PhreeqcMatrix::getMatrixEquilibrium() const {
+
+    std::vector<std::string> names;
+    std::vector<std::string> mat_names = this->get().names;
+    auto n = this->getIds().size();
+    for (auto i = 0; i<n; ++i) {
+	auto pqc_eqnames = this->getEquilibriumNames(i);
+	for (auto nam : pqc_eqnames ) {
+	    for (auto mat_name : mat_names){
+		if (mat_name.starts_with(nam)) {
+		    if (std::find(names.begin(), names.end(), mat_name) == names.end()) {
+			names.push_back(mat_name);
+		    }
+		}
+	    }
+	}
+    }
+    // std::sort(names.begin(), names.end());
+    // std::vector<std::string>::iterator it;
+    // it = std::unique(names.begin(), names.end());
+    // names.resize(std::distance(names.begin(),it) );
+    return names;
+}
+
+std::vector<std::string> PhreeqcMatrix::getMatrixMustTransport() const {
+    std::vector<std::string> names;
+    
+    const std::vector<std::string> to_remove = {
+	"tc", "patm", "SolVol", "pH", "pe"
+    };
+
+    // sols contains all solutes; we must remove { tc, patm, SolVol, pH, pe }
+    auto sols = this->getSolutionNames();
+    for (auto name : sols) {
+	if (std::find(to_remove.begin(), to_remove.end(), name) == to_remove.end()) {
+	    names.push_back(name);
+	}
+    }
+    
+    return names;
+}
+
+std::vector<std::string> PhreeqcMatrix::getMatrixOutOnly() const {
+    // MDL we must append here selected_output / user_punch
+    std::vector<std::string> defaultnames = {
+	"tc", "patm", "SolVol", "pH", "pe"
+    };
+    std::vector<std::string> ret;
+    for (auto nm : defaultnames) {
+	ret.push_back(nm);
+    }
+    return ret;
+}
+
+// std::vector<std::string> PhreeqcMatrix::getTransported() const {
+//   std::vector<std::string> names;
+
+//   const auto &first_element = _m_map.begin()->second;
+
+//   for (const auto &element : _m_map.begin()->second) {
+//     // assuming the element vector always starts with the solution components
+//     if (element.type != PhreeqcComponent::SOLUTION) {
+//       break;
+//     }
+
+//     names.push_back(element.name);
+//   }
+
+//   return names;
+// }
+

poet/test/testGetters.cpp

@@ -0,0 +1,103 @@
+//  Time-stamp: "Last modified 2025-07-28 19:54:57 delucia"
+#include <iostream>
+#include <iomanip>
+#include <linux/limits.h>
+#include <sstream>
+#include <fstream>
+#include <algorithm>
+#include <iterator>
+#include <vector>
+#include <string>
+#include <memory>
+#include <cmath>
+
+#include "PhreeqcMatrix.hpp"
+#include "PhreeqcEngine.hpp"
+#include "PhreeqcRunner.hpp"
+
+
+std::string readFile(const std::string &path) {
+  std::string string_rpath(PATH_MAX, '\0');
+
+  if (realpath(path.c_str(), string_rpath.data()) == nullptr) {
+    throw std::runtime_error(":: Failed to resolve the realpath to file " + path);
+  }
+
+  std::ifstream file(string_rpath);
+
+  if (!file.is_open()) {
+    throw std::runtime_error(":: Failed to open file: " + path);
+  }
+
+  std::stringstream buffer;
+  buffer << file.rdbuf();
+
+  return buffer.str();
+}
+
+// pretty print a vector, standard implementation from stackoverflow 
+template<typename T>
+std::ostream & operator<<(std::ostream & os, std::vector<T> vec)
+{
+    os << "{ ";
+    std::copy(vec.begin(), vec.end(), std::ostream_iterator<T>(os, ", "));
+    os << " }";
+    return os;
+}
+
+
+int main(int argc, char *argv[]) {
+    
+    if (argc < 3) {
+	std::cout << "::" << argv[0] << ": two args needed, script and database\n";
+	return 1;
+    }
+
+
+    ////// INITIALISATION
+    // read Script and Database and put it in a std::string
+    auto script = readFile(argv[1]);
+    auto db     = readFile(argv[2]);
+
+    // Create the matrix directly from database and init script
+    std::cout << ":: Creating a PhreeqcMatrix with valence states and H(0)/O(0) \n";
+    PhreeqcMatrix pqc_mat1(db, script, true, true);
+    
+    // How many different SOLUTIONS ("CELLS") are defined in the script?
+    const auto ids = pqc_mat1.getIds();
+    
+    int n = ids.size();
+    
+    std::cout << ":: Found " << n << " distinct PHREEQC problems \n";
+    std::cout << ids << "\n";
+
+    std::cout << ":: getSolutionsNames(): the common solutes across all problems: \n";
+    const auto solutes1  = pqc_mat1.getSolutionNames();
+    std::cout << solutes1 << "\n";
+
+    // auto expmat = pqc_mat1.get();
+    auto allvars = pqc_mat1.get().names;
+    std::cout << ":: pqc_mat1.get().names: \n";
+    std::cout << allvars << "\n";
+
+    auto MatNamesKin = pqc_mat1.getMatrixKinetics(); 
+    std::cout << ":: pqc_mat1.getMatrixKinetics(): \n";
+    std::cout << MatNamesKin << "\n";
+    
+    auto MatNamesEqui = pqc_mat1.getMatrixEquilibrium(); 
+    std::cout << ":: pqc_mat1.getMatrixEquilibrium(): \n";
+    std::cout << MatNamesEqui << "\n";
+
+    auto transported = pqc_mat1.getMatrixMustTransport(); 
+    std::cout << ":: pqc_mat1.getMatrixMustTransport(): \n";
+    std::cout << transported << "\n";
+
+    auto outonly = pqc_mat1.getMatrixOutOnly(); 
+    std::cout << ":: pqc_mat1.getMatrixOutOnly(): \n";
+    std::cout << outonly << "\n";
+    return 0;
+}
+
+// Oneliner for rz-vm278 relative to iphreeqc/poet/test!!
+
+// g++ testGolemRunner.cpp -o testG -Wall -I../../poet/include -I../../src -I../../src/phreeqcpp -I../../src/phreeqcpp/common -I../../src/phreeqcpp/PhreeqcKeywords -lIPhreeqc -lIPhreeqcPOET -L../../bbuild/ -L../../bbuild/poet