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