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

SpalartAllmaras.C

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 \*---------------------------------------------------------------------------*/
00025 
00026 #include "SpalartAllmaras.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 
00029 #include <compressibleRASModels/backwardsCompatibilityWallFunctions.H>
00030 
00031 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00032 
00033 namespace Foam
00034 {
00035 namespace compressible
00036 {
00037 namespace RASModels
00038 {
00039 
00040 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00041 
00042 defineTypeNameAndDebug(SpalartAllmaras, 0);
00043 addToRunTimeSelectionTable(RASModel, SpalartAllmaras, dictionary);
00044 
00045 // * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * * //
00046 
00047 tmp<volScalarField> SpalartAllmaras::chi() const
00048 {
00049     return rho_*nuTilda_/mu();
00050 }
00051 
00052 
00053 tmp<volScalarField> SpalartAllmaras::fv1(const volScalarField& chi) const
00054 {
00055     volScalarField chi3 = pow3(chi);
00056     return chi3/(chi3 + pow3(Cv1_));
00057 }
00058 
00059 
00060 tmp<volScalarField> SpalartAllmaras::fv2
00061 (
00062     const volScalarField& chi,
00063     const volScalarField& fv1
00064 ) const
00065 {
00066     return 1.0/pow3(scalar(1) + chi/Cv2_);
00067 }
00068 
00069 
00070 tmp<volScalarField> SpalartAllmaras::fv3
00071 (
00072     const volScalarField& chi,
00073     const volScalarField& fv1
00074 ) const
00075 {
00076     volScalarField chiByCv2 = (1/Cv2_)*chi;
00077 
00078     return
00079         (scalar(1) + chi*fv1)
00080        *(1/Cv2_)
00081        *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
00082        /pow3(scalar(1) + chiByCv2);
00083 }
00084 
00085 
00086 tmp<volScalarField> SpalartAllmaras::fw(const volScalarField& Stilda) const
00087 {
00088     volScalarField r = min
00089     (
00090         nuTilda_
00091        /(
00092            max(Stilda, dimensionedScalar("SMALL", Stilda.dimensions(), SMALL))
00093            *sqr(kappa_*d_)
00094         ),
00095         scalar(10.0)
00096     );
00097     r.boundaryField() == 0.0;
00098 
00099     volScalarField g = r + Cw2_*(pow6(r) - r);
00100 
00101     return g*pow((1.0 + pow6(Cw3_))/(pow6(g) + pow6(Cw3_)), 1.0/6.0);
00102 }
00103 
00104 
00105 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00106 
00107 SpalartAllmaras::SpalartAllmaras
00108 (
00109     const volScalarField& rho,
00110     const volVectorField& U,
00111     const surfaceScalarField& phi,
00112     const basicThermo& thermophysicalModel
00113 )
00114 :
00115     RASModel(typeName, rho, U, phi, thermophysicalModel),
00116 
00117     sigmaNut_
00118     (
00119         dimensioned<scalar>::lookupOrAddToDict
00120         (
00121             "sigmaNut",
00122             coeffDict_,
00123             0.66666
00124         )
00125     ),
00126     kappa_
00127     (
00128         dimensioned<scalar>::lookupOrAddToDict
00129         (
00130             "kappa",
00131             coeffDict_,
00132             0.41
00133         )
00134     ),
00135     Prt_
00136     (
00137         dimensioned<scalar>::lookupOrAddToDict
00138         (
00139             "Prt",
00140             coeffDict_,
00141             1.0
00142         )
00143     ),
00144 
00145     Cb1_
00146     (
00147         dimensioned<scalar>::lookupOrAddToDict
00148         (
00149             "Cb1",
00150             coeffDict_,
00151             0.1355
00152         )
00153     ),
00154     Cb2_
00155     (
00156         dimensioned<scalar>::lookupOrAddToDict
00157         (
00158             "Cb2",
00159             coeffDict_,
00160             0.622
00161         )
00162     ),
00163     Cw1_(Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_),
00164     Cw2_
00165     (
00166         dimensioned<scalar>::lookupOrAddToDict
00167         (
00168             "Cw2",
00169             coeffDict_,
00170             0.3
00171         )
00172     ),
00173     Cw3_
00174     (
00175         dimensioned<scalar>::lookupOrAddToDict
00176         (
00177             "Cw3",
00178             coeffDict_,
00179             2.0
00180         )
00181     ),
00182     Cv1_
00183     (
00184         dimensioned<scalar>::lookupOrAddToDict
00185         (
00186             "Cv1",
00187             coeffDict_,
00188             7.1
00189         )
00190     ),
00191     Cv2_
00192     (
00193         dimensioned<scalar>::lookupOrAddToDict
00194         (
00195             "Cv2",
00196             coeffDict_,
00197             5.0
00198         )
00199     ),
00200 
00201     nuTilda_
00202     (
00203         IOobject
00204         (
00205             "nuTilda",
00206             runTime_.timeName(),
00207             mesh_,
00208             IOobject::MUST_READ,
00209             IOobject::AUTO_WRITE
00210         ),
00211         mesh_
00212     ),
00213 
00214     mut_
00215     (
00216         IOobject
00217         (
00218             "mut",
00219             runTime_.timeName(),
00220             mesh_,
00221             IOobject::MUST_READ,
00222             IOobject::AUTO_WRITE
00223         ),
00224         mesh_
00225     ),
00226 
00227     alphat_
00228     (
00229         IOobject
00230         (
00231             "alphat",
00232             runTime_.timeName(),
00233             mesh_,
00234             IOobject::NO_READ,
00235             IOobject::AUTO_WRITE
00236         ),
00237         autoCreateAlphat("alphat", mesh_)
00238     ),
00239 
00240     d_(mesh_)
00241 {
00242     alphat_ = mut_/Prt_;
00243     alphat_.correctBoundaryConditions();
00244 
00245     printCoeffs();
00246 }
00247 
00248 
00249 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00250 
00251 tmp<volSymmTensorField> SpalartAllmaras::R() const
00252 {
00253     return tmp<volSymmTensorField>
00254     (
00255         new volSymmTensorField
00256         (
00257             IOobject
00258             (
00259                 "R",
00260                 runTime_.timeName(),
00261                 mesh_,
00262                 IOobject::NO_READ,
00263                 IOobject::NO_WRITE
00264             ),
00265             ((2.0/3.0)*I)*k() - (mut_/rho_)*dev(twoSymm(fvc::grad(U_)))
00266         )
00267     );
00268 }
00269 
00270 
00271 tmp<volSymmTensorField> SpalartAllmaras::devRhoReff() const
00272 {
00273     return tmp<volSymmTensorField>
00274     (
00275         new volSymmTensorField
00276         (
00277             IOobject
00278             (
00279                 "devRhoReff",
00280                 runTime_.timeName(),
00281                 mesh_,
00282                 IOobject::NO_READ,
00283                 IOobject::NO_WRITE
00284             ),
00285            -muEff()*dev(twoSymm(fvc::grad(U_)))
00286         )
00287     );
00288 }
00289 
00290 
00291 tmp<fvVectorMatrix> SpalartAllmaras::divDevRhoReff(volVectorField& U) const
00292 {
00293     volScalarField muEff_ = muEff();
00294 
00295     return
00296     (
00297       - fvm::laplacian(muEff_, U)
00298       - fvc::div(muEff_*dev2(fvc::grad(U)().T()))
00299     );
00300 }
00301 
00302 
00303 bool SpalartAllmaras::read()
00304 {
00305     if (RASModel::read())
00306     {
00307         sigmaNut_.readIfPresent(coeffDict());
00308         kappa_.readIfPresent(coeffDict());
00309         Prt_.readIfPresent(coeffDict());
00310 
00311         Cb1_.readIfPresent(coeffDict());
00312         Cb2_.readIfPresent(coeffDict());
00313         Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
00314         Cw2_.readIfPresent(coeffDict());
00315         Cw3_.readIfPresent(coeffDict());
00316         Cv1_.readIfPresent(coeffDict());
00317         Cv2_.readIfPresent(coeffDict());
00318 
00319         return true;
00320     }
00321     else
00322     {
00323         return false;
00324     }
00325 }
00326 
00327 
00328 void SpalartAllmaras::correct()
00329 {
00330     if (!turbulence_)
00331     {
00332         // Re-calculate viscosity
00333         mut_ = rho_*nuTilda_*fv1(chi());
00334         mut_.correctBoundaryConditions();
00335 
00336         // Re-calculate thermal diffusivity
00337         alphat_ = mut_/Prt_;
00338         alphat_.correctBoundaryConditions();
00339 
00340         return;
00341     }
00342 
00343     RASModel::correct();
00344 
00345     if (mesh_.changing())
00346     {
00347         d_.correct();
00348     }
00349 
00350     volScalarField chi = this->chi();
00351     volScalarField fv1 = this->fv1(chi);
00352 
00353     volScalarField Stilda =
00354         fv3(chi, fv1)*::sqrt(2.0)*mag(skew(fvc::grad(U_)))
00355       + fv2(chi, fv1)*nuTilda_/sqr(kappa_*d_);
00356 
00357     tmp<fvScalarMatrix> nuTildaEqn
00358     (
00359         fvm::ddt(rho_, nuTilda_)
00360       + fvm::div(phi_, nuTilda_)
00361       - fvm::laplacian(DnuTildaEff(), nuTilda_)
00362       - Cb2_/sigmaNut_*rho_*magSqr(fvc::grad(nuTilda_))
00363      ==
00364         Cb1_*rho_*Stilda*nuTilda_
00365       - fvm::Sp(Cw1_*fw(Stilda)*nuTilda_*rho_/sqr(d_), nuTilda_)
00366     );
00367 
00368     nuTildaEqn().relax();
00369     solve(nuTildaEqn);
00370     bound(nuTilda_, dimensionedScalar("0", nuTilda_.dimensions(), 0.0));
00371     nuTilda_.correctBoundaryConditions();
00372 
00373     // Re-calculate viscosity
00374     mut_.internalField() = fv1*nuTilda_.internalField()*rho_.internalField();
00375     mut_.correctBoundaryConditions();
00376 
00377     // Re-calculate thermal diffusivity
00378     alphat_ = mut_/Prt_;
00379     alphat_.correctBoundaryConditions();
00380 }
00381 
00382 
00383 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00384 
00385 } // End namespace RASModels
00386 } // End namespace compressible
00387 } // End namespace Foam
00388 
00389 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines