adiabaticPerfectFluidI.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) 2013-2016 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 "adiabaticPerfectFluid.H"
27 
28 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
29 
30 template<class Specie>
32 (
33  const Specie& sp,
34  const scalar p0,
35  const scalar rho0,
36  const scalar gamma,
37  const scalar B
38 )
39 :
40  Specie(sp),
41  p0_(p0),
42  rho0_(rho0),
43  gamma_(gamma),
44  B_(B)
45 {}
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
50 template<class Specie>
52 (
53  const word& name,
55 )
56 :
57  Specie(name, pf),
58  p0_(pf.p0_),
59  rho0_(pf.rho0_),
60  gamma_(pf.gamma_),
61  B_(pf.B_)
62 {}
63 
64 
65 template<class Specie>
68 {
70  (
72  );
73 }
74 
75 
76 template<class Specie>
79 {
81  (
83  );
84 }
85 
86 
87 template<class Specie>
90 (
91  const dictionary& dict
92 )
93 {
95  (
97  );
98 }
99 
100 
101 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
102 
103 template<class Specie>
105 (
106  scalar p,
107  scalar T
108 ) const
109 {
110  return rho0_*pow((p + B_)/(p0_ + B_), 1.0/gamma_);
111 }
112 
113 
114 template<class Specie>
115 inline Foam::scalar Foam::adiabaticPerfectFluid<Specie>::h
116 (
117  scalar p,
118  scalar T
119 ) const
120 {
121  return 0;
122 }
123 
124 
125 template<class Specie>
126 inline Foam::scalar Foam::adiabaticPerfectFluid<Specie>::cp
127 (
128  scalar p,
129  scalar T
130 ) const
131 {
132  return 0;
133 }
134 
135 
136 template<class Specie>
137 inline Foam::scalar Foam::adiabaticPerfectFluid<Specie>::s
138 (
139  scalar p,
140  scalar T
141 ) const
142 {
143  scalar n = 1 - 1.0/gamma_;
144  return
145  -pow(p0_ + B_, 1.0/gamma_)*(pow((p + B_), n) - pow((Pstd + B_), n))
146  /(rho0_*T*n);
147 }
148 
149 
150 template<class Specie>
152 (
153  scalar p,
154  scalar T
155 ) const
156 {
157  return
158  (rho0_/(gamma_*(p0_ + B_)))
159  *pow((p + B_)/(p0_ + B_), 1.0/gamma_ - 1.0);
160 }
161 
162 
163 template<class Specie>
164 inline Foam::scalar Foam::adiabaticPerfectFluid<Specie>::Z(scalar, scalar) const
165 {
166  return 1;
167 }
168 
169 
170 template<class Specie>
172 (
173  scalar p,
174  scalar T
175 ) const
176 {
177  return 0;
178 }
179 
180 
181 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
182 
183 template<class Specie>
184 inline void Foam::adiabaticPerfectFluid<Specie>::operator+=
185 (
187 )
188 {
189  scalar molr1 = this->nMoles();
190 
191  Specie::operator+=(pf);
192 
193  molr1 /= this->nMoles();
194  scalar molr2 = pf.nMoles()/this->nMoles();
195 
196  p0_ = molr1*p0_ + molr2*pf.p0_;
197  rho0_ = molr1*rho0_ + molr2*pf.rho0_;
198  gamma_ = molr1*gamma_ + molr2*pf.gamma_;
199  B_ = molr1*B_ + molr2*pf.B_;
200 }
201 
202 
203 template<class Specie>
204 inline void Foam::adiabaticPerfectFluid<Specie>::operator-=
205 (
207 )
208 {
209  scalar molr1 = this->nMoles();
210 
211  Specie::operator-=(pf);
212 
213  molr1 /= this->nMoles();
214  scalar molr2 = pf.nMoles()/this->nMoles();
215 
216  p0_ = molr1*p0_ - molr2*pf.p0_;
217  rho0_ = molr1*rho0_ - molr2*pf.rho0_;
218  gamma_ = molr1*gamma_ - molr2*pf.gamma_;
219  B_ = molr1*B_ - molr2*pf.B_;
220 }
221 
222 
223 template<class Specie>
225 {
226  Specie::operator*=(s);
227 }
228 
229 
230 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
231 
232 template<class Specie>
233 inline Foam::adiabaticPerfectFluid<Specie> Foam::operator+
234 (
237 )
238 {
239  scalar nMoles = pf1.nMoles() + pf2.nMoles();
240  scalar molr1 = pf1.nMoles()/nMoles;
241  scalar molr2 = pf2.nMoles()/nMoles;
242 
243  return rhoConst<Specie>
244  (
245  static_cast<const Specie&>(pf1)
246  + static_cast<const Specie&>(pf2),
247  molr1*pf1.p0_ + molr2*pf2.p0_,
248  molr1*pf1.rho0_ + molr2*pf2.rho0_,
249  molr1*pf1.gamma_ + molr2*pf2.gamma_,
250  molr1*pf1.B_ + molr2*pf2.B_
251  );
252 }
253 
254 
255 template<class Specie>
256 inline Foam::adiabaticPerfectFluid<Specie> Foam::operator-
257 (
260 )
261 {
262  scalar nMoles = pf1.nMoles() + pf2.nMoles();
263  scalar molr1 = pf1.nMoles()/nMoles;
264  scalar molr2 = pf2.nMoles()/nMoles;
265 
266  return rhoConst<Specie>
267  (
268  static_cast<const Specie&>(pf1)
269  - static_cast<const Specie&>(pf2),
270  molr1*pf1.p0_ - molr2*pf2.p0_,
271  molr1*pf1.rho0_ - molr2*pf2.rho0_,
272  molr1*pf1.gamma_ - molr2*pf2.gamma_,
273  molr1*pf1.B_ - molr2*pf2.B_
274  );
275 }
276 
277 
278 template<class Specie>
279 inline Foam::adiabaticPerfectFluid<Specie> Foam::operator*
280 (
281  const scalar s,
283 )
284 {
286  (
287  s*static_cast<const Specie&>(pf),
288  pf.p0_,
289  pf.rho0_,
290  pf.gamma_,
291  pf.B_
292  );
293 }
294 
295 
296 template<class Specie>
297 inline Foam::adiabaticPerfectFluid<Specie> Foam::operator==
298 (
301 )
302 {
303  return pf2 - pf1;
304 }
305 
306 
307 // ************************************************************************* //
dictionary dict
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
adiabaticPerfectFluid(const Specie &sp, const scalar p0, const scalar rho0, const scalar gamma, const scalar B)
Construct from components.
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
scalar cpMcv(scalar p, scalar T) const
Return (cp - cv) [J/(kmol K].
scalar psi(scalar p, scalar T) const
Return compressibility rho/p [s^2/m^2].
scalar Z(scalar p, scalar T) const
Return compression factor [].
static autoPtr< adiabaticPerfectFluid > New(Istream &is)
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
Return entropy [J/(kmol K)].
autoPtr< adiabaticPerfectFluid > clone() const
Construct and return a clone.
const dimensionedScalar Pstd
Standard pressure.
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
scalar cp(scalar p, scalar T) const
Return cp departure [J/(kmol K].
AdiabaticPerfect gas equation of state.
scalar h(const scalar p, const scalar T) const
Return enthalpy departure [J/kmol].
label n
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
RhoConst (rho = const) of state.
Definition: rhoConst.H:47