BoussinesqI.H
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2015-2017 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "Boussinesq.H"
29
30
31// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32
33template<class Specie>
35(
36 const Specie& sp, const scalar rho0, const scalar T0, const scalar beta
37)
38:
39 Specie(sp),
40 rho0_(rho0),
41 T0_(T0),
42 beta_(beta)
43{}
44
45
46template<class Specie>
48(
49 const word& name,
51)
52:
53 Specie(name, b),
54 rho0_(b.rho0_),
55 T0_(b.T0_),
56 beta_(b.beta_)
57{}
58
59
60template<class Specie>
63{
65}
66
67
68template<class Specie>
71(
72 const dictionary& dict
73)
74{
76}
77
78
79// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
80
81template<class Specie>
83(
84 scalar p,
85 scalar T
86) const
87{
88 return rho0_*(1.0 - beta_*(T - T0_));
89}
90
91
92template<class Specie>
93inline Foam::scalar Foam::Boussinesq<Specie>::H(scalar p, scalar T) const
94{
95 return 0;
96}
97
98
99template<class Specie>
100inline Foam::scalar Foam::Boussinesq<Specie>::Cp(scalar p, scalar T) const
101{
102 return 0;
103}
104
105template<class Specie>
106inline Foam::scalar Foam::Boussinesq<Specie>::E(scalar p, scalar T) const
107{
108 return 0;
109}
110
111
112template<class Specie>
113inline Foam::scalar Foam::Boussinesq<Specie>::Cv(scalar p, scalar T) const
114{
115 return 0;
116}
117
118
119template<class Specie>
121(
122 scalar p,
123 scalar T
124) const
125{
126 return 0;
127}
128
129
130template<class Specie>
132(
133 scalar p,
134 scalar T
135) const
136{
137 return 0;
138}
139
140
141template<class Specie>
143(
144 scalar p,
145 scalar T
146) const
147{
148 return 0;
149}
150
151
152template<class Specie>
154(
155 scalar p,
156 scalar T
157) const
158{
159 return 0;
160}
161
162
163// * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
164
165template<class Specie>
167(
168 const Boussinesq<Specie>& b
169)
170{
171 Specie::operator=(b);
172
173 rho0_ = b.rho0_;
174 T0_ = b.T0_;
175 beta_ = b.beta_;
176}
177
178
179template<class Specie>
181(
182 const Boussinesq<Specie>& b
183)
184{
185 scalar Y1 = this->Y();
186 Specie::operator+=(b);
187
188 if (mag(this->Y()) > SMALL)
189 {
190 Y1 /= this->Y();
191 const scalar Y2 = b.Y()/this->Y();
192
193 rho0_ = Y1*rho0_ + Y2*b.rho0_;
194 T0_ = Y1*T0_ + Y2*b.T0_;
195 beta_ = Y1*beta_ + Y2*b.beta_;
196 }
197}
198
199
200template<class Specie>
202{
203 Specie::operator*=(s);
204}
205
206
207// * * * * * * * * * * * * * * * Friend Operators * * * * * * * * * * * * * //
208
209template<class Specie>
210inline Foam::Boussinesq<Specie> Foam::operator+
211(
212 const Boussinesq<Specie>& b1,
213 const Boussinesq<Specie>& b2
214)
215{
216 Specie sp(static_cast<const Specie&>(b1) + static_cast<const Specie&>(b2));
217
218 if (mag(sp.Y()) < SMALL)
219 {
220 return Boussinesq<Specie>
221 (
222 sp,
223 b1.rho0_,
224 b1.T0_,
225 b1.beta_
226 );
227 }
228 else
229 {
230 const scalar Y1 = b1.Y()/sp.Y();
231 const scalar Y2 = b2.Y()/sp.Y();
232
233 return Boussinesq<Specie>
234 (
235 sp,
236 Y1*b1.rho0_ + Y2*b2.rho0_,
237 Y1*b1.T0_ + Y2*b2.T0_,
238 Y1*b1.beta_ + Y2*b2.beta_
239 );
240 }
241}
242
243
244template<class Specie>
245inline Foam::Boussinesq<Specie> Foam::operator*
246(
247 const scalar s,
248 const Boussinesq<Specie>& b
249)
250{
251 return Boussinesq<Specie>
252 (
253 s*static_cast<const Specie&>(b),
254 b.rho0_,
255 b.T0_,
256 b.beta_
257 );
258}
259
260
261template<class Specie>
262inline Foam::Boussinesq<Specie> Foam::operator==
263(
264 const Boussinesq<Specie>& b1,
265 const Boussinesq<Specie>& b2
266)
267{
268 Specie sp(static_cast<const Specie&>(b1) == static_cast<const Specie&>(b2));
269
270 const scalar Y1 = b1.Y()/sp.Y();
271 const scalar Y2 = b2.Y()/sp.Y();
272
273 return Boussinesq<Specie>
274 (
275 sp,
276 Y2*b2.rho0_ - Y1*b1.rho0_,
277 Y2*b2.T0_ - Y1*b1.T0_,
278 Y2*b2.beta_ - Y1*b1.beta_
279 );
280}
281
282
283// ************************************************************************* //
Incompressible gas equation of state using the Boussinesq approximation for the density as a function...
Definition: Boussinesq.H:96
autoPtr< Boussinesq > clone() const
Construct and return a clone.
Definition: BoussinesqI.H:62
scalar CpMCv(scalar p, scalar T) const
Return (Cp - Cv) [J/(kg K].
Definition: BoussinesqI.H:154
void operator*=(const scalar)
Definition: BoussinesqI.H:201
compactSpatialTensor S
The joint motion sub-space (3-DoF)
Definition: joint.H:132
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
tmp< GeometricField< Type, faPatchField, areaMesh > > H() const
Return the H operation source.
Definition: faMatrix.C:633
A class for handling words, derived from Foam::string.
Definition: word.H:68
volScalarField & p
PtrList< volScalarField > & Y
const volScalarField & psi
const volScalarField & T
const volScalarField & Cv
Definition: EEqn.H:8
const volScalarField & Cp
Definition: EEqn.H:7
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.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))
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh > > &tdf1, const word &name, const dimensionSet &dimensions)
Global function forwards to reuseTmpDimensionedField::New.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
scalar rho0
dictionary dict
volScalarField & b
Definition: createFields.H:27
scalar T0
Definition: createFields.H:22
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)