LienCubicKE.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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9 Copyright (C) 2019-2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27Class
28 Foam::incompressible::RASModels::LienCubicKE
29
30Group
31 grpIcoRASTurbulence
32
33Description
34 Lien cubic non-linear low-Reynolds k-epsilon turbulence models for
35 incompressible flows.
36
37 This turbulence model is described in:
38 \verbatim
39 Lien, F.S., Chen, W.L. & Leschziner, M.A. (1996).
40 Low-Reynolds-number eddy-viscosity modeling based on non-linear
41 stress-strain/vorticity relations.
42 Engineering Turbulence Modelling and Experiments 3, 91-100.
43 \endverbatim
44
45 Implemented according to the specification in:
46 <a href=
47 "http://personalpages.manchester.ac.uk/staff/david.d.apsley/specturb.pdf"
48 >Apsley: Turbulence Models 2002</a>
49
50 In addition to the low-Reynolds number damping functions support for
51 wall-functions is also included to allow for low- and high-Reynolds number
52 operation.
53
54See also
55 Foam::incompressible::RASModels::ShihQuadraticKE
56
57SourceFiles
58 LienCubicKE.C
59
60\*---------------------------------------------------------------------------*/
61
62#ifndef LienCubicKE_H
63#define LienCubicKE_H
64
67
68// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
69
70namespace Foam
71{
72namespace incompressible
73{
74namespace RASModels
75{
76
77/*---------------------------------------------------------------------------*\
78 Class LienCubicKE Declaration
79\*---------------------------------------------------------------------------*/
81class LienCubicKE
82:
83 public nonlinearEddyViscosity<incompressible::RASModel>
84{
85
86protected:
87
88 // Protected data
89
90 // Model coefficients
111
112
113 // Fields
117
118 //- Wall distance
119 // Note: different to wall distance in parent RASModel
120 // which is for near-wall cells only
121 const volScalarField& y_;
122
123
124 // Protected Member Functions
125
126 tmp<volScalarField> fMu() const;
127 tmp<volScalarField> f2() const;
129
130 virtual void correctNut();
131 virtual void correctNonlinearStress(const volTensorField& gradU);
132
133
134public:
135
136 //- Runtime type information
137 TypeName("LienCubicKE");
138
139 // Constructors
140
141 //- Construct from components
143 (
145 const geometricOneField& rho,
146 const volVectorField& U,
147 const surfaceScalarField& alphaRhoPhi,
148 const surfaceScalarField& phi,
149 const transportModel& transport,
150 const word& propertiesName = turbulenceModel::propertiesName,
151 const word& type = typeName
152 );
153
154
155 //- Destructor
156 virtual ~LienCubicKE() = default;
157
158
159 // Member Functions
160
161 //- Re-read model coefficients if they have changed
162 virtual bool read();
163
164 //- Return the effective diffusivity for k
166 {
168 (
169 new volScalarField("DkEff", nut_/sigmak_ + nu())
170 );
171 }
172
173 //- Return the effective diffusivity for epsilon
175 {
177 (
178 new volScalarField("DepsilonEff", nut_/sigmaEps_ + nu())
179 );
180 }
181
182 //- Return the turbulence kinetic energy
183 virtual tmp<volScalarField> k() const
184 {
185 return k_;
186 }
187
188 //- Return the turbulence kinetic energy dissipation rate
189 virtual tmp<volScalarField> epsilon() const
190 {
191 return epsilon_;
192 }
193
194 //- Solve the turbulence equations and correct the turbulence viscosity
195 virtual void correct();
196};
197
198
199// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
200
201} // End namespace RASModels
202} // Edn namespace incompressible
203} // End namespace Foam
204
205// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
206
207#endif
208
209// ************************************************************************* //
surfaceScalarField & phi
volScalarField nut_
Definition: eddyViscosity.H:66
BasicTurbulenceModel::transportModel transportModel
Definition: eddyViscosity.H:78
A class representing the concept of a GeometricField of 1 used to avoid unnecessary manipulations for...
Lien cubic non-linear low-Reynolds k-epsilon turbulence models for incompressible flows.
Definition: LienCubicKE.H:83
tmp< volScalarField > f2() const
Definition: LienCubicKE.C:60
virtual ~LienCubicKE()=default
Destructor.
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: LienCubicKE.C:372
TypeName("LienCubicKE")
Runtime type information.
tmp< volScalarField > DepsilonEff() const
Return the effective diffusivity for epsilon.
Definition: LienCubicKE.H:173
const volScalarField & y_
Wall distance.
Definition: LienCubicKE.H:120
tmp< volScalarField > E(const volScalarField &f2) const
Definition: LienCubicKE.C:68
virtual void correctNonlinearStress(const volTensorField &gradU)
Definition: LienCubicKE.C:88
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
Definition: LienCubicKE.H:182
tmp< volScalarField > fMu() const
Definition: LienCubicKE.C:50
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
Definition: LienCubicKE.H:188
virtual bool read()
Re-read model coefficients if they have changed.
Definition: LienCubicKE.C:343
tmp< volScalarField > DkEff() const
Return the effective diffusivity for k.
Definition: LienCubicKE.H:164
Eddy viscosity turbulence model with non-linear correction base class.
static const word propertiesName
Default name of the phase properties dictionary.
Definition: phaseSystem.H:290
A class for managing temporary objects.
Definition: tmp.H:65
A class for handling words, derived from Foam::string.
Definition: word.H:68
U
Definition: pEqn.H:72
Namespace for OpenFOAM.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
volScalarField & nu
volScalarField & alpha
#define TypeName(TypeNameString)
Declare a ClassName() with extra virtual type info.
Definition: typeInfo.H:73