fix: re-enable DHT correctly

feat: store excess H and O instead of total values
This commit is contained in:
Max Lübke 2023-04-24 16:58:48 +02:00
parent ec47cbe8e7
commit f9da7830f2
15 changed files with 265 additions and 208 deletions

View File

@ -70,7 +70,7 @@ void writeFieldsToR(RInside &R, const Field &trans, const Field &chem) {
void set_chem_parameters(poet::ChemistryModule &chem, uint32_t wp_size,
const std::string &database_path) {
chem.SetErrorHandlerMode(1);
chem.SetComponentH2O(false);
chem.SetComponentH2O(true);
chem.SetRebalanceFraction(0.5);
chem.SetRebalanceByCell(true);
chem.UseSolutionDensityVolume(false);
@ -132,9 +132,10 @@ inline double RunMasterLoop(SimParams &params, RInside &R,
chem.initializeField(diffusion.getField());
if (params.getNumParams().dht_enabled) {
chem.SetDHTEnabled(true, params.getNumParams().dht_size_per_process);
if (!params.getDHTSignifVector().empty()) {
chem.SetDHTSignifVector(params.getDHTSignifVector());
chem.SetDHTEnabled(true, params.getNumParams().dht_size_per_process,
chem_params.dht_species);
if (!chem_params.dht_signif.empty()) {
chem.SetDHTSignifVector(chem_params.dht_signif);
}
if (!params.getDHTPropTypeVector().empty()) {
chem.SetDHTPropTypeVector(params.getDHTPropTypeVector());

View File

@ -1,4 +1,4 @@
## Time-stamp: "Last modified 2023-04-14 11:02:31 mluebke"
## Time-stamp: "Last modified 2023-04-24 16:51:23 mluebke"
database <- normalizePath("../share/poet/bench/barite/db_barite.dat")
input_script <- normalizePath("../share/poet/bench/barite/barite.pqi")
@ -44,9 +44,12 @@ grid <- list(
## initial conditions
init_diffu <- list(
"H" = 110.0124,
"O" = 55.5087,
"Charge" = -1.217E-09,
#"H" = 110.0124,
"H" = 0.00000028904,
#"O" = 55.5087,
"O" = 0.000000165205,
#"Charge" = -1.217E-09,
"Charge" = -3.337E-08,
"Ba" = 1.E-10,
"Cl" = 1.E-10,
"S(6)" = 6.205E-4,
@ -55,9 +58,12 @@ init_diffu <- list(
injection_diff <- list(
list(
"H" = 111.0124,
"O" = 55.50622,
"Charge" = -3.337E-08,
#"H" = 111.0124,
"H" = 0.0000002890408,
#"O" = 55.50622,
"O" = 0.00002014464,
#"Charge" = -3.337E-08,
"Charge" = -3.337000004885E-08,
"Ba" = 0.1,
"Cl" = 0.2,
"S(6)" = 0,
@ -109,13 +115,16 @@ diffusion <- list(
## # Needed when using DHT
signif_vector <- c(9, 9, 10, 5, 5, 5, 5, 5, 5)
prop_type <- c("", "", "", "act", "act", "act", "act", "", "")
prop <- names(init_cell)
dht_species <- names(init_diffu)
#dht_signif <- rep(15, length(dht_species))
dht_signif <- c(10, 10, 3, 5, 5, 5, 5)
chemistry <- list(
database = database,
input_script = input_script
input_script = input_script,
dht_species = dht_species,
dht_signif = dht_signif
)
#################################################################

View File

@ -1,4 +1,4 @@
## Time-stamp: "Last modified 2023-04-24 14:11:01 mluebke"
## Time-stamp: "Last modified 2023-04-24 15:12:02 mluebke"
database <- normalizePath("../share/poet/examples/phreeqc_kin.dat")
input_script <- normalizePath("../share/poet/bench/dolo/dolo_inner.pqi")
@ -44,8 +44,8 @@ grid <- list(
## initial conditions
init_diffu <- list(
"H" = 110.683,
"O" = 55.3413,
"H" = 0.000211313883539788,
"O" = 0.00398302904424952,
"Charge" = -5.0822e-19,
"C(4)" = 1.2279E-4,
"Ca" = 1.2279E-4,
@ -67,8 +67,8 @@ alpha_diffu <- c(
## list of boundary conditions/inner nodes
vecinj_diffu <- list(
list(
"H" = 110.683,
"O" = 55.3413,
"H" = 0.0001540445,
"O" = 0.002148006,
"Charge" = 1.90431e-16,
"C(4)" = 0,
"Ca" = 0,
@ -76,8 +76,8 @@ vecinj_diffu <- list(
"Mg" = 0.001
),
list(
"H" = 110.683,
"O" = 55.3413,
"H" = 0.0001610193,
"O" = 0.002386934,
"Charge" = 1.90431e-16,
"C(4)" = 0,
"Ca" = 0.0,
@ -120,13 +120,16 @@ diffusion <- list(
## # Needed when using DHT
signif_vector <- c(10, 10, 2, 5, 5, 5, 5, 0, 5, 5)
prop_type <- c("", "", "", "act", "act", "act", "act", "ignore", "", "")
prop <- names(init_cell)
dht_species <- c("H", "O", "Charge", "C(4)", "Ca", "Cl", "Mg", "Calcite",
"Dolomite")
#dht_signif <- rep(15, length(dht_species))
dht_signif <- c(10, 10, 3, 5, 5, 5, 5, 5, 5)
chemistry <- list(
database = database,
input_script = input_script
input_script = input_script,
dht_species = dht_species,
dht_signif = dht_signif
)
#################################################################
@ -135,7 +138,7 @@ chemistry <- list(
#################################################################
iterations <- 1000
iterations <- 10
dt <- 200
setup <- list(
@ -144,6 +147,5 @@ setup <- list(
chemistry = chemistry,
iterations = iterations,
timesteps = rep(dt, iterations),
store_result = TRUE,
out_save = c(5, iterations)
store_result = TRUE
)

View File

@ -1,7 +1,7 @@
## Time-stamp: "Last modified 2023-04-17 12:27:27 mluebke"
## Time-stamp: "Last modified 2023-04-24 15:43:26 mluebke"
database <- normalizePath("../share/poet/examples/phreeqc_kin.dat")
input_script <- normalizePath("../share/poet/bench/dolo_inner.pqi")
input_script <- normalizePath("../share/poet/bench/dolo/dolo_inner.pqi")
#################################################################
## Section 1 ##
@ -44,10 +44,10 @@ grid <- list(
## initial conditions
init_diffu <- list(
"H" = 110.683,
"O" = 55.3413,
"H" = 0.000211313883539788,
"O" = 0.00398302904424952,
"Charge" = -5.0822e-19,
"C" = 1.2279E-4,
"C(4)" = 1.2279E-4,
"Ca" = 1.2279E-4,
"Cl" = 0,
"Mg" = 0
@ -58,7 +58,7 @@ alpha_diffu <- c(
"H" = 1E-6,
"O" = 1E-6,
"Charge" = 1E-6,
"C" = 1E-6,
"C(4)" = 1E-6,
"Ca" = 1E-6,
"Cl" = 1E-6,
"Mg" = 1E-6
@ -67,19 +67,19 @@ alpha_diffu <- c(
## list of boundary conditions/inner nodes
vecinj_diffu <- list(
list(
"H" = 110.683,
"O" = 55.3413,
"H" = 0.0001540445,
"O" = 0.002148006,
"Charge" = 1.90431e-16,
"C" = 0,
"C(4)" = 0,
"Ca" = 0,
"Cl" = 0.002,
"Mg" = 0.001
),
list(
"H" = 110.683,
"O" = 55.3413,
"H" = 0.0001610193,
"O" = 0.002386934,
"Charge" = 1.90431e-16,
"C" = 0,
"C(4)" = 0,
"Ca" = 0.0,
"Cl" = 0.004,
"Mg" = 0.002
@ -120,13 +120,16 @@ diffusion <- list(
## # Needed when using DHT
signif_vector <- c(10, 10, 2, 5, 5, 5, 5, 0, 5, 5)
prop_type <- c("", "", "", "act", "act", "act", "act", "ignore", "", "")
prop <- names(init_cell)
dht_species <- c("H", "O", "Charge", "C(4)", "Ca", "Cl", "Mg", "Calcite",
"Dolomite")
#dht_signif <- rep(15, length(dht_species))
dht_signif <- c(10, 10, 3, 5, 5, 5, 5, 5, 5)
chemistry <- list(
database = database,
input_script = input_script
input_script = input_script,
dht_species = dht_species,
dht_signif = dht_signif
)
#################################################################
@ -145,5 +148,5 @@ setup <- list(
iterations = iterations,
timesteps = rep(dt, iterations),
store_result = TRUE,
out_save = c(5, iterations)
out_save = seq(5, iterations, by=5)
)

@ -1 +1 @@
Subproject commit 6a4c14c2857f0ad3558af3f1e1aa7555a805e8e2
Subproject commit 6db8938592a487a701b3bca19858ab946b42f1a0

View File

@ -1,3 +1,5 @@
// Time-stamp: "Last modified 2023-04-24 14:30:06 mluebke"
#ifndef CHEMISTRYMODULE_H_
#define CHEMISTRYMODULE_H_
@ -158,8 +160,10 @@ public:
*
* \param enable Enables or disables the usage of the DHT.
* \param size_mb Size in megabyte to allocate for the DHT if enabled.
* \param key_species Name of species to be used for key creation.
*/
void SetDHTEnabled(bool enable, uint32_t size_mb);
void SetDHTEnabled(bool enable, uint32_t size_mb,
const std::vector<std::string> &key_species);
/**
* **Master only** Set DHT snapshots to specific mode.
*
@ -174,7 +178,7 @@ public:
* \param signif_vec Vector defining significant digit for each species. Order
* is defined by prop_type vector (ChemistryModule::GetPropNames).
*/
void SetDHTSignifVector(std::vector<uint32_t> signif_vec);
void SetDHTSignifVector(std::vector<uint32_t> &signif_vec);
/**
* **Master only** Set the DHT rounding type of each species. See
* DHT_PROP_TYPES enumemartion for explanation.
@ -343,7 +347,7 @@ protected:
void WorkerMetricsToMaster(int type);
IRM_RESULT WorkerRunWorkPackage(std::vector<double> &vecWP,
std::vector<int32_t> &vecMapping,
std::vector<std::uint32_t> &vecMapping,
double dSimTime, double dTimestep);
std::vector<uint32_t> CalculateWPSizesVector(uint32_t n_cells,
@ -355,6 +359,8 @@ protected:
void unshuffleField(const std::vector<double> &in_buffer,
uint32_t size_per_prop, uint32_t prop_count,
uint32_t wp_count, std::vector<double> &out_field);
std::vector<std::uint32_t>
parseDHTSpeciesVec(const std::vector<std::string> &species_vec) const;
int comm_size, comm_rank;
MPI_Comm group_comm;
@ -362,12 +368,12 @@ protected:
bool is_sequential;
bool is_master;
uint32_t wp_size;
bool dht_enabled = false;
int dht_snaps_type = DHT_FILES_DISABLED;
uint32_t wp_size{CHEM_DEFAULT_WP_SIZE};
bool dht_enabled{false};
int dht_snaps_type{DHT_FILES_DISABLED};
std::string dht_file_out_dir;
poet::DHT_Wrapper *dht = std::nullptr_t();
poet::DHT_Wrapper *dht = nullptr;
static constexpr uint32_t BUFFER_OFFSET = 5;
@ -383,6 +389,8 @@ protected:
double seq_t = 0.;
double send_recv_t = 0.;
std::array<double, 2> base_totals{0};
#endif
double chem_t = 0.;

View File

@ -2,7 +2,7 @@
#define DHT_TYPES_H_
namespace poet {
enum DHT_PROP_TYPES { DHT_TYPE_DEFAULT, DHT_TYPE_ACT, DHT_TYPE_IGNORE };
enum DHT_PROP_TYPES { DHT_TYPE_DEFAULT, DHT_TYPE_CHARGE, DHT_TYPE_IGNORE, DHT_TYPE_TOTAL };
}
#endif // DHT_TYPES_H_

View File

@ -1,3 +1,5 @@
// Time-stamp: "Last modified 2023-04-24 16:23:42 mluebke"
/*
** Copyright (C) 2018-2021 Alexander Lindemann, Max Luebke (University of
** Potsdam)
@ -21,6 +23,8 @@
#ifndef DHT_WRAPPER_H
#define DHT_WRAPPER_H
#include "poet/DHT_Types.hpp"
#include <array>
#include <cstdint>
#include <string>
#include <vector>
@ -33,19 +37,21 @@ extern "C" {
namespace poet {
using DHT_Keyelement = struct keyelem {
struct DHT_SCNotation {
std::int8_t exp : 8;
std::int64_t significant : 56;
};
union DHT_Keyelement {
double fp_elemet;
DHT_SCNotation sc_notation;
};
using DHT_ResultObject = struct DHTResobj {
uint32_t length;
std::vector<std::vector<DHT_Keyelement>> keys;
std::vector<std::vector<double>> results;
std::vector<bool> needPhreeqc;
void ResultsToMapping(std::vector<int32_t> &curr_mapping);
void ResultsToWP(std::vector<double> &curr_wp);
};
/**
@ -79,11 +85,14 @@ public:
*
* @param dht_comm Communicator which addresses all participating DHT
* processes
* @param buckets_per_process Count of buckets to allocate for each process
* @param key_count Count of key entries
* @param buckets_per_process Count of buckets to allocate for each
* process
* @param key_indices Vector indexing elements of one grid cell used
* for key creation.
* @param data_count Count of data entries
*/
DHT_Wrapper(MPI_Comm dht_comm, uint32_t dht_size, uint32_t key_count,
DHT_Wrapper(MPI_Comm dht_comm, uint32_t dht_size,
const std::vector<std::uint32_t> &key_indices,
uint32_t data_count);
/**
* @brief Destroy the dht wrapper object
@ -112,8 +121,9 @@ public:
* @param[in,out] work_package Pointer to current work package
* @param dt Current timestep of simulation
*/
auto checkDHT(int length, double dt, const std::vector<double> &work_package)
-> poet::DHT_ResultObject;
auto checkDHT(int length, double dt, const std::vector<double> &work_package,
std::vector<std::uint32_t> &curr_mapping)
-> const poet::DHT_ResultObject &;
/**
* @brief Write simulated values into DHT
@ -131,8 +141,9 @@ public:
* outputs of the PHREEQC simulation
* @param dt Current timestep of simulation
*/
void fillDHT(int length, const DHT_ResultObject &dht_check_data,
const std::vector<double> &results);
void fillDHT(int length, const std::vector<double> &work_package);
void resultsToWP(std::vector<double> &work_package);
/**
* @brief Dump current DHT state into file.
@ -189,6 +200,10 @@ public:
void SetSignifVector(std::vector<uint32_t> signif_vec);
void SetPropTypeVector(std::vector<uint32_t> prop_type_vec);
void setBaseTotals(const std::array<double, 2> &bt) {
this->base_totals = bt;
}
private:
uint32_t key_count;
uint32_t data_count;
@ -202,11 +217,16 @@ private:
uint32_t dht_evictions = 0;
std::vector<uint32_t> dht_signif_vector;
std::vector<uint32_t> dht_prop_type_vector;
std::vector<std::uint32_t> dht_prop_type_vector;
std::vector<std::uint32_t> input_key_elements;
static constexpr int DHT_KEY_SIGNIF_DEFAULT = 5;
static constexpr int DHT_KEY_SIGNIF_TOTALS = 12;
static constexpr int DHT_KEY_SIGNIF_TOTALS = 10;
static constexpr int DHT_KEY_SIGNIF_CHARGE = 3;
DHT_ResultObject dht_results;
std::array<double, 2> base_totals{0};
};
} // namespace poet

View File

@ -101,6 +101,8 @@ using DiffusionParams = struct s_DiffusionParams {
using ChemistryParams = struct s_ChemistryParams {
std::string database_path;
std::string input_script;
std::vector<std::string> dht_species;
std::vector<std::uint32_t> dht_signif;
s_ChemistryParams(RInside &R);
};

View File

@ -3,7 +3,6 @@
#include "poet/DHT_Wrapper.hpp"
#include <algorithm>
#include <bits/stdint-uintn.h>
#include <cassert>
#include <cstdint>
#include <map>
@ -146,7 +145,7 @@ void poet::ChemistryModule::initializeField(const Field &trans_field) {
}
// now sort the new values
std::sort(new_solution_names.begin() + 3, new_solution_names.end());
std::sort(new_solution_names.begin() + 4, new_solution_names.end());
// and append other processes than solutions
std::vector<std::string> new_prop_names = new_solution_names;
@ -184,14 +183,28 @@ void poet::ChemistryModule::initializeField(const Field &trans_field) {
}
}
void poet::ChemistryModule::SetDHTEnabled(bool enable, uint32_t size_mb) {
void poet::ChemistryModule::SetDHTEnabled(
bool enable, uint32_t size_mb,
const std::vector<std::string> &key_species) {
constexpr uint32_t MB_FACTOR = 1E6;
std::vector<std::uint32_t> key_inidices;
if (this->is_master) {
int ftype = CHEM_DHT_ENABLE;
PropagateFunctionType(ftype);
ChemBCast(&enable, 1, MPI_CXX_BOOL);
ChemBCast(&size_mb, 1, MPI_UINT32_T);
key_inidices = parseDHTSpeciesVec(key_species);
int vec_size = key_inidices.size();
ChemBCast(&vec_size, 1, MPI_INT);
ChemBCast(key_inidices.data(), key_inidices.size(), MPI_UINT32_T);
} else {
int vec_size;
ChemBCast(&vec_size, 1, MPI_INT);
key_inidices.resize(vec_size);
ChemBCast(key_inidices.data(), vec_size, MPI_UINT32_T);
}
this->dht_enabled = enable;
@ -214,10 +227,39 @@ void poet::ChemistryModule::SetDHTEnabled(bool enable, uint32_t size_mb) {
const uint32_t dht_size = size_mb * MB_FACTOR;
this->dht =
new DHT_Wrapper(dht_comm, dht_size, this->prop_count, this->prop_count);
new DHT_Wrapper(dht_comm, dht_size, key_inidices, this->prop_count);
// this->dht->setBaseTotals(this->base_totals);
}
}
inline std::vector<std::uint32_t> poet::ChemistryModule::parseDHTSpeciesVec(
const std::vector<std::string> &species_vec) const {
std::vector<uint32_t> species_indices;
if (species_vec.empty()) {
species_indices.resize(this->prop_count);
int i = 0;
std::generate(species_indices.begin(), species_indices.end(),
[&] { return i++; });
return species_indices;
}
species_indices.reserve(species_vec.size());
for (const auto &name : species_vec) {
auto it = std::find(this->prop_names.begin(), this->prop_names.end(), name);
if (it == prop_names.end()) {
throw std::runtime_error(
"DHT species name was not found in prop name vector!");
}
const std::uint32_t index = it - prop_names.begin();
species_indices.push_back(index);
}
return species_indices;
}
void poet::ChemistryModule::SetDHTSnaps(int type, const std::string &out_dir) {
if (this->is_master) {
int ftype = CHEM_DHT_SNAPS;
@ -230,24 +272,28 @@ void poet::ChemistryModule::SetDHTSnaps(int type, const std::string &out_dir) {
ChemBCast(const_cast<char *>(out_dir.data()), str_size, MPI_CHAR);
}
this->dht_snaps_type = type;
this->dht_file_out_dir = out_dir;
}
void poet::ChemistryModule::SetDHTSignifVector(
std::vector<uint32_t> signif_vec) {
std::vector<uint32_t> &signif_vec) {
if (this->is_master) {
if (signif_vec.size() != this->prop_count) {
throw std::runtime_error(
"Significant vector sizes mismatches prop count.");
}
int ftype = CHEM_DHT_SIGNIF_VEC;
PropagateFunctionType(ftype);
int data_count = signif_vec.size();
ChemBCast(&data_count, 1, MPI_INT);
ChemBCast(signif_vec.data(), signif_vec.size(), MPI_UINT32_T);
return;
}
int data_count;
ChemBCast(&data_count, 1, MPI_INT);
signif_vec.resize(data_count);
ChemBCast(signif_vec.data(), data_count, MPI_UINT32_T);
this->dht->SetSignifVector(signif_vec);
}

View File

@ -35,41 +35,24 @@
using namespace poet;
using namespace std;
inline DHT_Keyelement round_key_element(double value, std::uint32_t signif) {
std::int8_t exp = (std::int8_t)std::floor(std::log10(std::fabs(value)));
std::int64_t significant = value * std::pow(10, signif - exp - 1);
inline DHT_SCNotation round_key_element(double value, std::uint32_t signif) {
std::int8_t exp =
static_cast<std::int8_t>(std::floor(std::log10(std::fabs(value))));
std::int64_t significant =
static_cast<std::int64_t>(value * std::pow(10, signif - exp - 1));
DHT_Keyelement elem;
DHT_SCNotation elem;
elem.exp = exp;
elem.significant = significant;
return elem;
}
void poet::DHT_ResultObject::ResultsToMapping(
std::vector<int32_t> &curr_mapping) {
uint32_t iMappingIndex = 0;
for (uint32_t i = 0; i < this->length; i++) {
if (curr_mapping[i] == -1) {
continue;
}
curr_mapping[i] = (this->needPhreeqc[i] ? iMappingIndex++ : -1);
}
}
void poet::DHT_ResultObject::ResultsToWP(std::vector<double> &curr_wp) {
for (uint32_t i = 0; i < this->length; i++) {
if (!this->needPhreeqc[i]) {
const uint32_t length = this->results[i].end() - this->results[i].begin();
std::copy(this->results[i].begin(), this->results[i].end(),
curr_wp.begin() + (length * i));
}
}
}
DHT_Wrapper::DHT_Wrapper(MPI_Comm dht_comm, uint32_t dht_size,
uint32_t key_count, uint32_t data_count)
: key_count(key_count), data_count(data_count) {
const std::vector<std::uint32_t> &key_indices,
uint32_t data_count)
: key_count(key_indices.size()), data_count(data_count),
input_key_elements(key_indices) {
// initialize DHT object
uint32_t key_size = (key_count + 1) * sizeof(DHT_Keyelement);
uint32_t data_size = data_count * sizeof(double);
@ -77,45 +60,39 @@ DHT_Wrapper::DHT_Wrapper(MPI_Comm dht_comm, uint32_t dht_size,
dht_object = DHT_create(dht_comm, buckets_per_process, data_size, key_size,
&poet::Murmur2_64A);
// extract needed values from sim_param struct
// t_simparams tmp = params.getNumParams();
// this->dt_differ = tmp.dt_differ;
// this->dht_log = tmp.dht_log;
// this->dht_signif_vector = params.getDHTSignifVector();
// this->dht_prop_type_vector = params.getDHTPropTypeVector();
this->dht_signif_vector.resize(key_size, DHT_KEY_SIGNIF_DEFAULT);
this->dht_signif_vector[0] = DHT_KEY_SIGNIF_TOTALS;
this->dht_signif_vector[1] = DHT_KEY_SIGNIF_TOTALS;
this->dht_signif_vector[2] = DHT_KEY_SIGNIF_CHARGE;
this->dht_prop_type_vector.resize(key_size, DHT_TYPE_DEFAULT);
this->dht_prop_type_vector.resize(key_count, DHT_TYPE_DEFAULT);
this->dht_prop_type_vector[0] = DHT_TYPE_TOTAL;
this->dht_prop_type_vector[1] = DHT_TYPE_TOTAL;
this->dht_prop_type_vector[2] = DHT_TYPE_CHARGE;
}
DHT_Wrapper::~DHT_Wrapper() {
// free DHT
DHT_free(dht_object, NULL, NULL);
}
auto DHT_Wrapper::checkDHT(int length, double dt,
const std::vector<double> &work_package)
-> poet::DHT_ResultObject {
const std::vector<double> &work_package,
std::vector<std::uint32_t> &curr_mapping)
-> const poet::DHT_ResultObject & {
DHT_ResultObject check_data;
dht_results.length = length;
dht_results.keys.resize(length);
dht_results.results.resize(length);
dht_results.needPhreeqc.resize(length);
check_data.length = length;
check_data.keys.resize(length);
check_data.results.resize(length);
check_data.needPhreeqc.resize(length);
std::vector<std::uint32_t> new_mapping;
// loop over every grid cell contained in work package
for (int i = 0; i < length; i++) {
// point to current grid cell
void *key = (void *)&(work_package[i * this->key_count]);
auto &data = check_data.results[i];
auto &key_vector = check_data.keys[i];
void *key = (void *)&(work_package[i * this->data_count]);
auto &data = dht_results.results[i];
auto &key_vector = dht_results.keys[i];
data.resize(this->data_count);
key_vector = fuzzForDHT(this->key_count, key, dt);
@ -125,27 +102,29 @@ auto DHT_Wrapper::checkDHT(int length, double dt,
switch (res) {
case DHT_SUCCESS:
check_data.needPhreeqc[i] = false;
dht_results.needPhreeqc[i] = false;
this->dht_hits++;
break;
case DHT_READ_MISS:
check_data.needPhreeqc[i] = true;
dht_results.needPhreeqc[i] = true;
new_mapping.push_back(curr_mapping[i]);
this->dht_miss++;
break;
}
}
return check_data;
curr_mapping = std::move(new_mapping);
return dht_results;
}
void DHT_Wrapper::fillDHT(int length, const DHT_ResultObject &dht_check_data,
const std::vector<double> &results) {
void DHT_Wrapper::fillDHT(int length, const std::vector<double> &work_package) {
// loop over every grid cell contained in work package
for (int i = 0; i < length; i++) {
// If true grid cell was simulated, needs to be inserted into dht
if (dht_check_data.needPhreeqc[i]) {
const auto &key = dht_check_data.keys[i];
void *data = (void *)&(results[i * this->data_count]);
if (dht_results.needPhreeqc[i]) {
const auto &key = dht_results.keys[i];
void *data = (void *)&(work_package[i * this->data_count]);
// fuzz data (round, logarithm etc.)
// insert simulated data with fuzzed key into DHT
@ -159,6 +138,15 @@ void DHT_Wrapper::fillDHT(int length, const DHT_ResultObject &dht_check_data,
}
}
void DHT_Wrapper::resultsToWP(std::vector<double> &work_package) {
for (int i = 0; i < dht_results.length; i++) {
if (!dht_results.needPhreeqc[i]) {
std::copy(dht_results.results[i].begin(), dht_results.results[i].end(),
work_package.begin() + (data_count * i));
}
}
}
int DHT_Wrapper::tableToFile(const char *filename) {
int res = DHT_to_file(dht_object, filename);
return res;
@ -194,25 +182,28 @@ std::vector<DHT_Keyelement> DHT_Wrapper::fuzzForDHT(int var_count, void *key,
std::vector<DHT_Keyelement> vecFuzz(var_count + 1);
std::memset(&vecFuzz[0], 0, sizeof(DHT_Keyelement) * var_count);
unsigned int i = 0;
int totals_i = 0;
// introduce fuzzing to allow more hits in DHT
// loop over every variable of grid cell
for (i = 0; i < (unsigned int)var_count; i++) {
double &curr_key = ((double *)key)[i];
for (std::uint32_t i = 0; i < input_key_elements.size(); i++) {
double curr_key = ((double *)key)[input_key_elements[i]];
if (curr_key != 0) {
if (curr_key < zero_val &&
this->dht_prop_type_vector[i] == DHT_TYPE_ACT) {
this->dht_prop_type_vector[i] == DHT_TYPE_DEFAULT) {
continue;
}
if (this->dht_prop_type_vector[i] == DHT_TYPE_IGNORE) {
continue;
}
vecFuzz[i] = round_key_element(curr_key, dht_signif_vector[i]);
if (this->dht_prop_type_vector[i] == DHT_TYPE_TOTAL) {
curr_key -= base_totals[totals_i++];
}
vecFuzz[i].sc_notation = round_key_element(curr_key, dht_signif_vector[i]);
}
}
// if timestep differs over iterations set current current time step at the
// end of fuzzing buffer
vecFuzz[var_count] = round_key_element(dt, 55);
vecFuzz[var_count].fp_elemet = dt;
return vecFuzz;
}

View File

@ -1,4 +1,4 @@
// Time-stamp: "Last modified 2023-03-31 15:00:11 mluebke"
// Time-stamp: "Last modified 2023-04-24 16:56:18 mluebke"
/*
**-----------------------------------------------------------------------------
** MurmurHash2 was written by Austin Appleby, and is placed in the public
@ -26,9 +26,6 @@
#include "poet/HashFunctions.hpp"
#include <cstddef>
#include <stdexcept>
#if defined(_MSC_VER)
#define BIG_CONSTANT(x) (x)

View File

@ -1,6 +1,10 @@
// Time-stamp: "Last modified 2023-04-24 16:55:34 mluebke"
#include "IrmResult.h"
#include "poet/ChemistryModule.hpp"
#include <algorithm>
#include <bits/stdint-uintn.h>
#include <cstddef>
#include <cstdint>
#include <iomanip>
@ -55,12 +59,13 @@ void poet::ChemistryModule::WorkerLoop() {
uint32_t size_mb;
ChemBCast(&size_mb, 1, MPI_UINT32_T);
SetDHTEnabled(enable, size_mb);
std::vector<std::string> name_dummy;
SetDHTEnabled(enable, size_mb, name_dummy);
break;
}
case CHEM_DHT_SIGNIF_VEC: {
std::vector<uint32_t> input_vec(this->prop_count);
ChemBCast(input_vec.data(), this->prop_count, MPI_UINT32_T);
std::vector<uint32_t> input_vec;
SetDHTSignifVector(input_vec);
break;
@ -185,46 +190,35 @@ void poet::ChemistryModule::WorkerDoWork(MPI_Status &probe_status,
/* 4th double value is currently a placeholder */
// placeholder = mpi_buffer[count+4];
// std::vector<double> vecCurrWP(
// mpi_buffer,
// mpi_buffer + (local_work_package_size * this->prop_names.size()));
vecCurrWP.resize(n_cells_times_props);
std::vector<int32_t> vecMappingWP(this->wp_size);
std::vector<std::uint32_t> vecMappingWP(local_work_package_size);
DHT_ResultObject DHT_Results;
for (uint32_t i = 0; i < local_work_package_size; i++) {
vecMappingWP[i] = i;
}
if (local_work_package_size != this->wp_size) {
// std::vector<double> vecFiller(
// (this->wp_size - local_work_package_size) * prop_count, 0);
// vecCurrWP.insert(vecCurrWP.end(), vecFiller.begin(), vecFiller.end());
// set all remaining cells to inactive
for (int i = local_work_package_size; i < this->wp_size; i++) {
vecMappingWP[i] = -1;
}
{
std::uint32_t i = 0;
std::generate(vecMappingWP.begin(), vecMappingWP.end(),
[&] { return i++; });
}
if (dht_enabled) {
/* check for values in DHT */
dht_get_start = MPI_Wtime();
DHT_Results = dht->checkDHT(local_work_package_size, dt, vecCurrWP);
dht->checkDHT(local_work_package_size, dt, vecCurrWP, vecMappingWP);
dht_get_end = MPI_Wtime();
DHT_Results.ResultsToMapping(vecMappingWP);
// DHT_Results.ResultsToMapping(vecMappingWP);
}
phreeqc_time_start = MPI_Wtime();
WorkerRunWorkPackage(vecCurrWP, vecMappingWP, current_sim_time, dt);
if (WorkerRunWorkPackage(vecCurrWP, vecMappingWP, current_sim_time, dt) !=
IRM_OK) {
throw std::runtime_error("Phreeqc threw an error!");
};
phreeqc_time_end = MPI_Wtime();
if (dht_enabled) {
DHT_Results.ResultsToWP(vecCurrWP);
dht->resultsToWP(vecCurrWP);
}
/* send results to master */
@ -235,7 +229,7 @@ void poet::ChemistryModule::WorkerDoWork(MPI_Status &probe_status,
if (dht_enabled) {
/* write results to DHT */
dht_fill_start = MPI_Wtime();
dht->fillDHT(local_work_package_size, DHT_Results, vecCurrWP);
dht->fillDHT(local_work_package_size, vecCurrWP);
dht_fill_end = MPI_Wtime();
timings.dht_get += dht_get_end - dht_get_start;
@ -260,37 +254,14 @@ void poet::ChemistryModule::WorkerPostIter(MPI_Status &prope_status,
}
}
void poet::ChemistryModule::WorkerPostSim(uint32_t iteration) {
/* before death, submit profiling/timings to master*/
// double timings_serialized[4];
// timings_serialized[0] = timings.phreeqc_t;
// timings_serialized[1] = timings.dht_get;
// timings_serialized[2] = timings.dht_fill;
// timings_serialized[3] = timings.idle_t;
// // timings
// MPI_Send(timings_serialized, 4, MPI_DOUBLE, 0, 0, this->group_comm);
// // MPI_Send(&phreeqc_count, 1, MPI_INT, 0, TAG_TIMING, MPI_COMM_WORLD);
// // MPI_Send(&idle_t, 1, MPI_DOUBLE, 0, TAG_TIMING, MPI_COMM_WORLD);
// if (this->dht_enabled) {
// // dht_perf
// int dht_perf[3];
// dht_perf[0] = dht->getHits();
// dht_perf[1] = dht->getMisses();
// dht_perf[2] = dht->getEvictions();
// MPI_Send(dht_perf, 3, MPI_INT, 0, 0, this->group_comm);
// }
if (this->dht_enabled && this->dht_snaps_type > DHT_FILES_SIMEND) {
if (this->dht_enabled && this->dht_snaps_type == DHT_FILES_SIMEND) {
WorkerWriteDHTDump(iteration);
}
}
void poet::ChemistryModule::WorkerWriteDHTDump(uint32_t iteration) {
std::stringstream out;
out << this->dht_file_out_dir << "/iter_" << std::setfill('0') << std::setw(4)
out << this->dht_file_out_dir << "/iter_" << std::setfill('0') << std::setw(3)
<< iteration << ".dht";
int res = dht->tableToFile(out.str().c_str());
if (res != DHT_SUCCESS && this->comm_rank == 2)
@ -321,13 +292,9 @@ void poet::ChemistryModule::WorkerReadDHTDump(
}
IRM_RESULT
poet::ChemistryModule::WorkerRunWorkPackage(std::vector<double> &vecWP,
std::vector<int32_t> &vecMapping,
double dSimTime, double dTimestep) {
if (this->wp_size != vecMapping.size()) {
return IRM_INVALIDARG;
}
poet::ChemistryModule::WorkerRunWorkPackage(
std::vector<double> &vecWP, std::vector<std::uint32_t> &vecMapping,
double dSimTime, double dTimestep) {
if ((this->wp_size * this->prop_count) != vecWP.size()) {
return IRM_INVALIDARG;
}
@ -346,7 +313,7 @@ poet::ChemistryModule::WorkerRunWorkPackage(std::vector<double> &vecWP,
}
IRM_RESULT result;
this->PhreeqcRM::CreateMapping(vecMapping);
this->PhreeqcRM::setPOETMapping(vecMapping);
this->setDumpedField(vecWP);
this->PhreeqcRM::SetTime(dSimTime);

View File

@ -21,6 +21,7 @@
#include "poet/DHT_Types.hpp"
#include <algorithm>
#include <cassert>
#include <cstdint>
#include <poet/SimParams.hpp>
#include <RInside.h>
@ -86,6 +87,16 @@ poet::ChemistryParams::s_ChemistryParams(RInside &R) {
Rcpp::as<std::string>(R.parseEval("mysetup$chemistry$database"));
this->input_script =
Rcpp::as<std::string>(R.parseEval("mysetup$chemistry$input_script"));
if (Rcpp::as<bool>(
R.parseEval("'dht_species' %in% names(mysetup$chemistry)"))) {
this->dht_species = Rcpp::as<std::vector<std::string>>(
R.parseEval("mysetup$chemistry$dht_species"));
}
if (Rcpp::as<bool>(
R.parseEval("'dht_signif' %in% names(mysetup$chemistry)"))) {
this->dht_signif = Rcpp::as<std::vector<std::uint32_t>>(
R.parseEval("mysetup$chemistry$dht_signif"));
}
}
SimParams::SimParams(int world_rank_, int world_size_) {
@ -225,7 +236,7 @@ void SimParams::initVectorParams(RInside &R) {
for (const auto &type : prop_type_R) {
if (type == "act") {
this->dht_prop_type_vector.push_back(DHT_TYPE_ACT);
this->dht_prop_type_vector.push_back(DHT_TYPE_CHARGE);
continue;
}

View File

@ -1,6 +1,6 @@
## Simple library of functions to assess and visualize the results of the coupled simulations
## Time-stamp: "Last modified 2023-01-18 16:02:58 delucia"
## Time-stamp: "Last modified 2023-04-24 16:09:55 mluebke"
require(RedModRphree)
require(Rmufits) ## essentially for PlotCartCellData
@ -59,10 +59,10 @@ ReadDHT <- function(file, new_scheme = TRUE) {
if (new_scheme) {
nkeys <- dims[1] / 8
keys <- res[, 1:nkeys]
keys <- res[, 1:nkeys - 1]
conv <- apply(keys, 2, ConvertDHTKey)
res[, 1:nkeys] <- conv
res[, 1:nkeys - 1] <- conv
}
return(res)