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