OpenFOAM logo
Open Source CFD Toolkit

LESmodel.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     LESmodel
00027 
00028 Description
00029     Class for all compressible flow LES SGS models.
00030 
00031     This class defines the basic interface for a compressible flow SGS model,
00032     and encapsulates data of value to all possible models. In particular
00033     this includes references to all the dependent fields (rho, U, phi),
00034     the physical viscosity mu, and the turbulenceProperties dictionary
00035     which contains the model selection and model coefficients.
00036 
00037 SourceFiles
00038     LESmodel.C
00039     newLESmodel.C
00040 
00041 \*---------------------------------------------------------------------------*/
00042 
00043 #ifndef compressibleLESmodel_H
00044 #define compressibleLESmodel_H
00045 
00046 #include "LESdelta.H"
00047 #include "fvm.H"
00048 #include "fvc.H"
00049 #include "fvMatrices.H"
00050 #include "basicThermo.H"
00051 #include "bound.H"
00052 #include "autoPtr.H"
00053 #include "runTimeSelectionTables.H"
00054 
00055 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00056 
00057 namespace Foam
00058 {
00059 namespace compressible
00060 {
00061 
00062 /*---------------------------------------------------------------------------*\
00063                            Class LESmodel Declaration
00064 \*---------------------------------------------------------------------------*/
00065 
00066 class LESmodel
00067 :
00068     public IOdictionary
00069 {
00070 
00071 protected:
00072 
00073     // Protected data
00074 
00075         const Time& runTime_;
00076         const fvMesh& mesh_;
00077 
00078 
00079 private:
00080 
00081     // Private data
00082 
00083         const volScalarField& rho_;
00084         const volVectorField& U_;
00085         const surfaceScalarField& phi_;
00086 
00087         const basicThermo& thermoPhysicalModel_;
00088 
00089         dictionary LESmodelProperties_;
00090 
00091         dimensionedScalar k0_;
00092 
00093         autoPtr<LESdelta> delta_;
00094 
00095 
00096     // Private Member Functions
00097 
00098         // Disallow default bitwise copy construct and assignment
00099         LESmodel(const LESmodel&);
00100         LESmodel& operator=(const LESmodel&);
00101 
00102 
00103 public:
00104 
00105     //- Runtime type information
00106     TypeName("LESmodel");
00107 
00108 
00109     // Declare run-time constructor selection table
00110 
00111         declareRunTimeSelectionTable
00112         (
00113             autoPtr,
00114             LESmodel,
00115             dictionary,
00116             (
00117                 const volScalarField& rho,
00118                 const volVectorField& U,
00119                 const surfaceScalarField& phi,
00120                 const basicThermo& thermoPhysicalModel
00121             ),
00122             (rho, U, phi, thermoPhysicalModel)
00123         );
00124 
00125 
00126     // Constructors
00127 
00128         //- Construct from components
00129         LESmodel
00130         (
00131             const word& type,
00132             const volScalarField& rho,
00133             const volVectorField& U,
00134             const surfaceScalarField& phi,
00135             const basicThermo& thermoPhysicalModel
00136         );
00137 
00138 
00139     // Selectors
00140 
00141         //- Return a reference to the selected LES model
00142         static autoPtr<LESmodel> New
00143         (
00144             const volScalarField& rho,
00145             const volVectorField& U,
00146             const surfaceScalarField& phi,
00147             const basicThermo& thermoPhysicalModel
00148         );
00149 
00150 
00151     // Destructor
00152 
00153         virtual ~LESmodel()
00154         {}
00155 
00156 
00157     // Member Functions
00158 
00159         // Access
00160 
00161             //- Access function to the density field
00162             inline const volScalarField& rho() const
00163             {
00164                 return rho_;
00165             }
00166 
00167             //- Access function to velocity field
00168             inline const volVectorField& U() const
00169             {
00170                 return U_;
00171             }
00172 
00173             //- Access function to flux field
00174             inline const surfaceScalarField& phi() const
00175             {
00176                 return phi_;
00177             }
00178 
00179             //- Access function to the thermophysical properties model
00180             inline const basicThermo& thermo() const
00181             {
00182                 return thermoPhysicalModel_;
00183             }
00184 
00185             //- Access the dictionary which provides info. about choice of
00186             //  models, and all related data (particularly model coefficients).
00187             inline const dictionary& LESmodelProperties()
00188             {
00189                 return LESmodelProperties_;
00190             }
00191 
00192             //- Access function to filter width
00193             inline const volScalarField& delta() const
00194             {
00195                 return delta_();
00196             }
00197 
00198             //- Return the value of k0 which k is not allowed to be less than
00199             const dimensionedScalar& k0() const
00200             {
00201                 return k0_;
00202             }
00203 
00204             //- Allow k0 to be changed
00205             dimensionedScalar& k0()
00206             {
00207                 return k0_;
00208             }
00209 
00210             //- Access function to laminar viscosity
00211             tmp<volScalarField> mu() const
00212             {
00213                 return thermoPhysicalModel_.mu();
00214             }
00215 
00216             //- Access function to laminar thermal conductivity
00217             tmp<volScalarField> alpha() const
00218             {
00219                 return thermoPhysicalModel_.alpha();
00220             }
00221 
00222 
00223         //- Return B.
00224         virtual tmp<volTensorField> B() const = 0;
00225 
00226         //- Return the SGS turbulent kinetic energy.
00227         virtual tmp<volScalarField> k() const = 0;
00228 
00229         //- Return the SGS turbulent dissipation.
00230         virtual tmp<volScalarField> epsilon() const = 0;
00231 
00232         //- Return the effective viscosity
00233         virtual tmp<volScalarField> muSgs() const = 0;
00234 
00235         //- Return the effective viscosity
00236         virtual tmp<volScalarField> muEff() const
00237         {
00238             return tmp<volScalarField>
00239             (
00240                 new volScalarField("muEff", muSgs() + mu())
00241             );
00242         }
00243 
00244         //- Return the SGS thermal conductivity.
00245         virtual tmp<volScalarField> alphaEff() const = 0;
00246 
00247         //- Returns div(B).
00248         // This is the additional term due to the filtering of the NSE.
00249         virtual tmp<fvVectorMatrix> divRhoB(volVectorField& U) const = 0;
00250 
00251         //- Correct Eddy-Viscosity and related properties.
00252         //  This calls correct(const tmp<volTensorField>& gradU) by supplying
00253         //  gradU calculated locally.
00254         void correct();
00255 
00256         //- Correct Eddy-Viscosity and related properties
00257         virtual void correct(const tmp<volTensorField>& gradU);
00258 
00259         //- Read turbulenceProperties dictionary
00260         virtual bool read() = 0;
00261 };
00262 
00263 
00264 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00265 
00266 } // End namespace compressible
00267 } // End namespace Foam
00268 
00269 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00270 
00271 #endif
00272 
00273 // ************************************************************************* //
For further information go to www.openfoam.org