src/turbulenceModels/incompressible/RAS/LamBremhorstKE/LamBremhorstKE.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | foam-extend: Open Source CFD
   4    \\    /   O peration     | Version:     3.2
   5     \\  /    A nd           | Web:         http://www.foam-extend.org
   6      \\/     M anipulation  | For copyright notice see file Copyright
   7 -------------------------------------------------------------------------------
   8 License
   9     This file is part of foam-extend.
  10
  11     foam-extend 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 3 of the License, or (at your
  14     option) any later version.
  15
  16     foam-extend is distributed in the hope that it will be useful, but
  17     WITHOUT ANY WARRANTY; without even the implied warranty of
  18     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19     General Public License for more details.
  20
  21     You should have received a copy of the GNU General Public License
  22     along with foam-extend.  If not, see <http://www.gnu.org/licenses/>.
  23
  24 \*---------------------------------------------------------------------------*/
  25
  26 #include "LamBremhorstKE.H"
  27 #include "addToRunTimeSelectionTable.H"
  28
  29 #include "backwardsCompatibilityWallFunctions.H"
  30
  31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  32
  33 namespace Foam
  34 {
  35 namespace incompressible
  36 {
  37 namespace RASModels
  38 {
  39
  40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  41
  42 defineTypeNameAndDebug(LamBremhorstKE, 0);
  43 addToRunTimeSelectionTable(RASModel, LamBremhorstKE, dictionary);
  44
  45 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  46
  47 LamBremhorstKE::LamBremhorstKE
  48 (
  49     const volVectorField& U,
  50     const surfaceScalarField& phi,
  51     transportModel& lamTransportModel
  52 )
  53 :
  54     RASModel(typeName, U, phi, lamTransportModel),
  55
  56     Cmu_
  57     (
  58         dimensionedScalar::lookupOrAddToDict
  59         (
  60             "Cmu",
  61             coeffDict_,
  62             0.09
  63         )
  64     ),
  65     C1_
  66     (
  67         dimensionedScalar::lookupOrAddToDict
  68         (
  69             "C1",
  70             coeffDict_,
  71             1.44
  72         )
  73     ),
  74     C2_
  75     (
  76         dimensionedScalar::lookupOrAddToDict
  77         (
  78             "C2",
  79             coeffDict_,
  80             1.92
  81         )
  82     ),
  83     sigmaEps_
  84     (
  85         dimensionedScalar::lookupOrAddToDict
  86         (
  87             "alphaEps",
  88             coeffDict_,
  89             1.3
  90         )
  91     ),
  92
  93     k_
  94     (
  95         IOobject
  96         (
  97             "k",
  98             runTime_.timeName(),
  99             U_.db(),
 100             IOobject::MUST_READ,
 101             IOobject::AUTO_WRITE
 102         ),
 103         mesh_
 104     ),
 105
 106     epsilon_
 107     (
 108         IOobject
 109         (
 110             "epsilon",
 111             runTime_.timeName(),
 112             U_.db(),
 113             IOobject::MUST_READ,
 114             IOobject::AUTO_WRITE
 115         ),
 116         mesh_
 117     ),
 118
 119     y_(mesh_),
 120
 121     Rt_(sqr(k_)/(nu()*epsilon_)),
 122
 123     fMu_
 124     (
 125         sqr(scalar(1) - exp(-0.0165*(sqrt(k_)*y_/nu())))
 126        *(scalar(1) + 20.5/(Rt_ + SMALL))
 127     ),
 128
 129     nut_
 130     (
 131         IOobject
 132         (
 133             "nut",
 134             runTime_.timeName(),
 135             U_.db(),
 136             IOobject::NO_READ,
 137             IOobject::AUTO_WRITE
 138         ),
 139         autoCreateLowReNut("nut", mesh_)
 140     )
 141 {
 142     nut_ = Cmu_*fMu_*sqr(k_)/(epsilon_ + epsilonSmall_);
 143     nut_ = min(nut_, nuRatio()*nu());
 144     nut_.correctBoundaryConditions();
 145
 146     printCoeffs();
 147 }
 148
 149
 150 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 151
 152 tmp<volSymmTensorField> LamBremhorstKE::R() const
 153 {
 154     return tmp<volSymmTensorField>
 155     (
 156         new volSymmTensorField
 157         (
 158             IOobject
 159             (
 160                 "R",
 161                 runTime_.timeName(),
 162                 U_.db(),
 163                 IOobject::NO_READ,
 164                 IOobject::NO_WRITE
 165             ),
 166             ((2.0/3.0)*I)*k_ - nut_*twoSymm(fvc::grad(U_)),
 167             k_.boundaryField().types()
 168         )
 169     );
 170 }
 171
 172
 173 tmp<volSymmTensorField> LamBremhorstKE::devReff() const
 174 {
 175     return tmp<volSymmTensorField>
 176     (
 177         new volSymmTensorField
 178         (
 179             IOobject
 180             (
 181                 "devRhoReff",
 182                 runTime_.timeName(),
 183                 U_.db(),
 184                 IOobject::NO_READ,
 185                 IOobject::NO_WRITE
 186             ),
 187            -nuEff()*dev(twoSymm(fvc::grad(U_)))
 188         )
 189     );
 190 }
 191
 192
 193 tmp<fvVectorMatrix> LamBremhorstKE::divDevReff(volVectorField& U) const
 194 {
 195     return
 196     (
 197       - fvm::laplacian(nuEff(), U)
 198       - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
 199     );
 200 }
 201
 202
 203 bool LamBremhorstKE::read()
 204 {
 205     if (RASModel::read())
 206     {
 207         Cmu_.readIfPresent(coeffDict());
 208         C1_.readIfPresent(coeffDict());
 209         C2_.readIfPresent(coeffDict());
 210         sigmaEps_.readIfPresent(coeffDict());
 211
 212         return true;
 213     }
 214     else
 215     {
 216         return false;
 217     }
 218 }
 219
 220
 221 void LamBremhorstKE::correct()
 222 {
 223     // Bound in case of topological change
 224     // HJ, 22/Aug/2007
 225     if (mesh_.changing())
 226     {
 227         bound(k_, k0_);
 228         bound(epsilon_, epsilon0_);
 229     }
 230
 231     RASModel::correct();
 232
 233     if (!turbulence_)
 234     {
 235         return;
 236     }
 237
 238     if (mesh_.changing())
 239     {
 240         y_.correct();
 241     }
 242
 243     volScalarField G("RASModel::G", nut_*2*magSqr(symm(fvc::grad(U_))));
 244
 245
 246     // Calculate parameters and coefficients for low-Reynolds number model
 247
 248     Rt_ = sqr(k_)/(nu()*epsilon_);
 249     volScalarField Ry = sqrt(k_)*y_/nu();
 250
 251     fMu_ = sqr(scalar(1) - exp(-0.0165*Ry))*(scalar(1) + 20.5/(Rt_ + SMALL));
 252
 253     volScalarField f1 = scalar(1) + pow(0.05/(fMu_ + SMALL), 3);
 254     volScalarField f2 = scalar(1) - exp(-sqr(Rt_));
 255
 256
 257     // Dissipation equation
 258
 259     tmp<fvScalarMatrix> epsEqn
 260     (
 261         fvm::ddt(epsilon_)
 262       + fvm::div(phi_, epsilon_)
 263       + fvm::SuSp(-fvc::div(phi_), epsilon_)
 264       - fvm::laplacian(DepsilonEff(), epsilon_)
 265      ==
 266         C1_*f1*G*epsilon_/k_
 267       - fvm::Sp(C2_*f2*epsilon_/k_, epsilon_)
 268     );
 269
 270     epsEqn().relax();
 271     solve(epsEqn);
 272     bound(epsilon_, epsilon0_);
 273
 274
 275     // Turbulent kinetic energy equation
 276
 277     tmp<fvScalarMatrix> kEqn
 278     (
 279         fvm::ddt(k_)
 280       + fvm::div(phi_, k_)
 281       + fvm::SuSp(-fvc::div(phi_), k_)
 282       - fvm::laplacian(DkEff(), k_)
 283      ==
 284         G - fvm::Sp(epsilon_/k_, k_)
 285     );
 286
 287     kEqn().relax();
 288     solve(kEqn);
 289     bound(k_, k0_);
 290
 291
 292     // Re-calculate viscosity
 293     nut_ == Cmu_*fMu_*sqr(k_)/epsilon_;
 294     nut_ == min(nut_, nuRatio()*nu());
 295 }
 296
 297
 298 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 299
 300 } // End namespace RASModels
 301 } // End namespace incompressible
 302 } // End namespace Foam
 303
 304 // ************************************************************************* //