RNGkEpsilon.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-2018 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 Description
28  Renormalization group k-epsilon turbulence model for incompressible and
29  compressible flows.
30 
31  Reference:
32  \verbatim
33  Yakhot, V., Orszag, S. A., Thangam, S.,
34  Gatski, T. B., & Speziale, C. G. (1992).
35  Development of turbulence models for shear flows
36  by a double expansion technique.
37  Physics of Fluids A: Fluid Dynamics (1989-1993), 4(7), 1510-1520.
38 
39  For the RDT-based compression term:
40  El Tahry, S. H. (1983).
41  k-epsilon equation for compressible reciprocating engine flows.
42  Journal of Energy, 7(4), 345-353.
43  \endverbatim
44 
45  The default model coefficients are
46  \verbatim
47  RNGkEpsilonCoeffs
48  {
49  Cmu 0.0845;
50  C1 1.42;
51  C2 1.68;
52  C3 -0.33;
53  sigmak 0.71942;
54  sigmaEps 0.71942;
55  eta0 4.38;
56  beta 0.012;
57  }
58  \endverbatim
59 
60 SourceFiles
61  RNGkEpsilon.C
62 
63 \*---------------------------------------------------------------------------*/
64 
65 #ifndef RNGkEpsilon_H
66 #define RNGkEpsilon_H
67 
68 #include "RASModel.H"
69 #include "eddyViscosity.H"
70 
71 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
72 
73 namespace Foam
74 {
75 namespace RASModels
76 {
77 
78 /*---------------------------------------------------------------------------*\
79  Class RNGkEpsilon Declaration
80 \*---------------------------------------------------------------------------*/
81 
82 template<class BasicTurbulenceModel>
83 class RNGkEpsilon
84 :
85  public eddyViscosity<RASModel<BasicTurbulenceModel>>
86 {
87  // Private Member Functions
88 
89  // Disallow default bitwise copy construct and assignment
90  RNGkEpsilon(const RNGkEpsilon&);
91  void operator=(const RNGkEpsilon&);
92 
93 
94 protected:
95 
96  // Protected data
97 
98  // Model coefficients
99 
108 
109 
110  // Fields
114 
115 
116  // Protected Member Functions
117 
118  virtual void correctNut();
119  virtual tmp<fvScalarMatrix> kSource() const;
120  virtual tmp<fvScalarMatrix> epsilonSource() const;
121 
122 
123 public:
125  typedef typename BasicTurbulenceModel::alphaField alphaField;
126  typedef typename BasicTurbulenceModel::rhoField rhoField;
127  typedef typename BasicTurbulenceModel::transportModel transportModel;
128 
129 
130  //- Runtime type information
131  TypeName("RNGkEpsilon");
132 
133 
134  // Constructors
135 
136  //- Construct from components
138  (
139  const alphaField& alpha,
140  const rhoField& rho,
141  const volVectorField& U,
142  const surfaceScalarField& alphaRhoPhi,
143  const surfaceScalarField& phi,
144  const transportModel& transport,
145  const word& propertiesName = turbulenceModel::propertiesName,
146  const word& type = typeName
147  );
148 
149 
150  //- Destructor
151  virtual ~RNGkEpsilon()
152  {}
153 
154 
155  // Member Functions
156 
157  //- Re-read model coefficients if they have changed
158  virtual bool read();
159 
160  //- Return the effective diffusivity for k
161  tmp<volScalarField> DkEff() const
162  {
163  return tmp<volScalarField>
164  (
165  new volScalarField
166  (
167  "DkEff",
168  (this->nut_/sigmak_ + this->nu())
169  )
170  );
171  }
172 
173  //- Return the effective diffusivity for epsilon
175  {
176  return tmp<volScalarField>
177  (
178  new volScalarField
179  (
180  "DepsilonEff",
181  (this->nut_/sigmaEps_ + this->nu())
182  )
183  );
184  }
185 
186  //- Return the turbulence kinetic energy
187  virtual tmp<volScalarField> k() const
188  {
189  return k_;
190  }
191 
192  //- Return the turbulence kinetic energy dissipation rate
193  virtual tmp<volScalarField> epsilon() const
194  {
195  return epsilon_;
196  }
197 
198  //- Solve the turbulence equations and correct the turbulence viscosity
199  virtual void correct();
200 };
201 
202 
203 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
204 
205 } // End namespace RASModels
206 } // End namespace Foam
207 
208 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
209 
210 #ifdef NoRepository
211  #include "RNGkEpsilon.C"
212 #endif
213 
214 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
215 
216 #endif
217 
218 // ************************************************************************* //
dimensionedScalar beta_
Definition: RNGkEpsilon.H:106
surfaceScalarField & phi
TypeName("RNGkEpsilon")
Runtime type information.
virtual tmp< fvScalarMatrix > kSource() const
Definition: RNGkEpsilon.C:51
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:49
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual tmp< fvScalarMatrix > epsilonSource() const
Definition: RNGkEpsilon.C:66
tmp< volScalarField > DkEff() const
Return the effective diffusivity for k.
Definition: RNGkEpsilon.H:160
dimensionedScalar sigmak_
Definition: RNGkEpsilon.H:103
static const word propertiesName
Default name of the turbulence properties dictionary.
BasicTurbulenceModel::alphaField alphaField
Definition: RNGkEpsilon.H:124
A class for handling words, derived from string.
Definition: word.H:59
tmp< volScalarField > DepsilonEff() const
Return the effective diffusivity for epsilon.
Definition: RNGkEpsilon.H:173
virtual ~RNGkEpsilon()
Destructor.
Definition: RNGkEpsilon.H:150
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
Definition: RNGkEpsilon.H:186
dimensionedScalar Cmu_
Definition: RNGkEpsilon.H:99
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
Definition: RNGkEpsilon.H:192
dimensionedScalar eta0_
Definition: RNGkEpsilon.H:105
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:481
Renormalization group k-epsilon turbulence model for incompressible and compressible flows...
Definition: RNGkEpsilon.H:82
virtual bool read()
Re-read model coefficients if they have changed.
Definition: RNGkEpsilon.C:218
U
Definition: pEqn.H:72
BasicTurbulenceModel::rhoField rhoField
Definition: RNGkEpsilon.H:125
BasicTurbulenceModel::transportModel transportModel
Definition: RNGkEpsilon.H:126
dimensionedScalar sigmaEps_
Definition: RNGkEpsilon.H:104
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
volScalarField & nu
Namespace for OpenFOAM.