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