applications/solvers/combustion/PDRFoam/XiModels/XiEqModels/SCOPEXiEq/SCOPEXiEq.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
   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 "SCOPEXiEq.H"
  27 #include "addToRunTimeSelectionTable.H"
  28
  29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  30
  31 namespace Foam
  32 {
  33 namespace XiEqModels
  34 {
  35     defineTypeNameAndDebug(SCOPEXiEq, 0);
  36     addToRunTimeSelectionTable(XiEqModel, SCOPEXiEq, dictionary);
  37 }
  38 }
  39
  40
  41 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  42
  43 Foam::XiEqModels::SCOPEXiEq::SCOPEXiEq
  44 (
  45     const dictionary& XiEqProperties,
  46     const hhuCombustionThermo& thermo,
  47     const compressible::RASModel& turbulence,
  48     const volScalarField& Su
  49 )
  50 :
  51     XiEqModel(XiEqProperties, thermo, turbulence, Su),
  52     XiEqCoef_(readScalar(XiEqModelCoeffs_.lookup("XiEqCoef"))),
  53     XiEqExp_(readScalar(XiEqModelCoeffs_.lookup("XiEqExp"))),
  54     lCoef_(readScalar(XiEqModelCoeffs_.lookup("lCoef"))),
  55     SuMin_(0.01*Su.average()),
  56     MaModel
  57     (
  58         IOdictionary
  59         (
  60             IOobject
  61             (
  62                 "combustionProperties",
  63                 Su.mesh().time().constant(),
  64                 Su.mesh(),
  65                 IOobject::MUST_READ
  66             )
  67         ),
  68         thermo
  69     )
  70 {}
  71
  72
  73 // * * * * * * * * * * * * * * * * Destructors * * * * * * * * * * * * * * * //
  74
  75 Foam::XiEqModels::SCOPEXiEq::~SCOPEXiEq()
  76 {}
  77
  78
  79 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
  80
  81 Foam::tmp<Foam::volScalarField> Foam::XiEqModels::SCOPEXiEq::XiEq() const
  82 {
  83     const volScalarField& k = turbulence_.k();
  84     const volScalarField& epsilon = turbulence_.epsilon();
  85
  86     volScalarField up(sqrt((2.0/3.0)*k));
  87     if (subGridSchelkin())
  88     {
  89         up.internalField() += calculateSchelkinEffect();
  90     }
  91
  92     volScalarField l(lCoef_*sqrt(3.0/2.0)*up*k/epsilon);
  93     volScalarField Rl(up*l*thermo_.rhou()/thermo_.muu());
  94
  95     volScalarField upBySu(up/(Su_ + SuMin_));
  96     volScalarField K(0.157*upBySu/sqrt(Rl));
  97     volScalarField Ma(MaModel.Ma());
  98
  99     tmp<volScalarField> tXiEq
 100     (
 101         new volScalarField
 102         (
 103             IOobject
 104             (
 105                 "XiEq",
 106                 epsilon.time().timeName(),
 107                 epsilon.db(),
 108                 IOobject::NO_READ,
 109                 IOobject::NO_WRITE
 110             ),
 111             epsilon.mesh(),
 112             dimensionedScalar("XiEq", dimless, 0.0)
 113         )
 114     );
 115     volScalarField& xieq = tXiEq();
 116
 117     forAll(xieq, celli)
 118     {
 119         if (Ma[celli] > 0.01)
 120         {
 121             xieq[celli] =
 122                 XiEqCoef_*pow(K[celli]*Ma[celli], -XiEqExp_)*upBySu[celli];
 123         }
 124     }
 125
 126     forAll(xieq.boundaryField(), patchi)
 127     {
 128         scalarField& xieqp = xieq.boundaryField()[patchi];
 129         const scalarField& Kp = K.boundaryField()[patchi];
 130         const scalarField& Map = Ma.boundaryField()[patchi];
 131         const scalarField& upBySup = upBySu.boundaryField()[patchi];
 132
 133         forAll(xieqp, facei)
 134         {
 135             if (Ma[facei] > 0.01)
 136             {
 137                 xieqp[facei] =
 138                     XiEqCoef_*pow(Kp[facei]*Map[facei], -XiEqExp_)
 139                    *upBySup[facei];
 140             }
 141         }
 142     }
 143
 144     return tXiEq;
 145 }
 146
 147
 148 bool Foam::XiEqModels::SCOPEXiEq::read(const dictionary& XiEqProperties)
 149 {
 150     XiEqModel::read(XiEqProperties);
 151
 152     XiEqModelCoeffs_.lookup("XiEqCoef") >> XiEqCoef_;
 153     XiEqModelCoeffs_.lookup("XiEqExp") >> XiEqExp_;
 154     XiEqModelCoeffs_.lookup("lCoef") >> lCoef_;
 155
 156     return true;
 157 }
 158
 159
 160 // ************************************************************************* //