Go to the documentation of this file.
49 { sunDirModel::mSunDirConstant,
"sunDirConstant" },
50 { sunDirModel::mSunDirTracking,
"sunDirTracking" },
60 { sunLModel::mSunLoadConstant,
"sunLoadConstant" },
62 sunLModel::mSunLoadFairWeatherConditions,
63 "sunLoadFairWeatherConditions"
65 { sunLModel::mSunLoadTheoreticalMaximum,
"sunLoadTheoreticalMaximum" },
71 void Foam::solarCalculator::calculateBetaTheta()
74 switch (sunDirectionModel_)
87 scalar LSM = 15.0*(dict_.get<scalar>(
"localStandardMeridian"));
89 scalar
D = dict_.get<scalar>(
"startDay") +
runTime/86400.0;
90 scalar
M = 6.24004 + 0.0172*
D;
91 scalar EOT = -7.659*
sin(
M) + 9.863*
sin(2*
M + 3.5932);
93 dict_.readEntry(
"startTime", startTime_);
95 scalar LST = startTime_ +
runTime/3600.0;
97 scalar LON = dict_.get<scalar>(
"longitude");
99 scalar AST = LST + EOT/60.0 + (LON - LSM)/15;
105 scalar
L =
degToRad(dict_.get<scalar>(
"latitude"));
125 void Foam::solarCalculator::calculateSunDirection()
132 new coordinateSystem(
"grid",
Zero, gridUp_, eastDir_)
136 direction_.z() = -
sin(beta_);
137 direction_.y() =
cos(beta_)*
cos(theta_);
138 direction_.x() =
cos(beta_)*
sin(theta_);
140 direction_.normalise();
143 <<
"Sun direction in absolute coordinates : " << direction_ <<
endl;
146 direction_ = coord_->transform(direction_);
149 <<
"Sun direction in the Grid coordinates : " << direction_ <<
endl;
153 void Foam::solarCalculator::init()
155 switch (sunDirectionModel_)
157 case mSunDirConstant:
159 if (dict_.readIfPresent(
"sunDirection", direction_))
161 direction_.normalise();
165 calculateBetaTheta();
166 calculateSunDirection();
171 case mSunDirTracking:
173 if (word(mesh_.ddtScheme(
"default")) ==
"steadyState")
176 <<
" Sun direction model can not be sunDirtracking if the "
182 "sunTrackingUpdateInterval",
183 sunTrackingUpdateInterval_
186 calculateBetaTheta();
187 calculateSunDirection();
192 switch (sunLoadModel_)
194 case mSunLoadConstant:
196 dict_.readEntry(
"directSolarRad", directSolarRad_);
197 dict_.readEntry(
"diffuseSolarRad", diffuseSolarRad_);
200 case mSunLoadFairWeatherConditions:
204 "skyCloudCoverFraction",
205 skyCloudCoverFraction_
208 dict_.readEntry(
"A", A_);
209 dict_.readEntry(
"B", B_);
211 if (!dict_.readIfPresent(
"beta", beta_))
213 calculateBetaTheta();
217 (1.0 - 0.75*
pow(skyCloudCoverFraction_, 3.0))
220 dict_.readEntry(
"groundReflectivity", groundReflectivity_);
223 case mSunLoadTheoreticalMaximum:
225 dict_.readEntry(
"Setrn", Setrn_);
226 dict_.readEntry(
"SunPrime", SunPrime_);
227 directSolarRad_ = Setrn_*SunPrime_;
229 dict_.readEntry(
"groundReflectivity", groundReflectivity_);
238 Foam::solarCalculator::solarCalculator
247 directSolarRad_(0.0),
248 diffuseSolarRad_(0.0),
249 groundReflectivity_(0.0),
254 skyCloudCoverFraction_(0.0),
257 C_(
dict.get<scalar>(
"C")),
260 sunDirectionModelTypeNames_.get(
"sunDirectionModel",
dict)
262 sunLoadModel_(sunLoadModelTypeNames_.get(
"sunLoadModel",
dict)),
273 switch (sunDirectionModel_)
275 case mSunDirConstant:
279 case mSunDirTracking:
281 calculateBetaTheta();
282 calculateSunDirection();
283 directSolarRad_ = A_/
exp(B_/
sin(
max(beta_, ROOTVSMALL)));
const vector L(dict.get< vector >("L"))
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
static constexpr const zero Zero
Global zero (0)
Different types of constants.
dimensionedScalar sin(const dimensionedScalar &ds)
volScalarField H(IOobject("H", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), mesh, dimensionedScalar(dimLength, Zero))
Unit conversion functions.
Ostream & endl(Ostream &os)
Add newline and flush stream.
dimensionedScalar exp(const dimensionedScalar &ds)
sunLModel
Direct sun load models.
constexpr scalar radToDeg(const scalar rad) noexcept
Conversion from radians to degrees.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Mesh data needed to do the Finite Volume discretisation.
Vector< scalar > vector
A scalar version of the templated Vector.
constexpr scalar degToRad(const scalar deg) noexcept
Conversion from degrees to radians.
errorManipArg< error, int > exit(error &err, const int errNo=1)
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
void correctSunDirection()
Recalculate.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
#define DebugInfo
Report an information message using Foam::Info.
constexpr scalar pi(M_PI)
static const Enum< sunDirModel > sunDirectionModelTypeNames_
Sun direction models.
dimensionedScalar acos(const dimensionedScalar &ds)
const dimensionedScalar e
Elementary charge.
const dimensionedScalar & D
defineTypeNameAndDebug(combustionModel, 0)
sunDirModel
Sun direction models.
static const Enum< sunLModel > sunLoadModelTypeNames_
Sun load models.
dimensionedScalar asin(const dimensionedScalar &ds)
dimensionedScalar cos(const dimensionedScalar &ds)