LamBremhorstKE.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::incompressible::RASModels::LamBremhorstKE
26 
27 Group
28  grpIcoRASTurbulence
29 
30 Description
31  Lam and Bremhorst low-Reynolds number k-epsilon turbulence model
32  for incompressible flows
33 
34  This turbulence model is described in:
35  \verbatim
36  Lam, C. K. G., & Bremhorst, K. (1981).
37  A modified form of the k-ε model for predicting wall turbulence.
38  Journal of Fluids Engineering, 103(3), 456-460.
39  \endverbatim
40 
41 SourceFiles
42  LamBremhorstKE.C
43 
44 \*---------------------------------------------------------------------------*/
45 
46 #ifndef LamBremhorstKE_H
47 #define LamBremhorstKE_H
48 
50 #include "eddyViscosity.H"
51 
52 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
53 
54 namespace Foam
55 {
56 namespace incompressible
57 {
58 namespace RASModels
59 {
60 
61 /*---------------------------------------------------------------------------*\
62  Class LamBremhorstKE Declaration
63 \*---------------------------------------------------------------------------*/
64 
65 class LamBremhorstKE
66 :
67  public eddyViscosity<incompressible::RASModel>
68 {
69  // Private Member Functions
70 
71  // Disallow default bitwise copy construct and assignment
73  void operator=(const LamBremhorstKE&);
74 
75  tmp<volScalarField> Rt() const;
76  tmp<volScalarField> fMu(const volScalarField& Rt) const;
77  tmp<volScalarField> f1(const volScalarField& fMu) const;
78  tmp<volScalarField> f2(const volScalarField& Rt) const;
79  void correctNut(const volScalarField& fMu);
80 
81 
82 protected:
83 
84  // Protected data
85 
90 
93 
94  //- Wall distance
95  // Note: different to wall distance in parent RASModel
96  // which is for near-wall cells only
97  const volScalarField& y_;
98 
99 
100  // Protected Member Functions
101 
102  virtual void correctNut();
103 
104 
105 public:
106 
107  //- Runtime type information
108  TypeName("LamBremhorstKE");
109 
110 
111  // Constructors
112 
113  //- Construct from components
115  (
116  const geometricOneField& alpha,
117  const geometricOneField& rho,
118  const volVectorField& U,
119  const surfaceScalarField& alphaRhoPhi,
120  const surfaceScalarField& phi,
121  const transportModel& transport,
122  const word& propertiesName = turbulenceModel::propertiesName,
123  const word& type = typeName
124  );
125 
126 
127  //- Destructor
128  virtual ~LamBremhorstKE()
129  {}
130 
131 
132  // Member Functions
133 
134  //- Read RASProperties dictionary
135  virtual bool read();
136 
137  //- Return the effective diffusivity for k
138  tmp<volScalarField> DkEff() const
139  {
140  return tmp<volScalarField>
141  (
142  new volScalarField("DkEff", nut_ + nu())
143  );
144  }
145 
146  //- Return the effective diffusivity for epsilon
148  {
149  return tmp<volScalarField>
150  (
151  new volScalarField("DepsilonEff", nut_/sigmaEps_ + nu())
152  );
153  }
154 
155  //- Return the turbulence kinetic energy
156  virtual tmp<volScalarField> k() const
157  {
158  return k_;
159  }
160 
161  //- Return the turbulence kinetic energy dissipation rate
162  virtual tmp<volScalarField> epsilon() const
163  {
164  return epsilon_;
165  }
166 
167  //- Solve the turbulence equations and correct the turbulence viscosity
168  virtual void correct();
169 };
170 
171 
172 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
173 
174 } // End namespace RASModels
175 } // End namespace incompressible
176 } // End namespace Foam
177 
178 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
179 
180 #endif
181 
182 // ************************************************************************* //
surfaceScalarField & phi
virtual bool read()
Read RASProperties dictionary.
U
Definition: pEqn.H:83
Lam and Bremhorst low-Reynolds number k-epsilon turbulence model for incompressible flows...
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:52
TypeName("LamBremhorstKE")
Runtime type information.
incompressible::RASModel::transportModel transportModel
Definition: eddyViscosity.H:75
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
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 > DkEff() const
Return the effective diffusivity for k.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
const volScalarField & y_
Wall distance.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:461
A class for managing temporary objects.
Definition: PtrList.H:54
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
volScalarField & nu
Namespace for OpenFOAM.
tmp< volScalarField > DepsilonEff() const
Return the effective diffusivity for epsilon.