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