FreeFOAM The Cross-Platform CFD Toolkit
radiativeIntensityRay.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::radiativeIntensityRay
26 
27 Description
28  Radiation intensity for a ray in a given direction
29 
30 SourceFiles
31  radiativeIntensityRay.C
32 
33 \*---------------------------------------------------------------------------*/
34 
35 #ifndef radiativeIntensityRay_H
36 #define radiativeIntensityRay_H
37 
40 
41 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42 
43 namespace Foam
44 {
45 namespace radiation
46 {
47 
48 // Forward declaration of classes
49 class fvDOM;
50 
51 /*---------------------------------------------------------------------------*\
52  Class radiativeIntensityRay Declaration
53 \*---------------------------------------------------------------------------*/
54 
56 {
57 public:
58 
59  static const word intensityPrefix;
60 
61 
62 private:
63 
64  // Private data
65 
66  //- Refence to the owner fvDOM object
67  const fvDOM& dom_;
68 
69  //- Reference to the mesh
70  const fvMesh& mesh_;
71 
72  //- Absorption/emission model
73  const absorptionEmissionModel& absorptionEmission_;
74 
75  //- Black body
76  const blackBodyEmission& blackBody_;
77 
78  //- Total radiative intensity / [W/m2]
79  volScalarField I_;
80 
81  //- Total radiative heat flux on boundary
82  volScalarField Qr_;
83 
84  //- Incident radiative heat flux on boundary
85  volScalarField Qin_;
86 
87  //- Emitted radiative heat flux on boundary
88  volScalarField Qem_;
89 
90  //- Direction
91  vector d_;
92 
93  //- Average direction vector inside the solid angle
94  vector dAve_;
95 
96  //- Theta angle
97  scalar theta_;
98 
99  //- Phi angle
100  scalar phi_;
101 
102  //- Solid angle
103  scalar omega_;
104 
105  //- Number of wavelengths/bands
106  label nLambda_;
107 
108  //- List of pointers to radiative intensity fields for given wavelengths
109  PtrList<volScalarField> ILambda_;
110 
111 
112  // Private member functions
113 
114  //- Disallow default bitwise copy construct
116 
117  //- Disallow default bitwise assignment
118  void operator=(const radiativeIntensityRay&);
119 
120 
121 public:
122 
123  // Constructors
124 
125  //- Construct form components
127  (
128  const fvDOM& dom,
129  const fvMesh& mesh,
130  const scalar phi,
131  const scalar theta,
132  const scalar deltaPhi,
133  const scalar deltaTheta,
134  const label lambda,
135  const absorptionEmissionModel& absEmmModel_,
136  const blackBodyEmission& blackBody,
137  const label rayId
138  );
139 
140 
141  // Destructor
143 
144 
145  // Member functions
146 
147  // Edit
148 
149  //- Update radiative intensity on i direction
150  scalar correct();
151 
152  //- Initialise the ray in i direction
153  void init
154  (
155  const scalar phi,
156  const scalar theta,
157  const scalar deltaPhi,
158  const scalar deltaTheta,
159  const scalar lambda
160  );
161 
162  //- Add radiative intensities from all the bands
163  void addIntensity();
164 
165 
166  // Access
167 
168  //- Return intensity
169  inline const volScalarField& I() const;
170 
171  //- Return const access to the boundary heat flux
172  inline const volScalarField& Qr() const;
173 
174  //- Return non-const access to the boundary heat flux
175  inline volScalarField& Qr();
176 
177  //- Return non-const access to the boundary incident heat flux
178  inline volScalarField& Qin();
179 
180  //- Return non-const access to the boundary emmited heat flux
181  inline volScalarField& Qem();
182 
183  //- Return const access to the boundary incident heat flux
184  inline const volScalarField& Qin() const;
185 
186  //- Return const access to the boundary emmited heat flux
187  inline const volScalarField& Qem() const;
188 
189  //- Return direction
190  inline const vector& d() const;
191 
192  //- Return the average vector inside the solid angle
193  inline const vector& dAve() const;
194 
195  //- Return the number of bands
196  inline scalar nLambda() const;
197 
198  //- Return the phi angle
199  inline scalar phi() const;
200 
201  //- Return the theta angle
202  inline scalar theta() const;
203 
204  //- Return the solid angle
205  inline scalar omega() const;
206 
207  //- Return the radiative intensity for a given wavelength
208  inline const volScalarField& ILambda(const label lambdaI) const;
209 };
210 
211 
212 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
213 
214 } // End namespace radiation
215 } // End namespace Foam
216 
217 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
218 
219 #include "radiativeIntensityRayI.H"
220 
221 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
222 
223 #endif
224 
225 // ************************ vim: set sw=4 sts=4 et: ************************ //