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