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