OpenFOAM: API Guide
v2112
The open source CFD toolbox
compressibleInterFilmFoam.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
-------------------------------------------------------------------------------
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
compressibleInterFoam
28
29
Description
30
Solver for two compressible, non-isothermal immiscible fluids using a VOF
31
(volume of fluid) phase-fraction based interface capturing approach.
32
33
The momentum and other fluid properties are of the "mixture" and a single
34
momentum equation is solved.
35
36
Turbulence modelling is generic, i.e. laminar, RAS or LES may be selected.
37
38
\*---------------------------------------------------------------------------*/
39
40
#include "
fvCFD.H
"
41
#include "
CMULES.H
"
42
#include "
EulerDdtScheme.H
"
43
#include "
localEulerDdtScheme.H
"
44
#include "
CrankNicolsonDdtScheme.H
"
45
#include "
subCycle.H
"
46
#include "
compressibleInterPhaseTransportModel.H
"
47
#include "
pimpleControl.H
"
48
#include "
SLGThermo.H
"
49
#include "
surfaceFilmModel.H
"
50
#include "
pimpleControl.H
"
51
#include "
fvOptions.H
"
52
#include "
fvcSmooth.H
"
53
54
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55
56
int
main(
int
argc,
char
*argv[])
57
{
58
argList::addNote
59
(
60
"Solver for two compressible, non-isothermal immiscible fluids"
61
" using VOF phase-fraction based interface capturing."
62
);
63
64
#include "
postProcess.H
"
65
66
#include "
addCheckCaseOptions.H
"
67
#include "
setRootCaseLists.H
"
68
#include "
createTime.H
"
69
#include "
createMesh.H
"
70
#include "createControl.H"
71
#include "
createTimeControls.H
"
72
#include "createFields.H"
73
#include "createSurfaceFilmModel.H"
74
75
volScalarField
&
p
=
mixture
.p();
76
volScalarField
&
T
=
mixture
.T();
77
const
volScalarField
&
psi1
=
mixture
.thermo1().psi();
78
const
volScalarField
&
psi2
=
mixture
.thermo2().psi();
79
80
regionModels::surfaceFilmModel&
surfaceFilm
=
tsurfaceFilm
();
81
82
if
(!
LTS
)
83
{
84
#include "
readTimeControls.H
"
85
#include "CourantNo.H"
86
#include "setInitialDeltaT.H"
87
}
88
89
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
90
91
Info
<<
"\nStarting time loop\n"
<<
endl
;
92
93
while
(
runTime
.run())
94
{
95
#include "
readTimeControls.H
"
96
97
if
(
LTS
)
98
{
99
#include "setRDeltaT.H"
100
}
101
else
102
{
103
#include "CourantNo.H"
104
#include "alphaCourantNo.H"
105
#include "setDeltaT.H"
106
}
107
108
++
runTime
;
109
110
Info
<<
"Time = "
<<
runTime
.timeName() <<
nl
<<
endl
;
111
112
surfaceFilm
.evolve();
113
114
// --- Pressure-velocity PIMPLE corrector loop
115
while
(
pimple
.loop())
116
{
117
#include "alphaControls.H"
118
#include "compressibleAlphaEqnSubCycle.H"
119
120
turbulence
.correctPhasePhi();
121
122
volScalarField::Internal Srho(
surfaceFilm
.Srho());
123
contErr
-=
posPart
(Srho);
124
125
#include "UEqn.H"
126
#include "TEqn.H"
127
128
// --- Pressure corrector loop
129
while
(
pimple
.correct())
130
{
131
#include "pEqn.H"
132
}
133
134
if
(
pimple
.turbCorr())
135
{
136
turbulence
.correct();
137
}
138
}
139
140
runTime
.write();
141
142
runTime
.printExecutionTime(Info);
143
}
144
145
Info
<<
"End\n"
<<
endl
;
146
147
return
0;
148
}
149
150
151
// ************************************************************************* //
CMULES.H
CMULES: Multidimensional universal limiter for explicit corrected implicit solution.
CrankNicolsonDdtScheme.H
EulerDdtScheme.H
SLGThermo.H
addCheckCaseOptions.H
Required Classes.
pimple
pimpleControl & pimple
Definition:
setRegionFluidFields.H:56
psi2
const volScalarField & psi2
Definition:
setRegionFluidFields.H:31
psi1
const volScalarField & psi1
Definition:
setRegionFluidFields.H:28
p
volScalarField & p
Definition:
createFieldRefs.H:8
T
const volScalarField & T
Definition:
createFieldRefs.H:2
surfaceFilm
regionModels::surfaceFilmModel & surfaceFilm
Definition:
createFieldRefs.H:3
tsurfaceFilm
Info<< "\nConstructing surface film model"<< endl;autoPtr< regionModels::surfaceFilmModel > tsurfaceFilm(regionModels::surfaceFilmModel::New(mesh, g))
contErr
volScalarField::Internal contErr((fvc::ddt(rho)+fvc::div(rhoPhi) -(fvOptions(alpha1, mixture.thermo1().rho())&rho1) -(fvOptions(alpha2, mixture.thermo2().rho())&rho2))())
compressibleInterPhaseTransportModel.H
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
createMesh.H
Required Variables.
LTS
bool LTS
Definition:
createRDeltaT.H:1
createTimeControls.H
Read the control parameters used by setDeltaT.
createTime.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...
localEulerDdtScheme.H
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::posPart
dimensionedScalar posPart(const dimensionedScalar &ds)
Definition:
dimensionedScalar.C:221
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.
readTimeControls.H
Read the control parameters used by setDeltaT.
setRootCaseLists.H
mixture
Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n"<< endl;autoPtr< temperaturePhaseChangeTwoPhaseMixture > mixture
Definition:
createFields.H:39
subCycle.H
surfaceFilmModel.H
applications
solvers
multiphase
compressibleInterFoam
compressibleInterFilmFoam
compressibleInterFilmFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.