buoyantBoussinesqPimpleFoam.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) 2021 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
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
27Application
28 buoyantBoussinesqPimpleFoam
29
30Group
31 grpHeatTransferSolvers
32
33Description
34 Transient solver for buoyant, turbulent flow of incompressible fluids,
35 with optional mesh motion and mesh topology changes.
36
37 Uses the Boussinesq approximation:
38 \f[
39 rho_{k} = 1 - beta(T - T_{ref})
40 \f]
41
42 where:
43 \f$ rho_{k} \f$ = the effective (driving) kinematic density
44 beta = thermal expansion coefficient [1/K]
45 T = temperature [K]
46 \f$ T_{ref} \f$ = reference temperature [K]
47
48 Valid when:
49 \f[
50 \frac{beta(T - T_{ref})}{rho_{ref}} << 1
51 \f]
52
53\*---------------------------------------------------------------------------*/
54
55#include "fvCFD.H"
56#include "dynamicFvMesh.H"
59#include "radiationModel.H"
60#include "CorrectPhi.H"
61#include "fvOptions.H"
62#include "pimpleControl.H"
63
64// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
65
66int main(int argc, char *argv[])
67{
68 argList::addNote
69 (
70 "Transient solver for buoyant, turbulent flow"
71 " of incompressible fluids, with optional mesh"
72 " motion and mesh topology changes.\n"
73 "Uses the Boussinesq approximation."
74 );
75
76 #include "postProcess.H"
77
78 #include "addCheckCaseOptions.H"
79 #include "setRootCaseLists.H"
80 #include "createTime.H"
81 #include "createDynamicFvMesh.H"
82 #include "createDyMControls.H"
83 #include "createFields.H"
84 #include "createUfIfPresent.H"
85 #include "CourantNo.H"
86 #include "setInitialDeltaT.H"
87 #include "initContinuityErrs.H"
88
89 turbulence->validate();
90
91 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
92
93 Info<< "\nStarting time loop\n" << endl;
94
95 while (runTime.run())
96 {
97 #include "readDyMControls.H"
98 #include "CourantNo.H"
99 #include "setDeltaT.H"
100
101 ++runTime;
102
103 Info<< "Time = " << runTime.timeName() << nl << endl;
104
105 // --- Pressure-velocity PIMPLE corrector loop
106 while (pimple.loop())
107 {
108 if (pimple.firstIter() || moveMeshOuterCorrectors)
109 {
110 // Do any mesh changes
111 mesh.controlledUpdate();
112
113 if (mesh.changing())
114 {
115 MRF.update();
116
117 if (correctPhi)
118 {
119 // Calculate absolute flux
120 // from the mapped surface velocity
121 phi = mesh.Sf() & Uf();
122
123 #include "correctPhi.H"
124
125 // Make the flux relative to the mesh motion
126 fvc::makeRelative(phi, U);
127 }
128
130 {
131 #include "meshCourantNo.H"
132 }
133 }
134 }
135
136 #include "UEqn.H"
137 #include "TEqn.H"
138
139 // --- Pressure corrector loop
140 while (pimple.correct())
141 {
142 #include "pEqn.H"
143 }
144
145 if (pimple.turbCorr())
146 {
147 laminarTransport.correct();
148 turbulence->correct();
149 }
150 }
151
152 runTime.write();
153
154 runTime.printExecutionTime(Info);
155 }
156
157 Info<< "End\n" << endl;
158
159 return 0;
160}
161
162
163// ************************************************************************* //
Required Classes.
surfaceScalarField & phi
IOMRFZoneList & MRF
pimpleControl & pimple
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
engineTime & runTime
Creates and initialises the velocity field Uf if required.
autoPtr< surfaceVectorField > Uf
compressible::turbulenceModel & turbulence
Calculates and outputs the mean and maximum Courant Numbers.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Execute application functionObjects to post-process existing results.
checkMeshCourantNo
correctPhi
moveMeshOuterCorrectors
singlePhaseTransportModel laminarTransport(U, phi)