src/turbulenceModels/incompressible/RAS/LamBremhorstKE/LamBremhorstKE.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 "LamBremhorstKE.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(LamBremhorstKE, 0);
  44 addToRunTimeSelectionTable(RASModel, LamBremhorstKE, dictionary);
  45
  46 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  47
  48 LamBremhorstKE::LamBremhorstKE
  49 (
  50     const volVectorField& U,
  51     const surfaceScalarField& phi,
  52     transportModel& lamTransportModel
  53 )
  54 :
  55     RASModel(typeName, U, phi, lamTransportModel),
  56
  57     Cmu_
  58     (
  59         dimensioned<scalar>::lookupOrAddToDict
  60         (
  61             "Cmu",
  62             coeffDict_,
  63             0.09
  64         )
  65     ),
  66     C1_
  67     (
  68         dimensioned<scalar>::lookupOrAddToDict
  69         (
  70             "C1",
  71             coeffDict_,
  72             1.44
  73         )
  74     ),
  75     C2_
  76     (
  77         dimensioned<scalar>::lookupOrAddToDict
  78         (
  79             "C2",
  80             coeffDict_,
  81             1.92
  82         )
  83     ),
  84     sigmaEps_
  85     (
  86         dimensioned<scalar>::lookupOrAddToDict
  87         (
  88             "alphaEps",
  89             coeffDict_,
  90             1.3
  91         )
  92     ),
  93
  94     k_
  95     (
  96         IOobject
  97         (
  98             "k",
  99             runTime_.timeName(),
 100             mesh_,
 101             IOobject::MUST_READ,
 102             IOobject::AUTO_WRITE
 103         ),
 104         mesh_
 105     ),
 106
 107     epsilon_
 108     (
 109         IOobject
 110         (
 111             "epsilon",
 112             runTime_.timeName(),
 113             mesh_,
 114             IOobject::MUST_READ,
 115             IOobject::AUTO_WRITE
 116         ),
 117         mesh_
 118     ),
 119
 120     y_(mesh_),
 121
 122     Rt_(sqr(k_)/(nu()*epsilon_)),
 123
 124     fMu_
 125     (
 126         sqr(scalar(1) - exp(-0.0165*(sqrt(k_)*y_/nu())))
 127        *(scalar(1) + 20.5/(Rt_ + SMALL))
 128     ),
 129
 130     nut_(Cmu_*fMu_*sqr(k_)/(epsilon_ + epsilonSmall_))
 131 {
 132     printCoeffs();
 133 }
 134
 135
 136 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 137
 138 tmp<volSymmTensorField> LamBremhorstKE::R() const
 139 {
 140     return tmp<volSymmTensorField>
 141     (
 142         new volSymmTensorField
 143         (
 144             IOobject
 145             (
 146                 "R",
 147                 runTime_.timeName(),
 148                 mesh_,
 149                 IOobject::NO_READ,
 150                 IOobject::NO_WRITE
 151             ),
 152             ((2.0/3.0)*I)*k_ - nut_*twoSymm(fvc::grad(U_)),
 153             k_.boundaryField().types()
 154         )
 155     );
 156 }
 157
 158
 159 tmp<volSymmTensorField> LamBremhorstKE::devReff() const
 160 {
 161     return tmp<volSymmTensorField>
 162     (
 163         new volSymmTensorField
 164         (
 165             IOobject
 166             (
 167                 "devRhoReff",
 168                 runTime_.timeName(),
 169                 mesh_,
 170                 IOobject::NO_READ,
 171                 IOobject::NO_WRITE
 172             ),
 173            -nuEff()*dev(twoSymm(fvc::grad(U_)))
 174         )
 175     );
 176 }
 177
 178
 179 tmp<fvVectorMatrix> LamBremhorstKE::divDevReff(volVectorField& U) const
 180 {
 181     return
 182     (
 183       - fvm::laplacian(nuEff(), U)
 184       - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
 185     );
 186 }
 187
 188
 189 bool LamBremhorstKE::read()
 190 {
 191     if (RASModel::read())
 192     {
 193         Cmu_.readIfPresent(coeffDict());
 194         C1_.readIfPresent(coeffDict());
 195         C2_.readIfPresent(coeffDict());
 196         sigmaEps_.readIfPresent(coeffDict());
 197
 198         return true;
 199     }
 200     else
 201     {
 202         return false;
 203     }
 204 }
 205
 206
 207 void LamBremhorstKE::correct()
 208 {
 209     RASModel::correct();
 210
 211     if (!turbulence_)
 212     {
 213         return;
 214     }
 215
 216     if (mesh_.changing())
 217     {
 218         y_.correct();
 219     }
 220
 221     volScalarField G("RASModel::G", nut_*2*magSqr(symm(fvc::grad(U_))));
 222
 223
 224     // Calculate parameters and coefficients for low-Reynolds number model
 225
 226     Rt_ = sqr(k_)/(nu()*epsilon_);
 227     volScalarField Ry = sqrt(k_)*y_/nu();
 228
 229     fMu_ = sqr(scalar(1) - exp(-0.0165*Ry))*(scalar(1) + 20.5/(Rt_ + SMALL));
 230
 231     volScalarField f1 = scalar(1) + pow(0.05/(fMu_ + SMALL), 3);
 232     volScalarField f2 = scalar(1) - exp(-sqr(Rt_));
 233
 234
 235     // Dissipation equation
 236
 237     tmp<fvScalarMatrix> epsEqn
 238     (
 239         fvm::ddt(epsilon_)
 240       + fvm::div(phi_, epsilon_)
 241       - fvm::laplacian(DepsilonEff(), epsilon_)
 242      ==
 243         C1_*f1*G*epsilon_/k_
 244       - fvm::Sp(C2_*f2*epsilon_/k_, epsilon_)
 245     );
 246
 247     epsEqn().relax();
 248     solve(epsEqn);
 249     bound(epsilon_, epsilon0_);
 250
 251
 252     // Turbulent kinetic energy equation
 253
 254     tmp<fvScalarMatrix> kEqn
 255     (
 256         fvm::ddt(k_)
 257       + fvm::div(phi_, k_)
 258       - fvm::laplacian(DkEff(), k_)
 259      ==
 260         G - fvm::Sp(epsilon_/k_, k_)
 261     );
 262
 263     kEqn().relax();
 264     solve(kEqn);
 265     bound(k_, k0_);
 266
 267
 268     // Re-calculate viscosity
 269     nut_ == Cmu_*fMu_*sqr(k_)/epsilon_;
 270 }
 271
 272
 273 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 274
 275 } // End namespace RASModels
 276 } // End namespace incompressible
 277 } // End namespace Foam
 278
 279 // ************************************************************************* //