OpenFOAM logo
Open Source CFD Toolkit

GenSGSStress.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     GenSGSStress
00027 
00028 Description
00029     General base class for all models that directly solve for the SGS stress
00030     tensor B.  Contains tensor fields B (the SGS stress tensor) as well as
00031     scalar fields for k (SGS turbulent energy) gamma (SGS viscosity) and
00032     epsilon (SGS dissipation).
00033 
00034 SourceFiles
00035     GenSGSStress.C
00036 
00037 \*---------------------------------------------------------------------------*/
00038 
00039 #ifndef compressibleGenSGSStress_H
00040 #define compressibleGenSGSStress_H
00041 
00042 #include "LESmodel.H"
00043 
00044 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00045 
00046 namespace Foam
00047 {
00048 namespace compressible
00049 {
00050 namespace LESmodels
00051 {
00052 
00053 /*---------------------------------------------------------------------------*\
00054                            Class GenSGSStress Declaration
00055 \*---------------------------------------------------------------------------*/
00056 
00057 class GenSGSStress
00058 :
00059     virtual public LESmodel
00060 {
00061     // Private Member Functions
00062 
00063         // Disallow default bitwise copy construct and assignment
00064         GenSGSStress(const GenSGSStress&);
00065         GenSGSStress& operator=(const GenSGSStress&);
00066 
00067 
00068 protected:
00069 
00070         dimensionedScalar ce_;
00071 
00072         volTensorField B_;
00073         volScalarField muSgs_;
00074 
00075 
00076 public:
00077 
00078     // Constructors
00079 
00080         //- Constructor from components
00081         GenSGSStress
00082         (
00083             const volScalarField& rho,
00084             const volVectorField& U,
00085             const surfaceScalarField& phi,
00086             const basicThermo& thermoPhysicalModel
00087         );
00088 
00089 
00090     // Destructor
00091 
00092         virtual ~GenSGSStress()
00093         {}
00094 
00095 
00096     // Member Functions
00097 
00098         //- Return the SGS stress.
00099         virtual tmp<volTensorField> B() const
00100         {
00101             return B_;
00102         }
00103 
00104         //- Return the SGS turbulent kinetic energy.
00105         virtual tmp<volScalarField> k() const
00106         {
00107             return 0.5*tr(B_);
00108         }
00109 
00110         //- Return the SGS turbulent dissipation.
00111         virtual tmp<volScalarField> epsilon() const
00112         {
00113             volScalarField K = k();
00114             return ce_*K*sqrt(K)/delta();
00115         }
00116 
00117         //- Return the SGS viscosity.
00118         virtual tmp<volScalarField> muSgs() const
00119         {
00120             return muSgs_;
00121         }
00122 
00123         //- Return thermal conductivity
00124         virtual tmp<volScalarField> alphaEff() const
00125         {
00126             return tmp<volScalarField>
00127             (
00128                 new volScalarField("alphaEff", muSgs_ + alpha())
00129             );
00130         }
00131 
00132         //- Returns divergence of B : i.e. the additional term in the
00133         //  filtered NSE.
00134         virtual tmp<fvVectorMatrix> divRhoB(volVectorField& U) const;
00135 
00136         //- Correct Eddy-Viscosity and related properties
00137         virtual void correct(const tmp<volTensorField>& gradU);
00138 
00139         //- Read turbulenceProperties dictionary
00140         bool read();
00141 };
00142 
00143 
00144 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00145 
00146 } // End namespace LESmodels
00147 } // End namespace compressible
00148 } // End namespace Foam
00149 
00150 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00151 
00152 #endif
00153 
00154 // ************************************************************************* //
For further information go to www.openfoam.org