rhoPimpleFoam.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-------------------------------------------------------------------------------
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 rhoPimpleFoam
29
30Group
31 grpCompressibleSolvers
32
33Description
34 Transient solver for turbulent flow of compressible fluids for HVAC and
35 similar applications, with optional mesh motion and mesh topology changes.
36
37 Uses the flexible PIMPLE (PISO-SIMPLE) solution for time-resolved and
38 pseudo-transient simulations.
39
40Note
41 The motion frequency of this solver can be influenced by the presence
42 of "updateControl" and "updateInterval" in the dynamicMeshDict.
43
44\*---------------------------------------------------------------------------*/
45
46#include "fvCFD.H"
47#include "dynamicFvMesh.H"
48#include "fluidThermo.H"
50#include "bound.H"
51#include "pimpleControl.H"
52#include "pressureControl.H"
53#include "CorrectPhi.H"
54#include "fvOptions.H"
55#include "localEulerDdtScheme.H"
56#include "fvcSmooth.H"
57
58// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
59
60int main(int argc, char *argv[])
61{
62 argList::addNote
63 (
64 "Transient solver for compressible turbulent flow.\n"
65 "With optional mesh motion and mesh topology changes."
66 );
67
68 #include "postProcess.H"
69
70 #include "addCheckCaseOptions.H"
71 #include "setRootCaseLists.H"
72 #include "createTime.H"
73 #include "createDynamicFvMesh.H"
74 #include "createDyMControls.H"
75 #include "initContinuityErrs.H"
76 #include "createFields.H"
77 #include "createFieldRefs.H"
78 #include "createRhoUfIfPresent.H"
79
80 turbulence->validate();
81
82 if (!LTS)
83 {
84 #include "compressibleCourantNo.H"
85 #include "setInitialDeltaT.H"
86 }
87
88 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
89
90 Info<< "\nStarting time loop\n" << endl;
91
92 while (runTime.run())
93 {
94 #include "readDyMControls.H"
95
96 // Store divrhoU from the previous mesh so that it can be mapped
97 // and used in correctPhi to ensure the corrected phi has the
98 // same divergence
99 autoPtr<volScalarField> divrhoU;
100 if (correctPhi)
101 {
102 divrhoU.reset
103 (
104 new volScalarField
105 (
106 "divrhoU",
107 fvc::div(fvc::absolute(phi, rho, U))
108 )
109 );
110 }
111
112 if (LTS)
113 {
114 #include "setRDeltaT.H"
115 }
116 else
117 {
118 #include "compressibleCourantNo.H"
119 #include "setDeltaT.H"
120 }
121
122 ++runTime;
123
124 Info<< "Time = " << runTime.timeName() << nl << endl;
125
126 // --- Pressure-velocity PIMPLE corrector loop
127 while (pimple.loop())
128 {
129 if (pimple.firstIter() || moveMeshOuterCorrectors)
130 {
131 // Store momentum to set rhoUf for introduced faces.
132 autoPtr<volVectorField> rhoU;
133 if (rhoUf.valid())
134 {
135 rhoU.reset(new volVectorField("rhoU", rho*U));
136 }
137
138 // Do any mesh changes
139 mesh.controlledUpdate();
140
141 if (mesh.changing())
142 {
143 MRF.update();
144
145 if (correctPhi)
146 {
147 // Calculate absolute flux
148 // from the mapped surface velocity
149 phi = mesh.Sf() & rhoUf();
150
151 #include "correctPhi.H"
152
153 // Make the fluxes relative to the mesh-motion
154 fvc::makeRelative(phi, rho, U);
155 }
156
158 {
159 #include "meshCourantNo.H"
160 }
161 }
162 }
163
164 if (pimple.firstIter() && !pimple.SIMPLErho())
165 {
166 #include "rhoEqn.H"
167 }
168
169 #include "UEqn.H"
170 #include "EEqn.H"
171
172 // --- Pressure corrector loop
173 while (pimple.correct())
174 {
175 if (pimple.consistent())
176 {
177 #include "pcEqn.H"
178 }
179 else
180 {
181 #include "pEqn.H"
182 }
183 }
184
185 if (pimple.turbCorr())
186 {
187 turbulence->correct();
188 }
189 }
190
191 rho = thermo.rho();
192
193 runTime.write();
194
195 runTime.printExecutionTime(Info);
196 }
197
198 Info<< "End\n" << endl;
199
200 return 0;
201}
202
203
204// ************************************************************************* //
Required Classes.
Bound the given scalar field if it has gone unbounded.
surfaceScalarField & phi
IOMRFZoneList & MRF
pimpleControl & pimple
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
engineTime & runTime
bool LTS
Definition: createRDeltaT.H:1
Creates and initialises the velocity field rhoUf if required.
autoPtr< surfaceVectorField > rhoUf
compressible::turbulenceModel & turbulence
Provides functions smooth spread and sweep which use the FaceCellWave algorithm to smooth and redistr...
Calculates and outputs the mean and maximum Courant Numbers.
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:83
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