hPowerThermoI.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) 2012-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 "hPowerThermo.H"
27 #include "specie.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 template<class EquationOfState>
33 (
34  const scalar T
35 ) const
36 {
37  if (T < 0)
38  {
40  (
41  "hPowerThermo<EquationOfState>::checkT(const scalar T) const"
42  ) << "attempt to evaluate hPowerThermo<EquationOfState>"
43  " for negative temperature " << T
44  << abort(FatalError);
45  }
46 }
47 
48 
49 template<class EquationOfState>
51 (
52  const word& name,
53  const hPowerThermo& jt
54 )
55 :
56  EquationOfState(name, jt),
57  c0_(jt.c0_),
58  n0_(jt.n0_),
59  Tref_(jt.Tref_),
60  Hf_(jt.Hf_)
61 {}
62 
63 
64 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
65 
66 template<class EquationOfState>
68 (
69  const EquationOfState& st,
70  const scalar c0,
71  const scalar n0,
72  const scalar Tref,
73  const scalar Hf
74 )
75 :
76  EquationOfState(st),
77  c0_(c0),
78  n0_(n0),
79  Tref_(Tref),
80  Hf_(Hf)
81 {}
82 
83 
84 template<class EquationOfState>
87 {
89  (
91  );
92 }
93 
94 
95 template<class EquationOfState>
98 {
100  (
102  );
103 }
104 
105 
106 template<class EquationOfState>
109 {
111  (
113  );
114 }
115 
116 
117 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
118 
119 template<class EquationOfState>
121 (
122  const scalar T
123 ) const
124 {
125  return T;
126 }
127 
128 
129 template<class EquationOfState>
130 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::cp
131 (
132  const scalar p, const scalar T
133 ) const
134 {
135  return c0_*pow(T/Tref_, n0_);
136 }
137 
138 
139 template<class EquationOfState>
140 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::ha
141 (
142  const scalar p, const scalar T
143 ) const
144 {
145  return hs(p, T) + hc();
146 }
147 
148 
149 template<class EquationOfState>
150 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::hs
151 (
152  const scalar p, const scalar T
153 ) const
154 {
155  return
156  c0_*(pow(T, n0_ + 1) - pow(Tstd, n0_ + 1))/(pow(Tref_, n0_)*(n0_ + 1));
157 }
158 
159 
160 template<class EquationOfState>
161 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::hc() const
162 {
163  return Hf_;
164 }
165 
166 
167 template<class EquationOfState>
168 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::s
169 (
170  const scalar p, const scalar T
171 ) const
172 {
173  return
174  c0_*(pow(T, n0_) - pow(Tstd, n0_))/(pow(Tref_, n0_)*n0_)
175  + EquationOfState::s(p, T);
176 }
177 
178 
179 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
180 
181 template<class EquationOfState>
182 inline void Foam::hPowerThermo<EquationOfState>::operator+=
183 (
185 )
186 {
187  scalar molr1 = this->nMoles();
188 
189  EquationOfState::operator+=(ct);
190  molr1 /= this->nMoles();
191  scalar molr2 = ct.nMoles()/this->nMoles();
192 
193  Hf_ = molr1*Hf_ + molr2*ct.Hf_;
194  c0_ = molr1*c0_ + molr2*ct.c0_;
195  n0_ = molr1*n0_ + molr2*ct.n0_;
196  Tref_ = molr1*Tref_ + molr2*ct.Tref_;
197 }
198 
199 
200 template<class EquationOfState>
201 inline void Foam::hPowerThermo<EquationOfState>::operator-=
202 (
204 )
205 {
206  scalar molr1 = this->nMoles();
207 
208  EquationOfState::operator-=(ct);
209 
210  molr1 /= this->nMoles();
211  scalar molr2 = ct.nMoles()/this->nMoles();
212 
213  Hf_ = molr1*Hf_ - molr2*ct.Hf_;
214  c0_ = (molr1*c0_ - molr2*ct.c0_);
215  n0_ = (molr1*n0_ - molr2*ct.n0_);
216  Tref_ = (molr1*Tref_ - molr2*ct.Tref_);
217 }
218 
219 
220 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
221 
222 template<class EquationOfState>
223 inline Foam::hPowerThermo<EquationOfState> Foam::operator+
224 (
227 )
228 {
229  EquationOfState eofs
230  (
231  static_cast<const EquationOfState&>(ct1)
232  + static_cast<const EquationOfState&>(ct2)
233  );
234 
236  (
237  eofs,
238  ct1.nMoles()/eofs.nMoles()*ct1.c0_
239  + ct2.nMoles()/eofs.nMoles()*ct2.c0_,
240  ct1.nMoles()/eofs.nMoles()*ct1.n0_
241  + ct2.nMoles()/eofs.nMoles()*ct2.n0_,
242  ct1.nMoles()/eofs.nMoles()*ct1.Tref_
243  + ct2.nMoles()/eofs.nMoles()*ct2.Tref_,
244  ct1.nMoles()/eofs.nMoles()*ct1.Hf_
245  + ct2.nMoles()/eofs.nMoles()*ct2.Hf_
246  );
247 }
248 
249 
250 template<class EquationOfState>
251 inline Foam::hPowerThermo<EquationOfState> Foam::operator-
252 (
255 )
256 {
257  EquationOfState eofs
258  (
259  static_cast<const EquationOfState&>(ct1)
260  + static_cast<const EquationOfState&>(ct2)
261  );
262 
264  (
265  eofs,
266  ct1.nMoles()/eofs.nMoles()*ct1.c0_
267  - ct2.nMoles()/eofs.nMoles()*ct2.c0_,
268  ct1.nMoles()/eofs.nMoles()*ct1.n0_
269  - ct2.nMoles()/eofs.nMoles()*ct2.n0_,
270  ct1.nMoles()/eofs.nMoles()*ct1.Tref_
271  - ct2.nMoles()/eofs.nMoles()*ct2.Tref_,
272  ct1.nMoles()/eofs.nMoles()*ct1.Hf_
273  - ct2.nMoles()/eofs.nMoles()*ct2.Hf_
274  );
275 }
276 
277 
278 template<class EquationOfState>
279 inline Foam::hPowerThermo<EquationOfState> Foam::operator*
280 (
281  const scalar s,
283 )
284 {
286  (
287  s*static_cast<const EquationOfState&>(ct),
288  ct.c0_,
289  ct.n0_,
290  ct.Tref_,
291  ct.Hf_
292  );
293 }
294 
295 
296 template<class EquationOfState>
297 inline Foam::hPowerThermo<EquationOfState> Foam::operator==
298 (
301 )
302 {
303  return ct2 - ct1;
304 }
305 
306 
307 // ************************************************************************* //
scalar hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
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 ))
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
A class for handling words, derived from string.
Definition: word.H:59
scalar s(const scalar p, const scalar T) const
Entropy [J/(kmol K)].
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
const volScalarField & T
Definition: createFields.H:25
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
dictionary dict
Power-function based thermodynamics package templated on EquationOfState.
Definition: hPowerThermo.H:54
const dimensionedScalar Tstd
Standard temperature.
autoPtr< hPowerThermo > clone() const
Construct and return a clone.
Definition: hPowerThermoI.H:86
errorManip< error > abort(error &err)
Definition: errorManip.H:131
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
error FatalError
scalar ha(const scalar p, const scalar T) const
Absolute enthalpy [J/kmol].
scalar cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kg K)].
static autoPtr< hPowerThermo > New(Istream &is)
Selector from Istream.
Definition: hPowerThermoI.H:97
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117
scalar hc() const
Chemical enthalpy [J/kg].