00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
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
00040
00041 defineTypeNameAndDebug(SpalartAllmaras, 0);
00042 addToRunTimeSelectionTable(LESModel, SpalartAllmaras, dictionary);
00043
00044
00045
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
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
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 }
00356 }
00357 }
00358
00359