Skip to content

Commit

Permalink
Implement getParameters and getSpeciesParameters
Browse files Browse the repository at this point in the history
  • Loading branch information
@jongyoonbae @speth
jongyoonbae authored and speth committed Mar 18, 2023
1 parent 2b7b54e commit 4d3c83b
Show file tree
Hide file tree
Showing 3 changed files with 189 additions and 13 deletions.
41 changes: 36 additions & 5 deletions include/cantera/thermo/CoverageDependentSurfPhase.h
View file
Expand Up @@ -233,20 +233,25 @@ class CoverageDependentSurfPhase : public SurfPhase

virtual void initThermo();
virtual bool addSpecies(shared_ptr<Species> spec);
virtual void getParameters(AnyMap& phaseNode) const;
virtual void getSpeciesParameters(const std::string& name,
AnyMap& speciesNode) const;

//! @name Methods calculating reference state thermodynamic properties
//! Reference state properties are evaluated at \f$ T \text{ and }
//! \theta^{ref} \f$. With coverage fixed at a reference value,
//! reference state properties are effectively only dependent on temperature.

//! @{
virtual void getEnthalpy_RT_ref(double* hrt) const;
virtual void getEntropy_R_ref(double* sr) const;
virtual void getCp_R_ref(double* cpr) const;
virtual void getGibbs_RT_ref(double* grt) const;
//! @}

//! @name Methods calculating standard state thermodynamic properties
//! Standard state properties are evaluated at \f$ T \text{ and } \theta \f$,
//! and thus are dependent both on temperature and coverage.
//! @{

//! Get the nondimensionalized standard state enthalpy vector.
/*!
Expand Down Expand Up @@ -302,10 +307,12 @@ class CoverageDependentSurfPhase : public SurfPhase
* \f]
*/
virtual void getStandardChemPotentials(double* mu0) const;
//! @}

//! @name Methods calculating partial molar thermodynamic properties
//! Partial molar properties are evaluated at \f$ T \text{ and } \theta \f$,
//! and thus are dependent both on temperature and coverage.
//! @{

//! Get the partial molar enthalpy vector. Units: J/kmol.
/*!
Expand Down Expand Up @@ -338,11 +345,14 @@ class CoverageDependentSurfPhase : public SurfPhase
* \f]
*/
virtual void getChemPotentials(double* mu) const;
//! @}

//! @name Methods calculating Phase thermodynamic properties
//! Phase properties are evaluated at \f$ T \text{ and } \theta \f$,
//! and thus are dependent both on temperature and coverage.

//! @{

//! Return the solution's molar enthalpy. Units: J/kmol
/*!
* \f[
Expand All @@ -366,6 +376,7 @@ class CoverageDependentSurfPhase : public SurfPhase
* \f]
*/
virtual double cp_mole() const;
//! @}

protected:
//! Temporary storage for the coverages.
Expand Down Expand Up @@ -410,13 +421,33 @@ class CoverageDependentSurfPhase : public SurfPhase
//! Array of heat capacity coverage dependency parameters.
std::vector<HeatCapacityDependency> m_HeatCapacityDependency;

private:
//! Last value of the state number processed.
mutable int m_stateNumlast;
//! Explicitly-specified k-j interaction parameters, to enable serialization.
//! For linear model. <name_k>: <name_j: PolynomialDependency_index>>.
std::map<std::string, std::map<std::string, size_t>> m_indexmap_lin;

//! Explicitly-specified k-j interaction parameters, to enable serialization.
//! For polynomial model. <name_k>: <name_j: PolynomialDependency_index>>.
std::map<std::string, std::map<std::string, size_t>> m_indexmap_poly;

//! Explicitly-specified k-j interaction parameters, to enable serialization.
//! For piecewise-linear model. <name_k>: <name_j: InterpolativeDependency_index>>.
std::map<std::string, std::map<std::string, size_t>> m_indexmap_pwlin;

//! Explicitly-specified k-j interaction parameters, to enable serialization.
//! For interpolative model. <name_k>: <name_j: InterpolativeDependency_index>>.
std::map<std::string, std::map<std::string, size_t>> m_indexmap_int;

//! Explicitly-specified k-j interaction parameters, to enable serialization.
//! <name_k>: <name_j: HeatCapacityDependency_index>>.
std::map<std::string, std::map<std::string, size_t>> m_indexmap_cp;

private:
//! Storage for the user-defined reference state coverage which has to be
//! greater than 0.0 and less than or equal to 1.0. default = 1.0.
mutable double m_theta_ref;
double m_theta_ref;

//! Last value of the state number processed.
mutable int m_stateNumlast;

//! Update the species coverage-dependent thermodynamic functions.
/*!
Expand Down
130 changes: 122 additions & 8 deletions src/thermo/CoverageDependentSurfPhase.cpp
View file
Expand Up @@ -24,8 +24,8 @@ namespace Cantera

CoverageDependentSurfPhase::CoverageDependentSurfPhase(const std::string& infile,
const std::string& id_):
m_stateNumlast(-2),
m_theta_ref(1.0)
m_theta_ref(1.0),
m_stateNumlast(-2)
{
setNDim(2);
initThermoFile(infile, id_);
Expand Down Expand Up @@ -76,11 +76,7 @@ void CoverageDependentSurfPhase::addInterpolativeDependency(const

void CoverageDependentSurfPhase::initThermo()
{
if (m_input.hasKey("site-density")) {
// Units are kmol/m^2
setSiteDensity(m_input.convert("site-density",
Units(1.0, 0, -static_cast<double>(m_ndim), 0, 0, 0, 1)));
}
SurfPhase::initThermo();
if (m_input.hasKey("reference-state-coverage")) {
m_theta_ref = m_input["reference-state-coverage"].as<double>();
if (m_theta_ref <= 0.0 || m_theta_ref > 1.0) {
Expand All @@ -90,6 +86,8 @@ void CoverageDependentSurfPhase::initThermo()
}
}
for (auto& item : m_species) {
// Dependency index trackers to be used in serialization.
size_t poly_ind, int_ind, cp_ind;
// Read enthalpy and entropy dependencies from species 'input' information
// (i.e. as specified in a YAML input file) for both self- and cross-
// interactions.
Expand Down Expand Up @@ -119,6 +117,8 @@ void CoverageDependentSurfPhase::initThermo()
"J/kmol/K");
}

poly_ind = m_PolynomialDependency.size();
m_indexmap_lin[item.first][item2.first] = poly_ind;
m_PolynomialDependency.push_back(poly_deps);
// For polynomial(4th) model
} else if (cov_map2["model"] == "polynomial") {
Expand All @@ -136,8 +136,10 @@ void CoverageDependentSurfPhase::initThermo()
poly_deps.entropy_coeffs.begin(), 0.0);
}

poly_ind = m_PolynomialDependency.size();
m_indexmap_poly[item.first][item2.first] = poly_ind;
m_PolynomialDependency.push_back(poly_deps);
// For piecewise linear model
// For piecewise-linear model
} else if (cov_map2["model"] == "piecewise-linear") {
InterpolativeDependency int_deps(k, j);
if (cov_map2.hasKey("enthalpy-low") ||
Expand All @@ -161,6 +163,8 @@ void CoverageDependentSurfPhase::initThermo()
+ int_deps.entropy_map[cov_change];
}

int_ind = m_InterpolativeDependency.size();
m_indexmap_pwlin[item.first][item2.first] = int_ind;
addInterpolativeDependency(int_deps);
// For interpolative model
} else if (cov_map2["model"] == "interpolative") {
Expand Down Expand Up @@ -196,6 +200,8 @@ void CoverageDependentSurfPhase::initThermo()
}
}

int_ind = m_InterpolativeDependency.size();
m_indexmap_int[item.first][item2.first] = int_ind;
addInterpolativeDependency(int_deps);
} else {
throw InputFileError("CoverageDependentSurfPhase::initThermo",
Expand All @@ -209,6 +215,8 @@ void CoverageDependentSurfPhase::initThermo()
cpcov_deps.coeff_a = cov_map2.convert("heat-capacity-a", "J/kmol/K");
cpcov_deps.coeff_b = cov_map2.convert("heat-capacity-b", "J/kmol/K");

cp_ind = m_HeatCapacityDependency.size();
m_indexmap_cp[item.first][item2.first] = cp_ind;
m_HeatCapacityDependency.push_back(cpcov_deps);
}
}
Expand All @@ -233,6 +241,112 @@ bool CoverageDependentSurfPhase::addSpecies(shared_ptr<Species> spec)
return added;
}

void CoverageDependentSurfPhase::getParameters(AnyMap& phaseNode) const
{
SurfPhase::getParameters(phaseNode);
phaseNode["reference-state-coverage"] = m_theta_ref;
}

void CoverageDependentSurfPhase::getSpeciesParameters(const std::string& name,
AnyMap& speciesNode) const
{
SurfPhase::getSpeciesParameters(name, speciesNode);
// Get linear model parameters from PolynomialDependency vector
if (m_indexmap_lin.count(name)) {
for (const auto& item: m_indexmap_lin.at(name)) {
auto& covdepNode =
speciesNode["coverage-dependencies"][item.first].getMapWhere(
"model", "linear", true);
covdepNode["enthalpy"].setQuantity(
m_PolynomialDependency[item.second].enthalpy_coeffs[1], "J/kmol");
covdepNode["entropy"].setQuantity(
m_PolynomialDependency[item.second].entropy_coeffs[1], "J/kmol/K");
}
}
// Get polynomial model parameters from PolynomialDependency vector
if (m_indexmap_poly.count(name)) {
for (const auto& item: m_indexmap_poly.at(name)) {
auto& covdepNode =
speciesNode["coverage-dependencies"][item.first].getMapWhere(
"model", "polynomial", true);
vector_fp hvec =
vector_fp(
m_PolynomialDependency[item.second].enthalpy_coeffs.begin() + 1,
m_PolynomialDependency[item.second].enthalpy_coeffs.end());
covdepNode["enthalpy-coefficients"].setQuantity(hvec, "J/kmol");
vector_fp svec =
vector_fp(
m_PolynomialDependency[item.second].entropy_coeffs.begin() + 1,
m_PolynomialDependency[item.second].entropy_coeffs.end());
covdepNode["entropy-coefficients"].setQuantity(svec, "J/kmol/K");
}
}
// Get piecewise-linear model parameters from InterpolativeDependency vector
if (m_indexmap_pwlin.count(name)) {
for (const auto& item: m_indexmap_pwlin.at(name)) {
auto& covdepNode =
speciesNode["coverage-dependencies"][item.first].getMapWhere(
"model", "piecewise-linear", true);
vector_fp hcovs, enthalpies, scovs, entropies;
for (const auto& hmap:
m_InterpolativeDependency[item.second].enthalpy_map) {
hcovs.push_back(hmap.first);
enthalpies.push_back(hmap.second);
}
for (const auto& smap:
m_InterpolativeDependency[item.second].entropy_map) {
scovs.push_back(smap.first);
entropies.push_back(smap.second);
}
covdepNode["enthalpy-change"] = hcovs[1];
covdepNode["enthalpy-low"].setQuantity(
(enthalpies[1] - enthalpies[0]) / (hcovs[1] - hcovs[0]), "J/kmol");
covdepNode["enthalpy-high"].setQuantity(
(enthalpies[2] - enthalpies[1]) / (hcovs[2] - hcovs[1]), "J/kmol");
covdepNode["entropy-change"] = scovs[1];
covdepNode["entropy-low"].setQuantity(
(entropies[1] - entropies[0]) / (scovs[1] - scovs[0]), "J/kmol/K");
covdepNode["entropy-high"].setQuantity(
(entropies[2] - entropies[1]) / (scovs[2] - scovs[1]), "J/kmol/K");
}
}
// Get interpolative model parameters from InterpolativeDependency vector
if (m_indexmap_int.count(name)) {
for (const auto& item: m_indexmap_int.at(name)) {
auto& covdepNode =
speciesNode["coverage-dependencies"][item.first].getMapWhere(
"model", "interpolative", true);
vector_fp hcovs, enthalpies, scovs, entropies;
for (const auto& hmap:
m_InterpolativeDependency[item.second].enthalpy_map) {
hcovs.push_back(hmap.first);
enthalpies.push_back(hmap.second);
}
for (const auto& smap:
m_InterpolativeDependency[item.second].entropy_map) {
scovs.push_back(smap.first);
entropies.push_back(smap.second);
}
covdepNode["enthalpy-coverages"] = std::move(hcovs);
covdepNode["enthalpies"].setQuantity(enthalpies, "J/kmol");
covdepNode["entropy-coverages"] = std::move(scovs);
covdepNode["entropies"].setQuantity(entropies, "J/kmol/K");
}
}
// Get heat capacity model parameters from InterpolativeDependency vector
if (m_indexmap_cp.count(name)) {
for (const auto& item: m_indexmap_cp.at(name)) {
auto& covdepNode =
speciesNode["coverage-dependencies"][item.first].getMapWhere(
"heat-capacity-a", "", true);
covdepNode["heat-capacity-a"].setQuantity(
m_HeatCapacityDependency[item.second].coeff_a, "J/kmol/K");
covdepNode["heat-capacity-b"].setQuantity(
m_HeatCapacityDependency[item.second].coeff_b, "J/kmol/K");
}
}
}

void CoverageDependentSurfPhase::getGibbs_RT_ref(double* grt) const
{
SurfPhase::_updateThermo();
Expand Down
31 changes: 31 additions & 0 deletions test/thermo/thermoToYaml.cpp
View file
Expand Up @@ -141,6 +141,37 @@ TEST_F(ThermoToYaml, Edge)
EXPECT_DOUBLE_EQ(data["site-density"].asDouble(), 5e-18);
}

TEST_F(ThermoToYaml, CoverageDependentSurface)
{
setup("copt_covdepsurf_example.yaml", "covdep");
EXPECT_EQ(data["thermo"], "ideal-surface");
EXPECT_DOUBLE_EQ(data["site-density"].asDouble(), 2.72e-8);
EXPECT_DOUBLE_EQ(data["reference-state-coverage"].asDouble(), 0.22);
UnitSystem us;
EXPECT_EQ(speciesData[1]["coverage-dependencies"]["OC_Pt"]["model"], "linear");
EXPECT_DOUBLE_EQ(
speciesData[1]["coverage-dependencies"]["OC_Pt"]["enthalpy"]
.asDouble(), us.convertFrom(0.48, "eV/molec"));
EXPECT_EQ(
speciesData[2]["coverage-dependencies"]["OC_Pt"]["model"], "piecewise-linear");
EXPECT_DOUBLE_EQ(
speciesData[2]["coverage-dependencies"]["OC_Pt"]["entropy-high"]
.asDouble(), us.convertFrom(-0.1, "eV/molec"));
EXPECT_EQ(
speciesData[1]["coverage-dependencies"]["O_Pt"]["model"], "interpolative");
EXPECT_DOUBLE_EQ(
speciesData[1]["coverage-dependencies"]["O_Pt"]["enthalpies"]
.asVector<double>()[3], us.convertFrom(2.7, "kcal/mol"));
EXPECT_EQ(
speciesData[4]["coverage-dependencies"]["O_Pt"]["model"], "polynomial");
EXPECT_DOUBLE_EQ(
speciesData[4]["coverage-dependencies"]["O_Pt"]["enthalpy-coefficients"]
.asVector<double>()[3], us.convertFrom(2.11, "eV/molec"));
EXPECT_DOUBLE_EQ(
speciesData[3]["coverage-dependencies"]["C_Pt"]["heat-capacity-a"]
.asDouble(), us.convertFrom(0.07e-3, "eV/molec"));
}

TEST_F(ThermoToYaml, IonsFromNeutral)
{
suppress_deprecation_warnings();
Expand Down

0 comments on commit 4d3c83b

Please sign in to comment.