Welcome to Greenhouses’s documentation!¶

Organization:	University of Liège
Date:	Aug 02, 2019

Greenhouses is an open-source library for dynamic modelling of greenhouse climate developed in the Modelica language. The library aims at providing a modeling framework capable of simulating not only the energy flows of a greenhouse climate but also the complex interactions and energy flows relative to systems coupling the greenhouse e.g. generation and storage units. A number of platforms are currently available for greenhouse climate simulation, but none of them is open-source. Moreover, no platform is currently available for simulating the energy integration of greenhouses with other systems. The goal of Greenhouses is to fill this gap.

For greenhouse-scale simulation, the proposed library includes the modeling of greenhouse construction elements, indoor climate, heating systems, ventilation systems and crop yield. To the end of performing system-scale simulations, the library includes robust performance-based models of several generation units (e.g. CHP and heat pumps) as well as a model of thermal energy storage. The library also includes several pre-programmed control systems for climate control (heating circuit, window’s aperture, supplementary lighting, etc.) and operation control of the HVAC units.

Explore the Examples package to get a feeling of the modeling possibility that the Greenhouses library offers!

Downloading Greenhouses¶

Greenhouses can be downloaded from its github repository (using the Clone or Download button on the right side of the screen): https://github.com/queraltab/Greenhouses-Library

License¶

The Greenhouses package is free software licensed under the Modelica License 2. It can be redistributed and/or modified under the terms of this license.

Main developers¶

Queralt Altes-Buch (University of Liège, Belgium)

Contents¶

Overview¶

A greenhouse climate model is a model that describes the indoor climate of a greenhouse resulting from the greenhouse design, the outdoor climate and a specific control. In greenhouses, the indoor climate is characterized by the temperature, the vapor pressure of water (i.e. the relative humidity) and the CO₂ concentration of the air. Together with the temperature of the heating pipes, the indoor climate constitutes the climate controller feed- back quantities. However, in order to attain the desired climate, the variables with an indirect influence on the climate also need to be modeled. These are mainly the characteristics relative to the canopy and the envelope (i.e. the cover, the floor and the thermal screen). The canopy temperature has an impact on its photosynthesis and transpiration, which decrease the CO₂ concentration and increase the vapor content of the air, respectively. Evaporation or condensation at surfaces may occur depending on the water vapor pressure difference with respect to the air. The temperature of the envelope influences the vapor pressure of water of the air, which is decreased by condensation at the cover and at the thermal screen. The thermal screen is a membrane used to reduce the energy requirement to heat the greenhouse. Given the porous nature of the screen, air and moisture is exchanged through its fabric. Air exchange with the outside decreases the vapor pressure of water and the CO₂ concentration of the air, which can be increased by supplementary CO₂ supply by an external source.

Greenhouse climate models have been the object of a substantial literature. While many models have been developed ([6][10][24][32][19]), most of them can only be used for a single location and for a specific greenhouse structure and climate. Recently, a more generic greenhouse climate model [34] combining the work of [6] and [10] was developed. The model was validated for a range of climates and greenhouse designs. For the purpose of this work, the model [34] has been implemented in the Modelica language. It should be noted that, when the screen is drawn, the air of the greenhouse is divided in two zones, i.e. below and above the screen. These zones are modeled separately and their climate is assumed to be homogeneous. In all the models of the library, the air zones below and above the thermal screen are going to be referred as main and top air zones, respectively.

Getting started¶

In this tutorial, the necessary steps to lead the correct packages to use Greenhouses with Dymola are reported.

Download the Greenhouses library from the Greenhouses github repository.
Launch Dymola on your computer.
From the Dymola window, use the Open icon to browse to the Greenhouses folder and open the Greenhouses package.mo file.

Dymola loads the Greenhouses library and the following will be displayed:

Congratulations! You have set the environment needed to use the Greenhouses library in Dymola. You are now ready to simulate a wide range of greenhouse climates and thermal systems. Enjoy it!

Library structure¶

The Greenhouses library is hierarchically structured into different packages, listed below:

Components is the central part of the library. It is organized in three sub-packages: Greenhouse, HVAC and CropYield. It contains all the models available in the library from the simple greenhouse components (e.g. cover) to already-build greenhouse models ready to use; generation units models and a yield model for tomato crop.
Flows contains models of the flows that are encountered in a greenhouse system. It is organized in eight sub-packages that model the heat, vapour mass and CO₂ mass transfer, as well as fluid flow. It also contains a sub-package of interfaces, which defined the type of connectors used in the library.
ControlSystems contains the control units. It is organized in two sub-packages: climate control (i.e. the control of the thermal screen, artificial lighting and window’s aperture) and HVAC control (i.e. the control of the operation of the generation units).
Examples contains models where the components of the library can be tested. It includes the simulation of a greenhouse and two system-scale models that simulate the greenhouse connected to thermal energy storage, a CHP and a heat pump.
Interfaces contains all the type of connectors used in the library.
Functions are the empirical correlations used to characterize some of the models presents in the library.
Icons defines the graphical interface for some of the models in the library.

The main packages are further divided into sub-packages and models. In Greenhouses, each model has a relative documentation explaining the main features. To access the documentation layer of each model you can click on the INFO icon on the Dymola menu.

Graphical user interface¶

The developed modeling framework, being object-oriented, is made of independent sub-models for each greenhouse component and exchanged flow, that are interconnected to build a greenhouse model. The sub-models interact together through standard interfaces called ports. For modeling heat transfer and fluid flow, the heat transfer and the thermo-fluid flow connectors from the Modelica Standard Library are used. For moisture and CO₂ mass transfer, two connectors have been developed. An extra connector is developed for the short-wave radiation heat inputs. In total, five connectors are distinguished:

Heat port: from the Modelica Standard Library. Graphically represented by a red square, it is a thermal port for 1-dim heat transfer. Temperature (\(T\) [K]) and heat flow rate (\(\dot{Q}\) [W]) are the potential and flow variables, respectively.
Water mass port: graphically represented by a blue square, it is a mass port for 1-dim moisture transfer. The vapor pressure of water (\(P_v\) [Pa]) and the vapor mass flow rate (\(\dot{M}_v\) [kg s⁻¹]) are the potential and flow variables, respectively.
CO₂ mass port: graphically represented by a grey square, it is a mass port for 1-dim CO₂ mass transfer. The CO₂ concentration (\(CO2\) [mg m⁻³]) and the CO₂ mass flow rate (\(\dot{m}_c\) [mg m⁻² s⁻¹]) are the potential and flow variables, respectively.
Short-wave radiation connector: graphically represented by a yellow triangle (single input) or circle (vector input), it is an input/output connector for forced radiation flows from the sun and supplementary lighting. The radiation flow (W m⁻²) is the potential variable.
Thermo-fluid port: from the Modelica Standard Library. Graphically represented by a blue circle, it is an interface for quasi one-dimensional fluid flow in a piping network (incompressible or compressible, one or more phases, one or more substances). The connector is defined by the pressure and the mass flow rate as the potential and flow variables. Specific enthalpy and mass fractions are stream variables.

For more information on the definition of connectors, check the Modelica users guide.

An example of greenhouse model is shown in the figure below. As it can be distinguished, the greenhouse modeled in this example consists of two levels of heating circuits, roof windows (but not side vents), natural ventilation (no forced ventilation) and a movable thermal screen. The majority of the flows distinguished in a greenhouse originate from convection at surfaces, ventilation processes, conduction at the soil and long-wave infrared radiation. Forced flows such as the short-wave radiation from the sun, latent heat flows or the sensible heat from supplementary lighting are also considered.

Components¶

The Components package is the central part of the library. It is organized in three sub-packages:

Greenhouse¶

The Greenhouse sub-package include models of the different greenhouse elements at which the energy and mass balances are performed to compute all the state variables (temperature, vapor pressure of water and CO₂ concentration) of these elements. These are the main and top air zones, the canopy, the cover, the floor and the thermal screen. For a given component, the energy balance is done by taking into account all the heat flows that are connected to the heat port (i.e. conduction, convection and/or long-wave radiation), the latent heat of the vapor flows connected to the vapor mass port and the forced short-wave radiation input connected to the short-wave radiation flow connector. The computation of the heat and vapor flows is done in individual models, described in the Flows section. Since no spatial differences in temperature, vapor pressure of water and CO₂ concentration are considered, all the model flows are described per square meter of greenhouse floor.

Except for the cover model, since the short-wave radiation can origin from two sources (the sun and supplementary illumination), the short-wave radiation input has the form of a vector. The parameter N_rad defines the dimension of the vector and depending on whether or not the greenhouse has supplementary lighting, the user can set it to two or one, respectively. The computation of the absorbed short-wave radiation flows is done in the Solar_model and the Illumination model, according to the provenance of the flow.

The Greenhouse sub-package also includes an already-build greenhouse model ready to use.

General nomenclature¶

Subscripts¶

Name	Description
Air	Greenhouse main air zone
Can	Canopy
cnd	Conduction
cnv	Convection
Cov	Cover
Ext	External source of CO₂
Flr	Floor
Glob	Global radiation
Ilu	Supplementary lighting
lat	Latent heat
Out	Outside air
Pipe	Pipe heating system
rad	Long-wave infrared radiation
Scr	Thermal screen
So(j)	The j^th soil layer
Top	Greenhouse top air zone
v	Vapor

Variables¶

Name	Units	Description
\(\alpha\)	–	Absorption coefficient
\(\rho\)	kg m⁻³, –	Density or reflection coefficient
\(\tau\)	–	Transmission coefficient
\(\phi\)	º	Inclination of the greenhouse cover
c_p	J kg⁻¹ K⁻¹	Specific heat capacity
LAI	m²{leaf} m⁻²{floor}	Ratio of leaf area over greenhouse floor area
h	m	Thickness or vertical dimension
\(\dot{m}\)	kg m⁻² s⁻¹	Mass flow rate averaged per square meter of greenhouse floor
\(\dot{q}\)	W m⁻²	Heat flow averaged per square meter of greenhouse floor
T	K	Temperature

Solar model¶

The solar radiation incident in a greenhouse can be split in three spectral parts: ultra violet (UV, from 0.3 to 0.4 μm), visible light (from 0.4 to 0.7 μm) and near infrared light (NIR, from 0.7 to 3 μm). The visible light has an interest for biological growth and is referred as photosynthetically active radiation (PAR) in greenhouse modeling. The fraction of UV is 6- 10% and of PAR is 45-60% of the global radiation [8]. However, for plant growth it is common to assign 50% to PAR, neglect the UV and assign the other 50% to NIR [10]. Besides the spectral division, the solar radiation can be divided in direct and diffuse radiation. The solar model of this work is simplified by making no distinction between diffuse and direct solar radiation and by assuming that the transmission coefficient of the greenhouse cover does not depend on the solar angle.

This model computes the short-wave radiation originated from the sun absorbed by the different components of a greenhouse. To that end, the optical properties of the greenhouse components must be known. Therefore, the main parameters of the model are:

The floor reflection coefficients for PAR and NIR (\(\rho_{Flr,PAR}\), \(\rho_{Flr,NIR}\))
The cover transmission and reflection coefficients for PAR and NIR (\(\rho_{Cov,PAR}\), \(\rho_{Cov,NIR}\), \(\tau_{Cov,PAR}\), \(\tau_{Cov,NIR}\))
The screen transmission and reflection coefficients for PAR and NIR (\(\rho_{Scr,PAR}\), \(\rho_{Scr,NIR}\), \(\tau_{Scr,PAR}\), \(\tau_{Scr,NIR}\))
The canopy reflection coefficients for PAR and NIR (\(\rho_{Can,PAR}\), \(\rho_{Can,NIR}\))
The canopy extinction coefficient for PAR and NIR (\(K_{Can,PAR}\), \(K_{Can,NIR}\))

The user can adapt the parameters to the characteristics of the modeled greenhouse. Default values of the parameters are given for tomato crop, single-glass cover, concrete floor and aluminised screen.

The radiation from the sun is partially absorbed by the cover and partially transmitted inside the greenhouse. The global radiation partially absorbed by the cover is described by:

\[\dot{q}_{SunCov} = (\alpha_{Cov,PAR} \cdot \eta_{Glob,PAR} + \alpha_{Cov,NIR} \cdot \eta_{Glob,NIR})\cdot I_{Glob}\]

where \(\alpha\)_Cov,PAR and \(\alpha\)_Cov,NIR are the absorption coefficients for photosynthetically active radiation (PAR) and near-infrarred radiation (NIR) of the cover, and \(\eta\)_Glob,PAR and \(\eta\)_Glob,NIR are the ratio of PAR and NIR of the global radiation.

The radiation that is not reflected or absorbed by the cover, is transmitted inside the greenhouse. The transmitted PAR, i.e. the PAR above the canopy, can be defined by:

\[\dot{q}_{PAR,\tau} = (1-\eta_{Glob,Air})\cdot \tau_{Cov,PAR} \cdot \eta_{Glob,PAR} \cdot I_{Glob}\]

where \(\eta\)_Glob,Air is the ratio of the radiation that is absorbed by the greenhouse elements and is later released to the air as long-wave radiation. When the thermal screen is drawn, it influences the transmission, reflection and absorption coefficients of the greenhouse. Thus, \(\tau\)_Cov,PAR is the coefficient of a lumped transmission of the greenhouse cover and the movable thermal screen.

The PAR absorbed by the canopy is the sum of the PAR transmitted by the cover and directly absorbed by the canopy (\(\dot{q}_{PAR,SunCan↓}\)) and the PAR reflected by the greenhouse floor and then absorbed by the canopy (\(\dot{q}_{PAR,FlrCan↑}\)):

\[\dot{q}_{PAR,SunCan} = \dot{q}_{PAR,SunCan↓} + \dot{q}_{PAR,FlrCan↑}\]

According to [28], \(\dot{q}_{PAR,SunCan}\) in a homogenous crop is described by an exponential decomposition of light with the LAI:

\[\dot{q}_{PAR,SunCan↓} = \dot{q}_{PAR,\tau} \cdot (1-\rho_{Can,PAR}) \cdot \left( 1-e^{-K_{PAR} \cdot LAI}\right)\]

The LAI, computed in the crop yield model, is an input of the solar model. In a similar way:

\[\dot{q}_{PAR,FlrCan↑} = \dot{q}_{PAR,\tau} \cdot e^{-K_{PAR} \cdot LAI} \cdot \rho_{Flr,PAR} \cdot (1-\rho_{Can,PAR}) \cdot \left( 1-e^{-K_{PAR} \cdot LAI}\right)\]

The PAR absorbed by the greenhouse floor is described by:

\[\dot{q}_{PAR,SunFlr} = \dot{q}_{PAR,\tau} \cdot (1-\rho_{Flr,PAR}) \cdot \left( 1-e^{-K_{PAR} \cdot LAI}\right)\]

The reflections coefficients of the canopy and floor are higher for NIR than for PAR. Thus, a larger amount of NIR is reflected by the canopy and the floor back into the greenhouse, which can be reflected again by the greenhouse cover. This leads to a considerable scattering of NIR in the greenhouse. The NIR absorbed by the canopy and the floor is thus determined by considering the lumped cover, the canopy and the floor as a multiple layer model [34]. To that end, virtual transmission coefficients of the lumped cover and floor are applied. The tranmission and reflection coefficients of the canopy, which depend on the LAI, are also determined using the following relationship:

\[\widehat{\tau}_{Can,NIR} = e^{-K_{NIR} LAI}\]

Afterwards, the NIR transmission, reflection and absorption coefficients of the multiple-layers are determined by using the NIR transmission and reflection coefficients of the individual layers in the multiple-layer model. Finally, the calculated absorption coefficient of the multiple-layers equals the overall NIR absorption coefficient of the canopy (\(\alpha_{Can,NIR}\)), and the calculated transmission coefficient of the multiple-layers equals the overall NIR absorption coefficient of the floor (\(\alpha_{Flr,NIR}\)). The NIR absorbed by the canopy and the floor is described by:

\[\dot{q}_{NIR,SunCan} = \alpha_{Can,NIR} \cdot (1-\eta_{Glob,Air}) \cdot \eta_{Glob,NIR} \cdot I_{Glob}\]

\[\dot{q}_{NIR,SunFlr} = \alpha_{Flr,NIR} \cdot (1-\eta_{Glob,Air}) \cdot \eta_{Glob,NIR} \cdot I_{Glob}\]

Part of the radiation transmitted in the greenhouse (that is not absorbed by the canopy or the floor) is absorbed by the greenhouse construction elements and later released as long-wave radiation to the greenhouse air:

\[\dot{q}_{SunAir} = \eta_{Glob,Air} \cdot I_{Glob} \cdot (\tau_{Cov,PAR} \cdot \eta_{Glob,PAR} +(\alpha_{Can,NIR}+\alpha_{Flr,NIR}) \cdot \eta_{Glob,NIR})\]

As shown in the figure below, the solar model has one input connector, for the global irradiation, and four output connectors, whose output values are equal to the absorbed global radiation by (from top to bottom):

The cover (\(\dot{q}_{SunCov}\))
The air (\(\dot{q}_{SunAir}\))
The canopy (\(\dot{q}_{SunCan}=\dot{q}_{PAR,SunCan}+\dot{q}_{NIR,SunCan}\))
The floor (\(\dot{q}_{SunFlr}=\dot{q}_{PAR,SunFlr}+\dot{q}_{NIR,SunFlr}\))

Each output connector must be connected to the short-wave input connector of its corresponding component. The model outputs the radiation values in Wm⁻². However, the PAR absorbed values are also available in the model in μmol{photons}m⁻²s⁻¹, for which a conversion factor from global radiation to PAR equal to 2.3 μmol{photons} J⁻¹ is used.

Cover¶

This model computes the energy balance on the cover. Because the properties of the cover are parameters of the model, the user has the possibility to adapt the model for any type of cover (single-glass, double-glass, plastic, etc.). The main parameters of the model are:

\(\rho\)_Cov: the cover density
c_p,Cov : the cover specific heat capacity
h_Cov: the cover thickness
\(\phi\): the cover inclination (angle with the horizontal)

At the cover, many exhanges occur. Sensible heat caused by convection is mainly exhanged with the air (top air zone when the screen is drawn, main air zone in the absence of a thermal screen) and the outside air. Latent flows caused by condensation may appear in the inner side of the cover. Long-wave radiation is exchanged with the heating pipes, the canopy, the floor, the thermal screen and the sky. Moreover, the cover absorbs part of the incident solar irradiation. The energy balance on the cover is therefore described by:

\[\rho_{Cov} c_{p,Cov} \dfrac{h_{Cov}}{cos \phi} \dfrac{d{T}_{Cov}}{dt} = \dot{q}_{SunCov} + \dot{q}_{cnv,TopCov} + \dot{q}_{lat,TopCov} + \dot{q}_{cnv,AirCov} + \dot{q}_{lat,AirCov} + \dot{q}_{rad,PipeCov} + \dot{q}_{rad,CanCov} + \dot{q}_{rad,FlrCov} + \dot{q}_{rad,ScrCov} - \dot{q}_{cnv,CovOut} - \dot{q}_{rad,CovSky}\]

The vapor pressure of water at the cover is defined by the saturated vapor pressure for its temperature.

Canopy¶

This model computes the energy balance on the canopy. The main parameter of the model is the heat capacity of a square meter of canopy leaves (c_Leaf), which was estimated at 1200 J K⁻¹ m⁻² [30]. The leaf area index (LAI), computed in the crop yield model, is an input of the model.

The canopy exchanges long-wave radiation with the heating pipes, the cover, the floor and the thermal screen. Also, it exchanges sensible and latent heat with the main air zone. As computed in the solar model, part of the transmitted radiation is absorbed by the canopy. The energy balance is described by:

\[LAI c_{Leaf} \dfrac{d{T}_{Can}}{dt} = \dot{q}_{SunCan} + \dot{q}_{IluCan} + \dot{q}_{rad,PipeCan} - \dot{q}_{cnv,CanAir} - \dot{q}_{lat,CanAir} - \dot{q}_{rad,CanCov} - \dot{q}_{rad,CanFlr} + \dot{q}_{rad,CanScr}\]

Air¶

This model computes the energy and vapor mass balance on the air of the main zone of the greenhouse. The main parameters of the model are:

\(\rho\)_Air: the air density
c_p,Air : the air specific heat capacity
h_Air: the height of the main air zone (i.e. height of the thermal screen when this is drawn, mean height of the greenhouse when the screen is open)

Sensible heat is exchanged between the air and the heating pipes, the floor, the canopy, the cover, the thermal screen, the top air zone and the outdoor air. The exchange between the two air zones through the thermal screen occurs because of the porosity material, nature of the latter. The exchange with the outside air accounts for infiltration/exfiltration and natural ventilation through the greenhouse windows. The energy balance on the air of the main zone is dscribed by:

\[\rho_{Air} c_{p,Air} h_{Air} \dfrac{d{T}_{Air}}{dt} = \dot{q}_{SunAir} + \dot{q}_{IluAir} + \dot{q}_{cnv,PipeAir} + \dot{q}_{cnv,CanAir} - \dot{q}_{cnv,AirFlr} - \dot{q}_{cnv,AirCov} - \dot{q}_{cnv,AirScr} - \dot{q}_{cnv,AirTop} - \dot{q}_{cnv,AirOut}\]

The main air zone exchanges air through ventilation processes with the top air compartment and the outside air. The transpiration process of the canopy increases the vapor content of the air. Given the conditions, condensation may occur at the cover and the thermal screen. The vapor pressure of water in the air, defined by the mass balance on the air, is therefore described by:

\[M_H \dfrac{h_{Air}}{RT} \dfrac{d{P}_{v,Air}}{dt} = \dot{m}_{v,CanAir} - \dot{m}_{v,AirCov} - \dot{m}_{v,AirScr} - \dot{m}_{v,AirTop} - \dot{m}_{v,AirOut}\]

where M_H is the molar mass of vapor. Although the vapor capacity is a function of temperature (T), this model applies a constant capacity, holding for a temperature of 291 K.

Air Top¶

This model is a simplified version of the air model designed for the top air zone in which the heat inputs from short-wave radiation are neglected. This is because the top zone is modeled when the screen is drawn (i.e. mostly at night) and because the supplementary lighting lamps are installed below the thermal screen (i.e. below the zone). The model has the same parameters than the Air model, considering that the height of the top air zone is equal to the mean greenhouse height minus the height of the thermal screen.

Floor¶

This model computes the energy balance on the floor. The main parameters of the model are:

\(\rho\)_Flr: the floor density
c_p,Flr : the floor specific heat capacity
h_Flr: the thickness of the floor

The model can thus be applied to a wide range of floor materials (e.g. soil, concrete). Long-wave radiation is exchanged between the floor and the heating pipes, the canopy, the thermal screen and the cover. Sensible heat is exchanged with the air by convection and with the first soil layer by conduction. Therefore, the energy balance on the floor is described by:

\[\rho_{Flr} c_{p,Flr} h_{Flr} \dfrac{d{T}_{Flr}}{dt} = \dot{q}_{SunFlr} + \dot{q}_{IluFlr} + \dot{q}_{cnv,AirFlr} + \dot{q}_{rad,PipeFlr} + \dot{q}_{rad,CanFlr} - \dot{q}_{cnd,FlrSo1} - \dot{q}_{rad,FlrScr} - \dot{q}_{rad,FlrCov}\]

Thermal screen¶

This model computes the energy balance on the thermal screen. The main parameters of the model are:

\(\rho\)_Scr: the screen density
c_p,Scr : the screen specific heat capacity
h_Scr: the thickness of the screen
\(\tau\)_FIR: long-wave radiation transmission coefficient of the screen

The main input of the model is the screen closure (u_SC)

Long-wave radiation is exchanged with the heating pipes, the canopy, the floor and the cover. Sensible heat is exchanged with the air. Latent heat flows may be present if condensation beneath the screen or evaporation above the screen occurs. The energy balance on the screen is therefore described by:

\[\rho_{Scr} c_{p,Scr} h_{Scr} \dfrac{d{T}_{Scr}}{dt} = \dot{q}_{rad,PipeScr} + \dot{q}_{rad,CanScr} + \dot{q}_{rad,FlrScr} + \dot{q}_{cnv,AirScr} + \dot{q}_{lat,AirScr} - \dot{q}_{cnv,ScrTop} - \dot{q}_{lat,ScrTop} - \dot{q}_{rad,ScrCov}\]

The present model assumes that the thermal screen is capable of transporting water from the lower side to the upper side through the fabric. However, the storage of moisture in the screen is neglected. This implies that the vapor that condensates at the screen is either evaporated at the upper side of drips from the screen. Therefore, the rate of evaporation is lower or equal to the rate of condensation.

Illumination¶

Although the contrbution of supplementary lighting is very small during summer, in winter it can double the sun input through a day and thus, have an important impact on crop growth. This model computes the short-wave radiation originated from supplementary lighting absorbed by the different components of a greenhouse (i.e. canopy, floor and air). The main parameters of the model are:

The installed electrical capacity per square meter of greenhouse floor
The floor reflection coefficients for PAR and NIR (\(\rho_{Flr,PAR}\), \(\rho_{Flr,NIR}\))
The canopy reflection coefficients for PAR and NIR (\(\rho_{Can,PAR}\), \(\rho_{Can,NIR}\))
The canopy extinction coefficient for PAR and NIR (\(K_{Can,PAR}\), \(K_{Can,NIR}\))

Only part of the electric consumption of the supplementary lighting is converted to short-wave radiation. This model is based on high pressure sodium (HPS) lamps, which are a type of high intensity discharge lamps. For HPS, it is common that 17% of the electrical power is converted to NIR and 25% to visible light [31]. Thus, the fraction of radiation absorbed by the greenhouse elements and later released to the air in the form of long-wave radiation is assumed to be 58%.

The PAR and NIR absorbed by the canopy and the floor are computed similarly than in the Solar model. As shown in the figure below, the illumination model has one input connector, for the ON-OFF control of the lamps, and three output connectors, whose output values are equal to the absorbed global radiation by (from left to right):

The floor (\(\dot{q}_{IluFlr}=\dot{q}_{PAR,IluFlr}+\dot{q}_{NIR,IluFlr}\))
The canopy (\(\dot{q}_{IluCan}=\dot{q}_{PAR,IluCan}+\dot{q}_{NIR,IluCan}\))
The air (\(\dot{q}_{IluAir}\))

Each output connector must be connected to the short-wave input connector of its corresponding component. The model outputs the radiation values in Wm⁻². However, the PAR absorbed values are also available in the model in μmol{photons}m⁻²s⁻¹, for which a conversion factor from global radiation to PAR equal to 1.8 μmol{photons} J⁻¹, characteristic of HPS lamps, is used.

Heating pipes¶

The fluid in the heating pipes from the greenhouse heating ciruit is modeled by means of the discretized model for incompressible flow explained in the Flows section, in which a dynamic energy balance and static mass and momentum balances are applied on the fluid cells. Heat is transferred by long-wave radiation to the canopy, floor and cover, and by convection to the air. Since the resistance to heat transport to the air from the outer pipe surface is about 100 times greater than the resistance from the inner surface to the outer one [10], it can be assumed that the temperature of the pipe surface is equal to the water temperature. The heat flow from the fluid to the air is thus computed by an ideal heat transfer model with constant heat transfer coefficient.

Greenhouse heating circuits are commonly made of several parallel heating loops. The main parameters of the model are:

The pipe’s diameter (\(d\))
The installed length per loop (\(l\))
The ground floor area (\(A\))
The number of parallel loops (\(N_p\))
The number of nodes in which a loop is discretized (\(N\))
The nominal mass flow rate (\(\dot{M}_{nom}\))

HVAC¶

CHP¶

The CHP model is a performance-based model that assumes constant natural gas consumption and constant total efficiency. The electrical efficiency is described by the second-law efficiency and the efficiency of Carnot:

\[\eta_{el} = \eta_{II} \eta_{Carnot}\]

The second law efficiency is assumed constant and is defined by the previous equation at the nominal conditions. The electrical and thermal power are computed by:

\[\dot{W}_{CHP} = \eta_{el} \dot{Q}_{nom,gas} U_{ONOFF}\]

\[\dot{Q}_{CHP} = \eta_{th} \dot{Q}_{nom,gas} U_{ONOFF}\]

where the thermal efficiency is described by:

\[\eta_{th} = \eta_{tot} - \eta_{el}\]

The heat source is assumed to be at a constant temperature. The properties of the primary fluid are computed using the incompressible Cell1Dim model. The model includes a boolean input connector, which defines the operational state of the CHP, and a real output connector for the electrical power. In the equations, the boolean input is translated to the variable U_ONOFF.

Heat Pump¶

This heat pump model does not consider part-load operation (ON/OFF regulation is assumed). The nominal heat flow is an input of the model, characteristic of the size of the heat pump. The COP is described by:

\[COP = \dfrac{\eta_{II}}{\eta_{Carnot}}\]

The second-law efficiency is assumed to remain unchanged in part-load operation. The electrical power consumed by the heat pump is described by:

\[\dot{W} = \dfrac{\dot{Q}}{COP}\]

In order to account for the lower heat capacity of the heat pump at lower evaporating temperature, the heating power and the heat soure temperature are computed by assuming a linear correlation between their actual and nominal values:

\[\dfrac{\dot{Q}}{\dot{Q}_{nom}} = \dfrac{T_c}{T_{c,nom}}\]

The primary fluid is modeled by means of 1-D incompressible fluid flow model (Cell1DimInc), in which a dynamic energy balance and static mass and momentum balances are applied on the fluid. The heat transfer in the primary fluid is modeled with a constant heat transfer coefficient. However, it can be changed to other heat transfer models through the HeatTransfer parameter in the fluid model.

Furthermore, the model includes a Boolean input connector on_off, which defines the operational state of the heat pump. In the equations, the Boolean input is translated to a variable value together with a first order block, which can take into account a start-up time constraint.

Heat Pump Consoclim¶

This model is used to determine the performances of a heat pump for different operating conditions. The ConsoClim model developed by the Ecole de Mines (Paris) is used. The model predicts the performances of the system with three polynomial laws. The parameters of the model are identified with manufacturer data.

The first and the second law (EIRFT and CAPFT) are used respectively to determine the COP and the heating capacity of the machine at full load. These two polynomial laws depend on the outside air temperature (T_a_out) and the temperature of the water at the exhaust of the condenser (T_w). The third law is used to determine the performances of the system at part load.

The primary side fluid is modeled by means of 1-D incompressible fluid flow model , in which a dynamic energy balance and static mass and momentum balances are applied on the fluid. The heat transfer in the primary fluid is modeled with a constant heat transfer coefficient. However, it can be changed to other heat transfer models through the HeatTransfer parameter in the fluid model.

Furthermore, the model includes a Boolean input connector on_off, which defines the operational state of the heat pump. In the equations, the Boolean input is translated to a variable value together with a first order block, which can take into account a start-up time constraint. The model also includes a real input connector W_dot_set to define the power that is injected to the heat pump. This connector is useful, for example, in the case where a heat pump is powered by the CHP electrical generation.

Thermal energy storage¶

Nodal model of a stratified tank with an internal heat exchanger and ambient heat losses. The following hypothesis are applied:

No heat transfer between the different nodes.
The internal heat exchanger is discretized in the same way as the tank: each cell of the heat exchanger corresponds to one cell of the tank and exchanges heat with that cell only.
Incompressible fluid in bboth the tank and the heat exchanger.

The tank is discretized using a modified version of the incompressible Cell1Dim model adding an additional heat port, which acts as a heating resistance. The model also includes a temperature sensor, whose height can be adjusted. The energy balance of the fluid in the tank is described by:

\[\dfrac{V_{tank}}{N} \rho \dfrac{dh}{dt} + \dot{M} (h_{ex}-h_{su}) - A_{hx} \dot{q}_{hx} = \dfrac{A_{amb}}{N} \dot{q}_{amb} + \dot{Q}_{res}\]

The heat exchanger is modeled using the Flow1Dim component and a wall component. The energy balance of the fluid in the heat exchanger is described by:

\[V_i \rho \dfrac{dh}{dt} + \dot{M} (h_{ex}-h_{su}) = A_i \dot{q}\]

Crop Yield¶

Model variables and parameters¶

Parameters¶

Name	Units	Description
\(c_{\Gamma}\)	umol{CO₂} mol⁻¹{air} K⁻¹	Determines the effect of canopy temperature on the CO₂ compensation point
\(C_{Buf}^{max}\)	mg{CH₂O} m⁻²	Maximum buffer capacity
\(C_{Buf}^{min}\)	mg{CH₂O} m⁻²	Minimum amount of carbohydrates in the buffer
\(c_{BufFruit_1}^{max}\)	fruits plant⁻¹ s⁻¹	Regression coefficient
\(c_{BufFruit_2}^{max}\)	fruits plant⁻¹ s⁻¹ ºC⁻¹	Regression coefficient
\(c_{Dev_1}\)	s⁻¹	Regression coefficient
\(c_{Dev_2}\)	s⁻¹ ºC⁻¹	Regression coefficient
\(c_{Fruit,g}\)	–	Growth respiration coefficient of fruit
\(c_{Fruit,m}\)	s⁻¹	Maintenance respiration coefficient of fruit
\(c_{Leaf,g}\)	–	Growth respiration coefficient of leaf
\(c_{Leaf,m}\)	s⁻¹	Maintenance respiration coefficient of leaf
\(c_{RGR}\)	s	Regression coefficient for maintenance respiration
\(c_{Stem,g}\)	–	Growth respiration coefficient of stem
\(c_{Stem,m}\)	s⁻¹	Maintenance respiration coefficient of stem
\(E_j\)	J mol⁻¹	Activation energy for \(J^{POT}\)
\(G^{max}\)	mg{CH₂O} fruit⁻¹	Potential fruit weight
\(H\)	J mol⁻¹	Deactivation energy
\(J_{25Leaf}^{max}\)	umol{e-} m⁻²{leaf} s⁻¹	Maximum rate of electron transport at 25ºC for the leaf
\(LAI^{max}\)	m²{leaf} m⁻²	Maximum leaf area index
\(M_{CH2O}\)	mg umol⁻¹	Molar mass of CH₂O
\(M_{CO2}\)	mg umol⁻¹	Molar mass of CO₂
\(n_{Dev}\)	–	Number of development stages
\(n_{Plants}\)	plants m⁻²	Plants density in the greenhouse
\(Q_{10,m}\)	–	\(Q_{10}\) value for temperature effect on maintenance respiration
\(r_{BufFruit}^{max,FrtSet}\)	mg{CH₂O} m⁻² s⁻¹	Carbohydrate flow from buffer to the fruits above which fruit set is maximal
\(R_g\)	J mol⁻¹ K⁻¹	Molar gas constant
\(rg_{Fruit}\)	mg{CH₂O} m⁻² s⁻¹	Potential fruit growth rate coefficient at 20ºC
\(rg_{Leaf}\)	mg{CH₂O} m⁻² s⁻¹	Potential leaf growth rate coefficient at 20ºC
\(rg_{Stem}\)	mg{CH₂O} m⁻² s⁻¹	Potential stem growth rate coefficient at 20ºC
\(S\)	J mol⁻¹ K⁻¹	Entropy term
\(SLA\)	m²{leaf} mg⁻¹{CH₂O}	Specific leaf area index
\(T_{25,K}\)	K	Reference temperature at 25ºC
\(\alpha\)	umol{e⁻} umol⁻¹{photons}	Conversion factor from photons to electrons
\(\eta_{CO2,Air-Stom}\)	–	Conversion factor from the CO₂ concentration in the greenhouse air to the stomata
\(\eta_{C-DM}\)	mg{DM} mg⁻¹{CH₂O}	Conversion factor from carbohydrate to dry matter
\(\eta_{Glob-PAR}\)	umol{photons} J⁻¹	Conversion factor from global radiation to PAR
\(\theta\)	–	Degree of curvature of the electron transport rate
\(\tau\)	s	Time constant to calculate the 24-h mean temperature

Variables¶

Name	Units	Description
\(B\)	d⁻¹	Steepness of the curve
\(C_{Buf}\)	mg m⁻²	Carbohydrates in the buffer
\(C_{Fruit[j]}\)	mg m⁻²	Amount of fruit carbohydrates in the development stages
\(C_{Leaf}\)	mg m⁻²	Carbohydrate weight of the leaves
\(C_{Leaf}^{max}\)	mg m⁻²	Maximum allowed carbohydrates stored in the leaves
\(C_{Stem}\)	mg m⁻²	Carbohydrate weight of stems and roots
\(CO2_{Stom}\)	umol{CO₂} mol⁻¹{æir}	CO₂ concentration in the stomata
\(DM_{Har}\)	mg m⁻²	Accumulated harvested fruit dry matter
\(FGP\)	d	Fruit growth period
\(GR[j]\)	mg u⁻¹ d⁻¹	Daily potential growth rate per fruit in a development stage
\(g_{T_{Can24}}\)	–	Growth rate dependency to temperature
\(h_{C_{Buf}}^{m_{C,AirBuf}}\)	–	Inhibition coefficient of the photosynthesis rate by saturation of the leaves with carbohydrates
\(h_{C_{Buf}}^{m_{C,BufOrg}}\)	–	Inhibition coefficient of insufficient carbohydrates in the buffer
\(h_{T_{Can}}\)	–	Inhibition coefficient of non-optimal instantaneous canopy temperature
\(h_{T_{Can24}}\)	–	Inhibition coefficient of non-optimal 24-h mean canopy temperature
\(h_{T_{CanSum}}\)	–	Inhibition coefficient of crop development stage
\(h_{T_{CanSum}}^{m_{N,Fruit}}\)	–	Inhibition coefficient of fruit flow
\(J\)	umol{e-} m⁻² s⁻¹	Electron transportation rate
\(J_{POT}\)	umol{e-} m⁻² s⁻¹	Potential rate of electron transport at 25ºC for the canopy
\(J_{25Can}^{max}\)	umol{e-} m⁻² s⁻¹	Maximum rate of electron transport at 25ºC for the canopy
\(LAI\)	m²{leaf} m⁻²	Leaf area index
\(M\)	d	Fruit development time in days where GR is maximal
\(\dot{m}_{C,AirBuf}\)	mg m⁻² s⁻¹	Net photosynthesis rate
\(\dot{m}_{C,AirCan}\)	mg m⁻² s⁻¹	Carbohydrate flow between the air and the canopy
\(\dot{m}_{C,BufAir}\)	mg m⁻² s⁻¹	Growth respiration of the total plant
\(\dot{m}_{C,BufFruit}\)	mg m⁻² s⁻¹	Carbohydrate flow from the buffer to the fruits
\(\dot{m}_{C,BufFruit[1]}\)	mg m⁻² s⁻¹	Carbohydrate flow to the first development stage
\(\dot{m}_{C,BufFruit[j]}\)	mg m⁻² s⁻¹	Carbohydrate flow from the buffer to the remaining development stages
\(\dot{m}_{C,BufLeaf}\)	mg m⁻² s⁻¹	Carbohydrate flow from the buffer to the leaves
\(\dot{m}_{C,BufStem}\)	mg m⁻² s⁻¹	Carbohydrate flow from the buffer to the stems and roots
\(\dot{m}_{C,Fruit[j]Fruit[j+1]}\)	mg m⁻² s⁻¹	Carbohydrate flow through the fruit development stages
\(\dot{m}_{C,FruitAir}\)	mg m⁻² s⁻¹	Total maintenance respiration of the fruits
\(\dot{m}_{C,FruitAir,g}\)	mg m⁻² s⁻¹	Growth respiration of the fruits
\(\dot{m}_{C,FruitAir[j]}\)	mg m⁻² s⁻¹	Maintenance respiration of the fruits at the development stages
\(\dot{m}_{C,FruitHar}\)	mg m⁻² s⁻¹	Carbohydrate outflow of the last fruit development stage
\(\dot{m}_{C,LeafAir}\)	mg m⁻² s⁻¹	Maintenance respiration of the leaves
\(\dot{m}_{C,LeafAir,g}\)	mg m⁻² s⁻¹	Growth respiration of the leaves
\(\dot{m}_{C,LeafHar}\)	mg m⁻² s⁻¹	Leaf haverst rate
\(\dot{m}_{C,StemAir}\)	mg m⁻² s⁻¹	Maintenance respiration of the stems and roots
\(\dot{m}_{C,StemAir,g}\)	mg m⁻² s⁻¹	Growth respiration of the stems and roots
\(\dot{m}_{N,BufFruit[1]}\)	u m⁻² s⁻¹	Fruit set in the first development stage
\(\dot{m}_{N,BufFruit[1]}^{max}\)	u m⁻² s⁻¹	Maximum fruit set in the first development stage
\(\dot{m}_{N,Fruit[j]Fruit[j+1]}\)	u m⁻² s⁻¹	Fruit flow through the fruit development stages
\(N_{Fruit[j]}\)	u m⁻²	Number of fruits in the development stage j
\(P\)	umol m⁻² s⁻¹	Gross photosynthesis rate at canopy level
\(PAR_{Can}\)	umol{photon} m⁻² s⁻¹	Total PAR absorbed by the canopy computed in the solar model
\(R\)	umol m⁻² s⁻¹	Photorespiration during the photosynthesis process
\(r_{dev}\)	s⁻¹	Fruit development rate
\(RGR_{Fruit[j]}\)	s⁻¹	Net relative growth rate of fruits
\(RGR_{Leaf}\)	s⁻¹	Net relative growth rate of leaves
\(RGR_{Stem}\)	s⁻¹	Net relative growth rate of stems and roots
\(T_{Can}\) (\(T_{Can,K}\))	ºC (K)	Canopy temperature
\(T_{Can}^{24}\) (\(T_{Can,K}^{24}\))	ºC (K)	Mean 24-h canopy temperature
\(T_{Can}^{Sum}\) (\(T_{Can,K}^{Sum}\))	ºC (K)	Temperature sum
\(t_{j,FGP}\)	d	Number of days after fruit set for development stage j
\(W_{Fruit_1}^{Pot}\)	mg u⁻¹	Potential dry matter per fruit in first fruit development stage
\(\Gamma\)	umol{CO₂} mol⁻¹{air}	CO₂ compensation point
\(\eta\) _BufFruit	d m² mg⁻¹	Conversion factor to ensure that \(\dot{m}_{C,BufFruit}\) equals the sum of the carbohydrates that flow to the different fruit development stages

Model description¶

A dynamic tomato crop yield model was implemented to account for the effects of the indoor climate on crop growth and thereby on the harvested dry matter. Although crop growth is related to photosynthesis, most of the existent crop models directly relate these two with the absence of a carbohydrate buffer. The function of the buffer is to store the carbohydrates from the photosynthesis (inflow) and to distribute them to the plant organs (outflow). It has a maximum capacity, above which carbohydrates cannot be stored anymore, and a lower limit, below which the carbohydrate outflow stops. Thus, the in- and out-flows depend on the level of carbohydrates in the buffer and thereby, may not be simultaneous. An approach based on not considering the buffer neglects the non-simultaneous character of the flows. For example, it can neglect the crop growth after dusk, when photosynthesis stops but there may still be carbohydrate distribution if the buffer level is higher than its lower limit. The presence of a carbohydrate buffer is thus important when modeling crop growth.

Models with a common carbohydrate buffer are available in the current literature (e.g. [9], [18], [23], [25], [29]). In this work, a recent yield model [33] developed and validated for a variety of temperatures has been implemented. The model computes the carbohydrates distribution flows in the presence of a buffer, as shown in the figure below.

Schematic representation of the crop yield model adapted from [33]. Boxes define state variables (blocks), semi-state variables (dotted blocks) and carbohydrate flows (valves). Arrows define mass fluxes (solid lines) and information fluxes (dotted lines). For the purpose of readability, the grey box is a simplified scheme of the mass flows in the fruit development stages. A more detailed scheme of the latter can be found in the next figure.

To that end, the model applies carbohydrates mass balances on the buffer, fruits, stems and leaves. The mass balance on the buffer is defined by:

\[\dot{C}_{Buf} = \dot{m}_{C,AirBuf} - \dot{m}_{C,BufFruit} - \dot{m}_{C,BufLeaf} - \dot{m}_{C,BufStem} - \dot{m}_{C,BufAir}\]

The carbohydrates produced by the photosynthesis are stored in the buffer. Whenever carbohydrates are available in the buffer, carbohydrates flow to the fruits, leaves and stems. This flow stops when the buffer reaches its capacity lower limit. In a similar way, the carbohydrate inflow from photosynthesis stops if the buffer reaches its maximum capacity.

The available carbohydrates for fruit growth are distributed to the different fruit development stages. The fruit growth period, defined as the time between fruit set and fruit harvest, is modeled using the “fixed boxcar train” method [22]. With this method, fruit development is modeled by distinguishing several successive fruit development stages, at which the number of fruits and the amount of carbohydrates are described. As shown in the figure below, the number of fruits and carbohydrates are considered to flow from one stage to the next with a specific development rate.

Schematic representation of the fruit development model adapted from [33]. Boxes define state variables (blocks) and carbohydrate flows (valves). Arrows define mass fluxes (solid lines) and information fluxes (dotted lines).

The carbohydrates stored in a fruit development stage j are described by:

\[\dot{C}_{Fruit[j]} = \dot{m}_{C,BufFruit[j]} + \dot{m}_{C,Fruit[j-1]Fruit[j]} - \dot{m}_{C,Fruit[j]Fruit[j+1]} - \dot{m}_{C,FruitAir[j]}\]

for j = 1,2 … n_dev, where n_dev is the total number of fruit development stages. For the first development stage, the carbohydrate inflow from the previous stage \(\dot{m}_{C,Fruit[j-1]Fruit[j]}\) is zero. For the last development stage, the carbohydrate outflow to the next stage \(\dot{m}_{C,Fruit[j]Fruit[j+1]}\) is \(\dot{m}_{C,FruitHar}\).

The evolution of the number of fruits at a fruit development stage j is described by:

\[\dot{N}_{Fruit[j]} = \dot{m}_{N,Fruit[j-1]Fruit[j]} - \dot{m}_{N,Fruit[j]Fruit[j+1]}\]

for j = 1,2 … n_dev.

The evolution of the carbohydrates stored in the leaves is described by:

\[\dot{C}_{Leaf} = \dot{m}_{C,BufLeaf} - \dot{m}_{C,LeafAir} - \dot{m}_{C,LeafHar}\]

The Leaf Area Index (LAI), defined as the leaf area per unit of ground area, i.e. of greenhouse floor, is a semi-state variable of the model and determined by:

\[LAI = SLA \cdot C_{Leaf}\]

where SLA is the specific leaf area, whose value can be found in the literature.

The evolution of carbohydrates stored in the stems and roots is defined by:

\[\dot{C}_{Stem} = \dot{m}_{C,BufStem} - \dot{m}_{C,StemAir}\]

The harvested tomato dry matter (DM) is assumed to evolve with a continuous harvest rate. Thus, the accumulated DM equals the carbohydrate outflow from the last fruit development stage:

\[\dot{DM}_{Har} = \eta_{C-DM} \dot{m}_{C,FruitHar}\]

The development stages of the crop are defined by the evolution of the canopy temperature:

\[\dot{T}_{Can}^{Sum} = 86400^{-1} T_{Can}\]

Finally, the 24 hour mean canopy temperature is determined by a 1 st order approach:

\[\dot{T}_{Can}^{24} = \tau^{-1} (T_{Can}-T_{Can}^{24})\]

Model flows¶

Canopy photosynthesis¶

The net photosynthesis rate, which is equal to the gross photosynthesis rate minus the photorespiration, is described by:

\[\dot{m}_{C,AirBuf} = M_{CH2O} h_{C_{Buf}}^{m_{C,AirBuf}} (P-R)\]

The inhibition of the photosynthesis rate by saturation of the leaves with carbohydrates occurs when the carbohydrate amount in the buffer exceeds its maximum storage capacity. The inhibition is described by:

\[ h_{C_{Buf}}^{m_{C,AirBuf}} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 0,& C_{Buf} &> C_{Buf}^{max} \\ 1,& C_{Buf} &\leq C_{Buf}^{max} \end{array} \right. \]

Photosynthesis rate at the canopy level is described by:

\[P = \dfrac{J(CO2_{Stom}-\Gamma)}{4(CO2_{Stom}+2\Gamma)}\]

The photorespiration is described by:

\[R = P \dfrac{\Gamma}{CO2_{Stom}}\]

The electron transport rate is function of the potential rate and the PAR absorbed by the canopy, which is computed in the solar model and used in xx umol{photons} m⁻² s⁻¹.

\[J = \dfrac{J^{POT} + \alpha PAR_{Can} - \sqrt{(J^{POT}+\alpha PAR_{Can})^2 - 4 \theta J^{POT} \alpha PAR_{Can}}}{2\theta}\]

The potential electron transport rate depends on temperature:

\[J^{POT} = J_{25Can}^{max} e^{E_j \dfrac{T_{Can,K}-T_{25,K}}{R_g T_{Can,K} T_{25,K}}} \dfrac{1+e^\dfrac{S T_{25,K}-H}{R_g T_{25,K}}}{1+e^\dfrac{S T_{Can,K}-H}{R_g T_{Can,K}}}\]

where

\[J_{25Can}^{max} = LAI \cdot J_{25Leaf}^{max}\]

The CO₂ concentration in the stomata is expressed as a fraction of the CO₂ concentration in the greenhouse air [15]:

\[CO2_{Stom} = \eta_{CO2,Air-Stom} CO2_{Air}\]

Finally, the CO₂ compensation point (\(\Gamma\)) is described by:

\[\Gamma = \dfrac{J_{25Leaf}^{max}}{J_{25Can}^{max}} c_{\Gamma} T_{Can} + 20 c_{\Gamma} \left( 1 - \dfrac{J_{25Leaf}^{max}}{J_{25Can}^{max}} \right)\]

The carbohydrate flow to the fruits, leaves and stems¶

The carbohydrate flow from the buffer to the fruits is function of the potential of fruit growth coefficient, the effect of temperature on the flow and the inhibition factors (\(0<h<1\)):

\[\dot{m}_{C,BufFruit} = h_{C_{Buf}}^{m_{C,BufOrg}} h_{T_{Can}} h_{T_{Can24}} h_{T_{CanSum}} g_{T_{Can24}} rg_{Fruit}\]

The carbohydrate flows from the buffer to the leaves and stem are not influenced by the instantaneous temperature and are therefore described by:

\[\dot{m}_{C,BufLeaf} = h_{C_{Buf}}^{m_{C,BufOrg}} h_{T_{Can24}} g_{T_{Can24}} rg_{Leaf}\]

\[\dot{m}_{C,BufStem} = h_{C_{Buf}}^{m_{C,BufOrg}} h_{T_{Can24}} g_{T_{Can24}} rg_{Stem}\]

The inhibition of the carbohydrate flow to the fruits, leaves or stems caused by insufficient carbohydrates in the buffer is defined by its lower limit, which is equal to 5% of the buffer’s maximum capacity:

\[ h_{C_{Buf}}^{m_{C,BufOrg}} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 0,& C_{Buf} &\leq C_{Buf}^{min} \\ 1,& C_{Buf} &> C_{Buf}^{min} \end{array} \right. \]

Crop growth is also inhibited by non-optimal levels of the instantaneous and 24-hour mean temperature. The inhibitions coefficients for these temperatures (\(h_{T_{Can}}\) and \(h_{T_{Can24}}\)) can be described by two trapezoid functions (solid lines in figure below). Since the functions are non-differentiable, they have been smoothed (dotted lines in figure below) to make them suitable for dynamic simulation.

At the first development stage of the plant, all carbohydrates are used for leaf and stem growth. The first development stage is thus a vegetative stage. When a given temperature sum is reached, the generative stage starts and the carbohydrates are divided over the fruits, leaves and stems. The fruit growth rate is assumed to start at zero and increase linearly to full potential with increasing temperature sum:

\[ h_{T_{CanSum}} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 0,& T_{Can}^{Sum} &\leq T_{Start}^{Sum} & \\ \dfrac{T_{Can}^{Sum}}{T_{End}^{Sum}},& T_{Start}^{Sum} &< T_{Can}^{Sum} &\leq T_{End}^{Sum} \\ 1,& T_{Can}^{Sum} &> T_{End}^{Sum} & \end{array} \right. \]

Since at the start of the generative stage \(T_{Can}^{Sum}\) is zero, the temperature sum when the generative stage starts \(T_{Start}^{Sum}\) is zero. Moreover, it is assumed that the fruit growth rate is maximal after one fruit growth period. Thus, the temperature sum when the fruit growth rate is at full potential \(T_{End}^{Sum}\) is 1035 ºC.

The temperature effect on the growth rate coefficient is considered to be linear, as described by [21]:

\[g_{T_{Can24}} = 0.047 T_{Can24} + 0.060\]

The coefficient is 1 at 20ºC beacause the growth rate coefficients were defined at 20ºC.

The fruit flow to fruit development stages¶

The number of fruits in the development stages (\(N_{Fruit[j]}\)) depend on the fruit set of the first development stage (\(\dot{m}_{N,BufFruit[1]}\)) and the fruit flow to the remaining development stages (\(\dot{m}_{N,Fruit[j]Fruit[j+1]}\)). The fruit set of the first development stage depends on the carbohydrate flow from the buffer to the fruits and on the maximum fruit set:

\[ \dot{m}_{N,BufFruit[1]} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } \dfrac{\dot{m}_{C,BufFruit}}{r_{BufFruit}^{max,FrtSet}} \dot{m}_{N,BufFruit[1]}^{max},& \dot{m}_{C,BufFruit} &\leq r_{BufFruit}^{max,FrtSet} \\ \dot{m}_{N,BufFruit[1]}^{max},& \dot{m}_{C,BufFruit} &> r_{BufFruit}^{max,FrtSet} \end{array} \right. \]

where the maximum fruit set is depends on the plant density and the canopy temperature:

\[\dot{m}_{N,BufFruit[1]}^{max} = n_{Plants} \left( c_{BufFruit_1}^{max} + c_{BufFruit_2}^{max} T_{Can24} \right)\]

The fruit flow through the fruit development stages is computed based on the fixed boxcar train mechanism of [22]:

\[\dot{m}_{N,Fruit[j]Fruit[j+1]} = r_{Dev} n_{Dev} h_{T_{CanSum}}^{m_{N,Fruit}} N_{Fruit[j]}\]

for \(j=1,2...n_{Dev}\). The fruit development rate is defined by [21]:

\[r_{Dev} = c_{Dev_1} + c_{Dev_2} T_{Can24}\]

The fruit flow inhibition coefficient, used to ensure that the fruits stay in vegetative stage at the first development stage, is described by:

\[ h_{T_{CanSum}}^{m_{N,Fruit}} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 0,& T_{Can}^{Sum} &\leq T_{Start}^{Sum} \\ 1,& T_{Can}^{Sum} &> T_{Start}^{Sum} \end{array} \right. \]

The carbohydrate flow to fruit development stages¶

The amount of fruit carbohydrates in the development stages (\(C_{Fruit[j]}\)) depend on the carbohydrate flow from the buffer to a development stage (\(\dot{m}_{C,BufFruit[j]}\)) and the carbohydrate flow through the development stages (\(\dot{m}_{C,Fruit[j]Fruit[j+1]}\)).

The carbohydrate flow through the development stages is described by:

\[\dot{m}_{C,Fruit[j]Fruit[j+1]} = r_{Dev} n_{Dev} C_{Fruit[j]}\]

The carbohydrate flow to the first fruit development stage depends on the fruit set and the potential dry matter per fruit in the development stage one, as described by:

\[\dot{m}_{C,BufFruit[1]} = W_{Fruit[1]}^{POT} \dot{m}_{N,BufFruit[1]}\]

The carbohydrate flow from the buffer to the remaining fruit development stages, depends on the number of fruits, the fruit growth rate and the available carbohydrates for fruit growth (i.e. the total amount from the buffer minus the amount used in the first stage):

\[\dot{m}_{C,BufFruit[j]} = \eta_{BufFruit} N_{Fruit[j]} GR[j] \left( \dot{m}_{C,BufFruit} - \dot{m}_{C,BufFruit[1]} \right)\]

for \(j=2,3...n_{Dev}\), where the conversion factor \(\eta_{BufFruit}\) is described by:

\[\eta_{BufFruit} = \dfrac{1}{\sum\limits_{j=1}^{j=n_{Dev}} N_{Fruit[j]} GR[j]}\]

The fruit growth rate depends on the fruit development stage and is described by [21]:

\[GR[j] = G^{max} e^{-e^{-B(t_{[j]}^{FGP}-M)}} B e^{-B(t_{[j]}^{FGP}-M)}\]

where

\[FGP = \dfrac{1}{r_{Dev} 86400}\]

\[M = -4.93 + 0.548 FGP\]

\[B = \dfrac{1}{2.44 + 0.403 M}\]

\[t_{[j]}^{FGP} = \dfrac{(j-1)+0.5}{n_{Dev}} FGP\]

Growth and maintenance respiration¶

The growth respiration of the total plant (\(\dot{m}_{C,BufAir}\)) is equal to the sum of the growth respiration of the fruits, leaves and stems. For the leaves, the growth respiration is defined by:

\[\dot{m}_{C,LeafAir,g} = c_{Leaf,g} \dot{m}_{C,BufLeaf}\]

The maintenance respiration of the leaves is described by:

\[\dot{m}_{C,LeafAir} = c_{Leaf,m} Q_{10,m}^{0.1(T_{Can24} - 25)} C_{Leaf} \left( 1-e^{-c_{RGR} RGR} \right)\]

The growth and maintenance respirations of the fruits and stems are described analogously to the previous equations.

Leaf pruning¶

It is assumed that leaves are pruned if the simulated LAI exceeds the maximum allowed value. The maximum allowed amount of stored carbohydrate in the leaves is described by:

\[C_{Leaf}^{max} = \dfrac{LAI^{max}}{SLA}\]

The leaf harvest rate is determined by:

\[ \dot{m}_{C,LeafHar} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 0,& C_{Leaf} &< C_{Leaf}^{max} \\ C_{Leaf}-C_{Leaf}^{max},& C_{Leaf} &\geq C_{Leaf}^{max} \end{array} \right. \]

Flows¶

The Flows package contains models of the flows that are encountered in a greenhouse system. It is organized in eight sub-packages that model the heat, vapor mass and CO₂ mass transfer, as well as fluid flow. It also contains a sub-package of interfaces, which defined the type of connectors used in the library.

Heat flows¶

The HeatTransfer sub-package models the heat flows in a greenhouse, which can originate from convection at surfaces, ventilation processes, conduction at the soil and long-wave infrared radiation (FIR).

Convective and conductive heat flows are function of the heat exchange coefficient and are described by:

\[\dot{q}_{cnv,ij} = U_{ij}(T_i-T_j)\]

FreeConvection¶

Typically, convective processes in greenhouses are governed by free convection. In this case, the Nusselt (Nu) number describing the convective exchange process can be defined as a function of the Rayleigh (Ra) number [5]. The heat exchange coefficients are therefore modeled based on the Nu-Ra relation, as presented in [10].

Upward or downward heat exchange by free convection from an horizontal or tilted surface. In the greenhouse, this comprises the convection at the cover, the thermal screen and at the floor.

\[ U_{AirFlr} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 1.7(T_{Flr} - T_{Air})^{0.33},& T_{Flr} &> T_{Air} \\ 1.3(T_{Air} - T_{Flr})^{0.25},& T_{Flr} &\leq T_{Air} \end{array} \right. \]

\[U_{AirScr} = 1.7 u_{Scr} |T_{Air}-T_{Scr}|^{0.33}\]

\[U_{AirCov} = 1.7 (T_{Air}-T_{Cov})^{0.33} cos(\phi)^{-0.66}\]

PipeConvection¶

Free convection at the pipes with the greenhouse air. For the pipes situated close to the canopy and the floor, the heat exchange is considered to be hindered, compared to a pipe in free air. The heat exchange coefficients of these forced processes are modeled by experimental results [6].

\[U_{FreePipeAir} = 1.28 \pi \psi_{Pipe}^{0.75} l_{Pipe} |T_{Pipe}-T_{Air}|^{0.25}\]

\[U_{HinderedPipeAir} = 1.99 \pi \psi_{Pipe} l_{Pipe} |T_{Pipe}-T_{Air}|^{0.32}\]

PipeConvection_N¶

This model is a variant of the PipeConvection model, in which a single-port is replaced by a multi-port, thus enabling the computation of the heat flow when using discretized pipes models.

CanopyFreeConvection¶

The leaves heat exchange by free convection with the air is function of the Leaf Area Index (LAI) and the leaves heat transfer coefficient.

\[U_{CanAir} = 2 \alpha_{LeafAir} LAI\]

OutsideAirConvection¶

Convective heat exchange at the cover with the outside air. Convection, driven by wind speed, is considered to be forced.

\[ U_{CovOut} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } (2.8+1.2 v_w) \dfrac{1}{cos(\phi)},& v_w &< 4 \\ 2.5 v_w^{0.8} \dfrac{1}{cos(\phi)},& v_w &\geq 4 \end{array} \right. \]

SoilConduction¶

The only conductive flow considered in greenhouse modeling is the conduction through the greenhouse soil. The soil under the greenhouse floor represents a big thermal capacity with a poor thermal conductivity. The floor surface can show temperature variations of 10 K during a day. To be able to describe the temperature gradient, the soil is modeled in several layers, using the following heat exchange coefficient.

\[U_{So(j-1)So(j)} = \dfrac{2}{h_{So(j-1)}/\lambda_{So(j-1)}+h_{So(j)}/\lambda_{So(j)}}\]

FreeVentilation¶

Convective flows caused by ventilation processes are modeled based on the air exchange rate f_ij between two air volumes i and j, as described by:

\[U_{vent,AirOut} = \rho_{Air} c_{p,Air} (f_{AirOut}+f_{leakage})\]

The library offers two models (NaturalVentilationRate_1 and NaturalVentilationRate_2) to compute the air ventilation rate caused by natural ventilation with the outside air. The models are based on two different models from the literature. By default, f_AirOut is described by the model NaturalVentilationRate_2, which is based on Boulard and Baille (1993).

The leakage rate through the greenhouse structure is dependent on the wind speed and the leakage coefficient of the greenhouse, characteristic of its structure. It can be described by:

\[ f_{leakage} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 0.25 c_{leakage},& v_w &< 0.25 \\ v_w c_{leakage},& v_w &\geq 0.25 \end{array} \right. \]

AirThroughScreen¶

The air ventilation between the main and top air zones is caused by two mechanisms: the air through the openings in the fabric of the screen and the air through a gap when the screen is opened. Balemans, 1989studied the temperature driven air exchange through fully closed screens (u_Scr =1) and derived a fitted function through experimental data. When the screen is open (u_Scr <1), the air exchanged through the gap, caused by density difference, will dominate the exchange through the screen. This exchange was theoretically modeled by Miguel, 1998 using the Navier-Stokes equation. Combining the air flow through the screen and through the gap, the total air ventilation rate between the air and top zones is described by:

\[U_{vent,AirTop} = \rho_{Air} c_{p,Air} f_{AirTop}\]

\[f_{AirTop} = u_{Scr} K_{Scr} |T_{Air}-T_{Top}|^{0.66} + \dfrac{1-u_{Scr}}{\overline{\rho}_{Air}} \sqrt{0.5 \overline{\rho}_{Air} W (1-u_{Scr}) g |\rho_{Air}-\rho_{Top}|}\]

Radiation_T4¶

The thermal radiation, i.e. the electromagnetic radiation emitted between two bodies i and j as a result of their temperatures, is described by the Stefan-Boltzman equation:

\[\dot{q}_{rad,ij} = \epsilon_{i} \epsilon_{j} F_{ij} \sigma (T_i^4-T_j^4)\]

The view factors of the greenhouse elements are computed according to [10] in each component model (i.e. the components described in the Greenhouse section). The exchange with the sky, whose temperature is estimated from meteorological data by an approach proposed in [10], is also considered.The emission coefficients are characteristic of the surfaces. For the greenhouse elements, the following values are proposed:

Component Value

Glass cover 0.84

Pipes 0.88

Canopy leaves 1.00

Concrete floor 0.89

Thermal screen 1.00

Radiation_N¶

This model is a variant of the Radiation_T4 model, in which a single-port is replaced by a multi-port, thus enabling the computation of the radiative flow when using discretized pipes models.

Vapor flows¶

MV_cnv_condensation¶

The vapor mass transfer caused by condensation at a surface is linearly related to its convective heat exchange coefficient by a conversion factor. In the greenhouse, condensation may occur at the lower side of the cover and the thermal screen. The model excludes evaporation at these surfaces.

\[ \dot{m}_{v,ij} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 0,& P_{v,i} &< P_{v,j} \\ 6.4·10^{-9} U_{ij} (P_{v,i}-P_{v,j}),& P_{v,i} &\geq P_{v,j} \end{array} \right. \]

Because of the direction nature of this flow, the model is not reversible and must be connected as following: air (filled port) - surface (empty port).

MV_cnv_evaporation¶

The vapor mass transfer caused by evaporation at a surface is linearly related to its convective heat exchange coefficient by a conversion factor. In the greenhouse, evaporation may occur at the upper side of the thermal screen. The model excludes condensation at this surface. By allowing a mass flow rate from the upper surface of the screen to the top air compartment, the model assumes that the screen is capable of transporting water through its fabric. Water is transported from the lower side to the upper and storage of water in the screen is neglected. Therefore, evaporation from the upper side is only possible when condensation takes place at the lower side. Moreover, the evaporation rate must be lower or equal than the condensation rate.

\[ \dot{m}_{v,ScrTop} = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 0,& P_{v,Scr} &< P_{v,Top} \\ min \left( 6.4·10^{-9} U_{ScrTop}, 6.4·10^{-9} U_{AirTop} \dfrac{P_{v,Air}-P_{v,Scr}}{P_{v,Scr}-P_{v,Top}} \right) (P_{v,Scr}-P_{v,Top}),& P_{v,Scr} &\geq P_{v,Top} \end{array} \right. \]

Because of the direction nature of this flow, the model is not reversible and must be connected as following: surface (filled port) - air (empty port).

MV_Ventilation¶

Mass transfer also occurs in ventilation processes. The computation of the vapor flow exchanged by ventilation from the indoor to outdoor air is modeled by:

\[\dot{m}_{v,ij} = \dfrac{M_{Water} f_{ij}}{R} \left( \dfrac{P_{v,i}}{T_i}-\dfrac{P_{v,j}}{T_j} \right)\]

where f_ij can be f_AirOut or f_TopOut.

MV_AirThroughScreen¶

The computation of the vapor flow exchanged by ventilation between the main and top air zones is described similarly than in MV_Ventilation, but applying the air exchange coefficient f_AirTop.

MV_CanopyTranspiration¶

The vapor flow from the canopy to the greenhouse air originates from a phase interface somewhere inside the cavities of a leaf. The resistance to vapor transport from the canopy leaves to the greenhouse air is made of an internal resistance and a boundary layer resistance [30]. According to the latter, the canopy transpiration can be defined by:

\[\dot{m}_{v,CanAir} = \dfrac{2 \rho_{Air} c_{p,Air} LAI}{ \Delta H \gamma (r_b+r_s)} (P_{v,Can}-P_{v,Air})\]

Heat and Vapor Flows¶

In the vapor model, all flows result from convective exchange processes and in order to compute them, the heat exchange coefficient of these convective processes is used. Therefore, in order to reduce the number of connections and inputs when building a greenhouse model, the heat and vapor models of convective processes are lumped into single models in which both computations are performed simultaneously. The HeatAndVaporTransfer sub-package includes the lumped models.

Convection_Condensation¶

Combines the equations of FreeConvection and MV_cnv_condensation.

Convection_Evaporation¶

Combines the equations of FreeConvection and MV_cnv_evaporation.

Ventilation¶

Combines the equations of FreeVentilation and MV_ventilation.

AirThroughScreen¶

Combines the equations of AirThroughScreen from the HeatTransfer sub-package and MV_AirThroughScreen.

VentilationRates¶

This sub-package contains two different models for computing the air exchange rate in convective processes, and a model for computing the air rate due to a forced ventilation system.

NaturalVentilationRate_1: based on [20]. The air exchange rate is modeled in function of the wind and temperature. The contribution of the temperature driven ventilation in the total ventilation is small but can be important during nighttime and winter. The wind speed driven ventilation is computed differently for vents in the windward side and the leeside side. The air exchange is related to the wind speed and the opening of a window.

\[f_{AirOut} = 0.5 fr_{window} \sqrt{ \Phi_{wind}^2 + \Phi_{temp}^2}\]

\[\Phi_{wind} = \left( 2.29·10^{-2} (1- exp(-\theta /21.1) + 1.2·10^{-3} \theta exp(\theta /211) \right) A_{window} u_{wind}\]

\[\Phi_{temp} = C_f l/3 \sqrt{|g \beta \Delta T|} h^{1.5} \left[ \left( sin(\psi)-sin(\psi-\theta_l) \right)^{1.5} + \left( sin(\psi)-sin(\psi-\theta_w) \right)^{1.5} \right]\]

NaturalVentilationRate_2: based on [7]. The air ventilation ratedepends mainly on the windows opening (u_vent) and is influenced by the wind pressure coefficient and the coefficient of energy discharge caused by friction at the windows.

\[f_{AirOut} = \dfrac{u_{vent} A_{Roof} C_{d}}{2 A_{Flr}} \sqrt{g \dfrac{h_{vent}}{2} \dfrac{T_{Air}-T_{Out}}{\overline{T}} + C_w v_w^2 }\]

ForcedVentilationRate: The air exchange rate caused by a mechanical ventilation system is function of the speficic air flow capacity of the ventilation system and its control.

\[f_{ventForced} = U_{VentForced} \dfrac{\phi_{VentForced}}{A_{Floor}}\]

CO₂ flows¶

In the greenhouse, there are three CO₂ flows associated to the ventilation processes and two forced flows, i.e. the canopy consumption and the CO₂ enrichment.

MC_ventilation¶

The CO₂ flow accompanying an air flow is function of the air flow rate and can be described by:

\[\dot{m}_{c,ij} = f_{ij}(CO_{2,i}-CO_{2,j})\]

MC_AirCan¶

The greenhouse CO₂ net assimilation rate by the canopy is computed in the yield model and used as an input in this model.

Fluid flows¶

Cell1DimInc¶

Fluid flows are modeled using the finite volume approach by means of a discretized model for incompressible flow, adapted from [26]. The model distinguishes between two types of variables: cell and node variables. The basic fluid flow component is a cell in which the dynamic energy balance and static mass and momentum balances are applied. Node variables correspond to the inlet and outlet nodes of each cell. The relation between the cell and node values depend on the selected discretization scheme (upwind or central differences). In the cell, uniform velocity through the cross section and constant pressure are assumed. Axial thermal energy transfer is neglected.

Flow1DimInc¶

The overall flow model can be build by connecting several cells in series. The model is compatible with the textit{Media} package of the Modelica Standard Library, at the condition that the considered fluid is incompressible.

Control Systems¶

Greenhouses have high requirements on indoor climate control. The climate controller adjusts heating, ventilation and CO₂ supply to attain the desired climate. In this work, several control systems are developed, based on the control strategies proposed in the literature. In the library, depending on the nature of the strategies, two implementation approaches are distinguished:

State graphs: using the StateGraph package from the Modelica Standard Library, and
Proportional-integral-derivative (PID) controllers: using an ISA PID controller with anti-windup (PID model in the ControlSystems package of the Greenhouses Lirary), adapted from [26]).

The presented control strategies have been tested in [2], where several climate simulation results are presented and compared for three different cases.

Set-points definition¶

The determination of set-points is at the top of the functionality of the climate controller. Temperature set-points differ from day-time to night-time and are sometimes adapted to the level of radiation. CO₂ is supplied during daylight to enhance photosynthesis. As measured in xx Nilsen et al., 1983, different combinations of CO₂ concentration and air temperature lead to different photosynthesis rates. Although a sharp reduction in photosynthesis is measured at non-optimal temperatures, similar values are measured for close-to-optimal temperatures (i.e. optimal temperature ±5ºC). Therefore, temperature and CO₂ set-points can be optimized not only in terms of crop growth but also in terms of energy use. In fact, the definition of temperature set-points for optimal crop growth and energy use has been the subject of a substantial literature (e.g. [14]; [1]; [13]; [11]). However, since this work does not focus on climate set-points optimization, no innovative control is proposed. Instead, the strategy proposed in [1] is implemented in Python and the set-points are inputted as a time-series “.txt” file in the model. This strategy consists in minimizing energy consumption while maintaining a crop growth close to the maximal growth rate by:

computing a 2-D array of photosynthesis rates for a range of CO₂ and temperature values at a given PAR,
selecting the pairs of CO₂ and temperature that ensure at least 80% of the photosynthesis rate (being 100% the maximum value of the 2-D array), and
defining the temperature and CO₂ set-point (\(T_{Air,SP}\) and \(CO2_{Air,SP,th}\)) as the pair in the previous selection with the lowest temperature.

In this work, a PI controller is responsible for adjusting the heating power output by varying the mass flow rate of the heating pipes according to the difference between the air temperature set-point and the actual value. The control strategy for CO₂ is based on a maximal supply rate, defined by the capacity of the CO₂ enrichment system. This capacity is function of the CO₂ source, which is commonly a combination of fossil fuel combustion gases and CO₂ stored in liquid phase. While respecting the enrichment capacity, the supply rate is adapted to attain the CO₂ set-point. However, in high ventilation conditions, CO₂ enrichment is commonly reduced due to the high exchange rate to the outside air. To take this into account, the theoretical CO₂ set-point proposed by the control strategy is modified so that it decreases proportionally with the increase in the ventilation rate. This is done as defined by:

\[CO2_{Air,SP} = f(u_{vent}) \cdot \left( CO2_{Air,SP,th} - CO2_{Ext}^{Min} \right) + CO2_{Ext}^{Min}\]

where

\[ f(u_{vent}) = \left\{ \begin{array}{ @{}% no padding l@{\quad}% some padding r@{}% no padding >{{}}r@{}% no padding >{{}}l@{}% no padding } 1- \dfrac{u_{vent}}{u_{vent}^{Max}},& u_{vent} &< u_{vent}^{Max} \\ 0,& u_{vent} &\geq u_{vent}^{Max} \end{array} \right. \]

Supplementary lighting¶

The most popular lamp type for commercial supplementary lighting in horticulture is high pressure sodium (HPS) lamps. HPS lamps are the most efficient in the PAR spectrum range, with an emission highly concentrated between 500 and 650 nm. HPS lighting is not designed for frequent cycling because it dramatically reduces lamp lifespan. Thus, regardless of the control method, it is best to set up constraints to operate lighting for extended periods. The implemented control strategy for the lighting is based on the following:

Lighting window: allow lights to be turned on between \(h_{illu,ON}^{min}\) and \(h_{illu,ON}^{max}\) (e.g. 5 AM and 10 PM).
Lighting set-point: allow lights to be turned on during the lighting window if light levels decrease below \(I_{illu,ON}\) (e.g. 40 Wm-2) and to be turned off when light levels increase above \(I_{illu,OFF}\) (e.g. 120 Wm-2).
Light accumulation: turn off lights or do not allow turning them on if the daily accumulated light exceeds \(I_{acc}^{max}\) (e.g. 5 kWh).
Proving time: light levels must be below the set-point for at least \(t_{illu,proving}\) (e.g. 30 minutes).
Minimum on time: to prevent cycling, lights must remain on for minimum \(t_{illu,ON}^{min}\) (e.g. 2 hours) once they are turned on, regardless of other conditions.

The strategy sets up a time window for lighting, during which a lighting set-point condition is applied. The proving time and minimum on time strategies are implemented to prevent cycling.

Windows aperture¶

Windows in the greenhouse can be opened either for dehumidification or for cooling the greenhouse. Excessive humidity can cause fungal diseases or physiological disorders [16]. Humidity in greenhouses is controlled by means of a strategy related to a constraint rather than a specific set-point. The constraint is based on allowing a maximum value of relative humidity in the air, commonly set at 85%. The most common technique for dehumidification is the combination of ventilation and heating. Although this technique is energy consuming and thus expensive, dehumidifying systems based on refrigerant cycles, e.g. heat pumps, have not proved to be economically feasible [31]. Windows are also used for cooling the air in the case of excessive temperatures, since they have a negative impact on the harvest rate. For example, in [33] the harvest rate at daylight temperatures of 40ºC was 54.5% of that at 25ºC. Moreover, temperatures above 25ºC can penalize fruit quality e.g. size and color [31]. In this work, a proportional (P) controller is used to select the opening of the windows according to the following:

Air sanitation: A maximum value \(RH_{vent,ON}\) is allowed for humidity.
Air cooling: A maximum value \(T_{vent,ON}\) is allowed for air temperature.

Thermal screen closure¶

As previously mentioned, thermal losses to the outside can be reduced from 38% to 60% by using a thermal screen [4]. This capability of reducing thermal losses is defined by the screen material, which is selected according to the climate of the region. In fact, depending on the nature of the screen, the light transmission coefficient can vary from 15% to 88%. Thus, when drawn, the screen reduces considerably the transmitted light above the canopy. The most conventional method to operate the screen is therefore to deploy it at sunset, when heating demand becomes significant, and remove it at sunrise, to profit from the available sun light. The removal of the screen must be operated progressively to avoid a thermal shock. A way of further reducing energy consumption is to deploy the screen before sunset or to delay the removal until after sunrise. However, this implies a loss of crop production caused by a reduction on the available light. A good approach would be to study the threshold between energy saving and production loss in order to define the optimal deployment and removing times. However, estimating the reduction of plant growth is a complex task that, although it has been the object of some studies (e.g. [1]; [4]), it commonly has many uncertainties and thereby requires many assumptions. A simpler approach is to define the deployment of the screen in function of the outside irradiation. In fact, the photosynthetic activity of the plant achieves its maximal potential about one hour after sunrise and diminishes just before sunset [17]. In Dutch-conditions, deploying the screen after 50 Wm-2 (instead of 5 Wm-2 usually practiced) allows to decrease energy consumption by an extra 3% without penalizing crop growth [12]. Depending on the night, a small temporary opening of the screen may be necessary to regulate humidity or temperature. As defined by the ventilation rate through the screen xx , screen gaps increase the air exchange between the main and top air zones and therefore decrease temperature and humidity.

In this work, the developed screen control strategy is based on the following:

Opening/closing set-point: the screen is opened (closed) if irradiation increases (decreases) above (below) a certain value \(I_{Scr,ON}\) (\(I_{Scr,OFF}\)) (e.g. 35 Wm-2).
Opening/closing time: the screen is opened progressively by 1% per minute (with an interval pause of 3 minutes) followed by a full opening after 30%. This approach has proven not to generate cold air flows on top of the canopy [17]. The opening percentage and pause time is adapted to the outside weather. The time to fully open the screen is about 45 min to 60 min in cold days and 30 min in mild days.
Humidity gap: the screen is opened in steps of 1% (with an interval of 3 min) up to a maximum of 4% if relative humidity exceeds its set-point, as recommended in [12].

Examples¶

Greenhouse_1¶

This example intends to illustrate the simulation of the climate in a greenhouse. The greenhouse model is built by connecting all the greenhouse components to the energy and mass flows presents in a greenhouse. As it can be distinguished, the greenhouse modeled in this example consists of two levels of heating circuits, roof windows (but not side vents), natural ventilation (no forced ventilation) and a movable thermal screen. It should be noted that, when the screen is drawn, the air of the greenhouse is divided in two zones, i.e. below and above the screen. These zones are modeled separately (models air and air_Top) and their climate is assumed to be homogeneous. The models parameters have been set to typical values for Venlo-type greenhouse construction design dedicated to tomato crop cultivation. The greenhouse floor area and the mean greenhouse height are set in two individual block sources.

The simulated greenhouse is located in Belgium and the simulation period is from the 10th of December to the 22nd of November. Two input data files are required:

Weather data: The input weather data for the simulation period is extracted from a TMY for Brussels and can be found in ‘Greenhouses/Resources/Data/10Dec-22Nov.txt’. The file contains data for the outside air temperature, air pressure, wind speed and global irradiation. The sky temperature, previously computed in a Python script, is also included in this file.
Climate control set-points: The temperature and CO₂ set-points for the simulation period are calculated according to the strategy presented in Control Systems and can be found in ‘Greenhouses/Resources/Data/SP_10Dec-22Nov.txt’.

These .txt files are accessed by means of TMY_and_control and SP_new, which are two CombiTimeTables models from the Modelica Standard Library.

The goal of this example is to show the energy flows interacting in a greenhouse. Thus, no generation units are included. Instead, the heating pipes are connected to a water source and sink model. The model includes the following controls:

PID_Mdot: A PI controller adjusts the output mass flow rate of the water source connected to the heating pipes by compaing the air temperature set-point and present value.
PID_CO2: A PI controller adjusts the output of the CO₂ external source by comparing the actual CO₂ concentration of the air to its set-point.
Ctrl_SC: A state graph adjusts the screen closure (SC) according to the strategy presented in Control Systems. The real inputs must be connected to the air relative humidity, the outdoor temperature, the indoor air temperature set-point and the usable hours of the screen. The usable hours are 1h30 before dusk, 1h30 after dawn and during night. In the global he global outside irradiation
vents: A PI controller adjusts the opening of the windows according to the strategy presented in Control Systems. The opening depends mainly on the indoor air relative humidity and temperature.
OnOff: controls the ON/OFF operation of the supplementary lighting according to the strategy presented in Control Systems. The control output, previously computed in a Python script, is input as a .txt file by means of the TMY_and_control model.

GlobalSystem_1¶

This example illustrates the energy flows of a greenhouse connected to a CHP unit and a thermal storage tank. The storage tank is connected as an open buffer, i.e. the CHP is not directly connected to the greenhouse but only to the buffer. The greenhouse is modeled by the ready-to-use greenhouse model included in Greenhouse/Unit sub-package, in which the parameters are already set to typical values for Venlo-type greenhouse construction design dedicated to tomato crop cultivation. The CHP unit is sized to provide the peak load of the thermal demand of the greenhouse, which can be extracted from the example Greenhouse_1. The thermal energy storage unit is sized considering a storage period of 8h and given a nominal temperature difference (\(\Delta T\)) between the inlet and the outlet.

The CHP consumes natural gas and is providing the heating and electrical demand of the greenhouse. A heat-driven control decides when to run the CHP unit depending on the greenhouse demand and the storage level. The limited electrical load of greenhouses implies an excess in electricity generation when using CHP units, which have a power-to-heat ratio usually close to one. In this example, the excess of electricity is sold back to the grid. Since the control is heat-driven and there is no electrical storage, it is possible that at some points there is electrical demand but the CHP is not running. At this moment, electricity is bought from the grid. The buying and selling prices of electricity, as well as the buying price of gas, are considered constant during the whole simulation. The difference between the buying and selling prices of electricity occurs because, in the absence of subsidies, the electricity excess fed to the grid is remunerated at a price close to the wholesale price of electricity, which is much lower than the retail price (i.e. buying price).

The following sankey diagram illustrates the total production and consumption of the CHP and the greenhouse over the simulation period. The values are in \(MWh m^{-2}\), and they include: the total CHP consumption of gas (\(E_{gas}\)), the CHP thermal production (\(E_{th,CHP}\)), electrical production (\(E_{el,CHP}\)) and losses (\(E_{loss,CHP}\)); the greenhouse thermal demand (\(E_{th,gh}\)) and electrical demand (\(E_{el,gh}\)); and the electricity bought (\(E_{b,grid}\)) and sold back (\(E_{s,grid}\)) to the grid.

GlobalSystem_2¶

In the previous example GlobalSystem_1, a considerable part of the produced electricity is sold back to the grid. As explained, this electricity, in the absence of subsidies, is remunerated at a price close to the wholesale price of electricity. Because the retail price of electricity is significantly higher than wholesale price, prosumers have a clear advantage at maximizing their level of self-consumption [27]. In order to evaluate the potential of such activity and to illustrate the use of the library, a third simulation is developed. In this example, the particular case of maximizing the self-consumption level of the system simulated in GlobalSystem_1 is evaluated. To that end, a heat pump is added to the system. The heating demand of the greenhouse is produced by both the heat pump and the CHP. Therefore, the CHP capacity is reduced, and the heat pump is sized so that its nominal electrical capacity is equal to the excess of electricity production of the CHP in nominal conditions. Similarly than in GlobalSystem_1, a heat-driven control decides when to run the CHP and the heat pump. Electricity excess not consumed by the heat pump is sold to the grid. The greenhouse electrical demand not covered by the CHP is covered by the grid. The electricity and gas prices are the same than in GlobalSystem_1.

Similarly than for GlobalSystem_1, the total consumption and production of the system over the simulation period is illustrated in the following sankey diagram. In this diagram the electrical consumption (\(E_{el,HP}\)) and thermal production (\(E_{th,HP}\)) of the heat pump are also indicated.

A detailed comparison of the results obtained by GlobalSystem_1 and GlobalSystem_2 is presented in [3]. Because of the lower CHP capacity, the gas consumption is lower in GlobalSystem_2 than in case GlobalSystem_1. In a similar way, the total net electrical generation of case GlobalSystem_2 is also lower. As expected, the heat pump reduces the amount the electricity sold back to the grid.

References¶

[1]	J. M. Aaslyng, J. B. Lund, N. Ehler, and E. Rosenqvist. IntelliGrow: a greenhouse component-based climate control system. Environmental Modelling & Software, 18(7):657–666, September 2003. doi:10.1016/S1364-8152(03)00052-5.

[2]	Queralt Altes-Buch and Vincent Lemort. Modeling framework for the simulation and control of greenhouse climate. In Proceedings of the 10th International Conference on System Simulation in Buildings. Liege, December 2018.

[3]

Queralt Altes-Buch, Sylvain Quoilin, and Vincent Lemort. Modeling and control of CHP generation for greenhouse cultivation including thermal energy storage. In Proceedings of the 31st international conference on efficiency, cost, optimization, simulation and environmental impact of energy systems, 13. Guimaraes, Portugal, June 2018. URL: https://orbi.uliege.be/handle/2268/226034.

[4]	B. J. Bailey. Control strategies to enhance the performance of greenhouse thermal screens. Journal of Agricultural Engineering Research, 40(3):187–198, July 1988. doi:10.1016/0021-8634(88)90206-5.

[5]	Luc Balemans. Assessment of criteria for energetic effectiveness of greenhouse screens. PhD thesis, Agricultural University, Ghent, 1989.

[6]	G.P.A Bot. Greenhouse climate : from physical processes to a dynamic model. PhD thesis, Wageningen University, 1983.

[7]	T. Boulard and A. Baille. A simple greenhouse climate control model incorporating effects of ventilation and evaporative cooling. Agricultural and Forest Meteorology, 65(3):145–157, August 1993. doi:10.1016/0168-1923(93)90001-X.

[8]	Kinsell L. Coulson. Solar and Terrestrial Radiation. Elsevier, 1975. ISBN 978-0-12-192950-3. doi:10.1016/B978-0-12-192950-3.X5001-3.

[9]	E. Dayan, H. van Keulen, J. W. Jones, I. Zipori, D. Shmuel, and H. Challa. Development, calibration and validation of a greenhouse tomato growth model: I. Description of the model. Agricultural Systems, 43(2):145–163, January 1993. doi:10.1016/0308-521X(93)90024-V.

[10]	H.F. De Zwart. Analyzing energy-saving options in greenhouse cultivation using a simulation model. PhD thesis, Wageningen University, 1996.

[11]	J Dieleman, E Meinen, L.F.M. Marcelis, Zwart , and E.J. Van Henten. Optimisation of CO2 and Temperature in Terms of Crop Growth and Energy Use. Acta Horticulturae, October 2005. doi:10.17660/ActaHortic.2005.691.16.

[12]	J.A. Dieleman and F.L.K. Kempkes. Energy screens in tomato: determining the optimal opening strategy. Acta Horticulturae, pages 599–606, October 2006. doi:10.17660/ActaHortic.2006.718.70.

[13]	J.A. Dieleman, L.F.M. Marcelis, A. Elings, T.A. Dueck, and E. Meinen. Energy Saving in Greenhouses: Optimal Use of Climate Conditions and Crop Management. Acta Horticulturae, pages 203–210, October 2006. doi:10.17660/ActaHortic.2006.718.22.

[14]	A. Elings, H. F. de Zwart, J. Janse, L. F. M. Marcelis, and F. Buwalda. Multiple-day temperature settings on the basis of the assimilate balance: a simulation study. Acta Horticulturae, 2006.

[15]	J. R. Evans and G. D. Farquhar. Modelling canopy photosynthesis from the biochemistry of the C3 chloroplast. In Modeling crop photosynthesis - from biochemistry to canopy (Eds. K.J. Boote, R.S. Loomis), pages 1–15. CSSA Madison, Wisconsin, USA, 1991.

[16]	R. I. Grange and D. W. Hand. A review of the effects of atmospheric humidity on the growth of horticultural crops. Journal of Horticultural Science, 62(2):125–134, January 1987. doi:10.1080/14620316.1987.11515760.

[17]	Ariane Grisey and Eric Brajeul. Serres chauffées: réduire ses dépenses énergétiques. Centre technique interprofessionnel des fruits et légumes (CTIFL), 2007.

[18]	E. Heuvelink. Tomato growth and yield : quantitative analysis and synthesis. PhD thesis, Wageningen University, 1996.

[19]	I. Impron, S. Hemming, and G. P. A. Bot. Simple greenhouse climate model as a design tool for greenhouses in tropical lowland. Biosystems Engineering, 98(1):79–89, September 2007. doi:10.1016/j.biosystemseng.2007.03.028.

[20]	T. de Jong. Natural ventilation of large multispan greenhouses. PhD thesis, Wageningen University, 1991.

[21]	A. N. M. de Koning. Development and dry matter distribution in glasshouse tomato : a quantitative approach. PhD thesis, Wageningen University, Wageningen, 1994. URL: http://library.wur.nl/WebQuery/wurpubs/24498.

[22]	P. A. Leffelaar and Th J. Ferrari. Some elements of dynamic simulation. In Simulation and systems management in crop protection, pages 19–45. Pudoc, Wageningen, 1989. URL: http://library.wur.nl/WebQuery/wurpubs/8869.

[23]	R. Linker, I. Seginer, and F. Buwalda. Description and calibration of a dynamic model for lettuce grown in a nitrate-limiting environment. Mathematical and Computer Modelling, 40(9):1009–1024, November 2004. doi:10.1016/j.mcm.2004.12.001.

[24]

Weihong Luo, Hendrik Feije de Zwart, Jianfeng DaiI, Xiaohan Wang, Cecilia Stanghellini, and Chongxing Bu. Simulation of Greenhouse Management in the Subtropics, Part I: Model Validation and Scenario Study for the Winter Season. Biosystems Engineering, 90(3):307–318, March 2005. doi:10.1016/j.biosystemseng.2004.11.008.

[25]	L. F. M Marcelis, E Heuvelink, and J Goudriaan. Modelling biomass production and yield of horticultural crops: a review. Scientia Horticulturae, 74(1):83–111, April 1998. doi:10.1016/S0304-4238(98)00083-1.

[26]	Sylvain Quoilin, Adriano Desideri, Jorrit Wronski, Ian Bell, and Vincent Lemort. ThermoCycle: A Modelica library for the simulation of thermodynamic systems. In Proceedings of the 10th International Modelica Conference 2014. 2014.

[27]

Sylvain Quoilin, Konstantinos Kavvadias, Arnaud Mercier, Irene Pappone, and Andreas Zucker. Quantifying self-consumption linked to solar home battery systems: Statistical analysis and economic assessment. Applied Energy, 182:58–67, November 2016. URL: http://www.sciencedirect.com/science/article/pii/S0306261916311643, doi:10.1016/j.apenergy.2016.08.077.

[28]	J Ross. Radiative transfer in plant communities. In Vegetation and Atmosphere (Ed. J. L. Monteith), pages 13–55. Academic Press, London, UK, 1975.

[29]	Ido Seginer, Christian Gary, and Marc Tchamitchian. Optimal temperature regimes for a greenhouse crop with a carbohydrate pool: A modelling study. Scientia Horticulturae, 60(1):55–80, December 1994. doi:10.1016/0304-4238(94)90062-0.

[30]	C. Stanghellini. Transpiration of greenhouse crops : an aid to climate management. PhD thesis, Wageningen University, 1987.

[31]	Laurent Urban and Isabelle Urban. La production sous serre: La gestion du climat. Volume 1. Lavoisier, 2nd edition, 2010.

[32]	R. J. C. van Ooteghem. Optimal Control Design for a Solar Greenhouse. IFAC Proceedings Volumes, 43(26):304–309, January 2010. doi:10.3182/20101206-3-JP-3009.00054.

[33]	B. H. E. Vanthoor, P. H. B. de Visser, C. Stanghellini, and E. J. van Henten. A methodology for model-based greenhouse design: Part 2, description and validation of a tomato yield model. Biosystems Engineering, 110(4):378–395, December 2011. doi:10.1016/j.biosystemseng.2011.08.005.

[34]

B. H. E. Vanthoor, C. Stanghellini, E. J. van Henten, and P. H. B. de Visser. A methodology for model-based greenhouse design: Part 1, a greenhouse climate model for a broad range of designs and climates. Biosystems Engineering, 110(4):363–377, December 2011. doi:10.1016/j.biosystemseng.2011.06.001.

Welcome to Greenhouses’s documentation!¶

Downloading Greenhouses¶

License¶

Main developers¶

Contents¶

Overview¶

Getting started¶

Library structure¶

Graphical user interface¶

Components¶

Greenhouse¶

General nomenclature¶

Subscripts¶

Variables¶

Solar model¶

Cover¶

Canopy¶

Air¶

Air Top¶

Floor¶

Thermal screen¶

Illumination¶

Heating pipes¶

HVAC¶

CHP¶

Heat Pump¶

Heat Pump Consoclim¶

Thermal energy storage¶

Crop Yield¶

Model variables and parameters¶

Parameters¶

Variables¶

Model description¶

Model flows¶

Canopy photosynthesis¶

The carbohydrate flow to the fruits, leaves and stems¶

The fruit flow to fruit development stages¶

The carbohydrate flow to fruit development stages¶

Growth and maintenance respiration¶

Leaf pruning¶

Flows¶

Heat flows¶

FreeConvection¶

PipeConvection¶

PipeConvection_N¶

CanopyFreeConvection¶

OutsideAirConvection¶

SoilConduction¶

FreeVentilation¶

AirThroughScreen¶

Radiation_T4¶

Radiation_N¶

Vapor flows¶

MV_cnv_condensation¶

MV_cnv_evaporation¶

MV_Ventilation¶

MV_AirThroughScreen¶

MV_CanopyTranspiration¶

Heat and Vapor Flows¶

Convection_Condensation¶

Convection_Evaporation¶

Ventilation¶

AirThroughScreen¶

VentilationRates¶

CO2 flows¶

MC_ventilation¶

MC_AirCan¶

Fluid flows¶

Cell1DimInc¶

Flow1DimInc¶

Control Systems¶

Set-points definition¶

Supplementary lighting¶

Windows aperture¶

Thermal screen closure¶

Examples¶

Greenhouse_1¶

GlobalSystem_1¶

GlobalSystem_2¶

References¶

CO₂ flows¶