chtMultiRegionFoam.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-2016 OpenFOAM Foundation
9 Copyright (C) 2017-2019,2022 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 chtMultiRegionFoam
29
30Group
31 grpHeatTransferSolvers
32
33Description
34 Transient solver for buoyant, turbulent fluid flow and solid heat
35 conduction with conjugate heat transfer between solid and fluid regions.
36
37 It handles secondary fluid or solid circuits which can be coupled
38 thermally with the main fluid region. i.e radiators, etc.
39
40\*---------------------------------------------------------------------------*/
41
42#include "fvCFD.H"
44#include "rhoReactionThermo.H"
45#include "CombustionModel.H"
47#include "regionProperties.H"
48#include "compressibleCourantNo.H"
49#include "solidRegionDiffNo.H"
50#include "solidThermo.H"
51#include "radiationModel.H"
52#include "fvOptions.H"
53#include "coordinateSystem.H"
54#include "loopControl.H"
55#include "pressureControl.H"
56
57
58// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
59
60int main(int argc, char *argv[])
61{
62 argList::addNote
63 (
64 "Transient solver for buoyant, turbulent fluid flow and solid heat"
65 " conduction with conjugate heat transfer"
66 " between solid and fluid regions."
67 );
68
69 #define NO_CONTROL
70 #define CREATE_MESH createMeshesPostProcess.H
71 #include "postProcess.H"
72
73 #include "addCheckCaseOptions.H"
74 #include "setRootCaseLists.H"
75 #include "createTime.H"
76 #include "createMeshes.H"
77 #include "createFields.H"
78 #include "initContinuityErrs.H"
79 #include "createTimeControls.H"
81 #include "compressibleMultiRegionCourantNo.H"
82 #include "solidRegionDiffusionNo.H"
84
85 #include "createCoupledRegions.H"
86
87 while (runTime.run())
88 {
89 #include "readTimeControls.H"
91 #include "readPIMPLEControls.H"
92
93 #include "compressibleMultiRegionCourantNo.H"
94 #include "solidRegionDiffusionNo.H"
95 #include "setMultiRegionDeltaT.H"
96
97 ++runTime;
98
99 Info<< "Time = " << runTime.timeName() << nl << endl;
100
101 if (nOuterCorr != 1)
102 {
104 {
105 #include "storeOldFluidFields.H"
106 }
107 }
108
109 // --- PIMPLE loop
110 for (int oCorr=0; oCorr<nOuterCorr; ++oCorr)
111 {
112 const bool finalIter = (oCorr == nOuterCorr-1);
113
115 {
116 #include "setRegionFluidFields.H"
117 #include "readFluidMultiRegionPIMPLEControls.H"
118 #include "solveFluid.H"
119 }
120
122 {
123 #include "setRegionSolidFields.H"
125 #include "solveSolid.H"
126 }
127
128 if (coupled)
129 {
130 Info<< "\nSolving energy coupled regions " << endl;
131 fvMatrixAssemblyPtr->solve();
132 #include "correctThermos.H"
133
135 {
136 #include "setRegionFluidFields.H"
137 #include "readFluidMultiRegionPIMPLEControls.H"
138 if (!frozenFlow)
139 {
140 Info<< "\nSolving for fluid region "
141 << fluidRegions[i].name() << endl;
142 // --- PISO loop
143 for (int corr=0; corr<nCorr; corr++)
144 {
145 #include "pEqn.H"
146 }
147 turbulence.correct();
148 }
149
150 rho = thermo.rho();
151 Info<< "Min/max T:" << min(thermo.T()).value() << ' '
152 << max(thermo.T()).value() << endl;
153 }
154
155 fvMatrixAssemblyPtr->clear();
156 }
157
158 // Additional loops for energy solution only
159 if (!oCorr && nOuterCorr > 1)
160 {
161 loopControl looping(runTime, pimple, "energyCoupling");
162
163 while (looping.loop())
164 {
165 Info<< nl << looping << nl;
166
168 {
169 Info<< "\nSolving for fluid region "
170 << fluidRegions[i].name() << endl;
171 #include "setRegionFluidFields.H"
172 #include "readFluidMultiRegionPIMPLEControls.H"
173 frozenFlow = true;
174 #include "solveFluid.H"
175 }
176
178 {
179 Info<< "\nSolving for solid region "
180 << solidRegions[i].name() << endl;
181 #include "setRegionSolidFields.H"
183 #include "solveSolid.H"
184 }
185
186 if (coupled)
187 {
188 Info<< "\nSolving energy coupled regions " << endl;
189 fvMatrixAssemblyPtr->solve();
190 #include "correctThermos.H"
191
193 {
194 #include "setRegionFluidFields.H"
195 rho = thermo.rho();
196 }
197
198 fvMatrixAssemblyPtr->clear();
199 }
200 }
201 }
202 }
203
204 runTime.write();
205
206 runTime.printExecutionTime(Info);
207 }
208
209 Info<< "End\n" << endl;
210
211 return 0;
212}
213
214
215// ************************************************************************* //
Y[inertIndex] max(0.0)
Required Classes.
PtrList< fvMesh > fluidRegions(fluidNames.size())
bool frozenFlow
pimpleControl & pimple
Basic thermodynamics type based on the use of fitting functions for cp, h, s obtained from the templa...
autoPtr< fvMatrix< scalar > > fvMatrixAssemblyPtr
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
engineTime & runTime
PtrList< fvMesh > solidRegions(solidNames.size())
Read the control parameters used by setDeltaT.
compressible::turbulenceModel & turbulence
messageStream Info
Information stream (stdout output on master, null elsewhere)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
constexpr char nl
The newline '\n' character (0x0a)
Definition: Ostream.H:53
Execute application functionObjects to post-process existing results.
const int nOuterCorr
Read the control parameters used in the solid.
Read the control parameters used by setDeltaT.
Set the initial timestep for the CHT MultiRegion solver.
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333