OpenFOAM: API Guide
v2112
The open source CFD toolbox
XiEngineFoam.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
Copyright (C) 2020 OpenCFD Ltd.
10
-------------------------------------------------------------------------------
11
License
12
This file is part of OpenFOAM.
13
14
OpenFOAM is free software: you can redistribute it and/or modify it
15
under the terms of the GNU General Public License as published by
16
the Free Software Foundation, either version 3 of the License, or
17
(at your option) any later version.
18
19
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22
for more details.
23
24
You should have received a copy of the GNU General Public License
25
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27
Application
28
XiEngineFoam
29
30
Description
31
Solver for internal combustion engines.
32
33
Combusting RANS code using the b-Xi two-equation model.
34
Xi may be obtained by either the solution of the Xi transport
35
equation or from an algebraic expression. Both approaches are
36
based on Gulder's flame speed correlation which has been shown
37
to be appropriate by comparison with the results from the
38
spectral model.
39
40
Strain effects are encorporated directly into the Xi equation
41
but not in the algebraic approximation. Further work need to be
42
done on this issue, particularly regarding the enhanced removal rate
43
caused by flame compression. Analysis using results of the spectral
44
model will be required.
45
46
For cases involving very lean Propane flames or other flames which are
47
very strain-sensitive, a transport equation for the laminar flame
48
speed is present. This equation is derived using heuristic arguments
49
involving the strain time scale and the strain-rate at extinction.
50
the transport velocity is the same as that for the Xi equation.
51
52
\*---------------------------------------------------------------------------*/
53
54
#include "
fvCFD.H
"
55
#include "
engineTime.H
"
56
#include "
engineMesh.H
"
57
#include "
psiuReactionThermo.H
"
58
#include "
turbulentFluidThermoModel.H
"
59
#include "
laminarFlameSpeed.H
"
60
#include "
ignition.H
"
61
#include "
Switch.H
"
62
#include "
OFstream.H
"
63
#include "
mathematicalConstants.H
"
64
#include "
pimpleControl.H
"
65
#include "
fvOptions.H
"
66
67
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
68
69
int
main(
int
argc,
char
*argv[])
70
{
71
argList::addNote
72
(
73
"Solver for compressible premixed/partially-premixed combustion with"
74
" turbulence modelling in internal combustion engines."
75
);
76
77
#define CREATE_TIME createEngineTime.H
78
#define CREATE_MESH createEngineMesh.H
79
#include "
postProcess.H
"
80
81
#include "
setRootCaseLists.H
"
82
#include "
createEngineTime.H
"
83
#include "
createEngineMesh.H
"
84
#include "createControl.H"
85
#include "readCombustionProperties.H"
86
#include "createFields.H"
87
#include "createFieldRefs.H"
88
#include "
createRhoUf.H
"
89
#include "initContinuityErrs.H"
90
#include "
readEngineTimeControls.H
"
91
#include "compressibleCourantNo.H"
92
#include "setInitialDeltaT.H"
93
#include "startSummary.H"
94
95
turbulence
->validate();
96
97
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
98
99
Info
<<
"\nStarting time loop\n"
<<
endl
;
100
101
while
(
runTime
.run())
102
{
103
#include "
readEngineTimeControls.H
"
104
#include "compressibleCourantNo.H"
105
#include "setDeltaT.H"
106
107
++
runTime
;
108
109
Info
<<
"Crank angle = "
<<
runTime
.theta() <<
" CA-deg"
<<
endl
;
110
111
mesh
.move();
112
113
#include "rhoEqn.H"
114
115
// --- Pressure-velocity PIMPLE corrector loop
116
while
(
pimple
.loop())
117
{
118
#include "UEqn.H"
119
120
#include "ftEqn.H"
121
#include "bEqn.H"
122
#include "EauEqn.H"
123
#include "EaEqn.H"
124
125
if
(!ign.ignited())
126
{
127
thermo
.heu() ==
thermo
.he();
128
}
129
130
// --- Pressure corrector loop
131
while
(
pimple
.correct())
132
{
133
#include "pEqn.H"
134
}
135
136
if
(
pimple
.turbCorr())
137
{
138
turbulence
->correct();
139
}
140
}
141
142
#include "logSummary.H"
143
144
rho
=
thermo
.rho();
145
146
runTime
.write();
147
148
runTime
.printExecutionTime(Info);
149
}
150
151
Info
<<
"End\n"
<<
endl
;
152
153
return
0;
154
}
155
156
157
// ************************************************************************* //
OFstream.H
Switch.H
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...
mesh
dynamicFvMesh & mesh
Definition:
createDynamicFvMesh.H:6
createEngineMesh.H
createEngineTime.H
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
createRhoUf.H
Creates and initialises the velocity velocity field rhoUf.
engineMesh.H
engineTime.H
turbulence
compressible::turbulenceModel & turbulence
Definition:
setRegionFluidFields.H:30
fvCFD.H
fvOptions.H
ignition.H
laminarFlameSpeed.H
mathematicalConstants.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
readEngineTimeControls.H
rho
rho
Definition:
readInitialConditions.H:88
setRootCaseLists.H
turbulentFluidThermoModel.H
applications
solvers
combustion
XiFoam
XiEngineFoam
XiEngineFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.