src/turbulenceModels/incompressible/RAS/kOmega/kOmega.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 "kOmega.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(kOmega, 0);
  44 addToRunTimeSelectionTable(RASModel, kOmega, dictionary);
  45
  46 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  47
  48 kOmega::kOmega
  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             "betaStar",
  62             coeffDict_,
  63             0.09
  64         )
  65     ),
  66     beta_
  67     (
  68         dimensioned<scalar>::lookupOrAddToDict
  69         (
  70             "beta",
  71             coeffDict_,
  72             0.072
  73         )
  74     ),
  75     alpha_
  76     (
  77         dimensioned<scalar>::lookupOrAddToDict
  78         (
  79             "alpha",
  80             coeffDict_,
  81             0.52
  82         )
  83     ),
  84     alphaK_
  85     (
  86         dimensioned<scalar>::lookupOrAddToDict
  87         (
  88             "alphaK",
  89             coeffDict_,
  90             0.5
  91         )
  92     ),
  93     alphaOmega_
  94     (
  95         dimensioned<scalar>::lookupOrAddToDict
  96         (
  97             "alphaOmega",
  98             coeffDict_,
  99             0.5
 100         )
 101     ),
 102
 103     k_
 104     (
 105         IOobject
 106         (
 107             "k",
 108             runTime_.timeName(),
 109             mesh_,
 110             IOobject::NO_READ,
 111             IOobject::AUTO_WRITE
 112         ),
 113         autoCreateK("k", mesh_)
 114     ),
 115     omega_
 116     (
 117         IOobject
 118         (
 119             "omega",
 120             runTime_.timeName(),
 121             mesh_,
 122             IOobject::NO_READ,
 123             IOobject::AUTO_WRITE
 124         ),
 125         autoCreateOmega("omega", mesh_)
 126     ),
 127     nut_
 128     (
 129         IOobject
 130         (
 131             "nut",
 132             runTime_.timeName(),
 133             mesh_,
 134             IOobject::NO_READ,
 135             IOobject::AUTO_WRITE
 136         ),
 137         autoCreateNut("nut", mesh_)
 138     )
 139 {
 140     nut_ = k_/(omega_ + omegaSmall_);
 141     nut_.correctBoundaryConditions();
 142
 143     printCoeffs();
 144 }
 145
 146
 147 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 148
 149 tmp<volSymmTensorField> kOmega::R() const
 150 {
 151     return tmp<volSymmTensorField>
 152     (
 153         new volSymmTensorField
 154         (
 155             IOobject
 156             (
 157                 "R",
 158                 runTime_.timeName(),
 159                 mesh_,
 160                 IOobject::NO_READ,
 161                 IOobject::NO_WRITE
 162             ),
 163             ((2.0/3.0)*I)*k_ - nut_*twoSymm(fvc::grad(U_)),
 164             k_.boundaryField().types()
 165         )
 166     );
 167 }
 168
 169
 170 tmp<volSymmTensorField> kOmega::devReff() const
 171 {
 172     return tmp<volSymmTensorField>
 173     (
 174         new volSymmTensorField
 175         (
 176             IOobject
 177             (
 178                 "devRhoReff",
 179                 runTime_.timeName(),
 180                 mesh_,
 181                 IOobject::NO_READ,
 182                 IOobject::NO_WRITE
 183             ),
 184            -nuEff()*dev(twoSymm(fvc::grad(U_)))
 185         )
 186     );
 187 }
 188
 189
 190 tmp<fvVectorMatrix> kOmega::divDevReff(volVectorField& U) const
 191 {
 192     return
 193     (
 194       - fvm::laplacian(nuEff(), U)
 195       - fvc::div(nuEff()*dev(fvc::grad(U)().T()))
 196     );
 197 }
 198
 199
 200 bool kOmega::read()
 201 {
 202     if (RASModel::read())
 203     {
 204         Cmu_.readIfPresent(coeffDict());
 205         beta_.readIfPresent(coeffDict());
 206         alphaK_.readIfPresent(coeffDict());
 207         alphaOmega_.readIfPresent(coeffDict());
 208
 209         return true;
 210     }
 211     else
 212     {
 213         return false;
 214     }
 215 }
 216
 217
 218 void kOmega::correct()
 219 {
 220     RASModel::correct();
 221
 222     if (!turbulence_)
 223     {
 224         return;
 225     }
 226
 227     volScalarField G("RASModel::G", nut_*2*magSqr(symm(fvc::grad(U_))));
 228
 229     // Update omega and G at the wall
 230     omega_.boundaryField().updateCoeffs();
 231
 232     // Turbulence specific dissipation rate equation
 233     tmp<fvScalarMatrix> omegaEqn
 234     (
 235         fvm::ddt(omega_)
 236       + fvm::div(phi_, omega_)
 237       - fvm::Sp(fvc::div(phi_), omega_)
 238       - fvm::laplacian(DomegaEff(), omega_)
 239      ==
 240         alpha_*G*omega_/k_
 241       - fvm::Sp(beta_*omega_, omega_)
 242     );
 243
 244     omegaEqn().relax();
 245
 246     omegaEqn().boundaryManipulate(omega_.boundaryField());
 247
 248     solve(omegaEqn);
 249     bound(omega_, omega0_);
 250
 251
 252     // Turbulent kinetic energy equation
 253     tmp<fvScalarMatrix> kEqn
 254     (
 255         fvm::ddt(k_)
 256       + fvm::div(phi_, k_)
 257       - fvm::Sp(fvc::div(phi_), k_)
 258       - fvm::laplacian(DkEff(), k_)
 259      ==
 260         G
 261       - fvm::Sp(Cmu_*omega_, k_)
 262     );
 263
 264     kEqn().relax();
 265     solve(kEqn);
 266     bound(k_, k0_);
 267
 268
 269     // Re-calculate viscosity
 270     nut_ = k_/(omega_ + omegaSmall_);
 271     nut_.correctBoundaryConditions();
 272 }
 273
 274
 275 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 276
 277 } // End namespace RASModels
 278 } // End namespace incompressible
 279 } // End namespace Foam
 280
 281 // ************************************************************************* //