OpenFOAM: API Guide
v2112
The open source CFD toolbox
interIsoFoam.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) 2016 DHI
10
Copyright (C) 2017 OpenCFD Ltd.
11
Copyright (C) 2018 Johan Roenby
12
Copyright (C) 2019-2020 DLR
13
Copyright (C) 2020 OpenCFD Ltd.
14
-------------------------------------------------------------------------------
15
License
16
This file is part of OpenFOAM.
17
18
OpenFOAM is free software: you can redistribute it and/or modify it
19
under the terms of the GNU General Public License as published by
20
the Free Software Foundation, either version 3 of the License, or
21
(at your option) any later version.
22
23
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
24
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
25
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
26
for more details.
27
28
You should have received a copy of the GNU General Public License
29
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
30
31
Application
32
interIsoFoam
33
34
Group
35
grpMultiphaseSolvers
36
37
Description
38
Solver derived from interFoam for two incompressible, isothermal immiscible
39
fluids using the isoAdvector phase-fraction based interface capturing
40
approach, with optional mesh motion and mesh topology changes including
41
adaptive re-meshing.
42
43
Reference:
44
\verbatim
45
Roenby, J., Bredmose, H. and Jasak, H. (2016).
46
A computational method for sharp interface advection
47
Royal Society Open Science, 3
48
doi 10.1098/rsos.160405
49
\endverbatim
50
51
isoAdvector code supplied by Johan Roenby, STROMNING (2018)
52
53
\*---------------------------------------------------------------------------*/
54
55
#include "
fvCFD.H
"
56
#include "
dynamicFvMesh.H
"
57
#include "
isoAdvection.H
"
58
#include "
EulerDdtScheme.H
"
59
#include "
localEulerDdtScheme.H
"
60
#include "
CrankNicolsonDdtScheme.H
"
61
#include "
subCycle.H
"
62
#include "
immiscibleIncompressibleTwoPhaseMixture.H
"
63
#include "
incompressibleInterPhaseTransportModel.H
"
64
#include "
pimpleControl.H
"
65
#include "
fvOptions.H
"
66
#include "
CorrectPhi.H
"
67
#include "
fvcSmooth.H
"
68
#include "
dynamicRefineFvMesh.H
"
69
70
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
71
72
int
main(
int
argc,
char
*argv[])
73
{
74
argList::addNote
75
(
76
"Solver for two incompressible, isothermal immiscible fluids"
77
" using isoAdvector phase-fraction based interface capturing.\n"
78
"With optional mesh motion and mesh topology changes including"
79
" adaptive re-meshing.\n"
80
"The solver is derived from interFoam"
81
);
82
83
#include "
postProcess.H
"
84
85
#include "
addCheckCaseOptions.H
"
86
#include "
setRootCaseLists.H
"
87
#include "
createTime.H
"
88
#include "
createDynamicFvMesh.H
"
89
#include "initContinuityErrs.H"
90
#include "createDyMControls.H"
91
#include "createFields.H"
92
#include "
initCorrectPhi.H
"
93
#include "
createUfIfPresent.H
"
94
95
#include "
porousCourantNo.H
"
96
#include "setInitialDeltaT.H"
97
98
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
99
Info
<<
"\nStarting time loop\n"
<<
endl
;
100
101
while
(
runTime
.run())
102
{
103
#include "
readDyMControls.H
"
104
#include "
porousCourantNo.H
"
105
#include "
porousAlphaCourantNo.H
"
106
#include "setDeltaT.H"
107
108
++
runTime
;
109
110
Info
<<
"Time = "
<<
runTime
.timeName() <<
nl
<<
endl
;
111
112
// --- Pressure-velocity PIMPLE corrector loop
113
while
(
pimple
.loop())
114
{
115
if
(
pimple
.firstIter() ||
moveMeshOuterCorrectors
)
116
{
117
if
(isA<dynamicRefineFvMesh>(
mesh
))
118
{
119
advector
.surf().reconstruct();
120
}
121
122
mesh
.update();
123
124
if
(
mesh
.changing())
125
{
126
gh
= (
g
&
mesh
.C()) - ghRef;
127
ghf
= (
g
&
mesh
.Cf()) - ghRef;
128
129
if
(isA<dynamicRefineFvMesh>(
mesh
))
130
{
131
advector
.surf().mapAlphaField();
132
alpha2
= 1.0 -
alpha1
;
133
alpha2
.correctBoundaryConditions();
134
rho
==
alpha1
*
rho1
+
alpha2
*
rho2
;
135
rho
.correctBoundaryConditions();
136
rho
.oldTime() =
rho
;
137
alpha2
.oldTime() =
alpha2
;
138
}
139
140
MRF
.update();
141
142
if
(
correctPhi
)
143
{
144
// Calculate absolute flux
145
// from the mapped surface velocity
146
phi
=
mesh
.Sf() &
Uf
();
147
148
#include "correctPhi.H"
149
150
// Make the flux relative to the mesh motion
151
fvc::makeRelative(
phi
,
U
);
152
153
mixture
.correct();
154
}
155
156
if
(
checkMeshCourantNo
)
157
{
158
#include "
meshCourantNo.H
"
159
}
160
}
161
}
162
163
#include "alphaControls.H"
164
#include "alphaEqnSubCycle.H"
165
166
mixture
.correct();
167
168
if
(
pimple
.frozenFlow())
169
{
170
continue
;
171
}
172
173
#include "UEqn.H"
174
175
// --- Pressure corrector loop
176
while
(
pimple
.correct())
177
{
178
#include "pEqn.H"
179
}
180
181
if
(
pimple
.turbCorr())
182
{
183
turbulence
->correct();
184
}
185
}
186
187
runTime
.write();
188
189
runTime
.printExecutionTime(Info);
190
}
191
192
Info
<<
"End\n"
<<
endl
;
193
194
return
0;
195
}
196
197
198
// ************************************************************************* //
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
alpha1
const volScalarField & alpha1
Definition:
setRegionFluidFields.H:8
rho2
volScalarField & rho2
Definition:
setRegionFluidFields.H:30
alpha2
const volScalarField & alpha2
Definition:
setRegionFluidFields.H:9
rho1
volScalarField & rho1
Definition:
setRegionFluidFields.H:27
U
U
Definition:
pEqn.H:72
createDynamicFvMesh.H
mesh
dynamicFvMesh & mesh
Definition:
createDynamicFvMesh.H:6
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
createTime.H
createUfIfPresent.H
Creates and initialises the velocity field Uf if required.
Uf
autoPtr< surfaceVectorField > Uf
Definition:
createUfIfPresent.H:33
dynamicFvMesh.H
dynamicRefineFvMesh.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
isoAdvection.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
porousAlphaCourantNo.H
Calculates and outputs the mean and maximum Courant Numbers.
porousCourantNo.H
Calculates and outputs the mean and maximum Courant Numbers.
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
rho
rho
Definition:
readInitialConditions.H:88
setRootCaseLists.H
mixture
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition:
createFields.H:39
advector
isoAdvection advector(alpha1, phi, U)
subCycle.H
applications
solvers
multiphase
interIsoFoam
interIsoFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.