OpenFOAM: API Guide
v2112
The open source CFD toolbox
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
-------------------------------------------------------------------------------
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
chtMultiRegionFoam
29
30
Group
31
grpHeatTransferSolvers
32
33
Description
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
"
43
#include "
turbulentFluidThermoModel.H
"
44
#include "
rhoReactionThermo.H
"
45
#include "
CombustionModel.H
"
46
#include "
fixedGradientFvPatchFields.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
60
int
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
"
80
#include "
readSolidTimeControls.H
"
81
#include "compressibleMultiRegionCourantNo.H"
82
#include "solidRegionDiffusionNo.H"
83
#include "
setInitialMultiRegionDeltaT.H
"
84
85
#include "
createCoupledRegions.H
"
86
87
while
(
runTime
.run())
88
{
89
#include "
readTimeControls.H
"
90
#include "
readSolidTimeControls.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
{
103
forAll
(
fluidRegions
, i)
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
114
forAll
(
fluidRegions
, i)
115
{
116
#include "setRegionFluidFields.H"
117
#include "readFluidMultiRegionPIMPLEControls.H"
118
#include "solveFluid.H"
119
}
120
121
forAll
(
solidRegions
, i)
122
{
123
#include "
setRegionSolidFields.H
"
124
#include "
readSolidMultiRegionPIMPLEControls.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
134
forAll
(
fluidRegions
, i)
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
167
forAll
(
fluidRegions
, i)
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
177
forAll
(
solidRegions
, i)
178
{
179
Info
<<
"\nSolving for solid region "
180
<<
solidRegions
[i].name() <<
endl
;
181
#include "
setRegionSolidFields.H
"
182
#include "
readSolidMultiRegionPIMPLEControls.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
192
forAll
(
fluidRegions
, i)
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
// ************************************************************************* //
CombustionModel.H
max
Y[inertIndex] max(0.0)
addCheckCaseOptions.H
Required Classes.
fluidRegions
PtrList< fvMesh > fluidRegions(fluidNames.size())
frozenFlow
bool frozenFlow
Definition:
setRegionFluidFields.H:34
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...
coordinateSystem.H
correctThermos.H
createCoupledRegions.H
fvMatrixAssemblyPtr
autoPtr< fvMatrix< scalar > > fvMatrixAssemblyPtr
Definition:
createCoupledRegions.H:5
coupled
bool coupled(solutionDict.getOrDefault("coupledEnergyField", false))
runTime
engineTime & runTime
Definition:
createEngineTime.H:13
createMeshes.H
solidRegions
PtrList< fvMesh > solidRegions(solidNames.size())
createTimeControls.H
Read the control parameters used by setDeltaT.
createTime.H
fixedGradientFvPatchFields.H
nCorr
const int nCorr
Definition:
readFluidMultiRegionPIMPLEControls.H:3
turbulence
compressible::turbulenceModel & turbulence
Definition:
setRegionFluidFields.H:30
fvCFD.H
fvOptions.H
loopControl.H
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::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition:
hashSets.C:33
Foam::nl
constexpr char nl
The newline '\n' character (0x0a)
Definition:
Ostream.H:53
postProcess.H
Execute application functionObjects to post-process existing results.
pressureControl.H
radiationModel.H
rho
rho
Definition:
readInitialConditions.H:88
readPIMPLEControls.H
nOuterCorr
const int nOuterCorr
Definition:
readPIMPLEControls.H:7
readSolidMultiRegionPIMPLEControls.H
readSolidTimeControls.H
Read the control parameters used in the solid.
readTimeControls.H
Read the control parameters used by setDeltaT.
regionProperties.H
rhoReactionThermo.H
setInitialMultiRegionDeltaT.H
Set the initial timestep for the CHT MultiRegion solver.
setRegionSolidFields.H
setRootCaseLists.H
solidThermo.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition:
stdFoam.H:333
turbulentFluidThermoModel.H
applications
solvers
heatTransfer
chtMultiRegionFoam
chtMultiRegionFoam.C
Generated by
1.9.5
OPENFOAM® is a registered
trademark
of OpenCFD Ltd.