LienCubicKE.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 "LienCubicKE.H"
30#include "wallDist.H"
31#include "bound.H"
33
34// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35
36namespace Foam
37{
38namespace incompressible
39{
40namespace RASModels
41{
42
43// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
44
47
48// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
49
51{
52 const volScalarField yStar(sqrt(k_)*y_/nu());
53
54 return
55 (scalar(1) - exp(-Anu_*yStar))
56 *(scalar(1) + (2*kappa_/(pow(Cmu_, 0.75))/(yStar + SMALL)));
57}
58
59
61{
63
64 return scalar(1) - 0.3*exp(-sqr(Rt));
65}
66
67
68tmp<volScalarField> LienCubicKE::E(const volScalarField& f2) const
69{
70 const volScalarField yStar(sqrt(k_)*y_/nu());
71 const volScalarField le
72 (
73 kappa_*y_/(scalar(1) + (2*kappa_/(pow(Cmu_, 0.75))/(yStar + SMALL)))
74 );
75
76 return
77 (Ceps2_*pow(Cmu_, 0.75))
78 *(f2*sqrt(k_)*epsilon_/le)*exp(-AE_*sqr(yStar));
79}
80
81
83{
85}
86
87
89{
90 volSymmTensorField S(symm(gradU));
91 volTensorField W(skew(gradU));
92
93 volScalarField sBar((k_/epsilon_)*sqrt(2.0)*mag(S));
94 volScalarField wBar((k_/epsilon_)*sqrt(2.0)*mag(W));
95
96 volScalarField Cmu((2.0/3.0)/(Cmu1_ + sBar + Cmu2_*wBar));
97 volScalarField fMu(this->fMu());
98
99 nut_ = Cmu*fMu*sqr(k_)/epsilon_;
101
103 fMu*k_
104 *(
105 // Quadratic terms
106 sqr(k_/epsilon_)/(Cbeta_ + pow3(sBar))
107 *(
109 + Cbeta2_*twoSymm(S&W)
110 + Cbeta3_*dev(symm(W&W))
111 )
112
113 // Cubic terms
114 - pow3(Cmu*k_/epsilon_)
115 *(
116 (Cgamma1_*magSqr(S) - Cgamma2_*magSqr(W))*S
117 + Cgamma4_*twoSymm((innerSqr(S)&W))
118 )
119 );
120}
121
122
123// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
124
126(
128 const geometricOneField& rho,
129 const volVectorField& U,
130 const surfaceScalarField& alphaRhoPhi,
131 const surfaceScalarField& phi,
132 const transportModel& transport,
133 const word& propertiesName,
134 const word& type
135)
136:
137 nonlinearEddyViscosity<incompressible::RASModel>
138 (
139 type,
140 alpha,
141 rho,
142 U,
143 alphaRhoPhi,
144 phi,
145 transport,
146 propertiesName
147 ),
148
149 Ceps1_
150 (
151 dimensioned<scalar>::getOrAddToDict
152 (
153 "Ceps1",
154 coeffDict_,
155 1.44
156 )
157 ),
158 Ceps2_
159 (
160 dimensioned<scalar>::getOrAddToDict
161 (
162 "Ceps2",
163 coeffDict_,
164 1.92
165 )
166 ),
167 sigmak_
168 (
169 dimensioned<scalar>::getOrAddToDict
170 (
171 "sigmak",
172 coeffDict_,
173 1.0
174 )
175 ),
176 sigmaEps_
177 (
178 dimensioned<scalar>::getOrAddToDict
179 (
180 "sigmaEps",
181 coeffDict_,
182 1.3
183 )
184 ),
185 Cmu1_
186 (
187 dimensioned<scalar>::getOrAddToDict
188 (
189 "Cmu1",
190 coeffDict_,
191 1.25
192 )
193 ),
194 Cmu2_
195 (
196 dimensioned<scalar>::getOrAddToDict
197 (
198 "Cmu2",
199 coeffDict_,
200 0.9
201 )
202 ),
203 Cbeta_
204 (
205 dimensioned<scalar>::getOrAddToDict
206 (
207 "Cbeta",
208 coeffDict_,
209 1000.0
210 )
211 ),
212 Cbeta1_
213 (
214 dimensioned<scalar>::getOrAddToDict
215 (
216 "Cbeta1",
217 coeffDict_,
218 3.0
219 )
220 ),
221 Cbeta2_
222 (
223 dimensioned<scalar>::getOrAddToDict
224 (
225 "Cbeta2",
226 coeffDict_,
227 15.0
228 )
229 ),
230 Cbeta3_
231 (
232 dimensioned<scalar>::getOrAddToDict
233 (
234 "Cbeta3",
235 coeffDict_,
236 -19.0
237 )
238 ),
239 Cgamma1_
240 (
241 dimensioned<scalar>::getOrAddToDict
242 (
243 "Cgamma1",
244 coeffDict_,
245 16.0
246 )
247 ),
248 Cgamma2_
249 (
250 dimensioned<scalar>::getOrAddToDict
251 (
252 "Cgamma2",
253 coeffDict_,
254 16.0
255 )
256 ),
257 Cgamma4_
258 (
259 dimensioned<scalar>::getOrAddToDict
260 (
261 "Cgamma4",
262 coeffDict_,
263 -80.0
264 )
265 ),
266 Cmu_
267 (
268 dimensioned<scalar>::getOrAddToDict
269 (
270 "Cmu",
271 coeffDict_,
272 0.09
273 )
274 ),
275 kappa_
276 (
277 dimensioned<scalar>::getOrAddToDict
278 (
279 "kappa",
280 coeffDict_,
281 0.41
282 )
283 ),
284 Anu_
285 (
286 dimensioned<scalar>::getOrAddToDict
287 (
288 "Anu",
289 coeffDict_,
290 0.0198
291 )
292 ),
293 AE_
294 (
295 dimensioned<scalar>::getOrAddToDict
296 (
297 "AE",
298 coeffDict_,
299 0.00375
300 )
301 ),
302
303 k_
304 (
306 (
307 IOobject::groupName("k", alphaRhoPhi.group()),
308 runTime_.timeName(),
309 mesh_,
310 IOobject::MUST_READ,
311 IOobject::AUTO_WRITE
312 ),
313 mesh_
314 ),
315
316 epsilon_
317 (
319 (
320 IOobject::groupName("epsilon", alphaRhoPhi.group()),
321 runTime_.timeName(),
322 mesh_,
323 IOobject::MUST_READ,
324 IOobject::AUTO_WRITE
325 ),
326 mesh_
327 ),
328
329 y_(wallDist::New(mesh_).y())
330{
331 bound(k_, kMin_);
332 bound(epsilon_, epsilonMin_);
333
334 if (type == typeName)
335 {
336 printCoeffs(type);
337 }
338}
339
340
341// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
342
344{
346 {
347 Ceps1_.readIfPresent(coeffDict());
348 Ceps2_.readIfPresent(coeffDict());
349 sigmak_.readIfPresent(coeffDict());
350 sigmaEps_.readIfPresent(coeffDict());
351 Cmu1_.readIfPresent(coeffDict());
352 Cmu2_.readIfPresent(coeffDict());
353 Cbeta_.readIfPresent(coeffDict());
354 Cbeta1_.readIfPresent(coeffDict());
355 Cbeta2_.readIfPresent(coeffDict());
356 Cbeta3_.readIfPresent(coeffDict());
357 Cgamma1_.readIfPresent(coeffDict());
358 Cgamma2_.readIfPresent(coeffDict());
359 Cgamma4_.readIfPresent(coeffDict());
360 Cmu_.readIfPresent(coeffDict());
361 kappa_.readIfPresent(coeffDict());
362 Anu_.readIfPresent(coeffDict());
363 AE_.readIfPresent(coeffDict());
364
365 return true;
366 }
367
368 return false;
369}
370
371
373{
374 if (!turbulence_)
375 {
376 return;
377 }
378
380
381 tmp<volTensorField> tgradU = fvc::grad(U_);
382 const volTensorField& gradU = tgradU();
383
385 (
386 GName(),
387 (nut_*twoSymm(gradU) - nonlinearStress_) && gradU
388 );
389
390
391 // Update epsilon and G at the wall
393
394 const volScalarField f2(this->f2());
395
396 // Dissipation equation
398 (
400 + fvm::div(phi_, epsilon_)
402 ==
405 + E(f2)
406 );
407
408 epsEqn.ref().relax();
409 epsEqn.ref().boundaryManipulate(epsilon_.boundaryFieldRef());
410 solve(epsEqn);
411 bound(epsilon_, epsilonMin_);
412
413
414 // Turbulent kinetic energy equation
416 (
417 fvm::ddt(k_)
418 + fvm::div(phi_, k_)
420 ==
421 G
423 );
424
425 kEqn.ref().relax();
426 solve(kEqn);
427 bound(k_, kMin_);
428
429
430 // Re-calculate viscosity and non-linear stress
432}
433
434
435// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
436
437} // End namespace RASModels
438} // End namespace incompressible
439} // End namespace Foam
440
441// ************************************************************************* //
scalar y
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 updateCoeffs()
Update the boundary condition coefficients.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
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)
volScalarField nut_
Definition: eddyViscosity.H:66
BasicTurbulenceModel::transportModel transportModel
Definition: eddyViscosity.H:78
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
Lien cubic non-linear low-Reynolds k-epsilon turbulence models for incompressible flows.
Definition: LienCubicKE.H:83
tmp< volScalarField > f2() const
Definition: LienCubicKE.C:60
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: LienCubicKE.C:372
tmp< volScalarField > DepsilonEff() const
Return the effective diffusivity for epsilon.
Definition: LienCubicKE.H:173
const volScalarField & y_
Wall distance.
Definition: LienCubicKE.H:120
tmp< volScalarField > E(const volScalarField &f2) const
Definition: LienCubicKE.C:68
virtual void correctNonlinearStress(const volTensorField &gradU)
Definition: LienCubicKE.C:88
tmp< volScalarField > fMu() const
Definition: LienCubicKE.C:50
virtual bool read()
Re-read model coefficients if they have changed.
Definition: LienCubicKE.C:343
tmp< volScalarField > DkEff() const
Return the effective diffusivity for k.
Definition: LienCubicKE.H:164
Eddy viscosity turbulence model with non-linear correction base class.
A class for managing temporary objects.
Definition: tmp.H:65
T & ref() const
Definition: tmpI.H:227
Interface to run-time selectable methods to calculate the distance-to-wall and normal-to-wall fields.
Definition: wallDist.H:78
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< 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 Sp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
autoPtr< BasicCompressibleTurbulenceModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const typename BasicCompressibleTurbulenceModel::transportModel &transport, const word &propertiesName)
RASModel< turbulenceModel > RASModel
Namespace for OpenFOAM.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
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
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
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 pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dimensionedSymmTensor innerSqr(const dimensionedSymmTensor &dt)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
dimensionedTensor skew(const dimensionedTensor &dt)
volScalarField & nu
volScalarField & alpha
CEqn solve()