hPowerThermoI.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) 2012-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 \*---------------------------------------------------------------------------*/
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  << "attempt to evaluate hPowerThermo<EquationOfState>"
41  " for negative temperature " << T
42  << abort(FatalError);
43  }
44 }
45 
46 
47 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
48 
49 template<class EquationOfState>
51 (
52  const EquationOfState& st,
53  const scalar c0,
54  const scalar n0,
55  const scalar Tref,
56  const scalar hf
57 )
58 :
59  EquationOfState(st),
60  c0_(c0),
61  n0_(n0),
62  Tref_(Tref),
63  hf_(hf)
64 {}
65 
66 
67 template<class EquationOfState>
69 (
70  const word& name,
71  const hPowerThermo& jt
72 )
73 :
74  EquationOfState(name, jt),
75  c0_(jt.c0_),
76  n0_(jt.n0_),
77  Tref_(jt.Tref_),
78  hf_(jt.hf_)
79 {}
80 
81 
82 template<class EquationOfState>
85 {
87  (
89  );
90 }
91 
92 
93 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
94 
95 template<class EquationOfState>
97 (
98  const scalar T
99 ) const
100 {
101  return T;
102 }
103 
104 
105 template<class EquationOfState>
107 (
108  const scalar p,
109  const scalar T
110 ) const
111 {
112  return c0_*pow(T/Tref_, n0_) + EquationOfState::Cp(p, T);
113 }
114 
115 
116 template<class EquationOfState>
118 (
119  const scalar p,
120  const scalar T
121 ) const
122 {
123  return
124  c0_*(pow(T, n0_ + 1) - pow(Tstd, n0_ + 1))/(pow(Tref_, n0_)*(n0_ + 1))
125  + EquationOfState::h(p, T);
126 }
127 
128 
129 template<class EquationOfState>
131 (
132  const scalar p,
133  const scalar T
134 ) const
135 {
136  return hs(p, T) + hf();
137 }
138 
139 
140 template<class EquationOfState>
141 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::hf() const
142 {
143  return hf_;
144 }
145 
146 
147 template<class EquationOfState>
149 (
150  const scalar p,
151  const scalar T
152 ) const
153 {
154  return
155  c0_*(pow(T, n0_) - pow(Tstd, n0_))/(pow(Tref_, n0_)*n0_)
156  + EquationOfState::sp(p, T);
157 }
158 
159 
160 template<class EquationOfState>
162 (
163  const scalar T
164 ) const
165 {
166  return
167  c0_*(pow(T, n0_ + 1) - pow(Tstd, n0_ + 1))/(pow(Tref_, n0_)*(n0_ + 1))
168  + hf()
169  - c0_*(pow(T, n0_) - pow(Tstd, n0_))*T/(pow(Tref_, n0_)*n0_);
170 }
171 
172 
173 template<class EquationOfState>
175 (
176  const scalar p,
177  const scalar T
178 ) const
179 {
180  // To be implemented
182  return 0;
183 }
184 
185 
186 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
187 
188 template<class EquationOfState>
190 (
192 )
193 {
194  scalar Y1 = this->Y();
195 
196  EquationOfState::operator+=(ct);
197 
198  if (mag(this->Y()) > small)
199  {
200  Y1 /= this->Y();
201  const scalar Y2 = ct.Y()/this->Y();
202 
203  hf_ = Y1*hf_ + Y2*ct.hf_;
204  c0_ = Y1*c0_ + Y2*ct.c0_;
205  n0_ = Y1*n0_ + Y2*ct.n0_;
206  Tref_ = Y1*Tref_ + Y2*ct.Tref_;
207  }
208 }
209 
210 
211 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
212 
213 template<class EquationOfState>
214 inline Foam::hPowerThermo<EquationOfState> Foam::operator+
215 (
218 )
219 {
220  EquationOfState eofs
221  (
222  static_cast<const EquationOfState&>(ct1)
223  + static_cast<const EquationOfState&>(ct2)
224  );
225 
226  if (mag(eofs.Y()) < small)
227  {
229  (
230  eofs,
231  ct1.c0_,
232  ct1.n0_,
233  ct1.Tref_,
234  ct1.hf_
235  );
236  }
237  else
238  {
239  return hPowerThermo<EquationOfState>
240  (
241  eofs,
242  ct1.Y()/eofs.Y()*ct1.c0_
243  + ct2.Y()/eofs.Y()*ct2.c0_,
244  ct1.Y()/eofs.Y()*ct1.n0_
245  + ct2.Y()/eofs.Y()*ct2.n0_,
246  ct1.Y()/eofs.Y()*ct1.Tref_
247  + ct2.Y()/eofs.Y()*ct2.Tref_,
248  ct1.Y()/eofs.Y()*ct1.hf_
249  + ct2.Y()/eofs.Y()*ct2.hf_
250  );
251  }
252 }
253 
254 
255 template<class EquationOfState>
256 inline Foam::hPowerThermo<EquationOfState> Foam::operator*
257 (
258  const scalar s,
259  const hPowerThermo<EquationOfState>& ct
260 )
261 {
262  return hPowerThermo<EquationOfState>
263  (
264  s*static_cast<const EquationOfState&>(ct),
265  ct.c0_,
266  ct.n0_,
267  ct.Tref_,
268  ct.hf_
269  );
270 }
271 
272 
273 template<class EquationOfState>
274 inline Foam::hPowerThermo<EquationOfState> Foam::operator==
275 (
276  const hPowerThermo<EquationOfState>& ct1,
277  const hPowerThermo<EquationOfState>& ct2
278 )
279 {
280  EquationOfState eofs
281  (
282  static_cast<const EquationOfState&>(ct1)
283  == static_cast<const EquationOfState&>(ct2)
284  );
285 
286  return hPowerThermo<EquationOfState>
287  (
288  eofs,
289  ct2.Y()/eofs.Y()*ct2.c0_
290  - ct1.Y()/eofs.Y()*ct1.c0_,
291  ct2.Y()/eofs.Y()*ct2.n0_
292  - ct1.Y()/eofs.Y()*ct1.n0_,
293  ct2.Y()/eofs.Y()*ct2.Tref_
294  - ct1.Y()/eofs.Y()*ct1.Tref_,
295  ct2.Y()/eofs.Y()*ct2.hf_
296  - ct1.Y()/eofs.Y()*ct1.hf_
297  );
298 }
299 
300 
301 // ************************************************************************* //
scalar hs(const scalar p, const scalar T) const
Definition: EtoHthermo.H:11
scalar Cp(const scalar p, const scalar T) const
Definition: EtoHthermo.H:2
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
Enthalpy based thermodynamics package using a power function of temperature for the constant heat cap...
Definition: hPowerThermo.H:129
scalar gStd(const scalar T) const
Gibbs free energy of the mixture in the standard state [J/kg].
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
Definition: hPowerThermoI.H:97
autoPtr< hPowerThermo > clone() const
Construct and return a clone.
Definition: hPowerThermoI.H:84
scalar s(const scalar p, const scalar T) const
Entropy [J/kg/K].
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
scalar Cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/kg/K].
scalar ha(const scalar p, const scalar T) const
Absolute enthalpy [J/kg].
scalar hf() const
Enthalpy of formation [J/kg].
scalar hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
hPowerThermo(const EquationOfState &st, const scalar c0, const scalar n0, const scalar Tref, const scalar hf)
Construct from components.
Definition: hPowerThermoI.H:51
A class for handling words, derived from string.
Definition: word.H:62
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:381
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
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))
const dimensionedScalar Tstd
Standard temperature.
const dimensionedScalar h
Planck constant.
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
errorManip< error > abort(error &err)
Definition: errorManip.H:131
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensioned< scalar > mag(const dimensioned< Type > &)
error FatalError
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
volScalarField & p
PtrList< volScalarField > & Y