reactingMultiphaseEulerFoam.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-2018 OpenFOAM Foundation
9 Copyright (C) 2020 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 reactingMultiphaseEulerFoam
29
30Description
31 Solver for a system of any number of compressible fluid phases with a
32 common pressure, but otherwise separate properties. The type of phase model
33 is run time selectable and can optionally represent multiple species and
34 in-phase reactions. The phase system is also run time selectable and can
35 optionally represent different types of momentum, heat and mass transfer.
36
37\*---------------------------------------------------------------------------*/
38
39#include "fvCFD.H"
40#include "multiphaseSystem.H"
41#include "phaseCompressibleTurbulenceModel.H"
42#include "pimpleControl.H"
43#include "localEulerDdtScheme.H"
44#include "fvcSmooth.H"
45
46// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47
48int main(int argc, char *argv[])
49{
50 #include "postProcess.H"
51
52 #include "addCheckCaseOptions.H"
53 #include "setRootCase.H"
54 #include "createTime.H"
55 #include "createMesh.H"
56 #include "createControl.H"
57 #include "createTimeControls.H"
58 #include "createFields.H"
59 #include "createFieldRefs.H"
60
61 if (!LTS)
62 {
63 #include "CourantNo.H"
64 #include "setInitialDeltaT.H"
65 }
66
67 Switch faceMomentum
68 (
69 pimple.dict().getOrDefault<Switch>("faceMomentum", false)
70 );
71 Switch partialElimination
72 (
73 pimple.dict().getOrDefault<Switch>("partialElimination", false)
74 );
75
76 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
77
78 Info<< "\nStarting time loop\n" << endl;
79
80 while (runTime.run())
81 {
82 #include "readTimeControls.H"
83
85 (
86 pimple.dict().getOrDefault<int>("nEnergyCorrectors", 1)
87 );
88
89 if (LTS)
90 {
91 #include "setRDeltaT.H"
92 }
93 else
94 {
95 #include "CourantNo.H"
96 #include "setDeltaT.H"
97 }
98
99 runTime++;
100 Info<< "Time = " << runTime.timeName() << nl << endl;
101
102 // --- Pressure-velocity PIMPLE corrector loop
103 while (pimple.loop())
104 {
105 fluid.solve();
106 fluid.correct();
107
108 #include "YEqns.H"
109
110 if (faceMomentum)
111 {
112 #include "pUf/UEqns.H"
113 #include "EEqns.H"
114 #include "pUf/pEqn.H"
115 }
116 else
117 {
118 #include "pU/UEqns.H"
119 #include "EEqns.H"
120 #include "pU/pEqn.H"
121 }
122
123 fluid.correctKinematics();
124
125 if (pimple.turbCorr())
126 {
127 fluid.correctTurbulence();
128 }
129 }
130
131 runTime.write();
132
133 runTime.printExecutionTime(Info);
134 }
135
136 Info<< "End\n" << endl;
137
138 return 0;
139}
140
141
142// ************************************************************************* //
Required Classes.
Switch faceMomentum(pimpleDict.getOrDefault< Switch >("faceMomentum", false))
int nEnergyCorrectors(pimpleDict.getOrDefault< int >("nEnergyCorrectors", 1))
pimpleControl & pimple
twoPhaseSystem & fluid
engineTime & runTime
Required Variables.
bool LTS
Definition: createRDeltaT.H:1
Read the control parameters used by setDeltaT.
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Execute application functionObjects to post-process existing results.
Read the control parameters used by setDeltaT.