perfectFluidI.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) 2012-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 "perfectFluid.H"
27 #include "specie.H"
28 
29 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
30 
31 template<class Specie>
33 (
34  const Specie& sp,
35  const scalar R,
36  const scalar rho0
37 )
38 :
39  Specie(sp),
40  R_(R),
41  rho0_(rho0)
42 {}
43 
44 
45 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
46 
47 template<class Specie>
49 (
50  const word& name,
51  const perfectFluid<Specie>& pf
52 )
53 :
54  Specie(name, pf),
55  R_(pf.R_),
56  rho0_(pf.rho0_)
57 {}
58 
59 
60 template<class Specie>
63 {
65 }
66 
67 
68 template<class Specie>
71 {
73 }
74 
75 
76 template<class Specie>
79 (
80  const dictionary& dict
81 )
82 {
84 }
85 
86 
87 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
88 
89 template<class Specie>
90 inline Foam::scalar Foam::perfectFluid<Specie>::R() const
91 {
92  return R_;
93 }
94 
95 
96 template<class Specie>
97 inline Foam::scalar Foam::perfectFluid<Specie>::rho(scalar p, scalar T) const
98 {
99  return rho0_ + p/(this->R()*T);
100 }
101 
102 
103 template<class Specie>
104 inline Foam::scalar Foam::perfectFluid<Specie>::h(scalar p, scalar T) const
105 {
106  return 0;
107 }
108 
109 
110 template<class Specie>
111 inline Foam::scalar Foam::perfectFluid<Specie>::cp(scalar p, scalar T) const
112 {
113  return 0;
114 }
115 
116 
117 template<class Specie>
118 inline Foam::scalar Foam::perfectFluid<Specie>::s(scalar p, scalar T) const
119 {
120  return -RR*log(p/Pstd);
121 }
122 
123 
124 template<class Specie>
125 inline Foam::scalar Foam::perfectFluid<Specie>::psi(scalar p, scalar T) const
126 {
127  return 1.0/(this->R()*T);
128 }
129 
130 
131 template<class Specie>
132 inline Foam::scalar Foam::perfectFluid<Specie>::Z(scalar p, scalar T) const
133 {
134  return 1;
135 }
136 
137 
138 template<class Specie>
139 inline Foam::scalar Foam::perfectFluid<Specie>::cpMcv(scalar p, scalar T) const
140 {
141  return 0;
142 }
143 
144 
145 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
146 
147 template<class Specie>
148 inline void Foam::perfectFluid<Specie>::operator+=
149 (
150  const perfectFluid<Specie>& pf
151 )
152 {
153  scalar molr1 = this->nMoles();
154 
155  Specie::operator+=(pf);
156 
157  molr1 /= this->nMoles();
158  scalar molr2 = pf.nMoles()/this->nMoles();
159 
160  R_ = 1.0/(molr1/R_ + molr2/pf.R_);
161  rho0_ = molr1*rho0_ + molr2*pf.rho0_;
162 }
163 
164 
165 template<class Specie>
166 inline void Foam::perfectFluid<Specie>::operator-=
167 (
168  const perfectFluid<Specie>& pf
169 )
170 {
171  scalar molr1 = this->nMoles();
172 
173  Specie::operator-=(pf);
174 
175  molr1 /= this->nMoles();
176  scalar molr2 = pf.nMoles()/this->nMoles();
177 
178  R_ = 1.0/(molr1/R_ - molr2/pf.R_);
179  rho0_ = molr1*rho0_ - molr2*pf.rho0_;
180 }
181 
182 
183 template<class Specie>
184 inline void Foam::perfectFluid<Specie>::operator*=(const scalar s)
185 {
186  Specie::operator*=(s);
187 }
188 
189 
190 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
191 
192 template<class Specie>
193 inline Foam::perfectFluid<Specie> Foam::operator+
194 (
195  const perfectFluid<Specie>& pf1,
196  const perfectFluid<Specie>& pf2
197 )
198 {
199  scalar nMoles = pf1.nMoles() + pf2.nMoles();
200  scalar molr1 = pf1.nMoles()/nMoles;
201  scalar molr2 = pf2.nMoles()/nMoles;
202 
203  return perfectFluid<Specie>
204  (
205  static_cast<const Specie&>(pf1)
206  + static_cast<const Specie&>(pf2),
207  1.0/(molr1/pf1.R_ + molr2/pf2.R_),
208  molr1*pf1.rho0_ + molr2*pf2.rho0_
209  );
210 }
211 
212 
213 template<class Specie>
214 inline Foam::perfectFluid<Specie> Foam::operator-
215 (
216  const perfectFluid<Specie>& pf1,
217  const perfectFluid<Specie>& pf2
218 )
219 {
220  scalar nMoles = pf1.nMoles() + pf2.nMoles();
221  scalar molr1 = pf1.nMoles()/nMoles;
222  scalar molr2 = pf2.nMoles()/nMoles;
223 
224  return perfectFluid<Specie>
225  (
226  static_cast<const Specie&>(pf1)
227  - static_cast<const Specie&>(pf2),
228  1.0/(molr1/pf1.R_ - molr2/pf2.R_),
229  molr1*pf1.rho0_ - molr2*pf2.rho0_
230  );
231 }
232 
233 
234 template<class Specie>
235 inline Foam::perfectFluid<Specie> Foam::operator*
236 (
237  const scalar s,
238  const perfectFluid<Specie>& pf
239 )
240 {
241  return perfectFluid<Specie>
242  (
243  s*static_cast<const Specie&>(pf),
244  pf.R_,
245  pf.rho0_
246  );
247 }
248 
249 
250 template<class Specie>
251 inline Foam::perfectFluid<Specie> Foam::operator==
252 (
253  const perfectFluid<Specie>& pf1,
254  const perfectFluid<Specie>& pf2
255 )
256 {
257  return pf2 - pf1;
258 }
259 
260 
261 // ************************************************************************* //
dictionary dict
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 R() const
Return fluid constant [J/(kg K)].
Definition: perfectFluidI.H:90
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].
perfectFluid(const Specie &sp, const scalar R, const scalar rho0)
Construct from components.
Definition: perfectFluidI.H:33
scalar h(const scalar p, const scalar T) const
Return enthalpy departure [J/kmol].
static autoPtr< perfectFluid > New(Istream &is)
Definition: perfectFluidI.H:70
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
void operator*=(const scalar)
const dimensionedScalar Pstd
Standard pressure.
scalar Z(scalar p, scalar T) const
Return compression factor [].
const volScalarField & T
Perfect gas equation of state.
Definition: perfectFluid.H:47
scalar s(const scalar p, const scalar T) const
Return entropy [J/(kmol K)].
#define R(A, B, C, D, E, F, K, M)
scalar psi(scalar p, scalar T) const
Return compressibility rho/p [s^2/m^2].
const scalar RR
Universal gas constant (default in [J/(kmol K)])
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:53
volScalarField & p
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
Definition: perfectFluidI.H:97
autoPtr< perfectFluid > clone() const
Construct and return a clone.
Definition: perfectFluidI.H:62
scalar cp(scalar p, scalar T) const
Return cp departure [J/(kmol K].