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