All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Smagorinsky.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) 2011-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::Smagorinsky
26 
27 Description
28  The Smagorinsky SGS model.
29 
30  Reference:
31  \verbatim
32  Smagorinsky, J. (1963).
33  General circulation experiments with the primitive equations: I.
34  The basic experiment*.
35  Monthly weather review, 91(3), 99-164.
36  \endverbatim
37 
38  The form of the Smagorinsky model implemented is obtained from the
39  k-equation model assuming local equilibrium which provides estimates of both
40  k and epsilon separate from the sub-grid scale viscosity:
41 
42  \verbatim
43  B = 2/3*k*I - 2*nuSgs*dev(D)
44 
45  where
46 
47  D = symm(grad(U));
48  k from D:B + Ce*k^3/2/delta = 0
49  nuSgs = Ck*sqrt(k)*delta
50  \endverbatim
51 
52  The default model coefficients are
53  \verbatim
54  SmagorinskyCoeffs
55  {
56  Ck 0.094;
57  Ce 1.048;
58  }
59  \endverbatim
60 
61 See also
62  Foam::LESModels::kEqn
63 
64 SourceFiles
65  Smagorinsky.C
66 
67 \*---------------------------------------------------------------------------*/
68 
69 #ifndef Smagorinsky_H
70 #define Smagorinsky_H
71 
72 #include "LESModel.H"
73 #include "LESeddyViscosity.H"
74 
75 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
76 
77 namespace Foam
78 {
79 namespace LESModels
80 {
81 
82 /*---------------------------------------------------------------------------*\
83  Class Smagorinsky Declaration
84 \*---------------------------------------------------------------------------*/
85 
86 template<class BasicMomentumTransportModel>
87 class Smagorinsky
88 :
89  public LESeddyViscosity<BasicMomentumTransportModel>
90 {
91 protected:
92 
93  // Protected data
94 
96 
97 
98  // Protected Member Functions
99 
100  //- Return SGS kinetic energy
101  // calculated from the given velocity gradient
102  tmp<volScalarField> k(const tmp<volTensorField>& gradU) const;
103 
104  //- Update the SGS eddy viscosity
105  virtual void correctNut();
106 
107 
108 public:
110  typedef typename BasicMomentumTransportModel::alphaField alphaField;
111  typedef typename BasicMomentumTransportModel::rhoField rhoField;
112  typedef typename BasicMomentumTransportModel::transportModel transportModel;
113 
114 
115  //- Runtime type information
116  TypeName("Smagorinsky");
117 
118 
119  // Constructors
120 
121  //- Construct from components
123  (
124  const alphaField& alpha,
125  const rhoField& rho,
126  const volVectorField& U,
127  const surfaceScalarField& alphaRhoPhi,
128  const surfaceScalarField& phi,
129  const transportModel& transport,
130  const word& type = typeName
131  );
132 
133  //- Disallow default bitwise copy construction
134  Smagorinsky(const Smagorinsky&) = delete;
135 
136 
137  //- Destructor
138  virtual ~Smagorinsky()
139  {}
140 
141 
142  // Member Functions
143 
144  //- Read model coefficients if they have changed
145  virtual bool read();
146 
147  //- Return SGS kinetic energy
148  virtual tmp<volScalarField> k() const
149  {
150  return k(fvc::grad(this->U_));
151  }
152 
153  //- Return sub-grid disipation rate
154  virtual tmp<volScalarField> epsilon() const;
155 
156  //- Correct Eddy-Viscosity and related properties
157  virtual void correct();
158 
159 
160  // Member Operators
161 
162  //- Disallow default bitwise assignment
163  void operator=(const Smagorinsky&) = delete;
164 };
165 
166 
167 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
168 
169 } // End namespace LESModels
170 } // End namespace Foam
171 
172 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
173 
174 #ifdef NoRepository
175  #include "Smagorinsky.C"
176 #endif
177 
178 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
179 
180 #endif
181 
182 // ************************************************************************* //
Eddy viscosity LES SGS model base class.
TypeName("Smagorinsky")
Runtime type information.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
virtual tmp< volScalarField > epsilon() const
Return sub-grid disipation rate.
Definition: Smagorinsky.C:131
virtual bool read()
Read model coefficients if they have changed.
Definition: Smagorinsky.C:115
BasicMomentumTransportModel::transportModel transportModel
Definition: Smagorinsky.H:111
void operator=(const Smagorinsky &)=delete
Disallow default bitwise assignment.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
virtual void correctNut()
Update the SGS eddy viscosity.
Definition: Smagorinsky.C:60
The Smagorinsky SGS model.
Definition: Smagorinsky.H:86
Smagorinsky(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: Smagorinsky.C:74
A class for handling words, derived from string.
Definition: word.H:59
phi
Definition: correctPhi.H:3
BasicMomentumTransportModel::rhoField rhoField
Definition: Smagorinsky.H:110
virtual void correct()
Correct Eddy-Viscosity and related properties.
Definition: Smagorinsky.C:144
BasicMomentumTransportModel::alphaField alphaField
Definition: Smagorinsky.H:109
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.
Definition: Smagorinsky.H:147
virtual ~Smagorinsky()
Destructor.
Definition: Smagorinsky.H:137
U
Definition: pEqn.H:72
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dimensionedScalar Ck_
Definition: Smagorinsky.H:94
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.