From 280c055ef6322cd5ca0851b370268a496aa7746f Mon Sep 17 00:00:00 2001 From: Henry Weller Date: Sat, 19 Oct 2019 23:08:34 +0100 Subject: [PATCH] functionObjects::comfort: New functionObject to calculate fields relating to thermal comfort Description Calculates the thermal comfort quantities predicted mean vote (PMV) and predicted percentage of dissatisfaction (PPD) based on DIN ISO EN 7730:2005. Usage \table Property | Description | Required | Default value clothing | The insulation value of the cloth | no | 0 metabolicRate | The metabolic rate | no | 0.8 extWork | The external work | no | 0 Trad | Radiation temperature | no | -1 relHumidity | Relative humidity of the air | no | 50 pSat | Saturation pressure of water | no | -1 tolerance | Residual control for the cloth temperature | no | 1e-5 maxClothIter | Maximum number of iterations | no | 0 meanVelocity | Use a constant mean velocity in the whole domain | no |\ false \endtable \table Predicted Mean Vote (PMV) | evaluation + 3 | hot + 2 | warm + 1 | slightly warm + 0 | neutral - 1 | slightly cool - 2 | cool - 3 | cold \endtable \verbatim comfortAnalysis { type comfort; libs ("libfieldFunctionObjects.so"); executeControl writeTime; writeControl writeTime; } \endverbatim The new tutorial case heatTransfer/buoyantSimpleFoam/comfortHotRoom is provided to demonstrate the calculation of PMV and PPD using the comfort functionObject. This work is based on code and case contributed by Tobias Holzmann. --- src/functionObjects/field/Make/files | 2 + src/functionObjects/field/age/age.C | 4 +- src/functionObjects/field/comfort/comfort.C | 430 ++++++++++++++++++ src/functionObjects/field/comfort/comfort.H | 206 +++++++++ .../buoyantSimpleFoam/comfortHotRoom/0/T | 42 ++ .../buoyantSimpleFoam/comfortHotRoom/0/U | 42 ++ .../buoyantSimpleFoam/comfortHotRoom/0/alphat | 44 ++ .../comfortHotRoom/0/epsilon | 45 ++ .../buoyantSimpleFoam/comfortHotRoom/0/k | 45 ++ .../buoyantSimpleFoam/comfortHotRoom/0/nut | 43 ++ .../buoyantSimpleFoam/comfortHotRoom/0/p | 42 ++ .../buoyantSimpleFoam/comfortHotRoom/0/p_rgh | 43 ++ .../buoyantSimpleFoam/comfortHotRoom/Allrun | 15 + .../comfortHotRoom/constant/g | 21 + .../constant/thermophysicalProperties | 54 +++ .../constant/turbulenceProperties | 30 ++ .../comfortHotRoom/system/blockMeshDict | 49 ++ .../comfortHotRoom/system/controlDict | 76 ++++ .../comfortHotRoom/system/createPatchDict | 47 ++ .../comfortHotRoom/system/fvSchemes | 61 +++ .../comfortHotRoom/system/fvSolution | 74 +++ .../comfortHotRoom/system/topoSetDict | 41 ++ 22 files changed, 1453 insertions(+), 3 deletions(-) create mode 100644 src/functionObjects/field/comfort/comfort.C create mode 100644 src/functionObjects/field/comfort/comfort.H create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/T create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/U create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/alphat create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/epsilon create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/k create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/nut create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/p create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/p_rgh create mode 100755 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/Allrun create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/g create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/thermophysicalProperties create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/turbulenceProperties create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/blockMeshDict create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/controlDict create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/createPatchDict create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/fvSchemes create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/fvSolution create mode 100644 tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/topoSetDict diff --git a/src/functionObjects/field/Make/files b/src/functionObjects/field/Make/files index fbba9b2287..2b950187e2 100644 --- a/src/functionObjects/field/Make/files +++ b/src/functionObjects/field/Make/files @@ -69,4 +69,6 @@ interfaceHeight/interfaceHeight.C age/age.C +comfort/comfort.C + LIB = $(FOAM_LIBBIN)/libfieldFunctionObjects diff --git a/src/functionObjects/field/age/age.C b/src/functionObjects/field/age/age.C index 65f87f08ee..75953dd601 100644 --- a/src/functionObjects/field/age/age.C +++ b/src/functionObjects/field/age/age.C @@ -243,9 +243,7 @@ bool Foam::functionObjects::age::execute() Info<< "Min/max age:" << min(age).value() << ' ' << max(age).value() << endl; - store(tage); - - return true; + return store(tage); } diff --git a/src/functionObjects/field/comfort/comfort.C b/src/functionObjects/field/comfort/comfort.C new file mode 100644 index 0000000000..5aad71231c --- /dev/null +++ b/src/functionObjects/field/comfort/comfort.C @@ -0,0 +1,430 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Copyright (C) 2019 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +\*---------------------------------------------------------------------------*/ + +#include "comfort.H" +#include "wallFvPatch.H" +#include "addToRunTimeSelectionTable.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +namespace Foam +{ +namespace functionObjects +{ + defineTypeNameAndDebug(comfort, 0); + addToRunTimeSelectionTable(functionObject, comfort, dictionary); +} +} + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +Foam::tmp Foam::functionObjects::comfort::magU() const +{ + tmp tmagU = mag(lookupObject("U")); + + // Switch to use the averaged velocity field in the domain. + // Consistent with EN ISO 7730 but does not make physical sense + if (meanVelocity_) + { + tmagU.ref() = tmagU->weightedAverage(mesh_.V()); + } + + return tmagU; +} + + +Foam::dimensionedScalar Foam::functionObjects::comfort::Trad() const +{ + dimensionedScalar Trad(Trad_); + + // The mean radiation is calculated by the mean wall temperatures + // which are summed and divided by the area | only walls are taken into + // account. This approach might be correct for a squared room but will + // defintely be inconsistent for complex room geometries. The norm does + // not provide any information about the calculation of this quantity. + if (!TradSet_) + { + const volScalarField::Boundary& TBf = + lookupObject("T").boundaryField(); + + scalar areaIntegral = 0; + scalar TareaIntegral = 0; + + forAll(TBf, patchi) + { + const fvPatchScalarField& pT = TBf[patchi]; + const fvPatch& pTBf = TBf[patchi].patch(); + const scalarField& pSf = pTBf.magSf(); + + if (isType(pTBf)) + { + areaIntegral += gSum(pSf); + TareaIntegral += gSum(pSf*pT); + } + } + + Trad.value() = TareaIntegral/areaIntegral; + } + + // Bounds based on EN ISO 7730 + if ((Trad.value() < 283.15) || (Trad.value() > 313.15)) + { + WarningInFunction + << "The calculated mean wall radiation temperature is out of the\n" + << "bounds specified in EN ISO 7730:2006\n" + << "Valid range is 10 degC < T < 40 degC\n" + << "The actual value is: " << Trad - 273.15 << nl << endl; + } + + return Trad; +} + + +Foam::tmp Foam::functionObjects::comfort::pSat() const +{ + static const dimensionedScalar kPaToPa(dimPressure, 1000); + static const dimensionedScalar A(dimless, 16.6563); + static const dimensionedScalar B(dimTemperature, 4030.183); + static const dimensionedScalar C(dimTemperature, -38.15); + + tmp tpSat(volScalarField::New("pSat", mesh_, pSat_)); + + // Calculate the saturation pressure if no user input is given + if (pSat_.value() == 0) + { + const volScalarField& T = lookupObject("T"); + + // Equation based on ISO 7730:2006 + tpSat = kPaToPa*exp(A - B/(T + C)); + } + + return tpSat; +} + + +Foam::tmp Foam::functionObjects::comfort::Tcloth +( + volScalarField& hc, + const dimensionedScalar& metabolicRateSI, + const dimensionedScalar& extWorkSI, + const volScalarField& T, + const dimensionedScalar& Trad +) +{ + const dimensionedScalar factor1(dimTemperature, 308.85); + + const dimensionedScalar factor2 + ( + dimTemperature/metabolicRateSI.dimensions(), + 0.028 + ); + + const dimensionedScalar factor3 + ( + dimensionSet(1, 0, -3, -4, 0, 0, 0), + 3.96e-8 + ); + + // Heat transfer coefficient based on forced convection [W/m^2/K] + const volScalarField hcForced + ( + dimensionedScalar(hc.dimensions()/sqrt(dimVelocity), 12.1) + *sqrt(magU()) + ); + + // Tcl [K] (surface cloth temperature) + tmp tTcl + ( + volScalarField::New + ( + "Tcl", + T.mesh(), + dimTemperature + ) + ); + + volScalarField& Tcl = tTcl.ref(); + + // Initial guess + Tcl = T; + + label i = 0; + + Tcl.storePrevIter(); + + // Iterative solving of equation (2) + do + { + Tcl = (Tcl + Tcl.prevIter())/2; + Tcl.storePrevIter(); + + // Heat transfer coefficient based on natural convection + volScalarField hcNatural + ( + dimensionedScalar(hc.dimensions()/pow025(dimTemperature), 2.38) + *pow025(mag(Tcl - T)) + ); + + // Set heat transfer coefficient based on equation (3) + hc = + pos(hcForced - hcNatural)*hcForced + + neg0(hcForced - hcNatural)*hcNatural; + + // Calculate surface temperature based on equation (2) + Tcl = + factor1 + - factor2*(metabolicRateSI - extWorkSI) + - Icl_*factor3*fcl_*(pow4(Tcl) - pow4(Trad)) + - Icl_*fcl_*hc*(Tcl - T); + + } while (!converged(Tcl) && i++ < maxClothIter_); + + if (i == maxClothIter_) + { + WarningInFunction + << "The surface cloth temperature did not converge within " << i + << " iterations\n"; + } + + return tTcl; +} + + +bool Foam::functionObjects::comfort::converged +( + const volScalarField& phi +) const +{ + return + max(mag(phi.primitiveField() - phi.prevIter().primitiveField())) + < tolerance_; +} + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::functionObjects::comfort::comfort +( + const word& name, + const Time& runTime, + const dictionary& dict +) +: + fvMeshFunctionObject(name, runTime, dict), + clothing_("clothing", dimless, 0), + metabolicRate_("metabolicRate", dimMass/pow3(dimTime), 0.8), + extWork_("extWork", dimMass/pow3(dimTime), 0), + TradSet_(false), + Trad_("Trad", dimTemperature, 0), + relHumidity_("relHumidity", dimless, 0.5), + pSat_("pSat", dimPressure, 0), + Icl_("Icl", dimensionSet(-1, 0, 3, 1, 0, 0, 0), 0), + fcl_("fcl", dimless, 0), + tolerance_(1e-4), + maxClothIter_(100), + meanVelocity_(false) +{ + read(dict); +} + + +// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * // + +Foam::functionObjects::comfort::~comfort() +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +bool Foam::functionObjects::comfort::read(const dictionary& dict) +{ + clothing_.readIfPresent(dict); + metabolicRate_.readIfPresent(dict); + extWork_.readIfPresent(dict); + pSat_.readIfPresent(dict); + tolerance_ = dict.lookupOrDefault("tolerance", 1e-4); + maxClothIter_ = dict.lookupOrDefault("maxClothIter", 100); + meanVelocity_ = dict.lookupOrDefault("meanVelocity", false); + + // Read relative humidity if provided and convert from % to fraction + if (dict.found(relHumidity_.name())) + { + relHumidity_.read(dict); + relHumidity_ /= 100; + } + + // Read radiation temperature if provided + if (dict.found(Trad_.name())) + { + TradSet_ = true; + Trad_.read(dict); + } + else + { + TradSet_ = false; + } + + Icl_ = dimensionedScalar(Icl_.dimensions(), 0.155)*clothing_; + + fcl_.value() = + Icl_.value() <= 0.078 + ? 1.0 + 1.290*Icl_.value() + : 1.05 + 0.645*Icl_.value(); + + return true; +} + + +bool Foam::functionObjects::comfort::execute() +{ + const dimensionedScalar Trad(this->Trad()); + const volScalarField pSat(this->pSat()); + + const dimensionedScalar metabolicRateSI(58.15*metabolicRate_); + const dimensionedScalar extWorkSI(58.15*extWork_); + + const volScalarField& T = lookupObject("T"); + + // Heat transfer coefficient [W/m^2/K] + // This field is updated in Tcloth() + volScalarField hc + ( + IOobject + ( + "hc", + mesh_.time().timeName(), + mesh_ + ), + mesh_, + dimensionedScalar(dimensionSet(1, 0, -3, -1, 0, 0, 0), 0) + ); + + // Calculate the surface temperature of the cloth by an iterative + // process using equation (2) from DIN EN ISO 7730 [degC] + const volScalarField Tcloth + ( + this->Tcloth + ( + hc, + metabolicRateSI, + extWorkSI, + T, + Trad + ) + ); + + // Calculate the PMV quantity + const dimensionedScalar factor1(dimensionSet(-1, 0, 3, 0, 0, 0, 0), 0.303); + const dimensionedScalar factor2 + ( + dimless/metabolicRateSI.dimensions(), + -0.036 + ); + const dimensionedScalar factor3(factor1.dimensions(), 0.028); + const dimensionedScalar factor4(dimLength/dimTime, 3.05e-3); + const dimensionedScalar factor5(dimPressure, 5733); + const dimensionedScalar factor6(dimTime/dimLength, 6.99); + const dimensionedScalar factor8(metabolicRateSI.dimensions(), 58.15); + const dimensionedScalar factor9(dimless/dimPressure, 1.7e-5); + const dimensionedScalar factor10(dimPressure, 5867); + const dimensionedScalar factor11(dimless/dimTemperature, 0.0014); + const dimensionedScalar factor12(dimTemperature, 307.15); + const dimensionedScalar factor13 + ( + dimensionSet(1, 0, -3, -4, 0, 0, 0), + 3.96e-8 + ); + + const scalar factor7 + ( + // Special treatment of Term4 + // if metaRate - extWork < factor8, set to zero + (metabolicRateSI - extWorkSI).value() < factor8.value() ? 0 : 0.42 + ); + + Info<< "Calculating the predicted mean vote (PMV)\n"; + + // Equation (1) + tmp PMV + ( + volScalarField::New + ( + "PMV", + + // Term1: Thermal sensation transfer coefficient + (factor1*exp(factor2*metabolicRateSI) + factor3) + *( + (metabolicRateSI - extWorkSI) + + // Term2: Heat loss difference through skin + - factor4 + *( + factor5 + - factor6*(metabolicRateSI - extWorkSI) + - pSat*relHumidity_ + ) + + // Term3: Heat loss through sweating + - factor7*(metabolicRateSI - extWorkSI - factor8) + + // Term4: Heat loss through latent respiration + - factor9*metabolicRateSI*(factor10 - pSat*relHumidity_) + + // Term5: Heat loss through dry respiration + - factor11*metabolicRateSI*(factor12 - T) + + // Term6: Heat loss through radiation + - factor13*fcl_*(pow4(Tcloth) - pow4(Trad)) + + // Term7: Heat loss through convection + - fcl_*hc*(Tcloth - T) + ) + ) + ); + + Info<< "Calculating the predicted percentage of dissatisfaction (PPD)\n"; + + // Equation (5) in EN ISO + tmp PPD + ( + volScalarField::New + ( + "PPD", + 100 - 95*exp(-0.03353*pow4(PMV()) - 0.21790*sqr(PMV())) + ) + ); + + return store(PMV) && store(PPD); +} + + +bool Foam::functionObjects::comfort::write() +{ + return writeObject("PMV") && writeObject("PPD"); +} + + +// ************************************************************************* // diff --git a/src/functionObjects/field/comfort/comfort.H b/src/functionObjects/field/comfort/comfort.H new file mode 100644 index 0000000000..e6d41a3047 --- /dev/null +++ b/src/functionObjects/field/comfort/comfort.H @@ -0,0 +1,206 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Copyright (C) 2019 OpenFOAM Foundation + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software: you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM. If not, see . + +Class + Foam::functionObjects::comfort + +Description + Calculates the thermal comfort quantities predicted mean vote (PMV) and + predicted percentage of dissatisfaction (PPD) based on DIN ISO EN 7730:2005. + +Usage + \table + Property | Description | Required | Default value + clothing | The insulation value of the cloth | no | 0 + metabolicRate | The metabolic rate | no | 0.8 + extWork | The external work | no | 0 + Trad | Radiation temperature | no | -1 + relHumidity | Relative humidity of the air | no | 50 + pSat | Saturation pressure of water | no | -1 + tolerance | Residual control for the cloth temperature | no | 1e-5 + maxClothIter | Maximum number of iterations | no | 0 + meanVelocity | Use a constant mean velocity in the whole domain | no |\ + false + \endtable + + \table + Predicted Mean Vote (PMV) | evaluation + + 3 | hot + + 2 | warm + + 1 | slightly warm + + 0 | neutral + - 1 | slightly cool + - 2 | cool + - 3 | cold + \endtable + + \verbatim + comfortAnalysis + { + type comfort; + libs ("libfieldFunctionObjects.so"); + + executeControl writeTime; + writeControl writeTime; + } + \endverbatim + +SourceFiles + comfort.C + +\*---------------------------------------------------------------------------*/ + +#ifndef comfort_H +#define comfort_H + +#include "fvMeshFunctionObject.H" +#include "volFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ +namespace functionObjects +{ + +/*---------------------------------------------------------------------------*\ + Class comfort Declaration +\*---------------------------------------------------------------------------*/ + +class comfort +: + public fvMeshFunctionObject +{ + // Private Data + + //- Clothing [-] + dimensionedScalar clothing_; + + //- Metabolic rate [kg/s^3] + dimensionedScalar metabolicRate_; + + //- External work [kg/s^3] + dimensionedScalar extWork_; + + //- Switch set to true if the radiation temperature is provided + Switch TradSet_; + + //- Mean radiation temperature [K] + dimensionedScalar Trad_; + + //- Relative humidity [percentage] + dimensionedScalar relHumidity_; + + //- Saturation pressure of water [Pa] + dimensionedScalar pSat_; + + //- Thermal insulation of clothing [W/m^2/K] + dimensionedScalar Icl_; + + //- Prefactor of cloth area [-] + dimensionedScalar fcl_; + + //- Tolerance criteria for iterative process to find Tcl + scalar tolerance_; + + //- Maximum number of correctors for cloth temperature + int maxClothIter_; + + //- Switch to use volume weighted velocity field for caluclation + Switch meanVelocity_; + + + // Private Member Functions + + // Calculate the magnitude of the velocity [m/s] + tmp magU() const; + + // Calculate the radiation temperature in the domain using a simple + // approach [K] + dimensionedScalar Trad() const; + + // Calculate the saturation pressure based on 7730:2006 + // Possible options: adding different calculation methods such as + // the formulation based on Magnus or others [Pa] + tmp pSat() const; + + // Calculate and return the surface temperature of the cloth [K] + // and the heat transfer coefficient hc [W/m^2/K] + tmp Tcloth + ( + volScalarField& hc, + const dimensionedScalar& metabolicRateSI, + const dimensionedScalar& extWorkSI, + const volScalarField& TdegC, + const dimensionedScalar& Trad + ); + + // Return true if the cloth temperature iteration has converged + bool converged(const volScalarField&) const; + + +public: + + //- Runtime type information + TypeName("comfort"); + + + // Constructors + + //- Construct from Time and dictionary + comfort + ( + const word& name, + const Time& runTime, + const dictionary& dict + ); + + + //- Destructor + virtual ~comfort(); + + + // Member Functions + + //- Read the data needed for the comfort calculation + virtual bool read(const dictionary&); + + //- Calculate the predicted mean vote (PMV) + // and predicted percentage dissatisfaction (PPD) fields + virtual bool execute(); + + //- Write the PPD and PMV fields + virtual bool write(); +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace functionObjects +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/T b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/T new file mode 100644 index 0000000000..bade247ec2 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/T @@ -0,0 +1,42 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volScalarField; + object T; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 0 0 1 0 0 0]; + +internalField uniform 302; + +boundaryField +{ + walls + { + type fixedValue; + value $internalField; + } + + inlet + { + type fixedValue; + value uniform 290; + } + + outlet + { + type zeroGradient; + } +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/U b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/U new file mode 100644 index 0000000000..d90ce79d2d --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/U @@ -0,0 +1,42 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volVectorField; + object U; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 1 -1 0 0 0 0]; + +internalField uniform (0 0 0); + +boundaryField +{ + walls + { + type noSlip; + } + + inlet + { + type fixedValue; + value uniform (0.2 0 0); + } + + outlet + { + type pressureInletOutletVelocity; + value $internalField; + } +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/alphat b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/alphat new file mode 100644 index 0000000000..d41b25109f --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/alphat @@ -0,0 +1,44 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volScalarField; + object alphat; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [1 -1 -1 0 0 0 0]; + +internalField uniform 0; + +boundaryField +{ + walls + { + type compressible::alphatJayatillekeWallFunction; + Prt 0.85; + value $internalField; + } + + inlet + { + type calculated; + value $internalField; + } + + outlet + { + type calculated; + value $internalField; + } +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/epsilon b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/epsilon new file mode 100644 index 0000000000..54ccea40e8 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/epsilon @@ -0,0 +1,45 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volScalarField; + location "0"; + object epsilon; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 2 -3 0 0 0 0]; + +internalField uniform 0.23; + +boundaryField +{ + walls + { + type epsilonWallFunction; + value $internalField; + } + + inlet + { + type turbulentMixingLengthDissipationRateInlet; + mixingLength 0.0168; + value $internalField; + } + + outlet + { + type inletOutlet; + inletValue $internalField; + } +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/k b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/k new file mode 100644 index 0000000000..341479b846 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/k @@ -0,0 +1,45 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volScalarField; + location "0"; + object k; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 2 -2 0 0 0 0]; + +internalField uniform 8e-2; + +boundaryField +{ + walls + { + type kqRWallFunction; + value $internalField; + } + + inlet + { + type turbulentIntensityKineticEnergyInlet; + intensity 0.14; + value $internalField; + } + + outlet + { + type inletOutlet; + inletValue $internalField; + } +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/nut b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/nut new file mode 100644 index 0000000000..e9d2679dbb --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/nut @@ -0,0 +1,43 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volScalarField; + object nut; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 2 -1 0 0 0 0]; + +internalField uniform 0; + +boundaryField +{ + walls + { + type nutkWallFunction; + value $internalField; + } + + inlet + { + type calculated; + value $internalField; + } + + outlet + { + type calculated; + value $internalField; + } +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/p b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/p new file mode 100644 index 0000000000..decdf102ea --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/p @@ -0,0 +1,42 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volScalarField; + object p; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [1 -1 -2 0 0 0 0]; + +internalField uniform 1e5; + +boundaryField +{ + walls + { + type calculated; + value $internalField; + } + + inlet + { + type calculated; + value $internalField; + } + + outlet + { + type calculated; + value $internalField; + } +} + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/p_rgh b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/p_rgh new file mode 100644 index 0000000000..dad2174915 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/0/p_rgh @@ -0,0 +1,43 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class volScalarField; + object p_rgh; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [1 -1 -2 0 0 0 0]; + +internalField uniform 1e5; + +boundaryField +{ + walls + { + type fixedFluxPressure; + value $internalField; + } + + outlet + { + type prghPressure; + p $internalField; + value $internalField; + } + + inlet + { + type fixedFluxPressure; + value $internalField; + } +} + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/Allrun b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/Allrun new file mode 100755 index 0000000000..edee003c99 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/Allrun @@ -0,0 +1,15 @@ +#!/bin/sh +cd ${0%/*} || exit 1 # Run from this directory + +# Source tutorial run functions +. $WM_PROJECT_DIR/bin/tools/RunFunctions + +# Get application name +application=$(getApplication) + +runApplication blockMesh +runApplication topoSet +runApplication createPatch -overwrite +runApplication $application + +#------------------------------------------------------------------------------ diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/g b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/g new file mode 100644 index 0000000000..f7279e4c5c --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/g @@ -0,0 +1,21 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class uniformDimensionedVectorField; + location "constant"; + object g; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +dimensions [0 1 -2 0 0 0 0]; +value (0 0 -9.81); + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/thermophysicalProperties b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/thermophysicalProperties new file mode 100644 index 0000000000..5a81d7afe8 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/thermophysicalProperties @@ -0,0 +1,54 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + location "constant"; + object thermophysicalProperties; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +thermoType +{ + type heRhoThermo; + mixture pureMixture; + transport const; + thermo eConst; + equationOfState Boussinesq; + specie specie; + energy sensibleInternalEnergy; +} + +mixture +{ + specie + { + molWeight 28.9; + } + equationOfState + { + rho0 1; + T0 300; + beta 3e-03; + } + thermodynamics + { + Cv 712; + Hf 0; + } + transport + { + mu 1e-05; + Pr 0.7; + } +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/turbulenceProperties b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/turbulenceProperties new file mode 100644 index 0000000000..a652d98461 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/constant/turbulenceProperties @@ -0,0 +1,30 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + location "constant"; + object turbulenceProperties; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +simulationType RAS; + +RAS +{ + RASModel kEpsilon; + + turbulence on; + + printCoeffs on; +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/blockMeshDict b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/blockMeshDict new file mode 100644 index 0000000000..136132ac3f --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/blockMeshDict @@ -0,0 +1,49 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object blockMeshDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +convertToMeters 1; + +vertices +( + (0 0 0) + (0 0 1.6) + (0 3 1.6) + (0 3 0) + + (4 0 0) + (4 0 1.6) + (4 3 1.6) + (4 3 0) +); + +blocks +( + hex (0 3 2 1 4 7 6 5) + (40 20 60) + simpleGrading (1 1 1) +); + +defaultPatch +{ + name walls; + type wall; +} + +boundary +(); + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/controlDict b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/controlDict new file mode 100644 index 0000000000..9d2661b2cf --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/controlDict @@ -0,0 +1,76 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + location "system"; + object controlDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +application buoyantSimpleFoam; + +startFrom startTime; + +startTime 0; + +stopAt endTime; + +endTime 3000; + +deltaT 1; + +writeControl timeStep; + +writeInterval 100; + +purgeWrite 0; + +writeFormat binary; + +writePrecision 6; + +writeCompression off; + +timeFormat general; + +timePrecision 6; + +runTimeModifiable true; + +functions +{ + age + { + libs ("libfieldFunctionObjects.so"); + type age; + + diffusion on; + + writeControl writeTime; + executeControl writeTime; + } + + comfort + { + libs ("libfieldFunctionObjects.so"); + type comfort; + + clothing 0.5; + metabolicRate 1.2; + extWork 0; + relHumidity 60; + + writeControl writeTime; + executeControl writeTime; + } +} + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/createPatchDict b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/createPatchDict new file mode 100644 index 0000000000..5a81c48388 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/createPatchDict @@ -0,0 +1,47 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object createPatchDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +pointSync false; + +writeCyclicMatch false; + +patches +( + { + name inlet; + + patchInfo + { + type patch; + } + + constructFrom set; + set inlet; + } + { + name outlet; + + patchInfo + { + type patch; + } + + constructFrom set; + set outlet; + } +); + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/fvSchemes b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/fvSchemes new file mode 100644 index 0000000000..7a7c724f17 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/fvSchemes @@ -0,0 +1,61 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + location "system"; + object fvSchemes; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +ddtSchemes +{ + default steadyState; +} + +gradSchemes +{ + default Gauss linear; +} + +divSchemes +{ + default none; + + div(phi,U) bounded Gauss upwind; + div(phi,e) bounded Gauss upwind; + + div(phi,k) bounded Gauss upwind; + div(phi,epsilon) bounded Gauss upwind; + + div(phi,Ekp) bounded Gauss linear; + + div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear; + + div(phi,age) bounded Gauss upwind; +} + +laplacianSchemes +{ + default Gauss linear orthogonal; +} + +interpolationSchemes +{ + default linear; +} + +snGradSchemes +{ + default orthogonal; +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/fvSolution b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/fvSolution new file mode 100644 index 0000000000..44cf98e9c5 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/fvSolution @@ -0,0 +1,74 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + location "system"; + object fvSolution; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +solvers +{ + p_rgh + { + solver PCG; + preconditioner DIC; + tolerance 1e-8; + relTol 0.01; + } + + "(U|e|k|epsilon)" + { + solver PBiCGStab; + preconditioner DILU; + tolerance 1e-7; + relTol 0.1; + } + + age + { + $U; + relTol 0.001; + } +} + +SIMPLE +{ + nNonOrthogonalCorrectors 0; + + residualControl + { + p_rgh 1e-2; + U 1e-4; + e 1e-2; + + "(k|epsilon|omega)" 1e-3; + } +} + +relaxationFactors +{ + fields + { + p_rgh 0.7; + } + + equations + { + U 0.2; + e 0.1; + "(k|epsilon|R)" 0.7; + age 1; + } +} + + +// ************************************************************************* // diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/topoSetDict b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/topoSetDict new file mode 100644 index 0000000000..2a438d4c09 --- /dev/null +++ b/tutorials/heatTransfer/buoyantSimpleFoam/comfortHotRoom/system/topoSetDict @@ -0,0 +1,41 @@ +/*--------------------------------*- C++ -*----------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | Website: https://openfoam.org + \\ / A nd | Version: dev + \\/ M anipulation | +\*---------------------------------------------------------------------------*/ +FoamFile +{ + version 2.0; + format ascii; + class dictionary; + object topoSetDict; +} +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +actions +( + { + name inlet; + type faceSet; + action new; + source boxToFace; + sourceInfo + { + box (-0.001 0.25 1.1)(0.001 0.75 1.3); + } + } + { + name outlet; + type faceSet; + action new; + source boxToFace; + sourceInfo + { + box (1.75 2.999 0.3)(2.25 3.001 0.5); + } + } +); + +// ************************************************************************* //