DeardorffDiffStress.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::DeardorffDiffStress
26 
27 Description
28  Differential SGS Stress Equation Model for incompressible and
29  compressible flows
30 
31  Reference:
32  \verbatim
33  Deardorff, J. W. (1973).
34  The use of subgrid transport equations in a three-dimensional model
35  of atmospheric turbulence.
36  Journal of Fluids Engineering, 95(3), 429-438.
37  \endverbatim
38 
39  This SGS model uses a full balance equation for the SGS stress tensor to
40  simulate the behaviour of B.
41 
42  This implementation is as described in the above paper except that the
43  triple correlation model of Donaldson is replaced with the generalized
44  gradient diffusion model of Daly and Harlow:
45  \verbatim
46  Daly, B. J., & Harlow, F. H. (1970).
47  Transport equations in turbulence.
48  Physics of Fluids (1958-1988), 13(11), 2634-2649.
49  \endverbatim
50  with the default value for the coefficient Cs of 0.25 from
51  \verbatim
52  Launder, B. E., Reece, G. J., & Rodi, W. (1975).
53  Progress in the development of a Reynolds-stress turbulence closure.
54  Journal of fluid mechanics, 68(03), 537-566.
55  \endverbatim
56 
57 SourceFiles
58  DeardorffDiffStress.C
59 
60 \*---------------------------------------------------------------------------*/
61 
62 #ifndef DeardorffDiffStress_H
63 #define DeardorffDiffStress_H
64 
65 #include "LESModel.H"
66 #include "ReynoldsStress.H"
67 
68 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
69 
70 namespace Foam
71 {
72 namespace LESModels
73 {
74 
75 /*---------------------------------------------------------------------------*\
76  Class DeardorffDiffStress Declaration
77 \*---------------------------------------------------------------------------*/
78 
79 template<class BasicMomentumTransportModel>
81 :
82  public ReynoldsStress<LESModel<BasicMomentumTransportModel>>
83 {
84 protected:
85 
86  // Protected data
87 
88  // Model constants
89 
94 
95 
96  // Protected Member Functions
97 
98  //- Update the eddy-viscosity
99  virtual void correctNut();
100 
101 
102 public:
104  typedef typename BasicMomentumTransportModel::alphaField alphaField;
105  typedef typename BasicMomentumTransportModel::rhoField rhoField;
106  typedef typename BasicMomentumTransportModel::transportModel transportModel;
107 
108 
109  //- Runtime type information
110  TypeName("DeardorffDiffStress");
111 
112 
113  // Constructors
114 
115  //- Constructor from components
117  (
118  const alphaField& alpha,
119  const rhoField& rho,
120  const volVectorField& U,
121  const surfaceScalarField& alphaRhoPhi,
122  const surfaceScalarField& phi,
123  const transportModel& transport,
124  const word& type = typeName
125  );
126 
127  //- Disallow default bitwise copy construction
128  DeardorffDiffStress(const DeardorffDiffStress&) = delete;
129 
130 
131  //- Destructor
132  virtual ~DeardorffDiffStress()
133  {}
134 
135 
136  // Member Functions
137 
138  //- Read model coefficients if they have changed
139  virtual bool read();
140 
141  //- Return the turbulence kinetic energy dissipation rate
142  virtual tmp<volScalarField> epsilon() const;
143 
144  //- Correct sub-grid stress, eddy-Viscosity and related properties
145  virtual void correct();
146 
147 
148  // Member Operators
149 
150  //- Disallow default bitwise assignment
151  void operator=(const DeardorffDiffStress&) = delete;
152 };
153 
154 
155 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
156 
157 } // End namespace LESModels
158 } // End namespace Foam
159 
160 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
161 
162 #ifdef NoRepository
163  #include "DeardorffDiffStress.C"
164 #endif
165 
166 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
167 
168 #endif
169 
170 // ************************************************************************* //
Differential SGS Stress Equation Model for incompressible and compressible flows. ...
BasicMomentumTransportModel::rhoField rhoField
BasicMomentumTransportModel::alphaField alphaField
phi
Definition: pEqn.H:104
DeardorffDiffStress(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &type=typeName)
Constructor from components.
TypeName("DeardorffDiffStress")
Runtime type information.
virtual void correct()
Correct sub-grid stress, eddy-Viscosity and related properties.
A class for handling words, derived from string.
Definition: word.H:59
void operator=(const DeardorffDiffStress &)=delete
Disallow default bitwise assignment.
U
Definition: pEqn.H:72
BasicMomentumTransportModel::transportModel transportModel
virtual void correctNut()
Update the eddy-viscosity.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
virtual bool read()
Read model coefficients if they have changed.
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.
Reynolds-stress turbulence model base class.