incompressiblePerfectGasI.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-2015 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 
27 
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class Specie>
33 (
34  const Specie& sp, const scalar pRef
35 )
36 :
37  Specie(sp),
38  pRef_(pRef)
39 {}
40 
41 
42 template<class Specie>
44 (
45  const incompressiblePerfectGas& ipg
46 )
47 :
48  Specie(ipg),
49  pRef_(ipg.pRef_)
50 {}
51 
52 
53 template<class Specie>
55 (
56  const word& name,
58 )
59 :
60  Specie(name, ipg),
61  pRef_(ipg.pRef_)
62 {}
63 
64 
65 template<class Specie>
68 {
70  (
72  );
73 }
74 
75 
76 template<class Specie>
79 (
80  Istream& is
81 )
82 {
84  (
86  );
87 }
88 
89 
90 template<class Specie>
93 (
94  const dictionary& dict
95 )
96 {
98  (
100  );
101 }
102 
103 
104 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
105 
106 template<class Specie>
108 (
109  scalar p,
110  scalar T
111 ) const
112 {
113  return pRef_/(this->R()*T);
114 }
115 
116 
117 template<class Specie>
119 (
120  scalar p,
121  scalar T
122 ) const
123 {
124  return 0;
125 }
126 
127 
128 template<class Specie>
130 (
131  scalar p,
132  scalar T
133 ) const
134 {
135  return 0;
136 }
137 
138 
139 template<class Specie>
141 (
142  scalar p,
143  scalar T
144 ) const
145 {
146  return 0;
147 }
148 
149 
150 template<class Specie>
152 (
153  scalar p,
154  scalar T
155 ) const
156 {
157  return RR;
158 }
159 
160 
161 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
162 
163 template<class Specie>
165 Foam::incompressiblePerfectGas<Specie>::operator=
166 (
168 )
169 {
170  Specie::operator=(ipg);
171 
172  pRef_ = ipg.pRef_;
173 
174  return *this;
175 }
176 
177 template<class Specie>
178 inline void Foam::incompressiblePerfectGas<Specie>::operator+=
179 (
181 )
182 {
183  scalar molr1 = this->nMoles();
184  Specie::operator+=(ipg);
185  molr1 /= this->nMoles();
186  scalar molr2 = ipg.nMoles()/this->nMoles();
187 
188  pRef_ = molr1*pRef_ + molr2*ipg.pRef_;
189 }
190 
191 
192 template<class Specie>
193 inline void Foam::incompressiblePerfectGas<Specie>::operator-=
194 (
196 )
197 {
198  Specie::operator-=(ipg);
199  pRef_ = ipg.pRef_;
200 }
201 
202 
203 template<class Specie>
205 {
206  Specie::operator*=(s);
207 }
208 
209 
210 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
211 
212 template<class Specie>
213 inline Foam::incompressiblePerfectGas<Specie> Foam::operator+
214 (
217 )
218 {
219  scalar nMoles = ipg1.nMoles() + ipg2.nMoles();
220  scalar molr1 = ipg1.nMoles()/nMoles;
221  scalar molr2 = ipg2.nMoles()/nMoles;
222 
224  (
225  static_cast<const Specie&>(ipg1)
226  + static_cast<const Specie&>(ipg2),
227  molr1*ipg1.pRef_ + molr2*ipg2.pRef_
228  );
229 }
230 
231 
232 template<class Specie>
233 inline Foam::incompressiblePerfectGas<Specie> Foam::operator-
234 (
237 )
238 {
240  (
241  static_cast<const Specie&>(ipg1)
242  - static_cast<const Specie&>(ipg2),
243  ipg1.pRef_
244  );
245 }
246 
247 
248 template<class Specie>
249 inline Foam::incompressiblePerfectGas<Specie> Foam::operator*
250 (
251  const scalar s,
253 )
254 {
256  (
257  s*static_cast<const Specie&>(ipg),
258  ipg.pRef_
259  );
260 }
261 
262 
263 template<class Specie>
264 inline Foam::incompressiblePerfectGas<Specie> Foam::operator==
265 (
268 )
269 {
270  return pg2 - pg1;
271 }
272 
273 
274 // ************************************************************************* //
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 ))
scalar cpMcv(scalar p, scalar T) const
Return (cp - cv) [J/(kmol K].
#define R(A, B, C, D, E, F, K, M)
An Istream is an abstract base class for all input systems (streams, files, token lists etc)...
Definition: Istream.H:57
A class for handling words, derived from string.
Definition: word.H:59
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
scalar psi(scalar p, scalar T) const
Return compressibility rho/p [s^2/m^2].
const volScalarField & T
Definition: createFields.H:25
scalar Z(scalar p, scalar T) const
Return compression factor [].
dictionary dict
incompressiblePerfectGas(const Specie &sp, const scalar pRef)
Construct from components.
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
scalar s(const scalar p, const scalar T) const
Return entropy [J/(kmol K)].
static autoPtr< incompressiblePerfectGas > New(Istream &is)
const scalar RR
Universal gas constant (default in [J/(kmol K)])
Incompressible gas equation of state using a constant reference pressure in the perfect gas equation ...
autoPtr< incompressiblePerfectGas > clone() const
Construct and return a clone.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:117