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-2022 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 generalised
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 
107 
108  //- Runtime type information
109  TypeName("DeardorffDiffStress");
110 
111 
112  // Constructors
113 
114  //- Constructor from components
116  (
117  const alphaField& alpha,
118  const rhoField& rho,
119  const volVectorField& U,
120  const surfaceScalarField& alphaRhoPhi,
121  const surfaceScalarField& phi,
122  const viscosity& viscosity,
123  const word& type = typeName
124  );
125 
126  //- Disallow default bitwise copy construction
127  DeardorffDiffStress(const DeardorffDiffStress&) = delete;
128 
129 
130  //- Destructor
131  virtual ~DeardorffDiffStress()
132  {}
133 
134 
135  // Member Functions
136 
137  //- Read model coefficients if they have changed
138  virtual bool read();
139 
140  //- Return the turbulence kinetic energy dissipation rate
141  virtual tmp<volScalarField> epsilon() const;
142 
143  //- Return the turbulence specific dissipation rate
144  virtual tmp<volScalarField> omega() const;
145 
146  //- Correct sub-grid stress, eddy-Viscosity and related properties
147  virtual void correct();
148 
149 
150  // Member Operators
151 
152  //- Disallow default bitwise assignment
153  void operator=(const DeardorffDiffStress&) = delete;
154 };
155 
156 
157 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
158 
159 } // End namespace LESModels
160 } // End namespace Foam
161 
162 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
163 
164 #ifdef NoRepository
165  #include "DeardorffDiffStress.C"
166 #endif
167 
168 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
169 
170 #endif
171 
172 // ************************************************************************* //
U
Definition: pEqn.H:72
Differential SGS Stress Equation Model for incompressible and compressible flows. ...
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
BasicMomentumTransportModel::rhoField rhoField
BasicMomentumTransportModel::alphaField alphaField
TypeName("DeardorffDiffStress")
Runtime type information.
virtual void correct()
Correct sub-grid stress, eddy-Viscosity and related properties.
DeardorffDiffStress(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity, const word &type=typeName)
Constructor from components.
virtual tmp< volScalarField > omega() const
Return the turbulence specific dissipation rate.
A class for handling words, derived from string.
Definition: word.H:59
phi
Definition: correctPhi.H:3
void operator=(const DeardorffDiffStress &)=delete
Disallow default bitwise assignment.
Abstract base class for all fluid physical properties.
Definition: viscosity.H:49
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.
A class for managing temporary objects.
Definition: PtrList.H:53
Namespace for OpenFOAM.
Reynolds-stress turbulence model base class.