research/cantera
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

[Input] Read local 'units' definitions in reaction and thermo entries

Browse source 
Move 'convert' functions that work with missing keys to AnyMap class.
This commit is contained in:
Ray Speth 2018-12-20 14:45:37 -05:00
parent 6e3053dd1c
commit 84fa231ea2
11 changed files with 269 additions and 157 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

74
include/cantera/base/AnyMap.h
View file

@ -5,6 +5,7 @@
#include "cantera/base/ct_defs.h"
#include "cantera/base/global.h"
#include "cantera/base/Units.h"
#include "cantera/base/ctexceptions.h"
#include <string>
@ -102,6 +103,9 @@ public:
template<class T>
std::map<std::string, T> asMap() const;
//! @see AnyMap::applyUnits
void applyUnits(const UnitSystem& units);
private:
std::string demangle(const std::type_info& type) const;
@ -229,6 +233,53 @@ public:
const std::string& getString(const std::string& key,
const std::string& default_) const;
//! Convert the item stored by the given `key` to the units specified in
//! `units`. If the stored value is a double, convert it using the default
//! units. If the input is a string, treat this as a dimensioned value, e.g.
//! '988 kg/m^3' and convert from the specified units.
double convert(const std::string& key, const std::string& units) const;
//! Convert the item stored by the given `key` to the units specified in
//! `units`. If the stored value is a double, convert it using the default
//! units. If the input is a string, treat this as a dimensioned value, e.g.
//! '988 kg/m^3' and convert from the specified units. If the key is
//! missing, the `default_` value is returned.
double convert(const std::string& key, const std::string& units,
double default_) const;
//! Convert a vector of dimensional values
/*!
* For each item in the vector, if the stored value is a double, convert it
* using the default units. If the value is a string, treat it as a
* dimensioned value, e.g. '988 kg/m^3', and convert from the specified
* units.
*
* @param key Location of the vector in this AnyMap
* @param units Units to convert to
* @param nMin Minimum allowed length of the vector. If #nMax is not
* specified, this is also taken to be the maximum length. An exception
* is thrown if this condition is not met.
* @param nMax Maximum allowed length of the vector. An exception is
* thrown if this condition is not met.
*/
vector_fp convertVector(const std::string& key, const std::string& units,
size_t nMin=npos, size_t nMax=npos) const;
//! Convert the item stored by the given `key` to the units specified in
//! `units`. If the stored value is a double, convert it using the default
//! units. If the input is a string, treat this as a dimensioned value, e.g.
//! '2.7e4 J/kmol' and convert from the specified units.
double convertMolarEnergy(const std::string& key,
const std::string& units) const;
//! Convert the item stored by the given `key` to the units specified in
//! `units`. If the stored value is a double, convert it using the default
//! units. If the stored value is a string, treat it as a dimensioned value,
//! e.g. '2.7e4 J/kmol' and convert from the specified units. If the key is
//! missing, the `default_` value is returned.
double convertMolarEnergy(const std::string& key, const std::string& units,
double default_) const;
// Define begin() and end() to allow use with range-based for loops
using const_iterator = std::unordered_map<std::string, AnyValue>::const_iterator;
const_iterator begin() const {
@ -239,12 +290,35 @@ public:
return m_data.end();
}
size_t size() {
return m_data.size();
};
//! Return the default units that should be used to convert stored values
const UnitSystem& units() const { return m_units; }
//! Use the supplied UnitSystem to set the default units, and recursively
//! process overrides from nodes named `units`.
/*!
* If a `units` node is present in a map that contains other keys, the
* specified units are taken to be the defaults for that map. If the map
* contains only a `units` node, and is the first item in a list of maps,
* then the specified units are taken to be the defaults for all the maps in
* the list.
*
* After being processed, the `units` nodes are removed, so this function
* should be called only once, on the root AnyMap. This function is called
* automatically by the fromYamlFile() and fromYamlString() constructors.
*/
void applyUnits(const UnitSystem& units);
private:
template <class T>
const T& get(const std::string& key, const T& default_,
std::function<const T&(const AnyValue*)> getter) const;
std::unordered_map<std::string, AnyValue> m_data;
UnitSystem m_units;
friend class AnyValue;
};

17
include/cantera/base/Units.h
View file

@ -146,16 +146,6 @@ public:
double convert(const AnyValue& val, const std::string& dest) const;
double convert(const AnyValue& val, const Units& dest) const;
//! Convert the value at `node[key]` to the units specified in `dest`. If
//! the input is a double, convert it using the default units. If the input
//! is a string, treat this as a dimensioned value, e.g. '988 kg/m^3' and
//! convert from the specified units. If the key is missing, the `default_`
//! value is returned.
double convert(const AnyMap& node, const std::string key,
const std::string& dest, double default_) const;
double convert(const AnyMap& node, const std::string key,
const Units& dest, double default_) const;
//! Convert an array of AnyValue nodes to the units specified in `dest`. For
//! each node, if the value is a double, convert it using the default units,
//! and if it is a string, treat it as a value with the given dimensions.
@ -179,13 +169,6 @@ public:
//! convert from the specified units.
double convertMolarEnergy(const AnyValue& val, const std::string& dest) const;
//! Convert the value at `node[key]` to the molar energy units specified in
//! `dest`. If the input is a double, convert it using the default units. If
//! the input is a string, treat this as a dimensioned value, e.g. '988
//! cal/mol' and convert from the specified units. If the key is missing,
//! the `default_` value is returned.
double convertMolarEnergy(const AnyMap& node, const std::string& key,
const std::string& dest, double default_) const;
private:
//! Factor to convert mass from this unit system to kg

4
include/cantera/kinetics/Reaction.h
View file

@ -17,7 +17,6 @@ namespace Cantera
class Kinetics;
class AnyMap;
class UnitSystem;
//! Intermediate class which stores data about a reaction and its rate
//! parameterization so that it can be added to a Kinetics object.
@ -260,8 +259,7 @@ public:
shared_ptr<Reaction> newReaction(const XML_Node& rxn_node);
//! Create a new Reaction object using the specified parameters
unique_ptr<Reaction> newReaction(const AnyMap& rxn_node, const Kinetics& kin,
const UnitSystem& units);
unique_ptr<Reaction> newReaction(const AnyMap& rxn_node, const Kinetics& kin);
//! Create Reaction objects for all `<reaction>` nodes in an XML document.
//!

5
include/cantera/thermo/SpeciesThermoFactory.h
View file

@ -18,7 +18,6 @@ namespace Cantera
class XML_Node;
class AnyMap;
class UnitSystem;
//! Create a new SpeciesThermoInterpType object given a corresponding constant.
/*!
@ -56,10 +55,8 @@ SpeciesThermoInterpType* newSpeciesThermoInterpType(const XML_Node& thermoNode);
/*!
* @param thermo_node An AnyMap specifying the model type (e.g. "NASA") and any
* model parameters necessary to instantiate the object
* @param units Specification for the unit system to convert from
*/
unique_ptr<SpeciesThermoInterpType> newSpeciesThermo(
const AnyMap& thermo_node, const UnitSystem& units);
unique_ptr<SpeciesThermoInterpType> newSpeciesThermo(const AnyMap& thermo_node);
}

91
src/base/AnyMap.cpp
View file

@ -361,6 +361,39 @@ AnyValue& AnyValue::operator=(AnyMap&& value) {
return *this;
}
void AnyValue::applyUnits(const UnitSystem& units)
{
if (is<AnyMap>()) {
// Units declaration applicable to this map
as<AnyMap>().applyUnits(units);
} else if (is<std::vector<AnyMap>>()) {
auto& list = as<std::vector<AnyMap>>();
if (list.size() && list[0].hasKey("units") && list[0].size() == 1) {
// First item in the list is a units declaration, which applies to
// the items in the list
UnitSystem newUnits = units;
newUnits.setDefaults(list[0].at("units").asMap<std::string>());
list[0].m_data.erase("units");
for (auto& item : list) {
// Any additional units declarations are errors
if (item.size() == 1 && item.hasKey("units")) {
throw CanteraError("AnyValue::applyUnits",
"Found units entry as not the first item in a list.");
}
item.applyUnits(newUnits);
}
// Remove the "units" map after it has been applied
list.erase(list.begin());
} else {
// Simple downward propagation of the current units
for (auto& item : list) {
item.applyUnits(units);
}
}
}
}
std::string AnyValue::demangle(const std::type_info& type) const
{
if (s_typenames.find(type.name()) != s_typenames.end()) {
@ -596,14 +629,68 @@ const std::string& AnyMap::getString(const std::string& key,
return get<std::string>(key, default_, &AnyValue::asString);
}
double AnyMap::convert(const std::string& key, const std::string& dest) const
{
return units().convert(at(key), dest);
}
double AnyMap::convert(const std::string& key, const std::string& dest,
double default_) const
{
if (hasKey(key)) {
return units().convert(at(key), dest);
} else {
return default_;
}
}
vector_fp AnyMap::convertVector(const std::string& key, const std::string& dest,
size_t nMin, size_t nMax) const
{
return units().convert(at(key).asVector<AnyValue>(nMin, nMax), dest);
}
double AnyMap::convertMolarEnergy(const std::string& key,
const std::string& dest) const
{
return units().convertMolarEnergy(at(key), dest);
}
double AnyMap::convertMolarEnergy(const std::string& key,
const std::string& dest,
double default_) const
{
if (hasKey(key)) {
return units().convertMolarEnergy(at(key), dest);
} else {
return default_;
}
}
void AnyMap::applyUnits(const UnitSystem& units) {
m_units = units;
if (hasKey("units")) {
m_units.setDefaults(at("units").asMap<std::string>());
m_data.erase("units");
}
for (auto& item : m_data) {
item.second.applyUnits(m_units);
}
}
AnyMap AnyMap::fromYamlString(const std::string& yaml) {
YAML::Node node = YAML::Load(yaml);
return node.as<AnyMap>();
AnyMap amap = node.as<AnyMap>();
amap.applyUnits(UnitSystem());
return amap;
}
AnyMap AnyMap::fromYamlFile(const std::string& name) {
YAML::Node node = YAML::LoadFile(findInputFile(name));
return node.as<AnyMap>();
AnyMap amap = node.as<AnyMap>();
amap.applyUnits(UnitSystem());
return amap;
}
AnyMap::const_iterator begin(const AnyValue& v) {

33
src/base/Units.cpp
View file

@ -212,7 +212,6 @@ std::string Units::str() const {
m_temperature_dim, m_current_dim, m_quantity_dim);
}
UnitSystem::UnitSystem(std::initializer_list<std::string> units)
: m_mass_factor(1.0)
, m_length_factor(1.0)
@ -358,27 +357,6 @@ double UnitSystem::convert(const AnyValue& v, const Units& dest) const
}
}
double UnitSystem::convert(const AnyMap& node, const std::string key,
const std::string& dest, double default_) const
{
if (node.hasKey(key)) {
return convert(node.at(key), Units(dest));
} else {
return default_;
}
}
double UnitSystem::convert(const AnyMap& node, const std::string key,
const Units& dest, double default_) const
{
if (node.hasKey(key)) {
return convert(node.at(key), dest);
} else {
return default_;
}
}
vector_fp UnitSystem::convert(const std::vector<AnyValue>& vals,
const std::string& dest) const
{
@ -455,15 +433,4 @@ double UnitSystem::convertMolarEnergy(const AnyValue& v,
}
}
double UnitSystem::convertMolarEnergy(
const AnyMap& node, const std::string& key,
const std::string& dest, double default_) const
{
if (node.hasKey(key)) {
return convertMolarEnergy(node.at(key), dest);
} else {
return default_;
}
}
}

48
src/kinetics/Reaction.cpp
View file

@ -11,7 +11,6 @@
#include "cantera/base/ctml.h"
#include "cantera/base/Array.h"
#include "cantera/base/AnyMap.h"
#include "cantera/base/Units.h"
#include <sstream>
#include <boost/algorithm/string.hpp>
@ -528,11 +527,11 @@ void setupElementaryReaction(ElementaryReaction& R, const XML_Node& rxn_node)
}
void setupElementaryReaction(ElementaryReaction& R, const AnyMap& node,
const Kinetics& kin, const UnitSystem& units)
const Kinetics& kin)
{
setupReaction(R, node);
R.allow_negative_pre_exponential_factor = node.getBool("negative-A", false);
R.rate = readArrhenius(R, node.at("rate-constant"), kin, units);
R.rate = readArrhenius(R, node.at("rate-constant"), kin, node.units());
}
void setupThreeBodyReaction(ThreeBodyReaction& R, const XML_Node& rxn_node)
@ -542,9 +541,9 @@ void setupThreeBodyReaction(ThreeBodyReaction& R, const XML_Node& rxn_node)
}
void setupThreeBodyReaction(ThreeBodyReaction& R, const AnyMap& node,
const Kinetics& kin, const UnitSystem& units)
const Kinetics& kin)
{
setupElementaryReaction(R, node, kin, units);
setupElementaryReaction(R, node, kin);
if (R.reactants.count("M") != 1 || R.products.count("M") != 1) {
throw CanteraError("setupThreeBodyReaction",
"Reaction equation '{}' does not contain third body 'M'",
@ -584,7 +583,7 @@ void setupFalloffReaction(FalloffReaction& R, const XML_Node& rxn_node)
}
void setupFalloffReaction(FalloffReaction& R, const AnyMap& node,
const Kinetics& kin, const UnitSystem& units)
const Kinetics& kin)
{
setupReaction(R, node);
// setupReaction sets the stoichiometric coefficient for the falloff third
@ -621,11 +620,15 @@ void setupFalloffReaction(FalloffReaction& R, const AnyMap& node,
}
if (node.at("type").asString() == "falloff") {
R.low_rate = readArrhenius(R, node.at("low-P-rate-constant"), kin, units, 1);
R.high_rate = readArrhenius(R, node.at("high-P-rate-constant"), kin, units);
R.low_rate = readArrhenius(R, node.at("low-P-rate-constant"), kin,
node.units(), 1);
R.high_rate = readArrhenius(R, node.at("high-P-rate-constant"), kin,
node.units());
} else { // type == "chemically-activated"
R.low_rate = readArrhenius(R, node.at("low-P-rate-constant"), kin, units);
R.high_rate = readArrhenius(R, node.at("high-P-rate-constant"), kin, units, -1);
R.low_rate = readArrhenius(R, node.at("low-P-rate-constant"), kin,
node.units());
R.high_rate = readArrhenius(R, node.at("high-P-rate-constant"), kin,
node.units(), -1);
}
readFalloff(R, node);
@ -672,15 +675,14 @@ void setupPlogReaction(PlogReaction& R, const XML_Node& rxn_node)
setupReaction(R, rxn_node);
}
void setupPlogReaction(PlogReaction& R, const AnyMap& node,
const Kinetics& kin, const UnitSystem& units)
void setupPlogReaction(PlogReaction& R, const AnyMap& node, const Kinetics& kin)
{
setupReaction(R, node);
std::multimap<double, Arrhenius> rates;
for (const auto& rate : node.at("rate-constants").asVector<AnyValue>()) {
const auto& p_rate = rate.asVector<AnyValue>(2);
rates.insert({units.convert(p_rate[0], "Pa"),
readArrhenius(R, p_rate[1], kin, units)});
rates.insert({node.units().convert(p_rate[0], "Pa"),
readArrhenius(R, p_rate[1], kin, node.units())});
}
R.rate = Plog(rates);
}
@ -715,7 +717,7 @@ void setupChebyshevReaction(ChebyshevReaction& R, const XML_Node& rxn_node)
}
void setupChebyshevReaction(ChebyshevReaction&R, const AnyMap& node,
const Kinetics& kin, const UnitSystem& units)
const Kinetics& kin)
{
setupReaction(R, node);
R.reactants.erase("(+M)"); // remove optional third body notation
@ -729,6 +731,7 @@ void setupChebyshevReaction(ChebyshevReaction&R, const AnyMap& node,
coeffs(i, j) = vcoeffs[i][j];
}
}
const UnitSystem& units = node.units();
Units rcUnits = rateCoeffUnits(R, kin);
coeffs(0, 0) += std::log10(units.convert(1.0, rcUnits));
R.rate = ChebyshevRate(units.convert(T_range[0], "K"),
@ -916,8 +919,7 @@ shared_ptr<Reaction> newReaction(const XML_Node& rxn_node)
}
}
unique_ptr<Reaction> newReaction(const AnyMap& node, const Kinetics& kin,
const UnitSystem& units)
unique_ptr<Reaction> newReaction(const AnyMap& node, const Kinetics& kin)
{
std::string type = "elementary";
if (node.hasKey("type")) {
@ -926,27 +928,27 @@ unique_ptr<Reaction> newReaction(const AnyMap& node, const Kinetics& kin,
if (type == "elementary") {
unique_ptr<ElementaryReaction> R(new ElementaryReaction());
setupElementaryReaction(*R, node, kin, units);
setupElementaryReaction(*R, node, kin);
return unique_ptr<Reaction>(move(R));
} else if (type == "three-body") {
unique_ptr<ThreeBodyReaction> R(new ThreeBodyReaction());
setupThreeBodyReaction(*R, node, kin, units);
setupThreeBodyReaction(*R, node, kin);
return unique_ptr<Reaction>(move(R));
} else if (type == "falloff") {
unique_ptr<FalloffReaction> R(new FalloffReaction());
setupFalloffReaction(*R, node, kin, units);
setupFalloffReaction(*R, node, kin);
return unique_ptr<Reaction>(move(R));
} else if (type == "chemically-activated") {
unique_ptr<ChemicallyActivatedReaction> R(new ChemicallyActivatedReaction());
setupFalloffReaction(*R, node, kin, units);
setupFalloffReaction(*R, node, kin);
return unique_ptr<Reaction>(move(R));
} else if (type == "pressure-dependent-Arrhenius") {
unique_ptr<PlogReaction> R(new PlogReaction());
setupPlogReaction(*R, node, kin, units);
setupPlogReaction(*R, node, kin);
return unique_ptr<Reaction>(move(R));
} else if (type == "Chebyshev") {
unique_ptr<ChebyshevReaction> R(new ChebyshevReaction());
setupChebyshevReaction(*R, node, kin, units);
setupChebyshevReaction(*R, node, kin);
return unique_ptr<Reaction>(move(R));
} else {
throw CanteraError("newReaction", "Unknown reaction type '{}'", type);

62
src/thermo/SpeciesThermoFactory.cpp
View file

@ -143,18 +143,16 @@ static SpeciesThermoInterpType* newNasaThermoFromXML(vector<XML_Node*> nodes)
}
void setupSpeciesThermo(SpeciesThermoInterpType& thermo,
const AnyMap& node, const UnitSystem& units)
const AnyMap& node)
{
double Pref = units.convert(node, "reference-pressure", "Pa", OneAtm);
double Pref = node.convert("reference-pressure", "Pa", OneAtm);
thermo.setRefPressure(Pref);
}
void setupNasaPoly(NasaPoly2& thermo, const AnyMap& node,
const UnitSystem& units)
void setupNasaPoly(NasaPoly2& thermo, const AnyMap& node)
{
setupSpeciesThermo(thermo, node, units);
vector_fp Tranges = units.convert(
node.at("temperature-ranges").asVector<AnyValue>(2, 3), "K");
setupSpeciesThermo(thermo, node);
vector_fp Tranges = node.convertVector("temperature-ranges", "K", 2, 3);
const auto& data = node.at("data").asVector<vector_fp>(Tranges.size()-1);
for (const auto& poly : data) {
if (poly.size() != 7) {
@ -247,12 +245,10 @@ static SpeciesThermoInterpType* newShomateThermoFromXML(
}
void setupShomatePoly(ShomatePoly2& thermo, const AnyMap& node,
const UnitSystem& units)
void setupShomatePoly(ShomatePoly2& thermo, const AnyMap& node)
{
setupSpeciesThermo(thermo, node, units);
vector_fp Tranges = units.convert(
node.at("temperature-ranges").asVector<AnyValue>(2, 3), "K");
setupSpeciesThermo(thermo, node);
vector_fp Tranges = node.convertVector("temperature-ranges", "K", 2, 3);
const auto& data = node.at("data").asVector<vector_fp>(Tranges.size()-1);
for (const auto& poly : data) {
if (poly.size() != 7) {
@ -294,13 +290,12 @@ static SpeciesThermoInterpType* newConstCpThermoFromXML(XML_Node& f)
return newSpeciesThermoInterpType(CONSTANT_CP, tmin, tmax, p0, &c[0]);
}
void setupConstCp(ConstCpPoly& thermo, const AnyMap& node,
const UnitSystem& units)
void setupConstCp(ConstCpPoly& thermo, const AnyMap& node)
{
double T0 = units.convert(node.at("T0"), "K");
double h0 = units.convert(node, "h0", "J/kmol", 0.0);
double s0 = units.convert(node, "s0", "J/kmol/K", 0.0);
double cp0 = units.convert(node, "cp0", "J/kmol/K", 0.0);
double T0 = node.convert("T0", "K");
double h0 = node.convert("h0", "J/kmol", 0.0);
double s0 = node.convert("s0", "J/kmol/K", 0.0);
double cp0 = node.convert("cp0", "J/kmol/K", 0.0);
thermo.setParameters(T0, h0, s0, cp0);
}
@ -350,12 +345,10 @@ static SpeciesThermoInterpType* newNasa9ThermoFromXML(
}
void setupNasa9Poly(Nasa9PolyMultiTempRegion& thermo, const AnyMap& node,
const UnitSystem& units)
void setupNasa9Poly(Nasa9PolyMultiTempRegion& thermo, const AnyMap& node)
{
setupSpeciesThermo(thermo, node, units);
vector_fp Tranges = units.convert(
node.at("temperature-ranges").asVector<AnyValue>(2, 999), "K");
setupSpeciesThermo(thermo, node);
vector_fp Tranges = node.convertVector("temperature-ranges", "K", 2, 999);
const auto& data = node.at("data").asVector<vector_fp>(Tranges.size()-1);
map<double, vector_fp> regions;
for (size_t i = 0; i < data.size(); i++) {
@ -371,18 +364,18 @@ void setupNasa9Poly(Nasa9PolyMultiTempRegion& thermo, const AnyMap& node,
}
void setupMu0(Mu0Poly& thermo, const AnyMap& node, const UnitSystem& units)
void setupMu0(Mu0Poly& thermo, const AnyMap& node)
{
setupSpeciesThermo(thermo, node, units);
setupSpeciesThermo(thermo, node);
bool dimensionless = node.getBool("dimensionless", false);
double h0 = units.convertMolarEnergy(node, "h0", "J/kmol", 0.0);
double h0 = node.convertMolarEnergy("h0", "J/kmol", 0.0);
map<double, double> T_mu;
for (const auto& item : node.at("data")) {
double T = units.convert(fpValueCheck(item.first), "K");
double T = node.units().convert(fpValueCheck(item.first), "K");
if (dimensionless) {
T_mu[T] = item.second.asDouble() * GasConstant * T;
} else {
T_mu[T] = units.convertMolarEnergy(item.second, "J/kmol");
T_mu[T] = node.units().convertMolarEnergy(item.second, "J/kmol");
}
}
thermo.setParameters(h0, T_mu);
@ -443,29 +436,28 @@ SpeciesThermoInterpType* newSpeciesThermoInterpType(const XML_Node& thermo)
}
unique_ptr<SpeciesThermoInterpType> newSpeciesThermo(
const AnyMap& node, const UnitSystem& units)
unique_ptr<SpeciesThermoInterpType> newSpeciesThermo(const AnyMap& node)
{
std::string model = node.at("model").asString();
if (model == "NASA7") {
unique_ptr<NasaPoly2> thermo(new NasaPoly2());
setupNasaPoly(*thermo, node, units);
setupNasaPoly(*thermo, node);
return unique_ptr<SpeciesThermoInterpType>(move(thermo));
} else if (model == "Shomate") {
unique_ptr<ShomatePoly2> thermo(new ShomatePoly2());
setupShomatePoly(*thermo, node, units);
setupShomatePoly(*thermo, node);
return unique_ptr<SpeciesThermoInterpType>(move(thermo));
} else if (model == "NASA9") {
unique_ptr<Nasa9PolyMultiTempRegion> thermo(new Nasa9PolyMultiTempRegion());
setupNasa9Poly(*thermo, node, units);
setupNasa9Poly(*thermo, node);
return unique_ptr<SpeciesThermoInterpType>(move(thermo));
} else if (model == "constant-cp") {
unique_ptr<ConstCpPoly> thermo(new ConstCpPoly());
setupConstCp(*thermo, node, units);
setupConstCp(*thermo, node);
return unique_ptr<SpeciesThermoInterpType>(move(thermo));
} else if (model == "piecewise-Gibbs") {
unique_ptr<Mu0Poly> thermo(new Mu0Poly());
setupMu0(*thermo, node, units);
setupMu0(*thermo, node);
return unique_ptr<SpeciesThermoInterpType>(move(thermo));
} else {
throw CanteraError("newSpeciesThermo",

42
test/general/test_units.cpp
View file

@ -100,10 +100,11 @@ TEST(Units, from_anymap) {
"{p: 12 bar, v: 10, A: 1 cm^2, V: 1,"
" k1: [5e2, 2, 29000], k2: [1e14, -1, 1300 cal/kmol]}");
UnitSystem U({"mm", "min", "atm"});
EXPECT_DOUBLE_EQ(U.convert(m["p"], "Pa"), 12e5);
EXPECT_DOUBLE_EQ(U.convert(m["v"], "cm/min"), 1.0);
EXPECT_DOUBLE_EQ(U.convert(m["A"], "mm^2"), 100);
EXPECT_DOUBLE_EQ(U.convert(m["V"], "m^3"), 1e-9);
m.applyUnits(U);
EXPECT_DOUBLE_EQ(m.convert("p", "Pa"), 12e5);
EXPECT_DOUBLE_EQ(m.convert("v", "cm/min"), 1.0);
EXPECT_DOUBLE_EQ(m.convert("A", "mm^2"), 100);
EXPECT_DOUBLE_EQ(m.convert("V", "m^3"), 1e-9);
auto k1 = m["k1"].asVector<AnyValue>();
EXPECT_DOUBLE_EQ(U.convert(k1[0], "m^3/kmol"), 1e-9*5e2);
EXPECT_DOUBLE_EQ(U.convertMolarEnergy(k1[2], "J/kmol"), 29000);
@ -111,9 +112,32 @@ TEST(Units, from_anymap) {
TEST(Units, from_anymap_default) {
AnyMap m = AnyMap::fromYamlString("{p0: 10 atm, h0: 10 cal/kmol}");
UnitSystem U;
EXPECT_DOUBLE_EQ(U.convert(m, "p0", "Pa", 999), 10*OneAtm);
EXPECT_DOUBLE_EQ(U.convert(m, "p1", "Pa", 999), 999);
EXPECT_DOUBLE_EQ(U.convert(m, "h0", "J/kmol", 999), 41.84);
EXPECT_DOUBLE_EQ(U.convert(m, "h1", "J/kmol", 999), 999);
EXPECT_DOUBLE_EQ(m.convert("p0", "Pa", 999), 10*OneAtm);
EXPECT_DOUBLE_EQ(m.convert("p1", "Pa", 999), 999);
EXPECT_DOUBLE_EQ(m.convert("h0", "J/kmol", 999), 41.84);
EXPECT_DOUBLE_EQ(m.convert("h1", "J/kmol", 999), 999);
}
TEST(Units, from_yaml) {
AnyMap m = AnyMap::fromYamlString(
"units: {length: km}\n"
"foo:\n"
"- units: {length: cm}\n" // applies to items in foo
"- bar: 0.6\n"
"- baz: 0.2\n"
" units: {length: mm}\n" // applies to just this entry (with "baz")
"spam:\n"
"- eggs: 3\n"
"- ham: [0.1, 0.3, 0.5]\n"
);
EXPECT_FALSE(m.hasKey("units"));
EXPECT_DOUBLE_EQ(m.units().convert(1, "m"), 1000);
auto& foo = m["foo"].asVector<AnyMap>();
EXPECT_DOUBLE_EQ(foo[0].units().convert(1, "m"), 0.01);
EXPECT_DOUBLE_EQ(foo[1].units().convert(1, "m"), 0.001);
EXPECT_DOUBLE_EQ(foo[0].convert("bar", "m"), 0.006);
auto& spam = m["spam"].asVector<AnyMap>();
EXPECT_DOUBLE_EQ(spam[0].convert("eggs", "m"), 3000);
EXPECT_DOUBLE_EQ(spam[1].convertVector("ham", "m")[2], 500);
}

27
test/kinetics/kineticsFromYaml.cpp
View file

@ -14,8 +14,7 @@ TEST(Reaction, ElementaryFromYaml)
" rate-constant: [-2.70000E+13 cm^3/mol/s, 0, 355 cal/mol],"
" negative-A: true}");
UnitSystem U;
auto R = newReaction(rxn, gas, U);
auto R = newReaction(rxn, gas);
EXPECT_EQ(R->reactants.at("NO"), 1);
EXPECT_EQ(R->products.at("N2"), 1);
EXPECT_EQ(R->reaction_type, ELEMENTARY_RXN);
@ -36,8 +35,7 @@ TEST(Reaction, ThreeBodyFromYaml1)
" rate-constant: [1.20000E+17 cm^6/mol^2/s, -1, 0],"
" efficiencies: {AR: 0.83, H2O: 5}}");
UnitSystem U;
auto R = newReaction(rxn, gas, U);
auto R = newReaction(rxn, gas);
EXPECT_EQ(R->reactants.count("M"), (size_t) 0);
auto TBR = dynamic_cast<ThreeBodyReaction&>(*R);
@ -54,8 +52,7 @@ TEST(Reaction, ThreeBodyFromYaml2)
" type: three-body,"
" rate-constant: [1.20000E+17, -1, 0]}");
UnitSystem U;
EXPECT_THROW(newReaction(rxn, gas, U), CanteraError);
EXPECT_THROW(newReaction(rxn, gas), CanteraError);
}
TEST(Reaction, FalloffFromYaml1)
@ -69,8 +66,7 @@ TEST(Reaction, FalloffFromYaml1)
" SRI: {A: 1.1, B: 700.0, C: 1234.0, D: 56.0, E: 0.7},"
" efficiencies: {AR: 0.625}}");
UnitSystem U;
auto R = newReaction(rxn, gas, U);
auto R = newReaction(rxn, gas);
auto FR = dynamic_cast<FalloffReaction&>(*R);
EXPECT_DOUBLE_EQ(FR.third_body.efficiency("AR"), 0.625);
EXPECT_DOUBLE_EQ(FR.third_body.efficiency("N2"), 1.0);
@ -86,8 +82,7 @@ TEST(Reaction, FalloffFromYaml2)
" low-P-rate-constant: [1.04000E+26 cm^6/mol^2/s, -2.76, 1600],"
" Troe: {A: 0.562, T3: 91, T1: 5836}}");
UnitSystem U;
auto R = newReaction(rxn, gas, U);
auto R = newReaction(rxn, gas);
auto FR = dynamic_cast<FalloffReaction&>(*R);
EXPECT_DOUBLE_EQ(FR.third_body.efficiency("N2"), 1.0);
EXPECT_DOUBLE_EQ(FR.third_body.efficiency("H2O"), 0.0);
@ -105,13 +100,12 @@ TEST(Reaction, ChemicallyActivatedFromYaml)
IdealGasMix gas("gri30.xml");
AnyMap rxn = AnyMap::fromYamlString(
"{equation: CH3 + OH (+M) <=> CH2O + H2 (+M),"
" units: {length: cm, quantity: mol},"
" type: chemically-activated,"
" high-P-rate-constant: [5.88E-14, 6.721, -3022.227],"
" low-P-rate-constant: [282320.078, 1.46878, -3270.56495]}");
UnitSystem U;
U.setDefaults({"cm", "mol"});
auto R = newReaction(rxn, gas, U);
auto R = newReaction(rxn, gas);
auto CAR = dynamic_cast<ChemicallyActivatedReaction&>(*R);
EXPECT_DOUBLE_EQ(CAR.high_rate.preExponentialFactor(), 5.88e-14);
EXPECT_DOUBLE_EQ(CAR.low_rate.preExponentialFactor(), 2.82320078e2);
@ -123,6 +117,7 @@ TEST(Reaction, PlogFromYaml)
IdealGasMix gas("gri30.xml");
AnyMap rxn = AnyMap::fromYamlString(
"equation: 'H + CH4 <=> H2 + CH3'\n"
"units: {pressure: atm}\n"
"type: pressure-dependent-Arrhenius\n"
"rate-constants:\n"
"- [0.039474, [2.720000e+09 cm^3/mol/s, 1.2, 6834.0]]\n"
@ -130,8 +125,7 @@ TEST(Reaction, PlogFromYaml)
"- [1.0 atm, [1.230000e+04, 2.68, 6335.0]]\n"
"- [1.01325 MPa, [1.680000e+16, -0.6, 14754.0]]");
UnitSystem U({"atm"});
auto R = newReaction(rxn, gas, U);
auto R = newReaction(rxn, gas);
auto PR = dynamic_cast<PlogReaction&>(*R);
const auto& rates = PR.rate.rates();
EXPECT_EQ(rates.size(), (size_t) 4);
@ -157,8 +151,7 @@ TEST(Reaction, ChebyshevFromYaml)
" [-2.26210e-01, 1.69190e-01, 4.85810e-03, -2.38030e-03],\n"
" [-1.43220e-01, 7.71110e-02, 1.27080e-02, -6.41540e-04]]\n");
UnitSystem U;
auto R = newReaction(rxn, gas, U);
auto R = newReaction(rxn, gas);
EXPECT_EQ(R->reactants.size(), (size_t) 1);
auto CR = dynamic_cast<ChebyshevReaction&>(*R);
double logP = std::log10(2e6);

23
test/thermo/thermoParameterizations.cpp
View file

@ -7,7 +7,6 @@
#include "cantera/thermo/ShomatePoly.h"
#include "cantera/thermo/PDSS_HKFT.h"
#include "cantera/base/stringUtils.h"
#include "cantera/base/Units.h"
#include "thermo_data.h"
#include <sstream>
@ -134,9 +133,8 @@ TEST(SpeciesThermo, NasaPoly2FromYaml1) {
" 7.835056400E-12, 2.896617900E+03, 6.311991700E+00]\n"
"- [4.884754200E+00, 2.172395600E-03, -8.280690600E-07, 1.574751000E-10,\n"
" -1.051089500E-14, 2.316498300E+03, -1.174169500E-01]\n");
UnitSystem U;
double cp_R, h_RT, s_R;
auto st = newSpeciesThermo(data, U);
auto st = newSpeciesThermo(data);
st->validate("NO2");
st->updatePropertiesTemp(300, &cp_R, &h_RT, &s_R);
EXPECT_DOUBLE_EQ(st->refPressure(), OneAtm);
@ -148,15 +146,16 @@ TEST(SpeciesThermo, NasaPoly2FromYaml1) {
TEST(SpeciesThermo, NasaPoly2FromYaml2) {
AnyMap data = AnyMap::fromYamlString(
"model: NASA7\n"
"reference-pressure: 1 atm\n"
"units: {pressure: atm}\n"
"reference-pressure: 1\n"
"temperature-ranges: [200 K, 1000 K]\n"
"data:\n"
"- [3.944031200E+00, -1.585429000E-03, 1.665781200E-05, -2.047542600E-08,\n"
" 7.835056400E-12, 2.896617900E+03, 6.311991700E+00]\n");
UnitSystem U;
double cp_R, h_RT, s_R;
auto st = newSpeciesThermo(data, U);
auto st = newSpeciesThermo(data);
st->validate("NO2");
EXPECT_DOUBLE_EQ(st->refPressure(), OneAtm);
st->updatePropertiesTemp(300, &cp_R, &h_RT, &s_R);
EXPECT_DOUBLE_EQ(st->maxTemp(), 1000);
EXPECT_DOUBLE_EQ(cp_R, 4.47823303484);
@ -171,9 +170,8 @@ TEST(SpeciesThermo, Shomate2FromYaml1) {
"data:\n"
"- [25.56759, 6.096130, 4.054656, -2.671301, 0.131021, -118.0089, 227.3665]\n"
"- [35.15070, 1.300095, -0.205921, 0.013550, -3.282780, -127.8375, 231.7120]\n");
UnitSystem U;
double cp_R, h_RT, s_R;
auto st = newSpeciesThermo(data, U);
auto st = newSpeciesThermo(data);
st->validate("CO");
st->updatePropertiesTemp(1500, &cp_R, &h_RT, &s_R);
EXPECT_DOUBLE_EQ(st->refPressure(), OneAtm);
@ -197,9 +195,8 @@ TEST(SpeciesThermo, Nasa9PolyFromYaml) {
"- [8.310139160E+08, -6.420733540E+05, 2.020264635E+02, -3.065092046E-02,\n"
" 2.486903333E-06, -9.705954110E-11, 1.437538881E-15, 4.938707040E+06,\n"
" -1.672099740E+03]");
UnitSystem U;
double cp_R, h_RT, s_R;
auto st = newSpeciesThermo(data, U);
auto st = newSpeciesThermo(data);
EXPECT_DOUBLE_EQ(st->refPressure(), 1e5);
st->updatePropertiesTemp(2000, &cp_R, &h_RT, &s_R);
EXPECT_DOUBLE_EQ(cp_R, 4.326181187976);
@ -214,9 +211,8 @@ TEST(SpeciesThermo, ConstCpPolyFromYaml) {
"h0: 9.22 kcal/mol\n"
"s0: -3.02 cal/mol/K\n"
"cp0: 5.95 cal/mol/K\n");
UnitSystem U;
double cp_R, h_RT, s_R;
auto st = newSpeciesThermo(data, U);
auto st = newSpeciesThermo(data);
st->updatePropertiesTemp(1100, &cp_R, &h_RT, &s_R);
EXPECT_DOUBLE_EQ(cp_R * GasConst_cal_mol_K, 5.95);
EXPECT_DOUBLE_EQ(h_RT * GasConst_cal_mol_K * 1100, 9.22e3 + 100 * 5.95);
@ -229,8 +225,7 @@ TEST(SpeciesThermo, Mu0PolyFromYaml) {
" h0: -890 kJ/mol,"
" dimensionless: true,"
" data: {298.15: -363.2104, 323.15: -300}}");
UnitSystem U;
auto st = newSpeciesThermo(data, U);
auto st = newSpeciesThermo(data);
double cp_R, h_RT, s_R;
st->updatePropertiesTemp(310, &cp_R, &h_RT, &s_R);
EXPECT_DOUBLE_EQ(cp_R, -11226.315195362145);