OpenFOAM logo
Open Source CFD Toolkit

LaunderSharmaKE.H

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2005 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software; you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by the
00013     Free Software Foundation; either version 2 of the License, or (at your
00014     option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM; if not, write to the Free Software Foundation,
00023     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
00024 
00025 Class
00026     LaunderSharmaKE
00027 
00028 Description
00029     Launder and Sharma low-Reynolds number k-epsilon turbulence model.
00030 
00031 SourceFiles
00032     LaunderSharmaKE.C
00033     LaunderSharmaKECorrect.C
00034 
00035 \*---------------------------------------------------------------------------*/
00036 
00037 #ifndef LaunderSharmaKE_H
00038 #define LaunderSharmaKE_H
00039 
00040 #include "turbulenceModel.H"
00041 
00042 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00043 
00044 namespace Foam
00045 {
00046 namespace turbulenceModels
00047 {
00048 
00049 /*---------------------------------------------------------------------------*\
00050                            Class LaunderSharmaKE Declaration
00051 \*---------------------------------------------------------------------------*/
00052 
00053 class LaunderSharmaKE
00054 :
00055     public turbulenceModel
00056 {
00057     // Private data
00058 
00059         dimensionedScalar Cmu;
00060         dimensionedScalar C1;
00061         dimensionedScalar C2;
00062         dimensionedScalar alphaEps;
00063 
00064         volScalarField k_;
00065         volScalarField epsilonTilda_;
00066 
00067         volScalarField nut_;
00068 
00069 
00070     // Private member functions
00071 
00072         tmp<volScalarField> fMu() const;
00073         tmp<volScalarField> f2() const;
00074 
00075 
00076 public:
00077 
00078     //- Runtime type information
00079     TypeName("LaunderSharmaKE");
00080 
00081     // Constructors
00082 
00083         //- from components
00084         LaunderSharmaKE
00085         (
00086             const volVectorField& U,
00087             const surfaceScalarField& phi,
00088             transportModel& transport
00089         );
00090 
00091 
00092     // Destructor
00093 
00094         ~LaunderSharmaKE(){}
00095 
00096 
00097     // Member Functions
00098 
00099         tmp<volScalarField> nut() const
00100         {
00101             return nut_;
00102         }
00103 
00104         //- Return the effective diffusivity for k
00105         tmp<volScalarField> DkEff() const
00106         {
00107             return tmp<volScalarField> 
00108             (
00109                 new volScalarField("DkEff", nut_ + nu())
00110             );
00111         }
00112 
00113         //- Return the effective diffusivity for epsilon
00114         tmp<volScalarField> DepsilonEff() const
00115         {
00116             return tmp<volScalarField> 
00117             (
00118                 new volScalarField("DepsilonEff", alphaEps*nut_ + nu())
00119             );
00120         }
00121 
00122         tmp<volScalarField> k() const
00123         {
00124             return k_;
00125         }
00126 
00127         //- Note that epsilonTilda is returned as epsilon.
00128         //  This is the appropriate variable for most purposes.
00129         tmp<volScalarField> epsilon() const
00130         {
00131             return epsilonTilda_;
00132         }
00133 
00134         tmp<volTensorField> R() const;
00135 
00136         tmp<fvVectorMatrix> divR(volVectorField& U) const;
00137 
00138         void correct();
00139 
00140         //- Read turbulenceProperties dictionary
00141         bool read();
00142 };
00143 
00144 
00145 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00146 
00147 } // End namespace turbulenceModels
00148 } // End namespace Foam
00149 
00150 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00151 
00152 #endif
00153 
00154 // ************************************************************************* //
For further information go to www.openfoam.org