Cantera  3.0.0
Loading...
Searching...
No Matches
SurfPhase.cpp
Go to the documentation of this file.
1/**
2 * @file SurfPhase.cpp
3 * Definitions for a simple thermodynamic model of a surface phase
4 * derived from ThermoPhase, assuming an ideal solution model
5 * (see @ref thermoprops and class
6 * @link Cantera::SurfPhase SurfPhase@endlink).
7 */
8
9// This file is part of Cantera. See License.txt in the top-level directory or
10// at https://cantera.org/license.txt for license and copyright information.
11
18
19namespace Cantera
20{
21
22SurfPhase::SurfPhase(const string& infile, const string& id_)
23{
24 setNDim(2);
25 initThermoFile(infile, id_);
26}
27
29{
30 if (m_n0 <= 0.0) {
31 return 0.0;
32 }
34 return mean_X(m_h0);
35}
36
38{
39 return enthalpy_mole();
40}
41
43{
45 double s = 0.0;
46 for (size_t k = 0; k < m_kk; k++) {
47 s += moleFraction(k) * (m_s0[k] -
48 GasConstant * log(std::max(concentration(k) * size(k)/m_n0, SmallNumber)));
49 }
50 return s;
51}
52
53double SurfPhase::cp_mole() const
54{
56 return mean_X(m_cp0);
57}
58
59double SurfPhase::cv_mole() const
60{
61 return cp_mole();
62}
63
65{
66 getEnthalpy_RT(hbar);
67 for (size_t k = 0; k < m_kk; k++) {
68 hbar[k] *= RT();
69 }
70}
71
73{
74 getEntropy_R(sbar);
76 for (size_t k = 0; k < m_kk; k++) {
77 sbar[k] -= log(std::max(m_work[k], SmallNumber)) - logStandardConc(k);
78 sbar[k] *= GasConstant;
79 }
80}
81
82void SurfPhase::getPartialMolarCp(double* cpbar) const
83{
84 getCp_R(cpbar);
85 for (size_t k = 0; k < m_kk; k++) {
86 cpbar[k] *= GasConstant;
87 }
88}
89
90// HKM 9/1/11 The partial molar volumes returned here are really partial molar areas.
91// Partial molar volumes for this phase should actually be equal to zero.
92void SurfPhase::getPartialMolarVolumes(double* vbar) const
93{
95}
96
98{
100 copy(m_mu0.begin(), m_mu0.end(), mu0);
101}
102
103void SurfPhase::getChemPotentials(double* mu) const
104{
106 copy(m_mu0.begin(), m_mu0.end(), mu);
108 for (size_t k = 0; k < m_kk; k++) {
109 mu[k] += RT() * (log(std::max(m_work[k], SmallNumber)) - logStandardConc(k));
110 }
111}
112
114{
116}
117
119{
120 return m_n0/size(k);
121}
122
123double SurfPhase::logStandardConc(size_t k) const
124{
125 return m_logn0 - m_logsize[k];
126}
127
128void SurfPhase::getPureGibbs(double* g) const
129{
131 copy(m_mu0.begin(), m_mu0.end(), g);
132}
133
134void SurfPhase::getGibbs_RT(double* grt) const
135{
137 scale(m_mu0.begin(), m_mu0.end(), grt, 1.0/RT());
138}
139
140void SurfPhase::getEnthalpy_RT(double* hrt) const
141{
143 scale(m_h0.begin(), m_h0.end(), hrt, 1.0/RT());
144}
145
146void SurfPhase::getEntropy_R(double* sr) const
147{
149 scale(m_s0.begin(), m_s0.end(), sr, 1.0/GasConstant);
150}
151
152void SurfPhase::getCp_R(double* cpr) const
153{
155 scale(m_cp0.begin(), m_cp0.end(), cpr, 1.0/GasConstant);
156}
157
158void SurfPhase::getStandardVolumes(double* vol) const
159{
160 for (size_t k = 0; k < m_kk; k++) {
161 vol[k] = 0.0;
162 }
163}
164
165void SurfPhase::getGibbs_RT_ref(double* grt) const
166{
167 getGibbs_RT(grt);
168}
169
170void SurfPhase::getEnthalpy_RT_ref(double* hrt) const
171{
172 getEnthalpy_RT(hrt);
173}
174
175void SurfPhase::getEntropy_R_ref(double* sr) const
176{
177 getEntropy_R(sr);
178}
179
180void SurfPhase::getCp_R_ref(double* cprt) const
181{
182 getCp_R(cprt);
183}
184
185bool SurfPhase::addSpecies(shared_ptr<Species> spec)
186{
187 bool added = ThermoPhase::addSpecies(spec);
188 if (added) {
189 m_h0.push_back(0.0);
190 m_s0.push_back(0.0);
191 m_cp0.push_back(0.0);
192 m_mu0.push_back(0.0);
193 m_work.push_back(0.0);
194 m_speciesSize.push_back(spec->size);
195 m_logsize.push_back(log(spec->size));
196 if (m_kk == 1) {
197 vector<double> cov{1.0};
198 setCoverages(cov.data());
199 }
200 }
201 return added;
202}
203
204void SurfPhase::setMolarDensity(const double vm) {
205 warn_deprecated("SurfPhase::setMolarDensity", "To be removed after Cantera 3.0");
206 if (vm != 0.0) {
207 throw CanteraError("SurfPhase::setMolarDensity",
208 "The volume of an interface is zero");
209 }
210}
211
213{
214 if (n0 <= 0.0) {
215 throw CanteraError("SurfPhase::setSiteDensity",
216 "Site density must be positive. Got {}", n0);
217 }
218 m_n0 = n0;
220 m_logn0 = log(m_n0);
221}
222
223void SurfPhase::setCoverages(const double* theta)
224{
225 double sum = 0.0;
226 for (size_t k = 0; k < m_kk; k++) {
227 sum += theta[k] / size(k);
228 }
229 if (sum <= 0.0) {
230 throw CanteraError("SurfPhase::setCoverages",
231 "Sum of Coverage fractions is zero or negative");
232 }
233 for (size_t k = 0; k < m_kk; k++) {
234 m_work[k] = theta[k] / (sum * size(k));
235 }
236 setMoleFractions(m_work.data());
237}
238
239void SurfPhase::setCoveragesNoNorm(const double* theta)
240{
241 double sum = 0.0;
242 double sum2 = 0.0;
243 for (size_t k = 0; k < m_kk; k++) {
244 sum += theta[k] / size(k);
245 sum2 += theta[k];
246 }
247 if (sum <= 0.0) {
248 throw CanteraError("SurfPhase::setCoverages",
249 "Sum of Coverage fractions is zero or negative");
250 }
251 for (size_t k = 0; k < m_kk; k++) {
252 m_work[k] = theta[k] * sum2 / (sum * size(k));
253 }
255}
256
257void SurfPhase::getCoverages(double* theta) const
258{
259 double sum_X = 0.0;
260 double sum_X_s = 0.0;
261 getMoleFractions(theta);
262 for (size_t k = 0; k < m_kk; k++) {
263 sum_X += theta[k];
264 sum_X_s += theta[k] * size(k);
265 }
266 for (size_t k = 0; k < m_kk; k++) {
267 theta[k] *= size(k) * sum_X / sum_X_s;
268 }
269}
270
271void SurfPhase::setCoveragesByName(const string& cov)
272{
274}
275
277{
278 vector<double> cv(m_kk, 0.0);
279 bool ifound = false;
280 for (size_t k = 0; k < m_kk; k++) {
281 double c = getValue(cov, speciesName(k), 0.0);
282 if (c > 0.0) {
283 ifound = true;
284 cv[k] = c;
285 }
286 }
287 if (!ifound) {
288 throw CanteraError("SurfPhase::setCoveragesByName",
289 "Input coverages are all zero or negative");
290 }
291 setCoverages(cv.data());
292}
293
294void SurfPhase::setState(const AnyMap& state) {
295 if (state.hasKey("coverages")) {
296 if (state["coverages"].is<string>()) {
297 setCoveragesByName(state["coverages"].asString());
298 } else {
299 setCoveragesByName(state["coverages"].asMap<double>());
300 }
301 }
303}
304
306{
309}
310
311void SurfPhase::_updateThermo(bool force) const
312{
313 double tnow = temperature();
314 if (m_tlast != tnow || force) {
315 m_spthermo.update(tnow, m_cp0.data(), m_h0.data(), m_s0.data());
316 m_tlast = tnow;
317 for (size_t k = 0; k < m_kk; k++) {
318 m_h0[k] *= GasConstant * tnow;
319 m_s0[k] *= GasConstant;
320 m_cp0[k] *= GasConstant;
321 m_mu0[k] = m_h0[k] - tnow*m_s0[k];
322 }
323 m_tlast = tnow;
324 }
325}
326
328{
329 if (m_input.hasKey("site-density")) {
330 // Units are kmol/m^2 for surface phases or kmol/m for edge phases
331 setSiteDensity(m_input.convert("site-density",
332 Units(1.0, 0, -static_cast<double>(m_ndim), 0, 0, 0, 1)));
333 }
334}
335
336void SurfPhase::getParameters(AnyMap& phaseNode) const
337{
339 phaseNode["site-density"].setQuantity(
340 m_n0, Units(1.0, 0, -static_cast<double>(m_ndim), 0, 0, 0, 1));
341}
342
343EdgePhase::EdgePhase(const string& infile, const string& id_)
344{
345 setNDim(1);
346 initThermoFile(infile, id_);
347}
348
349}
Declarations for the EdgePhase ThermoPhase object, which models the interface between two surfaces (s...
Declaration for class Cantera::Species.
Header for a simple thermodynamics model of a surface phase derived from ThermoPhase,...
Headers for the factory class that can create known ThermoPhase objects (see 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
double convert(const string &key, const string &units) const
Convert the item stored by the given key to the units specified in units.
Definition AnyMap.cpp:1535
Base class for exceptions thrown by Cantera classes.
EdgePhase(const string &infile="", const string &id="")
Construct and initialize an EdgePhase directly from an input file.
virtual void update(double T, double *cp_R, double *h_RT, double *s_R) const
Compute the reference-state properties for all species.
virtual void getConcentrations(double *const c) const
Get the species concentrations (kmol/m^3).
Definition Phase.cpp:595
void assignDensity(const double density_)
Set the internally stored constant density (kg/m^3) of the phase.
Definition Phase.cpp:718
virtual void setMoleFractions(const double *const x)
Set the mole fractions to the specified values.
Definition Phase.cpp:304
void setNDim(size_t ndim)
Set the number of spatial dimensions (1, 2, or 3).
Definition Phase.h:653
size_t m_kk
Number of species in the phase.
Definition Phase.h:947
size_t m_ndim
Dimensionality of the phase.
Definition Phase.h:951
double temperature() const
Temperature (K).
Definition Phase.h:662
double meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
Definition Phase.h:760
virtual double concentration(const size_t k) const
Concentration of species k.
Definition Phase.cpp:589
string speciesName(size_t k) const
Name of the species with index k.
Definition Phase.cpp:151
void getMoleFractions(double *const x) const
Get the species mole fraction vector.
Definition Phase.cpp:540
double moleFraction(size_t k) const
Return the mole fraction of a single species.
Definition Phase.cpp:545
virtual void compositionChanged()
Apply changes to the state which are needed after the composition changes.
Definition Phase.cpp:1026
double mean_X(const double *const Q) const
Evaluate the mole-fraction-weighted mean of an array Q.
Definition Phase.cpp:737
const vector< string > & speciesNames() const
Return a const reference to the vector of species names.
Definition Phase.cpp:157
virtual void setMoleFractions_NoNorm(const double *const x)
Set the mole fractions to the specified values without normalizing.
Definition Phase.cpp:336
SurfPhase(const string &infile="", const string &id="")
Construct and initialize a SurfPhase ThermoPhase object directly from an input file.
Definition SurfPhase.cpp:22
void setMolarDensity(const double vm) override
Since interface phases have no volume, setting this to a value other than 0.0 raises an exception.
void getStandardChemPotentials(double *mu0) const override
Get the array of chemical potentials at unit activity for the species at their standard states at the...
Definition SurfPhase.cpp:97
void getPureGibbs(double *g) const override
Get the Gibbs functions for the standard state of the species at the current T and P of the solution.
double enthalpy_mole() const override
Return the Molar Enthalpy. Units: J/kmol.
Definition SurfPhase.cpp:28
void setSiteDensity(double n0)
Set the site density of the surface phase (kmol m-2)
double logStandardConc(size_t k=0) const override
Natural logarithm of the standard concentration of the kth species.
void setState(const AnyMap &state) override
Set the state using an AnyMap containing any combination of properties supported by the thermodynamic...
void getPartialMolarEnthalpies(double *hbar) const override
Returns an array of partial molar enthalpies for the species in the mixture.
Definition SurfPhase.cpp:64
void getChemPotentials(double *mu) const override
Get the species chemical potentials. Units: J/kmol.
void getEntropy_R(double *sr) const override
Get the array of nondimensional Entropy functions for the standard state species at the current T and...
vector< double > m_logsize
vector storing the log of the size of each species.
Definition SurfPhase.h:347
vector< double > m_work
Temporary work array.
Definition SurfPhase.h:340
void getCp_R(double *cpr) const override
Get the nondimensional Heat Capacities at constant pressure for the species standard states at the cu...
void getParameters(AnyMap &phaseNode) const override
Store the parameters of a ThermoPhase object such that an identical one could be reconstructed using ...
void initThermo() override
Initialize the ThermoPhase object after all species have been set up.
void getActivityConcentrations(double *c) const override
Return a vector of activity concentrations for each species.
double m_n0
Surface site density (kmol m-2)
Definition SurfPhase.h:316
double size(size_t k) const
Returns the number of sites occupied by one molecule of species k.
Definition SurfPhase.h:226
vector< double > m_h0
Temporary storage for the reference state enthalpies.
Definition SurfPhase.h:328
void getPartialMolarVolumes(double *vbar) const override
Return an array of partial molar volumes for the species in the mixture.
Definition SurfPhase.cpp:92
double cv_mole() const override
Molar heat capacity at constant volume. Units: J/kmol/K.
Definition SurfPhase.cpp:59
vector< double > m_s0
Temporary storage for the reference state entropies.
Definition SurfPhase.h:331
void setCoverages(const double *theta)
Set the surface site fractions to a specified state.
vector< double > m_cp0
Temporary storage for the reference state heat capacities.
Definition SurfPhase.h:334
vector< double > m_speciesSize
Vector of species sizes (number of sites occupied). length m_kk.
Definition SurfPhase.h:319
void getEnthalpy_RT(double *hrt) const override
Get the nondimensional Enthalpy functions for the species at their standard states at the current T a...
void getEntropy_R_ref(double *er) const override
Returns the vector of nondimensional entropies of the reference state at the current temperature of t...
void getGibbs_RT(double *grt) const override
Get the nondimensional Gibbs functions for the species in their standard states at the current T and ...
double intEnergy_mole() const override
Return the Molar Internal Energy. Units: J/kmol.
Definition SurfPhase.cpp:37
double entropy_mole() const override
Return the Molar Entropy. Units: J/kmol-K.
Definition SurfPhase.cpp:42
void _updateThermo(bool force=false) const
Update the species reference state thermodynamic functions.
void getCp_R_ref(double *cprt) const override
Returns the vector of nondimensional constant pressure heat capacities of the reference state at the ...
void getStandardVolumes(double *vol) const override
Get the molar volumes of the species standard states at the current T and P of the solution.
void setCoveragesNoNorm(const double *theta)
Set the surface site fractions to a specified state.
void getCoverages(double *theta) const
Return a vector of surface coverages.
double cp_mole() const override
Molar heat capacity at constant pressure. Units: J/kmol/K.
Definition SurfPhase.cpp:53
void getPartialMolarCp(double *cpbar) const override
Return an array of partial molar heat capacities for the species in the mixture.
Definition SurfPhase.cpp:82
void compositionChanged() override
Apply changes to the state which are needed after the composition changes.
double standardConcentration(size_t k=0) const override
Return the standard concentration for the kth species.
vector< double > m_mu0
Temporary storage for the reference state Gibbs energies.
Definition SurfPhase.h:337
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
void setCoveragesByName(const string &cov)
Set the coverages from a string of colon-separated name:value pairs.
void getGibbs_RT_ref(double *grt) const override
Returns the vector of nondimensional Gibbs Free Energies of the reference state at the current temper...
double m_logn0
log of the surface site density
Definition SurfPhase.h:322
void getPartialMolarEntropies(double *sbar) const override
Returns an array of partial molar entropies of the species in the solution.
Definition SurfPhase.cpp:72
void getEnthalpy_RT_ref(double *hrt) const override
Returns the vector of nondimensional enthalpies of the reference state at the current temperature of ...
virtual void getParameters(AnyMap &phaseNode) const
Store the parameters of a ThermoPhase object such that an identical one could be reconstructed using ...
virtual void setState(const AnyMap &state)
Set the state using an AnyMap containing any combination of properties supported by the thermodynamic...
double RT() const
Return the Gas Constant multiplied by the current temperature.
double m_tlast
last value of the temperature processed by reference state
void initThermoFile(const string &inputFile, const string &id)
Initialize a ThermoPhase object using an input file.
MultiSpeciesThermo m_spthermo
Pointer to the calculation manager for species reference-state thermodynamic properties.
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
AnyMap m_input
Data supplied via setParameters.
A representation of the units associated with a dimensional quantity.
Definition Units.h:35
Composition parseCompString(const string &ss, const vector< string > &names)
Parse a composition string into a map consisting of individual key:composition pairs.
void scale(InputIter begin, InputIter end, OutputIter out, S scale_factor)
Multiply elements of an array by a scale factor.
Definition utilities.h:104
const double GasConstant
Universal Gas Constant [J/kmol/K].
Definition ct_defs.h:120
Namespace for the Cantera kernel.
Definition AnyMap.cpp:564
const double SmallNumber
smallest number to compare to zero.
Definition ct_defs.h:158
void warn_deprecated(const string &source, const AnyBase &node, const string &message)
A deprecation warning for syntax in an input file.
Definition AnyMap.cpp:1926
const U & getValue(const map< T, U > &m, const T &key, const U &default_val)
Const accessor for a value in a map.
Definition utilities.h:190
map< string, double > Composition
Map from string names to doubles.
Definition ct_defs.h:184
Contains declarations for string manipulation functions within Cantera.
Various templated functions that carry out common vector and polynomial operations (see Templated Arr...