janafThermoI.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 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  const bool convertCoeffs
41 )
42 :
43  EquationOfState(st),
44  Tlow_(Tlow),
45  Thigh_(Thigh),
46  Tcommon_(Tcommon)
47 {
48  if (convertCoeffs)
49  {
50  for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
51  {
52  highCpCoeffs_[coefLabel] = highCpCoeffs[coefLabel]*this->R();
53  lowCpCoeffs_[coefLabel] = lowCpCoeffs[coefLabel]*this->R();
54  }
55  }
56  else
57  {
58  for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
59  {
60  highCpCoeffs_[coefLabel] = highCpCoeffs[coefLabel];
61  lowCpCoeffs_[coefLabel] = lowCpCoeffs[coefLabel];
62  }
63  }
64 }
65 
66 
67 template<class EquationOfState>
70 (
71  const scalar T
72 ) const
73 {
74  if (T < Tcommon_)
75  {
76  return lowCpCoeffs_;
77  }
78  else
79  {
80  return highCpCoeffs_;
81  }
82 }
83 
84 
85 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
86 
87 template<class EquationOfState>
89 (
90  const word& name,
91  const janafThermo& jt
92 )
93 :
94  EquationOfState(name, jt),
95  Tlow_(jt.Tlow_),
96  Thigh_(jt.Thigh_),
97  Tcommon_(jt.Tcommon_)
98 {
99  for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
100  {
101  highCpCoeffs_[coefLabel] = jt.highCpCoeffs_[coefLabel];
102  lowCpCoeffs_[coefLabel] = jt.lowCpCoeffs_[coefLabel];
103  }
104 }
105 
106 
107 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
108 
109 template<class EquationOfState>
111 (
112  const scalar T
113 ) const
114 {
115  if (T < Tlow_ || T > Thigh_)
116  {
118  << "attempt to use janafThermo<EquationOfState>"
119  " out of temperature range "
120  << Tlow_ << " -> " << Thigh_ << "; T = " << T
121  << endl;
122 
123  return min(max(T, Tlow_), Thigh_);
124  }
125  else
126  {
127  return T;
128  }
129 }
130 
131 
132 template<class EquationOfState>
133 inline Foam::scalar Foam::janafThermo<EquationOfState>::Tlow() const
134 {
135  return Tlow_;
136 }
137 
138 
139 template<class EquationOfState>
140 inline Foam::scalar Foam::janafThermo<EquationOfState>::Thigh() const
141 {
142  return Thigh_;
143 }
144 
145 
146 template<class EquationOfState>
148 {
149  return Tcommon_;
150 }
151 
152 
153 template<class EquationOfState>
156 {
157  return highCpCoeffs_;
158 }
159 
160 
161 template<class EquationOfState>
164 {
165  return lowCpCoeffs_;
166 }
167 
168 
169 template<class EquationOfState>
170 inline Foam::scalar Foam::janafThermo<EquationOfState>::Cp
171 (
172  const scalar p,
173  const scalar T
174 ) const
175 {
176  const coeffArray& a = coeffs(T);
177  return
178  ((((a[4]*T + a[3])*T + a[2])*T + a[1])*T + a[0])
179  + EquationOfState::Cp(p, T);
180 }
181 
182 
183 template<class EquationOfState>
184 inline Foam::scalar Foam::janafThermo<EquationOfState>::Hs
185 (
186  const scalar p,
187  const scalar T
188 ) const
189 {
190  return Ha(p, T) - Hf();
191 }
192 
193 
194 template<class EquationOfState>
195 inline Foam::scalar Foam::janafThermo<EquationOfState>::Ha
196 (
197  const scalar p,
198  const scalar T
199 ) const
200 {
201  const coeffArray& a = coeffs(T);
202  return
203  (
204  ((((a[4]/5.0*T + a[3]/4.0)*T + a[2]/3.0)*T + a[1]/2.0)*T + a[0])*T
205  + a[5]
206  ) + EquationOfState::H(p, T);
207 }
208 
209 
210 template<class EquationOfState>
211 inline Foam::scalar Foam::janafThermo<EquationOfState>::Hf() const
212 {
213  const coeffArray& a = lowCpCoeffs_;
214  return
215  (
216  (
217  (((a[4]/5.0*Tstd + a[3]/4.0)*Tstd + a[2]/3.0)*Tstd + a[1]/2.0)*Tstd
218  + a[0]
219  )*Tstd + a[5]
220  );
221 }
222 
223 
224 template<class EquationOfState>
225 inline Foam::scalar Foam::janafThermo<EquationOfState>::S
226 (
227  const scalar p,
228  const scalar T
229 ) const
230 {
231  const coeffArray& a = coeffs(T);
232  return
233  (
234  (((a[4]/4.0*T + a[3]/3.0)*T + a[2]/2.0)*T + a[1])*T + a[0]*log(T)
235  + a[6]
236  ) + EquationOfState::Sp(p, T);
237 }
238 
239 
240 template<class EquationOfState>
242 (
243  const scalar T
244 ) const
245 {
246  const coeffArray& a = coeffs(T);
247  return
248  (
249  (
250  a[0]*(1 - log(T))
251  - (((a[4]/20.0*T + a[3]/12.0)*T + a[2]/6.0)*T + a[1]/2.0)*T
252  - a[6]
253  )*T
254  + a[5]
255  );
256 }
257 
258 
259 template<class EquationOfState>
261 (
262  const scalar p,
263  const scalar T
264 ) const
265 {
266  const coeffArray& a = coeffs(T);
267  return (((4*a[4]*T + 3*a[3])*T + 2*a[2])*T + a[1]);
268 }
269 
270 
271 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
272 
273 template<class EquationOfState>
274 inline void Foam::janafThermo<EquationOfState>::operator+=
275 (
277 )
278 {
279  scalar Y1 = this->Y();
280 
281  EquationOfState::operator+=(jt);
282 
283  if (mag(this->Y()) > small)
284  {
285  Y1 /= this->Y();
286  const scalar Y2 = jt.Y()/this->Y();
287 
288  Tlow_ = max(Tlow_, jt.Tlow_);
289  Thigh_ = min(Thigh_, jt.Thigh_);
290 
291  if
292  (
294  && notEqual(Tcommon_, jt.Tcommon_)
295  )
296  {
298  << "Tcommon " << Tcommon_ << " for "
299  << (this->name().size() ? this->name() : "others")
300  << " != " << jt.Tcommon_ << " for "
301  << (jt.name().size() ? jt.name() : "others")
302  << exit(FatalError);
303  }
304 
305  for
306  (
307  label coefLabel=0;
308  coefLabel<janafThermo<EquationOfState>::nCoeffs_;
309  coefLabel++
310  )
311  {
312  highCpCoeffs_[coefLabel] =
313  Y1*highCpCoeffs_[coefLabel]
314  + Y2*jt.highCpCoeffs_[coefLabel];
315 
316  lowCpCoeffs_[coefLabel] =
317  Y1*lowCpCoeffs_[coefLabel]
318  + Y2*jt.lowCpCoeffs_[coefLabel];
319  }
320  }
321 }
322 
323 
324 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
325 
326 template<class EquationOfState>
327 inline Foam::janafThermo<EquationOfState> Foam::operator+
328 (
329  const janafThermo<EquationOfState>& jt1,
331 )
332 {
333  EquationOfState eofs = jt1;
334  eofs += jt2;
335 
336  if (mag(eofs.Y()) < small)
337  {
339  (
340  eofs,
341  jt1.Tlow_,
342  jt1.Thigh_,
343  jt1.Tcommon_,
344  jt1.highCpCoeffs_,
345  jt1.lowCpCoeffs_
346  );
347  }
348  else
349  {
350  const scalar Y1 = jt1.Y()/eofs.Y();
351  const scalar Y2 = jt2.Y()/eofs.Y();
352 
355 
356  for
357  (
358  label coefLabel=0;
359  coefLabel<janafThermo<EquationOfState>::nCoeffs_;
360  coefLabel++
361  )
362  {
363  highCpCoeffs[coefLabel] =
364  Y1*jt1.highCpCoeffs_[coefLabel]
365  + Y2*jt2.highCpCoeffs_[coefLabel];
366 
367  lowCpCoeffs[coefLabel] =
368  Y1*jt1.lowCpCoeffs_[coefLabel]
369  + Y2*jt2.lowCpCoeffs_[coefLabel];
370  }
371 
372  if
373  (
375  && notEqual(jt1.Tcommon_, jt2.Tcommon_)
376  )
377  {
379  << "Tcommon " << jt1.Tcommon_ << " for "
380  << (jt1.name().size() ? jt1.name() : "others")
381  << " != " << jt2.Tcommon_ << " for "
382  << (jt2.name().size() ? jt2.name() : "others")
383  << exit(FatalError);
384  }
385 
387  (
388  eofs,
389  max(jt1.Tlow_, jt2.Tlow_),
390  min(jt1.Thigh_, jt2.Thigh_),
391  jt1.Tcommon_,
392  highCpCoeffs,
394  );
395  }
396 }
397 
398 
399 template<class EquationOfState>
400 inline Foam::janafThermo<EquationOfState> Foam::operator*
401 (
402  const scalar s,
404 )
405 {
407  (
408  s*static_cast<const EquationOfState&>(jt),
409  jt.Tlow_,
410  jt.Thigh_,
411  jt.Tcommon_,
412  jt.highCpCoeffs_,
413  jt.lowCpCoeffs_
414  );
415 }
416 
417 
418 template<class EquationOfState>
419 inline Foam::janafThermo<EquationOfState> Foam::operator==
420 (
421  const janafThermo<EquationOfState>& jt1,
423 )
424 {
425  EquationOfState eofs
426  (
427  static_cast<const EquationOfState&>(jt1)
428  == static_cast<const EquationOfState&>(jt2)
429  );
430 
431  const scalar Y1 = jt2.Y()/eofs.Y();
432  const scalar Y2 = jt1.Y()/eofs.Y();
433 
436 
437  for
438  (
439  label coefLabel=0;
440  coefLabel<janafThermo<EquationOfState>::nCoeffs_;
441  coefLabel++
442  )
443  {
444  highCpCoeffs[coefLabel] =
445  Y1*jt2.highCpCoeffs_[coefLabel]
446  - Y2*jt1.highCpCoeffs_[coefLabel];
447 
448  lowCpCoeffs[coefLabel] =
449  Y1*jt2.lowCpCoeffs_[coefLabel]
450  - Y2*jt1.lowCpCoeffs_[coefLabel];
451  }
452 
453  if
454  (
456  && notEqual(jt2.Tcommon_, jt1.Tcommon_)
457  )
458  {
460  << "Tcommon " << jt2.Tcommon_ << " for "
461  << (jt2.name().size() ? jt2.name() : "others")
462  << " != " << jt1.Tcommon_ << " for "
463  << (jt1.name().size() ? jt1.name() : "others")
464  << exit(FatalError);
465  }
466 
468  (
469  eofs,
470  max(jt2.Tlow_, jt1.Tlow_),
471  min(jt2.Thigh_, jt1.Thigh_),
472  jt2.Tcommon_,
473  highCpCoeffs,
475  );
476 }
477 
478 
479 // ************************************************************************* //
bool notEqual(const Scalar s1, const Scalar s2)
Definition: Scalar.H:215
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
FvWallInfoData< WallInfo, label > label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
static const int nCoeffs_
Definition: janafThermo.H:158
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
const dimensionedScalar Tstd
Standard temperature.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:306
scalar Thigh() const
Return const access to the high temperature limit.
Definition: janafThermoI.H:140
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
const coeffArray & highCpCoeffs() const
Return const access to the high temperature poly coefficients.
Definition: janafThermoI.H:155
scalar Gstd(const scalar T) const
Gibbs free energy of the mixture in the standard state [J/kg].
Definition: janafThermoI.H:242
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))
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
Definition: janafThermoI.H:261
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
Definition: janafThermoI.H:111
A class for handling words, derived from string.
Definition: word.H:59
scalar S(const scalar p, const scalar T) const
Entropy [J/kg/K].
Definition: janafThermoI.H:226
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
scalar Tlow() const
Return const access to the low temperature limit.
Definition: janafThermoI.H:133
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
Definition: janafThermoI.H:185
janafThermo(const EquationOfState &st, const scalar Tlow, const scalar Thigh, const scalar Tcommon, const coeffArray &highCpCoeffs, const coeffArray &lowCpCoeffs, const bool convertCoeffs=false)
Construct from components.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
const tmp< volScalarField::Internal > & Sp
Definition: alphaSuSp.H:7
scalar Hf() const
Enthalpy of formation [J/kg].
Definition: janafThermoI.H:211
#define R(A, B, C, D, E, F, K, M)
#define WarningInFunction
Report a warning using Foam::Warning.
PtrList< volScalarField > & Y
const coeffArray & lowCpCoeffs() const
Return const access to the low temperature poly coefficients.
Definition: janafThermoI.H:163
Enthalpy based thermodynamics package using JANAF tables:
Definition: janafThermo.H:113
scalar Cp(const scalar p, const scalar T) const
Definition: EtoHthermo.H:2
scalar Ha(const scalar p, const scalar T) const
Absolute enthalpy [J/kg].
Definition: janafThermoI.H:196
dimensioned< scalar > mag(const dimensioned< Type > &)
scalar Cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/kg/K].
Definition: janafThermoI.H:171
scalar Tcommon() const
Return const access to the common temperature.
Definition: janafThermoI.H:147