Cantera  3.1.0a1
Loading...
Searching...
No Matches
PlogRate.cpp
Go to the documentation of this file.
1//! @file PlogRate.cpp
2
3// This file is part of Cantera. See License.txt in the top-level directory or
4// at https://cantera.org/license.txt for license and copyright information.
5
8
9namespace Cantera
10{
11
12void PlogData::update(double T)
13{
14 throw CanteraError("PlogData::update",
15 "Missing state information: 'PlogData' requires pressure.");
16}
17
18bool PlogData::update(const ThermoPhase& phase, const Kinetics& kin)
19{
20 double T = phase.temperature();
21 double P = phase.pressure();
22 if (P != pressure || T != temperature) {
23 update(T, P);
24 return true;
25 }
26 return false;
27}
28
29void PlogData::perturbPressure(double deltaP)
30{
31 if (m_pressure_buf > 0.) {
32 throw CanteraError("PlogData::perturbPressure",
33 "Cannot apply another perturbation as state is already perturbed.");
34 }
36 update(temperature, pressure * (1. + deltaP));
37}
38
40{
42 // only restore if there is a valid buffered value
43 if (m_pressure_buf < 0.) {
44 return;
45 }
47 m_pressure_buf = -1.;
48}
49
50// Methods of class PlogRate
51
52PlogRate::PlogRate(const std::multimap<double, ArrheniusRate>& rates)
53{
54 setRates(rates);
55}
56
57PlogRate::PlogRate(const AnyMap& node, const UnitStack& rate_units)
58{
59 setParameters(node, rate_units);
60}
61
62void PlogRate::setParameters(const AnyMap& node, const UnitStack& rate_units)
63{
64 ReactionRate::setParameters(node, rate_units);
65 std::multimap<double, ArrheniusRate> multi_rates;
66 if (node.hasKey("rate-constants")) {
67 auto& rates = node["rate-constants"].asVector<AnyMap>();
68 for (const auto& rate : rates) {
69 multi_rates.insert({rate.convert("P", "Pa"),
70 ArrheniusRate(AnyValue(rate), node.units(), rate_units)});
71 }
72 }
73 setRates(multi_rates);
74}
75
76void PlogRate::getParameters(AnyMap& rateNode, const Units& rate_units) const
77{
78 vector<AnyMap> rateList;
79 if (!valid()) {
80 // object not fully set up
81 return;
82 }
83 for (const auto& [pressure, rate] : getRates()) {
84 AnyMap rateNode_;
85 rateNode_["P"].setQuantity(pressure, "Pa");
86 rate.getRateParameters(rateNode_);
87 rateList.push_back(std::move(rateNode_));
88 }
89 rateNode["rate-constants"] = std::move(rateList);
90}
91
92void PlogRate::setRates(const std::multimap<double, ArrheniusRate>& rates)
93{
94 size_t j = 0;
95 rates_.clear();
96 pressures_.clear();
97 m_valid = !rates.empty();
98 rates_.reserve(rates.size());
99 // Insert intermediate pressures
100 for (const auto& [pressure, rate] : rates) {
101 double logp = std::log(pressure);
102 if (pressures_.empty() || pressures_.rbegin()->first != logp) {
103 // starting a new group
104 pressures_[logp] = {j, j+1};
105 } else {
106 // another rate expression at the same pressure
107 pressures_[logp].second = j+1;
108 }
109
110 j++;
111 rates_.push_back(rate);
112 }
113 if (!m_valid) {
114 // ensure that reaction rate can be evaluated (but returns NaN)
115 rates_.reserve(1);
116 pressures_[std::log(OneBar)] = {0, 0};
117 rates_.push_back(ArrheniusRate());
118 }
119
120 // Duplicate the first and last groups to handle P < P_0 and P > P_N
121 pressures_.insert({-1000.0, pressures_.begin()->second});
122 pressures_.insert({1000.0, pressures_.rbegin()->second});
123}
124
125void PlogRate::validate(const string& equation, const Kinetics& kin)
126{
127 if (!valid()) {
128 throw InputFileError("PlogRate::validate", m_input,
129 "Rate object for reaction '{}' is not configured.", equation);
130 }
131
132 fmt::memory_buffer err_reactions;
133 double T[] = {300.0, 500.0, 1000.0, 2000.0, 5000.0, 10000.0};
134 PlogData data;
135
136 for (auto iter = ++pressures_.begin(); iter->first < 1000; iter++) {
137 data.update(T[0], exp(iter->first)); // iter->first contains log(p)
138 updateFromStruct(data);
139 for (size_t i=0; i < 6; i++) {
140 double k = 0;
141 for (size_t p = ilow1_; p < ilow2_; p++) {
142 k += rates_.at(p).evalRate(log(T[i]), 1.0 / T[i]);
143 }
144 if (!(k > 0)) {
145 fmt_append(err_reactions,
146 "at P = {:.5g}, T = {:.1f}\n", std::exp(iter->first), T[i]);
147 }
148 }
149 }
150 if (err_reactions.size()) {
151 throw InputFileError("PlogRate::validate", m_input,
152 "\nInvalid rate coefficient for reaction '{}'\n{}\n"
153 "To fix this error, remove this reaction or contact the author of the\n"
154 "reaction/mechanism in question, because the rate expression is\n"
155 "mathematically unsound at the temperatures and pressures noted above.\n",
156 equation, to_string(err_reactions));
157 }
158}
159
160std::multimap<double, ArrheniusRate> PlogRate::getRates() const
161{
162 std::multimap<double, ArrheniusRate> rateMap;
163 // initial preincrement to skip rate for P --> 0
164 for (auto iter = ++pressures_.begin();
165 iter->first < 1000; // skip rates for (P --> infinity)
166 ++iter) {
167 for (size_t i = iter->second.first; i < iter->second.second; i++) {
168 rateMap.insert({std::exp(iter->first), rates_[i]});
169 }
170 }
171 return rateMap;
172}
173
174}
Header file for class ThermoPhase, the base class for phases with thermodynamic properties,...
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:427
bool hasKey(const string &key) const
Returns true if the map contains an item named key.
Definition AnyMap.cpp:1423
const UnitSystem & units() const
Return the default units that should be used to convert stored values.
Definition AnyMap.h:630
A wrapper for a variable whose type is determined at runtime.
Definition AnyMap.h:86
Arrhenius reaction rate type depends only on temperature.
Definition Arrhenius.h:170
Base class for exceptions thrown by Cantera classes.
Error thrown for problems processing information contained in an AnyMap or AnyValue.
Definition AnyMap.h:738
Public interface for kinetics managers.
Definition Kinetics.h:125
double temperature() const
Temperature (K).
Definition Phase.h:562
virtual double pressure() const
Return the thermodynamic pressure (Pa).
Definition Phase.h:580
map< double, pair< size_t, size_t > > pressures_
log(p) to (index range) in the rates_ vector
Definition PlogRate.h:181
void updateFromStruct(const PlogData &shared_data)
Update information specific to reaction.
Definition PlogRate.h:110
PlogRate()=default
Default constructor.
std::multimap< double, ArrheniusRate > getRates() const
Return the pressures and Arrhenius expressions which comprise this reaction.
Definition PlogRate.cpp:160
void setParameters(const AnyMap &node, const UnitStack &rate_units) override
Perform object setup based on AnyMap node information.
Definition PlogRate.cpp:62
void setRates(const std::multimap< double, ArrheniusRate > &rates)
Set up Plog object.
Definition PlogRate.cpp:92
void validate(const string &equation, const Kinetics &kin) override
Check to make sure that the rate expression is finite over a range of temperatures at each interpolat...
Definition PlogRate.cpp:125
size_t ilow1_
Indices to the ranges within rates_ for the lower / upper pressure, such that rates_[ilow1_] through ...
Definition PlogRate.h:194
virtual void setParameters(const AnyMap &node, const UnitStack &units)
Set parameters.
bool valid() const
Get flag indicating whether reaction rate is set up correctly.
bool m_valid
Flag indicating whether reaction rate is set up correctly.
AnyMap m_input
Input data used for specific models.
Base class for a phase with thermodynamic properties.
A representation of the units associated with a dimensional quantity.
Definition Units.h:35
void fmt_append(fmt::memory_buffer &b, Args... args)
Versions 6.2.0 and 6.2.1 of fmtlib do not include this define before they include windows....
Definition fmt.h:29
const double OneBar
One bar [Pa].
Definition ct_defs.h:99
Namespace for the Cantera kernel.
Definition AnyMap.cpp:564
Data container holding shared data specific to PlogRate.
Definition PlogRate.h:20
void perturbPressure(double deltaP)
Perturb pressure of data container.
Definition PlogRate.cpp:29
double m_pressure_buf
buffered pressure
Definition PlogRate.h:53
void update(double T) override
Update data container based on temperature T
Definition PlogRate.cpp:12
void restore() override
Restore data container after a perturbation.
Definition PlogRate.cpp:39
double pressure
pressure
Definition PlogRate.h:49
double temperature
temperature
virtual void restore()
Restore data container after a perturbation.
Unit aggregation utility.
Definition Units.h:105