src/turbulenceModels/incompressible/RAS/LienLeschzinerLowRe/LienLeschzinerLowRe.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   6      \\/     M anipulation  |
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of OpenFOAM.
  10
  11     OpenFOAM is free software: you can redistribute it and/or modify it
  12     under the terms of the GNU General Public License as published by
  13     the Free Software Foundation, either version 3 of the License, or
  14     (at your option) any later version.
  15
  16     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
  17     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
  18     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
  19     for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 \*---------------------------------------------------------------------------*/
  25
  26 #include "LienLeschzinerLowRe.H"
  27 #include "wallFvPatch.H"
  28 #include "addToRunTimeSelectionTable.H"
  29
  30 #include "backwardsCompatibilityWallFunctions.H"
  31
  32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  33
  34 namespace Foam
  35 {
  36 namespace incompressible
  37 {
  38 namespace RASModels
  39 {
  40
  41 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  42
  43 defineTypeNameAndDebug(LienLeschzinerLowRe, 0);
  44 addToRunTimeSelectionTable(RASModel, LienLeschzinerLowRe, dictionary);
  45
  46 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  47
  48 LienLeschzinerLowRe::LienLeschzinerLowRe
  49 (
  50     const volVectorField& U,
  51     const surfaceScalarField& phi,
  52     transportModel& transport,
  53     const word& turbulenceModelName,
  54     const word& modelName
  55 )
  56 :
  57     RASModel(modelName, U, phi, transport, turbulenceModelName),
  58
  59     C1_
  60     (
  61         dimensioned<scalar>::lookupOrAddToDict
  62         (
  63             "C1",
  64             coeffDict_,
  65             1.44
  66         )
  67     ),
  68     C2_
  69     (
  70         dimensioned<scalar>::lookupOrAddToDict
  71         (
  72             "C2",
  73             coeffDict_,
  74             1.92
  75         )
  76     ),
  77     sigmak_
  78     (
  79         dimensioned<scalar>::lookupOrAddToDict
  80         (
  81             "sigmak",
  82             coeffDict_,
  83             1.0
  84         )
  85     ),
  86     sigmaEps_
  87     (
  88         dimensioned<scalar>::lookupOrAddToDict
  89         (
  90             "sigmaEps",
  91             coeffDict_,
  92             1.3
  93         )
  94     ),
  95     Cmu_
  96     (
  97         dimensioned<scalar>::lookupOrAddToDict
  98         (
  99             "Cmu",
 100             coeffDict_,
 101             0.09
 102         )
 103     ),
 104     kappa_
 105     (
 106         dimensioned<scalar>::lookupOrAddToDict
 107         (
 108             "kappa",
 109             coeffDict_,
 110             0.41
 111         )
 112     ),
 113     Am_
 114     (
 115         dimensioned<scalar>::lookupOrAddToDict
 116         (
 117             "Am",
 118             coeffDict_,
 119             0.016
 120         )
 121     ),
 122     Aepsilon_
 123     (
 124         dimensioned<scalar>::lookupOrAddToDict
 125         (
 126             "Aepsilon",
 127             coeffDict_,
 128             0.263
 129         )
 130     ),
 131     Amu_
 132     (
 133         dimensioned<scalar>::lookupOrAddToDict
 134         (
 135             "Amu",
 136             coeffDict_,
 137             0.00222
 138         )
 139     ),
 140
 141     k_
 142     (
 143         IOobject
 144         (
 145             "k",
 146             runTime_.timeName(),
 147             mesh_,
 148             IOobject::MUST_READ,
 149             IOobject::AUTO_WRITE
 150         ),
 151         mesh_
 152     ),
 153
 154     epsilon_
 155     (
 156         IOobject
 157         (
 158             "epsilon",
 159             runTime_.timeName(),
 160             mesh_,
 161             IOobject::MUST_READ,
 162             IOobject::AUTO_WRITE
 163         ),
 164         mesh_
 165     ),
 166
 167     y_(mesh_),
 168
 169     yStar_(sqrt(k_)*y_/nu() + SMALL),
 170
 171     nut_
 172     (
 173         IOobject
 174         (
 175             "nut",
 176             runTime_.timeName(),
 177             mesh_,
 178             IOobject::NO_READ,
 179             IOobject::AUTO_WRITE
 180         ),
 181         autoCreateLowReNut("nut", mesh_)
 182     )
 183 {
 184     bound(k_, kMin_);
 185     bound(epsilon_, epsilonMin_);
 186
 187     nut_ = Cmu_*(scalar(1) - exp(-Am_*yStar_))
 188        /(scalar(1) - exp(-Aepsilon_*yStar_) + SMALL)*sqr(k_)
 189        /(epsilon_);
 190     nut_.correctBoundaryConditions();
 191
 192     printCoeffs();
 193 }
 194
 195
 196 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 197
 198 tmp<volSymmTensorField> LienLeschzinerLowRe::R() const
 199 {
 200     tmp<volTensorField> gradU = fvc::grad(U_);
 201
 202     return tmp<volSymmTensorField>
 203     (
 204         new volSymmTensorField
 205         (
 206             IOobject
 207             (
 208                 "R",
 209                 runTime_.timeName(),
 210                 mesh_,
 211                 IOobject::NO_READ,
 212                 IOobject::NO_WRITE
 213             ),
 214             ((2.0/3.0)*I)*k_ - nut_*twoSymm(gradU),
 215             k_.boundaryField().types()
 216         )
 217     );
 218 }
 219
 220
 221 tmp<volSymmTensorField> LienLeschzinerLowRe::devReff() const
 222 {
 223     return tmp<volSymmTensorField>
 224     (
 225         new volSymmTensorField
 226         (
 227             IOobject
 228             (
 229                 "devRhoReff",
 230                 runTime_.timeName(),
 231                 mesh_,
 232                 IOobject::NO_READ,
 233                 IOobject::NO_WRITE
 234             ),
 235            -nuEff()*dev(twoSymm(fvc::grad(U_)))
 236         )
 237     );
 238 }
 239
 240
 241 tmp<fvVectorMatrix> LienLeschzinerLowRe::divDevReff(volVectorField& U) const
 242 {
 243     return
 244     (
 245       - fvm::laplacian(nuEff(), U)
 246     //- (fvc::grad(U) & fvc::grad(nuEff()))
 247       - fvc::div(nuEff()*T(fvc::grad(U)))
 248     );
 249 }
 250
 251
 252 bool LienLeschzinerLowRe::read()
 253 {
 254     if (RASModel::read())
 255     {
 256         C1_.readIfPresent(coeffDict());
 257         C2_.readIfPresent(coeffDict());
 258         sigmak_.readIfPresent(coeffDict());
 259         sigmaEps_.readIfPresent(coeffDict());
 260         Cmu_.readIfPresent(coeffDict());
 261         kappa_.readIfPresent(coeffDict());
 262         Am_.readIfPresent(coeffDict());
 263         Aepsilon_.readIfPresent(coeffDict());
 264         Amu_.readIfPresent(coeffDict());
 265
 266         return true;
 267     }
 268     else
 269     {
 270         return false;
 271     }
 272 }
 273
 274
 275 void LienLeschzinerLowRe::correct()
 276 {
 277     RASModel::correct();
 278
 279     if (!turbulence_)
 280     {
 281         return;
 282     }
 283
 284     if (mesh_.changing())
 285     {
 286         y_.correct();
 287     }
 288
 289     scalar Cmu75 = pow(Cmu_.value(), 0.75);
 290
 291     const volTensorField gradU(fvc::grad(U_));
 292
 293     // generation term
 294     tmp<volScalarField> S2 = symm(gradU) && gradU;
 295
 296     yStar_ = sqrt(k_)*y_/nu() + SMALL;
 297     tmp<volScalarField> Rt = sqr(k_)/(nu()*epsilon_);
 298
 299     volScalarField fMu
 300     (
 301         (scalar(1) - exp(-Am_*yStar_))
 302        /(scalar(1) - exp(-Aepsilon_*yStar_) + SMALL)
 303     );
 304
 305     const volScalarField f2(scalar(1) - 0.3*exp(-sqr(Rt)));
 306
 307     volScalarField G("RASModel::G", Cmu_*fMu*sqr(k_)/epsilon_*S2);
 308
 309
 310     // Dissipation equation
 311     tmp<fvScalarMatrix> epsEqn
 312     (
 313         fvm::ddt(epsilon_)
 314       + fvm::div(phi_, epsilon_)
 315       - fvm::laplacian(DepsilonEff(), epsilon_)
 316       ==
 317         C1_*G*epsilon_/k_
 318         // E-term
 319         + C2_*f2*Cmu75*sqrt(k_)
 320         /(kappa_*y_*(scalar(1) - exp(-Aepsilon_*yStar_)))
 321        *exp(-Amu_*sqr(yStar_))*epsilon_
 322       - fvm::Sp(C2_*f2*epsilon_/k_, epsilon_)
 323     );
 324
 325     epsEqn().relax();
 326
 327     #include "LienLeschzinerLowReSetWallDissipation.H"
 328     #include "wallDissipationI.H"
 329
 330     solve(epsEqn);
 331     bound(epsilon_, epsilonMin_);
 332
 333
 334     // Turbulent kinetic energy equation
 335
 336     tmp<fvScalarMatrix> kEqn
 337     (
 338         fvm::ddt(k_)
 339       + fvm::div(phi_, k_)
 340       - fvm::laplacian(DkEff(), k_)
 341       ==
 342         G
 343       - fvm::Sp(epsilon_/k_, k_)
 344     );
 345
 346     kEqn().relax();
 347     solve(kEqn);
 348     bound(k_, kMin_);
 349
 350
 351     // Re-calculate viscosity
 352     nut_ = Cmu_*fMu*sqr(k_)/epsilon_;
 353 }
 354
 355
 356 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 357
 358 } // End namespace RASModels
 359 } // End namespace incompressible
 360 } // End namespace Foam
 361
 362 // ************************************************************************* //