solarCalculator.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) 2015-2021 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26Class
27 Foam::solarCalculator
28
29Description
30 A solar calculator model providing models
31 for the solar direction and solar loads.
32
33 Available models for the solar direction:
34 - \c constant: Constant sunbeam direction.
35 - \c tracking: Transient model calculating sunbeam direction
36 based on a given set of parameters.
37
38 Available models for the solar load:
39 - \c constant: Constant solar load.
40 - \c timeDependent: Time-dependent solar load.
41 - \c fairWeather: Solar fluxes are calculated following
42 the "Fair Weather Conditions Method from the ASHRAE Handbook".
43 - \c theoreticalMaximum: Theoretically maximum solar load.
44
45Usage
46 Minimal examples by using \c constant/radiationProperties:
47
48 \c sunDirectionModel - Option-1:
49 \verbatim
50 solarLoadCoeffs
51 {
52 sunDirectionModel constant;
53 sunDirection (1 0 -1);
54 }
55 \endverbatim
56
57 where the entries mean:
58 \table
59 Property | Description | Type | Reqd | Deflt
60 sunDirection | Sunbeam direction | vector | no | calculated
61 \endtable
62
63
64 \c sunDirectionModel - Option-2:
65 \verbatim
66 solarLoadCoeffs
67 {
68 sunDirectionModel tracking;
69 sunTrackingUpdateInterval 800;
70 localStandardMeridian 9;
71 startDay 204;
72 startTime 15;
73 longitude 139.74;
74 latitude 35.658;
75 gridUp (0 0 1);
76 gridEast (1 0 0);
77 }
78 \endverbatim
79
80 where the entries mean:
81 \table
82 Property | Description | Type | Reqd | Deflt
83 sunTrackingUpdateInterval | Interval to update the Sun direction <!--
84 --> [decimal hours] | scalar | yes | -
85 localStandardMeridian | GMT (Local Zone Meridian) [decimal hours]<!--
86 --> | scalar | yes | -
87 startDay | Day from 1 to 365 | scalar | yes | -
88 startTime | Start time for the Sun position [decimal hours] <!--
89 --> | scalar | yes | -
90 longitude | Geographic coordinate specifying the east–west <!--
91 --> position of a point on the surface of a planetary <!--
92 --> body [degree] | scalar | yes | -
93 latitude | Geographic coordinate specifying the north–south <!--
94 --> position of a point on the surface of a planetary <!--
95 --> body [degree] | scalar | yes | -
96 gridUp | Grid orientation upwards | vector | yes | -
97 gridEast | Grid orientation eastwards | vector | yes | -
98 \endtable
99
100
101 \c sunLoadModel - Option-1:
102 \verbatim
103 solarLoadCoeffs
104 {
105 sunLoadModel constant;
106 directSolarRad 100;
107 diffuseSolarRad 0;
108 }
109 \endverbatim
110
111 where the entries mean:
112 \table
113 Property | Description | Type | Reqd | Deflt
114 directSolarRad | Direct solar irradiation [W/m2] | scalar | yes | -
115 diffuseSolarRad | Diffuse solar irradiation on vertical surfaces <!--
116 --> [W/m2] | scalar | yes | -
117 \endtable
118
119
120 \c sunLoadModel - Option-2:
121 \verbatim
122 solarLoadCoeffs
123 {
124 sunLoadModel timeDependent;
125 directSolarRad <Function1<scalar>>;
126 diffuseSolarRad <Function1<scalar>>;
127 }
128 \endverbatim
129
130 where the entries mean:
131 \table
132 Property | Description | Type | Reqd | Deflt
133 directSolarRad | Time-series of direct solar irradiation <!--
134 --> [W/m2] | Function1<scalar> | yes | -
135 diffuseSolarRad | Time-series of diffuse solar irradiation on <!--
136 --> vertical surfaces [W/m2] <!--
137 --> | Function1<scalar> | yes | -
138 \endtable
139
140
141 \c sunLoadModel - Option-3:
142 \verbatim
143 solarLoadCoeffs
144 {
145 sunLoadModel fairWeather;
146 skyCloudCoverFraction 0.25;
147 groundReflectivity 1.0;
148 A 0.1;
149 B 0.2;
150 C 0.058;
151 beta 0.15;
152 }
153 \endverbatim
154
155 where the entries mean:
156 \table
157 Property | Description | Type | Reqd | Deflt
158 A | Apparent solar irradiation at air mass m = 0 <!--
159 --> | scalar | yes | -
160 B | Atmospheric extinction coefficient <!--
161 --> | scalar | yes | -
162 C | Solar diffusivity constant | scalar | yes | -
163 groundReflectivity | Ground reflectivity | scalar | yes | -
164 skyCloudCoverFraction | Fraction of sky covered by clouds [0,1] <!--
165 --> | scalar | no | 0
166 beta | Solar altitude (in degrees) above the horizontal <!--
167 --> | scalar | no | calculated
168 \endtable
169
170 In this model the flux is calculated as:
171 \verbatim
172 directSolarRad = (1 - 0.75*skyCloudCoverFraction^3)*A/exp(B/sin(beta));
173 \endverbatim
174
175
176 \c sunLoadModel - Option-4:
177 \verbatim
178 solarLoadCoeffs
179 {
180 sunLoadModel theoreticalMaximum;
181 Setrn 1.0;
182 SunPrime 4.0;
183 groundReflectivity 1.0;
184 C 0.058;
185 }
186 \endverbatim
187
188 where the entries mean:
189 \table
190 Property | Description | Type | Reqd | Deflt
191 Setrn | Parameter in maximum theoretical direct solar <!--
192 --> model | scalar | yes | -
193 SunPrime | Parameter in maximum theoretical direct solar <!--
194 --> model | scalar | yes | -
195 groundReflectivity | Ground reflectivity | scalar | yes | -
196 C | Solar diffusivity constant | scalar | yes | -
197 \endtable
198
199 In this model the flux is calculated as:
200 \verbatim
201 directSolarRad = Setrn*SunPrime;
202 \endverbatim
203
204Note
205 - The \c sunDirectionModel:tracking can only be used
206 in transient calculations.
207 - The keyword \c sunTrackingUpdateInterval (in hours) specifies on which
208 interval is the Sun direction updated.
209 - The diffuse on vertical/horizontal walls and ground-reflected radiation are
210 calculated following the ASHRAE Handbook.
211 - The range of \c skyCloudCoverFraction is [0,1].
212
213SourceFiles
214 solarCalculator.C
215
216\*---------------------------------------------------------------------------*/
217
218#ifndef solarCalculator_H
219#define solarCalculator_H
220
221#include "fvMesh.H"
222#include "meshTools.H"
223#include "DynamicField.H"
224#include "HashSet.H"
225#include "coordinateSystem.H"
226#include "Function1.H"
227
228// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
229
230namespace Foam
231{
232
233/*---------------------------------------------------------------------------*\
234 Class solarCalculator Declaration
235\*---------------------------------------------------------------------------*/
236
237class solarCalculator
238{
239public:
240
241 // Public Enumeration
242
243 //- Options for the Sun direction models
244 enum sunDirModel
245 {
248 };
249
250 //- Names for sunDirModel
251 static const Enum<sunDirModel> sunDirectionModelTypeNames_;
252
253 //- Options for the Sun load models
254 enum sunLModel
255 {
260 };
261
262 //- Names for sunLModel
263 static const Enum<sunLModel> sunLModelTypeNames_;
264
265
266private:
267
268 // Private Data
269
270 //- Reference to mesh
271 const fvMesh& mesh_;
272
273 //- Dictionary
274 dictionary dict_;
275
276 //- Sun direction model
277 sunDirModel sunDirectionModel_;
278
279 //- Sun load model
280 sunLModel sunLoadModel_;
281
282
283 // sunDirectionModel = constant
284
285 //- Sunbeam direction
286 vector direction_;
287
288
289 // sunDirectionModel = tracking
290
291 //- Interval to update the Sun direction [decimal hours]
292 scalar sunTrackingUpdateInterval_;
293
294 //- Start time for the Sun position [decimal hours]
295 scalar startTime_;
296
297 //- Up grid orientation
298 vector gridUp_;
299
300 //- East grid orientation
301 vector eastDir_;
302
303 //- Grid coordinate system
304 autoPtr<coordinateSystem> coord_;
305
306
307 // sunLoadModel = constant
308
309 //- Direct solar irradiation
310 scalar directSolarRad_;
311
312 //- Diffuse solar irradiation on vertical surfaces
313 scalar diffuseSolarRad_;
314
315
316 // sunLoadModel = timeDependent
317
318 //- Time-series of direct solar irradiation
319 autoPtr<Function1<scalar>> directSolarRads_;
320
321 //- Time-series of diffuse solar irradiation on vertical surfaces
322 autoPtr<Function1<scalar>> diffuseSolarRads_;
323
324
325 // sunLoadModel = fairWeather
326
327 //- Sky cloud cover fraction [0-1]
328 scalar skyCloudCoverFraction_;
329
330 //- Ground reflectivity
331 scalar groundReflectivity_;
332
333 //- Fair weather direct solar load model parameters
334 scalar A_;
335 scalar B_;
336 scalar beta_;
337 scalar theta_;
338
339 //- Diffusive solar load model parameter
340 scalar C_;
341
342
343 // sunLoadModel = theoreticalMaximum
344
345 //- Maximum theoretical direct solar load model parameters
346 scalar Setrn_;
347 scalar SunPrime_;
348
349
350 //- No copy construct
351 solarCalculator(const solarCalculator&) = delete;
352
353 //- No copy assignment
354 void operator=(const solarCalculator&) = delete;
355
356
357 // Private Member Functions
358
359 //- Initialise model parameters
360 void initialise();
361
362 //- Calculate beta and theta angles
363 void calculateBetaTheta();
364
365 //- Calculate the Sun direction
366 void calculateSunDirection();
367
368
369public:
370
371 // Declare name of the class and its debug switch
372 ClassName("solarCalculator");
373
374
375 // Constructors
376
377 //- Construct from dictionary
378 solarCalculator(const dictionary&, const fvMesh&);
379
380
381 //- Destructor
382 ~solarCalculator() = default;
383
384
385 // Member Functions
386
387 // Evaluation
388
389 //- Correct the Sun direction
390 void correctSunDirection();
391
392 //- Correct direct solar irradiation
394
395 //- Correct diffuse solar irradiation
397
398
399 // Access
400
401 //- Return const access to the Sun direction model
403 {
404 return sunDirectionModel_;
405 }
406
407 //- Return const access to the Sun load model
408 const sunLModel& sunLoadModel() const noexcept
409 {
410 return sunLoadModel_;
411 }
412
413 //- Return non-const access to the Sun direction
415 {
416 return direction_;
417 }
418
419 //- Return const access to the Sun direction
420 const vector& direction() const noexcept
421 {
422 return direction_;
423 }
424
425 //- Return non-const access to the direct solar irradiation
426 scalar& directSolarRad()
427 {
428 return directSolarRad_;
429 }
430
431 //- Return const access to the direct solar irradiation
432 const scalar& directSolarRad() const noexcept
433 {
434 return directSolarRad_;
435 }
436
437 //- Return non-const access to the diffuse solar irradiation
438 scalar& diffuseSolarRad()
439 {
440 return diffuseSolarRad_;
441 }
442
443 //- Return const access to the diffuse solar irradiation
444 const scalar& diffuseSolarRad() const noexcept
445 {
446 return diffuseSolarRad_;
447 }
448
449 //- Return const access to the C constant
450 scalar C() const noexcept
452 return C_;
454
455 //- Return const access to beta
456 scalar beta() const noexcept
457 {
458 return beta_;
459 }
460
461 //- Return const access to theta
462 scalar theta() const noexcept
464 return theta_;
466
467 //- Return const access to the ground reflectivity
468 scalar groundReflectivity() const noexcept
469 {
470 return groundReflectivity_;
471 }
472
473 //- Return const access to the coordinate system
474 const coordinateSystem& coord() const noexcept
475 {
476 return *coord_;
477 }
478
479 //- Return const access to sunTrackingUpdateInterval
481 {
482 return sunTrackingUpdateInterval_;
483 }
484
485 //- Return const access to startTime
486 scalar startTime() const noexcept
487 {
488 return startTime_;
489 }
490};
491
492
493// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
494
495} // End namespace Foam
496
497// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
498
499#endif
500
501// ************************************************************************* //
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: Enum.H:61
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
Definition: autoPtr.H:66
Base class for coordinate system specification, the default coordinate system type is cartesian .
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
A solar calculator model providing models for the solar direction and solar loads.
void correctSunDirection()
Correct the Sun direction.
sunDirModel
Options for the Sun direction models.
scalar C() const noexcept
Return const access to the C constant.
vector & direction()
Return non-const access to the Sun direction.
ClassName("solarCalculator")
static const Enum< sunLModel > sunLModelTypeNames_
Names for sunLModel.
void correctDiffuseSolarRad()
Correct diffuse solar irradiation.
scalar sunTrackingUpdateInterval() const noexcept
Return const access to sunTrackingUpdateInterval.
scalar beta() const noexcept
Return const access to beta.
void correctDirectSolarRad()
Correct direct solar irradiation.
const sunLModel & sunLoadModel() const noexcept
Return const access to the Sun load model.
const vector & direction() const noexcept
Return const access to the Sun direction.
const scalar & diffuseSolarRad() const noexcept
Return const access to the diffuse solar irradiation.
scalar startTime() const noexcept
Return const access to startTime.
~solarCalculator()=default
Destructor.
scalar & diffuseSolarRad()
Return non-const access to the diffuse solar irradiation.
scalar groundReflectivity() const noexcept
Return const access to the ground reflectivity.
const sunDirModel & sunDirectionModel() const noexcept
Return const access to the Sun direction model.
scalar theta() const noexcept
Return const access to theta.
const coordinateSystem & coord() const noexcept
Return const access to the coordinate system.
scalar & directSolarRad()
Return non-const access to the direct solar irradiation.
const scalar & directSolarRad() const noexcept
Return const access to the direct solar irradiation.
static const Enum< sunDirModel > sunDirectionModelTypeNames_
Names for sunDirModel.
sunLModel
Options for the Sun load models.
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
#define ClassName(TypeNameString)
Add typeName information from argument TypeNameString to a class.
Definition: className.H:67
Namespace for OpenFOAM.
const direction noexcept
Definition: Scalar.H:223