All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
incompressiblePerfectGas.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) 2012-2020 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 Class
25  Foam::incompressiblePerfectGas
26 
27 Description
28  Incompressible gas equation of state using a constant reference pressure in
29  the perfect gas equation of state rather than the local pressure so that the
30  density only varies with temperature and composition.
31 
32 SourceFiles
33  incompressiblePerfectGasI.H
34  incompressiblePerfectGas.C
35 
36 \*---------------------------------------------------------------------------*/
37 
38 #ifndef incompressiblePerfectGas_H
39 #define incompressiblePerfectGas_H
40 
41 #include "autoPtr.H"
42 
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
44 
45 namespace Foam
46 {
47 
48 // Forward declaration of friend functions and operators
49 
50 template<class Specie> class incompressiblePerfectGas;
51 
52 template<class Specie>
53 inline incompressiblePerfectGas<Specie> operator+
54 (
56  const incompressiblePerfectGas<Specie>&
57 );
58 
59 template<class Specie>
60 inline incompressiblePerfectGas<Specie> operator*
61 (
62  const scalar,
63  const incompressiblePerfectGas<Specie>&
64 );
65 
66 template<class Specie>
67 inline incompressiblePerfectGas<Specie> operator==
68 (
69  const incompressiblePerfectGas<Specie>&,
70  const incompressiblePerfectGas<Specie>&
71 );
72 
73 template<class Specie>
74 Ostream& operator<<
75 (
76  Ostream&,
77  const incompressiblePerfectGas<Specie>&
78 );
79 
80 
81 /*---------------------------------------------------------------------------*\
82  Class incompressiblePerfectGas Declaration
83 \*---------------------------------------------------------------------------*/
84 
85 template<class Specie>
87 :
88  public Specie
89 {
90  // Private Data
91 
92  //- Reference pressure
93  scalar pRef_;
94 
95 
96 public:
97 
98  // Constructors
99 
100  //- Construct from components
101  inline incompressiblePerfectGas(const Specie& sp, const scalar pRef);
102 
103  //- Construct from dictionary
105 
106  //- Construct as named copy
108  (
109  const word& name,
111  );
112 
113  //- Construct and return a clone
115 
116  // Selector from dictionary
118  (
119  const dictionary& dict
120  );
121 
122 
123  // Member Functions
124 
125  //- Return the instantiated type name
126  static word typeName()
127  {
128  return
129  "incompressiblePerfectGas<"
130  + word(Specie::typeName_()) + '>';
131  }
132 
133 
134  // Fundamental properties
135 
136  //- Is the equation of state is incompressible i.e. rho != f(p)
137  static const bool incompressible = true;
138 
139  //- Is the equation of state is isochoric i.e. rho = const
140  static const bool isochoric = false;
141 
142  //- Return density [kg/m^3]
143  inline scalar rho(scalar p, scalar T) const;
144 
145  //- Return enthalpy contribution [J/kg]
146  inline scalar H(const scalar p, const scalar T) const;
147 
148  //- Return Cp contribution [J/(kg K]
149  inline scalar Cp(scalar p, scalar T) const;
150 
151  //- Return internal energy contribution [J/kg]
152  inline scalar E(const scalar p, const scalar T) const;
153 
154  //- Return Cv contribution [J/(kg K]
155  inline scalar Cv(scalar p, scalar T) const;
156 
157  //- Return entropy contribution to the integral of Cp/T [J/kg/K]
158  inline scalar Sp(const scalar p, const scalar T) const;
159 
160  //- Return entropy contribution to the integral of Cv/T [J/kg/K]
161  inline scalar Sv(const scalar p, const scalar T) const;
162 
163  //- Return compressibility [s^2/m^2]
164  inline scalar psi(scalar p, scalar T) const;
165 
166  //- Return compression factor []
167  inline scalar Z(scalar p, scalar T) const;
168 
169  //- Return (Cp - Cv) [J/(kg K]
170  inline scalar CpMCv(scalar p, scalar T) const;
171 
172 
173  // IO
174 
175  //- Write to Ostream
176  void write(Ostream& os) const;
177 
178 
179  // Member Operators
180 
181  inline void operator+=(const incompressiblePerfectGas&);
182  inline void operator*=(const scalar);
183 
184 
185  // Friend operators
186 
187  friend incompressiblePerfectGas operator+ <Specie>
188  (
190  const incompressiblePerfectGas&
191  );
192 
193  friend incompressiblePerfectGas operator* <Specie>
194  (
195  const scalar s,
196  const incompressiblePerfectGas&
197  );
198 
199  friend incompressiblePerfectGas operator== <Specie>
200  (
201  const incompressiblePerfectGas&,
202  const incompressiblePerfectGas&
203  );
204 
205 
206  // Ostream Operator
207 
208  friend Ostream& operator<< <Specie>
209  (
210  Ostream&,
211  const incompressiblePerfectGas&
212  );
213 };
214 
215 
216 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
217 
218 } // End namespace Foam
219 
220 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
221 
223 
224 #ifdef NoRepository
225  #include "incompressiblePerfectGas.C"
226 #endif
227 
228 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
229 
230 #endif
231 
232 // ************************************************************************* //
void write(Ostream &os) const
Write to Ostream.
dictionary dict
scalar Cp(scalar p, scalar T) const
Return Cp contribution [J/(kg K].
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
Incompressible gas equation of state using a constant reference pressure in the perfect gas equation ...
void operator+=(const incompressiblePerfectGas &)
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
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 Z(scalar p, scalar T) const
Return compression factor [].
A class for handling words, derived from string.
Definition: word.H:59
incompressiblePerfectGas(const Specie &sp, const scalar pRef)
Construct from components.
static const bool isochoric
Is the equation of state is isochoric i.e. rho = const.
scalar E(const scalar p, const scalar T) const
Return internal energy contribution [J/kg].
static const bool incompressible
Is the equation of state is incompressible i.e. rho != f(p)
scalar Sp(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cp/T [J/kg/K].
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
scalar Cv(scalar p, scalar T) const
Return Cv contribution [J/(kg K].
scalar pRef
Definition: createFields.H:19
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
scalar H(const scalar p, const scalar T) const
Return enthalpy contribution [J/kg].
scalar psi(scalar p, scalar T) const
Return compressibility [s^2/m^2].
scalar Sv(const scalar p, const scalar T) const
Return entropy contribution to the integral of Cv/T [J/kg/K].
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
volScalarField & p
autoPtr< incompressiblePerfectGas > clone() const
Construct and return a clone.
static autoPtr< incompressiblePerfectGas > New(const dictionary &dict)
static word typeName()
Return the instantiated type name.
Namespace for OpenFOAM.