src/turbulenceModels/incompressible/LES/dynOneEqEddy/dynOneEqEddy.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 "dynOneEqEddy.H"
  28 #include "addToRunTimeSelectionTable.H"
  29
  30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  31
  32 namespace Foam
  33 {
  34 namespace incompressible
  35 {
  36 namespace LESModels
  37 {
  38
  39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  40
  41 defineTypeNameAndDebug(dynOneEqEddy, 0);
  42 addToRunTimeSelectionTable(LESModel, dynOneEqEddy, dictionary);
  43
  44 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  45
  46 void dynOneEqEddy::updateSubGridScaleFields(const volSymmTensorField& D)
  47 {
  48     nuSgs_ = ck(D)*sqrt(k_)*delta();
  49     nuSgs_.correctBoundaryConditions();
  50 }
  51
  52
  53 dimensionedScalar dynOneEqEddy::ck(const volSymmTensorField& D) const
  54 {
  55     volScalarField KK = 0.5*(filter_(magSqr(U())) - magSqr(filter_(U())));
  56
  57     volSymmTensorField LL = dev(filter_(sqr(U())) - sqr(filter_(U())));
  58
  59     volSymmTensorField MM =
  60         delta()*(filter_(sqrt(k_)*D) - 2*sqrt(KK + filter_(k_))*filter_(D));
  61
  62     dimensionedScalar MMMM = average(magSqr(MM));
  63
  64     if (MMMM.value() > VSMALL)
  65     {
  66         return average(LL && MM)/MMMM;
  67     }
  68     else
  69     {
  70         return 0.0;
  71     }
  72 }
  73
  74
  75 dimensionedScalar dynOneEqEddy::ce(const volSymmTensorField& D) const
  76 {
  77     volScalarField KK = 0.5*(filter_(magSqr(U())) - magSqr(filter_(U())));
  78
  79     volScalarField mm =
  80         pow(KK + filter_(k_), 1.5)/(2*delta()) - filter_(pow(k_, 1.5))/delta();
  81
  82     volScalarField ee =
  83         2*delta()*ck(D)
  84        *(
  85            filter_(sqrt(k_)*magSqr(D))
  86          - 2*sqrt(KK + filter_(k_))*magSqr(filter_(D))
  87         );
  88
  89     dimensionedScalar mmmm = average(magSqr(mm));
  90
  91     if (mmmm.value() > VSMALL)
  92     {
  93         return average(ee*mm)/mmmm;
  94     }
  95     else
  96     {
  97         return 0.0;
  98     }
  99 }
 100
 101
 102 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 103
 104 dynOneEqEddy::dynOneEqEddy
 105 (
 106     const volVectorField& U,
 107     const surfaceScalarField& phi,
 108     transportModel& transport
 109 )
 110 :
 111     LESModel(typeName, U, phi, transport),
 112     GenEddyVisc(U, phi, transport),
 113
 114     k_
 115     (
 116         IOobject
 117         (
 118             "k",
 119             runTime_.timeName(),
 120             mesh_,
 121             IOobject::MUST_READ,
 122             IOobject::AUTO_WRITE
 123         ),
 124         mesh_
 125     ),
 126
 127     filterPtr_(LESfilter::New(U.mesh(), coeffDict())),
 128     filter_(filterPtr_())
 129 {
 130     updateSubGridScaleFields(symm(fvc::grad(U)));
 131
 132     printCoeffs();
 133 }
 134
 135
 136 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 137
 138 void dynOneEqEddy::correct(const tmp<volTensorField>& gradU)
 139 {
 140     GenEddyVisc::correct(gradU);
 141
 142     volSymmTensorField D = symm(gradU);
 143
 144     volScalarField P = 2.0*nuSgs_*magSqr(D);
 145
 146     fvScalarMatrix kEqn
 147     (
 148        fvm::ddt(k_)
 149      + fvm::div(phi(), k_)
 150      - fvm::laplacian(DkEff(), k_)
 151     ==
 152        P
 153      - fvm::Sp(ce(D)*sqrt(k_)/delta(), k_)
 154     );
 155
 156     kEqn.relax();
 157     kEqn.solve();
 158
 159     bound(k_, k0());
 160
 161     updateSubGridScaleFields(D);
 162 }
 163
 164
 165 bool dynOneEqEddy::read()
 166 {
 167     if (GenEddyVisc::read())
 168     {
 169         filter_.read(coeffDict());
 170
 171         return true;
 172     }
 173     else
 174     {
 175         return false;
 176     }
 177 }
 178
 179
 180 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 181
 182 } // End namespace LESModels
 183 } // End namespace incompressible
 184 } // End namespace Foam
 185
 186 // ************************************************************************* //