FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

kOmegaSST.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-2010 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
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your 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, see <http://www.gnu.org/licenses/>.
00023 
00024 Class
00025     Foam::compressible::RASModels::kOmegaSST
00026 
00027 Description
00028     Implementation of the k-omega-SST turbulence model for compressible flows.
00029 
00030     Turbulence model described in:
00031     @verbatim
00032         Menter, F., Esch, T.
00033         "Elements of Industrial Heat Transfer Prediction"
00034         16th Brazilian Congress of Mechanical Engineering (COBEM),
00035         Nov. 2001
00036     @endverbatim
00037 
00038     Note that this implementation is written in terms of alpha diffusion
00039     coefficients rather than the more traditional sigma (alpha = 1/sigma) so
00040     that the blending can be applied to all coefficuients in a consistent
00041     manner.  The paper suggests that sigma is blended but this would not be
00042     consistent with the blending of the k-epsilon and k-omega models.
00043 
00044     Also note that the error in the last term of equation (2) relating to
00045     sigma has been corrected.
00046 
00047     Wall-functions are applied in this implementation by using equations (14)
00048     to specify the near-wall omega as appropriate.
00049 
00050     The blending functions (15) and (16) are not currently used because of the
00051     uncertainty in their origin, range of applicability and that is y+ becomes
00052     sufficiently small blending u_tau in this manner clearly becomes nonsense.
00053 
00054     The default model coefficients correspond to the following:
00055     @verbatim
00056         kOmegaSSTCoeffs
00057         {
00058             alphaK1     0.85034;
00059             alphaK2     1.0;
00060             alphaOmega1 0.5;
00061             alphaOmega2 0.85616;
00062             Prt         1.0;    // only for compressible
00063             beta1       0.075;
00064             beta2       0.0828;
00065             betaStar    0.09;
00066             gamma1      0.5532;
00067             gamma2      0.4403;
00068             a1          0.31;
00069             c1          10.0;
00070         }
00071     @endverbatim
00072 
00073 SourceFiles
00074     kOmegaSST.C
00075     kOmegaWallFunctionsI.H
00076     kOmegaWallViscosityI.H
00077     wallOmegaI.H
00078 
00079 \*---------------------------------------------------------------------------*/
00080 
00081 #ifndef compressiblekOmegaSST_H
00082 #define compressiblekOmegaSST_H
00083 
00084 #include <compressibleRASModels/RASModel.H>
00085 #include <finiteVolume/wallDist.H>
00086 
00087 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00088 
00089 namespace Foam
00090 {
00091 namespace compressible
00092 {
00093 namespace RASModels
00094 {
00095 
00096 /*---------------------------------------------------------------------------*\
00097                           Class kOmegaSST Declaration
00098 \*---------------------------------------------------------------------------*/
00099 
00100 class kOmegaSST
00101 :
00102     public RASModel
00103 {
00104     // Private data
00105 
00106         // Model coefficients
00107 
00108             dimensionedScalar alphaK1_;
00109             dimensionedScalar alphaK2_;
00110 
00111             dimensionedScalar alphaOmega1_;
00112             dimensionedScalar alphaOmega2_;
00113 
00114             dimensionedScalar Prt_;
00115 
00116             dimensionedScalar gamma1_;
00117             dimensionedScalar gamma2_;
00118 
00119             dimensionedScalar beta1_;
00120             dimensionedScalar beta2_;
00121 
00122             dimensionedScalar betaStar_;
00123 
00124             dimensionedScalar a1_;
00125             dimensionedScalar c1_;
00126 
00127 
00128         //- Wall distance
00129         //  Note: different to wall distance in parent RASModel
00130         wallDist y_;
00131 
00132         // Fields
00133 
00134             volScalarField k_;
00135             volScalarField omega_;
00136             volScalarField mut_;
00137             volScalarField alphat_;
00138 
00139 
00140     // Private member functions
00141 
00142         tmp<volScalarField> F1(const volScalarField& CDkOmega) const;
00143         tmp<volScalarField> F2() const;
00144 
00145         tmp<volScalarField> blend
00146         (
00147             const volScalarField& F1,
00148             const dimensionedScalar& psi1,
00149             const dimensionedScalar& psi2
00150         ) const
00151         {
00152             return F1*(psi1 - psi2) + psi2;
00153         }
00154 
00155         tmp<volScalarField> alphaK(const volScalarField& F1) const
00156         {
00157             return blend(F1, alphaK1_, alphaK2_);
00158         }
00159 
00160         tmp<volScalarField> alphaOmega(const volScalarField& F1) const
00161         {
00162             return blend(F1, alphaOmega1_, alphaOmega2_);
00163         }
00164 
00165         tmp<volScalarField> beta(const volScalarField& F1) const
00166         {
00167             return blend(F1, beta1_, beta2_);
00168         }
00169 
00170         tmp<volScalarField> gamma(const volScalarField& F1) const
00171         {
00172             return blend(F1, gamma1_, gamma2_);
00173         }
00174 
00175 
00176 public:
00177 
00178     //- Runtime type information
00179     TypeName("kOmegaSST");
00180 
00181 
00182     // Constructors
00183 
00184         //- Construct from components
00185         kOmegaSST
00186         (
00187             const volScalarField& rho,
00188             const volVectorField& U,
00189             const surfaceScalarField& phi,
00190             const basicThermo& thermophysicalModel
00191         );
00192 
00193 
00194     //- Destructor
00195     virtual ~kOmegaSST()
00196     {}
00197 
00198 
00199     // Member Functions
00200 
00201         //- Return the effective diffusivity for k
00202         tmp<volScalarField> DkEff(const volScalarField& F1) const
00203         {
00204             return tmp<volScalarField>
00205             (
00206                 new volScalarField("DkEff", alphaK(F1)*mut_ + mu())
00207             );
00208         }
00209 
00210         //- Return the effective diffusivity for omega
00211         tmp<volScalarField> DomegaEff(const volScalarField& F1) const
00212         {
00213             return tmp<volScalarField>
00214             (
00215                 new volScalarField("DomegaEff", alphaOmega(F1)*mut_ + mu())
00216             );
00217         }
00218 
00219         virtual tmp<volScalarField> mut() const
00220         {
00221             return mut_;
00222         }
00223 
00224         //- Return the effective turbulent thermal diffusivity
00225         virtual tmp<volScalarField> alphaEff() const
00226         {
00227             return tmp<volScalarField>
00228             (
00229                 new volScalarField("alphaEff", alphat_ + alpha())
00230             );
00231         }
00232 
00233         //- Return the turbulence kinetic energy
00234         virtual tmp<volScalarField> k() const
00235         {
00236             return k_;
00237         }
00238 
00239         virtual tmp<volScalarField> omega() const
00240         {
00241             return omega_;
00242         }
00243 
00244         //- Return the turbulence kinetic energy dissipation rate
00245         virtual tmp<volScalarField> epsilon() const
00246         {
00247             return tmp<volScalarField>
00248             (
00249                 new volScalarField
00250                 (
00251                     IOobject
00252                     (
00253                         "epsilon",
00254                         mesh_.time().timeName(),
00255                         mesh_
00256                     ),
00257                     betaStar_*k_*omega_,
00258                     omega_.boundaryField().types()
00259                 )
00260             );
00261         }
00262 
00263         //- Return the Reynolds stress tensor
00264         virtual tmp<volSymmTensorField> R() const;
00265 
00266         //- Return the effective stress tensor including the laminar stress
00267         virtual tmp<volSymmTensorField> devRhoReff() const;
00268 
00269         //- Return the source term for the momentum equation
00270         virtual tmp<fvVectorMatrix> divDevRhoReff(volVectorField& U) const;
00271 
00272         //- Solve the turbulence equations and correct the turbulence viscosity
00273         virtual void correct();
00274 
00275         //- Read RASProperties dictionary
00276         virtual bool read();
00277 };
00278 
00279 
00280 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00281 
00282 } // End namespace RASModels
00283 } // End namespace compressible
00284 } // End namespace Foam
00285 
00286 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00287 
00288 #endif
00289 
00290 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines