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 | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2016 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
28#include "SCOPEXiEq.H"
30
31// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32
33namespace Foam
34{
35namespace XiEqModels
36{
37 defineTypeNameAndDebug(SCOPEXiEq, 0);
38 addToRunTimeSelectionTable(XiEqModel, SCOPEXiEq, dictionary);
39}
40}
41
42
43// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44
46(
47 const dictionary& XiEqProperties,
48 const psiuReactionThermo& thermo,
49 const compressible::RASModel& turbulence,
50 const volScalarField& Su
51)
52:
53 XiEqModel(XiEqProperties, thermo, turbulence, Su),
54 XiEqCoef_(XiEqModelCoeffs_.get<scalar>("XiEqCoef")),
55 XiEqExp_(XiEqModelCoeffs_.get<scalar>("XiEqExp")),
56 lCoef_(XiEqModelCoeffs_.get<scalar>("lCoef")),
57 SuMin_(0.01*Su.average()),
58 uPrimeCoef_(XiEqModelCoeffs_.get<scalar>("uPrimeCoef")),
59 subGridSchelkin_(XiEqModelCoeffs_.get<bool>("subGridSchelkin")),
60 MaModel
61 (
62 Su.mesh().lookupObject<IOdictionary>("combustionProperties"),
63 thermo
64 )
65{}
66
67
68// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
69
71{}
72
73
74// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
75
77{
78 const volScalarField& k = turbulence_.k();
79 const volScalarField& epsilon = turbulence_.epsilon();
80
81 volScalarField up(sqrt((2.0/3.0)*k));
82 if (subGridSchelkin_)
83 {
84 up.primitiveFieldRef() += calculateSchelkinEffect(uPrimeCoef_);
85 }
86
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 (
97 (
98 IOobject
99 (
100 "XiEq",
101 epsilon.time().timeName(),
102 epsilon.db(),
105 ),
106 epsilon.mesh(),
108 )
109 );
110 volScalarField& xieq = tXiEq.ref();
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 volScalarField::Boundary& xieqBf = xieq.boundaryFieldRef();
122
123 forAll(xieq.boundaryField(), patchi)
124 {
125 scalarField& xieqp = xieqBf[patchi];
126 const scalarField& Kp = K.boundaryField()[patchi];
127 const scalarField& Map = Ma.boundaryField()[patchi];
128 const scalarField& upBySup = upBySu.boundaryField()[patchi];
129
130 forAll(xieqp, facei)
131 {
132 if (Ma[facei] > 0.01)
133 {
134 xieqp[facei] =
135 XiEqCoef_*pow(Kp[facei]*Map[facei], -XiEqExp_)
136 *upBySup[facei];
137 }
138 }
139 }
140
141 return tXiEq;
142}
143
144
145bool Foam::XiEqModels::SCOPEXiEq::read(const dictionary& XiEqProperties)
146{
147 XiEqModel::read(XiEqProperties);
148
149 XiEqModelCoeffs_.readEntry("XiEqCoef", XiEqCoef_);
150 XiEqModelCoeffs_.readEntry("XiEqExp", XiEqExp_);
151 XiEqModelCoeffs_.readEntry("lCoef", lCoef_);
152 XiEqModelCoeffs_.readEntry("uPrimeCoef", uPrimeCoef_);
153 XiEqModelCoeffs_.readEntry("subGridSchelkin", subGridSchelkin_);
154
155 return true;
156}
157
158
159// ************************************************************************* //
CGAL::Exact_predicates_exact_constructions_kernel K
label k
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
GeometricBoundaryField< scalar, fvPatchField, volMesh > Boundary
Type of boundary fields.
virtual bool read()
Re-read model coefficients if they have changed.
Simple SCOPEXiEq model for XiEq based on SCOPEXiEqs correlation with a linear correction function to ...
Definition: SCOPEXiEq.H:59
virtual ~SCOPEXiEq()
Destructor.
virtual tmp< volScalarField > XiEq() const
Return the flame-wrinkling XiEq.
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
A class for managing temporary objects.
Definition: tmp.H:65
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
bool
Definition: EEqn.H:20
dynamicFvMesh & mesh
scalar epsilon
compressible::turbulenceModel & turbulence
zeroField Su
Definition: alphaSuSp.H:1
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
tmp< GeometricField< Type, faPatchField, areaMesh > > average(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
Area-weighted average a edgeField creating a areaField.
Definition: facAverage.C:47
Namespace for OpenFOAM.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimless
Dimensionless.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333