qZeta.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) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2019-2020 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
29#include "qZeta.H"
30#include "bound.H"
32
33// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34
35namespace Foam
36{
37namespace incompressible
38{
39namespace RASModels
40{
41
42// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
43
46
47// * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * * //
48
50{
51 const volScalarField Rt(q_*k_/(2.0*nu()*zeta_));
52
53 if (anisotropic_)
54 {
55 return exp((-scalar(2.5) + Rt/20.0)/pow3(scalar(1) + Rt/130.0));
56 }
57 else
58 {
59 return
60 exp(-6.0/sqr(scalar(1) + Rt/50.0))
61 *(scalar(1) + 3.0*exp(-Rt/10.0));
62 }
63}
64
65
67{
68 tmp<volScalarField> Rt = q_*k_/(2.0*nu()*zeta_);
69 return scalar(1) - 0.3*exp(-sqr(Rt));
70}
71
72
74{
77}
78
79
80// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
81
83(
86 const volVectorField& U,
87 const surfaceScalarField& alphaRhoPhi,
89 const transportModel& transport,
90 const word& propertiesName,
91 const word& type
92)
93:
94 eddyViscosity<incompressible::RASModel>
95 (
96 type,
97 alpha,
98 rho,
99 U,
100 alphaRhoPhi,
101 phi,
102 transport,
103 propertiesName
104 ),
105
106 Cmu_
107 (
108 dimensioned<scalar>::getOrAddToDict
109 (
110 "Cmu",
111 coeffDict_,
112 0.09
113 )
114 ),
115 C1_
116 (
117 dimensioned<scalar>::getOrAddToDict
118 (
119 "C1",
120 coeffDict_,
121 1.44
122 )
123 ),
124 C2_
125 (
126 dimensioned<scalar>::getOrAddToDict
127 (
128 "C2",
129 coeffDict_,
130 1.92
131 )
132 ),
133 sigmaZeta_
134 (
135 dimensioned<scalar>::getOrAddToDict
136 (
137 "sigmaZeta",
138 coeffDict_,
139 1.3
140 )
141 ),
142 anisotropic_
143 (
144 Switch::getOrAddToDict
145 (
146 "anisotropic",
147 coeffDict_,
148 false
149 )
150 ),
151
152 qMin_("qMin", sqrt(kMin_)),
153 zetaMin_("zetaMin", epsilonMin_/(2*qMin_)),
154
155 k_
156 (
158 (
159 IOobject::groupName("k", alphaRhoPhi.group()),
160 runTime_.timeName(),
161 mesh_,
162 IOobject::MUST_READ,
163 IOobject::AUTO_WRITE
164 ),
165 mesh_
166 ),
167
168 epsilon_
169 (
171 (
172 IOobject::groupName("epsilon", alphaRhoPhi.group()),
173 runTime_.timeName(),
174 mesh_,
175 IOobject::MUST_READ,
176 IOobject::AUTO_WRITE
177 ),
178 mesh_
179 ),
180
181 q_
182 (
184 (
185 IOobject::groupName("q", alphaRhoPhi.group()),
186 runTime_.timeName(),
187 mesh_,
188 IOobject::NO_READ,
189 IOobject::AUTO_WRITE
190 ),
191 sqrt(bound(k_, kMin_)),
192 k_.boundaryField().types()
193 ),
194
195 zeta_
196 (
198 (
199 IOobject::groupName("zeta", alphaRhoPhi.group()),
200 runTime_.timeName(),
201 mesh_,
202 IOobject::NO_READ,
203 IOobject::AUTO_WRITE
204 ),
205 bound(epsilon_, epsilonMin_)/(2.0*q_),
206 epsilon_.boundaryField().types()
207 )
208{
210
211 if (type == typeName)
212 {
213 printCoeffs(type);
214 }
215}
216
217
218// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
219
221{
223 {
224 Cmu_.readIfPresent(coeffDict());
225 C1_.readIfPresent(coeffDict());
226 C2_.readIfPresent(coeffDict());
227 sigmaZeta_.readIfPresent(coeffDict());
228 anisotropic_.readIfPresent("anisotropic", coeffDict());
229
230 qMin_.readIfPresent(*this);
231 zetaMin_.readIfPresent(*this);
232
233 return true;
234 }
235
236 return false;
237}
238
239
241{
242 if (!turbulence_)
243 {
244 return;
245 }
246
248
249 volScalarField G(GName(), nut_/(2.0*q_)*2*magSqr(symm(fvc::grad(U_))));
251
252 // Zeta equation
253 tmp<fvScalarMatrix> zetaEqn
254 (
256 + fvm::div(phi_, zeta_)
258 ==
259 (2.0*C1_ - 1)*G*zeta_/q_
260 - fvm::SuSp((2.0*C2_*f2() - dimensionedScalar(1.0))*zeta_/q_, zeta_)
261 + E
262 );
263
264 zetaEqn.ref().relax();
265 solve(zetaEqn);
267
268
269 // q equation
271 (
272 fvm::ddt(q_)
273 + fvm::div(phi_, q_)
275 ==
276 G - fvm::Sp(zeta_/q_, q_)
277 );
278
279 qEqn.ref().relax();
280 solve(qEqn);
281 bound(q_, qMin_);
282
283
284 // Re-calculate k and epsilon
285 k_ = sqr(q_);
287
288 epsilon_ = 2*q_*zeta_;
290
291 correctNut();
292}
293
294
295// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
296
297} // End namespace RASModels
298} // End namespace incompressible
299} // End namespace Foam
300
301// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Bound the given scalar field if it has gone unbounded.
surfaceScalarField & phi
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:78
bool readIfPresent(const word &key, const dictionary &dict)
Update the value of the Switch if it is found in the dictionary.
Definition: Switch.C:335
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Generic dimensioned Type class.
bool readIfPresent(const dictionary &dict)
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:58
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
Gibson and Dafa'Alla's q-zeta two-equation low-Re turbulence model for incompressible flows.
Definition: qZeta.H:79
tmp< volScalarField > f2() const
Definition: qZeta.C:66
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: qZeta.C:240
dimensionedScalar sigmaZeta_
Definition: qZeta.H:90
dimensionedScalar qMin_
Lower limit of q.
Definition: qZeta.H:94
tmp< volScalarField > DqEff() const
Return the effective diffusivity for q.
Definition: qZeta.H:170
tmp< volScalarField > DzetaEff() const
Return the effective diffusivity for epsilon.
Definition: qZeta.H:179
dimensionedScalar zetaMin_
Lower limit of zeta.
Definition: qZeta.H:97
tmp< volScalarField > fMu() const
Definition: qZeta.C:49
virtual bool read()
Re-read model coefficients if they have changed.
Definition: qZeta.C:220
BasicTurbulenceModel::transportModel transportModel
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
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
U
Definition: pEqn.H:72
word timeName
Definition: getTimeIndex.H:3
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:54
tmp< volScalarField > magSqrGradGrad(const GeometricField< Type, fvPatchField, volMesh > &vf)
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:48
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:48
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:48
zeroField SuSp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
zeroField Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
RASModel< turbulenceModel > RASModel
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:35
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
volScalarField & nu
volScalarField & alpha
CEqn solve()