49 { sunDirModel::mSunDirConstant,
"constant" },
50 { sunDirModel::mSunDirTracking,
"tracking" },
53 { sunDirModel::mSunDirConstant,
"sunDirConstant" },
54 { sunDirModel::mSunDirTracking,
"sunDirTracking" }
64 { sunLModel::mSunLoadConstant,
"constant" },
65 { sunLModel::mSunLoadTimeDependent,
"timeDependent" },
66 { sunLModel::mSunLoadFairWeatherConditions,
"fairWeather" },
67 { sunLModel::mSunLoadTheoreticalMaximum,
"theoreticalMaximum" },
70 { sunLModel::mSunLoadConstant,
"sunLoadConstant" },
72 sunLModel::mSunLoadFairWeatherConditions,
73 "sunLoadFairWeatherConditions"
75 { sunLModel::mSunLoadTheoreticalMaximum,
"sunLoadTheoreticalMaximum" }
81void Foam::solarCalculator::calculateBetaTheta()
90 const scalar LSM = 15.0*(dict_.
get<scalar>(
"localStandardMeridian"));
92 const scalar
D = dict_.
get<scalar>(
"startDay") +
runTime/86400.0;
93 const scalar
M = 6.24004 + 0.0172*
D;
94 const scalar EOT = -7.659*
sin(
M) + 9.863*
sin(2*
M + 3.5932);
98 const scalar LST = startTime_ +
runTime/3600.0;
100 const scalar LON = dict_.
get<scalar>(
"longitude");
102 const scalar AST = LST + EOT/60.0 + (LON - LSM)/15;
106 const scalar
H =
degToRad(15*(AST - 12));
108 const scalar
L =
degToRad(dict_.
get<scalar>(
"latitude"));
128void Foam::solarCalculator::calculateSunDirection()
139 direction_.z() = -
sin(beta_);
140 direction_.y() =
cos(beta_)*
cos(theta_);
141 direction_.x() =
cos(beta_)*
sin(theta_);
143 direction_.normalise();
146 <<
"Sun direction in absolute coordinates : " << direction_ <<
endl;
149 direction_ = coord_->transform(direction_);
152 <<
"Sun direction in the Grid coordinates : " << direction_ <<
endl;
156void Foam::solarCalculator::initialise()
158 switch (sunDirectionModel_)
160 case mSunDirConstant:
162 if (dict_.readIfPresent(
"sunDirection", direction_))
164 direction_.normalise();
168 calculateBetaTheta();
169 calculateSunDirection();
173 case mSunDirTracking:
175 if (
word(mesh_.ddtScheme(
"default")) ==
"steadyState")
178 <<
" Sun direction model can not be sunDirtracking if the "
184 "sunTrackingUpdateInterval",
185 sunTrackingUpdateInterval_
188 calculateBetaTheta();
189 calculateSunDirection();
194 switch (sunLoadModel_)
196 case mSunLoadConstant:
198 dict_.readEntry(
"directSolarRad", directSolarRad_);
199 dict_.readEntry(
"diffuseSolarRad", diffuseSolarRad_);
202 case mSunLoadTimeDependent:
204 directSolarRads_.reset
214 diffuseSolarRads_.reset
225 directSolarRads_->value(mesh_.time().timeOutputValue());
227 diffuseSolarRads_->value(mesh_.time().timeOutputValue());
230 case mSunLoadFairWeatherConditions:
234 "skyCloudCoverFraction",
235 skyCloudCoverFraction_
238 dict_.readEntry(
"A", A_);
239 dict_.readEntry(
"B", B_);
240 dict_.readEntry(
"C", C_);
241 dict_.readEntry(
"groundReflectivity", groundReflectivity_);
242 if (!dict_.readIfPresent(
"beta", beta_))
244 calculateBetaTheta();
248 (1.0 - 0.75*
pow(skyCloudCoverFraction_, 3.0))
252 case mSunLoadTheoreticalMaximum:
254 dict_.readEntry(
"Setrn", Setrn_);
255 dict_.readEntry(
"SunPrime", SunPrime_);
256 dict_.readEntry(
"groundReflectivity", groundReflectivity_);
257 dict_.readEntry(
"C", C_);
259 directSolarRad_ = Setrn_*SunPrime_;
278 sunDirectionModelTypeNames_.get(
"sunDirectionModel",
dict)
280 sunLoadModel_(sunLModelTypeNames_.get(
"sunLoadModel",
dict)),
282 sunTrackingUpdateInterval_(0),
291 skyCloudCoverFraction_(0),
292 groundReflectivity_(0),
309 if (sunDirectionModel_ == mSunDirTracking)
311 calculateBetaTheta();
312 calculateSunDirection();
313 directSolarRad_ = A_/
exp(B_/
sin(
max(beta_, ROOTVSMALL)));
320 if (sunLoadModel_ == mSunLoadTimeDependent)
322 directSolarRad_ = directSolarRads_->value(mesh_.time().value());
329 if (sunLoadModel_ == mSunLoadTimeDependent)
331 diffuseSolarRad_ = diffuseSolarRads_->value(mesh_.time().value());
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Top level data entry class for use in dictionaries. Provides a mechanism to specify a variable as a c...
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,...
T get(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
const Type & value() const
Return const reference to value.
Mesh data needed to do the Finite Volume discretisation.
const Time & time() const
Return the top-level database.
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.
static const Enum< sunLModel > sunLModelTypeNames_
Names for sunLModel.
void correctDiffuseSolarRad()
Correct diffuse solar irradiation.
void correctDirectSolarRad()
Correct 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...
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define DebugInfo
Report an information message using Foam::Info.
constexpr scalar pi(M_PI)
Different types of constants.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
dimensionedScalar asin(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
quaternion normalised(const quaternion &q)
Return the normalised (unit) quaternion of the given quaternion.
dimensionedScalar sin(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Ostream & endl(Ostream &os)
Add newline and flush stream.
constexpr scalar degToRad() noexcept
Multiplication factor for degrees to radians conversion.
constexpr scalar radToDeg() noexcept
Multiplication factor for radians to degrees conversion.
static constexpr const zero Zero
Global zero (0)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar acos(const dimensionedScalar &ds)
constexpr char nl
The newline '\n' character (0x0a)
constexpr char tab
The tab '\t' character(0x09)
const dimensionedScalar & D
Unit conversion functions.
const vector L(dict.get< vector >("L"))