XiDyMFoam.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-2015 OpenFOAM Foundation
9 Copyright (C) 2015-2016 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
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
27Application
28 XiFoam
29
30Group
31 grpCombustionSolvers grpMovingMeshSolvers
32
33Description
34 Solver for compressible premixed/partially-premixed combustion with
35 turbulence modelling.
36
37 Combusting RANS code using the b-Xi two-equation model.
38 Xi may be obtained by either the solution of the Xi transport
39 equation or from an algebraic expression. Both approaches are
40 based on Gulder's flame speed correlation which has been shown
41 to be appropriate by comparison with the results from the
42 spectral model.
43
44 Strain effects are encorporated directly into the Xi equation
45 but not in the algebraic approximation. Further work need to be
46 done on this issue, particularly regarding the enhanced removal rate
47 caused by flame compression. Analysis using results of the spectral
48 model will be required.
49
50 For cases involving very lean Propane flames or other flames which are
51 very strain-sensitive, a transport equation for the laminar flame
52 speed is present. This equation is derived using heuristic arguments
53 involving the strain time scale and the strain-rate at extinction.
54 the transport velocity is the same as that for the Xi equation.
55
56\*---------------------------------------------------------------------------*/
57
58#include "fvCFD.H"
59#include "dynamicFvMesh.H"
60#include "psiuReactionThermo.H"
62#include "laminarFlameSpeed.H"
63#include "ignition.H"
64#include "Switch.H"
65#include "pimpleControl.H"
66#include "CorrectPhi.H"
67#include "fvOptions.H"
68
69// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
70
71int main(int argc, char *argv[])
72{
73 argList::addNote
74 (
75 "Solver for compressible premixed/partially-premixed combustion with"
76 " turbulence modelling."
77 );
78
79 #include "postProcess.H"
80
81 #include "setRootCaseLists.H"
82 #include "createTime.H"
83 #include "createDynamicFvMesh.H"
84 #include "createControl.H"
85 #include "readCombustionProperties.H"
86 #include "readGravitationalAcceleration.H"
87 #include "createFields.H"
88 #include "createFieldRefs.H"
89 #include "initContinuityErrs.H"
90 #include "createRhoUf.H"
91 #include "createControls.H"
92 #include "initContinuityErrs.H"
93 #include "compressibleCourantNo.H"
94 #include "setInitialDeltaT.H"
95
96 turbulence->validate();
97
98 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
99
100 Info<< "\nStarting time loop\n" << endl;
101
102 while (runTime.run())
103 {
104 #include "createTimeControls.H"
105
106 {
107 // Store divrhoU from the previous mesh so that it can be mapped
108 // and used in correctPhi to ensure the corrected phi has the
109 // same divergence
110 volScalarField divrhoU
111 (
112 "divrhoU",
113 fvc::div(fvc::absolute(phi, rho, U))
114 );
115
116 #include "compressibleCourantNo.H"
117 #include "setDeltaT.H"
118
119 ++runTime;
120
121 Info<< "Time = " << runTime.timeName() << nl << endl;
122
123 // Store momentum to set rhoUf for introduced faces.
124 volVectorField rhoU("rhoU", rho*U);
125
126 // Do any mesh changes
127 mesh.update();
128
129 if (mesh.changing() && correctPhi)
130 {
131 // Calculate absolute flux from the mapped surface velocity
132 phi = mesh.Sf() & rhoUf;
133
134 #include "correctPhi.H"
135
136 // Make the fluxes relative to the mesh-motion
137 fvc::makeRelative(phi, rho, U);
138 }
139 }
140
141 if (mesh.changing() && checkMeshCourantNo)
142 {
143 #include "meshCourantNo.H"
144 }
145
146 #include "rhoEqn.H"
147 Info<< "rho min/max : " << min(rho).value() << " " << max(rho).value()
148 << endl;
149
150 // --- Pressure-velocity PIMPLE corrector loop
151 while (pimple.loop())
152 {
153 #include "UEqn.H"
154
155 #include "ftEqn.H"
156 #include "bEqn.H"
157 #include "EauEqn.H"
158 #include "EaEqn.H"
159
160 if (!ign.ignited())
161 {
162 thermo.heu() == thermo.he();
163 }
164
165 // --- Pressure corrector loop
166 while (pimple.correct())
167 {
168 #include "pEqn.H"
169 }
170
171 if (pimple.turbCorr())
172 {
173 turbulence->correct();
174 }
175 }
176
177 rho = thermo.rho();
178
179 runTime.write();
180
181 runTime.printExecutionTime(Info);
182 }
183
184 Info<< "End\n" << endl;
185
186 return 0;
187}
188
189
190// ************************************************************************* //
Y[inertIndex] max(0.0)
surfaceScalarField & phi
pimpleControl & pimple
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
engineTime & runTime
autoPtr< surfaceVectorField > rhoUf
Creates and initialises the velocity velocity field rhoUf.
Read the control parameters used by setDeltaT.
compressible::turbulenceModel & turbulence
Calculates and outputs the mean and maximum Courant Numbers.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:83
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:82
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Execute application functionObjects to post-process existing results.
checkMeshCourantNo
correctPhi