Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
applications
solvers
combustion
PDRFoam
XiModels
XiModel
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 | Copyright (C) 1991-2010 OpenCFD Ltd.
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 inestability 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 inestability.
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 inestability. 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 <
OpenFOAM/IOdictionary.H
>
94
#include <
reactionThermophysicalModels/hhuCombustionThermo.H
>
95
#include <
compressibleRASModels/RASModel.H
>
96
#include <
finiteVolume/multivariateSurfaceInterpolationScheme.H
>
97
#include <
OpenFOAM/runTimeSelectionTables.H
>
98
99
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
100
101
namespace
Foam
102
{
103
104
/*---------------------------------------------------------------------------*\
105
Class XiModel Declaration
106
\*---------------------------------------------------------------------------*/
107
108
class
XiModel
109
{
110
111
protected
:
112
113
// Protected data
114
115
dictionary
XiModelCoeffs_
;
116
117
const
hhuCombustionThermo
&
thermo_
;
118
const
compressible::RASModel
&
turbulence_
;
119
const
volScalarField
&
Su_
;
120
const
volScalarField
&
rho_
;
121
const
volScalarField
&
b_
;
122
const
surfaceScalarField
&
phi_
;
123
124
//- Flame wrinking field
125
volScalarField
Xi_
;
126
127
128
private
:
129
130
// Private Member Functions
131
132
//- Disallow copy construct
133
XiModel
(
const
XiModel
&);
134
135
//- Disallow default bitwise assignment
136
void
operator=(
const
XiModel
&);
137
138
139
public
:
140
141
//- Runtime type information
142
TypeName
(
"XiModel"
);
143
144
145
// Declare run-time constructor selection table
146
147
declareRunTimeSelectionTable
148
(
149
autoPtr
,
150
XiModel
,
151
dictionary
,
152
(
153
const
dictionary
& XiProperties,
154
const
hhuCombustionThermo
&
thermo
,
155
const
compressible::RASModel
&
turbulence
,
156
const
volScalarField
&
Su
,
157
const
volScalarField
&
rho
,
158
const
volScalarField
&
b
,
159
const
surfaceScalarField
&
phi
160
),
161
(
162
XiProperties,
163
thermo,
164
turbulence,
165
Su,
166
rho,
167
b,
168
phi
169
)
170
);
171
172
173
// Selectors
174
175
//- Return a reference to the selected Xi model
176
static
autoPtr<XiModel>
New
177
(
178
const
dictionary
& XiProperties,
179
const
hhuCombustionThermo
& thermo,
180
const
compressible::RASModel
& turbulence,
181
const
volScalarField
& Su,
182
const
volScalarField
& rho,
183
const
volScalarField
& b,
184
const
surfaceScalarField
& phi
185
);
186
187
188
// Constructors
189
190
//- Construct from components
191
XiModel
192
(
193
const
dictionary
& XiProperties,
194
const
hhuCombustionThermo
& thermo,
195
const
compressible::RASModel
& turbulence,
196
const
volScalarField
& Su,
197
const
volScalarField
& rho,
198
const
volScalarField
& b,
199
const
surfaceScalarField
& phi
200
);
201
202
203
// Destructor
204
205
virtual
~XiModel
();
206
207
208
// Member Functions
209
210
//- Return the flame-wrinking Xi
211
virtual
const
volScalarField
&
Xi
()
const
212
{
213
return
Xi_
;
214
}
215
216
//- Return the flame diffusivity
217
virtual
tmp<volScalarField>
Db
()
const
218
{
219
return
turbulence_
.muEff();
220
}
221
222
//- Add Xi to the multivariateSurfaceInterpolationScheme table
223
// if required
224
virtual
void
addXi
225
(
226
multivariateSurfaceInterpolationScheme<scalar>::fieldTable
&
227
)
228
{}
229
230
//- Correct the flame-wrinking Xi
231
virtual
void
correct
() = 0;
232
233
//- Correct the flame-wrinking Xi using the given convection scheme
234
virtual
void
correct
(
const
fv::convectionScheme<scalar>
&)
235
{
236
correct
();
237
}
238
239
//- Update properties from given dictionary
240
virtual
bool
read
(
const
dictionary
& XiProperties) = 0;
241
};
242
243
244
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
245
246
}
// End namespace Foam
247
248
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
249
250
#endif
251
252
// ************************ vim: set sw=4 sts=4 et: ************************ //