feat(eos): EOS now uses composition module
EOS code updated to make use of composition module. This is part of the larger change to change all composition handling to use the composition module. Note that the current implimentation is a bit hacky since it simply copies back and forth to the alredy used HELMEOSInput and HELMEOSOutput structs. In furture this can be more tightly connected to avoid extra copies.
This commit is contained in:
@@ -25,6 +25,7 @@
|
||||
#include <unordered_map>
|
||||
#include <vector>
|
||||
#include <array>
|
||||
#include <ranges>
|
||||
|
||||
#include <utility>
|
||||
|
||||
@@ -32,13 +33,15 @@
|
||||
|
||||
namespace serif::composition {
|
||||
|
||||
CompositionEntry::CompositionEntry() : m_symbol("H-1"), m_isotope(chemSpecies::species.at("H-1")), m_initialized(false) {}
|
||||
CompositionEntry::CompositionEntry() : m_symbol("H-1"), m_isotope(chemSpecies::species.at("H-1")),
|
||||
m_initialized(false) {
|
||||
}
|
||||
|
||||
CompositionEntry::CompositionEntry(const std::string& symbol, bool massFracMode) : m_symbol(symbol), m_isotope(chemSpecies::species.at(symbol)), m_massFracMode(massFracMode) {
|
||||
CompositionEntry::CompositionEntry(const std::string& symbol, const bool massFracMode) : m_symbol(symbol), m_isotope(chemSpecies::species.at(symbol)), m_massFracMode(massFracMode) {
|
||||
setSpecies(symbol);
|
||||
}
|
||||
}
|
||||
|
||||
CompositionEntry::CompositionEntry(const CompositionEntry& entry) :
|
||||
CompositionEntry::CompositionEntry(const CompositionEntry& entry) :
|
||||
m_symbol(entry.m_symbol),
|
||||
m_isotope(entry.m_isotope),
|
||||
m_massFracMode(entry.m_massFracMode),
|
||||
@@ -47,7 +50,7 @@ CompositionEntry::CompositionEntry(const CompositionEntry& entry) :
|
||||
m_relAbundance(entry.m_relAbundance),
|
||||
m_initialized(entry.m_initialized) {}
|
||||
|
||||
void CompositionEntry::setSpecies(const std::string& symbol) {
|
||||
void CompositionEntry::setSpecies(const std::string& symbol) {
|
||||
if (m_initialized) {
|
||||
throw std::runtime_error("Composition entry is already initialized.");
|
||||
}
|
||||
@@ -57,100 +60,100 @@ void CompositionEntry::setSpecies(const std::string& symbol) {
|
||||
m_symbol = symbol;
|
||||
m_isotope = chemSpecies::species.at(symbol);
|
||||
m_initialized = true;
|
||||
}
|
||||
}
|
||||
|
||||
std::string CompositionEntry::symbol() const {
|
||||
std::string CompositionEntry::symbol() const {
|
||||
return m_symbol;
|
||||
}
|
||||
}
|
||||
|
||||
double CompositionEntry::mass_fraction() const {
|
||||
double CompositionEntry::mass_fraction() const {
|
||||
if (!m_massFracMode) {
|
||||
throw std::runtime_error("Composition entry is in number fraction mode.");
|
||||
}
|
||||
return m_massFraction;
|
||||
}
|
||||
}
|
||||
|
||||
double CompositionEntry::mass_fraction(double meanMolarMass) const {
|
||||
double CompositionEntry::mass_fraction(double meanMolarMass) const {
|
||||
if (m_massFracMode) {
|
||||
return m_massFraction;
|
||||
}
|
||||
return m_relAbundance / meanMolarMass;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
double CompositionEntry::number_fraction() const {
|
||||
double CompositionEntry::number_fraction() const {
|
||||
if (m_massFracMode) {
|
||||
throw std::runtime_error("Composition entry is in mass fraction mode.");
|
||||
}
|
||||
return m_numberFraction;
|
||||
}
|
||||
}
|
||||
|
||||
double CompositionEntry::number_fraction(double totalMoles) const {
|
||||
double CompositionEntry::number_fraction(double totalMoles) const {
|
||||
if (m_massFracMode) {
|
||||
return m_relAbundance / totalMoles;
|
||||
}
|
||||
return m_numberFraction;
|
||||
}
|
||||
}
|
||||
|
||||
double CompositionEntry::rel_abundance() const {
|
||||
double CompositionEntry::rel_abundance() const {
|
||||
return m_relAbundance;
|
||||
}
|
||||
}
|
||||
|
||||
chemSpecies::Species CompositionEntry::isotope() const {
|
||||
chemSpecies::Species CompositionEntry::isotope() const {
|
||||
return m_isotope;
|
||||
}
|
||||
}
|
||||
|
||||
void CompositionEntry::setMassFraction(double mass_fraction) {
|
||||
void CompositionEntry::setMassFraction(double mass_fraction) {
|
||||
if (!m_massFracMode) {
|
||||
throw std::runtime_error("Composition entry is in number fraction mode.");
|
||||
}
|
||||
m_massFraction = mass_fraction;
|
||||
m_relAbundance = m_massFraction / m_isotope.mass();
|
||||
}
|
||||
}
|
||||
|
||||
void CompositionEntry::setNumberFraction(double number_fraction) {
|
||||
void CompositionEntry::setNumberFraction(double number_fraction) {
|
||||
if (m_massFracMode) {
|
||||
throw std::runtime_error("Composition entry is in mass fraction mode.");
|
||||
}
|
||||
m_numberFraction = number_fraction;
|
||||
m_relAbundance = m_numberFraction * m_isotope.mass();
|
||||
}
|
||||
}
|
||||
|
||||
bool CompositionEntry::setMassFracMode(double meanParticleMass) {
|
||||
bool CompositionEntry::setMassFracMode(double meanParticleMass) {
|
||||
if (m_massFracMode) {
|
||||
return false;
|
||||
}
|
||||
m_massFracMode = true;
|
||||
m_massFraction = m_relAbundance / meanParticleMass;
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
bool CompositionEntry::setNumberFracMode(double specificNumberDensity) {
|
||||
bool CompositionEntry::setNumberFracMode(double specificNumberDensity) {
|
||||
if (!m_massFracMode) {
|
||||
return false;
|
||||
}
|
||||
m_massFracMode = false;
|
||||
m_numberFraction = m_relAbundance / specificNumberDensity;
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
bool CompositionEntry::getMassFracMode() const {
|
||||
bool CompositionEntry::getMassFracMode() const {
|
||||
return m_massFracMode;
|
||||
}
|
||||
}
|
||||
|
||||
Composition::Composition(const std::vector<std::string>& symbols) {
|
||||
Composition::Composition(const std::vector<std::string>& symbols) {
|
||||
for (const auto& symbol : symbols) {
|
||||
registerSymbol(symbol);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
Composition::Composition(const std::set<std::string>& symbols) {
|
||||
Composition::Composition(const std::set<std::string>& symbols) {
|
||||
for (const auto& symbol : symbols) {
|
||||
registerSymbol(symbol);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
Composition::Composition(const std::vector<std::string>& symbols, const std::vector<double>& fractions, bool massFracMode) : m_massFracMode(massFracMode) {
|
||||
Composition::Composition(const std::vector<std::string>& symbols, const std::vector<double>& fractions, bool massFracMode) : m_massFracMode(massFracMode) {
|
||||
if (symbols.size() != fractions.size()) {
|
||||
LOG_ERROR(m_logger, "The number of symbols and fractions must be equal.");
|
||||
throw std::runtime_error("The number of symbols and fractions must be equal.");
|
||||
@@ -170,9 +173,32 @@ Composition::Composition(const std::vector<std::string>& symbols, const std::vec
|
||||
}
|
||||
}
|
||||
finalize();
|
||||
}
|
||||
}
|
||||
|
||||
void Composition::registerSymbol(const std::string& symbol, bool massFracMode) {
|
||||
Composition::Composition(const Composition &composition) {
|
||||
m_finalized = composition.m_finalized;
|
||||
m_specificNumberDensity = composition.m_specificNumberDensity;
|
||||
m_meanParticleMass = composition.m_meanParticleMass;
|
||||
m_massFracMode = composition.m_massFracMode;
|
||||
m_registeredSymbols = composition.m_registeredSymbols;
|
||||
m_compositions = composition.m_compositions;
|
||||
}
|
||||
|
||||
Composition& Composition::operator=(const Composition &other) {
|
||||
if (this != &other) {
|
||||
m_finalized = other.m_finalized;
|
||||
m_specificNumberDensity = other.m_specificNumberDensity;
|
||||
m_meanParticleMass = other.m_meanParticleMass;
|
||||
m_massFracMode = other.m_massFracMode;
|
||||
m_registeredSymbols = other.m_registeredSymbols;
|
||||
m_compositions = other.m_compositions;
|
||||
// note: m_config remains bound to the same singleton, so we skip it
|
||||
}
|
||||
return *this;
|
||||
|
||||
}
|
||||
|
||||
void Composition::registerSymbol(const std::string& symbol, bool massFracMode) {
|
||||
if (!isValidSymbol(symbol)) {
|
||||
LOG_ERROR(m_logger, "Invalid symbol: {}", symbol);
|
||||
throw std::runtime_error("Invalid symbol.");
|
||||
@@ -197,26 +223,26 @@ void Composition::registerSymbol(const std::string& symbol, bool massFracMode) {
|
||||
CompositionEntry entry(symbol, m_massFracMode);
|
||||
m_compositions[symbol] = entry;
|
||||
LOG_INFO(m_logger, "Registered symbol: {}", symbol);
|
||||
}
|
||||
}
|
||||
|
||||
void Composition::registerSymbol(const std::vector<std::string>& symbols, bool massFracMode) {
|
||||
void Composition::registerSymbol(const std::vector<std::string>& symbols, bool massFracMode) {
|
||||
for (const auto& symbol : symbols) {
|
||||
registerSymbol(symbol, massFracMode);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
std::set<std::string> Composition::getRegisteredSymbols() const {
|
||||
std::set<std::string> Composition::getRegisteredSymbols() const {
|
||||
return m_registeredSymbols;
|
||||
}
|
||||
}
|
||||
|
||||
void Composition::validateComposition(const std::vector<double>& fractions) const {
|
||||
void Composition::validateComposition(const std::vector<double>& fractions) const {
|
||||
if (!isValidComposition(fractions)) {
|
||||
LOG_ERROR(m_logger, "Invalid composition.");
|
||||
throw std::runtime_error("Invalid composition.");
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
bool Composition::isValidComposition(const std::vector<double>& fractions) const {
|
||||
bool Composition::isValidComposition(const std::vector<double>& fractions) const {
|
||||
double sum = 0.0;
|
||||
for (const auto& fraction : fractions) {
|
||||
sum += fraction;
|
||||
@@ -227,13 +253,13 @@ bool Composition::isValidComposition(const std::vector<double>& fractions) const
|
||||
}
|
||||
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
bool Composition::isValidSymbol(const std::string& symbol) const {
|
||||
bool Composition::isValidSymbol(const std::string& symbol) const {
|
||||
return chemSpecies::species.count(symbol) > 0;
|
||||
}
|
||||
}
|
||||
|
||||
double Composition::setMassFraction(const std::string& symbol, const double& mass_fraction) {
|
||||
double Composition::setMassFraction(const std::string& symbol, const double& mass_fraction) {
|
||||
if (m_registeredSymbols.find(symbol) == m_registeredSymbols.end()) {
|
||||
LOG_ERROR(m_logger, "Symbol {} is not registered.", symbol);
|
||||
throw std::runtime_error("Symbol is not registered.");
|
||||
@@ -254,9 +280,9 @@ double Composition::setMassFraction(const std::string& symbol, const double& mas
|
||||
m_compositions.at(symbol).setMassFraction(mass_fraction);
|
||||
|
||||
return old_mass_fraction;
|
||||
}
|
||||
}
|
||||
|
||||
std::vector<double> Composition::setMassFraction(const std::vector<std::string>& symbols, const std::vector<double>& mass_fractions) {
|
||||
std::vector<double> Composition::setMassFraction(const std::vector<std::string>& symbols, const std::vector<double>& mass_fractions) {
|
||||
if (symbols.size() != mass_fractions.size()) {
|
||||
LOG_ERROR(m_logger, "The number of symbols and mass fractions must be equal.");
|
||||
throw std::runtime_error("The number of symbols and mass fractions must be equal.");
|
||||
@@ -268,9 +294,9 @@ std::vector<double> Composition::setMassFraction(const std::vector<std::string>&
|
||||
old_mass_fractions.push_back(setMassFraction(symbols[i], mass_fractions[i]));
|
||||
}
|
||||
return old_mass_fractions;
|
||||
}
|
||||
}
|
||||
|
||||
double Composition::setNumberFraction(const std::string& symbol, const double& number_fraction) {
|
||||
double Composition::setNumberFraction(const std::string& symbol, const double& number_fraction) {
|
||||
if (m_registeredSymbols.find(symbol) == m_registeredSymbols.end()) {
|
||||
LOG_ERROR(m_logger, "Symbol {} is not registered.", symbol);
|
||||
throw std::runtime_error("Symbol is not registered.");
|
||||
@@ -291,9 +317,9 @@ double Composition::setNumberFraction(const std::string& symbol, const double& n
|
||||
m_compositions.at(symbol).setNumberFraction(number_fraction);
|
||||
|
||||
return old_number_fraction;
|
||||
}
|
||||
}
|
||||
|
||||
std::vector<double> Composition::setNumberFraction(const std::vector<std::string>& symbols, const std::vector<double>& number_fractions) {
|
||||
std::vector<double> Composition::setNumberFraction(const std::vector<std::string>& symbols, const std::vector<double>& number_fractions) {
|
||||
if (symbols.size() != number_fractions.size()) {
|
||||
LOG_ERROR(m_logger, "The number of symbols and number fractions must be equal.");
|
||||
throw std::runtime_error("The number of symbols and number fractions must be equal.");
|
||||
@@ -305,9 +331,9 @@ std::vector<double> Composition::setNumberFraction(const std::vector<std::string
|
||||
old_number_fractions.push_back(setNumberFraction(symbols[i], number_fractions[i]));
|
||||
}
|
||||
return old_number_fractions;
|
||||
}
|
||||
}
|
||||
|
||||
bool Composition::finalize(bool norm) {
|
||||
bool Composition::finalize(bool norm) {
|
||||
bool finalized = false;
|
||||
if (m_massFracMode) {
|
||||
finalized = finalizeMassFracMode(norm);
|
||||
@@ -318,9 +344,9 @@ bool Composition::finalize(bool norm) {
|
||||
m_finalized = true;
|
||||
}
|
||||
return finalized;
|
||||
}
|
||||
}
|
||||
|
||||
bool Composition::finalizeMassFracMode(bool norm) {
|
||||
bool Composition::finalizeMassFracMode(bool norm) {
|
||||
std::vector<double> mass_fractions;
|
||||
mass_fractions.reserve(m_compositions.size());
|
||||
for (const auto& [_, entry] : m_compositions) {
|
||||
@@ -354,9 +380,9 @@ bool Composition::finalizeMassFracMode(bool norm) {
|
||||
}
|
||||
m_meanParticleMass = 1.0/m_specificNumberDensity;
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
bool Composition::finalizeNumberFracMode(bool norm) {
|
||||
bool Composition::finalizeNumberFracMode(bool norm) {
|
||||
std::vector<double> number_fractions;
|
||||
number_fractions.reserve(m_compositions.size());
|
||||
for (const auto& [_, entry] : m_compositions) {
|
||||
@@ -387,9 +413,9 @@ bool Composition::finalizeNumberFracMode(bool norm) {
|
||||
}
|
||||
m_specificNumberDensity = 1.0/m_meanParticleMass;
|
||||
return true;
|
||||
}
|
||||
}
|
||||
|
||||
Composition Composition::mix(const Composition& other, double fraction) const {
|
||||
Composition Composition::mix(const Composition& other, double fraction) const {
|
||||
if (!m_finalized || !other.m_finalized) {
|
||||
LOG_ERROR(m_logger, "Compositions have not both been finalized.");
|
||||
throw std::runtime_error("Compositions have not been finalized (Consider running .finalize()).");
|
||||
@@ -416,9 +442,9 @@ Composition Composition::mix(const Composition& other, double fraction) const {
|
||||
}
|
||||
mixedComposition.finalize();
|
||||
return mixedComposition;
|
||||
}
|
||||
}
|
||||
|
||||
double Composition::getMassFraction(const std::string& symbol) const {
|
||||
double Composition::getMassFraction(const std::string& symbol) const {
|
||||
if (!m_finalized) {
|
||||
LOG_ERROR(m_logger, "Composition has not been finalized.");
|
||||
throw std::runtime_error("Composition has not been finalized (Consider running .finalize()).");
|
||||
@@ -432,18 +458,18 @@ double Composition::getMassFraction(const std::string& symbol) const {
|
||||
} else {
|
||||
return m_compositions.at(symbol).mass_fraction(m_meanParticleMass);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
std::unordered_map<std::string, double> Composition::getMassFraction() const {
|
||||
std::unordered_map<std::string, double> Composition::getMassFraction() const {
|
||||
std::unordered_map<std::string, double> mass_fractions;
|
||||
for (const auto& [symbol, entry] : m_compositions) {
|
||||
mass_fractions[symbol] = getMassFraction(symbol);
|
||||
}
|
||||
return mass_fractions;
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
double Composition::getNumberFraction(const std::string& symbol) const {
|
||||
double Composition::getNumberFraction(const std::string& symbol) const {
|
||||
if (!m_finalized) {
|
||||
LOG_ERROR(m_logger, "Composition has not been finalized.");
|
||||
throw std::runtime_error("Composition has not been finalized (Consider running .finalize()).");
|
||||
@@ -457,17 +483,17 @@ double Composition::getNumberFraction(const std::string& symbol) const {
|
||||
} else {
|
||||
return m_compositions.at(symbol).number_fraction(m_specificNumberDensity);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
std::unordered_map<std::string, double> Composition::getNumberFraction() const {
|
||||
std::unordered_map<std::string, double> Composition::getNumberFraction() const {
|
||||
std::unordered_map<std::string, double> number_fractions;
|
||||
for (const auto& [symbol, entry] : m_compositions) {
|
||||
number_fractions[symbol] = getNumberFraction(symbol);
|
||||
}
|
||||
return number_fractions;
|
||||
}
|
||||
}
|
||||
|
||||
std::pair<CompositionEntry, GlobalComposition> Composition::getComposition(const std::string& symbol) const {
|
||||
std::pair<CompositionEntry, GlobalComposition> Composition::getComposition(const std::string& symbol) const {
|
||||
if (!m_finalized) {
|
||||
LOG_ERROR(m_logger, "Composition has not been finalized.");
|
||||
throw std::runtime_error("Composition has not been finalized (Consider running .finalize()).");
|
||||
@@ -477,17 +503,53 @@ std::pair<CompositionEntry, GlobalComposition> Composition::getComposition(const
|
||||
throw std::runtime_error("Symbol is not in the composition.");
|
||||
}
|
||||
return {m_compositions.at(symbol), {m_specificNumberDensity, m_meanParticleMass}};
|
||||
}
|
||||
}
|
||||
|
||||
std::pair<std::unordered_map<std::string, CompositionEntry>, GlobalComposition> Composition::getComposition() const {
|
||||
std::pair<std::unordered_map<std::string, CompositionEntry>, GlobalComposition> Composition::getComposition() const {
|
||||
if (!m_finalized) {
|
||||
LOG_ERROR(m_logger, "Composition has not been finalized.");
|
||||
throw std::runtime_error("Composition has not been finalized (Consider running .finalize()).");
|
||||
}
|
||||
return {m_compositions, {m_specificNumberDensity, m_meanParticleMass}};
|
||||
}
|
||||
}
|
||||
|
||||
Composition Composition::subset(const std::vector<std::string>& symbols, std::string method) const {
|
||||
double Composition::getMeanParticleMass() const {
|
||||
if (!m_finalized) {
|
||||
LOG_ERROR(m_logger, "Composition has not been finalized.");
|
||||
throw std::runtime_error("Composition has not been finalized (Consider running .finalize()).");
|
||||
}
|
||||
return m_meanParticleMass;
|
||||
}
|
||||
|
||||
double Composition::getMeanAtomicNumber() const {
|
||||
if (!m_finalized) {
|
||||
LOG_ERROR(m_logger, "Composition must be finalized before getting the mean atomic mass number.");
|
||||
throw std::runtime_error("Composition not finalized. Cannot retrieve mean atomic mass number.");
|
||||
}
|
||||
|
||||
double mean_A = 0.0;
|
||||
|
||||
// Loop through all registered species in the composition.
|
||||
for (const auto &val: m_compositions | std::views::values) {
|
||||
const CompositionEntry& entry = val;
|
||||
const chemSpecies::Species& species = entry.isotope();
|
||||
|
||||
const double mass_fraction = entry.mass_fraction();
|
||||
const double particle_mass_g = species.mass();
|
||||
const int mass_number = species.a();
|
||||
|
||||
// Avoid division by zero, though a valid species should have a positive mass.
|
||||
if (particle_mass_g > 0) {
|
||||
// Calculate the number fraction for this species.
|
||||
const double number_fraction = (mass_fraction / particle_mass_g) * m_meanParticleMass;
|
||||
mean_A += number_fraction * mass_number;
|
||||
}
|
||||
}
|
||||
|
||||
return mean_A;
|
||||
}
|
||||
|
||||
Composition Composition::subset(const std::vector<std::string>& symbols, std::string method) const {
|
||||
std::array<std::string, 2> methods = {"norm", "none"};
|
||||
|
||||
if (std::find(methods.begin(), methods.end(), method) == methods.end()) {
|
||||
@@ -514,9 +576,9 @@ Composition Composition::subset(const std::vector<std::string>& symbols, std::st
|
||||
}
|
||||
}
|
||||
return subsetComposition;
|
||||
}
|
||||
}
|
||||
|
||||
void Composition::setCompositionMode(bool massFracMode) {
|
||||
void Composition::setCompositionMode(bool massFracMode) {
|
||||
if (!m_finalized) {
|
||||
LOG_ERROR(m_logger, "Composition has not been finalized. Mode cannot be set unless composition is finalized.");
|
||||
throw std::runtime_error("Composition has not been finalized (Consider running .finalize()). The mode cannot be set unless the composition is finalized.");
|
||||
@@ -535,36 +597,36 @@ void Composition::setCompositionMode(bool massFracMode) {
|
||||
}
|
||||
}
|
||||
m_massFracMode = massFracMode;
|
||||
}
|
||||
}
|
||||
|
||||
bool Composition::hasSymbol(const std::string& symbol) const {
|
||||
bool Composition::hasSymbol(const std::string& symbol) const {
|
||||
return m_compositions.count(symbol) > 0;
|
||||
}
|
||||
}
|
||||
|
||||
/// OVERLOADS
|
||||
/// OVERLOADS
|
||||
|
||||
Composition Composition::operator+(const Composition& other) const {
|
||||
Composition Composition::operator+(const Composition& other) const {
|
||||
return mix(other, 0.5);
|
||||
}
|
||||
}
|
||||
|
||||
std::ostream& operator<<(std::ostream& os, const GlobalComposition& comp) {
|
||||
std::ostream& operator<<(std::ostream& os, const GlobalComposition& comp) {
|
||||
os << "Global Composition: \n";
|
||||
os << "\tSpecific Number Density: " << comp.specificNumberDensity << "\n";
|
||||
os << "\tMean Particle Mass: " << comp.meanParticleMass << "\n";
|
||||
return os;
|
||||
}
|
||||
}
|
||||
|
||||
std::ostream& operator<<(std::ostream& os, const CompositionEntry& entry) {
|
||||
std::ostream& operator<<(std::ostream& os, const CompositionEntry& entry) {
|
||||
os << "<" << entry.m_symbol << " : m_frac = " << entry.mass_fraction() << ">";
|
||||
return os;
|
||||
}
|
||||
}
|
||||
|
||||
std::ostream& operator<<(std::ostream& os, const Composition& composition) {
|
||||
std::ostream& operator<<(std::ostream& os, const Composition& composition) {
|
||||
os << "Composition: \n";
|
||||
for (const auto& [symbol, entry] : composition.m_compositions) {
|
||||
os << entry << "\n";
|
||||
}
|
||||
return os;
|
||||
}
|
||||
}
|
||||
|
||||
} // namespace serif::composition
|
||||
@@ -277,7 +277,7 @@ namespace serif::composition {
|
||||
* Composition comp(symbols);
|
||||
* @endcode
|
||||
*/
|
||||
Composition(const std::vector<std::string>& symbols);
|
||||
explicit Composition(const std::vector<std::string>& symbols);
|
||||
|
||||
/**
|
||||
* @brief Constructs a Composition with the given symbols as a set.
|
||||
@@ -288,7 +288,7 @@ namespace serif::composition {
|
||||
* Composition comp(symbols);
|
||||
* @endcode
|
||||
*/
|
||||
Composition(const std::set<std::string>& symbols);
|
||||
explicit Composition(const std::set<std::string>& symbols);
|
||||
|
||||
/**
|
||||
* @brief Constructs a Composition with the given symbols and mass fractions.
|
||||
@@ -304,6 +304,14 @@ namespace serif::composition {
|
||||
*/
|
||||
Composition(const std::vector<std::string>& symbols, const std::vector<double>& mass_fractions, bool massFracMode=true);
|
||||
|
||||
/**
|
||||
* @brief Constructs a Composition from another Composition.
|
||||
* @param composition The Composition to copy.
|
||||
*/
|
||||
Composition(const Composition& composition);
|
||||
|
||||
Composition& operator=(Composition const& other);
|
||||
|
||||
/**
|
||||
* @brief Registers a new symbol.
|
||||
* @param symbol The symbol to register.
|
||||
@@ -333,7 +341,7 @@ namespace serif::composition {
|
||||
* @brief Gets the registered symbols.
|
||||
* @return A set of registered symbols.
|
||||
*/
|
||||
std::set<std::string> getRegisteredSymbols() const;
|
||||
[[nodiscard]] std::set<std::string> getRegisteredSymbols() const;
|
||||
|
||||
/**
|
||||
* @brief Sets the mass fraction for a given symbol.
|
||||
@@ -390,40 +398,52 @@ namespace serif::composition {
|
||||
* @brief Gets the mass fractions of all compositions.
|
||||
* @return An unordered map of compositions with their mass fractions.
|
||||
*/
|
||||
std::unordered_map<std::string, double> getMassFraction() const;
|
||||
[[nodiscard]] std::unordered_map<std::string, double> getMassFraction() const;
|
||||
|
||||
/**
|
||||
* @brief Gets the mass fraction for a given symbol.
|
||||
* @param symbol The symbol to get the mass fraction for.
|
||||
* @return The mass fraction for the given symbol.
|
||||
*/
|
||||
double getMassFraction(const std::string& symbol) const;
|
||||
[[nodiscard]] double getMassFraction(const std::string& symbol) const;
|
||||
|
||||
/**
|
||||
* @brief Gets the number fraction for a given symbol.
|
||||
* @param symbol The symbol to get the number fraction for.
|
||||
* @return The number fraction for the given symbol.
|
||||
*/
|
||||
double getNumberFraction(const std::string& symbol) const;
|
||||
[[nodiscard]] double getNumberFraction(const std::string& symbol) const;
|
||||
|
||||
/**
|
||||
* @brief Gets the number fractions of all compositions.
|
||||
* @return An unordered map of compositions with their number fractions.
|
||||
*/
|
||||
std::unordered_map<std::string, double> getNumberFraction() const;
|
||||
[[nodiscard]] std::unordered_map<std::string, double> getNumberFraction() const;
|
||||
|
||||
/**
|
||||
* @brief Gets the composition entry and global composition for a given symbol.
|
||||
* @param symbol The symbol to get the composition for.
|
||||
* @return A pair containing the CompositionEntry and GlobalComposition for the given symbol.
|
||||
*/
|
||||
std::pair<CompositionEntry, GlobalComposition> getComposition(const std::string& symbol) const;
|
||||
[[nodiscard]] std::pair<CompositionEntry, GlobalComposition> getComposition(const std::string& symbol) const;
|
||||
|
||||
/**
|
||||
* @brief Gets all composition entries and the global composition.
|
||||
* @return A pair containing an unordered map of CompositionEntries and the GlobalComposition.
|
||||
*/
|
||||
std::pair<std::unordered_map<std::string, CompositionEntry>, GlobalComposition> getComposition() const;
|
||||
[[nodiscard]] std::pair<std::unordered_map<std::string, CompositionEntry>, GlobalComposition> getComposition() const;
|
||||
|
||||
/**
|
||||
* @brief Compute the mean particle mass of the composition.
|
||||
* @return Mean particle mass in g.
|
||||
*/
|
||||
[[nodiscard]] double getMeanParticleMass() const;
|
||||
|
||||
/**
|
||||
* @brief Compute the mean atomic mass number of the composition.
|
||||
* @return Mean atomic mass number.
|
||||
*/
|
||||
[[nodiscard]] double getMeanAtomicNumber() const;
|
||||
|
||||
/**
|
||||
* @brief Gets a subset of the composition.
|
||||
|
||||
@@ -1,15 +1,18 @@
|
||||
# Define the library
|
||||
eos_sources = files(
|
||||
'private/helm.cpp',
|
||||
'private/EOSio.cpp'
|
||||
'private/EOSio.cpp',
|
||||
'private/EOS.cpp'
|
||||
)
|
||||
|
||||
eos_headers = files(
|
||||
'public/helm.h',
|
||||
'public/EOSio.h'
|
||||
'public/EOSio.h',
|
||||
'public/EOS.h'
|
||||
)
|
||||
|
||||
dependencies = [
|
||||
composition_dep,
|
||||
const_dep,
|
||||
quill_dep,
|
||||
probe_dep,
|
||||
|
||||
64
src/eos/private/EOS.cpp
Normal file
64
src/eos/private/EOS.cpp
Normal file
@@ -0,0 +1,64 @@
|
||||
#include "EOS.h"
|
||||
#include "EOSio.h"
|
||||
#include "helm.h"
|
||||
#include <string>
|
||||
|
||||
namespace serif::eos {
|
||||
EOS::EOS(const EOSio& reader) : m_reader(reader) {}
|
||||
EOS::EOS(const std::string& filename, const EOSFormat format) : m_reader(EOSio(filename, format)) {}
|
||||
|
||||
EOSOutput EOS::get(const EOSInput& in) {
|
||||
EOSOutput output;
|
||||
if (getFormat() == EOSFormat::HELM) {
|
||||
helmholtz::HELMEOSInput q;
|
||||
q.T = in.temperature; // Temperature in K
|
||||
q.rho = in.density; // Density in g/cm^3
|
||||
|
||||
q.abar = in.composition.getMeanParticleMass(); // Mean atomic mass in g
|
||||
q.zbar = in.composition.getMeanAtomicNumber(); // Mean atomic number (dimensionless)
|
||||
|
||||
helmholtz::HELMEOSOutput tempOutput;
|
||||
tempOutput = helmholtz::get_helm_EOS(q, *std::get<std::unique_ptr<helmholtz::HELMTable>>(m_reader.getTable()));
|
||||
|
||||
output.electronFraction = tempOutput.ye;
|
||||
output.electronChemicalPotential = tempOutput.etaele;
|
||||
output.neutronExcessFraction = tempOutput.xnefer;
|
||||
|
||||
// --- Pressure Variables ---
|
||||
output.pressure.total = tempOutput.ptot;
|
||||
output.pressure.gas = tempOutput.pgas;
|
||||
output.pressure.radiation = tempOutput.ptot;
|
||||
output.pressure.dDensity = tempOutput.dpresdd;
|
||||
output.pressure.dTemperature = tempOutput.dpresdt;
|
||||
output.pressure.dMeanAtomicMassNumber = tempOutput.dpresda;
|
||||
output.pressure.dMeanAtomicNumber = tempOutput.dpresdz;
|
||||
|
||||
// --- Energy Variables ---
|
||||
output.energy.total = tempOutput.etot;
|
||||
output.energy.gas = tempOutput.egas;
|
||||
output.energy.radiation = tempOutput.erad;
|
||||
output.energy.dDensity = tempOutput.denerdd;
|
||||
output.energy.dTemperature = tempOutput.denerdt;
|
||||
output.energy.dMeanAtomicMassNumber = tempOutput.denerda;
|
||||
output.energy.dMeanAtomicNumber = tempOutput.denerdz;
|
||||
|
||||
// --- Entropy Variables ---
|
||||
output.entropy.total = tempOutput.stot;
|
||||
output.entropy.gas = tempOutput.sgas;
|
||||
output.entropy.radiation = tempOutput.srad;
|
||||
output.entropy.dDensity = tempOutput.dentrdd;
|
||||
output.entropy.dTemperature = tempOutput.dentrdt;
|
||||
output.entropy.dMeanAtomicMassNumber = tempOutput.dentrda;
|
||||
output.entropy.dMeanAtomicNumber = tempOutput.dentrdz;
|
||||
}
|
||||
return output;
|
||||
}
|
||||
|
||||
EOSFormat EOS::getFormat() const {
|
||||
return m_reader.getFormat();
|
||||
}
|
||||
|
||||
const EOSio& EOS::getReader() const {
|
||||
return m_reader;
|
||||
}
|
||||
}
|
||||
@@ -28,25 +28,32 @@
|
||||
#include <string>
|
||||
|
||||
namespace serif::eos {
|
||||
EOSio::EOSio(const std::string &filename) : m_filename(filename) {
|
||||
EOSio::EOSio(const std::string &filename, const EOSFormat format) : m_filename(filename), m_format(format){
|
||||
load();
|
||||
}
|
||||
|
||||
std::string EOSio::getFormat() const {
|
||||
EOSio::EOSio(const EOSio &other) {
|
||||
m_filename = other.m_filename;
|
||||
m_format = other.m_format;
|
||||
load();
|
||||
}
|
||||
|
||||
EOSFormat EOSio::getFormat() const {
|
||||
return m_format;
|
||||
}
|
||||
|
||||
|
||||
std::string EOSio::getFormatName() const {
|
||||
return FormatStringLookup.at(m_format);
|
||||
}
|
||||
|
||||
|
||||
EOSTable& EOSio::getTable() {
|
||||
return m_table;
|
||||
}
|
||||
|
||||
void EOSio::load() {
|
||||
// Load the EOS table from the file
|
||||
// For now, just set the format to HELM
|
||||
|
||||
m_format = "helm";
|
||||
if (m_format == "helm") {
|
||||
if (m_format == EOSFormat::HELM) {
|
||||
loadHelm();
|
||||
}
|
||||
}
|
||||
|
||||
@@ -240,7 +240,7 @@ namespace serif::eos::helmholtz {
|
||||
ion, radiation, electron-positron and Coulomb interaction
|
||||
and returns the calculated quantities in the input
|
||||
***/
|
||||
serif::eos::helmholtz::EOS get_helm_EOS(serif::eos::helmholtz::EOSInput &q, const serif::eos::helmholtz::HELMTable &table) {
|
||||
serif::eos::helmholtz::HELMEOSOutput get_helm_EOS(serif::eos::helmholtz::HELMEOSInput &q, const serif::eos::helmholtz::HELMTable &table) {
|
||||
serif::config::Config& config = serif::config::Config::getInstance();
|
||||
auto logFile = config.get<std::string>("EOS:Helm:LogFile", "log");
|
||||
serif::probe::LogManager& logManager = serif::probe::LogManager::getInstance();
|
||||
@@ -829,7 +829,7 @@ namespace serif::eos::helmholtz {
|
||||
double csound = clight * sqrt(gamma1/z);
|
||||
|
||||
// package in q:
|
||||
serif::eos::helmholtz::EOS eos;
|
||||
serif::eos::helmholtz::HELMEOSOutput eos;
|
||||
eos.ptot = ptot; eos.etot = etot; eos.stot = stot;
|
||||
eos.pgas = pgas; eos.egas = egas; eos.sgas = sgas;
|
||||
eos.prad = prad; eos.erad = erad; eos.srad = srad;
|
||||
|
||||
230
src/eos/public/EOS.h
Normal file
230
src/eos/public/EOS.h
Normal file
@@ -0,0 +1,230 @@
|
||||
#pragma once
|
||||
|
||||
#include "EOSio.h"
|
||||
#include "helm.h"
|
||||
#include <string>
|
||||
#include "composition.h"
|
||||
|
||||
namespace serif::eos {
|
||||
|
||||
/**
|
||||
* @brief Input parameters for an EOS calculation.
|
||||
*
|
||||
* This struct holds the necessary physical conditions (composition, density, temperature)
|
||||
* required to query the Equation of State.
|
||||
*/
|
||||
struct EOSInput {
|
||||
serif::composition::Composition composition; ///< The composition of the system.
|
||||
double density; ///< The density of the system in cgs (g/cm^3).
|
||||
double temperature; ///< The temperature of the system in cgs (K).
|
||||
};
|
||||
|
||||
/**
|
||||
* @brief Represents a thermodynamic parameter and its derivatives.
|
||||
*
|
||||
* This struct stores a specific thermodynamic quantity (e.g., pressure, energy, entropy),
|
||||
* its breakdown into gas and radiation components, and its partial derivatives
|
||||
* with respect to density, temperature, mean atomic mass number, and mean atomic number.
|
||||
* All values are in cgs units unless otherwise specified.
|
||||
*/
|
||||
struct EOSParameter {
|
||||
double total; ///< Total value of the parameter (gas + radiation) (cgs).
|
||||
double gas; ///< Gas contribution to the parameter (cgs).
|
||||
double radiation; ///< Radiation contribution to the parameter (cgs).
|
||||
|
||||
double dDensity; ///< Derivative of the total parameter with respect to density (cgs units / (g/cm^3)).
|
||||
double dTemperature; ///< Derivative of the total parameter with respect to temperature (cgs units / K).
|
||||
double dMeanAtomicMassNumber; ///< Derivative of the total parameter with respect to mean atomic mass number (Abar) (cgs units / (g/mol)).
|
||||
double dMeanAtomicNumber; ///< Derivative of the total parameter with respect to mean atomic number (Zbar) (cgs units / dimensionless).
|
||||
|
||||
std::string name; ///< Name of the parameter (e.g., "Pressure", "Energy", "Entropy").
|
||||
};
|
||||
|
||||
/**
|
||||
* @brief Output from an EOS calculation.
|
||||
*
|
||||
* This struct contains various thermodynamic quantities and their derivatives
|
||||
* calculated by the EOS for a given set of input conditions.
|
||||
* It includes fundamental properties like electron fraction and chemical potential,
|
||||
* as well as detailed breakdowns of pressure, energy, and entropy.
|
||||
* Additionally, it provides methods to calculate derived quantities like
|
||||
* susceptibilities, sound speed, adiabatic gradients, and specific heats.
|
||||
*/
|
||||
struct EOSOutput {
|
||||
double electronFraction{}; ///< Electron fraction (ye), dimensionless.
|
||||
double electronChemicalPotential{}; ///< Electron chemical potential (eta_e) in cgs (erg/g).
|
||||
double neutronExcessFraction{}; ///< Neutron excess fraction (xnefer), dimensionless.
|
||||
|
||||
EOSParameter pressure; ///< Pressure output data, including total, gas, radiation, and derivatives.
|
||||
EOSParameter energy; ///< Internal energy output data, including total, gas, radiation, and derivatives.
|
||||
EOSParameter entropy; ///< Entropy output data, including total, gas, radiation, and derivatives.
|
||||
|
||||
/**
|
||||
* @brief Calculates the temperature susceptibility (chi_T).
|
||||
* @return Temperature susceptibility, dimensionless.
|
||||
* @note Placeholder: Actual calculation needs to be implemented based on available EOS derivatives.
|
||||
* Typically, chi_T = (d ln P / d ln T)_rho.
|
||||
*/
|
||||
double chiTemperature();
|
||||
/**
|
||||
* @brief Calculates the density susceptibility (chi_rho).
|
||||
* @return Density susceptibility, dimensionless.
|
||||
* @note Placeholder: Actual calculation needs to be implemented based on available EOS derivatives.
|
||||
* Typically, chi_rho = (d ln P / d ln rho)_T.
|
||||
*/
|
||||
double chiRho();
|
||||
/**
|
||||
* @brief Calculates the adiabatic sound speed.
|
||||
* @return Sound speed in cgs (cm/s).
|
||||
* @note Placeholder: Actual calculation needs to be implemented based on available EOS derivatives.
|
||||
* Typically, c_s^2 = gamma1 * P / rho.
|
||||
*/
|
||||
double soundSpeed();
|
||||
/**
|
||||
* @brief Calculates the adiabatic temperature gradient (nabla_ad).
|
||||
* @return Adiabatic gradient, dimensionless.
|
||||
* @note Placeholder: Actual calculation needs to be implemented based on available EOS derivatives.
|
||||
* Typically, nabla_ad = (P * chi_T) / (rho * T * c_p * chi_rho).
|
||||
*/
|
||||
double adiabaticGradient();
|
||||
/**
|
||||
* @brief Calculates the first adiabatic index (Gamma1).
|
||||
* @return First adiabatic index, dimensionless.
|
||||
* @note Placeholder: Actual calculation needs to be implemented based on available EOS derivatives.
|
||||
* Typically, Gamma1 = (d ln P / d ln rho)_S.
|
||||
*/
|
||||
double gamma1();
|
||||
/**
|
||||
* @brief Calculates the second adiabatic index (Gamma2).
|
||||
* @return Second adiabatic index, dimensionless.
|
||||
* @note Placeholder: Actual calculation needs to be implemented based on available EOS derivatives.
|
||||
* Typically, Gamma2 / (Gamma2 - 1) = (d ln P / d ln T)_S.
|
||||
*/
|
||||
double gamma2();
|
||||
/**
|
||||
* @brief Calculates the third adiabatic index (Gamma3).
|
||||
* @return Third adiabatic index, dimensionless.
|
||||
* @note Placeholder: Actual calculation needs to be implemented based on available EOS derivatives.
|
||||
* Typically, Gamma3 - 1 = (d ln T / d ln rho)_S.
|
||||
*/
|
||||
double gamma3();
|
||||
/**
|
||||
* @brief Calculates the specific heat capacity at constant volume (c_v).
|
||||
* @return Specific heat capacity at constant volume in cgs (erg/K/g).
|
||||
* @note Placeholder: Actual calculation needs to be implemented based on available EOS derivatives.
|
||||
* Typically, c_v = (dE/dT)_rho.
|
||||
*/
|
||||
double specificHeatCapacityAtConstantVolume();
|
||||
/**
|
||||
* @brief Calculates the specific heat capacity at constant pressure (c_p).
|
||||
* @return Specific heat capacity at constant pressure in cgs (erg/K/g).
|
||||
* @note Placeholder: Actual calculation needs to be implemented based on available EOS derivatives.
|
||||
* Typically, c_p = c_v + (T / rho^2) * ( (d P / d T)_rho^2 / (d P / d rho)_T ).
|
||||
*/
|
||||
double specificHeatCapacityAtConstantPressure();
|
||||
|
||||
/**
|
||||
* @brief Returns the format of the EOS data used to generate this output.
|
||||
* @return The EOSFormat enum value (currently only EOSFormat::HELM).
|
||||
*/
|
||||
EOSFormat EOSFormat() const;
|
||||
};
|
||||
|
||||
/**
|
||||
* @class EOS
|
||||
* @brief Main class for accessing Equation of State data.
|
||||
*
|
||||
* This class provides an interface to an underlying EOS table (e.g., Helmholtz EOS).
|
||||
* It handles loading the EOS data and provides a method to retrieve thermodynamic
|
||||
* properties for given physical conditions.
|
||||
*
|
||||
* @example
|
||||
* @code
|
||||
* #include "EOS.h"
|
||||
* #include "composition.h" // For serif::composition::Composition
|
||||
* #include <iostream>
|
||||
*
|
||||
* int main() {
|
||||
* try {
|
||||
* // Initialize EOS from a Helmholtz table file
|
||||
* serif::eos::EOS helmEOS("path/to/helm_table.dat", serif::eos::EOSFormat::HELM);
|
||||
*
|
||||
* // Define input conditions
|
||||
* serif::eos::EOSInput input;
|
||||
* input.density = 1.0e6; // g/cm^3
|
||||
* input.temperature = 1.0e7; // K
|
||||
* // Assuming a simple composition (e.g., pure Helium-4 for demonstration)
|
||||
* // In a real scenario, initialize Composition properly.
|
||||
* // For example, if Composition has a constructor like:
|
||||
* // Composition(const std::map<std::pair<int, int>, double>& mass_fractions);
|
||||
* // std::map<std::pair<int, int>, double> he4_mass_fraction = {{{2, 4}, 1.0}};
|
||||
* // input.composition = serif::composition::Composition(he4_mass_fraction);
|
||||
* // For now, let's assume Composition can be default constructed or set up simply:
|
||||
* input.composition.addSpecies(2, 4, 1.0); // Z=2, A=4 (He-4), mass fraction 1.0
|
||||
*
|
||||
* // Get EOS output
|
||||
* serif::eos::EOSOutput output = helmEOS.get(input);
|
||||
*
|
||||
* // Access results
|
||||
* std::cout << "Pressure (total): " << output.pressure.total << " dyne/cm^2" << std::endl;
|
||||
* std::cout << "Energy (total): " << output.energy.total << " erg/g" << std::endl;
|
||||
* std::cout << "Entropy (total): " << output.entropy.total << " erg/K/g" << std::endl;
|
||||
* std::cout << "Electron fraction: " << output.electronFraction << std::endl;
|
||||
*
|
||||
* // Example of accessing derivatives
|
||||
* std::cout << "dP/dRho: " << output.pressure.dDensity << std::endl;
|
||||
* std::cout << "dE/dT: " << output.energy.dTemperature << std::endl;
|
||||
*
|
||||
* } catch (const std::exception& e) {
|
||||
* std::cerr << "An error occurred: " << e.what() << std::endl;
|
||||
* return 1;
|
||||
* }
|
||||
* return 0;
|
||||
* }
|
||||
* @endcode
|
||||
*/
|
||||
class EOS {
|
||||
public:
|
||||
/**
|
||||
* @brief Constructs an EOS object by loading data from a file.
|
||||
* @param filename The path to the EOS data file.
|
||||
* @param format The format of the EOS data file (e.g., EOSFormat::HELM).
|
||||
* @throw std::runtime_error If the file cannot be opened or read, or if the format is unsupported.
|
||||
*/
|
||||
explicit EOS(const std::string& filename, EOSFormat format=EOSFormat::HELM);
|
||||
/**
|
||||
* @brief Constructs an EOS object from an existing EOSio reader.
|
||||
* @param reader An EOSio object that has already loaded the EOS data.
|
||||
*/
|
||||
explicit EOS(const EOSio& reader);
|
||||
/**
|
||||
* @brief Default destructor.
|
||||
*/
|
||||
~EOS() = default;
|
||||
|
||||
/**
|
||||
* @brief Retrieves thermodynamic properties for the given input conditions.
|
||||
* @param in An EOSInput struct containing the density, temperature, and composition.
|
||||
* @return An EOSOutput struct containing the calculated thermodynamic properties.
|
||||
* @throw std::runtime_error If the underlying EOS calculation fails (e.g., out of table bounds for Helmholtz).
|
||||
*
|
||||
* This method queries the loaded EOS table (e.g., Helmholtz) using the provided
|
||||
* density, temperature, and composition (mean atomic mass Abar, mean atomic number Zbar).
|
||||
* It populates and returns an EOSOutput struct with various thermodynamic quantities
|
||||
* such as pressure, energy, entropy, their derivatives, electron fraction, etc.
|
||||
*/
|
||||
[[nodiscard]] EOSOutput get(const EOSInput& in);
|
||||
/**
|
||||
* @brief Gets the format of the loaded EOS data.
|
||||
* @return The EOSFormat enum value.
|
||||
*/
|
||||
[[nodiscard]] EOSFormat getFormat() const;
|
||||
/**
|
||||
* @brief Gets a constant reference to the internal EOSio reader.
|
||||
* @return A const reference to the EOSio object.
|
||||
*/
|
||||
[[nodiscard]] const EOSio& getReader() const;
|
||||
private:
|
||||
EOSio m_reader; ///< The EOS I/O handler responsible for reading and storing EOS table data.
|
||||
};
|
||||
}
|
||||
@@ -22,15 +22,20 @@
|
||||
#include <string>
|
||||
#include <variant>
|
||||
#include <memory>
|
||||
#include <unordered_map>
|
||||
#include "helm.h"
|
||||
|
||||
namespace serif::eos {
|
||||
enum EOSFormat {
|
||||
HELM, ///< Helmholtz EOS format.
|
||||
};
|
||||
|
||||
static inline std::unordered_map<EOSFormat, std::string> FormatStringLookup = {
|
||||
{HELM, "Helmholtz"}
|
||||
};
|
||||
// EOS table format includes
|
||||
|
||||
using EOSTable = std::variant<
|
||||
std::unique_ptr<serif::eos::helmholtz::HELMTable>
|
||||
>;
|
||||
using EOSTable = std::variant<std::unique_ptr<serif::eos::helmholtz::HELMTable>>;
|
||||
|
||||
/**
|
||||
* @class EOSio
|
||||
@@ -41,16 +46,19 @@ namespace serif::eos {
|
||||
*
|
||||
* Example usage:
|
||||
* @code
|
||||
* EosIO eosIO("path/to/file");
|
||||
* std::string format = eosIO.getFormat();
|
||||
* EOSTable& table = eosIO.getTable();
|
||||
* EOSio eosReader("path/to/file");
|
||||
* std::string format = eosReader.getFormatName();
|
||||
* EOSTable& table = eosReader.getTable();
|
||||
* @endcode
|
||||
*
|
||||
* @note The default format used for reading tables is HELM
|
||||
* @note Currently only the HELM format is implemented
|
||||
*/
|
||||
class EOSio {
|
||||
private:
|
||||
std::string m_filename; ///< The filename of the EOS table.
|
||||
bool m_loaded = false; ///< Flag indicating if the table is loaded.
|
||||
std::string m_format; ///< The format of the EOS table.
|
||||
EOSFormat m_format;
|
||||
EOSTable m_table; ///< The EOS table data.
|
||||
|
||||
/**
|
||||
@@ -66,8 +74,15 @@ namespace serif::eos {
|
||||
/**
|
||||
* @brief Constructs an EosIO object with the given filename.
|
||||
* @param filename The filename of the EOS table.
|
||||
* @param format The EOS file format (currently only HELM)
|
||||
*/
|
||||
explicit EOSio(const std::string &filename);
|
||||
explicit EOSio(const std::string &filename, EOSFormat format = EOSFormat::HELM);
|
||||
|
||||
/**
|
||||
* @brief Explicit copy constructor
|
||||
* @param other The EOSio to be copied
|
||||
*/
|
||||
EOSio(const EOSio& other);
|
||||
|
||||
/**
|
||||
* @brief Default destructor.
|
||||
@@ -75,10 +90,12 @@ namespace serif::eos {
|
||||
~EOSio() = default;
|
||||
|
||||
/**
|
||||
* @brief Gets the format of the EOS table.
|
||||
* @brief Gets the format name (as a string) of the EOS table.
|
||||
* @return The format of the EOS table as a string.
|
||||
*/
|
||||
[[nodiscard]] std::string getFormat() const;
|
||||
[[nodiscard]] std::string getFormatName() const;
|
||||
|
||||
[[nodiscard]] EOSFormat getFormat() const;
|
||||
|
||||
/**
|
||||
* @brief Gets the EOS table.
|
||||
@@ -87,6 +104,8 @@ namespace serif::eos {
|
||||
[[nodiscard]] EOSTable& getTable();
|
||||
|
||||
[[nodiscard]] std::string getFilename() const { return m_filename; }
|
||||
|
||||
bool isLoaded() const { return m_loaded; }
|
||||
};
|
||||
}
|
||||
|
||||
|
||||
@@ -169,18 +169,18 @@ namespace serif::eos::helmholtz {
|
||||
};
|
||||
|
||||
/**
|
||||
* @struct EOSInput
|
||||
* @struct HELMEOSInput
|
||||
* @brief Structure to hold the input parameters for the EOS calculation.
|
||||
*/
|
||||
struct EOSInput
|
||||
struct HELMEOSInput
|
||||
{
|
||||
double T; ///< Temperature.
|
||||
double rho; ///< Density.
|
||||
double abar; ///< Mean atomic mass.
|
||||
double zbar; ///< Mean atomic number.
|
||||
|
||||
friend std::ostream& operator<<(std::ostream& os, const helmholtz::EOSInput& eosInput) {
|
||||
os << "EOSInput Data:\n";
|
||||
friend std::ostream& operator<<(std::ostream& os, const helmholtz::HELMEOSInput& eosInput) {
|
||||
os << "HELMEOSInput Data:\n";
|
||||
os << " Temperature: " << eosInput.T << "\n";
|
||||
os << " Density: " << eosInput.rho << "\n";
|
||||
os << " Mean Atomic Mass: " << eosInput.abar << "\n";
|
||||
@@ -194,7 +194,7 @@ namespace serif::eos::helmholtz {
|
||||
* @struct EOS
|
||||
* @brief Structure to hold the output parameters and derivatives of the EOS calculation.
|
||||
*/
|
||||
struct EOS
|
||||
struct HELMEOSOutput
|
||||
{
|
||||
// output
|
||||
double ye, etaele, xnefer; //
|
||||
@@ -212,7 +212,7 @@ namespace serif::eos::helmholtz {
|
||||
double gamma1, gamma2, gamma3, cV, cP; // derived quantities
|
||||
double dse, dpe, dsp; // Maxwell relations
|
||||
|
||||
friend std::ostream& operator<<(std::ostream& os, const helmholtz::EOS& eos) {
|
||||
friend std::ostream& operator<<(std::ostream& os, const helmholtz::HELMEOSOutput& eos) {
|
||||
os << "EOS Data:\n" << std::setw(20) << std::left;
|
||||
os << " Electron Fraction: " << std::format("{0:24.16e}",eos.ye) << "\n";
|
||||
os << " Electron Chemical Potential: " << std::format("{0:24.16e}",eos.etaele) << "\n";
|
||||
@@ -379,7 +379,7 @@ namespace serif::eos::helmholtz {
|
||||
* @param table HELMTable structure containing the table data.
|
||||
* @return EOS structure containing the calculated quantities.
|
||||
*/
|
||||
EOS get_helm_EOS(EOSInput &q, const HELMTable &table);
|
||||
HELMEOSOutput get_helm_EOS(HELMEOSInput &q, const HELMTable &table);
|
||||
|
||||
}
|
||||
|
||||
|
||||
@@ -34,7 +34,7 @@ void register_eos_bindings(pybind11::module &eos_submodule) {
|
||||
}, py::return_value_policy::reference_internal, // IMPORTANT: Keep this policy!
|
||||
"Get the EOS table data.")
|
||||
.def("__repr__", [](const serif::eos::EOSio &eos) {
|
||||
return "<EOSio(filename='" + eos.getFilename() + "', format='" + eos.getFormat() + "')>";
|
||||
return "<EOSio(filename='" + eos.getFilename() + "', format='" + eos.getFormatName() + "')>";
|
||||
});
|
||||
|
||||
py::class_<serif::eos::EOSTable>(eos_submodule, "EOSTable");
|
||||
@@ -88,64 +88,64 @@ void register_eos_bindings(pybind11::module &eos_submodule) {
|
||||
);
|
||||
}, py::return_value_policy::reference_internal); // Keep parent 'table' alive
|
||||
|
||||
py::class_<serif::eos::helmholtz::EOS>(eos_submodule, "EOS")
|
||||
py::class_<serif::eos::helmholtz::HELMEOSOutput>(eos_submodule, "EOS")
|
||||
.def(py::init<>())
|
||||
.def_readonly("ye", &serif::eos::helmholtz::EOS::ye)
|
||||
.def_readonly("etaele", &serif::eos::helmholtz::EOS::etaele)
|
||||
.def_readonly("xnefer", &serif::eos::helmholtz::EOS::xnefer)
|
||||
.def_readonly("ye", &serif::eos::helmholtz::HELMEOSOutput::ye)
|
||||
.def_readonly("etaele", &serif::eos::helmholtz::HELMEOSOutput::etaele)
|
||||
.def_readonly("xnefer", &serif::eos::helmholtz::HELMEOSOutput::xnefer)
|
||||
|
||||
.def_readonly("ptot", &serif::eos::helmholtz::EOS::ptot)
|
||||
.def_readonly("pgas", &serif::eos::helmholtz::EOS::pgas)
|
||||
.def_readonly("prad", &serif::eos::helmholtz::EOS::prad)
|
||||
.def_readonly("ptot", &serif::eos::helmholtz::HELMEOSOutput::ptot)
|
||||
.def_readonly("pgas", &serif::eos::helmholtz::HELMEOSOutput::pgas)
|
||||
.def_readonly("prad", &serif::eos::helmholtz::HELMEOSOutput::prad)
|
||||
|
||||
.def_readonly("etot", &serif::eos::helmholtz::EOS::etot)
|
||||
.def_readonly("egas", &serif::eos::helmholtz::EOS::egas)
|
||||
.def_readonly("erad", &serif::eos::helmholtz::EOS::erad)
|
||||
.def_readonly("etot", &serif::eos::helmholtz::HELMEOSOutput::etot)
|
||||
.def_readonly("egas", &serif::eos::helmholtz::HELMEOSOutput::egas)
|
||||
.def_readonly("erad", &serif::eos::helmholtz::HELMEOSOutput::erad)
|
||||
|
||||
.def_readonly("stot", &serif::eos::helmholtz::EOS::stot)
|
||||
.def_readonly("sgas", &serif::eos::helmholtz::EOS::sgas)
|
||||
.def_readonly("srad", &serif::eos::helmholtz::EOS::srad)
|
||||
.def_readonly("stot", &serif::eos::helmholtz::HELMEOSOutput::stot)
|
||||
.def_readonly("sgas", &serif::eos::helmholtz::HELMEOSOutput::sgas)
|
||||
.def_readonly("srad", &serif::eos::helmholtz::HELMEOSOutput::srad)
|
||||
|
||||
.def_readonly("dpresdd", &serif::eos::helmholtz::EOS::dpresdd)
|
||||
.def_readonly("dpresdt", &serif::eos::helmholtz::EOS::dpresdt)
|
||||
.def_readonly("dpresda", &serif::eos::helmholtz::EOS::dpresda)
|
||||
.def_readonly("dpresdz", &serif::eos::helmholtz::EOS::dpresdz)
|
||||
// TODO: Finish adding all the derivatives to the bound class
|
||||
.def_readonly("dentrdd", &serif::eos::helmholtz::EOS::dentrdd)
|
||||
.def_readonly("dentrdt", &serif::eos::helmholtz::EOS::dentrdt)
|
||||
.def_readonly("dentrda", &serif::eos::helmholtz::EOS::dentrda)
|
||||
.def_readonly("dentrdz", &serif::eos::helmholtz::EOS::dentrdz)
|
||||
.def_readonly("dpresdd", &serif::eos::helmholtz::HELMEOSOutput::dpresdd)
|
||||
.def_readonly("dpresdt", &serif::eos::helmholtz::HELMEOSOutput::dpresdt)
|
||||
.def_readonly("dpresda", &serif::eos::helmholtz::HELMEOSOutput::dpresda)
|
||||
.def_readonly("dpresdz", &serif::eos::helmholtz::HELMEOSOutput::dpresdz)
|
||||
|
||||
.def_readonly("denerdd", &serif::eos::helmholtz::EOS::denerdd)
|
||||
.def_readonly("denerdt", &serif::eos::helmholtz::EOS::denerdt)
|
||||
.def_readonly("denerda", &serif::eos::helmholtz::EOS::denerda)
|
||||
.def_readonly("denerdz", &serif::eos::helmholtz::EOS::denerdz)
|
||||
.def_readonly("dentrdd", &serif::eos::helmholtz::HELMEOSOutput::dentrdd)
|
||||
.def_readonly("dentrdt", &serif::eos::helmholtz::HELMEOSOutput::dentrdt)
|
||||
.def_readonly("dentrda", &serif::eos::helmholtz::HELMEOSOutput::dentrda)
|
||||
.def_readonly("dentrdz", &serif::eos::helmholtz::HELMEOSOutput::dentrdz)
|
||||
|
||||
.def_readonly("chiT", &serif::eos::helmholtz::EOS::chiT)
|
||||
.def_readonly("chiRho", &serif::eos::helmholtz::EOS::chiRho)
|
||||
.def_readonly("csound", &serif::eos::helmholtz::EOS::csound)
|
||||
.def_readonly("grad_ad", &serif::eos::helmholtz::EOS::grad_ad)
|
||||
.def_readonly("gamma1", &serif::eos::helmholtz::EOS::gamma1)
|
||||
.def_readonly("gamma2", &serif::eos::helmholtz::EOS::gamma2)
|
||||
.def_readonly("gamma3", &serif::eos::helmholtz::EOS::gamma3)
|
||||
.def_readonly("cV", &serif::eos::helmholtz::EOS::cV)
|
||||
.def_readonly("cP", &serif::eos::helmholtz::EOS::cP)
|
||||
.def_readonly("dse", &serif::eos::helmholtz::EOS::dse)
|
||||
.def_readonly("dpe", &serif::eos::helmholtz::EOS::dpe)
|
||||
.def_readonly("dsp", &serif::eos::helmholtz::EOS::dsp)
|
||||
.def_readonly("denerdd", &serif::eos::helmholtz::HELMEOSOutput::denerdd)
|
||||
.def_readonly("denerdt", &serif::eos::helmholtz::HELMEOSOutput::denerdt)
|
||||
.def_readonly("denerda", &serif::eos::helmholtz::HELMEOSOutput::denerda)
|
||||
.def_readonly("denerdz", &serif::eos::helmholtz::HELMEOSOutput::denerdz)
|
||||
|
||||
.def("__repr__", [](const serif::eos::helmholtz::EOS &eos) {
|
||||
.def_readonly("chiT", &serif::eos::helmholtz::HELMEOSOutput::chiT)
|
||||
.def_readonly("chiRho", &serif::eos::helmholtz::HELMEOSOutput::chiRho)
|
||||
.def_readonly("csound", &serif::eos::helmholtz::HELMEOSOutput::csound)
|
||||
.def_readonly("grad_ad", &serif::eos::helmholtz::HELMEOSOutput::grad_ad)
|
||||
.def_readonly("gamma1", &serif::eos::helmholtz::HELMEOSOutput::gamma1)
|
||||
.def_readonly("gamma2", &serif::eos::helmholtz::HELMEOSOutput::gamma2)
|
||||
.def_readonly("gamma3", &serif::eos::helmholtz::HELMEOSOutput::gamma3)
|
||||
.def_readonly("cV", &serif::eos::helmholtz::HELMEOSOutput::cV)
|
||||
.def_readonly("cP", &serif::eos::helmholtz::HELMEOSOutput::cP)
|
||||
.def_readonly("dse", &serif::eos::helmholtz::HELMEOSOutput::dse)
|
||||
.def_readonly("dpe", &serif::eos::helmholtz::HELMEOSOutput::dpe)
|
||||
.def_readonly("dsp", &serif::eos::helmholtz::HELMEOSOutput::dsp)
|
||||
|
||||
.def("__repr__", [](const serif::eos::helmholtz::HELMEOSOutput &eos) {
|
||||
return "<EOS (output from helmholtz eos)>";
|
||||
});
|
||||
|
||||
py::class_<serif::eos::helmholtz::EOSInput>(eos_submodule, "EOSInput")
|
||||
py::class_<serif::eos::helmholtz::HELMEOSInput>(eos_submodule, "HELMEOSInput")
|
||||
.def(py::init<>())
|
||||
.def_readwrite("T", &serif::eos::helmholtz::EOSInput::T)
|
||||
.def_readwrite("rho", &serif::eos::helmholtz::EOSInput::rho)
|
||||
.def_readwrite("abar", &serif::eos::helmholtz::EOSInput::abar)
|
||||
.def_readwrite("zbar", &serif::eos::helmholtz::EOSInput::zbar)
|
||||
.def("__repr__", [](const serif::eos::helmholtz::EOSInput &input) {
|
||||
return "<EOSInput(T=" + std::to_string(input.T) +
|
||||
.def_readwrite("T", &serif::eos::helmholtz::HELMEOSInput::T)
|
||||
.def_readwrite("rho", &serif::eos::helmholtz::HELMEOSInput::rho)
|
||||
.def_readwrite("abar", &serif::eos::helmholtz::HELMEOSInput::abar)
|
||||
.def_readwrite("zbar", &serif::eos::helmholtz::HELMEOSInput::zbar)
|
||||
.def("__repr__", [](const serif::eos::helmholtz::HELMEOSInput &input) {
|
||||
return "<HELMEOSInput(T=" + std::to_string(input.T) +
|
||||
", rho=" + std::to_string(input.rho) +
|
||||
", abar=" + std::to_string(input.abar) +
|
||||
", zbar=" + std::to_string(input.zbar) + ")>";
|
||||
|
||||
Reference in New Issue
Block a user