All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
hRefConstThermoI.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) 2015-2017 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 cp,
33  const scalar hf,
34  const scalar tref,
35  const scalar href
36 )
37 :
38  EquationOfState(st),
39  Cp_(cp),
40  Hf_(hf),
41  Tref_(tref),
42  Href_(href)
43 {}
44 
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 template<class EquationOfState>
50 (
51  const word& name,
52  const hRefConstThermo& ct
53 )
54 :
55  EquationOfState(name, ct),
56  Cp_(ct.Cp_),
57  Hf_(ct.Hf_),
58  Tref_(ct.Tref_),
59  Href_(ct.Href_)
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 Cp_ + EquationOfState::Cp(p, T);
105 }
106 
107 
108 template<class EquationOfState>
110 (
111  const scalar p, const scalar T
112 ) const
113 {
114  return Cp_*(T-Tref_) + Href_ + Hf_ + EquationOfState::H(p, T);
115 }
116 
117 
118 template<class EquationOfState>
120 (
121  const scalar p, const scalar T
122 ) const
123 {
124  return Cp_*(T-Tref_) + Href_ + EquationOfState::H(p, T);
125 }
126 
127 
128 template<class EquationOfState>
130 {
131  return Hf_;
132 }
133 
134 
135 template<class EquationOfState>
137 (
138  const scalar p, const scalar T
139 ) const
140 {
141  return Cp_*log(T/Tstd) + EquationOfState::S(p, T);
142 }
143 
144 
145 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
146 
147 template<class EquationOfState>
148 inline void Foam::hRefConstThermo<EquationOfState>::operator+=
149 (
151 )
152 {
153  scalar Y1 = this->Y();
154 
155  EquationOfState::operator+=(ct);
156 
157  if (mag(this->Y()) > SMALL)
158  {
159  Y1 /= this->Y();
160  const scalar Y2 = ct.Y()/this->Y();
161 
162  Cp_ = Y1*Cp_ + Y2*ct.Cp_;
163  Hf_ = Y1*Hf_ + Y2*ct.Hf_;
164  }
165 }
166 
167 
168 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
169 
170 template<class EquationOfState>
171 inline Foam::hRefConstThermo<EquationOfState> Foam::operator+
172 (
175 )
176 {
177  EquationOfState eofs
178  (
179  static_cast<const EquationOfState&>(ct1)
180  + static_cast<const EquationOfState&>(ct2)
181  );
182 
183  if (mag(eofs.Y()) < SMALL)
184  {
186  (
187  eofs,
188  ct1.Cp_,
189  ct1.Hf_,
190  ct1.Tref_,
191  ct1.Href_
192  );
193  }
194  else
195  {
197  (
198  eofs,
199  ct1.Y()/eofs.Y()*ct1.Cp_
200  + ct2.Y()/eofs.Y()*ct2.Cp_,
201  ct1.Y()/eofs.Y()*ct1.Hf_
202  + ct2.Y()/eofs.Y()*ct2.Hf_,
203  ct1.Y()/eofs.Y()*ct1.Tref_
204  + ct2.Y()/eofs.Y()*ct2.Tref_,
205  ct1.Y()/eofs.Y()*ct1.Href_
206  + ct2.Y()/eofs.Y()*ct2.Href_
207  );
208  }
209 }
210 
211 
212 template<class EquationOfState>
213 inline Foam::hRefConstThermo<EquationOfState> Foam::operator*
214 (
215  const scalar s,
217 )
218 {
220  (
221  s*static_cast<const EquationOfState&>(ct),
222  ct.Cp_,
223  ct.Hf_,
224  ct.Tref_,
225  ct.Href_
226  );
227 }
228 
229 
230 template<class EquationOfState>
231 inline Foam::hRefConstThermo<EquationOfState> Foam::operator==
232 (
235 )
236 {
237  EquationOfState eofs
238  (
239  static_cast<const EquationOfState&>(ct1)
240  == static_cast<const EquationOfState&>(ct2)
241  );
242 
244  (
245  eofs,
246  ct2.Y()/eofs.Y()*ct2.Cp_
247  - ct1.Y()/eofs.Y()*ct1.Cp_,
248  ct2.Y()/eofs.Y()*ct2.Hf_
249  - ct1.Y()/eofs.Y()*ct1.Hf_
250  );
251 }
252 
253 
254 // ************************************************************************* //
PtrList< volScalarField > & Y1
Definition: YEqns.H:8
dictionary dict
scalar S(const scalar p, const scalar T) const
Entropy [J/(kg K)].
dimensionedScalar log(const dimensionedScalar &ds)
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
scalar Hs(const scalar p, const scalar T) const
Sensible enthalpy [J/kg].
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].
PtrList< volScalarField > & Y2
Definition: YEqns.H:9
static autoPtr< hRefConstThermo > New(const dictionary &dict)
Selector from dictionary.
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 limit(const scalar T) const
Limit the temperature to be in the range Tlow_ to Thigh_.
const dimensionedScalar Tstd
Standard temperature.
const volScalarField & cp
scalar Hc() const
Chemical enthalpy [J/kg].
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
PtrList< volScalarField > & Y
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.
autoPtr< hRefConstThermo > clone() const
Construct and return a clone.