Cantera  2.1.2
ShomateThermo.h
Go to the documentation of this file.
1 /**
2  * @file ShomateThermo.h
3  * Header for the 2 regions Shomate polynomial
4  * for multiple species in a phase, derived from the
7  */
8 // Copyright 2001 California Institute of Technology
9
10 #ifndef CT_SHOMATETHERMO_H
11 #define CT_SHOMATETHERMO_H
12
14 #include "ShomatePoly.h"
16
17 namespace Cantera
18 {
19 //! A species thermodynamic property manager for the Shomate polynomial parameterization.
20 /*!
21  * This is the parameterization used
22  * in the NIST Chemistry WebBook (http://webbook.nist.gov/chemistry)
23  * The parameterization assumes there are two temperature regions
24  * each with its own Shomate polynomial representation, for each
25  * species in the phase.
26  *
27  * \f[
28  * \tilde{c}_p^0(T) = A + B t + C t^2 + D t^3 + \frac{E}{t^2}
29  * \f]
30  * \f[
31  * \tilde{h}^0(T) = A t + \frac{B t^2}{2} + \frac{C t^3}{3}
32  + \frac{D t^4}{4} - \frac{E}{t} + F.
33  * \f]
34  * \f[
35  * \tilde{s}^0(T) = A\ln t + B t + \frac{C t^2}{2}
36  + \frac{D t^3}{3} - \frac{E}{2t^2} + G.
37  * \f]
38  *
39  * In the above expressions, the thermodynamic polynomials are expressed
40  * in dimensional units, but the temperature,\f$t \f$, is divided by 1000. The
41  * following dimensions are assumed in the above expressions:
42  *
43  * - \f$\tilde{c}_p^0(T)\f$ = Heat Capacity (J/gmol*K)
44  * - \f$\tilde{h}^0(T) \f$ = standard Enthalpy (kJ/gmol)
45  * - \f$\tilde{s}^0(T) \f$= standard Entropy (J/gmol*K)
46  * - \f$t \f$= temperature (K) / 1000.
47  *
48  * Note, the polynomial data (i.e., A, ... , G) is entered in dimensional form.
49  *
50  * This is in contrast to the NASA database polynomials which are entered in
51  * nondimensional form (i.e., NASA parameterizes C_p/R, while Shomate
52  * parameterizes C_p assuming units of J/gmol*K - and kJ/gmol*K for H).
53  * Note, also that the H - H_298.15 equation has units of kJ/gmol, because of
54  * the implicit integration of (t = T 1000), which provides a
55  * multiplier of 1000 to the Enthalpy equation.
56  *
57  * @ingroup mgrsrefcalc
58  */
60 {
61 public:
62  //! Initialized to the type of parameterization
63  /*!
64  * Note, this value is used in some template functions
65  */
66  const int ID;
67
68  //! constructor
70  ID(SHOMATE),
71  m_tlow_max(0.0),
72  m_thigh_min(1.e30),
73  m_p0(-1.0),
74  m_ngroups(0) {
75  m_t.resize(7);
76  }
77
78  //! Copy Constructor
79  /*!
80  * @param right Object to be copied
81  */
82  ShomateThermo(const ShomateThermo& right) :
83  ID(SHOMATE),
84  m_tlow_max(0.0),
85  m_thigh_min(1.e30),
86  m_p0(-1.0),
87  m_ngroups(0) {
88  *this = operator=(right);
89  }
90
91  //! Assignment Operator
92  /*!
93  * @param right Object to be copied
94  */
96  if (&right == this) {
97  return *this;
98  }
99
100  m_high = right.m_high;
101  m_low = right.m_low;
102  m_index = right.m_index;
103  m_tmid = right.m_tmid;
104  m_tlow_max = right.m_tlow_max;
105  m_thigh_min = right.m_thigh_min;
106  m_tlow = right.m_tlow;
107  m_thigh = right.m_thigh;
108  m_p0 = right.m_p0;
109  m_ngroups = right.m_ngroups;
110  m_t = right.m_t;
111  m_group_map = right.m_group_map;
113
114  return *this;
115  }
116
118  ShomateThermo* st = new ShomateThermo(*this);
119  return (SpeciesThermo*) st;
120  }
121
122  //! Install a new species thermodynamic property
123  //! parameterization for one species using Shomate polynomials
124  /*!
125  * Two temperature regions are assumed.
126  *
127  * @param name Name of the species
128  * @param index Species index
129  * @param type int flag specifying the type of parameterization to be
130  * installed.
131  * @param c Vector of coefficients for the parameterization.
132  * There are 15 coefficients for the 2-zone Shomate polynomial.
133  * The first coefficient is the value of Tmid. The next 7
134  * coefficients are the low temperature range Shomate coefficients.
135  * The last 7 are the high temperature range Shomate coefficients.
136  *
137  * @param minTemp minimum temperature for which this parameterization
138  * is valid.
139  * @param maxTemp maximum temperature for which this parameterization
140  * is valid.
141  * @param refPressure standard-state pressure for this
142  * parameterization.
143  *
144  * @see ShomatePoly
145  * @see ShomatePoly2
146  */
147  virtual void install(const std::string& name, size_t index, int type,
148  const doublereal* c,
149  doublereal minTemp, doublereal maxTemp,
150  doublereal refPressure) {
151  int imid = int(c[0]); // midpoint temp converted to integer
152  int igrp = m_index[imid]; // has this value been seen before?
153  if (igrp == 0) { // if not, prepare new group
154  std::vector<ShomatePoly> v;
155  m_high.push_back(v);
156  m_low.push_back(v);
157  m_tmid.push_back(c[0]);
158  m_index[imid] = igrp = static_cast<int>(m_high.size());
159  m_ngroups++;
160  }
161  m_group_map[index] = igrp;
162  m_posInGroup_map[index] = (int) m_low[igrp-1].size();
163  doublereal tlow = minTemp;
164  doublereal tmid = c[0];
165  doublereal thigh = maxTemp;
166
167  const doublereal* clow = c + 1;
168  const doublereal* chigh = c + 8;
169  m_high[igrp-1].push_back(ShomatePoly(index, tmid, thigh,
170  refPressure, chigh));
171  m_low[igrp-1].push_back(ShomatePoly(index, tlow, tmid,
172  refPressure, clow));
173  if (tlow > m_tlow_max) {
174  m_tlow_max = tlow;
175  }
176  if (thigh < m_thigh_min) {
177  m_thigh_min = thigh;
178  }
179
180  if (m_tlow.size() < index + 1) {
181  m_tlow.resize(index + 1, tlow);
182  m_thigh.resize(index + 1, thigh);
183  }
184  m_tlow[index] = tlow;
185  m_thigh[index] = thigh;
186
187  if (m_p0 < 0.0) {
188  m_p0 = refPressure;
189  } else if (fabs(m_p0 - refPressure) > 0.1) {
190  std::string logmsg = " ERROR ShomateThermo: New Species, " + name
191  + ", has a different reference pressure, "
192  + fp2str(refPressure) + ", than existing reference pressure, " + fp2str(m_p0) + "\n";
193  writelog(logmsg);
194  logmsg = " This is now a fatal error\n";
195  writelog(logmsg);
196  throw CanteraError("install()", "Species have different reference pressures");
197  }
198  m_p0 = refPressure;
199
200  }
201
202  virtual void install_STIT(SpeciesThermoInterpType* stit_ptr) {
203  throw CanteraError("install_STIT", "not implemented");
204  }
205
206  //! Like update(), but only updates the single species k.
207  /*!
208  * @param k species index
209  * @param t Temperature (Kelvin)
210  * @param cp_R Vector of Dimensionless heat capacities. (length m_kk).
211  * @param h_RT Vector of Dimensionless enthalpies. (length m_kk).
212  * @param s_R Vector of Dimensionless entropies. (length m_kk).
213  */
214  virtual void update_one(size_t k, doublereal t, doublereal* cp_R,
215  doublereal* h_RT, doublereal* s_R) const {
216  doublereal tt = 1.e-3*t;
217  m_t[0] = tt;
218  m_t[1] = tt*tt;
219  m_t[2] = m_t[1]*tt;
220  m_t[3] = 1.0/m_t[1];
221  m_t[4] = log(tt);
222  m_t[5] = 1.0/GasConstant;
223  m_t[6] = 1.0/(GasConstant * t);
224
225  size_t grp = m_group_map[k];
226  size_t pos = m_posInGroup_map[k];
227  const std::vector<ShomatePoly> &mlg = m_low[grp-1];
228  const ShomatePoly* nlow = &(mlg[pos]);
229
230  doublereal tmid = nlow->maxTemp();
231  if (t < tmid) {
232  nlow->updateProperties(&m_t[0], cp_R, h_RT, s_R);
233  } else {
234  const std::vector<ShomatePoly> &mhg = m_high[grp-1];
235  const ShomatePoly* nhigh = &(mhg[pos]);
236  nhigh->updateProperties(&m_t[0], cp_R, h_RT, s_R);
237  }
238  }
239
240  virtual void update(doublereal t, doublereal* cp_R,
241  doublereal* h_RT, doublereal* s_R) const {
242  int i;
243
244  doublereal tt = 1.e-3*t;
245  m_t[0] = tt;
246  m_t[1] = tt*tt;
247  m_t[2] = m_t[1]*tt;
248  m_t[3] = 1.0/m_t[1];
249  m_t[4] = log(tt);
250  m_t[5] = 1.0/GasConstant;
251  m_t[6] = 1.0/(GasConstant * t);
252
253  std::vector<ShomatePoly>::const_iterator _begin, _end;
254  for (i = 0; i != m_ngroups; i++) {
255  if (t > m_tmid[i]) {
256  _begin = m_high[i].begin();
257  _end = m_high[i].end();
258  } else {
259  _begin = m_low[i].begin();
260  _end = m_low[i].end();
261  }
262  for (; _begin != _end; ++_begin) {
263  _begin->updateProperties(&m_t[0], cp_R, h_RT, s_R);
264  }
265  }
266  }
267
268  virtual doublereal minTemp(size_t k=npos) const {
269  if (k == npos) {
270  return m_tlow_max;
271  } else {
272  return m_tlow[k];
273  }
274  }
275
276  virtual doublereal maxTemp(size_t k=npos) const {
277  if (k == npos) {
278  return m_thigh_min;
279  } else {
280  return m_thigh[k];
281  }
282  }
283
284  virtual doublereal refPressure(size_t k=npos) const {
285  return m_p0;
286  }
287
288  virtual int reportType(size_t index) const {
289  return SHOMATE;
290  }
291
292  //! @deprecated
293  virtual void reportParams(size_t index, int& type,
294  doublereal* const c,
295  doublereal& minTemp,
296  doublereal& maxTemp,
297  doublereal& refPressure) const {
298  warn_deprecated("ShomateThermo::reportParams");
299  type = reportType(index);
300  if (type == SHOMATE) {
301  size_t grp = m_group_map[index];
302  size_t pos = m_posInGroup_map[index];
303  int itype = SHOMATE;
304  const std::vector<ShomatePoly> &mlg = m_low[grp-1];
305  const std::vector<ShomatePoly> &mhg = m_high[grp-1];
306  const ShomatePoly* lowPoly = &(mlg[pos]);
307  const ShomatePoly* highPoly = &(mhg[pos]);
308  doublereal tmid = lowPoly->maxTemp();
309  c[0] = tmid;
310  size_t n;
311  double ttemp;
312  lowPoly->reportParameters(n, itype, minTemp, ttemp, refPressure,
313  c + 1);
314  if (n != index) {
315  throw CanteraError(" ", "confused");
316  }
317  if (itype != SHOMATE && itype != SHOMATE1) {
318  throw CanteraError(" ", "confused");
319  }
320  highPoly->reportParameters(n, itype, ttemp, maxTemp,
321  refPressure, c + 8);
322  if (n != index) {
323  throw CanteraError(" ", "confused");
324  }
325  if (itype != SHOMATE && itype != SHOMATE1) {
326  throw CanteraError(" ", "confused");
327  }
328  } else {
329  throw CanteraError(" ", "confused");
330  }
331  }
332
333 #ifdef H298MODIFY_CAPABILITY
334
335  virtual doublereal reportOneHf298(int k) const {
336  doublereal h;
337  doublereal t = 298.15;
338
339  int grp = m_group_map[k];
340  int pos = m_posInGroup_map[k];
341  const std::vector<ShomatePoly> &mlg = m_low[grp-1];
342  const ShomatePoly* nlow = &(mlg[pos]);
343
344  doublereal tmid = nlow->maxTemp();
345  if (t <= tmid) {
346  h = nlow->reportHf298();
347  } else {
348  const std::vector<ShomatePoly> &mhg = m_high[grp-1];
349  const ShomatePoly* nhigh = &(mhg[pos]);
350  h = nhigh->reportHf298();
351  }
352  return h;
353  }
354
355  virtual void modifyOneHf298(const int k, const doublereal Hf298New) {
356
357  int grp = m_group_map[k];
358  int pos = m_posInGroup_map[k];
359  std::vector<ShomatePoly> &mlg = m_low[grp-1];
360  ShomatePoly* nlow = &(mlg[pos]);
361  std::vector<ShomatePoly> &mhg = m_high[grp-1];
362  ShomatePoly* nhigh = &(mhg[pos]);
363  doublereal tmid = nlow->maxTemp();
364
365  double hnow = reportOneHf298(k);
366  double delH = Hf298New - hnow;
367  if (298.15 <= tmid) {
368  nlow->modifyOneHf298(k, Hf298New);
369  double h = nhigh->reportHf298(0);
370  double hnew = h + delH;
371  nhigh->modifyOneHf298(k, hnew);
372  } else {
373  nhigh->modifyOneHf298(k, Hf298New);
374  double h = nlow->reportHf298(0);
375  double hnew = h + delH;
376  nlow->modifyOneHf298(k, hnew);
377  }
378
379  }
380
381 #endif
382 protected:
383  //! Vector of vector of NasaPoly1's for the high temp region.
384  /*!
385  * This is the high temp region representation. The first Length is equal
386  * to the number of groups. The second vector is equal to the number of
387  * species in that particular group.
388  */
389  std::vector<std::vector<ShomatePoly> > m_high;
390
391  //! Vector of vector of NasaPoly1's for the low temp region.
392  /*!
393  * This is the low temp region representation. The first Length is equal
394  * to the number of groups. The second vector is equal to the number of
395  * species in that particular group.
396  */
397  std::vector<std::vector<ShomatePoly> > m_low;
398
399  //! Map between the midpoint temperature, as an int, to the group number
400  /*!
401  * Length is equal to the number of groups. Only used in the setup.
402  */
403  std::map<int, int> m_index;
404
405  //! Vector of log temperature limits
406  /*!
407  * Length is equal to the number of groups.
408  */
410
411  //! Maximum value of the low temperature limit
412  doublereal m_tlow_max;
413
414  //! Minimum value of the high temperature limit
415  doublereal m_thigh_min;
416
417  //! Vector of low temperature limits (species index)
418  /*!
419  * Length is equal to number of species
420  */
422
423  //! Vector of low temperature limits (species index)
424  /*!
425  * Length is equal to number of species
426  */
428
429  //! Reference pressure (Pa)
430  /*!
431  * all species must have the same reference pressure.
432  */
433  doublereal m_p0;
434
435  //! number of groups
437
438  //! Vector of temperature polynomials
439  mutable vector_fp m_t;
440
441  /*!
442  * This map takes as its index, the species index in the phase.
443  * It returns the group index, where the temperature polynomials
444  * for that species are stored. group indices start at 1,
445  * so a decrement is always performed to access vectors.
446  */
447  mutable std::map<size_t, size_t> m_group_map;
448
449  /*!
450  * This map takes as its index, the species index in the phase.
451  * It returns the position index within the group, where the
452  * temperature polynomials for that species are stored.
453  */
454  mutable std::map<size_t, size_t> m_posInGroup_map;
455 };
456
457 }
458
459 #endif
vector_fp m_tlow
Vector of low temperature limits (species index)
virtual void install(const std::string &name, size_t index, int type, const doublereal *c, doublereal minTemp, doublereal maxTemp, doublereal refPressure)
Install a new species thermodynamic property parameterization for one species using Shomate polynomia...
virtual int reportType(size_t index) const
This utility function reports the type of parameterization used for the species with index number ind...
Pure Virtual Base class for the thermodynamic manager for an individual species' reference state...
const int ID
Initialized to the type of parameterization.
Definition: ShomateThermo.h:66
#define SHOMATE
Two regions of Shomate Polynomials.
doublereal m_p0
Reference pressure (Pa)
const size_t npos
index returned by functions to indicate "no position"
Definition: ct_defs.h:173
virtual void reportParameters(size_t &n, int &type, doublereal &tlow, doublereal &thigh, doublereal &pref, doublereal *const coeffs) const
Definition: ShomatePoly.h:196
virtual void updateProperties(const doublereal *tt, doublereal *cp_R, doublereal *h_RT, doublereal *s_R) const
Update the properties for this species, given a temperature polynomial.
Definition: ShomatePoly.h:149
std::vector< std::vector< ShomatePoly > > m_low
Vector of vector of NasaPoly1's for the low temp region.
virtual void install_STIT(SpeciesThermoInterpType *stit_ptr)
Install a new species thermodynamic property parameterization for one species.
vector_fp m_t
Vector of temperature polynomials.
doublereal m_tlow_max
Maximum value of the low temperature limit.
void warn_deprecated(const std::string &method, const std::string &extra)
Print a warning indicating that method is deprecated.
Definition: global.cpp:76
virtual void reportParams(size_t index, int &type, doublereal *const c, doublereal &minTemp, doublereal &maxTemp, doublereal &refPressure) const
Pure Virtual base class for the species thermo manager classes.
ShomateThermo(const ShomateThermo &right)
Copy Constructor.
Definition: ShomateThermo.h:82
std::map< size_t, size_t > m_posInGroup_map
std::map< size_t, size_t > m_group_map
Contains const definitions for types of species reference-state thermodynamics managers (see Species ...
virtual doublereal maxTemp() const
Returns the maximum temperature that the thermo parameterization is valid.
std::vector< std::vector< ShomatePoly > > m_high
Vector of vector of NasaPoly1's for the high temp region.
doublereal m_thigh_min
Minimum value of the high temperature limit.
std::map< int, int > m_index
Map between the midpoint temperature, as an int, to the group number.
The Shomate polynomial parameterization for one temperature range for one species.
Definition: ShomatePoly.h:56
Header for a single-species standard state object derived from SpeciesThermoInterpType based on the S...
ShomateThermo & operator=(const ShomateThermo &right)
Assignment Operator.
Definition: ShomateThermo.h:95
std::string fp2str(const double x, const std::string &fmt)
Convert a double into a c++ string.
Definition: stringUtils.cpp:29
virtual void update(doublereal t, doublereal *cp_R, doublereal *h_RT, doublereal *s_R) const
Compute the reference-state properties for all species.
Base class for exceptions thrown by Cantera classes.
Definition: ctexceptions.h:68
vector_fp m_thigh
Vector of low temperature limits (species index)
virtual doublereal minTemp(size_t k=npos) const
Minimum temperature.
virtual doublereal refPressure(size_t k=npos) const
The reference-state pressure for species k.
int m_ngroups
number of groups
This file contains descriptions of templated subclasses of the virtual base class, SpeciesThermo, which includes SpeciesThermoDuo (see Managers for Calculating Reference-State Thermodynamics and class SpeciesThermoDuo)
A species thermodynamic property manager for the Shomate polynomial parameterization.
Definition: ShomateThermo.h:59
std::vector< double > vector_fp
Turn on the use of stl vectors for the basic array type within cantera Vector of doubles.
Definition: ct_defs.h:165
virtual doublereal maxTemp(size_t k=npos) const
Maximum temperature.
virtual SpeciesThermo * duplMyselfAsSpeciesThermo() const
Duplication routine for objects derived from SpeciesThermo.
const doublereal GasConstant
Universal Gas Constant. [J/kmol/K].
Definition: ct_defs.h:66
vector_fp m_tmid
Vector of log temperature limits.
ShomateThermo()
constructor
Definition: ShomateThermo.h:69
void writelog(const std::string &msg)
Write a message to the screen.
Definition: global.cpp:43
virtual void update_one(size_t k, doublereal t, doublereal *cp_R, doublereal *h_RT, doublereal *s_R) const
Like update(), but only updates the single species k.
#define SHOMATE1
one region of Shomate Polynomials used in NIST database This is implemented in the NIST database...