SCOPEXiEq.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2023 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"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace XiEqModels
34 {
37 }
38 }
39 
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
44 (
45  const dictionary& XiEqProperties,
48  const volScalarField& Su
49 )
50 :
51  XiEqModel(XiEqProperties, thermo, turbulence, Su),
52  XiEqCoef_(XiEqModelCoeffs_.lookup<scalar>("XiEqCoef")),
53  XiEqExp_(XiEqModelCoeffs_.lookup<scalar>("XiEqExp")),
54  lCoef_(XiEqModelCoeffs_.lookup<scalar>("lCoef")),
55  SuMin_(0.01*Su.average()),
56  uPrimeCoef_(XiEqModelCoeffs_.lookup<scalar>("uPrimeCoef")),
57  subGridSchelkin_
58  (
59  readBool(XiEqModelCoeffs_.lookup("subGridSchelkin"))
60  ),
61  MaModel
62  (
63  Su.mesh().lookupObject<IOdictionary>("combustionProperties"),
64  thermo
65  )
66 {}
67 
68 
69 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
70 
72 {}
73 
74 
75 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
76 
78 {
79  tmp<volScalarField> tk = turbulence_.k();
80  const volScalarField& k = tk();
81 
82  tmp<volScalarField> tepsilon = turbulence_.epsilon();
83  const volScalarField& epsilon = tepsilon();
84 
85  volScalarField up(sqrt((2.0/3.0)*k));
86  if (subGridSchelkin_)
87  {
88  up.primitiveFieldRef() += calculateSchelkinEffect(uPrimeCoef_);
89  }
90 
91  volScalarField l(lCoef_*sqrt(3.0/2.0)*up*k/epsilon);
92  volScalarField Rl(up*l*thermo_.rhou()/thermo_.muu());
93 
94  volScalarField upBySu(up/(Su_ + SuMin_));
95  volScalarField K(0.157*upBySu/sqrt(Rl));
96  volScalarField Ma(MaModel.Ma());
97 
99  (
101  (
102  "XiEq",
103  epsilon.mesh(),
105  )
106  );
107  volScalarField& xieq = tXiEq.ref();
108 
109  forAll(xieq, celli)
110  {
111  if (Ma[celli] > 0.01)
112  {
113  xieq[celli] =
114  XiEqCoef_*pow(K[celli]*Ma[celli], -XiEqExp_)*upBySu[celli];
115  }
116  }
117 
119 
120  forAll(xieq.boundaryField(), patchi)
121  {
122  scalarField& xieqp = xieqBf[patchi];
123  const scalarField& Kp = K.boundaryField()[patchi];
124  const scalarField& Map = Ma.boundaryField()[patchi];
125  const scalarField& upBySup = upBySu.boundaryField()[patchi];
126 
127  forAll(xieqp, facei)
128  {
129  if (Ma[facei] > 0.01)
130  {
131  xieqp[facei] =
132  XiEqCoef_*pow(Kp[facei]*Map[facei], -XiEqExp_)
133  *upBySup[facei];
134  }
135  }
136  }
137 
138  return tXiEq;
139 }
140 
141 
143 {
144  XiEqModel::read(XiEqProperties);
145 
146  XiEqModelCoeffs_.lookup("XiEqCoef") >> XiEqCoef_;
147  XiEqModelCoeffs_.lookup("XiEqExp") >> XiEqExp_;
148  XiEqModelCoeffs_.lookup("lCoef") >> lCoef_;
149  XiEqModelCoeffs_.lookup("uPrimeCoef") >> uPrimeCoef_;
150  XiEqModelCoeffs_.lookup("subGridSchelkin") >> subGridSchelkin_;
151 
152  return true;
153 }
154 
155 
156 // ************************************************************************* //
label k
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Generic GeometricBoundaryField class.
Generic GeometricField class.
Internal::FieldType & primitiveFieldRef()
Return a reference to the primitive field.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:57
A HashTable to objects of type <T> with a label key.
Definition: Map.H:52
Base-class for all XiEq models used by the b-XiEq combustion model. The available models are : basicX...
Definition: XiEqModel.H:58
virtual bool read(const dictionary &XiEqProperties)=0
Update properties from given dictionary.
Definition: XiEqModel.C:68
Simple SCOPEXiEq model for XiEq based on SCOPEXiEqs correlation with a linear correction function to ...
Definition: SCOPEXiEq.H:57
virtual bool read(const dictionary &XiEqProperties)
Update properties from given dictionary.
Definition: SCOPEXiEq.C:142
virtual tmp< volScalarField > XiEq() const
Return the flame-wrinkling XiEq.
Definition: SCOPEXiEq.C:77
SCOPEXiEq(const dictionary &XiEqProperties, const psiuMulticomponentThermo &thermo, const compressible::RASModel &turbulence, const volScalarField &Su)
Construct from components.
Definition: SCOPEXiEq.C:44
virtual ~SCOPEXiEq()
Destructor.
Definition: SCOPEXiEq.C:71
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Base-class for combustion fluid thermodynamic properties based on compressibility.
A class for managing temporary objects.
Definition: tmp.H:55
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:181
const scalar epsilon
label patchi
K
Definition: pEqn.H:75
addToRunTimeSelectionTable(XiEqModel, basicSubGrid, dictionary)
defineTypeNameAndDebug(basicSubGrid, 0)
RASModel< momentumTransportModel > RASModel
Typedefs for turbulence, RAS and LES models for compressible flow based on the standard laminar trans...
tmp< VolField< Type > > average(const SurfaceField< Type > &ssf)
Area-weighted average a surfaceField creating a volField.
Definition: fvcAverage.C:46
tmp< VolField< Type > > Su(const VolField< Type > &su, const VolField< Type > &vf)
Definition: fvcSup.C:44
Namespace for OpenFOAM.
bool readBool(Istream &)
Definition: boolIO.C:66
const dimensionSet dimless
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
autoPtr< incompressible::momentumTransportModel > turbulence(incompressible::momentumTransportModel::New(U, phi, viscosity))
fluidMulticomponentThermo & thermo
Definition: createFields.H:31