twoLiquidMixingFoam.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-2016 OpenFOAM Foundation
9-------------------------------------------------------------------------------
10License
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
26Application
27 twoLiquidMixingFoam
28
29Group
30 grpMultiphaseSolvers
31
32Description
33 Solver for mixing two incompressible fluids.
34
35 Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
36
37\*---------------------------------------------------------------------------*/
38
39#include "fvCFD.H"
40#include "MULES.H"
41#include "subCycle.H"
44#include "pimpleControl.H"
45
46// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
47
48int main(int argc, char *argv[])
49{
50 argList::addNote
51 (
52 "Solver for mixing two incompressible fluids"
53 );
54
55 #include "postProcess.H"
56
57 #include "addCheckCaseOptions.H"
58 #include "setRootCaseLists.H"
59 #include "createTime.H"
60 #include "createMesh.H"
61 #include "createControl.H"
62 #include "initContinuityErrs.H"
63 #include "createFields.H"
64 #include "createTimeControls.H"
65 #include "CourantNo.H"
66 #include "setInitialDeltaT.H"
67
68 turbulence->validate();
69
70 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
71
72 Info<< "\nStarting time loop\n" << endl;
73
74 while (runTime.run())
75 {
76 #include "readTimeControls.H"
77 #include "CourantNo.H"
78 #include "alphaCourantNo.H"
79 #include "setDeltaT.H"
80
81 ++runTime;
82
83 Info<< "Time = " << runTime.timeName() << nl << endl;
84
85 mixture.correct();
86
87 #include "alphaEqnSubCycle.H"
88 #include "alphaDiffusionEqn.H"
89
90 // --- Pressure-velocity PIMPLE corrector loop
91 while (pimple.loop())
92 {
93 #include "UEqn.H"
94
95 // --- Pressure corrector loop
96 while (pimple.correct())
97 {
98 #include "pEqn.H"
99 }
100
101 if (pimple.turbCorr())
102 {
103 turbulence->correct();
104 }
105 }
106
107 runTime.write();
108
109 runTime.printExecutionTime(Info);
110 }
111
112 Info<< "End\n" << endl;
113
114 return 0;
115}
116
117
118// ************************************************************************* //
MULES: Multidimensional universal limiter for explicit solution.
Required Classes.
pimpleControl & pimple
engineTime & runTime
Required Variables.
Read the control parameters used by setDeltaT.
compressible::turbulenceModel & turbulence
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.
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition: createFields.H:39