OpenFOAM: API Guide
v2112
The open source CFD toolbox
interCondensatingEvaporatingFoam.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) 2016-2020 OpenCFD Ltd.
9
-------------------------------------------------------------------------------
10
License
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
26
Application
27
interCondensatingEvaporatingFoam
28
29
Group
30
grpMultiphaseSolvers
31
32
Description
33
Solver for two incompressible, non-isothermal immiscible fluids with
34
phase-change (evaporation-condensation) between a fluid and its vapour.
35
Uses a VOF (volume of fluid) phase-fraction based interface capturing
36
approach.
37
38
The momentum, energy and other fluid properties are of the "mixture" and a
39
single momentum equation is solved.
40
41
Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
42
43
\*---------------------------------------------------------------------------*/
44
45
#include "
fvCFD.H
"
46
#include "
dynamicFvMesh.H
"
47
#include "
CMULES.H
"
48
#include "
EulerDdtScheme.H
"
49
#include "
localEulerDdtScheme.H
"
50
#include "
CrankNicolsonDdtScheme.H
"
51
#include "
subCycle.H
"
52
#include "
interfaceProperties.H
"
53
#include "
twoPhaseMixtureEThermo.H
"
54
#include "
temperaturePhaseChangeTwoPhaseMixture.H
"
55
#include "
turbulentTransportModel.H
"
56
#include "
turbulenceModel.H
"
57
#include "
pimpleControl.H
"
58
#include "
fvOptions.H
"
59
#include "
CorrectPhi.H
"
60
61
62
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
63
64
int
main(
int
argc,
char
*argv[])
65
{
66
argList::addNote
67
(
68
"Solver for two incompressible, non-isothermal immiscible fluids with"
69
" phase-change,"
70
" using VOF phase-fraction based interface capturing.\n"
71
"With optional mesh motion and mesh topology changes including"
72
" adaptive re-meshing."
73
);
74
75
#include "
postProcess.H
"
76
77
#include "
addCheckCaseOptions.H
"
78
#include "
setRootCaseLists.H
"
79
#include "
createTime.H
"
80
#include "
createDynamicFvMesh.H
"
81
#include "initContinuityErrs.H"
82
#include "createDyMControls.H"
83
#include "createFields.H"
84
#include "
createAlphaFluxes.H
"
85
#include "
initCorrectPhi.H
"
86
#include "
createUfIfPresent.H
"
87
88
#include "CourantNo.H"
89
#include "setInitialDeltaT.H"
90
91
volScalarField
&
T
=
thermo
->T();
92
93
turbulence
->validate();
94
95
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
96
97
Info
<<
"\nStarting time loop\n"
<<
endl
;
98
99
while
(
runTime
.run())
100
{
101
#include "
readDyMControls.H
"
102
103
// Store divU from the previous mesh so that it can be mapped
104
// and used in correctPhi to ensure the corrected phi has the
105
// same divergence
106
107
volScalarField
divU
(
"divU"
, fvc::div(fvc::absolute(
phi
,
U
)));
108
109
{
110
#include "CourantNo.H"
111
#include "alphaCourantNo.H"
112
#include "setDeltaT.H"
113
}
114
115
++
runTime
;
116
117
Info
<<
"Time = "
<<
runTime
.timeName() <<
nl
<<
endl
;
118
119
// --- Pressure-velocity PIMPLE corrector loop
120
while
(
pimple
.loop())
121
{
122
if
(
pimple
.firstIter() ||
moveMeshOuterCorrectors
)
123
{
124
mesh
.update();
125
126
if
(
mesh
.changing())
127
{
128
// Do not apply previous time-step mesh compression flux
129
// if the mesh topology changed
130
if
(
mesh
.topoChanging())
131
{
132
talphaPhi1Corr0
.clear();
133
}
134
135
gh
= (
g
&
mesh
.C()) - ghRef;
136
ghf
= (
g
&
mesh
.Cf()) - ghRef;
137
138
MRF
.update();
139
140
if
(
correctPhi
)
141
{
142
// Calculate absolute flux
143
// from the mapped surface velocity
144
phi
=
mesh
.Sf() &
Uf
();
145
146
#include "correctPhi.H"
147
148
// Make the flux relative to the mesh motion
149
fvc::makeRelative(
phi
,
U
);
150
}
151
152
if
(
checkMeshCourantNo
)
153
{
154
#include "
meshCourantNo.H
"
155
}
156
}
157
}
158
159
mixture
->correct();
160
161
#include "alphaControls.H"
162
#include "alphaEqnSubCycle.H"
163
164
interface
.correct();
165
166
#include "UEqn.H"
167
#include "TEqn.H"
168
169
// --- Pressure corrector loop
170
while
(
pimple
.correct())
171
{
172
#include "pEqn.H"
173
}
174
175
if
(
pimple
.turbCorr())
176
{
177
turbulence
->correct();
178
}
179
}
180
181
rho
=
alpha1
*
rho1
+
alpha2
*
rho2
;
182
183
runTime
.write();
184
185
runTime
.printExecutionTime(Info);
186
}
187
188
Info
<<
"End\n"
<<
endl
;
189
190
return
0;
191
}
192
193
194
// ************************************************************************* //
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
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
thermo
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
U
U
Definition:
pEqn.H:72
T
const volScalarField & T
Definition:
createFieldRefs.H:2
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
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
initCorrectPhi.H
divU
zeroField divU
Definition:
alphaSuSp.H:3
interfaceProperties.H
localEulerDdtScheme.H
meshCourantNo.H
Calculates and outputs the mean and maximum Courant Numbers.
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition:
volFieldsFwd.H:82
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
rho
rho
Definition:
readInitialConditions.H:88
setRootCaseLists.H
interface
interfaceProperties interface(alpha1, U, thermo->transportPropertiesDict())
mixture
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition:
createFields.H:39
subCycle.H
temperaturePhaseChangeTwoPhaseMixture.H
turbulenceModel.H
turbulentTransportModel.H
twoPhaseMixtureEThermo.H
applications
solvers
multiphase
interCondensatingEvaporatingFoam
interCondensatingEvaporatingFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.