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) 1991-2009 OpenCFD Ltd.
   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 the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     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, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 \*---------------------------------------------------------------------------*/
  26
  27 #include "LienLeschzinerLowRe.H"
  28 #include "wallFvPatch.H"
  29 #include "addToRunTimeSelectionTable.H"
  30
  31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  32
  33 namespace Foam
  34 {
  35 namespace incompressible
  36 {
  37 namespace RASModels
  38 {
  39
  40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  41
  42 defineTypeNameAndDebug(LienLeschzinerLowRe, 0);
  43 addToRunTimeSelectionTable(RASModel, LienLeschzinerLowRe, dictionary);
  44
  45 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  46
  47 LienLeschzinerLowRe::LienLeschzinerLowRe
  48 (
  49     const volVectorField& U,
  50     const surfaceScalarField& phi,
  51     transportModel& lamTransportModel
  52 )
  53 :
  54     RASModel(typeName, U, phi, lamTransportModel),
  55
  56     C1_
  57     (
  58         dimensioned<scalar>::lookupOrAddToDict
  59         (
  60             "C1",
  61             coeffDict_,
  62             1.44
  63         )
  64     ),
  65     C2_
  66     (
  67         dimensioned<scalar>::lookupOrAddToDict
  68         (
  69             "C2",
  70             coeffDict_,
  71             1.92
  72         )
  73     ),
  74     sigmak_
  75     (
  76         dimensioned<scalar>::lookupOrAddToDict
  77         (
  78             "sigmak",
  79             coeffDict_,
  80             1.0
  81         )
  82     ),
  83     sigmaEps_
  84     (
  85         dimensioned<scalar>::lookupOrAddToDict
  86         (
  87             "sigmaEps",
  88             coeffDict_,
  89             1.3
  90         )
  91     ),
  92     Cmu_
  93     (
  94         dimensioned<scalar>::lookupOrAddToDict
  95         (
  96             "Cmu",
  97             coeffDict_,
  98             0.09
  99         )
 100     ),
 101     kappa_
 102     (
 103         dimensioned<scalar>::lookupOrAddToDict
 104         (
 105             "kappa",
 106             coeffDict_,
 107             0.41
 108         )
 109     ),
 110     Am_
 111     (
 112         dimensioned<scalar>::lookupOrAddToDict
 113         (
 114             "Am",
 115             coeffDict_,
 116             0.016
 117         )
 118     ),
 119     Aepsilon_
 120     (
 121         dimensioned<scalar>::lookupOrAddToDict
 122         (
 123             "Aepsilon",
 124             coeffDict_,
 125             0.263
 126         )
 127     ),
 128     Amu_
 129     (
 130         dimensioned<scalar>::lookupOrAddToDict
 131         (
 132             "Amu",
 133             coeffDict_,
 134             0.00222
 135         )
 136     ),
 137
 138     k_
 139     (
 140         IOobject
 141         (
 142             "k",
 143             runTime_.timeName(),
 144             mesh_,
 145             IOobject::MUST_READ,
 146             IOobject::AUTO_WRITE
 147         ),
 148         mesh_
 149     ),
 150
 151     epsilon_
 152     (
 153         IOobject
 154         (
 155             "epsilon",
 156             runTime_.timeName(),
 157             mesh_,
 158             IOobject::MUST_READ,
 159             IOobject::AUTO_WRITE
 160         ),
 161         mesh_
 162     ),
 163
 164     y_(mesh_),
 165
 166     yStar_(sqrt(k_)*y_/nu() + SMALL),
 167
 168     nut_
 169     (
 170         Cmu_*(scalar(1) - exp(-Am_*yStar_))
 171        /(scalar(1) - exp(-Aepsilon_*yStar_) + SMALL)*sqr(k_)
 172        /(epsilon_ + epsilonSmall_)
 173     )
 174 {
 175     printCoeffs();
 176 }
 177
 178
 179 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 180
 181 tmp<volSymmTensorField> LienLeschzinerLowRe::R() const
 182 {
 183     volTensorField gradU = fvc::grad(U_);
 184
 185     return tmp<volSymmTensorField>
 186     (
 187         new volSymmTensorField
 188         (
 189             IOobject
 190             (
 191                 "R",
 192                 runTime_.timeName(),
 193                 mesh_,
 194                 IOobject::NO_READ,
 195                 IOobject::NO_WRITE
 196             ),
 197             ((2.0/3.0)*I)*k_ - nut_*twoSymm(gradU),
 198             k_.boundaryField().types()
 199         )
 200     );
 201 }
 202
 203
 204 tmp<volSymmTensorField> LienLeschzinerLowRe::devReff() const
 205 {
 206     return tmp<volSymmTensorField>
 207     (
 208         new volSymmTensorField
 209         (
 210             IOobject
 211             (
 212                 "devRhoReff",
 213                 runTime_.timeName(),
 214                 mesh_,
 215                 IOobject::NO_READ,
 216                 IOobject::NO_WRITE
 217             ),
 218            -nuEff()*dev(twoSymm(fvc::grad(U_)))
 219         )
 220     );
 221 }
 222
 223
 224 tmp<fvVectorMatrix> LienLeschzinerLowRe::divDevReff(volVectorField& U) const
 225 {
 226     return
 227     (
 228       - fvm::laplacian(nuEff(), U)
 229     //- (fvc::grad(U) & fvc::grad(nuEff()))
 230       - fvc::div(nuEff()*fvc::grad(U)().T())
 231     );
 232 }
 233
 234
 235 bool LienLeschzinerLowRe::read()
 236 {
 237     if (RASModel::read())
 238     {
 239         C1_.readIfPresent(coeffDict());
 240         C2_.readIfPresent(coeffDict());
 241         sigmak_.readIfPresent(coeffDict());
 242         sigmaEps_.readIfPresent(coeffDict());
 243         Cmu_.readIfPresent(coeffDict());
 244         Am_.readIfPresent(coeffDict());
 245         Aepsilon_.readIfPresent(coeffDict());
 246         Amu_.readIfPresent(coeffDict());
 247
 248         return true;
 249     }
 250     else
 251     {
 252         return false;
 253     }
 254 }
 255
 256
 257 void LienLeschzinerLowRe::correct()
 258 {
 259     RASModel::correct();
 260
 261     if (!turbulence_)
 262     {
 263         return;
 264     }
 265
 266     if (mesh_.changing())
 267     {
 268         y_.correct();
 269     }
 270
 271     scalar Cmu75 = pow(Cmu_.value(), 0.75);
 272
 273     volTensorField gradU = fvc::grad(U_);
 274
 275     // generation term
 276     volScalarField S2 = symm(gradU) && gradU;
 277
 278     yStar_ = sqrt(k_)*y_/nu() + SMALL;
 279     volScalarField Rt = sqr(k_)/(nu()*epsilon_);
 280
 281     volScalarField fMu =
 282         (scalar(1) - exp(-Am_*yStar_))
 283        /(scalar(1) - exp(-Aepsilon_*yStar_) + SMALL);
 284
 285     volScalarField f2 = scalar(1) - 0.3*exp(-sqr(Rt));
 286
 287     volScalarField G = Cmu_*fMu*sqr(k_)/epsilon_*S2;
 288
 289
 290     // Dissipation equation
 291     tmp<fvScalarMatrix> epsEqn
 292     (
 293         fvm::ddt(epsilon_)
 294       + fvm::div(phi_, epsilon_)
 295       - fvm::laplacian(DepsilonEff(), epsilon_)
 296       ==
 297         C1_*G*epsilon_/k_
 298         // E-term
 299         + C2_*f2*Cmu75*pow(k_, scalar(0.5))
 300         /(kappa_*y_*(scalar(1) - exp(-Aepsilon_*yStar_)))
 301        *exp(-Amu_*sqr(yStar_))*epsilon_
 302       - fvm::Sp(C2_*f2*epsilon_/k_, epsilon_)
 303     );
 304
 305     epsEqn().relax();
 306
 307 #   include "LienLeschzinerLowReSetWallDissipation.H"
 308 #   include "wallDissipationI.H"
 309
 310     solve(epsEqn);
 311     bound(epsilon_, epsilon0_);
 312
 313
 314     // Turbulent kinetic energy equation
 315
 316     tmp<fvScalarMatrix> kEqn
 317     (
 318         fvm::ddt(k_)
 319       + fvm::div(phi_, k_)
 320       - fvm::laplacian(DkEff(), k_)
 321       ==
 322         G
 323       - fvm::Sp(epsilon_/k_, k_)
 324     );
 325
 326     kEqn().relax();
 327     solve(kEqn);
 328     bound(k_, k0_);
 329
 330
 331     // Re-calculate viscosity
 332     nut_ = Cmu_*fMu*sqr(k_)/epsilon_;
 333 }
 334
 335
 336 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 337
 338 } // End namespace RASModels
 339 } // End namespace incompressible
 340 } // End namespace Foam
 341
 342 // ************************************************************************* //