OpenFOAM: API Guide
v2112
The open source CFD toolbox
interFoam.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
interFoam
29
30
Group
31
grpMultiphaseSolvers
32
33
Description
34
Solver for two incompressible, isothermal immiscible fluids using a VOF
35
(volume of fluid) phase-fraction based interface capturing approach,
36
with optional mesh motion and mesh topology changes including adaptive
37
re-meshing.
38
39
\*---------------------------------------------------------------------------*/
40
41
#include "
fvCFD.H
"
42
#include "
dynamicFvMesh.H
"
43
#include "
CMULES.H
"
44
#include "
EulerDdtScheme.H
"
45
#include "
localEulerDdtScheme.H
"
46
#include "
CrankNicolsonDdtScheme.H
"
47
#include "
subCycle.H
"
48
#include "
immiscibleIncompressibleTwoPhaseMixture.H
"
49
#include "
incompressibleInterPhaseTransportModel.H
"
50
#include "
turbulentTransportModel.H
"
51
#include "
pimpleControl.H
"
52
#include "
fvOptions.H
"
53
#include "
CorrectPhi.H
"
54
#include "
fvcSmooth.H
"
55
56
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
57
58
int
main(
int
argc,
char
*argv[])
59
{
60
argList::addNote
61
(
62
"Solver for two incompressible, isothermal immiscible fluids"
63
" using VOF phase-fraction based interface capturing.\n"
64
"With optional mesh motion and mesh topology changes including"
65
" adaptive re-meshing."
66
);
67
68
#include "
postProcess.H
"
69
70
#include "
addCheckCaseOptions.H
"
71
#include "
setRootCaseLists.H
"
72
#include "
createTime.H
"
73
#include "
createDynamicFvMesh.H
"
74
#include "initContinuityErrs.H"
75
#include "createDyMControls.H"
76
#include "createFields.H"
77
#include "
createAlphaFluxes.H
"
78
#include "
initCorrectPhi.H
"
79
#include "
createUfIfPresent.H
"
80
81
if
(!
LTS
)
82
{
83
#include "CourantNo.H"
84
#include "setInitialDeltaT.H"
85
}
86
87
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
88
Info
<<
"\nStarting time loop\n"
<<
endl
;
89
90
while
(
runTime
.run())
91
{
92
#include "
readDyMControls.H
"
93
94
if
(
LTS
)
95
{
96
#include "setRDeltaT.H"
97
}
98
else
99
{
100
#include "CourantNo.H"
101
#include "alphaCourantNo.H"
102
#include "setDeltaT.H"
103
}
104
105
++
runTime
;
106
107
Info
<<
"Time = "
<<
runTime
.timeName() <<
nl
<<
endl
;
108
109
// --- Pressure-velocity PIMPLE corrector loop
110
while
(
pimple
.loop())
111
{
112
if
(
pimple
.firstIter() ||
moveMeshOuterCorrectors
)
113
{
114
mesh
.update();
115
116
if
(
mesh
.changing())
117
{
118
// Do not apply previous time-step mesh compression flux
119
// if the mesh topology changed
120
if
(
mesh
.topoChanging())
121
{
122
talphaPhi1Corr0
.clear();
123
}
124
125
gh
= (
g
&
mesh
.C()) - ghRef;
126
ghf
= (
g
&
mesh
.Cf()) - ghRef;
127
128
MRF
.update();
129
130
if
(
correctPhi
)
131
{
132
// Calculate absolute flux
133
// from the mapped surface velocity
134
phi
=
mesh
.Sf() &
Uf
();
135
136
#include "correctPhi.H"
137
138
// Make the flux relative to the mesh motion
139
fvc::makeRelative(
phi
,
U
);
140
141
mixture
.correct();
142
}
143
144
if
(
checkMeshCourantNo
)
145
{
146
#include "
meshCourantNo.H
"
147
}
148
}
149
}
150
151
#include "alphaControls.H"
152
#include "alphaEqnSubCycle.H"
153
154
mixture
.correct();
155
156
if
(
pimple
.frozenFlow())
157
{
158
continue
;
159
}
160
161
#include "UEqn.H"
162
163
// --- Pressure corrector loop
164
while
(
pimple
.correct())
165
{
166
#include "pEqn.H"
167
}
168
169
if
(
pimple
.turbCorr())
170
{
171
turbulence
->correct();
172
}
173
}
174
175
runTime
.write();
176
177
runTime
.printExecutionTime(Info);
178
}
179
180
Info
<<
"End\n"
<<
endl
;
181
182
return
0;
183
}
184
185
186
// ************************************************************************* //
CMULES.H
CMULES: Multidimensional universal limiter for explicit corrected implicit solution.
CorrectPhi.H
CrankNicolsonDdtScheme.H
EulerDdtScheme.H
addCheckCaseOptions.H
Required Classes.
g
const uniformDimensionedVectorField & g
Definition:
createFluidFields.H:26
phi
surfaceScalarField & phi
Definition:
setRegionFluidFields.H:8
ghf
const surfaceScalarField & ghf
Definition:
setRegionFluidFields.H:18
MRF
IOMRFZoneList & MRF
Definition:
setRegionFluidFields.H:22
gh
const volScalarField & gh
Definition:
setRegionFluidFields.H:17
pimple
pimpleControl & pimple
Definition:
setRegionFluidFields.H:56
U
U
Definition:
pEqn.H:72
createAlphaFluxes.H
talphaPhi1Corr0
tmp< surfaceScalarField > talphaPhi1Corr0
Definition:
createAlphaFluxes.H:26
createDynamicFvMesh.H
mesh
dynamicFvMesh & mesh
Definition:
createDynamicFvMesh.H:6
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
LTS
bool LTS
Definition:
createRDeltaT.H:1
createTime.H
createUfIfPresent.H
Creates and initialises the velocity field Uf if required.
Uf
autoPtr< surfaceVectorField > Uf
Definition:
createUfIfPresent.H:33
dynamicFvMesh.H
turbulence
compressible::turbulenceModel & turbulence
Definition:
setRegionFluidFields.H:30
fvCFD.H
fvOptions.H
fvcSmooth.H
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
immiscibleIncompressibleTwoPhaseMixture.H
incompressibleInterPhaseTransportModel.H
initCorrectPhi.H
localEulerDdtScheme.H
meshCourantNo.H
Calculates and outputs the mean and maximum Courant Numbers.
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
Foam::nl
constexpr char nl
The newline '\n' character (0x0a)
Definition:
Ostream.H:53
pimpleControl.H
postProcess.H
Execute application functionObjects to post-process existing results.
readDyMControls.H
checkMeshCourantNo
checkMeshCourantNo
Definition:
readDyMControls.H:9
correctPhi
correctPhi
Definition:
readDyMControls.H:3
moveMeshOuterCorrectors
moveMeshOuterCorrectors
Definition:
readDyMControls.H:15
setRootCaseLists.H
mixture
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition:
createFields.H:39
subCycle.H
turbulentTransportModel.H
applications
solvers
multiphase
interFoam
interFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.