1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
26 Foam::incompressible::RASModels::LamBremhorstKE
29 Lam and Bremhorst low-Reynolds number k-epsilon turbulence model
30 for incompressible flows
35 \*---------------------------------------------------------------------------*/
37 #ifndef LamBremhorstKE_H
38 #define LamBremhorstKE_H
40 #include <incompressibleRASModels/RASModel.H>
41 #include <finiteVolume/wallDist.H>
43 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47 namespace incompressible
52 /*---------------------------------------------------------------------------*\
53 Class LamBremhorstKE Declaration
54 \*---------------------------------------------------------------------------*/
62 dimensionedScalar Cmu_;
63 dimensionedScalar C1_;
64 dimensionedScalar C2_;
65 dimensionedScalar sigmaEps_;
68 volScalarField epsilon_;
79 //- Runtime type information
80 TypeName("LamBremhorstKE");
85 //- Construct from components
88 const volVectorField& U,
89 const surfaceScalarField& phi,
90 transportModel& transport
95 virtual ~LamBremhorstKE()
101 //- Return the turbulence viscosity
102 virtual tmp<volScalarField> nut() const
107 //- Return the effective diffusivity for k
108 tmp<volScalarField> DkEff() const
110 return tmp<volScalarField>
112 new volScalarField("DkEff", nut_ + nu())
116 //- Return the effective diffusivity for epsilon
117 tmp<volScalarField> DepsilonEff() const
119 return tmp<volScalarField>
121 new volScalarField("DepsilonEff", nut_/sigmaEps_ + nu())
125 //- Return the turbulence kinetic energy
126 virtual tmp<volScalarField> k() const
131 //- Return the turbulence kinetic energy dissipation rate
132 virtual tmp<volScalarField> epsilon() const
137 //- Return the Reynolds stress tensor
138 virtual tmp<volSymmTensorField> R() const;
140 //- Return the effective stress tensor including the laminar stress
141 virtual tmp<volSymmTensorField> devReff() const;
143 //- Return the source term for the momentum equation
144 virtual tmp<fvVectorMatrix> divDevReff(volVectorField& U) const;
146 //- Solve the turbulence equations and correct the turbulence viscosity
147 virtual void correct();
149 //- Read RASProperties dictionary
154 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
156 } // End namespace RASModels
157 } // End namespace incompressible
158 } // End namespace Foam
160 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
164 // ************************ vim: set sw=4 sts=4 et: ************************ //