basicXiSubXiEq.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 "basicXiSubXiEq.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 
53  B_
54  (
55  IOobject
56  (
57  "B",
58  Su.mesh().facesInstance(),
59  Su.mesh(),
60  IOobject::MUST_READ,
61  IOobject::NO_WRITE
62  ),
63  Su.mesh()
64  ),
65 
66  XiEqModel_(XiEqModel::New(XiEqModelCoeffs_, thermo, turbulence, Su))
67 {}
68 
69 
70 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
71 
73 {}
74 
75 
76 // * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
77 
79 {
80  const fvMesh& mesh = Su_.mesh();
81  const volVectorField& U = mesh.lookupObject<volVectorField>("U");
82 
83  const volScalarField& Nv = mesh.lookupObject<volScalarField>("Nv");
84  const volSymmTensorField& nsv =
85  mesh.lookupObject<volSymmTensorField>("nsv");
86 
87  volScalarField magU(mag(U));
88  volVectorField Uhat
89  (
90  U/(mag(U) + dimensionedScalar(U.dimensions(), 1e-4))
91  );
92 
93  const scalarField Cw = pow(mesh.V(), 2.0/3.0);
94 
96  (
97  IOobject
98  (
99  "N",
100  mesh.time().constant(),
101  mesh,
104  ),
105  mesh,
106  dimensionedScalar(Nv.dimensions(), 0)
107  );
108  N.primitiveFieldRef() = Nv.primitiveField()*Cw;
109 
111  (
112  IOobject
113  (
114  "ns",
115  U.mesh().time().name(),
116  U.mesh(),
119  ),
120  U.mesh(),
122  (
123  "zero",
124  nsv.dimensions(),
125  Zero
126  )
127  );
128  ns.primitiveFieldRef() = nsv.primitiveField()*Cw;
129 
130  volScalarField n(max(N - (Uhat & ns & Uhat), scalar(1e-4)));
131  volScalarField b((Uhat & B_ & Uhat)/sqrt(n));
132  volScalarField up(sqrt((2.0/3.0)*turbulence_.k()));
133 
134  volScalarField XiSubEq
135  (
136  scalar(1)
137  + max(2.2*sqrt(b), min(0.34*magU/up*sqrt(b), scalar(1.6)))
138  * min(n, scalar(1))
139  );
140 
141  return (XiSubEq*XiEqModel_->XiEq());
142 }
143 
144 
146 {
147  XiEqModel::read(XiEqProperties);
148 
149  return XiEqModel_->read(XiEqModelCoeffs_);
150 }
151 
152 
153 // ************************************************************************* //
label n
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
Internal::FieldType & primitiveFieldRef()
Return a reference to the primitive field.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
static const word & constant()
Return constant name.
Definition: TimePaths.H:122
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
Basic sub-grid obstacle flame-wrinkling enhancement factor model. Details supplied by J Puttock 2/7/0...
virtual ~basicSubGrid()
Destructor.
basicSubGrid(const dictionary &XiEqProperties, const psiuMulticomponentThermo &thermo, const compressible::RASModel &turbulence, const volScalarField &Su)
Construct from components.
virtual bool read(const dictionary &XiEqProperties)
Update properties from given dictionary.
virtual tmp< volScalarField > XiEq() const
Return the flame-wrinkling XiEq.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:418
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
const polyMesh & mesh() const
Return reference to polyMesh.
Definition: fvMesh.H:441
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type and name.
Base-class for combustion fluid thermodynamic properties based on compressibility.
A class for managing temporary objects.
Definition: tmp.H:55
volScalarField & b
Definition: createFields.H:25
U
Definition: pEqn.H:72
addToRunTimeSelectionTable(XiEqModel, basicSubGrid, dictionary)
defineTypeNameAndDebug(basicSubGrid, 0)
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
RASModel< momentumTransportModel > RASModel
Typedefs for turbulence, RAS and LES models for compressible flow based on the standard laminar trans...
tmp< VolField< Type > > Su(const VolField< Type > &su, const VolField< Type > &vf)
Definition: fvcSup.C:44
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
const doubleScalar e
Definition: doubleScalar.H:106
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
autoPtr< incompressible::momentumTransportModel > turbulence(incompressible::momentumTransportModel::New(U, phi, viscosity))
fluidMulticomponentThermo & thermo
Definition: createFields.H:31