atmBuoyancyTurbSource.C
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) 2020 ENERCON GmbH
9 Copyright (C) 2020-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
31
32// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
33
34namespace Foam
35{
36namespace fv
37{
40}
41}
42
43
44// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
45
46void Foam::fv::atmBuoyancyTurbSource::calcB()
47{
48 //- Temperature field [K]
50
51 //- Kinematic turbulent thermal conductivity field [m2/s]
52 const volScalarField& alphat =
54
55 // (ARAL:Eq. 7)
56 B_ = beta_*alphat()*(fvc::grad(T) & g_)();
57}
58
59
61Foam::fv::atmBuoyancyTurbSource::calcC3
62(
66) const
67{
68 // Gradient Richardson number (ARAL:p. 4)
70 (
71 -B_/(G + dimensionedScalar(G.dimensions(), SMALL))
72 );
73
74 // Mixing-length scale estimation (P:Eq. 10.37 & p. 374) normalised by Lmax_
75 const volScalarField::Internal LbyLmax
76 (
77 (pow(Cmu_, 0.75)/Lmax_)*pow(k, 1.5)/epsilon
78 );
79
80 // (ARAL:Eq. 10), with a typo of (C2_) instead of using (C2_ - 1.0)
81 volScalarField::Internal alphaB(1.0 - LbyLmax);
82
83 alphaB =
84 neg0(Rig)*(1.0 - (1.0 + (C2_ - 1.0)/(C2_ - C1_))*LbyLmax)
85 + pos(Rig)*(1.0 - LbyLmax);
86
87 // (SKL:Eq. 18, rhs-term:3); (ARAL:Eq. 5, rhs-term:3) has a typo
88 return (C1_ - C2_)*alphaB + 1.0;
89}
90
91
93Foam::fv::atmBuoyancyTurbSource::calcC3
94(
96 const volScalarField::Internal& omega,
100) const
101{
102 // Gradient Richardson number (ARAL:p. 4)
104 (
105 -B_/(G + dimensionedScalar(G.dimensions(), SMALL))
106 );
107
108 // Mixing-length scale estimation (L:Eq. 3.20) normalised by Lmax_
109 const volScalarField::Internal LbyLmax
110 (
111 (1.0/(pow025(Cmu_)*Lmax_))*sqrt(k)/omega
112 );
113
114 // (ARAL:Eq. 10)
115 volScalarField::Internal alphaB(1.0 - LbyLmax);
116
117 alphaB =
118 neg0(Rig)*(1.0 - (1.0 + beta/(beta - gamma))*LbyLmax)
119 + pos(Rig)*(1.0 - LbyLmax);
120
121 // (SKL:Eq. 19, rhs-term:3); (ARAL:Eq. 5, rhs-term:3) has a typo
122 return (gamma - beta)*alphaB;
123}
124
125
126// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
127
129(
130 const word& sourceName,
131 const word& modelType,
132 const dictionary& dict,
133 const fvMesh& mesh
134)
135:
136 fv::cellSetOption(sourceName, modelType, dict, mesh),
137 isEpsilon_(true),
138 rhoName_(coeffs_.getOrDefault<word>("rho", "rho")),
139 Lmax_
140 (
142 (
143 dimLength,
144 coeffs_.getCheckOrDefault<scalar>
145 (
146 "Lmax",
147 41.575,
148 [&](const scalar Lmax){ return Lmax > SMALL; }
149 )
150 )
151 ),
152 beta_
153 (
155 (
157 coeffs_.getCheckOrDefault<scalar>
158 (
159 "beta",
160 3.3e-3,
161 [&](const scalar x){ return x > SMALL; }
162 )
163 )
164 ),
165 Cmu_(Zero),
166 C1_(Zero),
167 C2_(Zero),
168 g_
169 (
170 "g",
172 meshObjects::gravity::New(mesh_.time()).value()
173 ),
174 B_
175 (
177 (
178 "B",
179 mesh.time().timeName(),
180 mesh,
183 ),
184 mesh,
186 )
187{
188 const auto* turbPtr =
189 mesh_.findObject<turbulenceModel>
190 (
192 );
193
194 if (!turbPtr)
195 {
197 << "Unable to find a turbulence model."
198 << abort(FatalError);
199 }
200
201 fieldNames_.resize(2);
202
203 tmp<volScalarField> tepsilon = turbPtr->epsilon();
204 tmp<volScalarField> tomega = turbPtr->omega();
205
206 if (!tepsilon.isTmp())
207 {
208 fieldNames_[0] = tepsilon().name();
209
210 const dictionary& turbDict = turbPtr->coeffDict();
211
212 Cmu_.read("Cmu", turbDict);
213 C1_.read("C1", turbDict);
214 C2_.read("C2", turbDict);
215 }
216 else if (!tomega.isTmp())
217 {
218 isEpsilon_ = false;
219 fieldNames_[0] = tomega().name();
220
221 const dictionary& turbDict = turbPtr->coeffDict();
222
223 Cmu_.read("betaStar", turbDict);
224 }
225 else
226 {
228 << "Unable to find neither epsilon nor omega field." << nl
229 << "atmBuoyancyTurbSource needs either epsilon or omega field."
230 << abort(FatalError);
231 }
232
233 fieldNames_[1] = turbPtr->k()().name();
234
236
237 Log << " Applying atmBuoyancyTurbSource to: "
238 << fieldNames_[0] << " and " << fieldNames_[1]
239 << endl;
240}
241
242
243// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
244
246(
247 fvMatrix<scalar>& eqn,
248 const label fieldi
249)
250{
251 if (fieldi == 1)
252 {
253 atmBuoyancyTurbSourceK
254 (
257 eqn,
258 fieldi
259 );
260 return;
261 }
262
263 calcB();
264
265 if (isEpsilon_)
266 {
267 atmBuoyancyTurbSourceEpsilon
268 (
271 eqn,
272 fieldi
273 );
274 }
275 else
276 {
277 atmBuoyancyTurbSourceOmega
278 (
281 eqn,
282 fieldi
283 );
284 }
285}
286
287
289(
290 const volScalarField& rho,
291 fvMatrix<scalar>& eqn,
292 const label fieldi
293)
294{
295 if (fieldi == 1)
296 {
297 atmBuoyancyTurbSourceK(geometricOneField(), rho, eqn, fieldi);
298 return;
299 }
300
301 calcB();
302
303 if (isEpsilon_)
304 {
305 atmBuoyancyTurbSourceEpsilon(geometricOneField(), rho, eqn, fieldi);
306 }
307 else
308 {
309 atmBuoyancyTurbSourceOmega(geometricOneField(), rho, eqn, fieldi);
310 }
311}
312
313
315(
316 const volScalarField& alpha,
317 const volScalarField& rho,
318 fvMatrix<scalar>& eqn,
319 const label fieldi
320)
321{
322 if (fieldi == 1)
323 {
324 atmBuoyancyTurbSourceK(alpha, rho, eqn, fieldi);
325 return;
326 }
327
328 calcB();
329
330 if (isEpsilon_)
331 {
332 atmBuoyancyTurbSourceEpsilon(alpha, rho, eqn, fieldi);
333 }
334 else
335 {
336 atmBuoyancyTurbSourceOmega(alpha, rho, eqn, fieldi);
337 }
338}
339
340
341// ************************************************************************* //
label k
#define Log
Definition: PDRblock.C:35
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
DimensionedField< scalar, volMesh > Internal
The internal field type from which this GeometricField is derived.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
const dimensionSet & dimensions() const
Return const reference to dimensions.
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Applies sources on k and either epsilon or omega to incorporate effects of buoyancy for atmospheric b...
virtual void addSup(fvMatrix< scalar > &eqn, const label fieldi)
Intermediate abstract class for handling cell-set options for the derived fvOptions.
Base abstract class for handling finite volume options (i.e. fvOption).
Definition: fvOption.H:127
const fvMesh & mesh_
Reference to the mesh database.
Definition: fvOption.H:139
void resetApplied()
Resize/reset applied flag list for all fieldNames_ entries.
Definition: fvOption.C:48
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
Type & lookupObjectRef(const word &name, const bool recursive=false) const
A class for managing temporary objects.
Definition: tmp.H:65
Abstract base class for turbulence models (RAS, LES and laminar).
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
const volScalarField & T
const scalar gamma
Definition: EEqn.H:9
dynamicFvMesh & mesh
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
scalar epsilon
const dimensionedScalar G
Newtonian constant of gravitation.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:54
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pos(const dimensionedScalar &ds)
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
Definition: dimensionSets.H:54
dimensionedScalar sqrt(const dimensionedScalar &ds)
errorManip< error > abort(error &err)
Definition: errorManip.H:144
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
error FatalError
dimensionedScalar neg0(const dimensionedScalar &ds)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
dimensionedScalar pow025(const dimensionedScalar &ds)
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
labelList fv(nPoints)
volScalarField & alpha
dictionary dict
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)