WALE.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2015-2018 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "WALE.H"
27 #include "fvOptions.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace LESModels
34 {
35 
36 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
37 
38 template<class BasicTurbulenceModel>
39 tmp<volSymmTensorField> WALE<BasicTurbulenceModel>::Sd
40 (
41  const volTensorField& gradU
42 ) const
43 {
44  return dev(symm(gradU & gradU));
45 }
46 
47 
48 template<class BasicTurbulenceModel>
50 (
51  const volTensorField& gradU
52 ) const
53 {
54  volScalarField magSqrSd(magSqr(Sd(gradU)));
55 
56  return tmp<volScalarField>
57  (
58  new volScalarField
59  (
60  IOobject
61  (
62  IOobject::groupName("k", this->alphaRhoPhi_.group()),
63  this->runTime_.timeName(),
64  this->mesh_
65  ),
66  sqr(sqr(Cw_)*this->delta()/Ck_)*
67  (
68  pow3(magSqrSd)
69  /(
70  sqr
71  (
72  pow(magSqr(symm(gradU)), 5.0/2.0)
73  + pow(magSqrSd, 5.0/4.0)
74  )
76  (
77  "small",
78  dimensionSet(0, 0, -10, 0, 0),
79  small
80  )
81  )
82  )
83  )
84  );
85 }
86 
87 
88 template<class BasicTurbulenceModel>
90 {
91  this->nut_ = Ck_*this->delta()*sqrt(this->k(fvc::grad(this->U_)));
92  this->nut_.correctBoundaryConditions();
93  fv::options::New(this->mesh_).correct(this->nut_);
94 
95  BasicTurbulenceModel::correctNut();
96 }
97 
98 
99 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
100 
101 template<class BasicTurbulenceModel>
103 (
104  const alphaField& alpha,
105  const rhoField& rho,
106  const volVectorField& U,
107  const surfaceScalarField& alphaRhoPhi,
108  const surfaceScalarField& phi,
109  const transportModel& transport,
110  const word& propertiesName,
111  const word& type
112 )
113 :
115  (
116  type,
117  alpha,
118  rho,
119  U,
120  alphaRhoPhi,
121  phi,
122  transport,
123  propertiesName
124  ),
125 
126  Ck_
127  (
129  (
130  "Ck",
131  this->coeffDict_,
132  0.094
133  )
134  ),
135 
136  Cw_
137  (
139  (
140  "Cw",
141  this->coeffDict_,
142  0.325
143  )
144  )
145 {
146  if (type == typeName)
147  {
148  this->printCoeffs(type);
149  }
150 }
151 
152 
153 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
154 
155 template<class BasicTurbulenceModel>
157 {
159  {
160  Ck_.readIfPresent(this->coeffDict());
161  Cw_.readIfPresent(this->coeffDict());
162 
163  return true;
164  }
165  else
166  {
167  return false;
168  }
169 }
170 
171 
172 template<class BasicTurbulenceModel>
174 {
175  volScalarField k(this->k(fvc::grad(this->U_)));
176 
177  return tmp<volScalarField>
178  (
179  new volScalarField
180  (
181  IOobject
182  (
183  IOobject::groupName("epsilon", this->alphaRhoPhi_.group()),
184  this->runTime_.timeName(),
185  this->mesh_,
188  ),
189  this->Ce_*k*sqrt(k)/this->delta()
190  )
191  );
192 }
193 
194 
195 template<class BasicTurbulenceModel>
197 {
199  correctNut();
200 }
201 
202 
203 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
204 
205 } // End namespace LESModels
206 } // End namespace Foam
207 
208 // ************************************************************************* //
virtual bool read()
Read model coefficients if they have changed.
Definition: WALE.C:156
Eddy viscosity LES SGS model base class.
scalar delta
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
surfaceScalarField & phi
BasicTurbulenceModel::transportModel transportModel
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
virtual tmp< volScalarField > epsilon() const
Return sub-grid disipation rate.
Definition: WALE.C:173
dimensionedScalar sqrt(const dimensionedScalar &ds)
label k
Boltzmann constant.
Generic dimensioned Type class.
BasicTurbulenceModel::alphaField alphaField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual void correctNut()
Update the SGS eddy-viscosity.
Definition: WALE.C:89
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
Dimension set for the base types.
Definition: dimensionSet.H:120
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
tmp< volSymmTensorField > Sd(const volTensorField &gradU) const
Return the deviatoric symmetric part of the square of the given.
Definition: WALE.C:40
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:481
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar pow3(const dimensionedScalar &ds)
U
Definition: pEqn.H:72
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
The Wall-adapting local eddy-viscosity (WALE) SGS model.
Definition: WALE.H:73
BasicTurbulenceModel::rhoField rhoField
virtual void correct()
Correct Eddy-Viscosity and related properties.
Definition: WALE.C:196
A class for managing temporary objects.
Definition: PtrList.H:53
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
Definition: fvOptions.C:101
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.
Definition: WALE.H:144
Namespace for OpenFOAM.