Surplus thermal energy model of greenhouses and coefficient analysis for effective utilization

If a greenhouse in the temperate and subtropical regions is maintained in a closed condition, the indoor temperature commonly exceeds that required for optimal plant growth, even in the cold season. This study considered this excess energy as surplus thermal energy (STE), which can be recovered, stored and used when heating is necessary. To use the STE economically and effectively, the amount of STE must be estimated before designing a utilization system. Therefore, this study proposed an STE model using energy balance equations for the three steps of the STE generation process. The coefficients in the model were determined by the results of previous research and experiments using the test greenhouse. The proposed STE model produced monthly errors of 17.9%, 10.4% and 7.4% for December, January and February, respectively. Furthermore, the effects of the coefficients on the model accuracy were revealed by the estimation error assessment and linear regression analysis through fixing dynamic coefficients. A sensitivity analysis of the model coefficients indicated that the coefficients have to be determined carefully. This study also provides effective ways to increase the amount of STE.


Introduction
Modern greenhouses control plant growth conditions, such as temperature, humidity, carbon dioxide concentration, light and nutrients, and produce highquality products irrespective of the outdoor climate. However, controlling the climate requires a large amount of energy, especially for heating during the cold season. The energy costs diminish the benefits of greenhouse farming, and the use of fossil fuels for the energy contributes to global warming. Therefore, many studies have been conducted to reduce the use of fossil 2 control conditions and characteristics of crops were considered. Using the determined STE model, the coefficients that mostly affect the model accuracy were investigated, and a sensitivity analysis on the coefficients were conducted. We propose effective methods for estimating the amount of STE and its utilization.

Test greenhouse equipped with the STE utilization system
The test greenhouse was single span, and the covering material was double-layer glass with 7 mm thickness. The floor area was 99.36 m 2 with a width of 6.9 m and 14.4 m in length. The wall and the ridge were 3.4 m and 5.1 m high, respectively. Accordingly, the covering area was 267.33 m 2 , and the interior volume was 422.28 m 3 . The greenhouse was at latitude 37.2° N, in an east-west orientation rotated 20° clockwise. The STE utilization system consists of a heat pump, low and high temperature heat storage tanks (LST and HST, respectively), and 10 fan-coil units (FCUs) as shown in Fig. 1. When STE is generated, it is recovered by circulating cold water in the LST to the FCUs in the greenhouse. STE absorbed into the cold water is transferred to the HST as hot water using the heat pump. When the greenhouse needs heating, the hot water in the HST circulates into the FCUs in the greenhouse. For more detail system information, refer to Yang et al. (2012a).

Greenhouse operation and data measurement
The test greenhouse system was operated during the winter (December 2010 -February 2011). The cultivated plants were roses (Rosa hybrida) in pots, with 1200 pots laid on three growth trays. The set-points for the temperature control were 26°C in the day time and 16°C at night time. Thus, STE was recovered when the greenhouse indoor temperature exceeded 26°C. The amount of STE was calculated using Eq. [1] by measuring the water temperature and flow rate.
where E gr is the energy amount that supplied for or recovered from the greenhouse; Q wt is the water flow rate (m 3 /min); c wt and ρ wt are the specific heat of water (J/ kg/K) and the water density (kg/m 3 ), respectively; t en and t lv are the temperature (°C) entering and leaving the greenhouse at time T (min), respectively; f m-s is the unit factor for converting from minute to second (60 s/min).
al., 2009). This STE energy can be efficiently utilized with closed greenhouses (Bakker et al., 2006;Vadiee & Martin, 2012, 2013a. Yang et al. (2012a) investigated the daily change of greenhouse indoor temperatures and designed a greenhouse system to utilize STE using a heat pump system, which recovers and stores STE and then supplies it when necessary. This system operated in a greenhouse having a 100 m 2 floor area and recovered 8.71 GJ for 3 months in a cold season (Yang & Rhee, 2013). The economic feasibility of STE depends on the amount of it compared with the cost of its utilization system. Thus, estimating the amount of STE is inevitably necessary when deciding on the installation of the STE utilization system. Many studies on greenhouse thermal modeling have been conducted considering environmental control systems (Sharma et al., 1999;Chou et al., 2004;Hepbasli, 2011), thermal curtains and earth-air heat exchangers (Shukla et al., 2006), neural network or computer-aided modeling (Ferreira et al., 2002;Han et al., 2009), the effect of greenhouse orientation (Sethi, 2009), solar energy (Hamdan et al., 1992;Abdel-Ghany & Al-Helal, 2011), and thermal storages (Gauthier et al., 1997;Najjar & Hasan, 2008;Lee et al., 2011;Vadiee & Martin, 2012). In this study, a thermal model for the greenhouse was adopted into a test greenhouse including an STE utilization system. To fit the model at high accuracy, several coefficients related to environmental status, greenhouse structural specification,  Figure 1. Schematic diagram of the surplus thermal energy utilization system and the sensor positions. Arrows show the thermal energy flow. Circles indicate the position of sensors of air temperature (T a ), floor temperature (T f ), water temperature (T w ), water flow rate (F w ) and solar radiation (I s ).
3 Surplus thermal energy model of greenhouses and its coefficient analysis ratios of direct and diffuse radiation to global radiation, respectively; I s is the intensity of global solar radiation; τ gh is the solar radiation transmissivity of greenhousecovering materials; and T is time. Because the area affected by diffuse radiation (A dff ) is equal to the total covering area (A cv In winter, when the solar altitude is low, direct radiation reaches approximately half of the total covering area, A drt = 0.5A cv . The ratio of diffuse radiation to global radiation (Φ dff ) is 0.15 on a clear day, and thus, Φ drt = 0.85 (Kim et al., 1997), and f r becomes 0.575. The greenhouse covers τ gh varies according to the type, orientation and covering material of the greenhouse. The solar position and outside weather conditions also contribute to the τ gh . A variety of studies on the τ gh have considered the incidence angle of solar radiation (Bowman, 1970), weather conditions (Leonidopoulos, 2000;Cabrera et al., 2009), greenhouse orientation (Papadakis et al., 1998;Li et al., 2000), structural members (Critten, 1987), and dust and dirt on the covering (Geoola et al., 1998). Because the results on τ gh vary from 25% to 59%, according to experimental conditions, applying a specific τ gh for the energy input equation was difficult. To resolve this difficulty in τ gh determination, factors affecting the τ gh were divided into the incidence angle of solar radiation and the dust/dirt of the covering. Thus, the overall τ gh of a greenhouse is expressed by where τ pd is the solar radiation transmissivity for the perpendicular beam on the surface, and f ang and f dd are the incidence angle factor and the dust/dirt factor for the solar radiation transmissivity, respectively.
The water flow rate and the water temperature were measured using a water flow meter (HMD40-1b, Shinil meter tech, Korea), accuracy ±2%, and temperature sensors (NTC-10kΩ), accuracy ±0.3°C. For the greenhouse indoor and outdoor temperatures, temperature sensors (PT-100Ω), accuracy ±0.3°C, with the fan-aspirated radiation shield (Yang et al., 2012b) were installed in the center 1.5 m above the floor in the greenhouse and 5.5 m above the ground beside the greenhouse. Solar radiation was measured using a pyranometer (SYE-420M2007PM4, Shinyoung, Korea). Specifications for the pyranometer are provided in Table 1. The measurement points are presented in Fig. 1. Data were recorded every minute.

STE modeling
The elements of the greenhouse energy balance for the STE model are shown in Fig. 2. Solar radiation (I s ) transmits through the greenhouse covers, and its energy reduces. Total solar energy transmitted into the greenhouse (E so ) is converted to sensible heat by the plants and floor (E input ). If this sensible heat is greater than thermal energy loss of the greenhouse (E loss ), the greenhouse contains the STE (E sp ) as their difference. The detailed calculation procedure is elaborated in the following section.

Solar energy input
The energy input from solar radiation was determined by the covering areas, the ratios of the direct and the diffuse radiation to the global radiation, and light transmissivity of the covers. The ratios of the direct and the diffuse radiation to global radiation display a stable value at a specific location on a clear day (Al- Mohamad, 2004;Ahmed et al., 2009). Therefore, the energy input from solar radiation, E so , is expressed as where A drt and A dff are the areas affected by direct and diffuse radiation, respectively; Φ drt and Φ dff are the where inside temperature (t in ) in the thermal energy loss element (dE loss ) becomes the heating temperature if the heating system normally operates; i is the time sequence in minutes. Because it took approximately one hour to convert solar energy input into thermal energy, thermal energy input was compared with thermal energy loss after one hour. This time difference was empirically found through temperature measurement.
After the first step is satisfied, the indoor temperature is raised by the increase of solar energy. This is the second step, expressed by where the indoor temperature (t in ) in the thermal energy loss element (dE loss ) equals t in (i -1); C gh is the heat capacity within the greenhouse. Eq.
[8] is calculated repeatedly until the indoor temperature is increased to the optimal growth temperature. The final i is the initial time that STE is generated. When the difference between the thermal energy input and loss is a positive value, STE exists. The amount of STE is calculated by Eq. [9], which is the final step of the STE model.
where T 1 starts with the final i at the second step. T 2 is the time before the amount of STE calculated between i, and i+1 is the negative value. During cloudy days, many time-sections of T 1 and T 2 can be observed.

Determination of coefficients in the STE model for the test greenhouse
Incidence angle factor for the solar radiation transmissivity The incidence angle factor for the solar radiation transmissivity was determined based on Bowman (1970). The solar radiation transmissivity was between 0.92 and 1.0 for incidence angles of less than 60°, whereas solar radiation was reduced exponentially as incidence angle increased above 60°.
The incidence angle of solar radiation on the greenhouse covering, α i , was calculated using

Conversion of solar energy input to thermal energy
Solar radiation transmitted into a greenhouse reaches plants or the floor and converts to thermal energy. In here, coefficients to convert solar energy into thermal energy are considered, and they are designated the sensible heat emission factors (SHEFs). In particular, this factor for plants is related to the Bowen ratio, which is the ratio of the sensible heat loss to the evaporative heat loss (Taiz & Zeiger, 1991). Consequently, the amount of thermal energy in the greenhouse is expressed as where f sn-plt and f sn-flr are SHEFs for the plant and the floor, respectively; and mPLAI is the modified projected leaf area index. The mPLAI is the ratio of total projected leaf area to the floor area, which was modified by the PLAI concept (Barclay, 1998). In this study, the test plant was rose in pots, which takes a small portion of the greenhouse inside volume. Thus, minor parameters such as latent heat, evaporation and plant respiration were intended to be integrated into the SHEF. If plants took a large portion of the greenhouse inside volume such as tomatoes and paprika, these parameters need to be considered.

Thermal energy loss of the greenhouse
For greenhouses with closed windows, thermal energy loss, E loss , occurs by convection, conduction and infiltration through the cover and the floor, which is expressed by where U gh is the overall heat transfer coefficient of the greenhouse including the infiltration effect, and f wind is the wind effect factor, which was set to 1.0 in this study. Although the wind effect factor is a linear function of wind speed, the effect of the wind speed is very small in the case of closed greenhouses with the double-layered glass covers except in strong winds (Kim et al., 1997); h flr is the heat transfer coefficient through the floor, and t in , t out and t flr are the temperatures of the greenhouse inside, outside and floor, respectively.

STE model
The STE model consists of three steps. The first step is to find the time that greenhouse heating is 5 Surplus thermal energy model of greenhouses and its coefficient analysis shown in Fig. 4. After 33 days, the leaves of adjacent pots overlapped. Because the projected leaf area is acquired using only the upper leaves, the mPLAI did not increase further after 33 days. The number of total test pots for three trays was 1200.

SHEF of the plants and the greenhouse floor
The SHEF coefficients were determined using Eq.
[15] which was established from Eq.

Evaluation of the STE model and analysis of the coefficients
After the STE model was determined for the test greenhouse, it was evaluated through the error between the measured and estimated amounts of STE and linear regression analysis. Moreover, the STE model includes dynamic coefficients that are changed by time or temperature, and they complicate the calculation of the model. Thus, sinα h = sinφ ⋅sinδ + cosφ ⋅cosδ ⋅cosα t [11] cosα i = sinα h ⋅cosθ + cosα h ⋅sinθ ⋅cos α b − α a ( ) [12] where δ is the declination of the sun and n day the day of the year, for example, n day on March 1 is 60; α h is the solar altitude angle; φ the greenhouse latitude; α t the hour angle (noon = 0); and θ the inclined angle of the surface; α a is the surface azimuth angle; and α b the solar azimuth angle.

Heat transfer coefficients of the greenhouse
The overall heat transfer coefficient through the cover, U gh , and the heat transfer coefficient through the floor, h flr , were determined using the energy balance equation, Eq. [13]. Areas and temperatures were determined through measurement, and the supplied thermal energy, E gr , was determined by Eq. [1]. Thus, the heat transfer coefficients could be determined by regression analysis of Eq.
[13]. The datasets were prepared at hour intervals.
where A flr is the floor area, and t flr the floor temperature.

Modified projected leaf area index (mPLAI)
The projected leaf area was measured by assuming that it was a circular disc with diameter equal to the average length between the ends of leaves. The projected leaf area in a pot fluctuated from 0.00785 m 2 to 0.0616 m 2 between pinching and blooming as shown in Fig. 3. As there were 45 days from pinching to blooming, mPLAI is expressed using Eq. [14].
where d is the number of days after pinching, and n pot the number of pots. The arrangement of the pots is  5b. Based on the solar radiation transmissivity changes (Bowman, 1970), the incidence angle factors were determined as shown in Fig. 5c and 5d for the wall and the roof, respectively. The dirt/dust factor for the solar radiation transmissivity has been reported as 0.98 for 18 months and 0.95 for 24 months for glass covers (Kim et al., 1997). The test greenhouse for this research has been operating more than 2 years, so the dirt/dust factor was set as 0.95. Solar radiation transmissivity for the perpendicular beam on the surface is 0.71 for doublelayer glass (Hanan, 1998). Therefore, the solar radiation transmissivity for the test greenhouse is expressed as τ gh = 0.95⋅0.71 f ang = 0.675 f ang [17] where the incidence angle factor, f ang , is the dynamic coefficient as shown in Figs. 5c and 5d.

Heat transfer coefficients of the test greenhouse
The overall heat transfer coefficient through covers, U cv , was determined to be 2.98 W/m 2 /K (p<0.001), but the floor coefficient was not significant (p=0.63). The determined U cv is similar to the coefficient (3 W/m 2 /K) presented by ASAE (1988) for a double-glass covering. The number of samples for the regression analysis was 25. through fixing these dynamic coefficients, the effects of the coefficients on the model accuracy were investigated.
A sensitivity analysis on the model coefficients was also conducted. When the coefficients in the model were considered as input and the estimated amounts of STE were considered as output, the effect of the coefficient on the estimation result was assessed. Among several methods for sensitivity analysis of building energy (Lam & Hui, 1996), a method which is useful for multiple sets of data was chosen, as shown in: where OP base is base output, equal to the estimated value, and IP base is base input, which is equal to the coefficient. ΔOP and ΔIP are variations of output and input, respectively. High sensitivity coefficients indicate that the model is more sensitive to them and their values must be chosen very carefully.

Incidence angle factor for the solar radiation transmissivity
The incidence angle of the solar radiation changes according to the date and time as shown in Figs. 5a and 7 Surplus thermal energy model of greenhouses and its coefficient analysis Using the SHEF, the relation of sensible and latent heat in the greenhouse indoor air, and the relation of heat transfer between the greenhouse indoor air and plants or floor were simply solved. Joudi & Farhan (2015) have analyzed the thermal energy between the indoor air and soil in the greenhouse through dividing soil layers into 0.1 m, 0.2 m and 0.5 m depth. Because all these relations depend on the temperature, it is explained that the SHEF depends on the temperature.

STE model for the test greenhouse and effects of coefficients on the estimation accuracy
The coefficients determined for the test greenhouse are listed in Table 2. Measured and estimated amounts of STE were compared as shown in Fig. 7. The amount SHEF A total of 198 data points of temperature and supplied thermal energy were collected and regression analysis was conducted using Eq.
[15] for every temperature interval. For example, data between 19.5°C and 20.4°C were used for analysis of 20°C. Thus, the number of data for one temperature condition was approximately 20. Figure 6 shows the SHEF according to temperature. The SHEF for the plants and the floor could be expressed by Eqs. [18] and [19] with R 2 of 0.502 and R 2 of 0.844, respectively.
These results show that plants and floor contribute to increase the greenhouse temperature. However, according to temperature increase, the plants retain the SHEF between 0.35 and 0.40, whereas the SHEF of the floor keeps decreasing, which means the floor works as a heat sink at high indoor temperature. In addition, if it is assumed that the floor works as a heat source at low indoor temperature, the result that the SHEF of the floor was maximum at 16°C could be explained. These results provide a hint that greenhouses have to be operated at low temperature to recover the STE effectively.  0.3-0.5 f sn-flr [2] 0.3-0.5 ρ air 1.2 kg/m 3 at 293K τ gh =f age ·f ang ·τ pd 0.675f ang τ pd 0.71 [1] Time-dependent coefficient.
[2] Temperature-dependent coefficient. . Moreover, even though these three coefficients were fixed at the same time (All fixing), the errors were not seriously increased. However, when the SHEFs were fixed at 25ºC or 27ºC, and the mPLAI was fixed as 0 or 1, the errors seriously increased. Similar results were obtained through the linear regression analysis (Fig. 10). The over/under-estimation ratios and the R 2 values were sensitively changed be-of STE decreased in January and sharply increased in February. Monthly errors were evaluated to be 17.9%, 10.4% and 7.4% for December, January and February, respectively. As a result of the linear regression analysis (Fig. 8), it was found that the STE was underestimated by 7.2% compared with the measured STE (Eq. [20]), and the R 2 value was 0.850. In this study, because the STE data higher than 200 MJ/day were not collected sufficiently, the proposed STE model has an applicable limitation for the high STE condition. [20] where y and x are the STE estimated and measured, respectively. As shown in Table 2, the incidence angle factor (f ang ), the SHEFs (f sn-plt and f sn-flr ) and mPLAI are time or temperature dependent. We investigated the effects of these coefficients on the estimation accuracy of the STE model. Because the incidence angle factor is generally one during the time period that STE exists, the model was re-evaluated with a fixed incidence angle factor of 1. In addition, the model was re-evaluated with fixing of the SHEFs as the values at 25ºC, 26ºC and 27ºC, respectively. The temperature set-point for recovering STE was 26ºC. The mPLAI was fixed as 0, 0.575 and 1 for ignoring plants, full of plants and ignoring floor, respectively. Figure 9 shows the monthly errors according to fixing the coefficients. When the incidence angle factor was fixed at 1, the SHEFs were fixed at 26ºC (SHEF at 26ºC), and plants grew at the normal plant density (mPLAI=0.575). The monthly errors were slightly changed compared with the errors when the all coef-   9 Surplus thermal energy model of greenhouses and its coefficient analysis improves the solar radiation transmissivity, has been suggested (Geoola et al., 1998), and reinforcing sealing such repairing tears in covers and preventing the leakage may be helpful. Cleaning the covers also increases the photosynthetic rate of crops, and reinforcing sealing is necessary to save heating energy. As conclusions, this study was conducted to propose a model to estimate the surplus thermal energy (STE). Energy balance equations were basically used, and the sensible heat emission factor (SHEF) from plants and floor was newly introduced to develop the STE model. Monthly errors of the model were evaluated to be 17.9%, 10.4% and 7.4% for December, January and February, respectively. The STE model shows that the solar radiation transmissivity and the overall heat transfer coefficient of greenhouse covers are sensitively influenced to the amount of STE. Because this model was developed and verified through a small glass greenhouse, it must be tested in a larger greenhouse. As future research, it is necessary to find the coefficients in the STE model for different conditions of greenhouses.
cause of the SHEFs and the mPLAI. Consequently, it was revealed that the error of the STE model mainly relies on these two coefficients.

Sensitivity analysis of the STE model
The results of the sensitivity analysis of the STE model are shown in Fig. 11. In January, the STE model was mostly sensitive to all coefficients, and this might be attributed to the small amount of STE. Small changes in STE result in a significant increase of the sensitivity coefficient. For the mPLAI and the SHEFs, the sensitivity was stable for coefficient changes. However, as the sensitivity of the STE model was highly influenced by the solar radiation transmissivity and the overall heat transfer coefficient in January, these coefficients must be carefully determined. This result also means that improvement of the solar radiation transmissivity and the heat loss through the greenhouse cover is vital to recovering more STE. As simple ways to achieve this, cleaning the greenhouse covers, which