Skip to content

Commit

Permalink
Extend BinarySolutionTabulatedThermo by tabulated molar volumes
Browse files Browse the repository at this point in the history 
- Edit method docstrings in header file
- Change minor formatting issues
  • Loading branch information
@dschmider-HSOG
dschmider-HSOG committed Apr 13, 2022
1 parent 1fc652f commit 22c2045
Show file tree
Hide file tree
Showing 4 changed files with 52 additions and 50 deletions.
19 changes: 13 additions & 6 deletions include/cantera/thermo/BinarySolutionTabulatedThermo.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -205,7 +205,7 @@ class BinarySolutionTabulatedThermo : public IdealSolidSolnPhase
*
* @param vbar Output vector of partial molar volumes. Length: m_kk.
*/
virtual void getPartialMolarVolumes(doublereal* vbar) const;
virtual void getPartialMolarVolumes(double* vbar) const;

/**
* Overloads the calcDensity() method of IdealSolidSoln to also consider non-ideal
Expand All @@ -231,18 +231,25 @@ class BinarySolutionTabulatedThermo : public IdealSolidSolnPhase
* Tabulated values are only interpolated within the limits of the provided mole
* fraction. If these limits are exceeded, the values are capped at the lower or
* the upper limit.
*
* @param x Current mole fraction at which to interpolate.
* @param inputData Input vector of the data to be interpolated.
* @returns Linear interpolation of tabulated data at the current
* mole fraction x.
*/
double interpolate(double x, const vector_fp& molefrac,
const vector_fp& inputData) const;
double interpolate(const double x, const vector_fp& inputData) const;

//! Numerical derivation of the molar volume table
//! Numerical derivative of the molar volume table
/*!
* Tabulated values are only interpolated within the limits of the provided mole
* fraction. If these limits are exceeded, the values are capped at the lower or
* the upper limit.
*
* @param inputData Input vector of tabulated data to be derived.
* @param derivedData Output vector of tabulated data that is numerically
* derived with respect to the mole fraction.
*/
void derive(vector_fp& inputData, vector_fp& derivedData,
vector_fp& moleFraction);
void diff(const vector_fp& inputData, vector_fp& derivedData) const;

//! Current tabulated species index
size_t m_kk_tab;
Expand Down
6 changes: 3 additions & 3 deletions include/cantera/thermo/ThermoPhase.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -563,9 +563,9 @@ class ThermoPhase : public Phase
throw NotImplementedError("ThermoPhase::getPartialMolarVolumes");
}

//@}
/// @name Properties of the Standard State of the Species in the Solution
//@{
//! @}
//! @name Properties of the Standard State of the Species in the Solution
//! @{

//! Get the array of chemical potentials at unit activity for the species at
//! their standard states at the current *T* and *P* of the solution.
Expand Down
75 changes: 35 additions & 40 deletions src/thermo/BinarySolutionTabulatedThermo.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -49,8 +49,8 @@ void BinarySolutionTabulatedThermo::_updateThermo() const
bool x_changed = (m_xlast != xnow);

if (x_changed) {
m_h0_tab = interpolate(xnow, m_molefrac_tab, m_enthalpy_tab);
m_s0_tab = interpolate(xnow, m_molefrac_tab, m_entropy_tab);
m_h0_tab = interpolate(xnow, m_enthalpy_tab);
m_s0_tab = interpolate(xnow, m_entropy_tab);
if (xnow == 0) {
m_s0_tab = -BigNumber;
} else if (xnow == 1) {
Expand Down Expand Up @@ -79,6 +79,10 @@ void BinarySolutionTabulatedThermo::_updateThermo() const

bool BinarySolutionTabulatedThermo::addSpecies(shared_ptr<Species> spec)
{
if (m_kk == 2) {
throw CanteraError("BinarySolutionTabulatedThermo::addSpecies",
"No. of species should be equal to 2");
}
bool added = ThermoPhase::addSpecies(spec);
if (added) {
if (m_kk == 1) {
Expand All @@ -104,7 +108,7 @@ bool BinarySolutionTabulatedThermo::addSpecies(shared_ptr<Species> spec)
} else if (eos.hasKey("molar-volume")) {
mv = eos.convert("molar-volume", "m^3/kmol");
} else {
throw CanteraError("IdealSolidSolnPhase::addSpecies",
throw CanteraError("BinarySolutionTabulatedThermo::addSpecies",
"equation-of-state entry for species '{}' is missing "
"'density', 'molar-volume', or 'molar-density' "
"specification", spec->name);
Expand All @@ -114,7 +118,7 @@ bool BinarySolutionTabulatedThermo::addSpecies(shared_ptr<Species> spec)
// @Deprecated - remove this case for Cantera 3.0 with removal of the XML format
m_speciesMolarVolume.push_back(spec->input["molar_volume"].asDouble());
} else {
throw CanteraError("IdealSolidSolnPhase::addSpecies",
throw CanteraError("BinarySolutionTabulatedThermo::addSpecies",
"Molar volume not specified for species '{}'", spec->name);
}
}
Expand Down Expand Up @@ -142,13 +146,11 @@ void BinarySolutionTabulatedThermo::initThermo()
vector_fp h = table.convertVector("enthalpy", "J/kmol", N);
vector_fp s = table.convertVector("entropy", "J/kmol/K", N);
vector_fp vmol(N);
// Check for molarVolume key in tabulatedThermo table,
// Check for molar-volume key in tabulatedThermo table,
// otherwise calculate molar volume from pure species molar volumes
if (table.hasKey("molarVolume")) {
vmol = table.convertVector("molarVolume", "m^3/kmol", N);
}

else {
if (table.hasKey("molar-volume")) {
vmol = table.convertVector("molar-volume", "m^3/kmol", N);
} else {
for(size_t i = 0; i < N; i++) {
vmol[i] = x[i] * m_speciesMolarVolume[m_kk_tab] + (1-x[i])
* m_speciesMolarVolume[1-m_kk_tab];
Expand Down Expand Up @@ -182,7 +184,7 @@ void BinarySolutionTabulatedThermo::initThermo()
m_molar_volume_tab[i] = x_vmol[i].second;
}

derive(m_molar_volume_tab, m_derived_molar_volume_tab, m_molefrac_tab);
diff(m_molar_volume_tab, m_derived_molar_volume_tab);

for (size_t i = 0; i < N; i++) {
m_partial_molar_volume_1_tab[i] = m_molar_volume_tab[i] +
Expand Down Expand Up @@ -279,16 +281,15 @@ void BinarySolutionTabulatedThermo::initThermoXML(XML_Node& phaseNode, const std
m_molar_volume_tab[i] = x_vmol_temp[i].second;
}

derive(m_molar_volume_tab, m_derived_molar_volume_tab, m_molefrac_tab);
diff(m_molar_volume_tab, m_derived_molar_volume_tab);

for (size_t i = 0; i < x_h_temp.size(); i++) {
m_partial_molar_volume_1_tab[i] = m_molar_volume_tab[i] +
(1-m_molefrac_tab[i]) * m_derived_molar_volume_tab[i];
m_partial_molar_volume_2_tab[i] = m_molar_volume_tab[i] -
m_molefrac_tab[i] * m_derived_molar_volume_tab[i];
}
}
else {
} else {
throw CanteraError("BinarySolutionTabulatedThermo::initThermoXML",
"Unspecified tabulated species or thermo");
}
Expand Down Expand Up @@ -316,66 +317,60 @@ void BinarySolutionTabulatedThermo::initThermoXML(XML_Node& phaseNode, const std
ThermoPhase::initThermoXML(phaseNode, id_);
}

double BinarySolutionTabulatedThermo::interpolate(double x, const vector_fp& molefrac,
double BinarySolutionTabulatedThermo::interpolate(const double x,
const vector_fp& inputData) const
{
double c;
// Check if x is out of bound
if (x > molefrac.back()) {
if (x > m_molefrac_tab.back()) {
c = inputData.back();
return c;
}
if (x < molefrac.front()) {
if (x < m_molefrac_tab.front()) {
c = inputData.front();
return c;
}
size_t i = std::distance(molefrac.begin(),
std::lower_bound(molefrac.begin(), molefrac.end(), x));
size_t i = std::distance(m_molefrac_tab.begin(),
std::lower_bound(m_molefrac_tab.begin(), m_molefrac_tab.end(), x));
c = inputData[i-1] + (inputData[i] - inputData[i-1])
* (x - molefrac[i-1]) / (molefrac[i] - molefrac[i-1]);
* (x - m_molefrac_tab[i-1]) / (m_molefrac_tab[i] - m_molefrac_tab[i-1]);
return c;
}

void BinarySolutionTabulatedThermo::derive(vector_fp& inputData, vector_fp& derivedData,
vector_fp& moleFraction)
void BinarySolutionTabulatedThermo::diff(const vector_fp& inputData,
vector_fp& derivedData) const
{
if (inputData.size() > 1) {
derivedData[0] = (inputData[1] - inputData[0]) /
(moleFraction[1] - moleFraction[0]);
(m_molefrac_tab[1] - m_molefrac_tab[0]);
derivedData.back() = (inputData.back() - inputData[inputData.size()-2]) /
(moleFraction.back() - moleFraction[moleFraction.size()-2]);
(m_molefrac_tab.back() - m_molefrac_tab[m_molefrac_tab.size()-2]);

if (inputData.size() > 2) {
for (size_t i = 1; i < inputData.size()-1; i++) {
derivedData[i] = (inputData[i+1] - inputData[i-1]) /
(moleFraction[i+1] - moleFraction[i-1]);
(m_molefrac_tab[i+1] - m_molefrac_tab[i-1]);
}
}
}
else {
} else {
derivedData.front() = 0;
}
}

void BinarySolutionTabulatedThermo::getPartialMolarVolumes(doublereal* vbar) const
void BinarySolutionTabulatedThermo::getPartialMolarVolumes(double* vbar) const
{
vbar[m_kk_tab] = interpolate(Phase::moleFraction(m_kk_tab), m_molefrac_tab,
m_partial_molar_volume_1_tab);
vbar[1-m_kk_tab] = interpolate(Phase::moleFraction(m_kk_tab), m_molefrac_tab,
vbar[m_kk_tab] = interpolate(moleFraction(m_kk_tab), m_partial_molar_volume_1_tab);
vbar[1-m_kk_tab] = interpolate(moleFraction(m_kk_tab),
m_partial_molar_volume_2_tab);
}

void BinarySolutionTabulatedThermo::calcDensity()
{
{
double vmol = interpolate(Phase::moleFraction(m_kk_tab), m_molefrac_tab,
m_molar_volume_tab);
double dens = Phase::meanMolecularWeight()/vmol;

// Set the density in the parent State object directly, by calling the
// Phase::assignDensity() function.
Phase::assignDensity(dens);
}
double vmol = interpolate(moleFraction(m_kk_tab), m_molar_volume_tab);
double dens = meanMolecularWeight()/vmol;

// Set the density in the parent State object directly, by calling the
// Phase::assignDensity() function.
Phase::assignDensity(dens);
}
}
2 changes: 1 addition & 1 deletion test/data/BinarySolutionTabulatedThermo.yaml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -43,7 +43,7 @@ phases:
19.2971, 19.2977, 16.2213, 19.298, 19.2978, 19.2945, 19.2899, 19.2877,
19.2882, 19.2882, 19.2882, 19.2882, 19.2882, 19.2885, 19.2876, 19.2837,
19.2769, 19.285, 19.31, 19.3514]
molarVolume: [34.582988, 34.677730, 34.776275, 34.875760, 34.973980, 35.069268,
molar-volume: [34.582988, 34.677730, 34.776275, 34.875760, 34.973980, 35.069268,
35.160419, 35.246595, 35.327263, 35.402124, 35.471079, 35.534167, 35.591540,
35.643424, 35.690099, 35.731877, 35.769088, 35.802055, 35.831117, 35.856582,
35.878752, 35.897912, 35.914326, 35.928243, 35.939896, 35.949505, 35.957279,
Expand Down

0 comments on commit 22c2045

Please sign in to comment.