diff --git a/src/atmosphericModels/Make/files b/src/atmosphericModels/Make/files
index 62808d2debd74b05c7fe23be758d18ba691df0a0..0b61ec6284890e6b1788ecda7e31f9cee02f152c 100644
--- a/src/atmosphericModels/Make/files
+++ b/src/atmosphericModels/Make/files
@@ -35,6 +35,13 @@ fvOptions/atmPlantCanopyTSource/atmPlantCanopyTSource.C
fvOptions/atmNutSource/atmNutSource.C
fvOptions/treeTurbulence/treeTurbulence.C
+etht = fvOptions/evapotranspirationHeatTransfer
+ethtModels = $(etht)/evapotranspirationHeatTransferModels
+$(etht)/evapotranspirationHeatTransfer.C
+$(ethtModels)/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModel.C
+$(ethtModels)/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModelNew.C
+$(ethtModels)/tree/tree.C
+$(ethtModels)/grass/grass.C
/* functionObjects */
functionObjects/ObukhovLength/ObukhovLength.C
diff --git a/src/atmosphericModels/Make/options b/src/atmosphericModels/Make/options
index d2af502b1a05d1f68710f13284c2247eb8283f05..033971d12007a7133181d53bf5b718aba68c6dc5 100644
--- a/src/atmosphericModels/Make/options
+++ b/src/atmosphericModels/Make/options
@@ -8,6 +8,7 @@ EXE_INC = \
-I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-I$(LIB_SRC)/thermophysicalModels/solidThermo/lnInclude \
+ -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
-I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
@@ -18,6 +19,7 @@ LIB_LIBS = \
-lfvOptions \
-lmeshTools \
-lsurfMesh \
+ -lradiationModels \
-lfluidThermophysicalModels \
-lsolidThermo \
-lturbulenceModels \
diff --git a/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransfer.C b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransfer.C
new file mode 100644
index 0000000000000000000000000000000000000000..3c76ebe6ac45d471ebff668e317f133caef6ce79
--- /dev/null
+++ b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransfer.C
@@ -0,0 +1,172 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2022 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "evapotranspirationHeatTransfer.H"
+#include "evapotranspirationHeatTransferModel.H"
+#include "fvMesh.H"
+#include "fvMatrix.H"
+#include "basicThermo.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace fv
+{
+ defineTypeNameAndDebug(evapotranspirationHeatTransfer, 0);
+ addToRunTimeSelectionTable
+ (
+ option,
+ evapotranspirationHeatTransfer,
+ dictionary
+ );
+}
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::fv::evapotranspirationHeatTransfer::evapotranspirationHeatTransfer
+(
+ const word& sourceName,
+ const word& modelType,
+ const dictionary& dict,
+ const fvMesh& mesh
+)
+:
+ fv::cellSetOption(sourceName, modelType, dict, mesh),
+ ethtModelPtr_()
+{
+ // Set the field name to that of the energy
+ // field from which the temperature is obtained
+
+ const auto& thermo = mesh_.lookupObject<basicThermo>(basicThermo::dictName);
+
+ fieldNames_.resize(1, thermo.he().name());
+
+ fv::option::resetApplied();
+
+ Info<< " Applying evapotranspirationHeatTransfer to: " << fieldNames_[0]
+ << endl;
+
+ read(dict);
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::fv::evapotranspirationHeatTransfer::~evapotranspirationHeatTransfer()
+{} // evapotranspirationHeatTransferModel was forward declared
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+void Foam::fv::evapotranspirationHeatTransfer::addSup
+(
+ fvScalarMatrix& eqn,
+ const label fieldi
+)
+{
+ if (this->V() < VSMALL)
+ {
+ return;
+ }
+
+ const scalar V = this->V();
+ const scalarField& Vcells = mesh_.V();
+
+ const volScalarField& he = eqn.psi();
+ scalarField& heSource = eqn.source();
+
+ // Calculate evapotranspiration heat transfer rate per volume [J/s/m^3]
+ const scalarField Q(ethtModelPtr_->Q(cells_)/V);
+
+ if (he.dimensions() == dimEnergy/dimMass)
+ {
+ forAll(cells_, i)
+ {
+ const label celli = cells_[i];
+
+ heSource[celli] += Q[i]*Vcells[celli];
+ }
+ }
+ else if (he.dimensions() == dimTemperature)
+ {
+ const auto& thermo =
+ mesh_.lookupObject<basicThermo>(basicThermo::dictName);
+
+ // Calculate density*heat capacity at constant pressure/volume
+ const volScalarField rhoCpv(thermo.rho()*thermo.Cpv());
+
+ // heSource [K m^3/s] = [J/s/m^3] * m^3 / [kg/m^3] / [J/kg/K]
+ forAll(cells_, i)
+ {
+ const label celli = cells_[i];
+
+ heSource[celli] += Q[i]*Vcells[celli]*rhoCpv[i];
+ }
+ }
+}
+
+
+void Foam::fv::evapotranspirationHeatTransfer::addSup
+(
+ const volScalarField& rho,
+ fvScalarMatrix& eqn,
+ const label fieldi
+)
+{
+ addSup(eqn, fieldi);
+}
+
+
+bool Foam::fv::evapotranspirationHeatTransfer::read(const dictionary& dict)
+{
+ if (!fv::cellSetOption::read(dict))
+ {
+ return false;
+ }
+
+ if (selectionMode_ != selectionModeType::smCellZone)
+ {
+ FatalIOErrorInFunction(dict)
+ << "evapotranspirationHeatTransfer requires "
+ << selectionModeTypeNames_[selectionModeType::smCellZone]
+ << exit(FatalIOError);
+ }
+
+ ethtModelPtr_.reset
+ (
+ evapotranspirationHeatTransferModel::New(dict, mesh_)
+ );
+
+ return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransfer.H b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransfer.H
new file mode 100644
index 0000000000000000000000000000000000000000..95740bb74b7341974b93ed9cc264b8c3fbf8df56
--- /dev/null
+++ b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransfer.H
@@ -0,0 +1,183 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2022 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::fv::evapotranspirationHeatTransfer
+
+Description
+ Applies sources on temperature (\c T - incompressible)
+ or energy (\c h/e - compressible) equation to incorporate
+ evapotranspiration heat-transfer effects from the specified
+ plant canopy. Heat transfer is usually calculated based on
+ empirical relations between plants and solar radiation.
+
+Usage
+ Example by using \c constant/fvOptions:
+ \verbatim
+ evapotranspirationHeatTransfer1
+ {
+ // Mandatory entries
+ type evapotranspirationHeatTransfer;
+ model <word>;
+
+ // Model-specific entries
+ ...
+
+ // Inherited entries
+ ...
+ }
+ \endverbatim
+
+ where the entries mean:
+ \table
+ Property | Description | Type | Reqd | Deflt
+ type | Type name: evapotranspirationHeatTransfer | word | yes | -
+ model | Name of heat-transfer model | word | yes | -
+ \endtable
+
+ Options for the \c model entry:
+ \verbatim
+ tree | Regression-based heat-transfer model for trees
+ grass | Heat-transfer model for grass
+ \endverbatim
+
+ The inherited entries are elaborated in:
+ - \link fvOption.H \endlink
+ - \link cellSetOption.H \endlink
+
+Note
+ - Evapotranspiration mass transfer is not modelled.
+
+SourceFiles
+ evapotranspirationHeatTransfer.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Foam_fv_evapotranspirationHeatTransfer_H
+#define Foam_fv_evapotranspirationHeatTransfer_H
+
+#include "cellSetOption.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward Declarations
+class evapotranspirationHeatTransferModel;
+
+namespace fv
+{
+
+/*---------------------------------------------------------------------------*\
+ Class evapotranspirationHeatTransfer Declaration
+\*---------------------------------------------------------------------------*/
+
+class evapotranspirationHeatTransfer
+:
+ public fv::cellSetOption
+{
+ // Private Data
+
+ //- Runtime-selectable evapotranspiration heat transfer model
+ autoPtr<evapotranspirationHeatTransferModel> ethtModelPtr_;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("evapotranspirationHeatTransfer");
+
+
+ // Constructors
+
+ //- Construct from explicit source name and mesh
+ evapotranspirationHeatTransfer
+ (
+ const word& sourceName,
+ const word& modelType,
+ const dictionary& dict,
+ const fvMesh& mesh
+ );
+
+ //- No copy construct
+ evapotranspirationHeatTransfer(const evapotranspirationHeatTransfer&) =
+ delete;
+
+ //- No copy assignment
+ void operator=(const evapotranspirationHeatTransfer&) = delete;
+
+
+ //- Destructor
+ virtual ~evapotranspirationHeatTransfer();
+
+
+ // Member Functions
+
+ //- Add explicit contribution to energy equation
+ //- for incompressible flow computations
+ virtual void addSup
+ (
+ fvScalarMatrix& eqn,
+ const label fieldi
+ );
+
+ //- Add explicit contribution to energy equation
+ //- for compressible flow computations
+ virtual void addSup
+ (
+ const volScalarField& rho,
+ fvScalarMatrix& eqn,
+ const label fieldi
+ );
+
+ //- Add explicit contribution to energy equation
+ //- for multiphase flow computations
+ virtual void addSup
+ (
+ const volScalarField& alpha,
+ const volScalarField& rho,
+ fvScalarMatrix& eqn,
+ const label fieldi
+ )
+ {
+ NotImplemented;
+ }
+
+ //- Read dictionary
+ virtual bool read(const dictionary& dict);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace fv
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModel.C b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModel.C
new file mode 100644
index 0000000000000000000000000000000000000000..5051e01f80e82b15ef61f8e34b8c3f3a2bc8546a
--- /dev/null
+++ b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModel.C
@@ -0,0 +1,196 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2022 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "evapotranspirationHeatTransferModel.H"
+#include "radiationModel.H"
+#include "solarLoadBase.H"
+#include "fvMesh.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+ defineTypeNameAndDebug(evapotranspirationHeatTransferModel, 0);
+ defineRunTimeSelectionTable
+ (
+ evapotranspirationHeatTransferModel,
+ dictionary
+ );
+}
+
+
+// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
+
+Foam::volScalarField& Foam::evapotranspirationHeatTransferModel::getOrReadField
+(
+ const word& fieldName
+) const
+{
+ auto* ptr = mesh_.getObjectPtr<volScalarField>(fieldName);
+
+ if (!ptr)
+ {
+ ptr = new volScalarField
+ (
+ IOobject
+ (
+ fieldName,
+ mesh_.time().timeName(),
+ mesh_,
+ IOobject::MUST_READ,
+ IOobject::AUTO_WRITE
+ ),
+ mesh_
+ );
+ mesh_.objectRegistry::store(ptr);
+ }
+
+ return *ptr;
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModel::q() const
+{
+ // Retrieve solar-load information
+ const auto& base =
+ mesh().lookupObject<radiation::solarLoadBase>("solarLoadBase");
+ const solarCalculator& solar = base.solarCalculatorRef();
+
+ // Retrieve internal & boundary faces directly hit by solar rays
+ const faceShading& shade = base.faceShadingRef();
+ const labelList& hitFacesId = shade.rayStartFaces();
+
+ // Retrieve face zone information
+ const faceZone& zone = mesh_.faceZones()[faceZoneId_];
+ const vectorField& faceNormals = zone().faceNormals();
+ const scalarField& faceAreas = zone().magFaceAreas();
+
+ // Retrieve direct solar load [W/m^2]
+ const vector directSolarRad(solar.directSolarRad()*solar.direction());
+
+ // Identify zone faces within mesh
+ bitSet isZoneFace(mesh_.nFaces());
+ isZoneFace.set(zone);
+
+ // Calculate area-averaged incident solar load
+ // Assuming hit faces are updated by solarLoad
+ scalar q = 0;
+ scalar totalFaceArea = 0;
+
+ for (const label facei : hitFacesId)
+ {
+ if (isZoneFace[facei])
+ {
+ const label localFacei = zone.whichFace(facei);
+
+ const vector& faceNormal = faceNormals[localFacei];
+ const scalar faceArea = faceAreas[localFacei];
+
+ const scalar qIncident = directSolarRad & faceNormal;
+
+ q += qIncident*faceArea;
+ totalFaceArea += faceArea;
+ }
+ }
+ reduce(q, sumOp<scalar>());
+ reduce(totalFaceArea, sumOp<scalar>());
+ q /= totalFaceArea;
+
+
+ // Sum diffusive solar loads
+ switch(solar.sunLoadModel())
+ {
+ case solarCalculator::mSunLoadFairWeatherConditions:
+ case solarCalculator::mSunLoadTheoreticalMaximum:
+ {
+ // Calculate diffusive radiance
+ // contribution from sky and ground
+ tmp<scalarField> tdiffuseSolarRad =
+ solar.diffuseSolarRad(faceNormals);
+ const scalarField& diffuseSolarRad = tdiffuseSolarRad.cref();
+
+ // Calculate area-averaged diffusive solar load
+ scalar meanDiffuseSolarRad = 0;
+ forAll(faceAreas, i)
+ {
+ meanDiffuseSolarRad += diffuseSolarRad[i]*faceAreas[i];
+ }
+ meanDiffuseSolarRad /= totalFaceArea;
+
+ q += meanDiffuseSolarRad;
+ }
+ break;
+
+ case solarCalculator::mSunLoadConstant:
+ case solarCalculator::mSunLoadTimeDependent:
+ {
+ q += solar.diffuseSolarRad();
+ }
+ break;
+ }
+
+ return q;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::evapotranspirationHeatTransferModel::evapotranspirationHeatTransferModel
+(
+ const dictionary& dict,
+ const fvMesh& mesh
+)
+:
+ mesh_(mesh)
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
+
+bool Foam::evapotranspirationHeatTransferModel::read(const dictionary& dict)
+{
+ faceZoneId_ = mesh_.faceZones().findZoneID(dict.get<word>("faceZone"));
+
+ if (faceZoneId_ < 0)
+ {
+ const word faceZoneName(mesh_.faceZones()[faceZoneId_].name());
+
+ FatalIOErrorInFunction(dict)
+ << type() << ' ' << typeName << ": "
+ << " No matching face zone: " << faceZoneName << nl
+ << " Known face zones: "
+ << flatOutput(mesh_.faceZones().names()) << nl
+ << exit(FatalIOError);
+
+ return false;
+ }
+
+ return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModel.H b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModel.H
new file mode 100644
index 0000000000000000000000000000000000000000..1f34f91c945575132ff7943f088988a90e6b5f50
--- /dev/null
+++ b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModel.H
@@ -0,0 +1,180 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2022 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::evapotranspirationHeatTransferModel
+
+Description
+ Base class for evapotranspiration models
+ to handle general evapotranspiration characteristics.
+
+SourceFiles
+ evapotranspirationHeatTransferModel.C
+ evapotranspirationHeatTransferModelNew.C
+ evapotranspirationHeatTransferModelTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Foam_evapotranspirationHeatTransferModel_H
+#define Foam_evapotranspirationHeatTransferModel_H
+
+#include "dictionary.H"
+#include "HashSet.H"
+#include "volFields.H"
+#include "runTimeSelectionTables.H"
+#include "OFstream.H"
+#include "coordinateSystem.H"
+#include "writeFile.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward Declarations
+class fvMesh;
+
+/*---------------------------------------------------------------------------*\
+ Class evapotranspirationHeatTransferModel Declaration
+\*---------------------------------------------------------------------------*/
+
+class evapotranspirationHeatTransferModel
+{
+ // Private Data
+
+ //- Reference to the mesh
+ const fvMesh& mesh_;
+
+ //- Identifier of specified face zone
+ label faceZoneId_;
+
+
+protected:
+
+ // Protected Member Functions
+
+ //- Return requested field from the object registry
+ //- or read+register the field to the object registry
+ volScalarField& getOrReadField(const word& fieldName) const;
+
+ //- Area-averaged heat flux due to short-wave
+ //- solar rays on face zone [W/m^2]
+ // Short-wave: direct and diffusive solar rays
+ scalar q() const;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("evapotranspirationHeatTransferModel");
+
+
+ // Declare runtime constructor selection table
+
+ declareRunTimeSelectionTable
+ (
+ autoPtr,
+ evapotranspirationHeatTransferModel,
+ dictionary,
+ (
+ const dictionary& dict,
+ const fvMesh& mesh
+ ),
+ (dict, mesh)
+ );
+
+
+ // Selectors
+
+ //- Return a reference to the selected evapotranspiration model
+ static autoPtr<evapotranspirationHeatTransferModel> New
+ (
+ const dictionary& dict,
+ const fvMesh& mesh
+ );
+
+
+ // Constructors
+
+ //- Construct from components
+ evapotranspirationHeatTransferModel
+ (
+ const dictionary& dict,
+ const fvMesh& mesh
+ );
+
+ //- No copy construct
+ evapotranspirationHeatTransferModel
+ (
+ const evapotranspirationHeatTransferModel&
+ ) = delete;
+
+ //- No copy assignment
+ void operator=(const evapotranspirationHeatTransferModel&) = delete;
+
+
+ //- Destructor
+ virtual ~evapotranspirationHeatTransferModel() = default;
+
+
+ // Member Functions
+
+ // Access
+
+ //- Return const reference to the mesh
+ const fvMesh& mesh() const noexcept
+ {
+ return mesh_;
+ }
+
+ //- Return the face-zone identifier
+ label faceZoneId() const noexcept
+ {
+ return faceZoneId_;
+ }
+
+
+ // Evaluation
+
+ //- Return heat-transfer rate for speficied cells [J/s]
+ virtual tmp<scalarField> Q(const labelList& cells) const = 0;
+
+
+ // I-O
+
+ //- Read the dictionary
+ virtual bool read(const dictionary& dict);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModelNew.C b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModelNew.C
new file mode 100644
index 0000000000000000000000000000000000000000..cf8f876d9e68412914413fd64f202d4b96a280d5
--- /dev/null
+++ b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/evapotranspirationHeatTransferModel/evapotranspirationHeatTransferModelNew.C
@@ -0,0 +1,62 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2022 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "evapotranspirationHeatTransferModel.H"
+#include "fvMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::evapotranspirationHeatTransferModel>
+Foam::evapotranspirationHeatTransferModel::New
+(
+ const dictionary& dict,
+ const fvMesh& mesh
+)
+{
+ word modelType(dict.get<word>("model"));
+
+ auto cstrIter = dictionaryConstructorTablePtr_->cfind(modelType);
+
+ if (!cstrIter.found())
+ {
+ FatalIOErrorInLookup
+ (
+ dict,
+ "model",
+ modelType,
+ *dictionaryConstructorTablePtr_
+ ) << exit(FatalIOError);
+ }
+
+ return autoPtr<evapotranspirationHeatTransferModel>
+ (
+ cstrIter()(dict, mesh)
+ );
+}
+
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/grass/grass.C b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/grass/grass.C
new file mode 100644
index 0000000000000000000000000000000000000000..c22c1d12978144c9559cea6260541bd4381eb905
--- /dev/null
+++ b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/grass/grass.C
@@ -0,0 +1,418 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2022 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "grass.H"
+#include "basicThermo.H"
+#include "fluidThermo.H"
+#include "turbulentTransportModel.H"
+#include "turbulentFluidThermoModel.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace evapotranspirationHeatTransferModels
+{
+ defineTypeNameAndDebug(grass, 0);
+ addToRunTimeSelectionTable
+ (
+ evapotranspirationHeatTransferModel,
+ grass,
+ dictionary
+ );
+}
+}
+
+
+const Foam::Enum
+<
+ Foam::evapotranspirationHeatTransferModels::grass::soilHeatFluxType
+>
+Foam::evapotranspirationHeatTransferModels::grass::soilHeatFluxTypeNames
+({
+ {
+ soilHeatFluxType::PROPORTIONAL_TO_SOLAR_RADIATION,
+ "proportionalToSolarRadiation"
+ },
+ { soilHeatFluxType::BOUNDARY, "boundary" }
+});
+
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+Foam::tmp<Foam::scalarField>
+Foam::evapotranspirationHeatTransferModels::grass::E(const labelList& cells)
+const
+{
+ const scalar q = this->q();
+ const scalar G = this->G(cells);
+ const scalar Delta = this->Delta();
+ const scalar D = this->D();
+ const scalar ra = this->ra();
+ const scalar rs = this->rs();
+ const scalar gamma = this->gamma();
+
+ // (BSG:Eq. 11)
+ const scalar E =
+ (Delta*(q - G) + rho_*Cp_*D/ra)/(Delta + gamma*(1 + rs/ra));
+
+ return tmp<scalarField>::New(cells.size(), E);
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::grass::Delta() const
+{
+ // (BSG:Eq. 15), (APR:Eq. 13) - note 0.6108 -> 610.8 due to kPa -> Pa
+ const scalar Ta = Tref_ + 237.3;
+
+ return 4098*(610.8*exp(17.27*Tref_/Ta))/sqr(Ta);
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::grass::D() const
+{
+ // (BSG:Eq. 16)
+ const scalar RHref = RHptr_->value(mesh().time().timeOutputValue());
+
+ return (1 - RHref)*pSat();
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::grass::pSat() const
+{
+ // (BSG:Eq. 17; Arden Buck equation)
+ const scalar p1 = (1.0007 + 3.46e-8*pAtm_)*611.21;
+
+ return p1*exp((18.678 - Tref_/234.5)*Tref_/(Tref_ + 257.14));
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::grass::gamma() const
+{
+ // (APR:Eq. 8)
+ return Cp_*pAtm_/(epsilon_*lambda());
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::grass::lambda() const
+{
+ // (BSG:Eq. 12)
+ const scalar T1 = Tref_ + 273.15;
+ const scalar T2 = Tref_ + 239.24;
+
+ return 1.91846e6*sqr(T1/T2);
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::grass::ra() const
+{
+ // (APR:Eq. 4), (BSG:Eq. 14)
+ const scalar log1 = log((zRefU_ - d_*h_)/(zom_*h_));
+ const scalar log2 = log((zRefH_ - d_*h_)/(zoh_*h_));
+
+ return log1*log2/(sqr(kappa_)*uRef_);
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::grass::rs() const
+{
+ // (BSG:Eq. 13), (APR:Eq. 5; reduced form in p. 4-5)
+ return ri_/(scalar(12)*h_);
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::grass::S
+(
+ const labelList& cells
+) const
+{
+ // Mark fvOption cells within mesh
+ bitSet isZoneCell(mesh().nCells());
+ isZoneCell.set(cells);
+
+ scalar S = 0;
+
+ // Select cells next to specified patches
+ // Sum patch area that is covered by fvOption cells
+ for (const label patchi : patchSet_)
+ {
+ const scalarField& s = mesh().magSf().boundaryField()[patchi];
+
+ const polyPatch& pp = mesh().boundaryMesh()[patchi];
+ const labelList& faceCells = pp.faceCells();
+
+ forAll(faceCells, i)
+ {
+ const bool isCovered = isZoneCell[faceCells[i]];
+
+ if (isCovered)
+ {
+ S += s[i];
+ }
+ }
+ }
+ reduce(S, sumOp<scalar>());
+
+ if (mag(S) < SMALL)
+ {
+ FatalErrorInFunction
+ << "Area coverage of grass cannot be zero. "
+ << "Check whether cellZone has any boundary faces."
+ << exit(FatalError);
+ }
+
+ return S;
+}
+
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::grass::G
+(
+ const labelList& cells
+) const
+{
+ switch (soilHeatFluxMethod_)
+ {
+ case soilHeatFluxType::PROPORTIONAL_TO_SOLAR_RADIATION:
+ {
+ return Csoil_*q();
+ }
+
+ case soilHeatFluxType::BOUNDARY:
+ {
+ // Retrieve heat flux through patches
+ tmp<FieldField<Field, scalar>> tqBf = this->qBf();
+ const FieldField<Field, scalar>& qBf = tqBf.cref();
+
+ // Mark fvOption cells within mesh
+ bitSet isZoneCell(mesh().nCells());
+ isZoneCell.set(cells);
+
+ scalar G = 0;
+
+ for (const label patchi : patchSet_)
+ {
+ const scalarField& s = mesh().magSf().boundaryField()[patchi];
+
+ const scalarField& qfld = qBf[patchi];
+
+ const polyPatch& pp = mesh().boundaryMesh()[patchi];
+ const labelList& faceCells = pp.faceCells();
+
+ forAll(faceCells, i)
+ {
+ const bool isCovered = isZoneCell[faceCells[i]];
+
+ if (isCovered)
+ {
+ G += qfld[i]*s[i];
+ }
+ }
+ }
+ reduce(G, sumOp<scalar>());
+
+ return G/area_;
+ }
+ }
+
+ return -1;
+}
+
+
+Foam::tmp<Foam::FieldField<Foam::Field, Foam::scalar>>
+Foam::evapotranspirationHeatTransferModels::grass::qBf() const
+{
+ const auto& T = mesh().lookupObject<volScalarField>(TName_);
+ const volScalarField::Boundary& Tbf = T.boundaryField();
+
+ auto tq = tmp<FieldField<Field, scalar>>::New(Tbf.size());
+ auto& q = tq.ref();
+
+ forAll(q, patchi)
+ {
+ q.set(patchi, new Field<scalar>(Tbf[patchi].size(), Zero));
+ }
+
+ typedef compressible::turbulenceModel cmpTurbModel;
+
+ if (mesh().foundObject<cmpTurbModel>(cmpTurbModel::propertiesName))
+ {
+ const auto& turb =
+ mesh().lookupObject<cmpTurbModel>(cmpTurbModel::propertiesName);
+
+ // Note: calling he(p,T) instead of he()
+ const volScalarField he(turb.transport().he(turb.transport().p(), T));
+ const volScalarField::Boundary& hebf = he.boundaryField();
+
+ const volScalarField alphaEff(turb.alphaEff());
+ const volScalarField::Boundary& alphaEffbf = alphaEff.boundaryField();
+
+ for (const label patchi : patchSet_)
+ {
+ q[patchi] = alphaEffbf[patchi]*hebf[patchi].snGrad();
+ }
+ }
+ else if (mesh().foundObject<fluidThermo>(fluidThermo::dictName))
+ {
+ const auto& thermo =
+ mesh().lookupObject<fluidThermo>(fluidThermo::dictName);
+
+ // Note: calling he(p,T) instead of he()
+ const volScalarField he(thermo.he(thermo.p(), T));
+ const volScalarField::Boundary& hebf = he.boundaryField();
+
+ const volScalarField& alpha(thermo.alpha());
+ const volScalarField::Boundary& alphabf = alpha.boundaryField();
+
+ for (const label patchi : patchSet_)
+ {
+ q[patchi] = alphabf[patchi]*hebf[patchi].snGrad();
+ }
+ }
+ else
+ {
+ FatalErrorInFunction
+ << "Unable to find a valid thermo model to evaluate q. " << nl
+ << "Database contents are: " << mesh().objectRegistry::sortedToc()
+ << exit(FatalError);
+ }
+
+ // No radiative heat flux contribution is present
+
+ return tq;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::evapotranspirationHeatTransferModels::grass::grass
+(
+ const dictionary& dict,
+ const fvMesh& mesh
+)
+:
+ evapotranspirationHeatTransferModel(dict, mesh),
+ soilHeatFluxMethod_
+ (
+ soilHeatFluxTypeNames.getOrDefault
+ (
+ "soilHeatFluxMethod",
+ dict,
+ soilHeatFluxType::PROPORTIONAL_TO_SOLAR_RADIATION
+ )
+ ),
+ Tptr_(nullptr),
+ RHptr_(nullptr),
+ area_(-1),
+ Csoil_(),
+ rho_(),
+ Cp_(),
+ epsilon_(),
+ h_(),
+ kappa_(),
+ uRef_(),
+ zRefU_(),
+ zRefH_(),
+ zom_(),
+ zoh_(),
+ d_(),
+ ri_(),
+ pAtm_(),
+ Tref_(0),
+ TName_(),
+ patchSet_()
+{
+ Info<< " Activating evapotranspiration heat transfer model: "
+ << typeName << endl;
+
+ read(dict);
+}
+
+
+// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
+
+Foam::tmp<Foam::scalarField>
+Foam::evapotranspirationHeatTransferModels::grass::Q
+(
+ const labelList& cells
+) const
+{
+ if (area_ == -1)
+ {
+ area_ = S(cells);
+ }
+
+ Tref_ = Tptr_->value(mesh().time().timeOutputValue());
+
+ return E(cells)*area_;
+}
+
+
+bool Foam::evapotranspirationHeatTransferModels::grass::read
+(
+ const dictionary& dict
+)
+{
+ if (!evapotranspirationHeatTransferModel::read(dict))
+ {
+ return false;
+ }
+
+ Tptr_.reset(Function1<scalar>::New("Tref", dict, &mesh()));
+ RHptr_.reset(Function1<scalar>::New("RHref", dict, &mesh()));
+
+ area_ = -1;
+
+ const auto range = scalarMinMax::ge(SMALL);
+
+ Csoil_ = dict.getCheckOrDefault<scalar>("Csoil", 0.1, range);
+ rho_ = dict.getCheckOrDefault<scalar>("rho", 1.225, range);
+ Cp_ = dict.getCheckOrDefault<scalar>("Cp", 1013.0, range);
+ epsilon_ = dict.getCheckOrDefault<scalar>("epsilon", 0.622, range);
+ h_ = dict.getCheckOrDefault<scalar>("h", 0.1, range);
+ kappa_ = dict.getCheckOrDefault<scalar>("kappa", 0.41, range);
+ uRef_ = dict.getCheckOrDefault<scalar>("uRef", 2, range);
+ zRefU_ = dict.getCheckOrDefault<scalar>("zRefU", 10, range);
+ zRefH_ = dict.getCheckOrDefault<scalar>("zRefH", 10, range);
+ zom_ = dict.getCheckOrDefault<scalar>("zom", 0.123, range);
+ zoh_ = dict.getCheckOrDefault<scalar>("zoh", 0.0123, range);
+ d_ = dict.getCheckOrDefault<scalar>("d", 2.0/3.0, range);
+ ri_ = dict.getCheckOrDefault<scalar>("ri", 100, range);
+ pAtm_ = dict.getCheckOrDefault<scalar>("pAtm", 101.325, range);
+ TName_ = dict.getOrDefault<word>("T", "T");
+
+ if (soilHeatFluxMethod_ == soilHeatFluxType::BOUNDARY)
+ {
+ patchSet_ =
+ mesh().boundaryMesh().patchSet(dict.get<wordRes>("patches"));
+ }
+
+ return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/grass/grass.H b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/grass/grass.H
new file mode 100644
index 0000000000000000000000000000000000000000..b679f4e3ae2b2f48bfb7c25b6d3f143e6134613b
--- /dev/null
+++ b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/grass/grass.H
@@ -0,0 +1,342 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2022 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::evapotranspirationHeatTransferModels::grass
+
+Description
+ Applies sources on temperature (\c T - incompressible) or energy
+ (\c h/e - compressible) equation to incorporate evapotranspiration
+ heat-transfer effects from the specified grass canopy.
+ The model is based on Pemnan-Monteith Equation model.
+
+ Sources applied to either of the below, if exist:
+ \verbatim
+ T | Temperature [K]
+ e | Internal energy [m^2/s^2]
+ h | Enthalphy [m^2/s^2]
+ \endverbatim
+
+ Required fields:
+ \verbatim
+ T | Temperature [K]
+ e | Internal energy [m^2/s^2]
+ h | Enthalphy [m^2/s^2]
+ LAD | Leaf area density [m^2/m^3]
+ \endverbatim
+
+ References:
+ \verbatim
+ Governing equations (tag:BSG):
+ Brozovsky, J., Simonsen, A., & Gaitani, N. (2021).
+ Validation of a CFD model for the evaluation of urban
+ microclimate at high latitudes: A case study in Trondheim, Norway.
+ Building and Environment, 205, 108175.
+ DOI:10.1016/j.buildenv.2021.108175
+
+ Governing equations (tag:APR):
+ Allen, R. G., Pereira, L. S., Raes, D., & Smith, M. (1998).
+ Crop evapotranspiration-Guidelines for computing crop
+ water requirements-FAO Irrigation and drainage paper 56.
+ Fao, Rome, 300(9), D05109.
+ \endverbatim
+
+Usage
+ Example by using \c constant/fvOptions:
+ \verbatim
+ evapotranspirationHeatTransfer1
+ {
+ // Inherited entries
+ ...
+
+ // Mandatory entries
+ model grass;
+ Tref <Function1<scalar>>;
+ RHref <Function1<scalar>>;
+
+ // Conditional entries
+
+ // when soilHeatFluxMethod == boundary
+ patches <wordRes>;
+
+ // Optional entries
+ soilHeatFluxMethod <word>;
+ Csoil <scalar>;
+ rho <scalar>;
+ Cp <scalar>;
+ epsilon <scalar>;
+ h <scalar>;
+ kappa <scalar>;
+ uRef <scalar>;
+ zRefU <scalar>;
+ zRefH <scalar>;
+ zom <scalar>;
+ zoh <scalar>;
+ d <scalar>;
+ ri <scalar>;
+ pAtm <scalar>;
+ T <word>;
+ }
+ \endverbatim
+
+ where the entries mean:
+ \table
+ Property | Description | Type | Reqd | Deflt
+ model | Model name: grass | word | yes | -
+ Tref | Reference weather station air temperature [Celsius] <!--
+ --> | Function1\<scalar\> | yes | -
+ RHref | Reference weather station relative humidity [%] <!--
+ --> | Function1\<scalar\> | yes | -
+ patches | Names of ground patches | wordRes | yes | -
+ soilHeatFluxMethod | Method to calculate soil heat flux - see below <!--
+ --> | word | no | -
+ Csoil | Proportionality constant of soil heat-flux wrt <!--
+ --> solar radiation | scalar | no | 0.1
+ rho | Mean air density at constant pressure [kg/m^3] <!--
+ --> | scalar | no | 1.225
+ Cp | Specific heat at constant pressure [J/kg/C] <!--
+ --> | scalar | no | 1013.0
+ epsilon | Molecular-weight ratio of water vapour/dry air [-] <!--
+ --> | scalar | no | 0.622
+ h | Height of grass layer [m] | scalar | no | 0.1
+ kappa | Von Karman constant [-] | scalar | no | 0.41
+ uRef | Reference velocity magnitude [m/s] | scalar | no | 2.0
+ zRefU | Height of wind speed measurements [m] | scalar | no | 10.0
+ zRefH | Height of humidity measurements [m] | scalar | no | 10.0
+ zom | Roughness length governing momentum transfer coefficient <!--
+ --> [-] | scalar | no | 0.123
+ zoh | Roughness length governing transfer of heat and vapour <!--
+ --> coefficient [-] | scalar | no | 0.0123
+ d | Zero plane displacement height coefficient [-] <!--
+ --> | scalar | no | 0.666
+ ri | Bulk stomata resistance of leaf [s/m] | scalar | no | 100.0
+ pAtm | Atmospheric pressure [kPa] | scalar | no | 101.325
+ T | Name of temperature field | word | no | T
+ \endtable
+
+ Options for the \c soilHeatFluxMethod entry:
+ \verbatim
+ proportionalToSolarRadiation | Estimate from solar load
+ boundary | Obtain soil heat flux from boundary
+ \endverbatim
+
+ The inherited entries are elaborated in:
+ - \link evapotranspirationHeatTransfer.H \endlink
+ - \link evapotranspirationHeatTransferModel.H \endlink
+
+SourceFiles
+ grass.C
+ grassTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Foam_evapotranspirationHeatTransferModels_grass_H
+#define Foam_evapotranspirationHeatTransferModels_grass_H
+
+#include "evapotranspirationHeatTransferModel.H"
+#include "Function1.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace evapotranspirationHeatTransferModels
+{
+
+/*---------------------------------------------------------------------------*\
+ Class grass Declaration
+\*---------------------------------------------------------------------------*/
+
+class grass
+:
+ public evapotranspirationHeatTransferModel
+{
+ // Private Enumerations
+
+ //- Options for the soil heat flux
+ enum soilHeatFluxType : char
+ {
+ PROPORTIONAL_TO_SOLAR_RADIATION = 0, //!< "Estimate from solar load"
+ BOUNDARY, //!< "Obtain soil heat flux from boundary"
+ };
+
+ //- Names for soilHeatFluxType
+ static const Enum<soilHeatFluxType> soilHeatFluxTypeNames;
+
+
+ // Private Data
+
+ //- Soil heat-flux calculation method
+ enum soilHeatFluxType soilHeatFluxMethod_;
+
+ //- Reference weather station air temperature [Celsius]
+ autoPtr<Function1<scalar>> Tptr_;
+
+ //- Reference weather station relative humidity [%]
+ autoPtr<Function1<scalar>> RHptr_;
+
+ //- Area coverage of grass [m^2]
+ mutable scalar area_;
+
+ //- Proportionality constant of soil heat-flux wrt solar radiation [-]
+ scalar Csoil_;
+
+ //- Mean air density at constant pressure [kg/m^3]
+ scalar rho_;
+
+ //- Specific heat at constant pressure [J/kg/C]
+ scalar Cp_;
+
+ //- Molecular-weight ratio of water vapour/dry air [-]
+ scalar epsilon_;
+
+ //- Height of grass layer [m]
+ scalar h_;
+
+ //- Von Karman constant [-]
+ scalar kappa_;
+
+ //- Reference velocity magnitude [m/s]
+ scalar uRef_;
+
+ //- Height of wind speed measurements [m]
+ scalar zRefU_;
+
+ //- Height of humidity measurements [m]
+ scalar zRefH_;
+
+ //- Roughness length governing momentum transfer coefficient [-]
+ scalar zom_;
+
+ //- Roughness length governing transfer of heat and vapour coeff [-]
+ scalar zoh_;
+
+ //- Zero plane displacement height coefficient [-]
+ scalar d_;
+
+ //- Bulk stomata resistance of leaf [s/m]
+ scalar ri_;
+
+ //- Atmospheric pressure [kPa]
+ scalar pAtm_;
+
+ //- Cached reference weather station air temperature [Celsius]
+ mutable scalar Tref_;
+
+ //- Name of temperature field
+ word TName_;
+
+ //- List of patches to calculate boundary heat flux
+ labelHashSet patchSet_;
+
+
+ // Private Member Functions
+
+ //- Return evapotranspiration heat flux through grass layer [W/m^2]
+ tmp<scalarField> E(const labelList& cells) const;
+
+ //- Return slope of the relation between
+ //- vapour pressure-temperature [Pa/K]
+ scalar Delta() const;
+
+ //- Return atmospheric vapour pressure deficit [Pa]
+ scalar D() const;
+
+ //- Return saturated vapour pressure over water [Pa]
+ scalar pSat() const;
+
+ //- Return psychrometric constant [kPa/K]
+ scalar gamma() const;
+
+ //- Return specific latent heat of vaporisation [MJ/kg]
+ scalar lambda() const;
+
+ //- Return bulk aerodynamic resistance [Pa]
+ scalar ra() const;
+
+ //- Return bulk surface resistance [Pa]
+ scalar rs() const;
+
+ //- Return area coverage of grass [m^2]
+ scalar S(const labelList& cells) const;
+
+ //- Return area-averaged heat flux through ground [W/m^2]
+ scalar G(const labelList& cells) const;
+
+ //- Return heat-flux boundary fields
+ tmp<FieldField<Field, scalar>> qBf() const;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("grass");
+
+
+ // Constructors
+
+ //- Construct from components
+ grass
+ (
+ const dictionary& dict,
+ const fvMesh& mesh
+ );
+
+ //- No copy construct
+ grass(const grass&) = delete;
+
+ //- No copy assignment
+ void operator=(const grass&) = delete;
+
+
+ //- Destructor
+ virtual ~grass() = default;
+
+
+ // Member Functions
+
+ // Evaluation
+
+ //- Return heat-transfer rate for speficied cells [J/s]
+ virtual tmp<scalarField> Q(const labelList& cells) const;
+
+
+ // I-O
+
+ //- Read the dictionary
+ virtual bool read(const dictionary& dict);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace evapotranspirationHeatTransferModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/tree/tree.C b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/tree/tree.C
new file mode 100644
index 0000000000000000000000000000000000000000..9c5cb5daee5683dd5187ef0d34dbcc9e6a0b615e
--- /dev/null
+++ b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/tree/tree.C
@@ -0,0 +1,139 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2022 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "tree.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace evapotranspirationHeatTransferModels
+{
+ defineTypeNameAndDebug(tree, 0);
+ addToRunTimeSelectionTable
+ (
+ evapotranspirationHeatTransferModel,
+ tree,
+ dictionary
+ );
+}
+}
+
+
+// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
+
+Foam::scalar Foam::evapotranspirationHeatTransferModels::tree::Et() const
+{
+ return a_*q() + b_;
+}
+
+
+Foam::tmp<Foam::scalarField>
+Foam::evapotranspirationHeatTransferModels::tree::leafArea
+(
+ const labelList& cells
+) const
+{
+ const volScalarField& LAD = getOrReadField(LADname_);
+
+ const scalarField& V = mesh().V();
+
+ auto tleafArea = tmp<scalarField>::New(cells.size(), Zero);
+ auto& leafArea = tleafArea.ref();
+
+ forAll(cells, i)
+ {
+ const label celli = cells[i];
+
+ leafArea[i] = LAD[celli]*V[celli];
+ }
+
+ return tleafArea;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::evapotranspirationHeatTransferModels::tree::tree
+(
+ const dictionary& dict,
+ const fvMesh& mesh
+)
+:
+ evapotranspirationHeatTransferModel(dict, mesh),
+ a_(),
+ b_(),
+ lambda_(),
+ LADname_()
+{
+ Info<< " Activating evapotranspiration heat transfer model: "
+ << typeName << endl;
+
+ read(dict);
+}
+
+
+// * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * * //
+
+Foam::tmp<Foam::scalarField>
+Foam::evapotranspirationHeatTransferModels::tree::Q
+(
+ const labelList& cells
+) const
+{
+ // Convert units from [MJ g/hr] to [J kg/s]
+ static const scalar unitConverter = scalar(1000)/scalar(3600);
+
+ return unitConverter*lambda_*Et()*leafArea(cells);
+}
+
+
+bool Foam::evapotranspirationHeatTransferModels::tree::read
+(
+ const dictionary& dict
+)
+{
+ if (!evapotranspirationHeatTransferModel::read(dict))
+ {
+ return false;
+ }
+
+ const auto range = scalarMinMax::ge(SMALL);
+
+ a_ = dict.getOrDefault<scalar>("a", 0.3622);
+ b_ = dict.getOrDefault<scalar>("b", 60.758);
+ lambda_ = dict.getCheckOrDefault<scalar>("lambda", 2.44, range);
+ LADname_ = dict.getOrDefault<word>("LAD", "LAD");
+
+ (void) getOrReadField(LADname_);
+
+ return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/tree/tree.H b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/tree/tree.H
new file mode 100644
index 0000000000000000000000000000000000000000..c8704a2cb300d8c4d31c8727b6595dce12a8f807
--- /dev/null
+++ b/src/atmosphericModels/fvOptions/evapotranspirationHeatTransfer/evapotranspirationHeatTransferModels/tree/tree.H
@@ -0,0 +1,219 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | www.openfoam.com
+ \\/ M anipulation |
+-------------------------------------------------------------------------------
+ Copyright (C) 2022 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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 <http://www.gnu.org/licenses/>.
+
+Class
+ Foam::evapotranspirationHeatTransferModels::tree
+
+Description
+ Applies sources on temperature (\c T - incompressible) or energy
+ (\c h/e - compressible) equation to incorporate evapotranspiration
+ heat-transfer effects from the specified tree canopy.
+
+ The heat-transfer is calculated based on the following relations.
+ The heat transfer is given by the first equation and the
+ evapotranspiration is calculated linearly proportional to absorbed
+ solar radiation (direct and diffusive) with the empirical relation
+ shown in the second equation below.
+
+ \f[
+ Q = \lambda E_t \frac{1000}{3600} \int_{cellZone} LAD \, dV
+ \f]
+
+ where
+ \vartable
+ Q | Heat transfer from tree crown [W]
+ \lambda | Specific latent heat of vaporisation [MJ/kg]
+ E_t | Area-averaged evapotranspiration [g/m^2/hr]
+ LAD | Leaf area density [m^2/m^3]
+ V | Volume [m^3]
+ \endvartable
+
+ with
+
+ \f[
+ E_t = a R_n + b
+ \f]
+
+ where
+ \vartable
+ R_n | Heat flux through tree crown [W/m^2]
+ a | Linear-regression coefficient - slope parameter [g/W/hr]
+ b | Linear-regression coefficient - intercept parameter [g/m^2/hr]
+ \endvartable
+
+ Sources applied to either of the below, if exist:
+ \verbatim
+ T | Temperature [K]
+ e | Internal energy [m^2/s^2]
+ h | Enthalphy [m^2/s^2]
+ \endverbatim
+
+ Required fields:
+ \verbatim
+ T | Temperature [K]
+ e | Internal energy [m^2/s^2]
+ h | Enthalphy [m^2/s^2]
+ LAD | Leaf area density [m^2/m^3]
+ \endverbatim
+
+Usage
+ Example by using \c constant/fvOptions:
+ \verbatim
+ evapotranspirationHeatTransfer1
+ {
+ // Inherited entries
+ ...
+
+ // Mandatory entries
+ model tree;
+
+ // Optional entries
+ a <scalar>;
+ b <scalar>;
+ lambda <scalar>;
+ LAD <word>;
+ }
+ \endverbatim
+
+ where the entries mean:
+ \table
+ Property | Description | Type | Reqd | Deflt
+ model | Model name: tree | word | yes | -
+ a | Linear-regression coefficient - slope parameter [g/W/hr] | scalar <!--
+ --> | no | 0.3622
+ b | Linear-regression coefficient - intercept parameter [g/m^2/hr] <!--
+ --> | scalar | no | 60.758
+ lambda | Specific latent heat of vaporisation [MJ/kg] | scalar <!--
+ --> | no | 2.44
+ LAD | Name of leaf-area density field | word | no | LAD
+ \endtable
+
+ The inherited entries are elaborated in:
+ - \link evapotranspirationHeatTransfer.H \endlink
+ - \link evapotranspirationHeatTransferModel.H \endlink
+
+SourceFiles
+ tree.C
+ treeTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Foam_evapotranspirationHeatTransferModels_tree_H
+#define Foam_evapotranspirationHeatTransferModels_tree_H
+
+#include "evapotranspirationHeatTransferModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace evapotranspirationHeatTransferModels
+{
+
+/*---------------------------------------------------------------------------*\
+ Class tree Declaration
+\*---------------------------------------------------------------------------*/
+
+class tree
+:
+ public evapotranspirationHeatTransferModel
+{
+ // Private Data
+
+ //- Linear-regression coefficient - slope parameter [g/W/hr]
+ scalar a_;
+
+ //- Linear-regression coefficient - intercept parameter [g/m^2/hr]
+ scalar b_;
+
+ //- Specific latent heat of vaporisation [MJ/kg]
+ scalar lambda_;
+
+ //- Name of leaf-area density field
+ word LADname_;
+
+
+ // Private Member Functions
+
+ //- Return area-averaged transpiration [g/m^2/hr]
+ // Calculated from empirical regressions with heat flux
+ scalar Et() const;
+
+ //- Return volume weighted leaf area density (LAD) [m^2]
+ // LAD [m^2/m^3]: one-sided area of leaves per unit volume
+ tmp<scalarField> leafArea(const labelList& cells) const;
+
+
+public:
+
+ //- Runtime type information
+ TypeName("tree");
+
+
+ // Constructors
+
+ //- Construct from components
+ tree
+ (
+ const dictionary& dict,
+ const fvMesh& mesh
+ );
+
+ //- No copy construct
+ tree(const tree&) = delete;
+
+ //- No copy assignment
+ void operator=(const tree&) = delete;
+
+
+ //- Destructor
+ virtual ~tree() = default;
+
+
+ // Member Functions
+
+ // Evaluation
+
+ //- Return heat-transfer rate for speficied cells [J/s]
+ virtual tmp<scalarField> Q(const labelList& cells) const;
+
+
+ // I-O
+
+ //- Read the dictionary
+ virtual bool read(const dictionary& dict);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace evapotranspirationHeatTransferModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //