FreeFOAM The Cross-Platform CFD Toolkit
specieThermo.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) 1991-2010 OpenCFD Ltd.
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::specieThermo
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  specieThermoI.H
34  specieThermo.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef specieThermo_H
39 #define specieThermo_H
40 
41 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45 
46 // Forward declaration of friend functions and operators
47 
48 template<class thermo> class specieThermo;
49 
50 template<class thermo>
51 inline specieThermo<thermo> operator+
52 (
53  const specieThermo<thermo>&,
54  const specieThermo<thermo>&
55 );
56 
57 template<class thermo>
58 inline specieThermo<thermo> operator-
59 (
60  const specieThermo<thermo>&,
61  const specieThermo<thermo>&
62 );
63 
64 template<class thermo>
65 inline specieThermo<thermo> operator*
66 (
67  const scalar,
68  const specieThermo<thermo>&
69 );
70 
71 template<class thermo>
72 inline specieThermo<thermo> operator==
73 (
74  const specieThermo<thermo>&,
75  const specieThermo<thermo>&
76 );
77 
78 template<class thermo>
79 Ostream& operator<<
80 (
81  Ostream&,
82  const specieThermo<thermo>&
83 );
84 
85 
86 /*---------------------------------------------------------------------------*\
87  Class specieThermo Declaration
88 \*---------------------------------------------------------------------------*/
89 
90 template<class thermo>
92 :
93  public thermo
94 {
95  // Private data
96 
97  //- Convergence tolerance of energy -> temperature inversion functions
98  static const scalar tol_;
99 
100  //- Max number of iterations in energy->temperature inversion functions
101  static const int maxIter_;
102 
103 
104  // Private member functions
105 
106  //- return the temperature corresponding to the value of the
107  // thermodynamic property f, given the function f = F(T) and dF(T)/dT
108  inline scalar T
109  (
110  scalar f,
111  scalar T0,
112  scalar (specieThermo::*F)(const scalar) const,
113  scalar (specieThermo::*dFdT)(const scalar) const
114  ) const;
115 
116 
117 public:
118 
119  // Constructors
120 
121  //- construct from components
122  inline specieThermo(const thermo& sp);
123 
124  //- Construct from Istream
126 
127  //- Construct as named copy
128  inline specieThermo(const word& name, const specieThermo&);
129 
130 
131  // Member Functions
132 
133  // Fundamaental properties
134  // (These functions must be provided in derived types)
135 
136  // Heat capacity at constant pressure [J/(kmol K)]
137  //scalar cp(const scalar) const;
138 
139  // Enthalpy [J/kmol]
140  //scalar h(const scalar) const;
141 
142  // Sensible enthalpy [J/kmol]
143  //scalar hs(const scalar) const;
144 
145  // Chemical enthalpy [J/kmol]
146  //scalar hc(const scalar) const;
147 
148  // Entropy [J/(kmol K)]
149  //scalar s(const scalar) const;
150 
151 
152  // Calculate and return derived properties
153  // (These functions need not provided in derived types)
154 
155  // Mole specific properties
156 
157  //- Heat capacity at constant volume [J/(kmol K)]
158  inline scalar cv(const scalar T) const;
159 
160  //- gamma = cp/cv []
161  inline scalar gamma(const scalar T) const;
162 
163  //- Internal energy [J/kmol]
164  inline scalar e(const scalar T) const;
165 
166  //- Sensible internal energy [J/kmol]
167  inline scalar es(const scalar T) const;
168 
169  //- Gibbs free energy [J/kmol]
170  inline scalar g(const scalar T) const;
171 
172  //- Helmholtz free energy [J/kmol]
173  inline scalar a(const scalar T) const;
174 
175 
176  // Mass specific properties
177 
178  //- Heat capacity at constant pressure [J/(kg K)]
179  inline scalar Cp(const scalar T) const;
180 
181  //- Heat capacity at constant volume [J/(kg K)]
182  inline scalar Cv(const scalar T) const;
183 
184  //- Enthalpy [J/kg]
185  inline scalar H(const scalar T) const;
186 
187  //- Sensible enthalpy [J/kg]
188  inline scalar Hs(const scalar T) const;
189 
190  //- Chemical enthalpy [J/kg]
191  inline scalar Hc() const;
192 
193  //- Entropy [J/(kg K)]
194  inline scalar S(const scalar T) const;
195 
196  //- Internal energy [J/kg]
197  inline scalar E(const scalar T) const;
198 
199  //- Gibbs free energy [J/kg]
200  inline scalar G(const scalar T) const;
201 
202  //- Helmholtz free energy [J/kg]
203  inline scalar A(const scalar T) const;
204 
205 
206  // Equilibrium reaction thermodynamics
207 
208  //- Equilibrium constant [] i.t.o fugacities
209  // = PIi(fi/Pstd)^nui
210  inline scalar K(const scalar T) const;
211 
212  //- Equilibrium constant [] i.t.o. partial pressures
213  // = PIi(pi/Pstd)^nui
214  // For low pressures (where the gas mixture is near perfect) Kp = K
215  inline scalar Kp(const scalar T) const;
216 
217  //- Equilibrium constant i.t.o. molar concentration
218  // = PIi(ci/cstd)^nui
219  // For low pressures (where the gas mixture is near perfect)
220  // Kc = Kp(pstd/(RR*T))^nu
221  inline scalar Kc(const scalar T) const;
222 
223  //- Equilibrium constant [] i.t.o. mole-fractions
224  // For low pressures (where the gas mixture is near perfect)
225  // Kx = Kp(pstd/p)^nui
226  inline scalar Kx(const scalar T, const scalar p) const;
227 
228  //- Equilibrium constant [] i.t.o. number of moles
229  // For low pressures (where the gas mixture is near perfect)
230  // Kn = Kp(n*pstd/p)^nui where n = number of moles in mixture
231  inline scalar Kn
232  (
233  const scalar T,
234  const scalar p,
235  const scalar n
236  ) const;
237 
238 
239  // Energy->temperature inversion functions
240 
241  //- Temperature from Enthalpy given an initial temperature T0
242  inline scalar TH(const scalar H, const scalar T0) const;
243 
244  //- Temperature from sensible Enthalpy given an initial T0
245  inline scalar THs(const scalar Hs, const scalar T0) const;
246 
247  //- Temperature from internal energy given an initial temperature T0
248  inline scalar TE(const scalar E, const scalar T0) const;
249 
250 
251  // Member operators
252 
253  inline void operator+=(const specieThermo&);
254  inline void operator-=(const specieThermo&);
255 
256  inline void operator*=(const scalar);
257 
258 
259  // Friend operators
260 
261  friend specieThermo operator+ <thermo>
262  (
263  const specieThermo&,
264  const specieThermo&
265  );
266 
267  friend specieThermo operator- <thermo>
268  (
269  const specieThermo&,
270  const specieThermo&
271  );
272 
273  friend specieThermo operator* <thermo>
274  (
275  const scalar s,
276  const specieThermo&
277  );
278 
279  friend specieThermo operator== <thermo>
280  (
281  const specieThermo&,
282  const specieThermo&
283  );
284 
285 
286  // Ostream Operator
287 
288  friend Ostream& operator<< <thermo>
289  (
290  Ostream&,
291  const specieThermo&
292  );
293 };
294 
295 
296 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
297 
298 } // End namespace Foam
299 
300 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
301 
302 #include <specie/specieThermoI.H>
303 
304 #ifdef NoRepository
305 # include <specie/specieThermo.C>
306 #endif
307 
308 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
309 
310 #endif
311 
312 // ************************ vim: set sw=4 sts=4 et: ************************ //