kOmegaSSTIDDES.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) 2015-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 "kOmegaSSTIDDES.H"
30
31// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32
33namespace Foam
34{
35namespace LESModels
36{
37
38// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
39
40template<class BasicTurbulenceModel>
41const IDDESDelta& kOmegaSSTIDDES<BasicTurbulenceModel>::setDelta() const
42{
43 if (!isA<IDDESDelta>(this->delta_()))
44 {
46 << "The delta function must be set to a " << IDDESDelta::typeName
47 << " -based model" << exit(FatalError);
48 }
49
50 return refCast<const IDDESDelta>(this->delta_());
51}
52
53
54template<class BasicTurbulenceModel>
55tmp<volScalarField> kOmegaSSTIDDES<BasicTurbulenceModel>::alpha() const
56{
57 return max(0.25 - this->y_/IDDESDelta_.hmax(), scalar(-5));
58}
59
60
61template<class BasicTurbulenceModel>
62tmp<volScalarField> kOmegaSSTIDDES<BasicTurbulenceModel>::ft
63(
64 const volScalarField& magGradU
65) const
66{
67 return tanh(pow3(sqr(Ct_)*rd(this->nut_, magGradU)));
68}
69
70
71template<class BasicTurbulenceModel>
72tmp<volScalarField> kOmegaSSTIDDES<BasicTurbulenceModel>::fl
73(
74 const volScalarField& magGradU
75) const
76{
77 return tanh(pow(sqr(Cl_)*rd(this->nu(), magGradU), 10));
78}
79
80
81template<class BasicTurbulenceModel>
82tmp<volScalarField> kOmegaSSTIDDES<BasicTurbulenceModel>::rd
83(
84 const volScalarField& nur,
85 const volScalarField& magGradU
86) const
87{
88 tmp<volScalarField> tr
89 (
90 min
91 (
92 nur
93 /(
94 max
95 (
96 magGradU,
97 dimensionedScalar("SMALL", magGradU.dimensions(), SMALL)
98 )
99 *sqr(this->kappa_*this->y_)
100 ),
101 scalar(10)
102 )
103 );
104 tr.ref().boundaryFieldRef() == 0.0;
105
106 return tr;
107}
108
109
110// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
111
112template<class BasicTurbulenceModel>
113tmp<volScalarField> kOmegaSSTIDDES<BasicTurbulenceModel>::fdt
114(
115 const volScalarField& magGradU
116) const
117{
118 return 1 - tanh(pow(Cdt1_*rd(this->nut_, magGradU), Cdt2_));
119}
120
121
122template<class BasicTurbulenceModel>
124(
125 const volScalarField& magGradU,
126 const volScalarField& CDES
127) const
128{
129 const volScalarField& k = this->k_;
130 const volScalarField& omega = this->omega_;
131
132 const volScalarField lRAS(sqrt(k)/(this->betaStar_*omega));
133 const volScalarField lLES(CDES*this->delta());
134
135 const volScalarField alpha(this->alpha());
136 const volScalarField expTerm(exp(sqr(alpha)));
137
138 tmp<volScalarField> fB = min(2*pow(expTerm, -9.0), scalar(1));
140 2*(pos0(alpha)*pow(expTerm, -11.09) + neg(alpha)*pow(expTerm, -9.0));
141 tmp<volScalarField> fe2 = 1 - max(ft(magGradU), fl(magGradU));
142 tmp<volScalarField> fe = max(fe1 - 1, scalar(0))*fe2;
143
144 const volScalarField fdTilda(max(1 - fdt(magGradU), fB));
145
146 // Simplified formulation from Gritskevich et al. paper (2011) where fe = 0
147 // return max
148 // (
149 // fdTilda*lRAS + (1 - fdTilda)*lLES,
150 // dimensionedScalar("SMALL", dimLength, SMALL)
151 // );
152
153 // Original formulation from Shur et al. paper (2008)
154 return max
155 (
156 fdTilda*(1 + fe)*lRAS + (1 - fdTilda)*lLES,
157 dimensionedScalar("SMALL", dimLength, SMALL)
158 );
159}
160
161
162// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
163
164template<class BasicTurbulenceModel>
166(
167 const alphaField& alpha,
168 const rhoField& rho,
169 const volVectorField& U,
170 const surfaceScalarField& alphaRhoPhi,
171 const surfaceScalarField& phi,
172 const transportModel& transport,
173 const word& propertiesName,
174 const word& type
175)
176:
177 kOmegaSSTDES<BasicTurbulenceModel>
178 (
179 alpha,
180 rho,
181 U,
182 alphaRhoPhi,
183 phi,
184 transport,
185 propertiesName,
186 type
187 ),
188
189 Cdt1_
190 (
191 dimensioned<scalar>::getOrAddToDict
192 (
193 "Cdt1",
194 this->coeffDict_,
195 20
196 )
197 ),
198 Cdt2_
199 (
200 dimensioned<scalar>::getOrAddToDict
201 (
202 "Cdt2",
203 this->coeffDict_,
204 3
205 )
206 ),
207 Cl_
208 (
209 dimensioned<scalar>::getOrAddToDict
210 (
211 "Cl",
212 this->coeffDict_,
213 5
214 )
215 ),
216 Ct_
217 (
218 dimensioned<scalar>::getOrAddToDict
219 (
220 "Ct",
221 this->coeffDict_,
222 1.87
223 )
224 ),
225 IDDESDelta_(setDelta())
226{
227 if (type == typeName)
228 {
229 this->printCoeffs(type);
230 }
231}
232
233
234// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
235
236template<class BasicTurbulenceModel>
238{
240 {
241 Cdt1_.readIfPresent(this->coeffDict());
242 Cdt2_.readIfPresent(this->coeffDict());
243 Cl_.readIfPresent(this->coeffDict());
244 Ct_.readIfPresent(this->coeffDict());
245
246 return true;
247 }
248
249 return false;
250}
251
252
253// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
254
255} // End namespace LESModels
256} // End namespace Foam
257
258// ************************************************************************* //
label k
scalar delta
surfaceScalarField & phi
IDDESDelta used by the IDDES (improved low Re Spalart-Allmaras DES model) The min and max delta are c...
k-omega-SST DES turbulence model for incompressible and compressible flows
Definition: kOmegaSSTDES.H:74
k-omega-SST IDDES turbulence model for incompressible and compressible flows
virtual tmp< volScalarField > dTilda(const volScalarField &magGradU, const volScalarField &CDES) const
Length scale.
BasicTurbulenceModel::alphaField alphaField
BasicTurbulenceModel::rhoField rhoField
BasicTurbulenceModel::transportModel transportModel
virtual bool read()
Re-read model coefficients if they have changed.
Generic dimensioned Type class.
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
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
dimensionedScalar pos0(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
dimensionedScalar pow3(const dimensionedScalar &ds)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
dimensionedScalar tanh(const dimensionedScalar &ds)
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)
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)
error FatalError
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
volScalarField & nu
volScalarField & alpha
static const char *const typeName
The type name used in ensight case files.