eConstThermoI.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 \*---------------------------------------------------------------------------*/
25 
26 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
27 
28 template<class EquationOfState>
30 (
31  const EquationOfState& st,
32  const scalar Cv,
33  const scalar hf,
34  const scalar Tref,
35  const scalar esRef
36 )
37 :
38  EquationOfState(st),
39  Cv_(Cv),
40  hf_(hf),
41  Tref_(Tref),
42  esRef_(esRef)
43 {}
44 
45 
46 template<class EquationOfState>
48 (
49  const word& name,
50  const eConstThermo& ct
51 )
52 :
53  EquationOfState(name, ct),
54  Cv_(ct.Cv_),
55  hf_(ct.hf_),
56  Tref_(ct.Tref_),
57  esRef_(ct.esRef_)
58 {}
59 
60 
61 template<class EquationOfState>
64 {
66  (
68  );
69 }
70 
71 
72 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
73 
74 template<class EquationOfState>
76 (
77  const scalar T
78 ) const
79 {
80  return T;
81 }
82 
83 
84 template<class EquationOfState>
86 (
87  const scalar p,
88  const scalar T
89 ) const
90 {
91  return Cv_ + EquationOfState::Cv(p, T);
92 }
93 
94 
95 template<class EquationOfState>
97 (
98  const scalar p,
99  const scalar T
100 ) const
101 {
102  return Cv_*(T - Tref_) + esRef_ + EquationOfState::e(p, T);
103 }
104 
105 
106 template<class EquationOfState>
108 (
109  const scalar p,
110  const scalar T
111 ) const
112 {
113  return es(p, T) + hf();
114 }
115 
116 
117 template<class EquationOfState>
118 inline Foam::scalar Foam::eConstThermo<EquationOfState>::hf() const
119 {
120  return hf_;
121 }
122 
123 
124 template<class EquationOfState>
126 (
127  const scalar p,
128  const scalar T
129 ) const
130 {
131  return Cp(p, T)*log(T/Tstd) + EquationOfState::sv(p, T);
132 }
133 
134 
135 template<class EquationOfState>
137 (
138  const scalar T
139 ) const
140 {
141  return
142  Cv_*(T - Tref_) + esRef_ + hf() + Pstd/EquationOfState::rho(Pstd, T)
143  - s(Pstd, T)*T;
144 }
145 
146 
147 template<class EquationOfState>
149 (
150  const scalar p,
151  const scalar T
152 ) const
153 {
155  return 0; // EquationOfState::dCpdT
156 }
157 
158 
159 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
160 
161 template<class EquationOfState>
163 (
165 )
166 {
167  scalar Y1 = this->Y();
168 
169  EquationOfState::operator+=(ct);
170 
171  if (mag(this->Y()) > small)
172  {
173  if
174  (
176  && notEqual(Tref_, ct.Tref_)
177  )
178  {
180  << "Tref " << Tref_ << " for "
181  << (this->name().size() ? this->name() : "others")
182  << " != " << ct.Tref_ << " for "
183  << (ct.name().size() ? ct.name() : "others")
184  << exit(FatalError);
185  }
186 
187  Y1 /= this->Y();
188  const scalar Y2 = ct.Y()/this->Y();
189 
190  Cv_ = Y1*Cv_ + Y2*ct.Cv_;
191  hf_ = Y1*hf_ + Y2*ct.hf_;
192  esRef_ = Y1*esRef_ + Y2*ct.esRef_;
193  }
194 }
195 
196 
197 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
198 
199 template<class EquationOfState>
200 inline Foam::eConstThermo<EquationOfState> Foam::operator+
201 (
204 )
205 {
206  EquationOfState eofs
207  (
208  static_cast<const EquationOfState&>(ct1)
209  + static_cast<const EquationOfState&>(ct2)
210  );
211 
212  if (mag(eofs.Y()) < small)
213  {
215  (
216  eofs,
217  ct1.Cv_,
218  ct1.hf_,
219  ct1.Tref_,
220  ct1.esRef_
221  );
222  }
223  else
224  {
225  if
226  (
227  eConstThermo<EquationOfState>::debug
228  && notEqual(ct1.Tref_, ct2.Tref_)
229  )
230  {
232  << "Tref " << ct1.Tref_ << " for "
233  << (ct1.name().size() ? ct1.name() : "others")
234  << " != " << ct2.Tref_ << " for "
235  << (ct2.name().size() ? ct2.name() : "others")
236  << exit(FatalError);
237  }
238 
239  return eConstThermo<EquationOfState>
240  (
241  eofs,
242  ct1.Y()/eofs.Y()*ct1.Cv_
243  + ct2.Y()/eofs.Y()*ct2.Cv_,
244  ct1.Y()/eofs.Y()*ct1.hf_
245  + ct2.Y()/eofs.Y()*ct2.hf_,
246  ct1.Tref_,
247  ct1.Y()/eofs.Y()*ct1.esRef_
248  + ct2.Y()/eofs.Y()*ct2.esRef_
249  );
250  }
251 }
252 
253 
254 template<class EquationOfState>
255 inline Foam::eConstThermo<EquationOfState> Foam::operator*
256 (
257  const scalar s,
258  const eConstThermo<EquationOfState>& ct
259 )
260 {
261  return eConstThermo<EquationOfState>
262  (
263  s*static_cast<const EquationOfState&>(ct),
264  ct.Cv_,
265  ct.hf_,
266  ct.Tref_,
267  ct.esRef_
268  );
269 }
270 
271 
272 template<class EquationOfState>
273 inline Foam::eConstThermo<EquationOfState> Foam::operator==
274 (
275  const eConstThermo<EquationOfState>& ct1,
276  const eConstThermo<EquationOfState>& ct2
277 )
278 {
279  EquationOfState eofs
280  (
281  static_cast<const EquationOfState&>(ct1)
282  == static_cast<const EquationOfState&>(ct2)
283  );
284 
285  if
286  (
287  eConstThermo<EquationOfState>::debug
288  && notEqual(ct1.Tref_, ct2.Tref_)
289  )
290  {
292  << "Tref " << ct1.Tref_ << " for "
293  << (ct1.name().size() ? ct1.name() : "others")
294  << " != " << ct2.Tref_ << " for "
295  << (ct2.name().size() ? ct2.name() : "others")
296  << exit(FatalError);
297  }
298 
299  return eConstThermo<EquationOfState>
300  (
301  eofs,
302  ct2.Y()/eofs.Y()*ct2.Cv_
303  - ct1.Y()/eofs.Y()*ct1.Cv_,
304  ct2.Y()/eofs.Y()*ct2.hf_
305  - ct1.Y()/eofs.Y()*ct1.hf_,
306  ct1.Tref_,
307  ct2.Y()/eofs.Y()*ct2.esRef_
308  - ct1.Y()/eofs.Y()*ct1.esRef_
309  );
310 }
311 
312 
313 // ************************************************************************* //
scalar Cp(const scalar p, const scalar T) const
Definition: EtoHthermo.H:2
scalar es(const scalar p, const scalar T) const
Definition: HtoEthermo.H:11
scalar Cv(const scalar p, const scalar T) const
Definition: HtoEthermo.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
Internal energy based thermodynamics package using a constant heat capacity at constant volume.
Definition: eConstThermo.H:123
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: eConstThermoI.H:76
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 ea(const scalar p, const scalar T) const
Absolute internal energy [J/kg].
scalar Cv(const scalar p, const scalar T) const
Heat capacity at constant volume [J/kg/K].
Definition: eConstThermoI.H:86
scalar hf() const
Enthalpy of formation [J/kg].
scalar es(const scalar p, const scalar T) const
Sensible internal energy [J/kg].
Definition: eConstThermoI.H:97
autoPtr< eConstThermo > clone() const
Construct and return a clone.
Definition: eConstThermoI.H:63
eConstThermo(const EquationOfState &st, const scalar Cv, const scalar hf, const scalar Tref, const scalar esRef)
Construct from components.
Definition: eConstThermoI.H:30
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 e
Elementary charge.
const dimensionedScalar Pstd
Standard pressure.
const dimensionedScalar Tstd
Standard temperature.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
bool notEqual(const Scalar s1, const Scalar s2)
Definition: Scalar.H:216
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
dimensionedScalar log(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
error FatalError
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
volScalarField & p
PtrList< volScalarField > & Y