Home
Downloads
Documentation
Installation
User Guide
man-pages
API Documentation
README
Release Notes
Changes
License
Support
SourceForge Project
Main Page
Related Pages
Namespaces
Classes
Files
File List
File Members
All
Classes
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
src
thermophysicalModels
radiation
radiationModel
fvDOM
fvDOM
fvDOM.H
Go to the documentation of this file.
1
/*---------------------------------------------------------------------------*\
2
========= |
3
\\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4
\\ / O peration |
5
\\ / A nd | Copyright (C) 2008-2010 OpenCFD Ltd.
6
\\/ M anipulation |
7
-------------------------------------------------------------------------------
8
License
9
This file is part of OpenFOAM.
10
11
OpenFOAM is free software: you can redistribute it and/or modify it
12
under the terms of the GNU General Public License as published by
13
the Free Software Foundation, either version 3 of the License, or
14
(at your option) any later version.
15
16
OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18
FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19
for more details.
20
21
You should have received a copy of the GNU General Public License
22
along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23
24
Class
25
Foam::radiation::fvDOM
26
27
Description
28
29
Finite Volume Discrete Ordinates Method. Solves the RTE equation for n
30
directions in a participating media, not including scatter.
31
32
Available absorption models:
33
greyMeanAbsoprtionEmission
34
wideBandAbsorptionEmission
35
36
i.e. dictionary
37
fvDOMCoeffs
38
{
39
nPhi 1; // azimuthal angles in PI/2 on X-Y.(from Y to X)
40
nTheta 2; // polar angles in PI (from Z to X-Y plane)
41
convergence 1e-4; // convergence criteria for radiation iteration
42
}
43
44
solverFreq 1; // Number of flow iterations per radiation iteration
45
46
The total number of solid angles is 4*nPhi*nTheta.
47
48
In 1D the direction of the rays is X (nPhi and nTheta are ignored)
49
In 2D the direction of the rays is on X-Y plane (only nPhi is considered)
50
In 3D (nPhi and nTheta are considered)
51
52
SourceFiles
53
fvDOM.C
54
55
\*---------------------------------------------------------------------------*/
56
57
#ifndef radiationModelfvDOM_H
58
#define radiationModelfvDOM_H
59
60
#include <
radiation/radiativeIntensityRay.H
>
61
#include <
radiation/radiationModel.H
>
62
63
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
64
65
namespace
Foam
66
{
67
namespace
radiation
68
{
69
70
/*---------------------------------------------------------------------------*\
71
Class fvDOM Declaration
72
\*---------------------------------------------------------------------------*/
73
74
class
fvDOM
75
:
76
public
radiationModel
77
{
78
// Private data
79
80
//- Incident radiation [W/m2]
81
volScalarField
G_;
82
83
//- Total radiative heat flux [W/m2]
84
volScalarField
Qr_;
85
86
//- Emmited radiative heat flux [W/m2]
87
volScalarField
Qem_;
88
89
//- Incidet radiative heat flux [W/m2]
90
volScalarField
Qin_;
91
92
//- Total absorption coefficient [1/m]
93
volScalarField
a_;
94
95
//- Total emission coefficient [1/m]
96
volScalarField
e_;
97
98
//- Emission contribution [Kg/m/s^3]
99
volScalarField
E_;
100
101
//- Number of solid angles in theta
102
label nTheta_;
103
104
//- Number of solid angles in phi
105
label nPhi_ ;
106
107
//- Total number of rays (1 per direction)
108
label nRay_;
109
110
//- Number of wavelength bands
111
label nLambda_;
112
113
//- Wavelength total absorption coefficient [1/m]
114
PtrList<volScalarField>
aLambda_;
115
116
//- Black body
117
blackBodyEmission
blackBody_;
118
119
//- List of pointers to radiative intensity rays
120
PtrList<radiativeIntensityRay>
IRay_;
121
122
//- Convergence criterion
123
scalar convergence_;
124
125
//- Maximum number of iterations
126
scalar maxIter_;
127
128
129
// Private member functions
130
131
//- Disallow default bitwise copy construct
132
fvDOM
(
const
fvDOM
&);
133
134
//- Disallow default bitwise assignment
135
void
operator=(
const
fvDOM
&);
136
137
//- Update nlack body emission
138
void
updateBlackBodyEmission();
139
140
141
public
:
142
143
//- Runtime type information
144
TypeName
(
"fvDOM"
);
145
146
147
// Constructors
148
149
//- Construct from components
150
fvDOM
(
const
volScalarField
&
T
);
151
152
153
//- Destructor
154
virtual
~fvDOM
();
155
156
157
// Member functions
158
159
// Edit
160
161
//- Solve radiation equation(s)
162
void
calculate
();
163
164
//- Read radiation properties dictionary
165
bool
read
();
166
167
//- Update G and calculate total heat flux on boundary
168
void
updateG
();
169
170
//- Set the rayId and lambdaId from by decomposing an intensity
171
// field name
172
void
setRayIdLambdaId
173
(
174
const
word
&
name
,
175
label& rayId,
176
label& lambdaId
177
)
const
;
178
179
//- Source term component (for power of T^4)
180
virtual
tmp<volScalarField>
Rp
()
const
;
181
182
//- Source term component (constant)
183
virtual
tmp<DimensionedField<scalar, volMesh>
>
Ru
()
const
;
184
185
186
// Access
187
188
//- Ray intensity for rayI
189
inline
const
radiativeIntensityRay
&
IRay
(
const
label rayI)
const
;
190
191
//- Ray intensity for rayI and lambda bandwidth
192
inline
const
volScalarField
&
IRayLambda
193
(
194
const
label rayI,
195
const
label lambdaI
196
)
const
;
197
198
//- Number of angles in theta
199
inline
label
nTheta
()
const
;
200
201
//- Number of angles in phi
202
inline
label
nPhi
()
const
;
203
204
//- Number of rays
205
inline
label
nRay
()
const
;
206
207
//- Number of wavelengths
208
inline
label
nLambda
()
const
;
209
210
//- Const access to total absorption coefficient
211
inline
const
volScalarField
&
a
()
const
;
212
213
//- Const access to wavelength total absorption coefficient
214
inline
const
volScalarField
&
aLambda
(
const
label lambdaI)
const
;
215
216
//- Const access to incident radiation field
217
inline
const
volScalarField
&
G
()
const
;
218
219
//- Const access to total radiative heat flux field
220
inline
const
volScalarField
&
Qr
()
const
;
221
222
//- Const access to incident radiative heat flux field
223
inline
const
volScalarField
&
Qin
()
const
;
224
225
//- Const access to emitted radiative heat flux field
226
inline
const
volScalarField
&
Qem
()
const
;
227
228
//- Const access to black body
229
inline
const
blackBodyEmission
&
blackBody
()
const
;
230
};
231
232
233
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
234
235
#include "
fvDOMI.H
"
236
237
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
238
239
}
// End namespace radiation
240
}
// End namespace Foam
241
242
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
243
244
#endif
245
246
// ************************ vim: set sw=4 sts=4 et: ************************ //