RNGkEpsilon.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::RNGkEpsilon
26 
27 Group
28  grpRASTurbulence
29 
30 Description
31  Renormalization group k-epsilon turbulence model for incompressible and
32  compressible flows.
33 
34  Reference:
35  \verbatim
36  Yakhot, V., Orszag, S. A., Thangam, S.,
37  Gatski, T. B., & Speziale, C. G. (1992).
38  Development of turbulence models for shear flows
39  by a double expansion technique.
40  Physics of Fluids A: Fluid Dynamics (1989-1993), 4(7), 1510-1520.
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  RNGkEpsilonCoeffs
51  {
52  Cmu 0.0845;
53  C1 1.42;
54  C2 1.68;
55  C3 -0.33;
56  sigmak 0.71942;
57  sigmaEps 0.71942;
58  eta0 4.38;
59  beta 0.012;
60  }
61  \endverbatim
62 
63 SourceFiles
64  RNGkEpsilon.C
65 
66 \*---------------------------------------------------------------------------*/
67 
68 #ifndef RNGkEpsilon_H
69 #define RNGkEpsilon_H
70 
71 #include "RASModel.H"
72 #include "eddyViscosity.H"
73 
74 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
75 
76 namespace Foam
77 {
78 namespace RASModels
79 {
80 
81 /*---------------------------------------------------------------------------*\
82  Class RNGkEpsilon Declaration
83 \*---------------------------------------------------------------------------*/
84 
85 template<class BasicTurbulenceModel>
86 class RNGkEpsilon
87 :
88  public eddyViscosity<RASModel<BasicTurbulenceModel>>
89 {
90  // Private Member Functions
91 
92  // Disallow default bitwise copy construct and assignment
93  RNGkEpsilon(const RNGkEpsilon&);
94  void operator=(const RNGkEpsilon&);
95 
96 
97 protected:
98 
99  // Protected data
100 
101  // Model coefficients
111 
112 
113  // Fields
117 
118 
119  // Protected Member Functions
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("RNGkEpsilon");
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 ~RNGkEpsilon()
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 "RNGkEpsilon.C"
215 #endif
216 
217 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
218 
219 #endif
220 
221 // ************************************************************************* //
surfaceScalarField & phi
dimensionedScalar beta_
Definition: RNGkEpsilon.H:109
U
Definition: pEqn.H:83
TypeName("RNGkEpsilon")
Runtime type information.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: RNGkEpsilon.C:241
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:52
tmp< volScalarField > DepsilonEff() const
Return the effective diffusivity for epsilon.
Definition: RNGkEpsilon.H:176
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dimensionedScalar sigmak_
Definition: RNGkEpsilon.H:106
tmp< volScalarField > DkEff() const
Return the effective diffusivity for k.
Definition: RNGkEpsilon.H:163
static const word propertiesName
Default name of the turbulence properties dictionary.
BasicTurbulenceModel::alphaField alphaField
Definition: RNGkEpsilon.H:127
A class for handling words, derived from string.
Definition: word.H:59
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
Definition: RNGkEpsilon.H:189
virtual ~RNGkEpsilon()
Destructor.
Definition: RNGkEpsilon.H:153
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
Definition: RNGkEpsilon.H:195
virtual tmp< fvScalarMatrix > kSource() const
Definition: RNGkEpsilon.C:51
dimensionedScalar Cmu_
Definition: RNGkEpsilon.H:102
dimensionedScalar eta0_
Definition: RNGkEpsilon.H:108
Renormalization group k-epsilon turbulence model for incompressible and compressible flows...
Definition: RNGkEpsilon.H:85
virtual tmp< fvScalarMatrix > epsilonSource() const
Definition: RNGkEpsilon.C:66
virtual bool read()
Re-read model coefficients if they have changed.
Definition: RNGkEpsilon.C:218
BasicTurbulenceModel::rhoField rhoField
Definition: RNGkEpsilon.H:128
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
BasicTurbulenceModel::transportModel transportModel
Definition: RNGkEpsilon.H:129
dimensionedScalar sigmaEps_
Definition: RNGkEpsilon.H:107
A class for managing temporary objects.
Definition: PtrList.H:54
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
volScalarField & nu
Namespace for OpenFOAM.