OpenFOAM: API Guide
v2112
The open source CFD toolbox
PDRFoam.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-2017 OpenFOAM Foundation
9
-------------------------------------------------------------------------------
10
License
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
Application
27
PDRFoam
28
29
Group
30
grpCombustionSolvers
31
32
Description
33
Solver for compressible premixed/partially-premixed combustion with
34
turbulence modelling.
35
36
Combusting RANS code using the b-Xi two-equation model.
37
Xi may be obtained by either the solution of the Xi transport
38
equation or from an algebraic expression. Both approaches are
39
based on Gulder's flame speed correlation which has been shown
40
to be appropriate by comparison with the results from the
41
spectral model.
42
43
Strain effects are incorporated directly into the Xi equation
44
but not in the algebraic approximation. Further work need to be
45
done on this issue, particularly regarding the enhanced removal rate
46
caused by flame compression. Analysis using results of the spectral
47
model will be required.
48
49
For cases involving very lean Propane flames or other flames which are
50
very strain-sensitive, a transport equation for the laminar flame
51
speed is present. This equation is derived using heuristic arguments
52
involving the strain time scale and the strain-rate at extinction.
53
the transport velocity is the same as that for the Xi equation.
54
55
For large flames e.g. explosions additional modelling for the flame
56
wrinkling due to surface instabilities may be applied.
57
58
PDR (porosity/distributed resistance) modelling is included to handle
59
regions containing blockages which cannot be resolved by the mesh.
60
61
The fields used by this solver are:
62
\plaintable
63
betav | Volume porosity
64
Lobs | Average diameter of obstacle in cell (m)
65
Aw | Obstacle surface area per unit volume (1/m)
66
CR | Drag tensor (1/m)
67
CT | Turbulence generation parameter (1/m)
68
Nv | Number of obstacles in cell per unit volume (m^-2)
69
nsv | Tensor whose diagonal indicates the number to subtract from
70
| Nv to get the number of obstacles crossing the flow in each
71
| direction.
72
\endplaintable
73
74
\*---------------------------------------------------------------------------*/
75
76
#include "
fvCFD.H
"
77
#include "
psiuReactionThermo.H
"
78
#include "
turbulentFluidThermoModel.H
"
79
#include "
laminarFlameSpeed.H
"
80
#include "
XiModel.H
"
81
#include "
PDRDragModel.H
"
82
#include "
ignition.H
"
83
#include "
bound.H
"
84
#include "
pimpleControl.H
"
85
#include "
fvOptions.H
"
86
87
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
88
89
int
main(
int
argc,
char
*argv[])
90
{
91
argList::addNote
92
(
93
"Solver for compressible premixed/partially-premixed combustion with"
94
" turbulence modelling."
95
);
96
97
#include "
postProcess.H
"
98
99
#include "
addCheckCaseOptions.H
"
100
#include "
setRootCaseLists.H
"
101
#include "
createTime.H
"
102
#include "
createMesh.H
"
103
#include "createControl.H"
104
#include "readCombustionProperties.H"
105
#include "readGravitationalAcceleration.H"
106
#include "createFields.H"
107
#include "createFieldRefs.H"
108
#include "initContinuityErrs.H"
109
#include "
createTimeControls.H
"
110
#include "compressibleCourantNo.H"
111
#include "setInitialDeltaT.H"
112
113
turbulence
->validate();
114
scalar StCoNum = 0.0;
115
116
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
117
118
Info
<<
"\nStarting time loop\n"
<<
endl
;
119
120
while
(
runTime
.run())
121
{
122
#include "
readTimeControls.H
"
123
#include "compressibleCourantNo.H"
124
#include "setDeltaT.H"
125
126
++
runTime
;
127
Info
<<
"\n\nTime = "
<<
runTime
.timeName() <<
endl
;
128
129
#include "rhoEqn.H"
130
131
// --- Pressure-velocity PIMPLE corrector loop
132
while
(
pimple
.loop())
133
{
134
#include "UEqn.H"
135
136
// --- Pressure corrector loop
137
while
(
pimple
.correct())
138
{
139
#include "bEqn.H"
140
#include "ftEqn.H"
141
#include "EauEqn.H"
142
#include "EaEqn.H"
143
144
if
(!ign.ignited())
145
{
146
thermo
.heu() ==
thermo
.he();
147
}
148
149
#include "pEqn.H"
150
}
151
152
if
(
pimple
.turbCorr())
153
{
154
turbulence
->correct();
155
}
156
}
157
158
runTime
.write();
159
160
runTime
.printExecutionTime(Info);
161
}
162
163
Info
<<
"End\n"
;
164
165
return
0;
166
}
167
168
169
// ************************************************************************* //
PDRDragModel.H
XiModel.H
addCheckCaseOptions.H
Required Classes.
bound.H
Bound the given scalar field if it has gone unbounded.
pimple
pimpleControl & pimple
Definition:
setRegionFluidFields.H:56
thermo
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
createMesh.H
Required Variables.
createTimeControls.H
Read the control parameters used by setDeltaT.
createTime.H
turbulence
compressible::turbulenceModel & turbulence
Definition:
setRegionFluidFields.H:30
fvCFD.H
fvOptions.H
ignition.H
laminarFlameSpeed.H
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition:
Ostream.H:372
pimpleControl.H
postProcess.H
Execute application functionObjects to post-process existing results.
psiuReactionThermo.H
readTimeControls.H
Read the control parameters used by setDeltaT.
setRootCaseLists.H
turbulentFluidThermoModel.H
applications
solvers
combustion
PDRFoam
PDRFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.