FreeFOAM The Cross-Platform CFD Toolkit
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) 1991-2010 OpenCFD Ltd.
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::compressible::RASModels::RNGkEpsilon
26 
27 Description
28  Renormalisation group k-epsilon turbulence model for compressible flows.
29 
30  The default model coefficients correspond to the following:
31  @verbatim
32  RNGkEpsilonCoeffs
33  {
34  Cmu 0.0845;
35  C1 1.42;
36  C2 1.68;
37  C3 -0.33; // only for compressible
38  Prt 1.0; // only for compressible
39  sigmak 0.71942;
40  sigmaEps 0.71942;
41  eta0 4.38;
42  beta 0.012;
43  }
44  @endverbatim
45 
46 SourceFiles
47  RNGkEpsilon.C
48  RNGkEpsilonCorrect.C
49 
50 \*---------------------------------------------------------------------------*/
51 
52 #ifndef compressibleRNGkEpsilon_H
53 #define compressibleRNGkEpsilon_H
54 
56 
57 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
58 
59 namespace Foam
60 {
61 namespace compressible
62 {
63 namespace RASModels
64 {
65 
66 /*---------------------------------------------------------------------------*\
67  Class RNGkEpsilon Declaration
68 \*---------------------------------------------------------------------------*/
69 
71 :
72  public RASModel
73 {
74  // Private data
75 
76  // Model coefficients
77 
78  dimensionedScalar Cmu_;
82  dimensionedScalar sigmak_;
83  dimensionedScalar sigmaEps_;
84  dimensionedScalar Prt_;
85  dimensionedScalar eta0_;
86  dimensionedScalar beta_;
87 
88  // Fields
89 
90  volScalarField k_;
91  volScalarField epsilon_;
92  volScalarField mut_;
93  volScalarField alphat_;
94 
95 
96 public:
97 
98  //- Runtime type information
99  TypeName("RNGkEpsilon");
100 
101  // Constructors
102 
103  //- Construct from components
105  (
106  const volScalarField& rho,
107  const volVectorField& U,
108  const surfaceScalarField& phi,
109  const basicThermo& thermophysicalModel
110  );
111 
112 
113  //- Destructor
114  virtual ~RNGkEpsilon()
115  {}
116 
117 
118  // Member Functions
119 
120  //- Return the effective diffusivity for k
122  {
123  return tmp<volScalarField>
124  (
125  new volScalarField("DkEff", mut_/sigmak_ + mu())
126  );
127  }
128 
129  //- Return the effective diffusivity for epsilon
131  {
132  return tmp<volScalarField>
133  (
134  new volScalarField("DepsilonEff", mut_/sigmaEps_ + mu())
135  );
136  }
137 
138  //- Return the turbulence viscosity
139  virtual tmp<volScalarField> mut() const
140  {
141  return mut_;
142  }
143 
144  //- Return the effective turbulent thermal diffusivity
146  {
147  return tmp<volScalarField>
148  (
149  new volScalarField("alphaEff", alphat_ + alpha())
150  );
151  }
152 
153  //- Return the turbulence kinetic energy
154  virtual tmp<volScalarField> k() const
155  {
156  return k_;
157  }
158 
159  //- Return the turbulence kinetic energy dissipation rate
160  virtual tmp<volScalarField> epsilon() const
161  {
162  return epsilon_;
163  }
164 
165  //- Return the Reynolds stress tensor
166  virtual tmp<volSymmTensorField> R() const;
167 
168  //- Return the effective stress tensor including the laminar stress
169  virtual tmp<volSymmTensorField> devRhoReff() const;
170 
171  //- Return the effective stress tensor including the laminar stress
172  virtual tmp<fvVectorMatrix> divDevRhoReff(volVectorField& U) const;
173 
174  //- Solve the turbulence equations and correct the turbulence viscosity
175  virtual void correct();
176 
177  //- Read RASProperties dictionary
178  virtual bool read();
179 };
180 
181 
182 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
183 
184 } // End namespace RASModels
185 } // End namespace compressible
186 } // End namespace Foam
187 
188 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
189 
190 #endif
191 
192 // ************************ vim: set sw=4 sts=4 et: ************************ //