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