LaunderSharmaKE.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::RASModels::LaunderSharmaKE
26 
27 Group
28  grpRASTurbulence
29 
30 Description
31  Launder and Sharma low-Reynolds k-epsilon turbulence model for
32  incompressible and compressible and combusting flows including
33  rapid distortion theory (RDT) based compression term.
34 
35  References:
36  \verbatim
37  Launder, B. E., & Sharma, B. I. (1974).
38  Application of the energy-dissipation model of turbulence to the
39  calculation of flow near a spinning disc.
40  Letters in heat and mass transfer, 1(2), 131-137.
41 
42  For the RDT-based compression term:
43  El Tahry, S. H. (1983).
44  k-epsilon equation for compressible reciprocating engine flows.
45  Journal of Energy, 7(4), 345-353.
46  \endverbatim
47 
48  The default model coefficients are
49  \verbatim
50  LaunderSharmaKECoeffs
51  {
52  Cmu 0.09;
53  C1 1.44;
54  C2 1.92;
55  C3 -0.33;
56  alphah 1.0; // only for compressible
57  alphahk 1.0; // only for compressible
58  alphaEps 0.76923;
59  }
60  \endverbatim
61 
62 SourceFiles
63  LaunderSharmaKE.C
64 
65 \*---------------------------------------------------------------------------*/
66 
67 #ifndef LaunderSharmaKE_H
68 #define LaunderSharmaKE_H
69 
70 #include "RASModel.H"
71 #include "eddyViscosity.H"
72 
73 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
74 
75 namespace Foam
76 {
77 namespace RASModels
78 {
79 
80 /*---------------------------------------------------------------------------*\
81  Class LaunderSharmaKE Declaration
82 \*---------------------------------------------------------------------------*/
83 
84 template<class BasicTurbulenceModel>
85 class LaunderSharmaKE
86 :
87  public eddyViscosity<RASModel<BasicTurbulenceModel>>
88 {
89  // Private Member Functions
90 
91  // Disallow default bitwise copy construct and assignment
93  void operator=(const LaunderSharmaKE&);
94 
95 
96 protected:
97 
98  // Protected data
99 
100  // Model coefficients
108 
109 
110  // Fields
114 
115 
116  // Private Member Functions
117 
118  tmp<volScalarField> fMu() const;
119  tmp<volScalarField> f2() const;
120 
121  virtual void correctNut();
122  virtual tmp<fvScalarMatrix> kSource() const;
123  virtual tmp<fvScalarMatrix> epsilonSource() const;
124 
125 
126 public:
128  typedef typename BasicTurbulenceModel::alphaField alphaField;
129  typedef typename BasicTurbulenceModel::rhoField rhoField;
130  typedef typename BasicTurbulenceModel::transportModel transportModel;
131 
132 
133  //- Runtime type information
134  TypeName("LaunderSharmaKE");
135 
136 
137  // Constructors
138 
139  //- Construct from components
141  (
142  const alphaField& alpha,
143  const rhoField& rho,
144  const volVectorField& U,
145  const surfaceScalarField& alphaRhoPhi,
146  const surfaceScalarField& phi,
147  const transportModel& transport,
148  const word& propertiesName = turbulenceModel::propertiesName,
149  const word& type = typeName
150  );
151 
152 
153  //- Destructor
154  virtual ~LaunderSharmaKE()
155  {}
156 
157 
158  // Member Functions
159 
160  //- Re-read model coefficients if they have changed
161  virtual bool read();
162 
163  //- Return the effective diffusivity for k
164  tmp<volScalarField> DkEff() const
165  {
166  return tmp<volScalarField>
167  (
168  new volScalarField
169  (
170  "DkEff",
171  (this->nut_/sigmak_ + this->nu())
172  )
173  );
174  }
175 
176  //- Return the effective diffusivity for epsilon
178  {
179  return tmp<volScalarField>
180  (
181  new volScalarField
182  (
183  "DepsilonEff",
184  (this->nut_/sigmaEps_ + this->nu())
185  )
186  );
187  }
188 
189  //- Return the turbulence kinetic energy
190  virtual tmp<volScalarField> k() const
191  {
192  return k_;
193  }
194 
195  //- Return the turbulence kinetic energy dissipation rate
196  virtual tmp<volScalarField> epsilon() const
197  {
198  return epsilon_;
199  }
200 
201  //- Solve the turbulence equations and correct the turbulence viscosity
202  virtual void correct();
203 };
204 
205 
206 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
207 
208 } // End namespace RASModels
209 } // End namespace Foam
210 
211 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
212 
213 #ifdef NoRepository
214  #include "LaunderSharmaKE.C"
215 #endif
216 
217 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
218 
219 #endif
220 
221 // ************************************************************************* //
surfaceScalarField & phi
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
TypeName("LaunderSharmaKE")
Runtime type information.
U
Definition: pEqn.H:83
virtual ~LaunderSharmaKE()
Destructor.
virtual tmp< fvScalarMatrix > epsilonSource() const
virtual bool read()
Re-read model coefficients if they have changed.
tmp< volScalarField > DepsilonEff() const
Return the effective diffusivity for epsilon.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:52
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
BasicTurbulenceModel::transportModel transportModel
tmp< volScalarField > DkEff() const
Return the effective diffusivity for k.
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
virtual tmp< fvScalarMatrix > kSource() const
BasicTurbulenceModel::rhoField rhoField
BasicTurbulenceModel::alphaField alphaField
tmp< volScalarField > f2() const
A class for managing temporary objects.
Definition: PtrList.H:54
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
tmp< volScalarField > fMu() const
Launder and Sharma low-Reynolds k-epsilon turbulence model for incompressible and compressible and co...
volScalarField & nu
Namespace for OpenFOAM.