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-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 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 BasicTurbulenceModel>
87 class Smagorinsky
88 :
89  public LESeddyViscosity<BasicTurbulenceModel>
90 {
91  // Private Member Functions
92 
93  // Disallow default bitwise copy construct and assignment
94  Smagorinsky(const Smagorinsky&);
95  void operator=(const Smagorinsky&);
96 
97 
98 protected:
99 
100  // Protected data
103 
104 
105  // Protected Member Functions
106 
107  //- Return SGS kinetic energy
108  // calculated from the given velocity gradient
109  tmp<volScalarField> k(const tmp<volTensorField>& gradU) const;
110 
111  //- Update the SGS eddy viscosity
112  virtual void correctNut();
113 
114 
115 public:
117  typedef typename BasicTurbulenceModel::alphaField alphaField;
118  typedef typename BasicTurbulenceModel::rhoField rhoField;
119  typedef typename BasicTurbulenceModel::transportModel transportModel;
120 
121 
122  //- Runtime type information
123  TypeName("Smagorinsky");
124 
125 
126  // Constructors
127 
128  //- Construct from components
130  (
131  const alphaField& alpha,
132  const rhoField& rho,
133  const volVectorField& U,
134  const surfaceScalarField& alphaRhoPhi,
135  const surfaceScalarField& phi,
136  const transportModel& transport,
137  const word& propertiesName = turbulenceModel::propertiesName,
138  const word& type = typeName
139  );
140 
141 
142  //- Destructor
143  virtual ~Smagorinsky()
144  {}
145 
146 
147  // Member Functions
148 
149  //- Read model coefficients if they have changed
150  virtual bool read();
151 
152  //- Return SGS kinetic energy
153  virtual tmp<volScalarField> k() const
154  {
155  return k(fvc::grad(this->U_));
156  }
157 
158  //- Return sub-grid disipation rate
159  virtual tmp<volScalarField> epsilon() const;
160 
161  //- Correct Eddy-Viscosity and related properties
162  virtual void correct();
163 };
164 
165 
166 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
167 
168 } // End namespace LESModels
169 } // End namespace Foam
170 
171 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
172 
173 #ifdef NoRepository
174  #include "Smagorinsky.C"
175 #endif
176 
177 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
178 
179 #endif
180 
181 // ************************************************************************* //
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:142
BasicTurbulenceModel::alphaField alphaField
Definition: Smagorinsky.H:116
surfaceScalarField & phi
virtual bool read()
Read model coefficients if they have changed.
Definition: Smagorinsky.C:126
BasicTurbulenceModel::transportModel transportModel
Definition: Smagorinsky.H:118
virtual void correctNut()
Update the SGS eddy viscosity.
Definition: Smagorinsky.C:67
The Smagorinsky SGS model.
Definition: Smagorinsky.H:86
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:117
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:481
virtual tmp< volScalarField > k() const
Return SGS kinetic energy.
Definition: Smagorinsky.H:152
virtual ~Smagorinsky()
Destructor.
Definition: Smagorinsky.H:142
U
Definition: pEqn.H:72
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Namespace for OpenFOAM.