All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
fvDOM.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 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::radiationModels::fvDOM
26 
27 Description
28  Finite Volume Discrete Ordinates Method. Solves the RTE equation for n
29  directions in a participating media, not including scatter.
30 
31  In 1-D the ray directions are bound to one of the X, Y or Z directions. The
32  total number of solid angles is 2. nPhi and nTheta are ignored.
33 
34  In 2-D the ray directions are within one of the X-Y, X-Z or Y-Z planes. The
35  total number of solid angles is 4*nPhi. nTheta is ignored.
36 
37  In 3D the rays span all directions. The total number of solid angles is
38  4*nPhi*nTheta.
39 
40 Usage
41  \verbatim
42  fvDOMCoeffs
43  {
44  nPhi 4; // azimuthal angles in PI/2 on X-Y (from Y to X)
45  nTheta 0; // polar angles in PI (from Z to X-Y plane)
46  convergence 1e-3; // convergence criteria for radiation iteration
47  maxIter 4; // maximum number of iterations
48  }
49  solverFreq 1; // Number of flow iterations per radiation iteration
50  \endverbatim
51 
52 SourceFiles
53  fvDOM.C
54 
55 \*---------------------------------------------------------------------------*/
56 
57 #ifndef radiation_fvDOM_H
58 #define radiation_fvDOM_H
59 
60 #include "radiativeIntensityRay.H"
61 #include "radiationModel.H"
62 #include "fvMatrices.H"
63 
64 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
65 
66 namespace Foam
67 {
68 namespace radiationModels
69 {
70 
71 /*---------------------------------------------------------------------------*\
72  Class fvDOM Declaration
73 \*---------------------------------------------------------------------------*/
74 
75 class fvDOM
76 :
77  public radiationModel
78 {
79  // Private Data
80 
81 
82  //- Incident radiation [W/m^2]
83  volScalarField G_;
84 
85  //- Total radiative heat flux [W/m^2]
86  volScalarField qr_;
87 
88  //- Emitted radiative heat flux [W/m^2]
89  volScalarField qem_;
90 
91  //- Incident radiative heat flux [W/m^2]
92  volScalarField qin_;
93 
94  //- Total absorption coefficient [1/m]
96 
97  //- Number of solid angles in theta
98  label nTheta_;
99 
100  //- Number of solid angles in phi
101  label nPhi_ ;
102 
103  //- Total number of rays (1 per direction)
104  label nRay_;
105 
106  //- Number of wavelength bands
107  label nLambda_;
109  //- Wavelength total absorption coefficient [1/m]
110  PtrList<volScalarField> aLambda_;
111 
112  //- Black body
113  blackBodyEmission blackBody_;
115  //- List of pointers to radiative intensity rays
117 
118  //- Convergence tolerance
119  scalar tolerance_;
120 
121  //- Maximum number of iterations
122  label maxIter_;
123 
124  //- Maximum omega weight
125  scalar omegaMax_;
126 
127 
128  // Private Member Functions
129 
130  //- Initialise
131  void initialise();
132 
133  //- Update black body emission
134  void updateBlackBodyEmission();
135 
136 
137 public:
138 
139  //- Runtime type information
140  TypeName("fvDOM");
141 
142 
143  // Constructors
144 
145  //- Construct from components
146  fvDOM(const volScalarField& T);
147 
148  //- Construct from components
149  fvDOM(const dictionary& dict, const volScalarField& T);
150 
151  //- Disallow default bitwise copy construction
152  fvDOM(const fvDOM&) = delete;
153 
154 
155  //- Destructor
156  virtual ~fvDOM();
157 
158 
159  // Member Functions
160 
161  // Edit
162 
163  //- Solve radiation equation(s)
164  void calculate();
165 
166  //- Read radiation properties dictionary
167  bool read();
168 
169  //- Update G and calculate total heat flux on boundary
170  void updateG();
171 
172  //- Set the rayId and lambdaId from by decomposing an intensity
173  // field name
174  void setRayIdLambdaId
175  (
176  const word& name,
177  label& rayId,
178  label& lambdaId
179  ) const;
180 
181  //- Source term component (for power of T^4)
182  virtual tmp<volScalarField> Rp() const;
183 
184  //- Source term component (constant)
185  virtual tmp<volScalarField::Internal> Ru() const;
186 
187 
188  // Access
189 
190  //- Ray intensity for rayI
191  inline const radiativeIntensityRay& IRay(const label rayI) const;
192 
193  //- Ray intensity for rayI and lambda bandwidth
194  inline const volScalarField& IRayLambda
195  (
196  const label rayI,
197  const label lambdaI
198  ) const;
199 
200  //- Number of angles in theta
201  inline label nTheta() const;
202 
203  //- Number of angles in phi
204  inline label nPhi() const;
205 
206  //- Number of rays
207  inline label nRay() const;
208 
209  //- Number of wavelengths
210  inline label nLambda() const;
211 
212  //- Const access to total absorption coefficient
213  inline const volScalarField& a() const;
214 
215  //- Const access to wavelength total absorption coefficient
216  inline const volScalarField& aLambda(const label lambdaI) const;
217 
218  //- Const access to incident radiation field
219  inline const volScalarField& G() const;
220 
221  //- Const access to total radiative heat flux field
222  inline const volScalarField& qr() const;
223 
224  //- Const access to incident radiative heat flux field
225  inline const volScalarField& qin() const;
226 
227  //- Const access to emitted radiative heat flux field
228  inline const volScalarField& qem() const;
229 
230  //- Const access to black body
231  inline const blackBodyEmission& blackBody() const;
232 
233  //- Return omegaMax
234  inline scalar omegaMax() const;
235 
236 
237  // Member Operators
238 
239  //- Disallow default bitwise assignment
240  void operator=(const fvDOM&) = delete;
241 };
242 
243 
244 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
245 
246 #include "fvDOMI.H"
247 
248 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
249 
250 } // End namespace radiationModels
251 } // End namespace Foam
252 
253 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
254 
255 #endif
256 
257 // ************************************************************************* //
label nLambda() const
Number of wavelengths.
Definition: fvDOM.H:63
const volScalarField & qem() const
Const access to emitted radiative heat flux field.
Definition: fvDOM.H:101
Finite Volume Discrete Ordinates Method. Solves the RTE equation for n directions in a participating ...
Definition: fvDOM.H:74
dictionary dict
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 calculate()
Solve radiation equation(s)
Definition: fvDOM.C:409
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:158
const volScalarField & qin() const
Const access to incident radiative heat flux field.
Definition: fvDOM.H:95
const volScalarField & a() const
Const access to total absorption coefficient.
Definition: fvDOM.H:69
TypeName("fvDOM")
Runtime type information.
const volScalarField & IRayLambda(const label rayI, const label lambdaI) const
Ray intensity for rayI and lambda bandwidth.
Definition: fvDOM.H:36
virtual tmp< volScalarField::Internal > Ru() const
Source term component (constant)
Definition: fvDOM.C:480
label nTheta() const
Number of angles in theta.
Definition: fvDOM.H:45
const volScalarField & G() const
Const access to incident radiation field.
Definition: fvDOM.H:84
fvDOM(const volScalarField &T)
Construct from components.
Definition: fvDOM.C:201
const radiativeIntensityRay & IRay(const label rayI) const
Ray intensity for rayI.
Definition: fvDOM.H:28
const volScalarField & aLambda(const label lambdaI) const
Const access to wavelength total absorption coefficient.
Definition: fvDOM.H:76
A class for handling words, derived from string.
Definition: word.H:59
Top level model for radiation modelling.
bool read()
Read radiation properties dictionary.
Definition: fvDOM.C:389
label nPhi() const
Number of angles in phi.
Definition: fvDOM.H:51
const word & name() const
Name function is needed to disambiguate those inherited.
void setRayIdLambdaId(const word &name, label &rayId, label &lambdaId) const
Set the rayId and lambdaId from by decomposing an intensity.
Definition: fvDOM.C:552
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
const volScalarField & qr() const
Const access to total radiative heat flux field.
Definition: fvDOM.H:90
Radiation intensity for a ray in a given direction.
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
void operator=(const fvDOM &)=delete
Disallow default bitwise assignment.
scalar omegaMax() const
Return omegaMax.
Definition: fvDOM.H:114
A special matrix type and solver, designed for finite volume solutions of scalar equations.
label nRay() const
Number of rays.
Definition: fvDOM.H:57
virtual ~fvDOM()
Destructor.
Definition: fvDOM.C:383
A class for managing temporary objects.
Definition: PtrList.H:53
virtual tmp< volScalarField > Rp() const
Source term component (for power of T^4)
Definition: fvDOM.C:446
const blackBodyEmission & blackBody() const
Const access to black body.
Definition: fvDOM.H:108
void updateG()
Update G and calculate total heat flux on boundary.
Definition: fvDOM.C:533
Namespace for OpenFOAM.