eRefConstThermoI.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) 2018 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 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
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 eref
36 )
37 :
38  EquationOfState(st),
39  Cv_(cv),
40  Hf_(hf),
41  Tref_(tref),
42  Eref_(eref)
43 {}
44 
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 template<class EquationOfState>
50 (
51  const word& name,
52  const eRefConstThermo& ct
53 )
54 :
55  EquationOfState(name, ct),
56  Cv_(ct.Cv_),
57  Hf_(ct.Hf_),
58  Tref_(ct.Tref_),
59  Eref_(ct.Eref_)
60 {}
61 
62 
63 template<class EquationOfState>
66 {
68  (
70  );
71 }
72 
73 
74 template<class EquationOfState>
77 {
79  (
81  );
82 }
83 
84 
85 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
86 
87 template<class EquationOfState>
89 (
90  const scalar T
91 ) const
92 {
93  return T;
94 }
95 
96 
97 template<class EquationOfState>
99 (
100  const scalar p,
101  const scalar T
102 ) const
103 {
104  return Cv_ + EquationOfState::Cv(p, T);
105 }
106 
107 
108 template<class EquationOfState>
110 (
111  const scalar p, const scalar T
112 ) const
113 {
114  return Cv_*(T - Tref_) + Eref_ + EquationOfState::E(p, T);
115 }
116 
117 
118 template<class EquationOfState>
120 {
121  return Hf_;
122 }
123 
124 
125 template<class EquationOfState>
127 (
128  const scalar p, const scalar T
129 ) const
130 {
131  return Es(p, T) + Hc();
132 }
133 
134 
135 template<class EquationOfState>
137 (
138  const scalar p, const scalar T
139 ) const
140 {
141  return Cp(p, T)*log(T/Tstd) + EquationOfState::S(p, T);
142 }
143 
144 
145 template<class EquationOfState>
147 (
148  const scalar p, const scalar T
149 ) const
150 {
151  return 0;
152 }
153 
154 
155 template<class EquationOfState>
157 (
158  const scalar p, const scalar T
159 ) const
160 {
161  return 0;
162 }
163 
164 
165 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
166 
167 template<class EquationOfState>
168 inline void Foam::eRefConstThermo<EquationOfState>::operator+=
169 (
171 )
172 {
173  scalar Y1 = this->Y();
174 
175  EquationOfState::operator+=(ct);
176 
177  if (mag(this->Y()) > small)
178  {
179  Y1 /= this->Y();
180  const scalar Y2 = ct.Y()/this->Y();
181 
182  Cv_ = Y1*Cv_ + Y2*ct.Cv_;
183  Hf_ = Y1*Hf_ + Y2*ct.Hf_;
184  }
185 }
186 
187 
188 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
189 
190 template<class EquationOfState>
191 inline Foam::eRefConstThermo<EquationOfState> Foam::operator+
192 (
195 )
196 {
197  EquationOfState eofs
198  (
199  static_cast<const EquationOfState&>(ct1)
200  + static_cast<const EquationOfState&>(ct2)
201  );
202 
203  if (mag(eofs.Y()) < small)
204  {
206  (
207  eofs,
208  ct1.Cv_,
209  ct1.Hf_,
210  ct1.Tref_,
211  ct1.Eref_
212  );
213  }
214  else
215  {
217  (
218  eofs,
219  ct1.Y()/eofs.Y()*ct1.Cv_
220  + ct2.Y()/eofs.Y()*ct2.Cv_,
221  ct1.Y()/eofs.Y()*ct1.Hf_
222  + ct2.Y()/eofs.Y()*ct2.Hf_,
223  ct1.Y()/eofs.Y()*ct1.Tref_
224  + ct2.Y()/eofs.Y()*ct2.Tref_,
225  ct1.Y()/eofs.Y()*ct1.Eref_
226  + ct2.Y()/eofs.Y()*ct2.Eref_
227  );
228  }
229 }
230 
231 
232 template<class EquationOfState>
233 inline Foam::eRefConstThermo<EquationOfState> Foam::operator*
234 (
235  const scalar s,
237 )
238 {
240  (
241  s*static_cast<const EquationOfState&>(ct),
242  ct.Cv_,
243  ct.Hf_,
244  ct.Tref_,
245  ct.Eref_
246  );
247 }
248 
249 
250 template<class EquationOfState>
251 inline Foam::eRefConstThermo<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  ct2.Y()/eofs.Y()*ct2.Cv_
267  - ct1.Y()/eofs.Y()*ct1.Cv_,
268  ct2.Y()/eofs.Y()*ct2.Hf_
269  - ct1.Y()/eofs.Y()*ct1.Hf_
270  );
271 }
272 
273 
274 // ************************************************************************* //
scalar Cv(const scalar p, const scalar T) const
Heat capacity at constant volume [J/kg/K].
scalar Ea(const scalar p, const scalar T) const
Absolute internal energy [J/kg].
scalar limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
scalar Cv(const scalar p, const scalar T) const
Definition: HtoEthermo.H:2
dictionary dict
scalar dGdT(const scalar p, const scalar T) const
Derivative of Gibbs free energy w.r.t. temperature.
dimensionedScalar log(const dimensionedScalar &ds)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
scalar dCpdT(const scalar p, const scalar T) const
Temperature derivative of heat capacity at constant pressure.
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))
A class for handling words, derived from string.
Definition: word.H:59
scalar S(const scalar p, const scalar T) const
Entropy [J/kg/K].
const dimensionedScalar Tstd
Standard temperature.
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
scalar Hc() const
Chemical enthalpy [J/kg].
PtrList< volScalarField > & Y
static autoPtr< eRefConstThermo > New(const dictionary &dict)
Selector from dictionary.
autoPtr< eRefConstThermo > clone() const
Construct and return a clone.
scalar Es(const scalar p, const scalar T) const
Sensible internal energy [J/kg].
scalar Cp(const scalar p, const scalar T) const
Definition: EtoHthermo.H:2
dimensioned< scalar > mag(const dimensioned< Type > &)
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
Constant properties thermodynamics package templated into the EquationOfState.