OpenFOAM logo
Open Source CFD Toolkit

LRR.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     LRR
00027 
00028 Description
00029     Launder, Reece and Rodi Reynolds-stress turbulence model.
00030 
00031 SourceFiles
00032     LRR.C
00033     LRRcorrect.C
00034 
00035 \*---------------------------------------------------------------------------*/
00036 
00037 #ifndef LRR_H
00038 #define LRR_H
00039 
00040 #include "turbulenceModel.H"
00041 
00042 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00043 
00044 namespace Foam
00045 {
00046 namespace turbulenceModels
00047 {
00048 
00049 /*---------------------------------------------------------------------------*\
00050                            Class LRR Declaration
00051 \*---------------------------------------------------------------------------*/
00052 
00053 class LRR
00054 :
00055     public turbulenceModel
00056 {
00057     // Private data
00058 
00059         dimensionedScalar Cmu;
00060 
00061         dimensionedScalar Clrr1;
00062         dimensionedScalar Clrr2;
00063 
00064         dimensionedScalar C1;
00065         dimensionedScalar C2;
00066         dimensionedScalar Cs;
00067         dimensionedScalar Ceps;
00068         dimensionedScalar alphaEps;
00069 
00070         scalar couplingFactor_;
00071 
00072         volTensorField R_;
00073         volScalarField k_;
00074         volScalarField epsilon_;
00075         volScalarField nut_;
00076 
00077 
00078 public:
00079 
00080     //- Runtime type information
00081     TypeName("LRR");
00082 
00083     // Constructors
00084 
00085         //- from components
00086         LRR
00087         (
00088             const volVectorField& U,
00089             const surfaceScalarField& phi,
00090             transportModel& transport
00091         );
00092 
00093 
00094     // Destructor
00095 
00096         ~LRR(){}
00097 
00098 
00099     // Member Functions
00100 
00101         tmp<volScalarField> nut() const
00102         {
00103             return nut_;
00104         }
00105 
00106         //- Return the effective diffusivity for R
00107         tmp<volScalarField> DREff() const
00108         {
00109             return tmp<volScalarField> 
00110             (
00111                 new volScalarField("DREff", nut_ + nu())
00112             );
00113         }
00114 
00115         //- Return the effective diffusivity for epsilon
00116         tmp<volScalarField> DepsilonEff() const
00117         {
00118             return tmp<volScalarField> 
00119             (
00120                 new volScalarField("DepsilonEff", alphaEps*nut_ + nu())
00121             );
00122         }
00123 
00124         tmp<volScalarField> k() const
00125         {
00126             return k_;
00127         }
00128 
00129         tmp<volScalarField> epsilon() const
00130         {
00131             return epsilon_;
00132         }
00133 
00134         tmp<volTensorField> R() const
00135         {
00136             return R_;
00137         }
00138 
00139         tmp<fvVectorMatrix> divR(volVectorField& U) const;
00140 
00141         void correct();
00142 
00143         //- Read turbulenceProperties dictionary
00144         bool read();
00145 };
00146 
00147 
00148 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00149 
00150 } // End namespace turbulenceModels
00151 } // End namespace Foam
00152 
00153 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00154 
00155 #endif
00156 
00157 // ************************************************************************* //
For further information go to www.openfoam.org