OpenFOAM: API Guide
v2112
The open source CFD toolbox
overBuoyantPimpleDyMFoam.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) 2019 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
overBuoyantPimpleDymFoam
28
29
Group
30
grpHeatTransferSolvers
31
32
Description
33
Transient solver for buoyant, turbulent flow of compressible fluids for
34
ventilation and heat-transfer with overset feature
35
36
Turbulence is modelled using a run-time selectable compressible RAS or
37
LES model.
38
39
\*---------------------------------------------------------------------------*/
40
41
#include "
fvCFD.H
"
42
#include "
dynamicFvMesh.H
"
43
#include "
rhoThermo.H
"
44
#include "
turbulentFluidThermoModel.H
"
45
#include "
radiationModel.H
"
46
#include "
fvOptions.H
"
47
#include "
pimpleControl.H
"
48
#include "
pressureControl.H
"
49
50
#include "
CorrectPhi.H
"
51
#include "
cellCellStencilObject.H
"
52
#include "
localMin.H
"
53
54
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
55
56
int
main(
int
argc,
char
*argv[])
57
{
58
argList::addNote
59
(
60
"Transient solver for buoyant, turbulent fluid flow"
61
" of compressible fluids, including radiation."
62
);
63
64
#include "
postProcess.H
"
65
66
#include "
addCheckCaseOptions.H
"
67
#include "
setRootCaseLists.H
"
68
#include "
createTime.H
"
69
#include "
createDynamicFvMesh.H
"
70
#include "createDyMControls.H"
71
#include "createFields.H"
72
#include "createFieldRefs.H"
73
#include "initContinuityErrs.H"
74
75
#include "
createRhoUfIfPresent.H
"
76
#include "createControls.H"
77
78
#include "compressibleCourantNo.H"
79
#include "setInitialDeltaT.H"
80
81
turbulence
->validate();
82
83
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
84
85
Info
<<
"\nStarting time loop\n"
<<
endl
;
86
87
while
(
runTime
.run())
88
{
89
#include "
readTimeControls.H
"
90
91
#include "readControls.H"
92
#include "
readDyMControls.H
"
93
94
#include "compressibleCourantNo.H"
95
#include "setDeltaT.H"
96
97
// Store divrhoU from the previous mesh so that it can be mapped
98
// and used in correctPhi to ensure the corrected phi has the
99
// same divergence
100
autoPtr<volScalarField> divrhoU;
101
if
(
correctPhi
)
102
{
103
divrhoU.reset
104
(
105
new
volScalarField
106
(
107
"divrhoU"
,
108
fvc::div(fvc::absolute(
phi
,
rho
,
U
))
109
)
110
);
111
}
112
113
114
++
runTime
;
115
116
Info
<<
"Time = "
<<
runTime
.timeName() <<
nl
<<
endl
;
117
118
// --- Pressure-velocity PIMPLE corrector loop
119
while
(
pimple
.loop())
120
{
121
if
(
pimple
.firstIter() ||
moveMeshOuterCorrectors
)
122
{
123
// Do any mesh changes
124
mesh
.update();
125
126
if
(
mesh
.changing())
127
{
128
MRF
.update();
129
130
#include "
setCellMask.H
"
131
132
const
surfaceScalarField
faceMaskOld
133
(
134
localMin<scalar>(
mesh
).
interpolate
(cellMask.oldTime())
135
);
136
137
// Zero Uf on old faceMask (H-I)
138
rhoUf
() *= faceMaskOld;
139
140
//fvc::correctRhoUf(rhoUfint, rho, U, phi);
141
surfaceVectorField
rhoUfint(fvc::interpolate(
rho
*
U
));
142
143
// Update Uf and phi on new C-I faces
144
rhoUf
() += (1-faceMaskOld)*rhoUfint;
145
146
// Update Uf boundary
147
forAll
(
rhoUf
().boundaryField(), patchI)
148
{
149
rhoUf
().boundaryFieldRef()[patchI] =
150
rhoUfint.boundaryField()[patchI];
151
}
152
153
// Calculate absolute flux from the mapped surface velocity
154
phi
=
mesh
.Sf() &
rhoUf
();
155
156
if
(
correctPhi
)
157
{
158
#include "correctPhi.H"
159
}
160
161
// Zero phi on current H-I
162
const
surfaceScalarField
faceMask
163
(
164
localMin<scalar>(
mesh
).
interpolate
(cellMask)
165
);
166
167
phi
*=
faceMask
;
168
U
*= cellMask;
169
170
// Make the fluxes relative to the mesh-motion
171
fvc::makeRelative(
phi
,
rho
,
U
);
172
}
173
174
if
(
checkMeshCourantNo
)
175
{
176
#include "
meshCourantNo.H
"
177
}
178
}
179
180
if
(
pimple
.firstIter())
181
{
182
#include "rhoEqn.H"
183
}
184
185
#include "UEqn.H"
186
#include "EEqn.H"
187
188
// --- Pressure corrector loop
189
while
(
pimple
.correct())
190
{
191
#include "pEqn.H"
192
}
193
194
if
(
pimple
.turbCorr())
195
{
196
turbulence
->correct();
197
}
198
}
199
200
rho
=
thermo
.rho();
201
202
runTime
.write();
203
204
runTime
.printExecutionTime(Info);
205
}
206
207
Info
<<
"End\n"
<<
endl
;
208
209
return
0;
210
}
211
212
213
// ************************************************************************* //
CorrectPhi.H
addCheckCaseOptions.H
Required Classes.
cellCellStencilObject.H
phi
surfaceScalarField & phi
Definition:
setRegionFluidFields.H:8
MRF
IOMRFZoneList & MRF
Definition:
setRegionFluidFields.H:22
pimple
pimpleControl & pimple
Definition:
setRegionFluidFields.H:56
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
faceMask
surfaceScalarField faceMask(localMin< scalar >(mesh).interpolate(cellMask))
interpolate
mesh interpolate(rAU)
createDynamicFvMesh.H
mesh
dynamicFvMesh & mesh
Definition:
createDynamicFvMesh.H:6
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
createRhoUfIfPresent.H
Creates and initialises the velocity field rhoUf if required.
rhoUf
autoPtr< surfaceVectorField > rhoUf
Definition:
createRhoUfIfPresent.H:33
createTime.H
dynamicFvMesh.H
turbulence
compressible::turbulenceModel & turbulence
Definition:
setRegionFluidFields.H:30
fvCFD.H
fvOptions.H
localMin.H
meshCourantNo.H
Calculates and outputs the mean and maximum Courant Numbers.
Foam::surfaceVectorField
GeometricField< vector, fvsPatchField, surfaceMesh > surfaceVectorField
Definition:
surfaceFieldsFwd.H:67
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::surfaceScalarField
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
Definition:
surfaceFieldsFwd.H:66
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.
pressureControl.H
radiationModel.H
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
readTimeControls.H
Read the control parameters used by setDeltaT.
rhoThermo.H
setCellMask.H
Sets blocked cells mask field.
setRootCaseLists.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition:
stdFoam.H:333
turbulentFluidThermoModel.H
applications
solvers
heatTransfer
buoyantPimpleFoam
overBuoyantPimpleDyMFoam
overBuoyantPimpleDyMFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.