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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2012-2017 OpenFOAM Foundation
9  Copyright (C) 2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
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 "hPowerThermo.H"
30 #include "specie.H"
31 
32 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
33 
34 template<class EquationOfState>
36 (
37  const scalar T
38 ) const
39 {
40  if (T < 0)
41  {
43  << "attempt to evaluate hPowerThermo<EquationOfState>"
44  " for negative temperature " << T
45  << abort(FatalError);
46  }
47 }
48 
49 
50 template<class EquationOfState>
52 (
53  const word& name,
54  const hPowerThermo& jt
55 )
56 :
57  EquationOfState(name, jt),
58  c0_(jt.c0_),
59  n0_(jt.n0_),
60  Tref_(jt.Tref_),
61  Hf_(jt.Hf_)
62 {}
63 
64 
65 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
66 
67 template<class EquationOfState>
69 (
70  const EquationOfState& st,
71  const scalar c0,
72  const scalar n0,
73  const scalar Tref,
74  const scalar Hf
75 )
76 :
77  EquationOfState(st),
78  c0_(c0),
79  n0_(n0),
80  Tref_(Tref),
81  Hf_(Hf)
82 {}
83 
84 
85 template<class EquationOfState>
88 {
90 }
91 
92 
93 template<class EquationOfState>
96 {
98 }
99 
100 
101 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
102 
103 template<class EquationOfState>
105 (
106  const scalar T
107 ) const
108 {
109  return T;
110 }
111 
112 
113 template<class EquationOfState>
114 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::Cp
115 (
116  const scalar p, const scalar T
117 ) const
118 {
119  return c0_*pow(T/Tref_, n0_) + EquationOfState::Cp(p, T);
120 }
121 
122 
123 template<class EquationOfState>
124 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::Ha
125 (
126  const scalar p, const scalar T
127 ) const
128 {
129  return Hs(p, T) + Hc();
130 }
131 
132 
133 template<class EquationOfState>
134 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::Hs
135 (
136  const scalar p, const scalar T
137 ) const
138 {
139  return
140  c0_*(pow(T, n0_ + 1) - pow(Tstd, n0_ + 1))/(pow(Tref_, n0_)*(n0_ + 1))
141  + EquationOfState::H(p, T);
142 }
143 
144 
145 template<class EquationOfState>
146 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::Hc() const
147 {
148  return Hf_;
149 }
150 
151 
152 template<class EquationOfState>
153 inline Foam::scalar Foam::hPowerThermo<EquationOfState>::S
154 (
155  const scalar p, const scalar T
156 ) const
157 {
158  return
159  c0_*(pow(T, n0_) - pow(Tstd, n0_))/(pow(Tref_, n0_)*n0_)
160  + EquationOfState::S(p, T);
161 }
162 
163 
164 template<class EquationOfState>
166 (
167  const scalar T
168 ) const
169 {
170  return
171  c0_*(pow(T, n0_ + 1) - pow(Tstd, n0_ + 1))/(pow(Tref_, n0_)*(n0_ + 1))
172  + Hc()
173  - c0_*(pow(T, n0_) - pow(Tstd, n0_))*T/(pow(Tref_, n0_)*n0_);
174 }
175 
176 
177 template<class EquationOfState>
179 (
180  const scalar p, const scalar T
181 ) const
182 {
183  // To be implemented
185  return 0;
186 }
187 
188 
189 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
190 
191 template<class EquationOfState>
193 (
195 )
196 {
197  scalar Y1 = this->Y();
198 
199  EquationOfState::operator+=(ct);
200 
201  if (mag(this->Y()) > SMALL)
202  {
203  Y1 /= this->Y();
204  const scalar Y2 = ct.Y()/this->Y();
205 
206  Hf_ = Y1*Hf_ + Y2*ct.Hf_;
207  c0_ = Y1*c0_ + Y2*ct.c0_;
208  n0_ = Y1*n0_ + Y2*ct.n0_;
209  Tref_ = Y1*Tref_ + Y2*ct.Tref_;
210  }
211 }
212 
213 
214 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
215 
216 template<class EquationOfState>
217 inline Foam::hPowerThermo<EquationOfState> Foam::operator+
218 (
221 )
222 {
223  EquationOfState eofs
224  (
225  static_cast<const EquationOfState&>(ct1)
226  + static_cast<const EquationOfState&>(ct2)
227  );
228 
229  if (mag(eofs.Y()) < SMALL)
230  {
232  (
233  eofs,
234  ct1.c0_,
235  ct1.n0_,
236  ct1.Tref_,
237  ct1.Hf_
238  );
239  }
240  else
241  {
242  return hPowerThermo<EquationOfState>
243  (
244  eofs,
245  ct1.Y()/eofs.Y()*ct1.c0_
246  + ct2.Y()/eofs.Y()*ct2.c0_,
247  ct1.Y()/eofs.Y()*ct1.n0_
248  + ct2.Y()/eofs.Y()*ct2.n0_,
249  ct1.Y()/eofs.Y()*ct1.Tref_
250  + ct2.Y()/eofs.Y()*ct2.Tref_,
251  ct1.Y()/eofs.Y()*ct1.Hf_
252  + ct2.Y()/eofs.Y()*ct2.Hf_
253  );
254  }
255 }
256 
257 
258 template<class EquationOfState>
259 inline Foam::hPowerThermo<EquationOfState> Foam::operator*
260 (
261  const scalar s,
262  const hPowerThermo<EquationOfState>& ct
263 )
264 {
265  return hPowerThermo<EquationOfState>
266  (
267  s*static_cast<const EquationOfState&>(ct),
268  ct.c0_,
269  ct.n0_,
270  ct.Tref_,
271  ct.Hf_
272  );
273 }
274 
275 
276 template<class EquationOfState>
277 inline Foam::hPowerThermo<EquationOfState> Foam::operator==
278 (
279  const hPowerThermo<EquationOfState>& ct1,
280  const hPowerThermo<EquationOfState>& ct2
281 )
282 {
283  EquationOfState eofs
284  (
285  static_cast<const EquationOfState&>(ct1)
286  == static_cast<const EquationOfState&>(ct2)
287  );
288 
289  return hPowerThermo<EquationOfState>
290  (
291  eofs,
292  ct2.Y()/eofs.Y()*ct2.c0_
293  - ct1.Y()/eofs.Y()*ct1.c0_,
294  ct2.Y()/eofs.Y()*ct2.n0_
295  - ct1.Y()/eofs.Y()*ct1.n0_,
296  ct2.Y()/eofs.Y()*ct2.Tref_
297  - ct1.Y()/eofs.Y()*ct1.Tref_,
298  ct2.Y()/eofs.Y()*ct2.Hf_
299  - ct1.Y()/eofs.Y()*ct1.Hf_
300  );
301 }
302 
303 
304 // ************************************************************************* //
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::hPowerThermo
Power-function based thermodynamics package templated on EquationOfState.
Definition: hPowerThermo.H:60
s
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))
Definition: gmvOutputSpray.H:25
hPowerThermo.H
Hs
scalar Hs(const scalar p, const scalar T) const
Definition: EtoHthermo.H:17
H
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
Foam::hPowerThermo::limit
scalar limit(const scalar T) const
Limit temperature to be within the range.
Definition: hPowerThermoI.H:105
Foam::constant::standard::Tstd
const dimensionedScalar Tstd
Standard temperature.
Definition: thermodynamicConstants.C:49
specie.H
Foam::hPowerThermo::Hs
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
Definition: hPowerThermoI.H:135
NotImplemented
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:517
Foam::hPowerThermo::clone
autoPtr< hPowerThermo > clone() const
Construct and return a clone.
Definition: hPowerThermoI.H:87
Foam::T
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
Definition: FieldFieldFunctions.C:58
Foam::hPowerThermo::Hc
scalar Hc() const
Chemical enthalpy [J/kg].
Definition: hPowerThermoI.H:146
Foam::hPowerThermo::Ha
scalar Ha(const scalar p, const scalar T) const
Absolute Enthalpy [J/kg].
Definition: hPowerThermoI.H:125
Foam::pow
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Definition: dimensionedScalar.C:75
Foam::hPowerThermo::S
scalar S(const scalar p, const scalar T) const
Entropy [J/(kg K)].
Definition: hPowerThermoI.H:154
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::hPowerThermo::dCpdT
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
Definition: hPowerThermoI.H:179
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
T
const volScalarField & T
Definition: createFieldRefs.H:2
Y
PtrList< volScalarField > & Y
Definition: createFieldRefs.H:7
Foam::hPowerThermo::New
static autoPtr< hPowerThermo > New(const dictionary &dict)
Selector from dictionary.
Definition: hPowerThermoI.H:95
Foam::New
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
Definition: DimensionedFieldReuseFunctions.H:105
Foam::autoPtr
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: HashPtrTable.H:53
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::hPowerThermo::Cp
scalar Cp(const scalar p, const scalar T) const
Heat capacity at constant pressure [J/(kg K)].
Definition: hPowerThermoI.H:115
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Cp
const volScalarField & Cp
Definition: EEqn.H:7
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::hPowerThermo::Gstd
scalar Gstd(const scalar T) const
Gibbs free energy of the mixture in the standard state [J/kg].
Definition: hPowerThermoI.H:166