janafThermoI.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2015 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "janafThermo.H"
27 #include "specie.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 template<class EquationOfState>
33 (
34  const EquationOfState& st,
35  const scalar Tlow,
36  const scalar Thigh,
37  const scalar Tcommon,
38  const typename janafThermo<EquationOfState>::coeffArray& highCpCoeffs,
39  const typename janafThermo<EquationOfState>::coeffArray& lowCpCoeffs
40 )
41 :
42  EquationOfState(st),
43  Tlow_(Tlow),
44  Thigh_(Thigh),
45  Tcommon_(Tcommon)
46 {
47  for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
48  {
49  highCpCoeffs_[coefLabel] = highCpCoeffs[coefLabel];
50  lowCpCoeffs_[coefLabel] = lowCpCoeffs[coefLabel];
51  }
52 }
53 
54 
55 template<class EquationOfState>
58 (
59  const scalar T
60 ) const
61 {
62  if (T < Tcommon_)
63  {
64  return lowCpCoeffs_;
65  }
66  else
67  {
68  return highCpCoeffs_;
69  }
70 }
71 
72 
73 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
74 
75 template<class EquationOfState>
77 (
78  const word& name,
79  const janafThermo& jt
80 )
81 :
82  EquationOfState(name, jt),
83  Tlow_(jt.Tlow_),
84  Thigh_(jt.Thigh_),
85  Tcommon_(jt.Tcommon_)
86 {
87  for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
88  {
89  highCpCoeffs_[coefLabel] = jt.highCpCoeffs_[coefLabel];
90  lowCpCoeffs_[coefLabel] = jt.lowCpCoeffs_[coefLabel];
91  }
92 }
93 
94 
95 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
96 
97 template<class EquationOfState>
99 (
100  const scalar T
101 ) const
102 {
103  if (T < Tlow_ || T > Thigh_)
104  {
105  WarningIn
106  (
107  "janafThermo<EquationOfState>::limit(const scalar T) const"
108  ) << "attempt to use janafThermo<EquationOfState>"
109  " out of temperature range "
110  << Tlow_ << " -> " << Thigh_ << "; T = " << T
111  << endl;
112 
113  return min(max(T, Tlow_), Thigh_);
114  }
115  else
116  {
117  return T;
118  }
119 }
120 
121 
122 template<class EquationOfState>
123 inline Foam::scalar Foam::janafThermo<EquationOfState>::Tlow() const
124 {
125  return Tlow_;
126 }
127 
128 
129 template<class EquationOfState>
130 inline Foam::scalar Foam::janafThermo<EquationOfState>::Thigh() const
131 {
132  return Thigh_;
133 }
134 
135 
136 template<class EquationOfState>
138 {
139  return Tcommon_;
140 }
141 
142 
143 template<class EquationOfState>
146 {
147  return highCpCoeffs_;
148 }
149 
150 
151 template<class EquationOfState>
154 {
155  return lowCpCoeffs_;
156 }
157 
158 
159 template<class EquationOfState>
160 inline Foam::scalar Foam::janafThermo<EquationOfState>::cp
161 (
162  const scalar p,
163  const scalar T
164 ) const
165 {
166  const coeffArray& a = coeffs(T);
167  return RR*((((a[4]*T + a[3])*T + a[2])*T + a[1])*T + a[0]);
168 }
169 
170 
171 template<class EquationOfState>
172 inline Foam::scalar Foam::janafThermo<EquationOfState>::ha
173 (
174  const scalar p,
175  const scalar T
176 ) const
177 {
178  const coeffArray& a = coeffs(T);
179  return RR*
180  (
181  ((((a[4]/5.0*T + a[3]/4.0)*T + a[2]/3.0)*T + a[1]/2.0)*T + a[0])*T
182  + a[5]
183  );
184 }
185 
186 
187 template<class EquationOfState>
188 inline Foam::scalar Foam::janafThermo<EquationOfState>::hs
189 (
190  const scalar p,
191  const scalar T
192 ) const
193 {
194  return ha(p, T) - hc();
195 }
196 
197 
198 template<class EquationOfState>
199 inline Foam::scalar Foam::janafThermo<EquationOfState>::hc() const
200 {
201  const coeffArray& a = lowCpCoeffs_;
202  return RR*
203  (
204  (
205  (((a[4]/5.0*Tstd + a[3]/4.0)*Tstd + a[2]/3.0)*Tstd + a[1]/2.0)*Tstd
206  + a[0]
207  )*Tstd + a[5]
208  );
209 }
210 
211 
212 template<class EquationOfState>
213 inline Foam::scalar Foam::janafThermo<EquationOfState>::s
214 (
215  const scalar p,
216  const scalar T
217 ) const
218 {
219  const coeffArray& a = coeffs(T);
220  return
221  RR*
222  (
223  (((a[4]/4.0*T + a[3]/3.0)*T + a[2]/2.0)*T + a[1])*T + a[0]*log(T)
224  + a[6]
225  )
226  + EquationOfState::s(p, T);
227 }
228 
229 
230 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
231 
232 template<class EquationOfState>
233 inline void Foam::janafThermo<EquationOfState>::operator+=
234 (
236 )
237 {
238  scalar molr1 = this->nMoles();
239 
240  EquationOfState::operator+=(jt);
241 
242  molr1 /= this->nMoles();
243  scalar molr2 = jt.nMoles()/this->nMoles();
244 
245  Tlow_ = max(Tlow_, jt.Tlow_);
246  Thigh_ = min(Thigh_, jt.Thigh_);
247 
248  if (janafThermo<EquationOfState>::debug && notEqual(Tcommon_, jt.Tcommon_))
249  {
251  (
252  "janafThermo<EquationOfState>::operator+="
253  "(const janafThermo<EquationOfState>& jt) const"
254  ) << "Tcommon " << Tcommon_ << " for "
255  << (this->name().size() ? this->name() : "others")
256  << " != " << jt.Tcommon_ << " for "
257  << (jt.name().size() ? jt.name() : "others")
258  << exit(FatalError);
259  }
260 
261  for
262  (
263  label coefLabel=0;
264  coefLabel<janafThermo<EquationOfState>::nCoeffs_;
265  coefLabel++
266  )
267  {
268  highCpCoeffs_[coefLabel] =
269  molr1*highCpCoeffs_[coefLabel]
270  + molr2*jt.highCpCoeffs_[coefLabel];
271 
272  lowCpCoeffs_[coefLabel] =
273  molr1*lowCpCoeffs_[coefLabel]
274  + molr2*jt.lowCpCoeffs_[coefLabel];
275  }
276 }
277 
278 
279 template<class EquationOfState>
280 inline void Foam::janafThermo<EquationOfState>::operator-=
281 (
283 )
284 {
285  scalar molr1 = this->nMoles();
286 
287  EquationOfState::operator-=(jt);
288 
289  molr1 /= this->nMoles();
290  scalar molr2 = jt.nMoles()/this->nMoles();
291 
292  Tlow_ = max(Tlow_, jt.Tlow_);
293  Thigh_ = min(Thigh_, jt.Thigh_);
294 
295  if (janafThermo<EquationOfState>::debug && notEqual(Tcommon_, jt.Tcommon_))
296  {
298  (
299  "janafThermo<EquationOfState>::operator-="
300  "(const janafThermo<EquationOfState>& jt) const"
301  ) << "Tcommon " << Tcommon_ << " for "
302  << (this->name().size() ? this->name() : "others")
303  << " != " << jt.Tcommon_ << " for "
304  << (jt.name().size() ? jt.name() : "others")
305  << exit(FatalError);
306  }
307 
308  for
309  (
310  label coefLabel=0;
311  coefLabel<janafThermo<EquationOfState>::nCoeffs_;
312  coefLabel++
313  )
314  {
315  highCpCoeffs_[coefLabel] =
316  molr1*highCpCoeffs_[coefLabel]
317  - molr2*jt.highCpCoeffs_[coefLabel];
318 
319  lowCpCoeffs_[coefLabel] =
320  molr1*lowCpCoeffs_[coefLabel]
321  - molr2*jt.lowCpCoeffs_[coefLabel];
322  }
323 }
324 
325 
326 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
327 
328 template<class EquationOfState>
329 inline Foam::janafThermo<EquationOfState> Foam::operator+
330 (
331  const janafThermo<EquationOfState>& jt1,
333 )
334 {
335  EquationOfState eofs = jt1;
336  eofs += jt2;
337 
338  scalar molr1 = jt1.nMoles()/eofs.nMoles();
339  scalar molr2 = jt2.nMoles()/eofs.nMoles();
340 
343 
344  for
345  (
346  label coefLabel=0;
347  coefLabel<janafThermo<EquationOfState>::nCoeffs_;
348  coefLabel++
349  )
350  {
351  highCpCoeffs[coefLabel] =
352  molr1*jt1.highCpCoeffs_[coefLabel]
353  + molr2*jt2.highCpCoeffs_[coefLabel];
354 
355  lowCpCoeffs[coefLabel] =
356  molr1*jt1.lowCpCoeffs_[coefLabel]
357  + molr2*jt2.lowCpCoeffs_[coefLabel];
358  }
359 
360  if
361  (
363  && notEqual(jt1.Tcommon_, jt2.Tcommon_)
364  )
365  {
367  (
368  "operator+"
369  "(const janafThermo<EquationOfState>& jt1,"
370  " const janafThermo<EquationOfState>& jt2)"
371  ) << "Tcommon " << jt1.Tcommon_ << " for "
372  << (jt1.name().size() ? jt1.name() : "others")
373  << " != " << jt2.Tcommon_ << " for "
374  << (jt2.name().size() ? jt2.name() : "others")
375  << exit(FatalError);
376  }
377 
379  (
380  eofs,
381  max(jt1.Tlow_, jt2.Tlow_),
382  min(jt1.Thigh_, jt2.Thigh_),
383  jt1.Tcommon_,
384  highCpCoeffs,
386  );
387 }
388 
389 
390 template<class EquationOfState>
391 inline Foam::janafThermo<EquationOfState> Foam::operator-
392 (
393  const janafThermo<EquationOfState>& jt1,
395 )
396 {
397  EquationOfState eofs = jt1;
398  eofs -= jt2;
399 
400  scalar molr1 = jt1.nMoles()/eofs.nMoles();
401  scalar molr2 = jt2.nMoles()/eofs.nMoles();
402 
405 
406  for
407  (
408  label coefLabel=0;
409  coefLabel<janafThermo<EquationOfState>::nCoeffs_;
410  coefLabel++
411  )
412  {
413  highCpCoeffs[coefLabel] =
414  molr1*jt1.highCpCoeffs_[coefLabel]
415  - molr2*jt2.highCpCoeffs_[coefLabel];
416 
417  lowCpCoeffs[coefLabel] =
418  molr1*jt1.lowCpCoeffs_[coefLabel]
419  - molr2*jt2.lowCpCoeffs_[coefLabel];
420  }
421 
422  if
423  (
425  && notEqual(jt1.Tcommon_, jt2.Tcommon_)
426  )
427  {
429  (
430  "operator-"
431  "(const janafThermo<EquationOfState>& jt1,"
432  " const janafThermo<EquationOfState>& jt2)"
433  ) << "Tcommon " << jt1.Tcommon_ << " for "
434  << (jt1.name().size() ? jt1.name() : "others")
435  << " != " << jt2.Tcommon_ << " for "
436  << (jt2.name().size() ? jt2.name() : "others")
437  << exit(FatalError);
438  }
439 
441  (
442  eofs,
443  max(jt1.Tlow_, jt2.Tlow_),
444  min(jt1.Thigh_, jt2.Thigh_),
445  jt1.Tcommon_,
446  highCpCoeffs,
448  );
449 }
450 
451 
452 template<class EquationOfState>
453 inline Foam::janafThermo<EquationOfState> Foam::operator*
454 (
455  const scalar s,
457 )
458 {
460  (
461  s*static_cast<const EquationOfState&>(jt),
462  jt.Tlow_,
463  jt.Thigh_,
464  jt.Tcommon_,
465  jt.highCpCoeffs_,
466  jt.lowCpCoeffs_
467  );
468 }
469 
470 
471 template<class EquationOfState>
472 inline Foam::janafThermo<EquationOfState> Foam::operator==
473 (
474  const janafThermo<EquationOfState>& jt1,
476 )
477 {
478  return jt2 - jt1;
479 }
480 
481 
482 // ************************************************************************* //
gmvFile<< "tracers "<< particles.size()<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().x()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().y()<< " ";}gmvFile<< nl;forAllConstIter(Cloud< passiveParticle >, particles, iter){gmvFile<< iter().position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject( name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE ))
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
A class for handling words, derived from string.
Definition: word.H:59
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
scalar Thigh() const
Return const access to the high temperature limit.
Definition: janafThermoI.H:130
scalar cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kmol K)].
Definition: janafThermoI.H:161
scalar Tcommon() const
Return const access to the common temperature.
Definition: janafThermoI.H:137
dimensionedScalar log(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
#define WarningIn(functionName)
Report a warning using Foam::Warning.
scalar hc() const
Chemical enthalpy [J/kmol].
Definition: janafThermoI.H:199
const dimensionedScalar Tstd
Standard temperature.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
scalar s(const scalar p, const scalar T) const
Entropy [J/(kmol K)].
Definition: janafThermoI.H:214
janafThermo(const EquationOfState &st, const scalar Tlow, const scalar Thigh, const scalar Tcommon, const coeffArray &highCpCoeffs, const coeffArray &lowCpCoeffs)
Construct from components.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
const coeffArray & lowCpCoeffs() const
Return const access to the low temperature poly coefficients.
Definition: janafThermoI.H:153
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
Definition: janafThermoI.H:99
const scalar RR
Universal gas constant (default in [J/(kmol K)])
bool notEqual(const Scalar s1, const Scalar s2)
Definition: Scalar.H:152
error FatalError
scalar Tlow() const
Return const access to the low temperature limit.
Definition: janafThermoI.H:123
const coeffArray & highCpCoeffs() const
Return const access to the high temperature poly coefficients.
Definition: janafThermoI.H:145
static const int nCoeffs_
Definition: janafThermo.H:101
scalar hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kmol].
Definition: janafThermoI.H:189
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
JANAF tables based thermodynamics package templated into the equation of state.
Definition: janafThermo.H:49
scalar ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kmol].
Definition: janafThermoI.H:173