LienLeschziner.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-2015 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::LienLeschziner
26 
27 Group
28  grpIcoRASTurbulence
29 
30 Description
31  Lien and Leschziner low-Reynolds number k-epsilon turbulence model for
32  incompressible flows.
33 
34  This turbulence model is described in:
35  \verbatim
36  Lien, F. S., & Leschziner, M. A. (1993).
37  A pressure-velocity solution strategy for compressible flow
38  and its application to shock/boundary-layer interaction
39  using second-moment turbulence closure.
40  Journal of fluids engineering, 115(4), 717-725.
41  \endverbatim
42 
43  Implemented according to the specification in:
44  <a href=
45  "http://personalpages.manchester.ac.uk/staff/david.d.apsley/specturb.pdf"
46  >Apsley: Turbulence Models 2002</a>
47 
48  In addition to the low-Reynolds number damping functions support for
49  wall-functions is also included to allow for low- and high-Reynolds number
50  operation.
51 
52 SourceFiles
53  LienLeschziner.C
54 
55 \*---------------------------------------------------------------------------*/
56 
57 #ifndef LienLeschziner_H
58 #define LienLeschziner_H
59 
61 #include "eddyViscosity.H"
62 
63 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
64 
65 namespace Foam
66 {
67 namespace incompressible
68 {
69 namespace RASModels
70 {
71 
72 /*---------------------------------------------------------------------------*\
73  Class LienLeschziner Declaration
74 \*---------------------------------------------------------------------------*/
75 
76 class LienLeschziner
77 :
78  public eddyViscosity<incompressible::RASModel>
79 {
80 
81 protected:
82 
83  // Protected data
84 
85  // Model coefficients
86 
91 
94 
98 
99 
100  // Fields
104 
105  //- Wall distance
106  // Note: different to wall distance in parent RASModel
107  // which is for near-wall cells only
108  const volScalarField& y_;
109 
110 
111  // Protected Member Functions
112 
113  tmp<volScalarField> fMu() const;
114  tmp<volScalarField> f2() const;
115  tmp<volScalarField> E(const volScalarField& f2) const;
116 
117  virtual void correctNut();
118 
119 
120 public:
121 
122  TypeName("LienLeschziner");
123 
124  // Constructors
125 
126  //- Construct from components
128  (
129  const geometricOneField& alpha,
130  const geometricOneField& rho,
131  const volVectorField& U,
132  const surfaceScalarField& alphaRhoPhi,
133  const surfaceScalarField& phi,
134  const transportModel& transport,
135  const word& propertiesName = turbulenceModel::propertiesName,
136  const word& type = typeName
137  );
138 
139 
140  //- Destructor
141  virtual ~LienLeschziner()
142  {}
143 
144 
145  // Member Functions
146 
147  //- Read RASProperties dictionary
148  virtual bool read();
149 
150  //- Return the effective diffusivity for k
151  tmp<volScalarField> DkEff() const
152  {
153  return tmp<volScalarField>
154  (
155  new volScalarField("DkEff", nut_/sigmak_ + nu())
156  );
157  }
158 
159  //- Return the effective diffusivity for epsilon
161  {
162  return tmp<volScalarField>
163  (
164  new volScalarField("DepsilonEff", nut_/sigmaEps_ + nu())
165  );
166  }
167 
168  //- Return the turbulence kinetic energy
169  virtual tmp<volScalarField> k() const
170  {
171  return k_;
172  }
173 
174  //- Return the turbulence kinetic energy dissipation rate
175  virtual tmp<volScalarField> epsilon() const
176  {
177  return epsilon_;
178  }
179 
180  //- Solve the turbulence equations and correct the turbulence viscosity
181  virtual void correct();
182 };
183 
184 
185 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
186 
187 } // End namespace RASModels
188 } // End namespace incompressible
189 } // End namespace Foam
190 
191 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
192 
193 #endif
194 
195 // ************************************************************************* //
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
surfaceScalarField & phi
U
Definition: pEqn.H:83
LienLeschziner(const geometricOneField &alpha, const geometricOneField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
tmp< volScalarField > E(const volScalarField &f2) const
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:52
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...
tmp< volScalarField > DepsilonEff() const
Return the effective diffusivity for epsilon.
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
virtual bool read()
Read RASProperties dictionary.
Lien and Leschziner low-Reynolds number k-epsilon turbulence model for incompressible flows...
const volScalarField & y_
Wall distance.
fileName::Type type(const fileName &, const bool followLink=true)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:485
tmp< volScalarField > DkEff() const
Return the effective diffusivity for k.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
volScalarField & nu
Namespace for OpenFOAM.