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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2020 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "janafThermo.H"
30#include "specie.H"
31
32// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33
34template<class EquationOfState>
36(
37 const EquationOfState& st,
38 const scalar Tlow,
39 const scalar Thigh,
40 const scalar Tcommon,
41 const typename janafThermo<EquationOfState>::coeffArray& highCpCoeffs,
42 const typename janafThermo<EquationOfState>::coeffArray& lowCpCoeffs,
43 const bool convertCoeffs
44)
45:
46 EquationOfState(st),
47 Tlow_(Tlow),
48 Thigh_(Thigh),
49 Tcommon_(Tcommon)
50{
51 if (convertCoeffs)
52 {
53 for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
54 {
55 highCpCoeffs_[coefLabel] = highCpCoeffs[coefLabel]*this->R();
56 lowCpCoeffs_[coefLabel] = lowCpCoeffs[coefLabel]*this->R();
57 }
58 }
59 else
60 {
61 for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
62 {
63 highCpCoeffs_[coefLabel] = highCpCoeffs[coefLabel];
64 lowCpCoeffs_[coefLabel] = lowCpCoeffs[coefLabel];
65 }
66 }
67}
68
69
70template<class EquationOfState>
73(
74 const scalar T
75) const
76{
77 if (T < Tcommon_)
78 {
79 return lowCpCoeffs_;
80 }
81 else
82 {
83 return highCpCoeffs_;
84 }
85}
86
87
88// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
89
90template<class EquationOfState>
92(
93 const word& name,
94 const janafThermo& jt
95)
96:
97 EquationOfState(name, jt),
98 Tlow_(jt.Tlow_),
99 Thigh_(jt.Thigh_),
100 Tcommon_(jt.Tcommon_)
101{
102 for (label coefLabel=0; coefLabel<nCoeffs_; coefLabel++)
103 {
104 highCpCoeffs_[coefLabel] = jt.highCpCoeffs_[coefLabel];
105 lowCpCoeffs_[coefLabel] = jt.lowCpCoeffs_[coefLabel];
106 }
107}
108
109
110// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
111
112template<class EquationOfState>
114(
115 const scalar T
116) const
117{
118 if (T < Tlow_ || T > Thigh_)
119 {
121 << "attempt to use janafThermo<EquationOfState>"
122 " out of temperature range "
123 << Tlow_ << " -> " << Thigh_ << "; T = " << T
124 << endl;
125
126 return min(max(T, Tlow_), Thigh_);
127 }
128 else
129 {
130 return T;
131 }
132}
133
134
135template<class EquationOfState>
137{
138 return Tlow_;
139}
140
141
142template<class EquationOfState>
144{
145 return Thigh_;
146}
147
148
149template<class EquationOfState>
151{
152 return Tcommon_;
153}
154
155
156template<class EquationOfState>
159{
160 return highCpCoeffs_;
161}
162
163
164template<class EquationOfState>
167{
168 return lowCpCoeffs_;
169}
170
171
172template<class EquationOfState>
174(
175 const scalar p,
176 const scalar T
177) const
178{
179 const coeffArray& a = coeffs(T);
180 return
181 ((((a[4]*T + a[3])*T + a[2])*T + a[1])*T + a[0])
182 + EquationOfState::Cp(p, T);
183}
184
185
186template<class EquationOfState>
188(
189 const scalar p,
190 const scalar T
191) const
192{
193 const coeffArray& a = coeffs(T);
194 return
195 (
196 ((((a[4]/5.0*T + a[3]/4.0)*T + a[2]/3.0)*T + a[1]/2.0)*T + a[0])*T
197 + a[5]
198 ) + EquationOfState::H(p, T);
199}
200
201
202template<class EquationOfState>
204(
205 const scalar p,
206 const scalar T
207) const
208{
209 return Ha(p, T) - Hc();
210}
211
212
213template<class EquationOfState>
214inline Foam::scalar Foam::janafThermo<EquationOfState>::Hc() const
215{
216 const coeffArray& a = lowCpCoeffs_;
217 return
218 (
219 (
220 (((a[4]/5.0*Tstd + a[3]/4.0)*Tstd + a[2]/3.0)*Tstd + a[1]/2.0)*Tstd
221 + a[0]
222 )*Tstd + a[5]
223 );
224}
225
226
227template<class EquationOfState>
229(
230 const scalar p,
231 const scalar T
232) const
233{
234 const coeffArray& a = coeffs(T);
235 return
236 (
237 (((a[4]/4.0*T + a[3]/3.0)*T + a[2]/2.0)*T + a[1])*T + a[0]*log(T)
238 + a[6]
239 ) + EquationOfState::S(p, T);
240}
241
242
243template<class EquationOfState>
245(
246 const scalar T
247) const
248{
249 const coeffArray& a = coeffs(T);
250 return
251 (
252 (
253 a[0]*(1 - log(T))
254 - (((a[4]/20.0*T + a[3]/12.0)*T + a[2]/6.0)*T + a[1]/2.0)*T
255 - a[6]
256 )*T
257 + a[5]
258 );
259}
260
261
262template<class EquationOfState>
264(
265 const scalar p,
266 const scalar T
267) const
268{
269 const coeffArray& a = coeffs(T);
270 return
271 (((4*a[4]*T + 3*a[3])*T + 2*a[2])*T + a[1]);
272}
273
274// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
275
276template<class EquationOfState>
278(
280)
281{
282 scalar Y1 = this->Y();
283
284 EquationOfState::operator+=(jt);
285
286 if (mag(this->Y()) > SMALL)
287 {
288 Y1 /= this->Y();
289 const scalar Y2 = jt.Y()/this->Y();
290
291 Tlow_ = max(Tlow_, jt.Tlow_);
292 Thigh_ = min(Thigh_, jt.Thigh_);
293
294 if
295 (
297 && notEqual(Tcommon_, jt.Tcommon_)
298 )
299 {
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 Y1*highCpCoeffs_[coefLabel]
317 + Y2*jt.highCpCoeffs_[coefLabel];
318
319 lowCpCoeffs_[coefLabel] =
320 Y1*lowCpCoeffs_[coefLabel]
321 + Y2*jt.lowCpCoeffs_[coefLabel];
322 }
323 }
324}
325
326
327// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
328
329template<class EquationOfState>
330inline Foam::janafThermo<EquationOfState> Foam::operator+
331(
334)
335{
336 EquationOfState eofs = jt1;
337 eofs += jt2;
338
339 if (mag(eofs.Y()) < SMALL)
340 {
342 (
343 eofs,
344 jt1.Tlow_,
345 jt1.Thigh_,
346 jt1.Tcommon_,
347 jt1.highCpCoeffs_,
348 jt1.lowCpCoeffs_
349 );
350 }
351 else
352 {
353 const scalar Y1 = jt1.Y()/eofs.Y();
354 const scalar Y2 = jt2.Y()/eofs.Y();
355
356 typename janafThermo<EquationOfState>::coeffArray highCpCoeffs;
357 typename janafThermo<EquationOfState>::coeffArray lowCpCoeffs;
358
359 for
360 (
361 label coefLabel=0;
362 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
363 coefLabel++
364 )
365 {
366 highCpCoeffs[coefLabel] =
367 Y1*jt1.highCpCoeffs_[coefLabel]
368 + Y2*jt2.highCpCoeffs_[coefLabel];
369
370 lowCpCoeffs[coefLabel] =
371 Y1*jt1.lowCpCoeffs_[coefLabel]
372 + Y2*jt2.lowCpCoeffs_[coefLabel];
373 }
374
375 if
376 (
377 janafThermo<EquationOfState>::debug
378 && notEqual(jt1.Tcommon_, jt2.Tcommon_)
379 )
380 {
382 << "Tcommon " << jt1.Tcommon_ << " for "
383 << (jt1.name().size() ? jt1.name() : "others")
384 << " != " << jt2.Tcommon_ << " for "
385 << (jt2.name().size() ? jt2.name() : "others")
386 << exit(FatalError);
387 }
388
389 return janafThermo<EquationOfState>
390 (
391 eofs,
392 max(jt1.Tlow_, jt2.Tlow_),
393 min(jt1.Thigh_, jt2.Thigh_),
394 jt1.Tcommon_,
395 highCpCoeffs,
396 lowCpCoeffs
397 );
398 }
399}
400
401
402template<class EquationOfState>
403inline Foam::janafThermo<EquationOfState> Foam::operator*
404(
405 const scalar s,
406 const janafThermo<EquationOfState>& jt
407)
408{
409 return janafThermo<EquationOfState>
410 (
411 s*static_cast<const EquationOfState&>(jt),
412 jt.Tlow_,
413 jt.Thigh_,
414 jt.Tcommon_,
415 jt.highCpCoeffs_,
416 jt.lowCpCoeffs_
417 );
418}
419
420
421template<class EquationOfState>
422inline Foam::janafThermo<EquationOfState> Foam::operator==
423(
424 const janafThermo<EquationOfState>& jt1,
425 const janafThermo<EquationOfState>& jt2
426)
427{
428 EquationOfState eofs
429 (
430 static_cast<const EquationOfState&>(jt1)
431 == static_cast<const EquationOfState&>(jt2)
432 );
433
434 const scalar Y1 = jt2.Y()/eofs.Y();
435 const scalar Y2 = jt1.Y()/eofs.Y();
436
437 typename janafThermo<EquationOfState>::coeffArray highCpCoeffs;
438 typename janafThermo<EquationOfState>::coeffArray lowCpCoeffs;
439
440 for
441 (
442 label coefLabel=0;
443 coefLabel<janafThermo<EquationOfState>::nCoeffs_;
444 coefLabel++
445 )
446 {
447 highCpCoeffs[coefLabel] =
448 Y1*jt2.highCpCoeffs_[coefLabel]
449 - Y2*jt1.highCpCoeffs_[coefLabel];
450
451 lowCpCoeffs[coefLabel] =
452 Y1*jt2.lowCpCoeffs_[coefLabel]
453 - Y2*jt1.lowCpCoeffs_[coefLabel];
454 }
455
456 if
457 (
458 janafThermo<EquationOfState>::debug
459 && notEqual(jt2.Tcommon_, jt1.Tcommon_)
460 )
461 {
463 << "Tcommon " << jt2.Tcommon_ << " for "
464 << (jt2.name().size() ? jt2.name() : "others")
465 << " != " << jt1.Tcommon_ << " for "
466 << (jt1.name().size() ? jt1.name() : "others")
467 << exit(FatalError);
468 }
469
470 return janafThermo<EquationOfState>
471 (
472 eofs,
473 max(jt2.Tlow_, jt1.Tlow_),
474 min(jt2.Thigh_, jt1.Thigh_),
475 jt2.Tcommon_,
476 highCpCoeffs,
477 lowCpCoeffs
478 );
479}
480
481
482// ************************************************************************* //
scalar Ha(const scalar p, const scalar T) const
Definition: EtoHthermo.H:32
#define R(A, B, C, D, E, F, K, M)
compactSpatialTensor S
The joint motion sub-space (3-DoF)
Definition: joint.H:132
const dictionary & coeffs() const
Return const dictionary of the model.
friend complex limit(const complex &c1, const complex &c2)
Definition: complexI.H:263
JANAF tables based thermodynamics package templated into the equation of state.
Definition: janafThermo.H:94
static constexpr int nCoeffs_
Definition: janafThermo.H:99
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
Definition: janafThermoI.H:204
scalar Tlow() const
Return const access to the low temperature limit.
Definition: janafThermoI.H:136
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
Definition: janafThermoI.H:264
FixedList< scalar, nCoeffs_ > coeffArray
Definition: janafThermo.H:100
scalar Hc() const
Chemical enthalpy [J/kg].
Definition: janafThermoI.H:214
scalar Tcommon() const
Return const access to the common temperature.
Definition: janafThermoI.H:150
scalar Ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kg].
Definition: janafThermoI.H:188
scalar Thigh() const
Return const access to the high temperature limit.
Definition: janafThermoI.H:143
const coeffArray & highCpCoeffs() const
Return const access to the high temperature poly coefficients.
Definition: janafThermoI.H:158
scalar Gstd(const scalar T) const
Gibbs free energy of the mixture in the standard state [J/kg].
Definition: janafThermoI.H:245
const coeffArray & lowCpCoeffs() const
Return const access to the low temperature poly coefficients.
Definition: janafThermoI.H:166
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
PtrList< volScalarField > & Y
const volScalarField & T
const volScalarField & Cp
Definition: EEqn.H:7
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.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))
#define WarningInFunction
Report a warning using Foam::Warning.
const scalar Tstd
Standard temperature: default in [K].
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
bool notEqual(const Scalar s1, const Scalar s2)
Definition: Scalar.H:291
dimensionedScalar log(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
error FatalError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130