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