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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2017 OpenFOAM Foundation
9  Copyright (C) 2019 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 Class
28  Foam::radiation::fvDOM
29 
30 Group
31  grpRadiationModels
32 
33 Description
34 
35  Finite Volume Discrete Ordinates Method. Solves the RTE equation for n
36  directions in a participating media, not including scatter and reflective
37  walls.
38 
39  Available absorption models:
40 
41  - constantAbsorptionEmission
42  - greyMeanAbsoprtionEmission
43  - wideBandAbsorptionEmission
44  - multiBandAbsorptionEmission
45 
46  This model can handle non-grey participating media using
47  multiBandAbsorptionEmission model. Accordingly the BC for rays should
48  be wideBandDiffussive type
49 
50 Usage
51  \verbatim
52  fvDOMCoeffs
53  {
54  nPhi 4; // azimuthal angles in PI/2 on X-Y.
55  //(from Y to X)
56  nTheta 0; // polar angles in PI (from Z to X-Y plane)
57  tolerance 1e-3; // convergence tolerance for radiation
58  // iteration
59  maxIter 4; // maximum number of iterations
60  meshOrientation (1 1 1); //Mesh orientation used for 2D and 1D
61 
62  useSolarLoad false;
63  useExternalBeam true;
64  spectralDistribution (2 1);
65  }
66 
67  solverFreq 1; // Number of flow iterations per radiation iteration
68  \endverbatim
69 
70  The total number of solid angles is 4*nPhi*nTheta in 3-D.
71 
72  Operating modes:
73  - 1-D:
74  - ray directions are on X, Y or Z
75  - \c nPhi and \c nTheta entries are ignored
76  - \c meshOrientation vector can be used for any other 1-D direction.
77  - 2-D:
78  - ray directions are on X-Y, X-Z or Y-Z planes
79  - only the \c nPhi entry is considered
80  - \c meshOrientation vector can be used for non-aligned planes
81  specifying the plane normal vector.
82  - 3-D:
83  - rays geberated in 3-D using the \c nPhi and \c nTheta entries
84  - \c meshOrientation vector is not applicable.
85 
86  useSolarLoad calculates the primary and diffusive Sun fluxes on walls in
87  addition to the RTE equations
88 
89  useExternalBeam add an external collimated beam to the domain. This option
90  is not available if useSolarLoad is true.
91 
92  spectralDistribution is the energy spectral distribution of the collimated
93  external beam.
94 
95 
96 SourceFiles
97  fvDOM.C
98 
99 \*---------------------------------------------------------------------------*/
100 
101 #ifndef radiation_fvDOM_H
102 #define radiation_fvDOM_H
103 
104 #include "radiativeIntensityRay.H"
105 #include "radiationModel.H"
106 #include "fvMatrices.H"
107 #include "solarLoad.H"
108 
109 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
110 
111 namespace Foam
112 {
113 namespace radiation
114 {
115 
116 /*---------------------------------------------------------------------------*\
117  Class fvDOM Declaration
118 \*---------------------------------------------------------------------------*/
119 
120 class fvDOM
121 :
122  public radiationModel
123 {
124  // Private data
125 
126 
127  //- Incident radiation [W/m2]
128  volScalarField G_;
129 
130  //- Total radiative heat flux [W/m2]
131  volScalarField qr_;
132 
133  //- Emitted radiative heat flux [W/m2]
134  volScalarField qem_;
135 
136  //- Incident radiative heat flux [W/m2]
137  volScalarField qin_;
138 
139  //- Total absorption coefficient [1/m]
140  volScalarField a_;
141 
142  //- Number of solid angles in theta
143  label nTheta_;
144 
145  //- Number of solid angles in phi
146  label nPhi_ ;
147 
148  //- Total number of rays (1 per direction)
149  label nRay_;
150 
151  //- Number of wavelength bands
152  label nLambda_;
153 
154  //- Wavelength total absorption coefficient [1/m]
155  PtrList<volScalarField> aLambda_;
156 
157  //- Black body
158  blackBodyEmission blackBody_;
159 
160  //- List of pointers to radiative intensity rays
162 
163  //- Convergence tolerance
164  scalar tolerance_;
165 
166  //- Maximum number of iterations
167  label maxIter_;
168 
169  //- Maximum omega weight
170  scalar omegaMax_;
171 
172  //- Use Solar Load model
173  bool useSolarLoad_;
174 
175  //- Solar load radiation model
176  autoPtr<solarLoad> solarLoad_;
177 
178  //- Mesh orientation vector
179  vector meshOrientation_;
180 
181  //- Use external parallel irradiation beam
182  bool useExternalBeam_;
183 
184  //- Spectral energy distribution for the external beam
185  scalarList spectralDistribution_;
186 
187  //- Solar calculator
188  autoPtr<solarCalculator> solarCalculator_;
189 
190  //- Update Sun position index
191  label updateTimeIndex_;
192 
193 
194  // Private Member Functions
195 
196  //- Initialise
197  void initialise();
198 
199  //- No copy construct
200  fvDOM(const fvDOM&) = delete;
201 
202  //- No copy assignment
203  void operator=(const fvDOM&) = delete;
204 
205  //- Update black body emission
206  void updateBlackBodyEmission();
207 
208 
209 public:
210 
211  //- Runtime type information
212  TypeName("fvDOM");
213 
214 
215  // Constructors
216 
217  //- Construct from components
218  fvDOM(const volScalarField& T);
219 
220  //- Construct from components
221  fvDOM(const dictionary& dict, const volScalarField& T);
222 
223 
224  //- Destructor
225  virtual ~fvDOM();
226 
227 
228  // Member functions
229 
230  // Edit
231 
232  //- Solve radiation equation(s)
233  void calculate();
234 
235  //- Read radiation properties dictionary
236  bool read();
237 
238  //- Update G and calculate total heat flux on boundary
239  void updateG();
240 
241  //- Set the rayId and lambdaId from by decomposing an intensity
242  // field name
243  void setRayIdLambdaId
244  (
245  const word& name,
246  label& rayId,
247  label& lambdaId
248  ) const;
249 
250 
251  //- Rotate rays according to Sun direction
252  void updateRaysDir();
253 
254  //- Rotate rays spheric equator to sunDir
255  void rotateInitialRays(const vector& sunDir);
256 
257  //- Align closest ray to sunDir
258  void alignClosestRayToSun(const vector& sunDir);
259 
260  //- Source term component (for power of T^4)
261  virtual tmp<volScalarField> Rp() const;
262 
263  //- Source term component (constant)
264  virtual tmp<volScalarField::Internal> Ru() const;
265 
266 
267  // Access
268 
269  //- Solar calculator
270  const solarCalculator& solarCalc() const;
271 
272  //- Ray intensity for rayI
273  inline const radiativeIntensityRay& IRay(const label rayI) const;
274 
275  //- Ray intensity for rayI and lambda bandwidth
276  inline const volScalarField& IRayLambda
277  (
278  const label rayI,
279  const label lambdaI
280  ) const;
281 
282  //- Number of angles in theta
283  inline label nTheta() const;
284 
285  //- Number of angles in phi
286  inline label nPhi() const;
287 
288  //- Number of rays
289  inline label nRay() const;
290 
291  //- Number of wavelengths
292  inline label nLambda() const;
293 
294  //- Number of bands
295  inline label nBands() const;
296 
297  //- Const access to total absorption coefficient
298  inline const volScalarField& a() const;
299 
300  //- Const access to wavelength total absorption coefficient
301  inline const volScalarField& aLambda(const label lambdaI) const;
302 
303  //- Const access to incident radiation field
304  inline const volScalarField& G() const;
305 
306  //- Const access to total radiative heat flux field
307  inline const volScalarField& qr() const;
308 
309  //- Const access to incident radiative heat flux field
310  inline const volScalarField& qin() const;
311 
312  //- Const access to emitted radiative heat flux field
313  inline const volScalarField& qem() const;
314 
315  //- Const access to black body
316  inline const blackBodyEmission& blackBody() const;
317 
318  //- Return omegaMax
319  inline scalar omegaMax() const;
320 
321  //- Return meshOrientation
322  inline vector meshOrientation() const;
323 
324  //- Use solar load
325  inline bool useSolarLoad() const;
326 
327  //- Use external beam
328  inline bool useExternalBeam() const;
329 
330  //- Energy spectral distribution for external beam
331  inline const scalarList& spectralDistribution() const;
332 };
333 
334 
335 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
336 
337 #include "fvDOMI.H"
338 
339 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
340 
341 } // End namespace radiation
342 } // End namespace Foam
343 
344 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
345 
346 #endif
347 
348 // ************************************************************************* //
Foam::radiation::fvDOM::rotateInitialRays
void rotateInitialRays(const vector &sunDir)
Rotate rays spheric equator to sunDir.
Definition: fvDOM.C:54
Foam::radiation::fvDOM::nRay
label nRay() const
Number of rays.
Definition: fvDOMI.H:59
solarLoad.H
Foam::radiation::fvDOM::nPhi
label nPhi() const
Number of angles in phi.
Definition: fvDOMI.H:53
fvDOMI.H
Foam::radiation::fvDOM::a
const volScalarField & a() const
Const access to total absorption coefficient.
Definition: fvDOMI.H:76
Foam::radiation::fvDOM::~fvDOM
virtual ~fvDOM()
Destructor.
Definition: fvDOM.C:545
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:62
Foam::radiation::fvDOM::omegaMax
scalar omegaMax() const
Return omegaMax.
Definition: fvDOMI.H:121
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
radiativeIntensityRay.H
Foam::radiation::fvDOM::updateG
void updateG()
Update G and calculate total heat flux on boundary.
Definition: fvDOM.C:735
Foam::radiation::fvDOM::alignClosestRayToSun
void alignClosestRayToSun(const vector &sunDir)
Align closest ray to sunDir.
Definition: fvDOM.C:68
Foam::radiation::radiativeIntensityRay
Radiation intensity for a ray in a given direction.
Definition: radiativeIntensityRay.H:57
Foam::radiation::fvDOM::IRay
const radiativeIntensityRay & IRay(const label rayI) const
Ray intensity for rayI.
Definition: fvDOMI.H:30
Foam::radiation::fvDOM::qem
const volScalarField & qem() const
Const access to emitted radiative heat flux field.
Definition: fvDOMI.H:108
Foam::radiation::fvDOM::qin
const volScalarField & qin() const
Const access to incident radiative heat flux field.
Definition: fvDOMI.H:102
Foam::baseIOdictionary::name
const word & name() const
Definition: baseIOdictionary.C:82
Foam::radiation::fvDOM::nLambda
label nLambda() const
Number of wavelengths.
Definition: fvDOMI.H:65
fvMatrices.H
A special matrix type and solver, designed for finite volume solutions of scalar equations.
Foam::radiation::fvDOM::IRayLambda
const volScalarField & IRayLambda(const label rayI, const label lambdaI) const
Ray intensity for rayI and lambda bandwidth.
Definition: fvDOMI.H:38
Foam::radiation::fvDOM::useExternalBeam
bool useExternalBeam() const
Use external beam.
Definition: fvDOMI.H:139
Foam::solarCalculator
The solar calculator model provides information about the Sun direction and Sun load model....
Definition: solarCalculator.H:105
Foam::radiation::fvDOM::G
const volScalarField & G() const
Const access to incident radiation field.
Definition: fvDOMI.H:91
Foam::radiation::fvDOM::aLambda
const volScalarField & aLambda(const label lambdaI) const
Const access to wavelength total absorption coefficient.
Definition: fvDOMI.H:83
Foam::radiation::fvDOM::TypeName
TypeName("fvDOM")
Runtime type information.
Foam::T
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
Definition: FieldFieldFunctions.C:58
Foam::radiation::fvDOM::Rp
virtual tmp< volScalarField > Rp() const
Source term component (for power of T^4)
Definition: fvDOM.C:638
Foam::PtrList
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
Definition: List.H:62
Foam::radiation::fvDOM::calculate
void calculate()
Solve radiation equation(s)
Definition: fvDOM.C:569
radiation
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::radiation::fvDOM::nTheta
label nTheta() const
Number of angles in theta.
Definition: fvDOMI.H:47
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:121
Foam::radiation::fvDOM::updateRaysDir
void updateRaysDir()
Rotate rays according to Sun direction.
Definition: fvDOM.C:99
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::radiation::fvDOM::solarCalc
const solarCalculator & solarCalc() const
Solar calculator.
Definition: fvDOM.C:769
Foam::radiation::fvDOM::meshOrientation
vector meshOrientation() const
Return meshOrientation.
Definition: fvDOMI.H:127
Foam::radiation::fvDOM::qr
const volScalarField & qr() const
Const access to total radiative heat flux field.
Definition: fvDOMI.H:97
Foam::autoPtr
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: HashPtrTable.H:53
Foam::radiation::fvDOM::read
bool read()
Read radiation properties dictionary.
Definition: fvDOM.C:551
Foam::radiation::fvDOM::useSolarLoad
bool useSolarLoad() const
Use solar load.
Definition: fvDOMI.H:133
Foam::Vector< scalar >
Foam::List< scalar >
Foam::radiation::radiationModel
Top level model for radiation modelling.
Definition: radiationModel.H:75
Foam::radiation::fvDOM::spectralDistribution
const scalarList & spectralDistribution() const
Energy spectral distribution for external beam.
Definition: fvDOMI.H:146
Foam::radiation::blackBodyEmission
Class black body emission.
Definition: blackBodyEmission.H:58
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::radiation::fvDOM::blackBody
const blackBodyEmission & blackBody() const
Const access to black body.
Definition: fvDOMI.H:115
radiationModel.H
Foam::radiation::fvDOM::Ru
virtual tmp< volScalarField::Internal > Ru() const
Source term component (constant)
Definition: fvDOM.C:682
Foam::radiation::fvDOM::setRayIdLambdaId
void setRayIdLambdaId(const word &name, label &rayId, label &lambdaId) const
Set the rayId and lambdaId from by decomposing an intensity.
Definition: fvDOM.C:754
Foam::radiation::fvDOM::nBands
label nBands() const
Number of bands.
Definition: fvDOMI.H:71
Foam::radiation::fvDOM
Finite Volume Discrete Ordinates Method. Solves the RTE equation for n directions in a participating ...
Definition: fvDOM.H:119