OpenFOAM logo
Open Source CFD Toolkit

LienCubicKELowRe.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     LienCubicKELowRe
00027 
00028 Description
00029     Lien cubic non-linear low-Reynolds k-epsilon turbulence models.
00030 
00031 SourceFiles
00032     LienCubicKELowRe.C
00033     LienCubicKELowReCorrect.C
00034 
00035 \*---------------------------------------------------------------------------*/
00036 
00037 #ifndef LienCubicKELowRe_H
00038 #define LienCubicKELowRe_H
00039 
00040 #include "turbulenceModel.H"
00041 #include "wallDist.H"
00042 
00043 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00044 
00045 namespace Foam
00046 {
00047 namespace turbulenceModels
00048 {
00049 
00050 /*---------------------------------------------------------------------------*\
00051                            Class LienCubicKELowRe Declaration
00052 \*---------------------------------------------------------------------------*/
00053 
00054 class LienCubicKELowRe
00055 :
00056     public turbulenceModel
00057 {
00058     // Private data
00059 
00060         dimensionedScalar C1;
00061         dimensionedScalar C2;
00062         dimensionedScalar alphak;
00063         dimensionedScalar alphaEps;
00064         dimensionedScalar A1;
00065         dimensionedScalar A2;
00066         dimensionedScalar Ctau1;
00067         dimensionedScalar Ctau2;
00068         dimensionedScalar Ctau3;
00069         dimensionedScalar alphaKsi;
00070 
00071         dimensionedScalar CmuWall;
00072 
00073         dimensionedScalar Am;
00074         dimensionedScalar Aepsilon;
00075         dimensionedScalar Amu;
00076 
00077         volScalarField k_;
00078         volScalarField epsilon_;
00079 
00080         wallDist y_;
00081 
00082         volTensorField gradU;
00083         volScalarField eta;
00084         volScalarField ksi;
00085         volScalarField Cmu;
00086         volScalarField fEta;
00087         volScalarField C5viscosity;
00088 
00089         volScalarField yStar;
00090 
00091         volScalarField nut_;
00092 
00093         volTensorField nonlinearStress;
00094 
00095 public:
00096 
00097     //- Runtime type information
00098     TypeName("LienCubicKELowRe");
00099 
00100     // Constructors
00101 
00102         //- from components
00103         LienCubicKELowRe
00104         (
00105             const volVectorField& U,
00106             const surfaceScalarField& phi,
00107             transportModel& transport
00108         );
00109 
00110 
00111     // Destructor
00112 
00113         ~LienCubicKELowRe(){}
00114 
00115 
00116     // Member Functions
00117 
00118         tmp<volScalarField> nut() const
00119         {
00120             return nut_;
00121         }
00122 
00123         //- Return the effective diffusivity for k
00124         tmp<volScalarField> DkEff() const
00125         {
00126             return tmp<volScalarField> 
00127             (
00128                 new volScalarField("DkEff", alphak*nut_ + nu())
00129             );
00130         }
00131 
00132         //- Return the effective diffusivity for epsilon
00133         tmp<volScalarField> DepsilonEff() const
00134         {
00135             return tmp<volScalarField> 
00136             (
00137                 new volScalarField("DepsilonEff", alphaEps*nut_ + nu())
00138             );
00139         }
00140 
00141         tmp<volScalarField> k() const
00142         {
00143             return k_;
00144         }
00145 
00146         tmp<volScalarField> epsilon() const
00147         {
00148             return epsilon_;
00149         }
00150 
00151         tmp<volTensorField> R() const;
00152 
00153         tmp<fvVectorMatrix> divR(volVectorField& U) const;
00154 
00155         void correct();
00156 
00157         //- Read turbulenceProperties dictionary
00158         bool read();
00159 };
00160 
00161 
00162 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00163 
00164 } // End namespace turbulenceModels
00165 } // End namespace Foam
00166 
00167 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00168 
00169 #endif
00170 
00171 // ************************************************************************* //
For further information go to www.openfoam.org