Skip to content

Commit

Permalink
Fix interpolative map check routine and getSpeciesParameters method
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 20f428c commit 8e9913d
Show file tree
Hide file tree
Showing 2 changed files with 89 additions and 125 deletions.
37 changes: 13 additions & 24 deletions include/cantera/thermo/CoverageDependentSurfPhase.h
View file
Expand Up @@ -139,12 +139,14 @@ class CoverageDependentSurfPhase : public SurfPhase
* coverage-dependent entropy [J/kmol/K] in order of
* 1st-order, 2nd-order, 3rd-order, and 4th-order
* coefficients
* @param isLinear boolean indicating whether the dependency is linear
*/
PolynomialDependency(size_t k, size_t j,
vector_fp enthalpy_coeffs={0.0, 0.0, 0.0, 0.0, 0.0},
vector_fp entropy_coeffs={0.0, 0.0, 0.0, 0.0, 0.0}):
vector_fp entropy_coeffs={0.0, 0.0, 0.0, 0.0, 0.0},
bool isLinear=false):
k(k), j(j), enthalpy_coeffs(enthalpy_coeffs),
entropy_coeffs(entropy_coeffs) {}
entropy_coeffs(entropy_coeffs), isLinear(isLinear) {}
//! index of a target species whose enthalpy and entropy is calculated
size_t k;
//! index of a species whose coverage affects enthalpy and entropy of
Expand All @@ -158,6 +160,8 @@ class CoverageDependentSurfPhase : public SurfPhase
//! [J/kmol/K] in order of 1st-order, 2nd-order, 3rd-order, and 4th-order
//! coefficients
vector_fp entropy_coeffs;
//! boolean indicating whether the dependency is linear
bool isLinear;
};

//! A struct to store sets of parameters used in coverage-dependent enthalpy
Expand All @@ -171,14 +175,17 @@ class CoverageDependentSurfPhase : public SurfPhase
* a target species
* @param enthalpy_map map of <coverage[dimensionless], enthalpy[J/kmol]> pairs
* @param entropy_map map of <coverage[dimensionless], entropy[J/kmol/K]> pairs
* @param isPiecewise boolean indicating whether the dependency is
* piecewise-linear
*/
InterpolativeDependency(size_t k, size_t j,
std::map<double, double> enthalpy_map={{0.0, 0.0},
{1.0, 0.0}},
std::map<double, double> entropy_map={{0.0, 0.0},
{1.0, 0.0}}):
{1.0, 0.0}},
bool isPiecewise=false):
k(k), j(j), enthalpy_map(enthalpy_map),
entropy_map(entropy_map) {}
entropy_map(entropy_map), isPiecewise(isPiecewise) {}
//! index of a target species whose enthalpy and entropy are calculated
size_t k;
//! index of a species whose coverage affects enthalpy and entropy of
Expand All @@ -188,6 +195,8 @@ class CoverageDependentSurfPhase : public SurfPhase
std::map<double, double> enthalpy_map;
//! map of <coverage[dimensionless], entropy[J/kmol/K]> pairs
std::map<double, double> entropy_map;
//! boolean indicating whether the dependency is piecewise-linear
bool isPiecewise;
};

//! A struct to store sets of parameters used in coverage-dependent heat capacity
Expand Down Expand Up @@ -425,26 +434,6 @@ class CoverageDependentSurfPhase : public SurfPhase
//! Array of heat capacity coverage dependency parameters.
std::vector<HeatCapacityDependency> m_HeatCapacityDependency;

//! 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.
Expand Down
177 changes: 76 additions & 101 deletions src/thermo/CoverageDependentSurfPhase.cpp
View file
Expand Up @@ -34,20 +34,20 @@ CoverageDependentSurfPhase::CoverageDependentSurfPhase(const std::string& infile
void CoverageDependentSurfPhase::addInterpolativeDependency(
const InterpolativeDependency& int_deps)
{
if (int_deps.enthalpy_map.count(0.0) == 0) {
if (int_deps.enthalpy_map.begin()->first != 0.0) {
throw CanteraError("CoverageDependentSurfPhase::addInterpolativeDependency",
"The first element of enthalpy-coverages array must be 0.0.");
}
if (int_deps.enthalpy_map.count(1.0) == 0) {
if (int_deps.enthalpy_map.rbegin()->first != 1.0) {
throw CanteraError("CoverageDependentSurfPhase::addInterpolativeDependency",
"The last element of enthalpy-coverages array must be 1.0.");
}

if (int_deps.entropy_map.count(0.0) == 0) {
if (int_deps.entropy_map.begin()->first != 0.0) {
throw CanteraError("CoverageDependentSurfPhase::addInterpolativeDependency",
"The first element of entropy-coverages array must be 0.0.");
}
if (int_deps.entropy_map.count(1.0) == 0) {
if (int_deps.entropy_map.rbegin()->first != 1.0) {
throw CanteraError("CoverageDependentSurfPhase::addInterpolativeDependency",
"The last element of entropy-coverages array must be 1.0.");
}
Expand All @@ -67,8 +67,6 @@ 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 @@ -98,8 +96,7 @@ void CoverageDependentSurfPhase::initThermo()
"J/kmol/K");
}

poly_ind = m_PolynomialDependency.size();
m_indexmap_lin[item.first][item2.first] = poly_ind;
poly_deps.isLinear = true;
m_PolynomialDependency.push_back(poly_deps);
// For polynomial(4th) model
} else if (cov_map2["model"] == "polynomial") {
Expand All @@ -117,8 +114,6 @@ 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
} else if (cov_map2["model"] == "piecewise-linear") {
Expand All @@ -144,8 +139,7 @@ void CoverageDependentSurfPhase::initThermo()
+ int_deps.entropy_map[cov_change];
}

int_ind = m_InterpolativeDependency.size();
m_indexmap_pwlin[item.first][item2.first] = int_ind;
int_deps.isPiecewise = true;
addInterpolativeDependency(int_deps);
// For interpolative model
} else if (cov_map2["model"] == "interpolative") {
Expand Down Expand Up @@ -181,8 +175,6 @@ 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 @@ -196,8 +188,6 @@ 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 Down Expand Up @@ -232,98 +222,83 @@ 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);
size_t k = speciesIndex(name);
// Get linear and polynomial model parameters from PolynomialDependency vector
for (auto& item : m_PolynomialDependency) {
if (item.k == k) {
if (item.isLinear) {
auto& covdepNode =
speciesNode["coverage-dependencies"][speciesName(item.j)]
.getMapWhere("model", "linear", true);
covdepNode["enthalpy"].setQuantity(item.enthalpy_coeffs[1], "J/kmol");
covdepNode["entropy"].setQuantity(item.entropy_coeffs[1], "J/kmol/K");
} else {
auto& covdepNode =
speciesNode["coverage-dependencies"][speciesName(item.j)]
.getMapWhere("model", "polynomial", true);
vector_fp hvec (
item.enthalpy_coeffs.begin() + 1, item.enthalpy_coeffs.end());
covdepNode["enthalpy-coefficients"].setQuantity(hvec, "J/kmol");
vector_fp svec (
item.entropy_coeffs.begin() + 1, item.entropy_coeffs.end());
covdepNode["entropy-coefficients"].setQuantity(svec, "J/kmol/K");
}
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);
// Get piecewise-linear model parameters from InterpolativeDependency vector
for (auto& item : m_InterpolativeDependency) {
if (item.k == k) {
if (item.isPiecewise) {
auto& covdepNode =
speciesNode["coverage-dependencies"][speciesName(item.j)]
.getMapWhere("model", "piecewise-linear", true);
vector_fp hcovs, enthalpies, scovs, entropies;
for (const auto& hmap : item.enthalpy_map) {
hcovs.push_back(hmap.first);
enthalpies.push_back(hmap.second);
}
for (const auto& smap : item.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");
} else {
auto& covdepNode =
speciesNode["coverage-dependencies"][speciesName(item.j)]
.getMapWhere("model", "interpolative", true);
vector_fp hcovs, enthalpies, scovs, entropies;
for (const auto& hmap : item.enthalpy_map) {
hcovs.push_back(hmap.first);
enthalpies.push_back(hmap.second);
}
for (const auto& smap : item.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");
}
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)) {
// Get heat capacity model parameters from HeatCapacityDependency vector
for (auto& item : m_HeatCapacityDependency) {
if (item.k == k) {
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");
speciesNode["coverage-dependencies"][speciesName(item.j)]
.getMapWhere("heat-capacity-a", "", true);
covdepNode["heat-capacity-a"].setQuantity(item.coeff_a, "J/kmol/K");
covdepNode["heat-capacity-b"].setQuantity(item.coeff_b, "J/kmol/K");
}
}
}
Expand Down

0 comments on commit 8e9913d

Please sign in to comment.