37namespace functionObjects
48 tmp<volScalarField> tmagU = mag(lookupObject<volVectorField>(
"U"));
74 lookupObject<volScalarField>(
"T").boundaryField();
76 scalar areaIntegral = 0;
77 scalar TareaIntegral = 0;
82 const fvPatch& pTBf = TBf[patchi].
patch();
85 if (isType<wallFvPatch>(pTBf))
87 areaIntegral +=
gSum(pSf);
88 TareaIntegral +=
gSum(pSf*pT);
92 Trad.value() = TareaIntegral/areaIntegral;
96 if ((Trad.value() < 283.15) || (Trad.value() > 313.15))
99 <<
"The calculated mean wall radiation temperature is out of the\n"
100 <<
"bounds specified in EN ISO 7730:2005\n"
101 <<
"Valid range is 10 degC < T < 40 degC\n"
102 <<
"The actual value is: " << (Trad.value() - 273.15) <<
nl <<
endl;
119 if (pSat_.value() == 0)
121 const auto&
T = lookupObject<volScalarField>(
"T");
124 tpSat = kPaToPa*
exp(
A -
B/(
T + C));
162 tmp<volScalarField> tTcl
187 Tcl = (Tcl + Tcl.prevIter())/2;
199 pos(hcForced - hcNatural)*hcForced
200 +
neg0(hcForced - hcNatural)*hcNatural;
205 - factor2*(metabolicRateSI - extWorkSI)
206 - Icl_*factor3*fcl_*(
pow4(Tcl) -
pow4(Trad))
207 - Icl_*fcl_*hc*(Tcl -
T);
211 Tcl.clip(Tmin, Tmax);
213 }
while (!converged(Tcl) && i++ < maxClothIter_);
215 if (i == maxClothIter_)
218 <<
"The surface cloth temperature did not converge within " << i
219 <<
" iterations" <<
nl;
232 max(
mag(
phi.primitiveField() -
phi.prevIter().primitiveField()))
247 clothing_(
"clothing",
dimless, 0),
251 relHumidity_(
"relHumidity",
dimless, 0.5),
270 clothing_.readIfPresent(
dict);
271 metabolicRate_.readIfPresent(
dict);
272 extWork_.readIfPresent(
dict);
273 pSat_.readIfPresent(
dict);
274 tolerance_ =
dict.getOrDefault(
"tolerance", 1
e-4);
275 maxClothIter_ =
dict.getOrDefault(
"maxClothIter", 100);
276 meanVelocity_ =
dict.getOrDefault<
bool>(
"meanVelocity",
false);
279 if (
dict.found(relHumidity_.name()))
281 relHumidity_.read(
dict);
286 if (
dict.found(Trad_.name()))
295 Icl_.value() <= 0.078
296 ? 1.0 + 1.290*Icl_.value()
297 : 1.05 + 0.645*Icl_.value();
315 const auto&
T = lookupObject<volScalarField>(
"T");
324 mesh_.time().timeName(),
371 (metabolicRateSI - extWorkSI).value() < factor8.
value() ? 0 : 0.42
374 Info<<
"Calculating the predicted mean vote (PMV)" <<
endl;
380 (factor1*
exp(factor2*metabolicRateSI) + factor3)
382 (metabolicRateSI - extWorkSI)
388 - factor6*(metabolicRateSI - extWorkSI)
393 - factor7*(metabolicRateSI - extWorkSI - factor8)
396 - factor9*metabolicRateSI*(factor10 - pSat*relHumidity_)
399 - factor11*metabolicRateSI*(factor12 -
T)
402 - factor13*fcl_*(
pow4(Tcloth) -
pow4(Trad))
405 - fcl_*hc*(Tcloth -
T)
409 Info<<
"Calculating the predicted percentage of dissatisfaction (PPD)"
414 100 - 95*
exp(-0.03353*
pow4(PMV()) - 0.21790*
sqr(PMV()));
416 Info<<
"Calculating the draught rating (DR)\n";
425 Umag.
clip(Umin, Umax);
434 mesh_.time().timeName(),
441 if (foundObject<volScalarField>(
"k"))
443 const auto&
k = lookupObject<volScalarField>(
"k");
444 TI =
sqrt(2/3*
k)/Umag;
456 correctUnit*(factor12 -
T)*
pow(Umag - Umin, 0.62)*(pre*Umag*TI + C1);
462 word fieldNamePMV =
"PMV";
463 word fieldNamePPD =
"PPD";
464 word fieldNameDR =
"DR";
465 word fieldNameTop =
"Top";
469 store(fieldNamePMV, PMV)
470 && store(fieldNamePPD, PPD)
471 && store(fieldNameDR, DR)
472 && store(fieldNameTop, Top)
482 && writeObject(
"PPD")
484 && writeObject(
"Top")
static const Foam::dimensionedScalar A("", Foam::dimPressure, 611.21)
static const Foam::dimensionedScalar B("", Foam::dimless, 18.678)
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
dimensioned< Type > weightedAverage(const DimensionedField< scalar, GeoMesh > &weightField) const
Calculate and return weighted average.
Internal & ref(const bool updateAccessTime=true)
Return a reference to the dimensioned internal field.
GeometricBoundaryField< scalar, fvPatchField, volMesh > Boundary
Type of boundary fields.
void clip(const dimensioned< MinMax< Type > > &range)
Clip the field to be bounded within the specified range.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
virtual bool read()
Re-read model coefficients if they have changed.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Dimension set for the base types, which can be used to implement rigorous dimension checking for alge...
const dimensionSet & dimensions() const
Return const reference to dimensions.
const Type & value() const
Return const reference to value.
Abstract base-class for Time/database function objects.
Calculates the thermal comfort quantities predicted mean vote (PMV), predicted percentage of dissatis...
virtual bool write()
Write the PPD and PMV fields.
virtual bool read(const dictionary &)
Read the data needed for the comfort calculation.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
const fvMesh & mesh_
Reference to the fvMesh.
Computes the magnitude of an input field.
Computes the power of an input volScalarField.
const DimensionedField< scalar, volMesh > & V() const
Return cell volumes.
const fvPatch & patch() const
Return patch.
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionSet dimPressure
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
dimensionedScalar pos(const dimensionedScalar &ds)
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Type gSum(const FieldField< Field, Type > &f)
dimensionedScalar exp(const dimensionedScalar &ds)
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
GeometricField< scalar, fvPatchField, volMesh > volScalarField
const dimensionSet dimVelocity
messageStream Info
Information stream (stdout output on master, null elsewhere)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const dimensionSet dimTemperature(0, 0, 0, 1, 0, 0, 0)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
dimensionedScalar pow4(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
dimensionedScalar neg0(const dimensionedScalar &ds)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
fvPatchField< scalar > fvPatchScalarField
dimensionedScalar pow025(const dimensionedScalar &ds)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.