FreeFOAM The Cross-Platform CFD Toolkit
SCOPEXiEq.C
Go to the documentation of this file.
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 "SCOPEXiEq.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 
44 (
45  const dictionary& XiEqProperties,
46  const hhuCombustionThermo& thermo,
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 
76 {}
77 
78 
79 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
80 
82 {
83  const volScalarField& k = turbulence_.k();
84  const volScalarField& epsilon = turbulence_.epsilon();
85 
86  volScalarField up = sqrt((2.0/3.0)*k);
87  volScalarField l = (lCoef*sqrt(3.0/2.0))*up*k/epsilon;
88  volScalarField Rl = up*l*thermo_.rhou()/thermo_.muu();
89 
90  volScalarField upBySu = up/(Su_ + SuMin);
91  volScalarField K = 0.157*upBySu/sqrt(Rl);
92  volScalarField Ma = MaModel.Ma();
93 
94  tmp<volScalarField> tXiEq
95  (
96  new volScalarField
97  (
98  IOobject
99  (
100  "XiEq",
101  epsilon.time().timeName(),
102  epsilon.db(),
105  ),
106  epsilon.mesh(),
107  dimensionedScalar("XiEq", dimless, 0.0)
108  )
109  );
110  volScalarField& xieq = tXiEq();
111 
112  forAll(xieq, celli)
113  {
114  if (Ma[celli] > 0.01)
115  {
116  xieq[celli] =
117  XiEqCoef*pow(K[celli]*Ma[celli], -XiEqExp)*upBySu[celli];
118  }
119  }
120 
121  forAll(xieq.boundaryField(), patchi)
122  {
123  scalarField& xieqp = xieq.boundaryField()[patchi];
124  const scalarField& Kp = K.boundaryField()[patchi];
125  const scalarField& Map = Ma.boundaryField()[patchi];
126  const scalarField& upBySup = upBySu.boundaryField()[patchi];
127 
128  forAll(xieqp, facei)
129  {
130  if (Ma[facei] > 0.01)
131  {
132  xieqp[facei] =
133  XiEqCoef*pow(Kp[facei]*Map[facei], -XiEqExp)*upBySup[facei];
134  }
135  }
136  }
137 
138  return tXiEq;
139 }
140 
141 
142 bool Foam::XiEqModels::SCOPEXiEq::read(const dictionary& XiEqProperties)
143 {
144  XiEqModel::read(XiEqProperties);
145 
146  XiEqModelCoeffs_.lookup("XiEqCoef") >> XiEqCoef;
147  XiEqModelCoeffs_.lookup("XiEqExp") >> XiEqExp;
148  XiEqModelCoeffs_.lookup("lCoef") >> lCoef;
149 
150  return true;
151 }
152 
153 
154 // ************************ vim: set sw=4 sts=4 et: ************************ //