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) 2011-2017 OpenFOAM Foundation
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 
39 #include "blackBodyEmission.H"
40 
41 
42 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
43 
44 namespace Foam
45 {
46 namespace radiation
47 {
48 
49 // Forward declaration of classes
50 class fvDOM;
51 
52 /*---------------------------------------------------------------------------*\
53  Class radiativeIntensityRay Declaration
54 \*---------------------------------------------------------------------------*/
55 
57 {
58 public:
59 
60  static const word intensityPrefix;
61 
62 
63 private:
64 
65  // Private data
66 
67  //- Refence to the owner fvDOM object
68  const fvDOM& dom_;
69 
70  //- Reference to the mesh
71  const fvMesh& mesh_;
72 
73  //- Absorption/emission model
74  const absorptionEmissionModel& absorptionEmission_;
75 
76  //- Black body
77  const blackBodyEmission& blackBody_;
78 
79  //- Total radiative intensity / [W/m2]
80  volScalarField I_;
81 
82  //- Total radiative heat flux on boundary
83  volScalarField qr_;
84 
85  //- Incident radiative heat flux on boundary
86  volScalarField qin_;
87 
88  //- Emitted radiative heat flux on boundary
89  volScalarField qem_;
90 
91  //- Direction
92  vector d_;
93 
94  //- Average direction vector inside the solid angle
95  vector dAve_;
96 
97  //- Theta angle
98  scalar theta_;
99 
100  //- Phi angle
101  scalar phi_;
102 
103  //- Solid angle
104  scalar omega_;
105 
106  //- Number of wavelengths/bands
107  label nLambda_;
108 
109  //- List of pointers to radiative intensity fields for given wavelengths
110  PtrList<volScalarField> ILambda_;
111 
112  //- Global ray id - incremented in constructor
113  static label rayId;
114 
115  //- My ray Id
116  label myRayId_;
117 
118 
119  // Private Member Functions
120 
121  //- Disallow default bitwise copy construct
123 
124  //- Disallow default bitwise assignment
125  void operator=(const radiativeIntensityRay&);
126 
127 
128 public:
129 
130  // Constructors
131 
132  //- Construct form components
134  (
135  const fvDOM& dom,
136  const fvMesh& mesh,
137  const scalar phi,
138  const scalar theta,
139  const scalar deltaPhi,
140  const scalar deltaTheta,
141  const label lambda,
142  const absorptionEmissionModel& absEmmModel_,
143  const blackBodyEmission& blackBody,
144  const label rayId
145  );
146 
147 
148  //- Destructor
150 
151 
152  // Member functions
153 
154  // Edit
155 
156  //- Update radiative intensity on i direction
157  scalar correct();
158 
159  //- Initialise the ray in i direction
160  void init
161  (
162  const scalar phi,
163  const scalar theta,
164  const scalar deltaPhi,
165  const scalar deltaTheta,
166  const scalar lambda
167  );
168 
169  //- Add radiative intensities from all the bands
170  void addIntensity();
171 
172 
173  // Access
174 
175  //- Return intensity
176  inline const volScalarField& I() const;
177 
178  //- Return const access to the boundary heat flux
179  inline const volScalarField& qr() const;
180 
181  //- Return non-const access to the boundary heat flux
182  inline volScalarField& qr();
183 
184  //- Return non-const access to the boundary incident heat flux
185  inline volScalarField& qin();
186 
187  //- Return non-const access to the boundary emmited heat flux
188  inline volScalarField& qem();
189 
190  //- Return const access to the boundary incident heat flux
191  inline const volScalarField& qin() const;
192 
193  //- Return const access to the boundary emmited heat flux
194  inline const volScalarField& qem() const;
195 
196  //- Return direction
197  inline const vector& d() const;
198 
199  //- Return the average vector inside the solid angle
200  inline const vector& dAve() const;
201 
202  //- Return the number of bands
203  inline scalar nLambda() const;
204 
205  //- Return the phi angle
206  inline scalar phi() const;
207 
208  //- Return the theta angle
209  inline scalar theta() const;
210 
211  //- Return the solid angle
212  inline scalar omega() const;
213 
214  //- Return the radiative intensity for a given wavelength
215  inline const volScalarField& ILambda(const label lambdaI) const;
216 
217 };
218 
219 
220 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
221 
222 } // End namespace radiation
223 } // End namespace Foam
224 
225 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
226 
227 #include "radiativeIntensityRayI.H"
228 
229 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
230 
231 #endif
232 
233 // ************************************************************************* //
Class black body emission.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void init(const scalar phi, const scalar theta, const scalar deltaPhi, const scalar deltaTheta, const scalar lambda)
Initialise the ray in i direction.
volScalarField & qin()
Return non-const access to the boundary incident heat flux.
Model to supply absorption and emission coefficients for radiation modelling.
scalar correct()
Update radiative intensity on i direction.
dynamicFvMesh & mesh
scalar phi() const
Return the phi angle.
A class for handling words, derived from string.
Definition: word.H:59
const volScalarField & I() const
Return intensity.
scalar omega() const
Return the solid angle.
void addIntensity()
Add radiative intensities from all the bands.
scalar nLambda() const
Return the number of bands.
Radiation intensity for a ray in a given direction.
const vector & d() const
Return direction.
const volScalarField & ILambda(const label lambdaI) const
Return the radiative intensity for a given wavelength.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:63
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
dimensionedScalar lambda(laminarTransport.lookup("lambda"))
volScalarField & qem()
Return non-const access to the boundary emmited heat flux.
const vector & dAve() const
Return the average vector inside the solid angle.
Finite Volume Discrete Ordinates Method. Solves the RTE equation for n directions in a participating ...
Definition: fvDOM.H:80
const volScalarField & qr() const
Return const access to the boundary heat flux.
scalar theta() const
Return the theta angle.
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
Namespace for OpenFOAM.