thermo.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-2023 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 Class
25  Foam::species::thermo
26 
27 Description
28  Basic thermodynamics type based on the use of fitting functions for
29  cp, h, s obtained from the template argument type thermo. All other
30  properties are derived from these primitive functions.
31 
32 SourceFiles
33  thermoI.H
34  thermo.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef thermo_H
39 #define thermo_H
40 
41 #include "thermodynamicConstants.H"
42 using namespace Foam::constant::thermodynamic;
43 
44 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
45 
46 namespace Foam
47 {
48 namespace species
49 {
50 
51 // Forward declaration of friend functions and operators
52 
53 template<class Thermo, template<class> class Type> class thermo;
54 
55 template<class Thermo, template<class> class Type>
56 inline thermo<Thermo, Type> operator+
57 (
58  const thermo<Thermo, Type>&,
60 );
61 
62 template<class Thermo, template<class> class Type>
63 inline thermo<Thermo, Type> operator*
64 (
65  const scalar,
67 );
68 
69 template<class Thermo, template<class> class Type>
70 inline thermo<Thermo, Type> operator==
71 (
72  const thermo<Thermo, Type>&,
74 );
75 
76 template<class Thermo, template<class> class Type>
77 Ostream& operator<<
78 (
79  Ostream&,
81 );
82 
83 
84 /*---------------------------------------------------------------------------*\
85  Class thermo Declaration
86 \*---------------------------------------------------------------------------*/
87 
88 template<class Thermo, template<class> class Type>
89 class thermo
90 :
91  public Thermo,
92  public Type<thermo<Thermo, Type>>
93 {
94  // Private Data
95 
96  //- Convergence tolerance of energy -> temperature inversion functions
97  static const scalar tol_;
98 
99  //- Max number of iterations in energy->temperature inversion functions
100  static const int maxIter_;
101 
102 
103 public:
104 
105  // Public Typedefs
106 
107  //- The thermodynamics of the individual species'
109 
110 
111  // Constructors
112 
113  //- Construct from components
114  inline thermo(const Thermo& sp);
115 
116  //- Construct from name and dictionary
117  thermo(const word& name, const dictionary& dict);
118 
119  //- Construct as named copy
120  inline thermo(const word& name, const thermo&);
121 
122 
123  // Member Functions
124 
125  //- Return the instantiated type name
126  static word typeName()
127  {
128  return
129  Thermo::typeName() + ','
130  + Type<thermo<Thermo, Type>>::typeName();
131  }
132 
133  //- Return true if energy type is enthalpy
134  static inline bool enthalpy();
135 
136  //- Name of Enthalpy/Internal energy
137  static inline word heName();
138 
139 
140  // Fundamental properties
141  // (These functions must be provided in derived types)
142 
143  // Heat capacity at constant pressure [J/kg/K]
144  // inline scalar Cp(const scalar p, const scalar T) const;
145 
146  // Sensible enthalpy [J/kg]
147  // inline scalar hs(const scalar p, const scalar T) const;
148 
149  // Enthalpy of formation [J/kg]
150  // inline scalar hf() const;
151 
152  // Absolute enthalpy [J/kg]
153  // inline scalar ha(const scalar p, const scalar T) const;
154 
155  // Heat capacity at constant volume [J/kg/K]
156  // inline scalar Cv(const scalar p, const scalar T) const;
157 
158  // Sensible internal energy [J/kg]
159  // inline scalar es(const scalar p, const scalar T) const;
160 
161  // Absolute internal energy [J/kg]
162  // inline scalar ea(const scalar p, const scalar T) const;
163 
164  // Entropy [J/kg/K]
165  // inline scalar s(const scalar p, const scalar T) const;
166 
167 
168  // Derived properties
169 
170  //- Heat capacity at constant pressure/volume [J/kg/K]
171  inline scalar Cpv(const scalar p, const scalar T) const;
172 
173  //- Gamma = Cp/Cv []
174  inline scalar gamma(const scalar p, const scalar T) const;
175 
176  //- Enthalpy/Internal energy [J/kg]
177  inline scalar he(const scalar p, const scalar T) const;
178 
179  //- Gibbs free energy [J/kg]
180  inline scalar g(const scalar p, const scalar T) const;
181 
182  //- Helmholtz free energy [J/kg]
183  inline scalar a(const scalar p, const scalar T) const;
184 
185 
186  // Equilibrium reaction thermodynamics
187 
188  //- Equilibrium constant [] i.t.o fugacities
189  // = PIi(fi/Pstd)^nui
190  inline scalar K(const scalar p, const scalar T) const;
191 
192  //- Equilibrium constant [] i.t.o. partial pressures
193  // = PIi(pi/Pstd)^nui
194  // For low pressures (where the gas mixture is near perfect) Kp = K
195  inline scalar Kp(const scalar p, const scalar T) const;
196 
197  //- Equilibrium constant i.t.o. molar concentration
198  // = PIi(ci/cStd)^nui
199  // For low pressures (where the gas mixture is near perfect)
200  // Kc = Kp(Pstd/(RR*T))^nu
201  inline scalar Kc(const scalar p, const scalar T) const;
202 
203  //- Equilibrium constant [] i.t.o. mole-fractions
204  // For low pressures (where the gas mixture is near perfect)
205  // Kx = Kp(Pstd/p)^nui
206  inline scalar Kx(const scalar p, const scalar T) const;
207 
208  //- Equilibrium constant [] i.t.o. number of moles
209  // For low pressures (where the gas mixture is near perfect)
210  // Kn = Kp(n*Pstd/p)^nui where n = number of moles in mixture
211  inline scalar Kn
212  (
213  const scalar p,
214  const scalar T,
215  const scalar n
216  ) const;
217 
218 
219  // Energy->temperature inversion functions
220 
221  //- Return the temperature corresponding to the value of the
222  // thermodynamic property f, given the function f = F(p, T)
223  // and dF(p, T)/dT
224  template
225  <
226  class ThermoType,
227  class FType,
228  class dFdTType,
229  class LimitType
230  >
231  inline static scalar T
232  (
233  const ThermoType& thermo,
234  const scalar f,
235  const scalar p,
236  const scalar T0,
237  FType F,
238  dFdTType dFdT,
239  LimitType limit,
240  const bool diagnostics = false
241  );
242 
243  //- Temperature from enthalpy or internal energy
244  // given an initial temperature T0
245  inline scalar The
246  (
247  const scalar H,
248  const scalar p,
249  const scalar T0
250  ) const;
251 
252  //- Temperature from sensible enthalpy given an initial T0
253  inline scalar THs
254  (
255  const scalar hs,
256  const scalar p,
257  const scalar T0
258  ) const;
259 
260  //- Temperature from absolute enthalpy
261  // given an initial temperature T0
262  inline scalar Tha
263  (
264  const scalar H,
265  const scalar p,
266  const scalar T0
267  ) const;
268 
269  //- Temperature from sensible internal energy
270  // given an initial temperature T0
271  inline scalar Tes
272  (
273  const scalar E,
274  const scalar p,
275  const scalar T0
276  ) const;
277 
278  //- Temperature from absolute internal energy
279  // given an initial temperature T0
280  inline scalar Tea
281  (
282  const scalar E,
283  const scalar p,
284  const scalar T0
285  ) const;
286 
287 
288  // Derivative term used for Jacobian
289 
290  //- Derivative of B (according to Niemeyer et al.)
291  // w.r.t. temperature
292  inline scalar dKcdTbyKc(const scalar p, const scalar T) const;
293 
294 
295  // I-O
296 
297  //- Write to Ostream
298  void write(Ostream& os) const;
299 
300 
301  // Member Operators
302 
303  inline void operator+=(const thermo&);
304  inline void operator*=(const scalar);
305 
306 
307  // Friend operators
308 
309  friend thermo operator+ <Thermo, Type>
310  (
311  const thermo&,
312  const thermo&
313  );
314 
315  friend thermo operator* <Thermo, Type>
316  (
317  const scalar s,
318  const thermo&
319  );
320 
321  friend thermo operator== <Thermo, Type>
322  (
323  const thermo&,
324  const thermo&
325  );
326 
327 
328  // Ostream Operator
329 
330  friend Ostream& operator<< <Thermo, Type>
331  (
332  Ostream&,
333  const thermo&
334  );
335 };
336 
337 
338 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
339 
340 } // End namespace species
341 } // End namespace Foam
342 
343 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
344 
345 #include "thermoI.H"
346 
347 #ifdef NoRepository
348  #include "thermo.C"
349 #endif
350 
351 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
352 
353 #endif
354 
355 // ************************************************************************* //
scalar hs(const scalar p, const scalar T) const
Definition: EtoHthermo.H:11
label n
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
Definition: thermo.H:92
A class for handling words, derived from string.
Definition: word.H:62
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.name(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
K
Definition: pEqn.H:75
const dimensionedScalar F
Faraday constant: default SI units: [C/kmol].
Thermodynamic scalar constants.
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
Namespace for OpenFOAM.
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
complex limit(const complex &, const complex &)
Definition: complexI.H:202
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
thermo he()
labelList f(nPoints)
dictionary dict
volScalarField & p
fluidMulticomponentThermo & thermo
Definition: createFields.H:31
scalar T0
Definition: createFields.H:22