OpenFOAM: API Guide
v2112
The open source CFD toolbox
pimpleFoam.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) 2019 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
pimpleFoam.C
29
30
Group
31
grpIncompressibleSolvers
32
33
Description
34
Transient solver for incompressible, turbulent flow of Newtonian fluids
35
on a moving mesh.
36
37
\heading Solver details
38
The solver uses the PIMPLE (merged PISO-SIMPLE) algorithm to solve the
39
continuity equation:
40
41
\f[
42
\div \vec{U} = 0
43
\f]
44
45
and momentum equation:
46
47
\f[
48
\ddt{\vec{U}} + \div \left( \vec{U} \vec{U} \right) - \div \gvec{R}
49
= - \grad p + \vec{S}_U
50
\f]
51
52
Where:
53
\vartable
54
\vec{U} | Velocity
55
p | Pressure
56
\vec{R} | Stress tensor
57
\vec{S}_U | Momentum source
58
\endvartable
59
60
Sub-models include:
61
- turbulence modelling, i.e. laminar, RAS or LES
62
- run-time selectable MRF and finite volume options, e.g. explicit porosity
63
64
\heading Required fields
65
\plaintable
66
U | Velocity [m/s]
67
p | Kinematic pressure, p/rho [m2/s2]
68
<turbulence fields> | As required by user selection
69
\endplaintable
70
71
Note
72
The motion frequency of this solver can be influenced by the presence
73
of "updateControl" and "updateInterval" in the dynamicMeshDict.
74
75
\*---------------------------------------------------------------------------*/
76
77
#include "
fvCFD.H
"
78
#include "
dynamicFvMesh.H
"
79
#include "
singlePhaseTransportModel.H
"
80
#include "
turbulentTransportModel.H
"
81
#include "
pimpleControl.H
"
82
#include "
CorrectPhi.H
"
83
#include "
fvOptions.H
"
84
#include "
localEulerDdtScheme.H
"
85
#include "
fvcSmooth.H
"
86
87
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
88
89
int
main(
int
argc,
char
*argv[])
90
{
91
argList::addNote
92
(
93
"Transient solver for incompressible, turbulent flow"
94
" of Newtonian fluids on a moving mesh."
95
);
96
97
#include "
postProcess.H
"
98
99
#include "
addCheckCaseOptions.H
"
100
#include "
setRootCaseLists.H
"
101
#include "
createTime.H
"
102
#include "
createDynamicFvMesh.H
"
103
#include "initContinuityErrs.H"
104
#include "createDyMControls.H"
105
#include "createFields.H"
106
#include "
createUfIfPresent.H
"
107
#include "CourantNo.H"
108
#include "setInitialDeltaT.H"
109
110
turbulence
->validate();
111
112
if
(!
LTS
)
113
{
114
#include "CourantNo.H"
115
#include "setInitialDeltaT.H"
116
}
117
118
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
119
120
Info
<<
"\nStarting time loop\n"
<<
endl
;
121
122
while
(
runTime
.run())
123
{
124
#include "
readDyMControls.H
"
125
126
if
(
LTS
)
127
{
128
#include "setRDeltaT.H"
129
}
130
else
131
{
132
#include "CourantNo.H"
133
#include "setDeltaT.H"
134
}
135
136
++
runTime
;
137
138
Info
<<
"Time = "
<<
runTime
.timeName() <<
nl
<<
endl
;
139
140
// --- Pressure-velocity PIMPLE corrector loop
141
while
(
pimple
.loop())
142
{
143
if
(
pimple
.firstIter() ||
moveMeshOuterCorrectors
)
144
{
145
// Do any mesh changes
146
mesh
.controlledUpdate();
147
148
if
(
mesh
.changing())
149
{
150
MRF
.update();
151
152
if
(
correctPhi
)
153
{
154
// Calculate absolute flux
155
// from the mapped surface velocity
156
phi
=
mesh
.Sf() &
Uf
();
157
158
#include "correctPhi.H"
159
160
// Make the flux relative to the mesh motion
161
fvc::makeRelative(
phi
,
U
);
162
}
163
164
if
(
checkMeshCourantNo
)
165
{
166
#include "
meshCourantNo.H
"
167
}
168
}
169
}
170
171
#include "UEqn.H"
172
173
// --- Pressure corrector loop
174
while
(
pimple
.correct())
175
{
176
#include "pEqn.H"
177
}
178
179
if
(
pimple
.turbCorr())
180
{
181
laminarTransport
.correct();
182
turbulence
->correct();
183
}
184
}
185
186
runTime
.write();
187
188
runTime
.printExecutionTime(Info);
189
}
190
191
Info
<<
"End\n"
<<
endl
;
192
193
return
0;
194
}
195
196
197
// ************************************************************************* //
CorrectPhi.H
addCheckCaseOptions.H
Required Classes.
phi
surfaceScalarField & phi
Definition:
setRegionFluidFields.H:8
MRF
IOMRFZoneList & MRF
Definition:
setRegionFluidFields.H:22
pimple
pimpleControl & pimple
Definition:
setRegionFluidFields.H:56
U
U
Definition:
pEqn.H:72
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...
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
singlePhaseTransportModel.H
laminarTransport
singlePhaseTransportModel laminarTransport(U, phi)
turbulentTransportModel.H
applications
solvers
incompressible
pimpleFoam
pimpleFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.