diff --git a/poet/CMakeLists.txt b/poet/CMakeLists.txt index fa2c901b..c67413bf 100644 --- a/poet/CMakeLists.txt +++ b/poet/CMakeLists.txt @@ -1,5 +1,6 @@ set(POET_SOURCE_FILES - src/Engine.cpp + src/Engine.cpp + src/Runner.cpp #Wrappers src/Wrapper/EquilibriumWrapper.cpp src/Wrapper/EquilibriumCompWrapper.cpp @@ -33,6 +34,7 @@ if (BUILD_TESTING AND STANDALONE_BUILD) set(POET_TEST_SOURCE_FILES test/testPhreeqcEngine.cpp test/testPhreeqcMatrix.cpp + test/testPhreeqcRunner.cpp test/utils.cpp ) diff --git a/poet/include/PhreeqcMatrix.hpp b/poet/include/PhreeqcMatrix.hpp index 3162fc85..0af656be 100644 --- a/poet/include/PhreeqcMatrix.hpp +++ b/poet/include/PhreeqcMatrix.hpp @@ -57,7 +57,7 @@ public: * object! * @param other */ - PhreeqcMatrix(const PhreeqcMatrix &other); + PhreeqcMatrix(const PhreeqcMatrix &other) = default; /** * @brief Construct a new Phreeqc Matrix object @@ -66,7 +66,7 @@ public: * object! * @param other */ - PhreeqcMatrix(PhreeqcMatrix &&other); + PhreeqcMatrix(PhreeqcMatrix &&other) = default; /** * @brief Assignment operator @@ -75,7 +75,7 @@ public: * @param other * @return PhreeqcMatrix& */ - PhreeqcMatrix &operator=(const PhreeqcMatrix &other); + PhreeqcMatrix &operator=(const PhreeqcMatrix &other) = default; /** * @brief Assignment operator @@ -84,7 +84,7 @@ public: * @param other * @return PhreeqcMatrix& */ - PhreeqcMatrix &operator=(PhreeqcMatrix &&other); + PhreeqcMatrix &operator=(PhreeqcMatrix &&other) = default; /** * @brief Access the value of a given cell by name. diff --git a/poet/include/PhreeqcRunner.hpp b/poet/include/PhreeqcRunner.hpp new file mode 100644 index 00000000..5c5230a3 --- /dev/null +++ b/poet/include/PhreeqcRunner.hpp @@ -0,0 +1,66 @@ +#pragma once + +#include "PhreeqcEngine.hpp" +#include "PhreeqcMatrix.hpp" +#include +#include +#include +#include + +/** + * @class PhreeqcRunner + * @brief Manages the execution of Phreeqc simulations. + * + * The PhreeqcRunner class is responsible for managing and running Phreeqc + * simulations. It initializes a set of PhreeqcEngine instances based on the + * provided PhreeqcMatrix and provides functionality to run simulations and + * retrieve the number of engines. + * + * @note Copy and move operations are deleted to prevent unintended behavior. + */ +class PhreeqcRunner { +public: + PhreeqcRunner(const PhreeqcRunner &) = delete; + PhreeqcRunner(PhreeqcRunner &&) = delete; + PhreeqcRunner &operator=(const PhreeqcRunner &) = delete; + PhreeqcRunner &operator=(PhreeqcRunner &&) = delete; + + /** + * @brief Constructs a PhreeqcRunner object with the given PhreeqcMatrix. + * + * Constructs a PhreeqcRunner object with the given PhreeqcMatrix. For each + * cell found in the PhreeqcMatrix, a PhreeqcEngine is created and stored in + * an internal map. + * + * @param matrix A reference to a PhreeqcMatrix object used to initialize the + * PhreeqcRunner. + */ + PhreeqcRunner(const PhreeqcMatrix &matrix); + ~PhreeqcRunner() = default; + + /** + * @brief Runs the simulation with the given input and output data for a + * specified time step. + * + * @param simulationInOut A reference to a 2D vector containing the input data + * for the simulation. The vector will be modified to contain the output data + * after the simulation. + * @param time_step The time step for the simulation. + */ + void run(std::vector> &simulationInOut, + const double time_step); + + /** + * @brief Returns the number of engines currently stored. + * + * This function provides the count of engines that are currently + * stored in the _engineStorage container. + * + * @return std::size_t The number of engines. + */ + std::size_t numEngines() const { return _engineStorage.size(); } + +private: + std::unordered_map> _engineStorage; + std::vector _buffer; +}; \ No newline at end of file diff --git a/poet/src/Engine.cpp b/poet/src/Engine.cpp index d36c46d0..22df7a0a 100644 --- a/poet/src/Engine.cpp +++ b/poet/src/Engine.cpp @@ -110,8 +110,8 @@ void PhreeqcEngine::runCell(std::vector &cell_values, throw std::invalid_argument("Time step must be positive"); } - // skip ID - std::span cell_data{cell_values.begin() + 1, cell_values.end()}; + // ID is already skipped by PhreeqcRunner, so no need to start ahead + std::span cell_data{cell_values.begin(), cell_values.end()}; this->impl->set_essential_values(cell_data); this->impl->run(time_step); diff --git a/poet/src/PhreeqcMatrix/Access.cpp b/poet/src/PhreeqcMatrix/Access.cpp index fb770712..59b6fc6c 100644 --- a/poet/src/PhreeqcMatrix/Access.cpp +++ b/poet/src/PhreeqcMatrix/Access.cpp @@ -7,9 +7,10 @@ #include PhreeqcMatrix PhreeqcMatrix::subset(const std::vector &indices) const { - PhreeqcMatrix result; + PhreeqcMatrix result(*this); - result._m_pqc = _m_pqc; + result._m_map.clear(); + result._m_internal_names.clear(); for (const auto &index : indices) { result._m_map[index] = _m_map.at(index); diff --git a/poet/src/PhreeqcMatrix/Ctor.cpp b/poet/src/PhreeqcMatrix/Ctor.cpp index b31b6632..efdec903 100644 --- a/poet/src/PhreeqcMatrix/Ctor.cpp +++ b/poet/src/PhreeqcMatrix/Ctor.cpp @@ -17,32 +17,32 @@ PhreeqcMatrix::PhreeqcMatrix(const std::string &database, this->initialize(); } -PhreeqcMatrix::PhreeqcMatrix(const PhreeqcMatrix &other) - : _m_map(other._m_map), _m_internal_names(other._m_internal_names), - _m_surface_primaries(other._m_surface_primaries), _m_pqc(other._m_pqc), - _m_database(other._m_database) {} +// PhreeqcMatrix::PhreeqcMatrix(const PhreeqcMatrix &other) +// : _m_map(other._m_map), _m_internal_names(other._m_internal_names), +// _m_surface_primaries(other._m_surface_primaries), _m_pqc(other._m_pqc), +// _m_database(other._m_database) {} -PhreeqcMatrix::PhreeqcMatrix(PhreeqcMatrix &&other) - : _m_map(other._m_map), _m_internal_names(other._m_internal_names), - _m_surface_primaries(other._m_surface_primaries), _m_pqc(other._m_pqc), - _m_database(other._m_database) {} +// PhreeqcMatrix::PhreeqcMatrix(PhreeqcMatrix &&other) +// : _m_map(other._m_map), _m_internal_names(other._m_internal_names), +// _m_surface_primaries(other._m_surface_primaries), _m_pqc(other._m_pqc), +// _m_database(other._m_database) {} -PhreeqcMatrix &PhreeqcMatrix::operator=(const PhreeqcMatrix &other) { - _m_map = other._m_map; - _m_internal_names = other._m_internal_names; - _m_surface_primaries = other._m_surface_primaries; - _m_pqc = other._m_pqc; - _m_database = other._m_database; +// PhreeqcMatrix &PhreeqcMatrix::operator=(const PhreeqcMatrix &other) { +// _m_map = other._m_map; +// _m_internal_names = other._m_internal_names; +// _m_surface_primaries = other._m_surface_primaries; +// _m_pqc = other._m_pqc; +// _m_database = other._m_database; - return *this; -} +// return *this; +// } -PhreeqcMatrix &PhreeqcMatrix::operator=(PhreeqcMatrix &&other) { - _m_map = other._m_map; - _m_internal_names = other._m_internal_names; - _m_surface_primaries = other._m_surface_primaries; - _m_pqc = other._m_pqc; - _m_database = other._m_database; +// PhreeqcMatrix &PhreeqcMatrix::operator=(PhreeqcMatrix &&other) { +// _m_map = other._m_map; +// _m_internal_names = other._m_internal_names; +// _m_surface_primaries = other._m_surface_primaries; +// _m_pqc = other._m_pqc; +// _m_database = other._m_database; - return *this; -} \ No newline at end of file +// return *this; +// } \ No newline at end of file diff --git a/poet/src/Runner.cpp b/poet/src/Runner.cpp new file mode 100644 index 00000000..44913b20 --- /dev/null +++ b/poet/src/Runner.cpp @@ -0,0 +1,57 @@ +#include "PhreeqcEngine.hpp" +#include "PhreeqcMatrix.hpp" +#include "PhreeqcRunner.hpp" +#include +#include +#include +#include + +PhreeqcRunner::PhreeqcRunner(const PhreeqcMatrix &matrix) { + // first make sure to have enough space in our buffer + this->_buffer.reserve(matrix.get().names.size()); + + // Create a PhreeqcEngine for each id + for (const auto &id : matrix.getIds()) { + this->_engineStorage[id] = std::make_unique(matrix, id); + } +} + +static void copy_to_buffer(std::vector &buffer, + const std::vector ¤t_inout) { + for (std::size_t j = 1; j < current_inout.size(); j++) { + if (std::isnan(current_inout[j])) { + continue; + } + buffer.push_back(current_inout[j]); + } +} + +static void copy_from_buffer(const std::vector &buffer, + std::vector ¤t_inout) { + std::size_t buffer_index = 0; + for (std::size_t j = 1; j < current_inout.size(); j++) { + if (std::isnan(current_inout[j])) { + continue; + } + current_inout[j] = buffer[buffer_index++]; + } +} + +void PhreeqcRunner::run(std::vector> &simulationInOut, + const double time_step) { + for (std::size_t i = 0; i < simulationInOut.size(); i++) { + const auto ¤t_inout = simulationInOut[i]; + const auto pqc_id = static_cast(current_inout[0]); + + this->_buffer.clear(); + + // Copy the input to the buffer while ignoring the first element and NaNs + copy_to_buffer(this->_buffer, current_inout); + + this->_engineStorage.at(pqc_id)->runCell(this->_buffer, time_step); + + // Copy the buffer back to the output while ignoring the first element and + // NaNs + copy_from_buffer(this->_buffer, simulationInOut[i]); + } +} \ No newline at end of file diff --git a/poet/test/testInput.hpp.in b/poet/test/testInput.hpp.in index 6e011b8b..664c1263 100644 --- a/poet/test/testInput.hpp.in +++ b/poet/test/testInput.hpp.in @@ -4,6 +4,8 @@ #include #include +#include + #define POET_TEST(name) TEST(TestPOET, name) namespace base_test { diff --git a/poet/test/testPhreeqcRunner.cpp b/poet/test/testPhreeqcRunner.cpp new file mode 100644 index 00000000..41f0d1d5 --- /dev/null +++ b/poet/test/testPhreeqcRunner.cpp @@ -0,0 +1,75 @@ +#include "PhreeqcRunner.hpp" +#include "utils.hpp" + +#include +#include +#include +#include +#include + +const std::string test_database = readFile(barite_test::database); +const std::string test_script = readFile(barite_test::script); + +constexpr std::size_t num_cells = 10; + +POET_TEST(PhreeqcRunnerConstructor) { + PhreeqcMatrix pqc_mat(test_database, test_script); + EXPECT_NO_THROW(PhreeqcRunner tmp(pqc_mat)); +} + +POET_TEST(PhreeqcRunnerSimulation) { + PhreeqcMatrix pqc_mat(test_database, test_script); + const auto subsetted_pqc_mat = pqc_mat.subset({2, 3}); + PhreeqcRunner runner(subsetted_pqc_mat); + + const auto stl_mat = subsetted_pqc_mat.get(); + const auto matrix_values = stl_mat.values; + const auto num_columns = stl_mat.names.size(); + + std::vector> simulationInOut; + + for (std::size_t index = 0; index < num_cells; ++index) { + if (index < static_cast(num_cells / 2)) { + simulationInOut.push_back(std::vector( + matrix_values.begin(), matrix_values.begin() + num_columns)); + } else { + simulationInOut.push_back(std::vector( + matrix_values.begin() + num_columns, matrix_values.end())); + } + } + EXPECT_NO_THROW(runner.run(simulationInOut, 100)); + + constexpr std::size_t half_cells = num_cells / 2; + constexpr int expected_value_first_half = 2; + constexpr int expected_value_second_half = 3; + + for (std::size_t cell_index = 0; cell_index < simulationInOut.size(); + ++cell_index) { + 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])); + } else { + EXPECT_EQ(simulationInOut[cell_index][0], expected_value_second_half); + EXPECT_FALSE(std::isnan(simulationInOut[cell_index][11])); + } + + EXPECT_NEAR(simulationInOut[cell_index][1], 111, 1); + EXPECT_NEAR(simulationInOut[cell_index][2], 55, 1); + } +} + +POET_TEST(PhreeqcRunnerUnknownID) { + PhreeqcMatrix pqc_mat(test_database, test_script); + PhreeqcRunner runner(pqc_mat); + + std::vector> simulationInOut; + + for (std::size_t index = 0; index < num_cells; ++index) { + simulationInOut.push_back(std::vector(pqc_mat.get().names.size())); + } + + simulationInOut[0][0] = 1000; + + EXPECT_THROW(runner.run(simulationInOut, 100), std::out_of_range); +} \ No newline at end of file