XiModel.H
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 Class
25  Foam::XiModel
26 
27 Description
28  Base-class for all Xi models used by the b-Xi combustion model.
29  See Technical Report SH/RE/01R for details on the PDR modelling.
30 
31  Xi is given through an algebraic expression (\link algebraic.H \endlink),
32  by solving a transport equation (\link transport.H \endlink) or a
33  fixed value (\link fixed.H \endlink).
34 
35  See report TR/HGW/10 for details on the Weller two equations model.
36 
37  In the algebraic and transport methods \f$\Xi_{eq}\f$ is calculated in
38  similar way. In the algebraic approach, \f$\Xi_{eq}\f$ is the value used in
39  the \f$ b \f$ transport equation.
40 
41  \f$\Xi_{eq}\f$ is calculated as follows:
42 
43  \f$\Xi_{eq} = 1 + (1 + 2\Xi_{coeff}(0.5 - \dwea{b}))(\Xi^* - 1)\f$
44 
45  where:
46 
47  \f$ \dwea{b} \f$ is the regress variable.
48 
49  \f$ \Xi_{coeff} \f$ is a model constant.
50 
51  \f$ \Xi^* \f$ is the total equilibrium wrinkling combining the effects
52  of the flame instability and turbulence interaction and is given by
53 
54  \f[
55  \Xi^* = \frac {R}{R - G_\eta - G_{in}}
56  \f]
57 
58  where:
59 
60  \f$ G_\eta \f$ is the generation rate of wrinkling due to turbulence
61  interaction.
62 
63  \f$ G_{in} = \kappa \rho_{u}/\rho_{b} \f$ is the generation
64  rate due to the flame instability.
65 
66  By adding the removal rates of the two effects:
67 
68  \f[
69  R = G_\eta \frac{\Xi_{\eta_{eq}}}{\Xi_{\eta_{eq}} - 1}
70  + G_{in} \frac{\Xi_{{in}_{eq}}}{\Xi_{{in}_{eq}} - 1}
71  \f]
72 
73  where:
74 
75  \f$ R \f$ is the total removal.
76 
77  \f$ G_\eta \f$ is a model constant.
78 
79  \f$ \Xi_{\eta_{eq}} \f$ is the flame wrinkling due to turbulence.
80 
81  \f$ \Xi_{{in}_{eq}} \f$ is the equilibrium level of the flame wrinkling
82  generated by instability. It is a constant (default 2.5).
83 
84 
85 SourceFiles
86  XiModel.C
87 
88 \*---------------------------------------------------------------------------*/
89 
90 #ifndef XiModel_H
91 #define XiModel_H
92 
93 #include "IOdictionary.H"
97 #include "fvcDiv.H"
98 #include "runTimeSelectionTables.H"
99 
100 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
101 
102 namespace Foam
103 {
104 
105 /*---------------------------------------------------------------------------*\
106  Class XiModel Declaration
107 \*---------------------------------------------------------------------------*/
108 
109 class XiModel
110 {
111 
112 protected:
113 
114  // Protected data
115 
117 
123  const surfaceScalarField& phi_;
124 
125  //- Flame wrinkling field
127 
128 
129 public:
130 
131  //- Runtime type information
132  TypeName("XiModel");
133 
134 
135  // Declare run-time constructor selection table
136 
138  (
139  autoPtr,
140  XiModel,
141  dictionary,
142  (
143  const dictionary& XiProperties,
146  const volScalarField& Su,
147  const volScalarField& rho,
148  const volScalarField& b,
149  const surfaceScalarField& phi
150  ),
151  (
152  XiProperties,
153  thermo,
154  turbulence,
155  Su,
156  rho,
157  b,
158  phi
159  )
160  );
161 
162 
163  // Constructors
164 
165  //- Construct from components
166  XiModel
167  (
168  const dictionary& XiProperties,
171  const volScalarField& Su,
172  const volScalarField& rho,
173  const volScalarField& b,
174  const surfaceScalarField& phi
175  );
176 
177  //- Disallow default bitwise copy construction
178  XiModel(const XiModel&) = delete;
179 
180 
181  // Selectors
182 
183  //- Return a reference to the selected Xi model
184  static autoPtr<XiModel> New
185  (
186  const dictionary& XiProperties,
189  const volScalarField& Su,
190  const volScalarField& rho,
191  const volScalarField& b,
192  const surfaceScalarField& phi
193  );
194 
195 
196  //- Destructor
197  virtual ~XiModel();
198 
199 
200  // Member Functions
201 
202  //- Return the flame-wrinkling Xi
203  virtual const volScalarField& Xi() const
204  {
205  return Xi_;
206  }
207 
208  //- Return the flame diffusivity
209  virtual tmp<volScalarField> Db() const
210  {
211  return volScalarField::New
212  (
213  "Db",
214  turbulence_.rho()*turbulence_.nuEff()
215  );
216  }
217 
218  //- Add Xi to the multivariateSurfaceInterpolationScheme table
219  // if required
220  virtual void addXi
221  (
223  )
224  {}
225 
226  //- Correct the flame-wrinkling Xi
227  virtual void correct() = 0;
228 
229  //- Correct the flame-wrinkling Xi using the given convection scheme
230  virtual void correct(const fv::convectionScheme<scalar>&)
231  {
232  correct();
233  }
234 
235  //- Update properties from given dictionary
236  virtual bool read(const dictionary& XiProperties) = 0;
237 
238  //- Write fields related to Xi model
239  virtual void writeFields() = 0;
240 
241 
242  // Member Operators
243 
244  //- Disallow default bitwise assignment
245  void operator=(const XiModel&) = delete;
246 };
247 
248 
249 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
250 
251 } // End namespace Foam
252 
253 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
254 
255 #endif
256 
257 // ************************************************************************* //
Generic GeometricField class.
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,.
Base-class for all Xi models used by the b-Xi combustion model. See Technical Report SH/RE/01R for de...
Definition: XiModel.H:109
dictionary XiModelCoeffs_
Definition: XiModel.H:115
virtual ~XiModel()
Destructor.
Definition: XiModel.C:80
const volScalarField & b_
Definition: XiModel.H:121
static autoPtr< XiModel > New(const dictionary &XiProperties, const psiuMulticomponentThermo &thermo, const compressible::RASModel &turbulence, const volScalarField &Su, const volScalarField &rho, const volScalarField &b, const surfaceScalarField &phi)
Return a reference to the selected Xi model.
Definition: XiModelNew.C:31
void operator=(const XiModel &)=delete
Disallow default bitwise assignment.
XiModel(const dictionary &XiProperties, const psiuMulticomponentThermo &thermo, const compressible::RASModel &turbulence, const volScalarField &Su, const volScalarField &rho, const volScalarField &b, const surfaceScalarField &phi)
Construct from components.
Definition: XiModel.C:40
const surfaceScalarField & phi_
Definition: XiModel.H:122
declareRunTimeSelectionTable(autoPtr, XiModel, dictionary,(const dictionary &XiProperties, const psiuMulticomponentThermo &thermo, const compressible::RASModel &turbulence, const volScalarField &Su, const volScalarField &rho, const volScalarField &b, const surfaceScalarField &phi),(XiProperties, thermo, turbulence, Su, rho, b, phi))
virtual bool read(const dictionary &XiProperties)=0
Update properties from given dictionary.
Definition: XiModel.C:86
virtual tmp< volScalarField > Db() const
Return the flame diffusivity.
Definition: XiModel.H:208
virtual void correct()=0
Correct the flame-wrinkling Xi.
const psiuMulticomponentThermo & thermo_
Definition: XiModel.H:117
virtual void writeFields()=0
Write fields related to Xi model.
virtual const volScalarField & Xi() const
Return the flame-wrinkling Xi.
Definition: XiModel.H:202
volScalarField Xi_
Flame wrinkling field.
Definition: XiModel.H:125
const volScalarField & Su_
Definition: XiModel.H:119
const volScalarField & rho_
Definition: XiModel.H:120
TypeName("XiModel")
Runtime type information.
virtual void addXi(multivariateSurfaceInterpolationScheme< scalar >::fieldTable &)
Add Xi to the multivariateSurfaceInterpolationScheme table.
Definition: XiModel.H:220
const compressible::RASModel & turbulence_
Definition: XiModel.H:118
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base class for convection schemes.
Abstract base class for multi-variate surface interpolation schemes.
Base-class for combustion fluid thermodynamic properties based on compressibility.
A class for managing temporary objects.
Definition: tmp.H:55
Calculate the divergence of the given field.
volScalarField & b
Definition: createFields.H:25
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.
Macros to ease declaration of run-time selection tables.
autoPtr< incompressible::momentumTransportModel > turbulence(incompressible::momentumTransportModel::New(U, phi, viscosity))
fluidMulticomponentThermo & thermo
Definition: createFields.H:31