kOmegaSSTDES.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) 2015 OpenFOAM Foundation
9 Copyright (C) 2016-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 "kOmegaSSTDES.H"
30
31// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32
33namespace Foam
34{
35namespace LESModels
36{
37
38// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
39
40template<class BasicTurbulenceModel>
42{
43 // Correct the turbulence viscosity
45
46 // Correct the turbulence thermal diffusivity
47 BasicTurbulenceModel::correctNut();
48}
49
50
51template<class BasicTurbulenceModel>
53{
54 correctNut(2*magSqr(symm(fvc::grad(this->U_))));
55}
56
57
58template<class BasicTurbulenceModel>
60(
61 const volScalarField& magGradU,
62 const volScalarField& CDES
63) const
64{
65 const volScalarField& k = this->k_;
66 const volScalarField& omega = this->omega_;
67
68 return min(CDES*this->delta(), sqrt(k)/(this->betaStar_*omega));
69}
70
71
72template<class BasicTurbulenceModel>
74(
75 const volScalarField& F1,
76 const volTensorField& gradU
77) const
78{
79 volScalarField CDES(this->CDES(F1));
80 return sqrt(this->k_())/dTilda(mag(gradU), CDES)()();
81}
82
83
84template<class BasicTurbulenceModel>
86(
87 const volScalarField::Internal& GbyNu0,
90) const
91{
92 return GbyNu0; // Unlimited
93}
94
95
96// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
97
98template<class BasicTurbulenceModel>
100(
101 const alphaField& alpha,
102 const rhoField& rho,
103 const volVectorField& U,
104 const surfaceScalarField& alphaRhoPhi,
105 const surfaceScalarField& phi,
106 const transportModel& transport,
107 const word& propertiesName,
108 const word& type
109)
110:
111 kOmegaSSTBase<DESModel<BasicTurbulenceModel>>
112 (
113 type,
114 alpha,
115 rho,
116 U,
117 alphaRhoPhi,
118 phi,
119 transport,
120 propertiesName
121 ),
122
123 kappa_
124 (
125 dimensioned<scalar>::getOrAddToDict
126 (
127 "kappa",
128 this->coeffDict_,
129 0.41
130 )
131 ),
132 CDESkom_
133 (
134 dimensioned<scalar>::getOrAddToDict
135 (
136 "CDESkom",
137 this->coeffDict_,
138 0.82
139 )
140 ),
141 CDESkeps_
142 (
143 dimensioned<scalar>::getOrAddToDict
144 (
145 "CDESkeps",
146 this->coeffDict_,
147 0.60
148 )
149 )
150{
151 if (type == typeName)
152 {
153 this->printCoeffs(type);
154 }
155}
156
157
158// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
159
160template<class BasicTurbulenceModel>
162{
164 {
165 kappa_.readIfPresent(this->coeffDict());
166 CDESkom_.readIfPresent(this->coeffDict());
167 CDESkeps_.readIfPresent(this->coeffDict());
168
169 return true;
170 }
171
172 return false;
173}
174
175
176template<class BasicTurbulenceModel>
178{
179 const volScalarField& k = this->k_;
180 const volScalarField& omega = this->omega_;
181 const volVectorField& U = this->U_;
182
183 const volScalarField CDkOmega
184 (
185 (2*this->alphaOmega2_)*(fvc::grad(k) & fvc::grad(omega))/omega
186 );
187
188 const volScalarField F1(this->F1(CDkOmega));
189
190 tmp<volScalarField> tLESRegion
191 (
193 (
195 (
196 "DES::LESRegion",
197 this->mesh_.time().timeName(),
198 this->mesh_,
201 ),
202 neg
203 (
204 dTilda
205 (
206 mag(fvc::grad(U)),
207 F1*CDESkom_ + (1 - F1)*CDESkeps_
208 )
209 - sqrt(k)/(this->betaStar_*omega)
210 )
211 )
212 );
213
214 return tLESRegion;
215}
216
217
218// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
219
220} // End namespace LESModels
221} // End namespace Foam
222
223// ************************************************************************* //
label k
scalar delta
#define F2(B, C, D)
Definition: SHA1.C:151
#define F1(B, C, D)
Definition: SHA1.C:150
surfaceScalarField & phi
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
Templated abstract base class for DES models.
Definition: DESModel.H:58
k-omega-SST DES turbulence model for incompressible and compressible flows
Definition: kOmegaSSTDES.H:74
virtual tmp< volScalarField > dTilda(const volScalarField &magGradU, const volScalarField &CDES) const
Length scale.
Definition: kOmegaSSTDES.C:60
BasicTurbulenceModel::alphaField alphaField
Definition: kOmegaSSTDES.H:131
BasicTurbulenceModel::rhoField rhoField
Definition: kOmegaSSTDES.H:132
virtual tmp< volScalarField > LESRegion() const
Return the LES field indicator.
Definition: kOmegaSSTDES.C:177
virtual tmp< volScalarField::Internal > GbyNu(const volScalarField::Internal &GbyNu0, const volScalarField::Internal &F2, const volScalarField::Internal &S2) const
Return G/nu.
Definition: kOmegaSSTDES.C:86
BasicTurbulenceModel::transportModel transportModel
Definition: kOmegaSSTDES.H:133
virtual tmp< volScalarField::Internal > epsilonByk(const volScalarField &F1, const volTensorField &gradU) const
Return epsilon/k.
Definition: kOmegaSSTDES.C:74
virtual bool read()
Re-read model coefficients if they have changed.
Definition: kOmegaSSTDES.C:161
Generic dimensioned Type class.
Base class implementation of the k-omega-SST turbulence model for incompressible and compressible flo...
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
U
Definition: pEqn.H:72
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:54
Namespace for OpenFOAM.
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
bool read(const char *buf, int32_t &val)
Same as readInt32.
Definition: int32.H:108
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 > mag(const dimensioned< Type > &dt)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
dimensionedScalar neg(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
volScalarField & alpha