Cantera  3.0.0
Loading...
Searching...
No Matches
VPStandardStateTP.cpp
Go to the documentation of this file.
1/**
2 * @file VPStandardStateTP.cpp
3 * Definition file for a derived class of ThermoPhase that handles
4 * variable pressure standard state methods for calculating
5 * thermodynamic properties (see @ref thermoprops and
6 * class @link Cantera::VPStandardStateTP VPStandardStateTP@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
13#include "cantera/thermo/PDSS.h"
16#include "cantera/base/global.h"
17
18namespace Cantera
19{
20
22{
23 // Defined in .cpp to limit dependence on PDSS.h via vector<unique_ptr<PDSS>>
24}
25
26VPStandardStateTP::~VPStandardStateTP()
27{
28 // Defined in .cpp to limit dependence on PDSS.h
29}
30
32{
34}
35
37{
38 warn_deprecated("VPStandardStateTP::getChemPotentials_RT",
39 "To be removed after Cantera 3.0. Use getChemPotentials instead.");
41 for (size_t k = 0; k < m_kk; k++) {
42 muRT[k] *= 1.0 / RT();
43 }
44}
45
46// ----- Thermodynamic Values for the Species Standard States States ----
47
49{
50 getGibbs_RT(g);
51 for (size_t k = 0; k < m_kk; k++) {
52 g[k] *= RT();
53 }
54}
55
56void VPStandardStateTP::getEnthalpy_RT(double* hrt) const
57{
59 std::copy(m_hss_RT.begin(), m_hss_RT.end(), hrt);
60}
61
62void VPStandardStateTP::getEntropy_R(double* sr) const
63{
65 std::copy(m_sss_R.begin(), m_sss_R.end(), sr);
66}
67
68void VPStandardStateTP::getGibbs_RT(double* grt) const
69{
71 std::copy(m_gss_RT.begin(), m_gss_RT.end(), grt);
72}
73
75{
77 std::copy(m_gss_RT.begin(), m_gss_RT.end(), g);
78 scale(g, g+m_kk, g, RT());
79}
80
82{
84 std::copy(m_hss_RT.begin(), m_hss_RT.end(), urt);
85 for (size_t k = 0; k < m_kk; k++) {
86 urt[k] -= m_Plast_ss / RT() * m_Vss[k];
87 }
88}
89
90void VPStandardStateTP::getCp_R(double* cpr) const
91{
93 std::copy(m_cpss_R.begin(), m_cpss_R.end(), cpr);
94}
95
96void VPStandardStateTP::getStandardVolumes(double* vol) const
97{
99 std::copy(m_Vss.begin(), m_Vss.end(), vol);
100}
101const vector<double>& VPStandardStateTP::getStandardVolumes() const
102{
104 return m_Vss;
105}
106
107// ----- Thermodynamic Values for the Species Reference States ----
108
110{
112 std::copy(m_h0_RT.begin(), m_h0_RT.end(), hrt);
113}
114
116{
118 std::copy(m_g0_RT.begin(), m_g0_RT.end(), grt);
119}
120
122{
124 std::copy(m_g0_RT.begin(), m_g0_RT.end(), g);
125 scale(g, g+m_kk, g, RT());
126}
127
128const vector<double>& VPStandardStateTP::Gibbs_RT_ref() const
129{
131 return m_g0_RT;
132}
133
135{
137 std::copy(m_s0_R.begin(), m_s0_R.end(), sr);
138}
139
140void VPStandardStateTP::getCp_R_ref(double* cpr) const
141{
143 std::copy(m_cp0_R.begin(), m_cp0_R.end(), cpr);
144}
145
147{
149 std::copy(m_Vss.begin(), m_Vss.end(), vol);
150}
151
153{
155 for (size_t k = 0; k < m_kk; k++) {
156 PDSS* kPDSS = m_PDSS_storage[k].get();
157 if (kPDSS == 0) {
158 throw CanteraError("VPStandardStateTP::initThermo",
159 "No PDSS object for species {}", k);
160 }
161 kPDSS->initThermo();
162 }
163}
164
166 AnyMap& speciesNode) const
167{
168 AnyMap eos;
169 providePDSS(speciesIndex(name))->getParameters(eos);
170 speciesNode["equation-of-state"].getMapWhere(
171 "model", eos.getString("model", ""), true) = std::move(eos);
172}
173
174bool VPStandardStateTP::addSpecies(shared_ptr<Species> spec)
175{
176 // Specifically skip ThermoPhase::addSpecies since the Species object
177 // doesn't have an associated SpeciesThermoInterpType object
178 bool added = Phase::addSpecies(spec);
179 if (!added) {
180 return false;
181 }
182
183 // VPStandardState does not use m_spthermo - install a dummy object
184 m_spthermo.install_STIT(m_kk-1, make_shared<SpeciesThermoInterpType>());
185 m_h0_RT.push_back(0.0);
186 m_cp0_R.push_back(0.0);
187 m_g0_RT.push_back(0.0);
188 m_s0_R.push_back(0.0);
189 m_V0.push_back(0.0);
190 m_hss_RT.push_back(0.0);
191 m_cpss_R.push_back(0.0);
192 m_gss_RT.push_back(0.0);
193 m_sss_R.push_back(0.0);
194 m_Vss.push_back(0.0);
195 return true;
196}
197
199{
200 setState_TP(temp, m_Pcurrent);
202}
203
205{
208}
209
211{
212 throw NotImplementedError("VPStandardStateTP::calcDensity");
213}
214
215void VPStandardStateTP::setState_TP(double t, double pres)
216{
217 // A pretty tricky algorithm is needed here, due to problems involving
218 // standard states of real fluids. For those cases you need to combine the T
219 // and P specification for the standard state, or else you may venture into
220 // the forbidden zone, especially when nearing the triple point. Therefore,
221 // we need to do the standard state thermo calc with the (t, pres) combo.
223 m_Pcurrent = pres;
225
226 // Now, we still need to do the calculations for general ThermoPhase
227 // objects. So, we switch back to a virtual function call, setTemperature,
228 // and setPressure to recalculate stuff for child ThermoPhase objects of the
229 // VPStandardStateTP object. At this point, we haven't touched m_tlast or
230 // m_plast, so some calculations may still need to be done at the
231 // ThermoPhase object level.
232 calcDensity();
233}
234
235void VPStandardStateTP::installPDSS(size_t k, unique_ptr<PDSS>&& pdss)
236{
237 pdss->setParent(this, k);
238 pdss->setMolecularWeight(molecularWeight(k));
239 Species& spec = *species(k);
240 if (spec.thermo) {
241 pdss->setReferenceThermo(spec.thermo);
242 spec.thermo->validate(spec.name);
243 }
244 m_minTemp = std::max(m_minTemp, pdss->minTemp());
245 m_maxTemp = std::min(m_maxTemp, pdss->maxTemp());
246
247 if (m_PDSS_storage.size() < k+1) {
248 m_PDSS_storage.resize(k+1);
249 }
250 m_PDSS_storage[k].swap(pdss);
251}
252
253PDSS* VPStandardStateTP::providePDSS(size_t k)
254{
255 return m_PDSS_storage[k].get();
256}
257
258const PDSS* VPStandardStateTP::providePDSS(size_t k) const
259{
260 return m_PDSS_storage[k].get();
261}
262
264{
266 m_Tlast_ss += 0.0001234;
267}
268
270{
271 double Tnow = temperature();
272 for (size_t k = 0; k < m_kk; k++) {
273 PDSS* kPDSS = m_PDSS_storage[k].get();
274 kPDSS->setState_TP(Tnow, m_Pcurrent);
275 // reference state thermo
276 if (Tnow != m_tlast) {
277 m_h0_RT[k] = kPDSS->enthalpy_RT_ref();
278 m_s0_R[k] = kPDSS->entropy_R_ref();
279 m_g0_RT[k] = m_h0_RT[k] - m_s0_R[k];
280 m_cp0_R[k] = kPDSS->cp_R_ref();
281 m_V0[k] = kPDSS->molarVolume_ref();
282 }
283 // standard state thermo
284 m_hss_RT[k] = kPDSS->enthalpy_RT();
285 m_sss_R[k] = kPDSS->entropy_R();
286 m_gss_RT[k] = m_hss_RT[k] - m_sss_R[k];
287 m_cpss_R[k] = kPDSS->cp_R();
288 m_Vss[k] = kPDSS->molarVolume();
289 }
291 m_Tlast_ss = Tnow;
292 m_tlast = Tnow;
293}
294
296{
297 double Tnow = temperature();
298 if (Tnow != m_Tlast_ss || Tnow != m_tlast || m_Pcurrent != m_Plast_ss) {
300 }
301}
302
303double VPStandardStateTP::minTemp(size_t k) const
304{
305 if (k == npos) {
306 return m_minTemp;
307 } else {
308 return m_PDSS_storage.at(k)->minTemp();
309 }
310}
311
312double VPStandardStateTP::maxTemp(size_t k) const
313{
314 if (k == npos) {
315 return m_maxTemp;
316 } else {
317 return m_PDSS_storage.at(k)->maxTemp();
318 }
319}
320
321}
Declarations for the virtual base class PDSS (pressure dependent standard state) which handles calcul...
Declaration for class Cantera::Species.
Header file for a derived class of ThermoPhase that handles variable pressure standard state methods ...
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:427
const string & getString(const string &key, const string &default_) const
If key exists, return it as a string, otherwise return default_.
Definition AnyMap.cpp:1530
Base class for exceptions thrown by Cantera classes.
virtual void install_STIT(size_t index, shared_ptr< SpeciesThermoInterpType > stit)
Install a new species thermodynamic property parameterization for one species.
An error indicating that an unimplemented function has been called.
Virtual base class for a species with a pressure dependent standard state.
Definition PDSS.h:147
virtual void initThermo()
Initialization routine.
Definition PDSS.h:418
virtual double cp_R_ref() const
Return the molar heat capacity divided by R at reference pressure.
Definition PDSS.cpp:93
virtual double entropy_R_ref() const
Return the molar entropy divided by R at reference pressure.
Definition PDSS.cpp:88
virtual double molarVolume_ref() const
Return the molar volume at reference pressure.
Definition PDSS.cpp:98
virtual double enthalpy_RT() const
Return the standard state molar enthalpy divided by RT.
Definition PDSS.cpp:23
virtual double enthalpy_RT_ref() const
Return the molar enthalpy divided by RT at reference pressure.
Definition PDSS.cpp:83
virtual double entropy_R() const
Return the standard state entropy divided by RT.
Definition PDSS.cpp:38
virtual double cp_R() const
Return the molar const pressure heat capacity divided by RT.
Definition PDSS.cpp:58
virtual double molarVolume() const
Return the molar volume at standard state.
Definition PDSS.cpp:63
virtual void getParameters(AnyMap &eosNode) const
Store the parameters needed to reconstruct a copy of this PDSS object.
Definition PDSS.h:427
virtual void setState_TP(double temp, double pres)
Set the internal temperature and pressure.
Definition PDSS.cpp:176
virtual bool addSpecies(shared_ptr< Species > spec)
Add a Species to this Phase.
Definition Phase.cpp:822
size_t m_kk
Number of species in the phase.
Definition Phase.h:947
double temperature() const
Temperature (K).
Definition Phase.h:662
size_t speciesIndex(const string &name) const
Returns the index of a species named 'name' within the Phase object.
Definition Phase.cpp:138
virtual void setTemperature(double temp)
Set the internally stored temperature of the phase (K).
Definition Phase.h:728
double molecularWeight(size_t k) const
Molecular weight of species k.
Definition Phase.cpp:482
shared_ptr< Species > species(const string &name) const
Return the Species object for the named species.
Definition Phase.cpp:977
string name() const
Return the name of the phase.
Definition Phase.cpp:20
Contains data about a single chemical species.
Definition Species.h:25
string name
The name of the species.
Definition Species.h:41
shared_ptr< SpeciesThermoInterpType > thermo
Thermodynamic data for the species.
Definition Species.h:80
double RT() const
Return the Gas Constant multiplied by the current temperature.
double m_tlast
last value of the temperature processed by reference state
virtual void initThermo()
Initialize the ThermoPhase object after all species have been set up.
void invalidateCache() override
Invalidate any cached values which are normally updated only when a change in state is detected.
virtual void getChemPotentials(double *mu) const
Get the species chemical potentials. Units: J/kmol.
MultiSpeciesThermo m_spthermo
Pointer to the calculation manager for species reference-state thermodynamic properties.
double m_Plast_ss
The last pressure at which the Standard State thermodynamic properties were calculated at.
int standardStateConvention() const override
This method returns the convention used in specification of the standard state, of which there are cu...
vector< double > m_g0_RT
Vector containing the species reference Gibbs functions at T = m_tlast and P = p_ref.
vector< double > m_sss_R
Vector containing the species Standard State entropies at T = m_tlast and P = m_plast.
void installPDSS(size_t k, unique_ptr< PDSS > &&pdss)
Install a PDSS object for species k
void getSpeciesParameters(const string &name, AnyMap &speciesNode) const override
Get phase-specific parameters of a Species object such that an identical one could be reconstructed a...
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_h0_RT
Vector containing the species reference enthalpies at T = m_tlast and P = p_ref.
virtual void _updateStandardStateThermo() const
Updates the standard state thermodynamic functions at the current T and P of the solution.
void getGibbs_ref(double *g) const override
Returns the vector of the Gibbs function of the reference state at the current temperature of the sol...
void getStandardChemPotentials(double *mu) const override
Get the array of chemical potentials at unit activity for the species at their standard states at the...
void getCp_R(double *cpr) const override
Get the nondimensional Heat Capacities at constant pressure for the species standard states at the cu...
void initThermo() override
Initialize the ThermoPhase object after all species have been set up.
void setPressure(double p) override
Set the internally stored pressure (Pa) at constant temperature and composition.
vector< unique_ptr< PDSS > > m_PDSS_storage
Storage for the PDSS objects for the species.
void getStandardVolumes_ref(double *vol) const override
Get the molar volumes of the species reference states at the current T and P_ref of the solution.
vector< double > m_gss_RT
Vector containing the species Standard State Gibbs functions at T = m_tlast and P = m_plast.
double m_Tlast_ss
The last temperature at which the standard state thermodynamic properties were calculated at.
void getPureGibbs(double *gpure) const override
Get the Gibbs functions for the standard state of the species at the current T and P of the solution.
void getEnthalpy_RT(double *hrt) const override
Get the nondimensional Enthalpy functions for the species at their standard states at the current T a...
vector< double > m_cpss_R
Vector containing the species Standard State constant pressure heat capacities at T = m_tlast and P =...
double m_maxTemp
The maximum temperature at which data for all species is valid.
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 setTemperature(const double temp) override
Set the temperature of the phase.
double minTemp(size_t k=npos) const override
Minimum temperature for which the thermodynamic data for the species or phase are valid.
vector< double > m_s0_R
Vector containing the species reference entropies at T = m_tlast and P = p_ref.
vector< double > m_Vss
Vector containing the species standard state volumes at T = m_tlast and P = m_plast.
double m_minTemp
The minimum temperature at which data for all species is valid.
void getGibbs_RT(double *grt) const override
Get the nondimensional Gibbs functions for the species in their standard states at the current T and ...
vector< double > m_V0
Vector containing the species reference molar volumes.
void invalidateCache() override
Invalidate any cached values which are normally updated only when a change in state is detected.
void getCp_R_ref(double *cprt) const override
Returns the vector of nondimensional constant pressure heat capacities of the reference state at the ...
void setState_TP(double T, double pres) override
Set the temperature and pressure at the same time.
void getIntEnergy_RT(double *urt) const override
Returns the vector of nondimensional Internal Energies of the standard state species at the current T...
vector< double > m_cp0_R
Vector containing the species reference constant pressure heat capacities at T = m_tlast and P = p_re...
virtual void updateStandardStateThermo() const
Updates the standard state thermodynamic functions at the current T and P of the solution.
vector< double > m_hss_RT
Vector containing the species Standard State enthalpies at T = m_tlast and P = m_plast.
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
void getChemPotentials_RT(double *mu) const override
Get the array of non-dimensional species chemical potentials.
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 maxTemp(size_t k=npos) const override
Maximum temperature for which the thermodynamic data for the species are valid.
double m_Pcurrent
Current value of the pressure - state variable.
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 calcDensity()
Calculate the density of the mixture using the partial molar volumes and mole fractions as input.
This file contains definitions for utility functions and text for modules, inputfiles and logging,...
void scale(InputIter begin, InputIter end, OutputIter out, S scale_factor)
Multiply elements of an array by a scale factor.
Definition utilities.h:104
Namespace for the Cantera kernel.
Definition AnyMap.cpp:564
const size_t npos
index returned by functions to indicate "no position"
Definition ct_defs.h:195
const int cSS_CONVENTION_VPSS
Standard state uses the molality convention.
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
Various templated functions that carry out common vector and polynomial operations (see Templated Arr...