kOmegaSSTDDES.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 "kOmegaSSTDDES.H"
30
31// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32
33namespace Foam
34{
35namespace LESModels
36{
37
38// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
39
40template<class BasicTurbulenceModel>
41tmp<volScalarField> kOmegaSSTDDES<BasicTurbulenceModel>::rd
42(
43 const volScalarField& magGradU
44) const
45{
46 tmp<volScalarField> tr
47 (
48 min
49 (
50 this->nuEff()
51 /(
52 max
53 (
54 magGradU,
55 dimensionedScalar("SMALL", magGradU.dimensions(), SMALL)
56 )
57 *sqr(this->kappa_*this->y_)
58 ),
59 scalar(10)
60 )
61 );
62 tr.ref().boundaryFieldRef() == 0.0;
63
64 return tr;
65}
66
67
68template<class BasicTurbulenceModel>
69tmp<volScalarField> kOmegaSSTDDES<BasicTurbulenceModel>::fd
70(
71 const volScalarField& magGradU
72) const
73{
74 return 1 - tanh(pow(Cd1_*rd(magGradU), Cd2_));
75}
76
77
78// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
79
80template<class BasicTurbulenceModel>
82(
83 const volScalarField& magGradU,
84 const volScalarField& CDES
85) const
86{
87 const volScalarField& k = this->k_;
88 const volScalarField& omega = this->omega_;
89
90 const volScalarField lRAS(sqrt(k)/(this->betaStar_*omega));
91 const volScalarField lLES(CDES*this->delta());
92
93 return max
94 (
95 lRAS
96 - fd(magGradU)
97 *max
98 (
99 lRAS - lLES,
101 ),
102 dimensionedScalar("small", dimLength, SMALL)
103 );
104}
105
106
107// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
108
109template<class BasicTurbulenceModel>
111(
112 const alphaField& alpha,
113 const rhoField& rho,
114 const volVectorField& U,
115 const surfaceScalarField& alphaRhoPhi,
116 const surfaceScalarField& phi,
117 const transportModel& transport,
118 const word& propertiesName,
119 const word& type
120)
121:
122 kOmegaSSTDES<BasicTurbulenceModel>
123 (
124 alpha,
125 rho,
126 U,
127 alphaRhoPhi,
128 phi,
129 transport,
130 propertiesName,
131 type
132 ),
133
134 Cd1_
135 (
136 dimensioned<scalar>::getOrAddToDict
137 (
138 "Cd1",
139 this->coeffDict_,
140 20
141 )
142 ),
143 Cd2_
144 (
145 dimensioned<scalar>::getOrAddToDict
146 (
147 "Cd2",
148 this->coeffDict_,
149 3
150 )
151 )
152{
153 if (type == typeName)
154 {
155 this->printCoeffs(type);
156 }
157}
158
159
160// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
161
162template<class BasicTurbulenceModel>
164{
166 {
167 Cd1_.readIfPresent(this->coeffDict());
168 Cd2_.readIfPresent(this->coeffDict());
169
170 return true;
171 }
172
173 return false;
174}
175
176
177// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
178
179} // End namespace LESModels
180} // End namespace Foam
181
182// ************************************************************************* //
label k
scalar delta
surfaceScalarField & phi
k-omega-SST DDES turbulence model for incompressible and compressible flows
Definition: kOmegaSSTDDES.H:68
virtual tmp< volScalarField > dTilda(const volScalarField &magGradU, const volScalarField &CDES) const
Length scale.
Definition: kOmegaSSTDDES.C:82
BasicTurbulenceModel::alphaField alphaField
BasicTurbulenceModel::rhoField rhoField
BasicTurbulenceModel::transportModel transportModel
virtual bool read()
Re-read model coefficients if they have changed.
k-omega-SST DES turbulence model for incompressible and compressible flows
Definition: kOmegaSSTDES.H:74
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
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 tr(const dimensionedSphericalTensor &dt)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:52
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
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
volScalarField & alpha