WALE.H
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-2020 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 Class
25  Foam::LESModels::WALE
26 
27 Description
28  The Wall-adapting local eddy-viscosity (WALE) SGS model.
29 
30  Reference:
31  \verbatim
32  Nicoud, F., & Ducros, F. (1999).
33  Subgrid-scale stress modelling based on the square of the velocity
34  gradient tensor.
35  Flow, Turbulence and Combustion, 62(3), 183-200.
36  \endverbatim
37 
38  The default model coefficients are
39  \verbatim
40  WALECoeffs
41  {
42  Ck 0.094;
43  Ce 1.048;e
44  Cw 0.325;
45  }
46  \endverbatim
47 
48 See also
49  Foam::LESModels::Smagorinsky
50 
51 SourceFiles
52  WALE.C
53 
54 \*---------------------------------------------------------------------------*/
55 
56 #ifndef WALE_H
57 #define WALE_H
58 
59 #include "LESModel.H"
60 #include "LESeddyViscosity.H"
61 
62 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
63 
64 namespace Foam
65 {
66 namespace LESModels
67 {
68 
69 /*---------------------------------------------------------------------------*\
70  Class WALE Declaration
71 \*---------------------------------------------------------------------------*/
72 
73 template<class BasicMomentumTransportModel>
74 class WALE
75 :
76  public LESeddyViscosity<BasicMomentumTransportModel>
77 {
78 protected:
79 
80  // Protected data
81 
84 
85 
86  // Protected Member Functions
87 
88  //- Return the deviatoric symmetric part of the square of the given
89  // velocity gradient field
90  tmp<volSymmTensorField> Sd(const volTensorField& gradU) const;
91 
92  //- Return SGS kinetic energy
93  // calculated from the given velocity gradient
94  tmp<volScalarField> k(const volTensorField& gradU) const;
95 
96  //- Update the SGS eddy-viscosity
97  virtual void correctNut();
98 
99 
100 public:
102  typedef typename BasicMomentumTransportModel::alphaField alphaField;
103  typedef typename BasicMomentumTransportModel::rhoField rhoField;
104  typedef typename BasicMomentumTransportModel::transportModel transportModel;
105 
106 
107  //- Runtime type information
108  TypeName("WALE");
109 
110 
111  // Constructors
112 
113  //- Construct from components
114  WALE
115  (
116  const alphaField& alpha,
117  const rhoField& rho,
118  const volVectorField& U,
119  const surfaceScalarField& alphaRhoPhi,
120  const surfaceScalarField& phi,
121  const transportModel& transport,
122  const word& type = typeName
123  );
124 
125  //- Disallow default bitwise copy construction
126  WALE(const WALE&) = delete;
127 
128 
129  //- Destructor
130  virtual ~WALE()
131  {}
132 
133 
134  // Member Functions
135 
136  //- Read model coefficients if they have changed
137  virtual bool read();
138 
139  //- Return SGS kinetic energy
140  virtual tmp<volScalarField> k() const
141  {
142  return k(fvc::grad(this->U_));
143  }
144 
145  //- Return sub-grid disipation rate
146  virtual tmp<volScalarField> epsilon() const;
147 
148  //- Correct Eddy-Viscosity and related properties
149  virtual void correct();
150 
151 
152  // Member Operators
153 
154  //- Disallow default bitwise assignment
155  void operator=(const WALE&) = delete;
156 };
157 
158 
159 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
160 
161 } // End namespace LESModels
162 } // End namespace Foam
163 
164 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
165 
166 #ifdef NoRepository
167  #include "WALE.C"
168 #endif
169 
170 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
171 
172 #endif
173 
174 // ************************************************************************* //
virtual bool read()
Read model coefficients if they have changed.
Definition: WALE.C:145
Eddy viscosity LES SGS model base class.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
virtual ~WALE()
Destructor.
Definition: WALE.H:129
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
virtual tmp< volScalarField > epsilon() const
Return sub-grid disipation rate.
Definition: WALE.C:162
WALE(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &type=typeName)
Construct from components.
Definition: WALE.C:94
virtual void correctNut()
Update the SGS eddy-viscosity.
Definition: WALE.C:82
tmp< volSymmTensorField > Sd(const volTensorField &gradU) const
Return the deviatoric symmetric part of the square of the given.
Definition: WALE.C:41
BasicMomentumTransportModel::rhoField rhoField
Definition: WALE.H:102
A class for handling words, derived from string.
Definition: word.H:59
TypeName("WALE")
Runtime type information.
void operator=(const WALE &)=delete
Disallow default bitwise assignment.
phi
Definition: correctPhi.H:3
dimensionedScalar Cw_
Definition: WALE.H:82
dimensionedScalar Ck_
Definition: WALE.H:81
U
Definition: pEqn.H:72
The Wall-adapting local eddy-viscosity (WALE) SGS model.
Definition: WALE.H:73
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
virtual void correct()
Correct Eddy-Viscosity and related properties.
Definition: WALE.C:175
A class for managing temporary objects.
Definition: PtrList.H:53
BasicMomentumTransportModel::transportModel transportModel
Definition: WALE.H:103
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.
Definition: WALE.H:139
Namespace for OpenFOAM.
BasicMomentumTransportModel::alphaField alphaField
Definition: WALE.H:101