BoussinesqI.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-2021 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 
26 #include "Boussinesq.H"
27 
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class Specie>
33 (
34  const Specie& sp, const scalar rho0, const scalar T0, const scalar beta
35 )
36 :
37  Specie(sp),
38  rho0_(rho0),
39  T0_(T0),
40  beta_(beta)
41 {}
42 
43 
44 template<class Specie>
46 (
47  const word& name,
48  const Boussinesq<Specie>& b
49 )
50 :
51  Specie(name, b),
52  rho0_(b.rho0_),
53  T0_(b.T0_),
54  beta_(b.beta_)
55 {}
56 
57 
58 template<class Specie>
61 {
63  (
64  new Boussinesq<Specie>(*this)
65  );
66 }
67 
68 
69 template<class Specie>
72 (
73  const dictionary& dict
74 )
75 {
77  (
79  );
80 }
81 
82 
83 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
84 
85 template<class Specie>
86 inline Foam::scalar Foam::Boussinesq<Specie>::rho
87 (
88  scalar p,
89  scalar T
90 ) const
91 {
92  return rho0_*(1.0 - beta_*(T - T0_));
93 }
94 
95 
96 template<class Specie>
97 inline Foam::scalar Foam::Boussinesq<Specie>::H(scalar p, scalar T) const
98 {
99  return p/this->rho(p, T);
100 }
101 
102 
103 template<class Specie>
104 inline Foam::scalar Foam::Boussinesq<Specie>::Cp(scalar p, scalar T) const
105 {
106  return 0;
107 }
108 
109 
110 template<class Specie>
111 inline Foam::scalar Foam::Boussinesq<Specie>::E(scalar p, scalar T) const
112 {
113  return 0;
114 }
115 
116 
117 template<class Specie>
118 inline Foam::scalar Foam::Boussinesq<Specie>::Cv(scalar p, scalar T) const
119 {
120  return 0;
121 }
122 
123 
124 template<class Specie>
125 inline Foam::scalar Foam::Boussinesq<Specie>::Sp
126 (
127  scalar p,
128  scalar T
129 ) const
130 {
131  return 0;
132 }
133 
134 
135 template<class Specie>
136 inline Foam::scalar Foam::Boussinesq<Specie>::Sv
137 (
138  scalar p,
139  scalar T
140 ) const
141 {
143  return 0;
144 }
145 
146 
147 template<class Specie>
148 inline Foam::scalar Foam::Boussinesq<Specie>::psi
149 (
150  scalar p,
151  scalar T
152 ) const
153 {
154  return 0;
155 }
156 
157 
158 template<class Specie>
159 inline Foam::scalar Foam::Boussinesq<Specie>::Z
160 (
161  scalar p,
162  scalar T
163 ) const
164 {
165  return p/(rho(p, T)*this->R()*T);
166 }
167 
168 
169 template<class Specie>
170 inline Foam::scalar Foam::Boussinesq<Specie>::CpMCv
171 (
172  scalar p,
173  scalar T
174 ) const
175 {
176  return 0;
177 }
178 
179 
180 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
181 
182 template<class Specie>
183 inline void Foam::Boussinesq<Specie>::operator+=
184 (
185  const Boussinesq<Specie>& b
186 )
187 {
189 }
190 
191 
192 template<class Specie>
193 inline void Foam::Boussinesq<Specie>::operator*=(const scalar s)
194 {
195  Specie::operator*=(s);
196 }
197 
198 
199 // * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
200 
201 template<class Specie>
202 inline Foam::Boussinesq<Specie> Foam::operator+
203 (
204  const Boussinesq<Specie>& b1,
205  const Boussinesq<Specie>& b2
206 )
207 {
209  return b1;
210 }
211 
212 
213 template<class Specie>
214 inline Foam::Boussinesq<Specie> Foam::operator*
215 (
216  const scalar s,
217  const Boussinesq<Specie>& b
218 )
219 {
220  return Boussinesq<Specie>
221  (
222  s*static_cast<const Specie&>(b),
223  b.rho0_,
224  b.T0_,
225  b.beta_
226  );
227 }
228 
229 
230 template<class Specie>
231 inline Foam::Boussinesq<Specie> Foam::operator==
232 (
233  const Boussinesq<Specie>& b1,
234  const Boussinesq<Specie>& b2
235 )
236 {
238  return b1;
239 }
240 
241 
242 // ************************************************************************* //
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
const dimensionedScalar b
Wien displacement law constant: default SI units: [m K].
Definition: createFields.H:27
scalar rho(scalar p, scalar T) const
Return density [kg/m^3].
Definition: BoussinesqI.H:87
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
#define noCoefficientMixing(Type)
Definition: specie.H:159
const volScalarField & T
#define R(A, B, C, D, E, F, K, M)
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
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:370
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