feat(eosIO): added EosIO class to handle arbitrary eos data

EosIO class wraps all eos tables (like helm) so that they can be used in a more standard fashion
This commit is contained in:
2025-03-20 14:25:22 -04:00
parent efa4bdadff
commit 171fbf7961
4 changed files with 232 additions and 46 deletions

View File

@@ -24,6 +24,7 @@
#include <iostream>
#include <fstream>
#include <memory>
#include <sstream>
#include <cmath>
#include <string>
@@ -121,21 +122,22 @@ namespace helmholtz {
}
// this function reads in the HELM table and stores in the above arrays
HELMTable read_helm_table(const std::string filename) {
std::unique_ptr<HELMTable> read_helm_table(const std::string filename) {
Config& config = Config::getInstance();
std::string logFile = config.get<std::string>("EOS:Helm:LogFile", "log");
Probe::LogManager& logManager = Probe::LogManager::getInstance();
quill::Logger* logger = logManager.getLogger(logFile);
LOG_INFO(logger, "read_helm_table : Reading HELM table from file {}", filename);
HELMTable table;
// Make a unique pointer to the HELMTable
std::unique_ptr<HELMTable> table = std::make_unique<HELMTable>();
string data;
int i, j;
//set T and Rho (d) arrays
for (j=0; j<table.jmax; j++) { table.t[j] = pow(10, tlo + tstp*j); }
for (j=0; j<table->jmax; j++) { table->t[j] = pow(10, tlo + tstp*j); }
for (i=0; i<table.imax; i++) { table.d[i] = pow(10, dlo + dstp*i); }
for (i=0; i<table->imax; i++) { table->d[i] = pow(10, dlo + dstp*i); }
ifstream helm_table(filename);
if (!helm_table) {
@@ -144,83 +146,83 @@ namespace helmholtz {
throw std::runtime_error("Error (" + std::to_string(errorCode) + ") opening file " + filename);
}
//read the Helmholtz free energy and its derivatives
for (j=0; j<table.jmax; j++) {
for (i=0; i<table.imax; i++){
for (j=0; j<table->jmax; j++) {
for (i=0; i<table->imax; i++){
getline(helm_table, data);
stringstream id(data);
id >> table.f[i][j];
id >> table.fd[i][j];
id >> table.ft[i][j];
id >> table.fdd[i][j];
id >> table.ftt[i][j];
id >> table.fdt[i][j];
id >> table.fddt[i][j];
id >> table.fdtt[i][j];
id >> table.fddtt[i][j];
id >> table->f[i][j];
id >> table->fd[i][j];
id >> table->ft[i][j];
id >> table->fdd[i][j];
id >> table->ftt[i][j];
id >> table->fdt[i][j];
id >> table->fddt[i][j];
id >> table->fdtt[i][j];
id >> table->fddtt[i][j];
}
}
//read the pressure derivative with density
for (j=0; j<table.jmax; j++) {
for (i=0; i<table.imax; i++){
for (j=0; j<table->jmax; j++) {
for (i=0; i<table->imax; i++){
getline(helm_table, data);
stringstream id(data);
id >> table.dpdf[i][j];
id >> table.dpdfd[i][j];
id >> table.dpdft[i][j];
id >> table.dpdfdt[i][j];
id >> table->dpdf[i][j];
id >> table->dpdfd[i][j];
id >> table->dpdft[i][j];
id >> table->dpdfdt[i][j];
}
}
//read the electron chemical potential
for (j=0; j<table.jmax; j++) {
for (i=0; i<table.imax; i++){
for (j=0; j<table->jmax; j++) {
for (i=0; i<table->imax; i++){
getline(helm_table, data);
stringstream id(data);
id >> table.ef[i][j];
id >> table.efd[i][j];
id >> table.eft[i][j];
id >> table.efdt[i][j];
id >> table->ef[i][j];
id >> table->efd[i][j];
id >> table->eft[i][j];
id >> table->efdt[i][j];
}
}
//read the number density
for (j=0; j<table.jmax; j++) {
for (i=0; i<table.imax; i++){
for (j=0; j<table->jmax; j++) {
for (i=0; i<table->imax; i++){
getline(helm_table, data);
stringstream id(data);
id >> table.xf[i][j];
id >> table.xfd[i][j];
id >> table.xft[i][j];
id >> table.xfdt[i][j];
id >> table->xf[i][j];
id >> table->xfd[i][j];
id >> table->xft[i][j];
id >> table->xfdt[i][j];
}
}
helm_table.close(); //done reading
// construct the temperature and density deltas and their inverses
for (j=0; j<table.jmax; j++) {
double dth = table.t[j+1] - table.t[j];
for (j=0; j<table->jmax; j++) {
double dth = table->t[j+1] - table->t[j];
double dt2 = dth * dth;
double dti = 1.0/dth;
double dt2i = 1.0/dt2;
table.dt_sav[j] = dth;
table.dt2_sav[j] = dt2;
table.dti_sav[j] = dti;
table.dt2i_sav[j] = dt2i;
table->dt_sav[j] = dth;
table->dt2_sav[j] = dt2;
table->dti_sav[j] = dti;
table->dt2i_sav[j] = dt2i;
}
for (i=0; i<table.imax; i++) {
double dd = table.d[i+1] - table.d[i];
for (i=0; i<table->imax; i++) {
double dd = table->d[i+1] - table->d[i];
double dd2 = dd * dd;
double ddi = 1.0/dd;
table.dd_sav[i] = dd;
table.dd2_sav[i] = dd2;
table.ddi_sav[i] = ddi;
table->dd_sav[i] = dd;
table->dd2_sav[i] = dd2;
table->ddi_sav[i] = ddi;
}
table.loaded = true;
table->loaded = true;
return table;
}