Implementation of ventilation towers in a greenhouse established in low altitude tropical climate conditions: numerical approach to the behavior of the natural ventilation

Submitted on April 16, 2020 and accepted on September 21, 2020. 1 Corporación Colombiana de Investigación Agropecuaria AGROSAVIA, Centro de Investigación Tibaitata, Bogotá, Colombia. evillagran@agrosavia.co *Corresponding Author: evillagran@agrosavia.co Implementation of ventilation towers in a greenhouse established in low altitude tropical climate conditions: numerical approach to the behavior of the natural ventilation


INTRODUCTION
The increase in the world's population, together with the deterioration of natural resources and current and future food demands, is one of the frequent problems that receives the most attention within the working agendas of the governmental and non-governmental entities of the nations (Akrami et al., 2020). One of the alternatives of intensification of the agriculture more promoted for the food production is the agriculture under cover, at present it is estimated that for the production of fresh vegetables there are about 3,220,000 has under this type of productive system and of this total approximately 95% of the area are medium and low technology structures, within which one of the most relevant are passive type greenhouse structures, which are used in more than 115 countries worldwide (Baeza et al., 2020;Ezzaeri et al., 2020).
The main functions of a greenhouse are; to protect the crops from extreme climatic factors, to generate an adequate microclimate for the growth and development of the plants and to protect them from the attack of pests and diseases (Syed & Hachem, 2019;Villagrán et al., 2019). The benefits sought with the implementation of the crop under cover are; increases in crop yields, increase in the final quality of the harvested product, reduction of water and fertilizer consumption, reduction of pesticides and ensure agricultural production throughout the year (Pakari & Ghani, 2019a). Passive greenhouses are characterized by the lack of climate control equipment for tasks such as cooling and heating, therefore, this type of greenhouse is considered low cost and low environmental impact by not requiring energy from non-renewable sources for its operation (Graamans et al., 2018).
The management of the microclimate in passive greenhouses is carried out by opening and closing ventilation areas on the sides and roof of the greenhouse, a process known as natural ventilation. From this management of the ventilation areas, air flow movements are generated from physical phenomena such as natural convection by buoyancy or thermal effect and by differences in pressure between the external and internal environment known as wind effect (Bournet & Boulard, 2010;Espinoza et al., 2017). Through the natural ventilation phenomenon, an air exchange is induced between the external and internal environment of the greenhouse. This exchange allows the regulation of the thermal and humidity excesses in the internal environment of the greenhouse and it is also the only source of carbon enrichment in passive greenhouses (Singh et al., 2018).
The efficiency of natural ventilation depends on many factors such as: the geometry and size of the greenhouse, the area of the ventilation surfaces, the speed and direction of the outside wind, the thermal gradient between the indoor and outdoor environment, as well as the presence or not of insect screens also play an important role in the behavior of the airflow patterns (Baeza et al., 2009). Deficient natural ventilation rates are mainly characterized by inadequate airflow patterns within the greenhouses, thus generating microclimatic conditions that cause stress conditions in the plants that affect their growth and development, these stress conditions are usually more critical in hot climate regions .
Currently there is a need to seek strategies to improve the airflow patterns of passive greenhouses and the microclimate generated inside them using different renewable energy sources. A little explored alternative is the use of wind energy by coupling aeration towers to the sides of the greenhouse, these towers have the ability to generate air flow movements by taking advantage of the pressure differences that are created between the windward and leeward sides (Pakari & Ghani, 2019b). A technique for the study and optimization of air flow patterns inside greenhouses is the numerical simulation by computational fluid dynamics (CFD), this is a mature and robust technique with which it is possible to evaluate different alternatives in a virtual scenario of engineering problems that include transport phenomena and fluid movements (Dhiman et al., 2019;Mesmoudi et al., 2017).
The use of wind towers is not widespread in agricultural structures, in a recent study developed by Pakari & Ghani (2019b), the authors analyzed by means of CFD and an experimental validation in a wind tunnel the effect of a wind capture tower on the flow patterns inside a greenhouse type tunnel of 81 m 2 . This study only considered the wind effect of natural ventilation and the authors reported air flows with speeds inside the tunnel that varied between 0.41 and 1.81 m s -1 , for simulation scenarios with outside wind speeds between 1 and 5 m s -1 . The objective of this research was to evaluate, using an experimentally validated 2D CFD model, the effect that the implementation of an air capture tower in each side wall of the greenhouse would have on air flow patterns and thermal distribution. The numerical model considered the thermal and wind effect of the natural ventilation phenomenon.

Description of the greenhouse
The experiment was conducted in a naturally ventilated flat roof chapel type greenhouse (PG), built of structural steel and covered with polyethylene, located in the Colombian Caribbean region in the municipality of Seville, department of Magdalena, at the Caribbean Research Center of the Colombian Agricultural Research Corporation -AGROSAVIA. The dimensions of the greenhouse were 20 m in cross section and 25 m in length, which is equivalent to a covered area of 500 m 2 . The greenhouse was equipped with side and roof vents, these vents have an insect proof mesh ( Figure 1). This greenhouse model already built at full scale was used to carry out the comparisons of air flow and thermal distribution with the possible scenario of a greenhouse equipped with air capture towers (MPG). The capture and recording of experimental data will also be carried out on PG for the purpose of validating the CFD model.
The region of study presents climatological conditions where the average of the multiannual temperature for a period of 30 years is 28.16 °C, with a maximum and minimum average of 32.8 and 23.2 °C, respectively, the relative humidity presented an average behavior with values over 65% and an annual precipitation with average values of 1245.1 mm, the dominant wind speed values in the region of study vary between 0.5 and 3 m s -1 , with a dominant direction of in the E-W direction, which is perpendicular to the longitudinal axis and side walls of the greenhouse.

Description of air capture towers
The air capture towers proposed in this research should be attached to each of the lateral sides of the greenhouse studied ( Figure 2B). The towers can be built in the same structural and covering material in which the greenhouse is manufactured, in this case structural steel and polyethylene covering. The overall dimensions of the towers are 10 m and 9.5 m for maximum and minimum height, the width of the tower is 2.5 m, the length should be similar to the side length of the greenhouse in this case 20 m and at the top of the tower an area should be generated for air inlet into the greenhouse with a width of 3.5 m ( Figure 2A).

Numerical model
The non-linear partial differential equations known as Navier Stokes (RANS) averaged by Reynolds that govern the amount of energy, momentum and mass of a steady Rev. Ceres, Viçosa, v. 68, n.1, p. 010-022, jan/feb, 2021 state fluid in a two-dimensional field can be described as follows: (1) Where is the air density (kg m -3 ), is the air pressure (Pa), is the force of gravity (m s -2 ), is the thermal conductivity of the air (W m -1 K -1 ), y are the components of speed (m s -1 ) y is Reynolds' stress tensor (N).
The turbulence was incorporated into the numerical model through the standard empirical k -ε model, a model based on the transport equations that solve kinetic turbulent energy k and the dispair of this energy per unit volume ε. This model has been the most used and widely validated in greenhouse airflow studies showing to be efficient with the use of computational resources and providing realistic solutions (Mesmoudi et al., 2017). The transport equations for k and ε can be modeled as: Where is viscosity and is turbulent viscosity (kg m -1 s -1 ), G b is turbulent kinetic energy generation due to buoyancy, G k is the generation of turbulent kinetic energy due to speed gradients, σ K y σ E are Prandtl's turbulent numbers for k and ε, v is the coefficient of kinematic viscosity, Y M is the fluctuating expansion in turbulence due to the overall dissipation rate and C 1E , C 2E , C µ , σ K y σ E are constant with empirically determined values and set by default in the simulation software. Buoyancy forces due to differences in air density are added as a source term to the impulse equation through the Boussinesq approach.  The presence of anti-insect screens was modelled by using equations derived from the flow of a free and forced fluid through porous materials, taking into account their main characteristics of porosity and permeability (Valera et al., 2005). These equations can be derived from the Forchheimer equation: Where u is the speed of air (m s -1 ); µ is the dynamic viscosity of the fluid (kg m -1 s -1 ), K is the permeability of the medium (m 2 ); Cf is the inertial factor of the porous mesh; ñ is the air density (kg m -3 ) and ϑx the thickness of the porous material (m). In our research we use a twodimensional modeling approach, this is valid and generally reliable results are obtained for cases where there are dominant wind directions and in regions where the greenhouses establish their orientation under this concept (Villagrán et al., 2019).

Construction of the computer domain and meshing process
This type of research involving transport phenomena, where the final objective is to analyze the phenomenon of natural ventilation in the middle section of a greenhouse structure equipped with side openings. It is possible to implement a two-dimensional simulation model that allows to obtain results in a shorter period of time than in a threedimensional simulation model (Baxevanou et al., 2010;Benni et al., 2016). The 2D approach also allows to obtain a saving of the computational resources that are consumed in the process of development of the model, meshing and numerical simulation (Mesmoudi et al., 2017). Finally, through the two-dimensional models, precise and reliable solutions can be obtained when the fluid flow is perpendicular to the ventilation sections under study (Bartzanas et al., 2004: Villagrán & Bojacá, 2019b. For pre-processing stage of the CFD-RANS simulation, it is recommended to follow guidelines for the construction of the computational domain in engineering problems involving air movement and its interaction with buildings, such as those proposed by Blocken (2015). Where it is established that the entrance and upper limit of the domain must have a minimum separation with the building of 5H, while the exit limit must have a minimum distance of 10H with the building, where H is the maximum height of the building under analysis. In our case, this computational domain has dimensions of 70 m (y-axis) and 243 m (x-axis) (Figure 3).
The meshing process must be of high quality in order to ensure high quality results, independent of mesh size and without convergence problems (Tominaga et al., 2008). The computational domain was divided into an unstructured mesh of 520,183 and 561,235 square elements for PG and MPG respectively (Figure 4).
This mesh resolution was defined after developing a sensitivity analysis of 8 squares with different sizes, which allowed to verify the independence of the numerical solution to the mesh size ( Figure 5). The quality of the mesh was evaluated by means of 2X2 determinant and the orthogonality criterion, where for the first one values were obtained for 99.8% of the cells between 0.95 and 1 and for the second criterion a minimum value of 0.94 results that are considered of high quality (Reynafarje et al., 2020).

Numerical simulation and boundary conditions
In the upper part of the computer domain a limit condition of symmetry was established, the floor, the air intake towers and the walls of the greenhouse the limit condition established was that of an enhanced wall. The ventilation areas of the greenhouse with porous meshes for protection against insects were treated as porous media, the air inlets and outlets in the intake towers were treated with the limit condition of the interior, the flow output of the computational domain was set to the pressure output condition.
While the air intake side was set by a user-defined function (UDF) to a logarithmic pattern profile for air velocity, where vertical air movement profile parameters and turbulence conditions are set, this profile was generated following the procedure established by Hargreaves & Wright (2007). For each material considered in the computational domain such as; air, polyethylene, soil, porous screen, physical, optical and aerodynamic properties were established in the case of the porous screen, the values of these properties were previously established by Flores-Velazquez et al. (2013) and Villagrán et al. (2019; (Table 1).
As the objective of the work is to determine the effect on the ventilation rates and the thermal behavior inside the greenhouse due to the possible implementation of air capture towers, the numerical model does not consider either the presence of plants, nor was the phenomenon of solar radiation simulated with any of the available simulation models. For this last case and as a valid simplification a condition of heat flow dependent on the level of solar radiation was established as it was done in the work recently developed by Flores-Velázquez & Vega-García (2019) and by Villagrán & Bojaca (2019c). The numerical simulations were developed in steady state with a pressure-based, segregated solver and an implicit formulation for the linearization of the discretized equations by second order schemes, the convergence criteria according to the absolute residues for the equations of momentum, turbulence, mass and energy were established in 10 -6 (Bouhoun Ali et al., 2017). Rev. Ceres, Viçosa, v. 68, n.1, p. 010-022, jan/feb, 2021

Validation of the CFD model
To validate the model, an experimental process was developed that included the recording of climatic variables in the greenhouse's external environment, such as temperature, air humidity, wind speed and direction, and solar radiation, by means of a weather station Davis-Vantage 2 plus 6162 (Davis Instruments, Hayward, CA, EE. UU.). Inside the greenhouse on the cross section (xaxis), at a height of 1.5 m above ground level (y-axis), 5 micro weather stations of the series WatchDog 1000 (Spectrum Technologies, Aurora, IL), which were used to record and store the temperature values during the days between 1 November and 30 December 2019, the frequency of measurement and recording was 10 minutes both inside and outside the greenhouse.
Regarding the variable recorded inside the greenhouse for the validation process, it should be mentioned that it was defined from the available monitoring equipment, establishing temperature prediction as a parameter for evaluating the simulation efficiency of the CFD model. Additionally, for this case where no crops are present the results obtained from the CFD are a combination of the relationship between flow distribution patterns and temperature distribution. Therefore, temperature is an indirect measurement variable that allows us to have an idea of how close the   Rev. Ceres, Viçosa, v. 68, n.1, p. 010-022, jan/feb, 2021 numerical solution of the CFD model is to the real aerodynamic behavior of the structure (Villagrán & Bojacá, 2019a).
The data obtained during the measurement period were analyzed, selecting those of the hour 12. Which presented average values for the outdoor environment; 33.6 ± 0.67°C for the outdoor temperature, 695.9 ± 23.5 W m -2 for solar radiation and 1.56 ± 0.26 m s -1 for wind speed with a predominant direction to the west (W). This information was used to carry out a simulation establishing these parameters as initial conditions.
Once the simulation was completed, the post-process was carried out by extracting the data on air temperature inside the greenhouse at the x-points, which coincided with the internal sampling points. The data obtained through simulation and experimentation were compared through the calculation of the mean absolute error (MAE) and the root mean square error (RMSE). (7) where Tmj and Tsj, are the value of the measured and simulated temperature respectively and n is the number of data compared.

Scenarios considered for simulation
Once the numerical model was validated and its adequate prediction of natural ventilation and thermal behavior in the PG greenhouse was verified, this model was used to develop the simulations for two specific scenarios, the PG and MPG greenhouses, evaluated with 4 outdoor wind speeds ( Table 2). The initial conditions corresponded to the maximum values recorded in the study area for both temperature and solar radiation of 34.8 °C and 860 W m -2 , respectively, and were kept constant in each simulated scenario.

Validation of the CFD model
The calculation of the evaluated goodness-of-fit parameters for the simulated and measured temperature at the observation points inside PG showed values of 0.31 and 0.34 °C for MAE and RMSE respectively. Likewise, the qualitative behavior observed when graphing the spatial trend of the behavior of the simulated and measured data in the transversal axis of the greenhouse (x-axis), allows confirming that the CFD model has a high capacity to satisfactorily predict the trend of the thermal behavior inside PG, therefore, its use for the objective of this investigation is valid (Figure 6).
These results, where the difference between experimental and simulated data is less than 1 °C, are similar to those found in several CFD studies of natural ventilation in passive greenhouses such as those developed by He et al. (2015He et al. ( , 2018 and Villagrán & Bojacá (2019b).

Air flow patterns
The airflow patterns inside PG and MPG will be influenced by the natural convection inside the greenhouse and the effect generated by the wind and the wind speed from the outside environment. The airflow patterns simulated for PG show two characteristic behaviors, the first one occurs for the S1 scenario which can be considered as a low outdoor wind speed scenario, under these conditions it is possible to observe an airflow entering the greenhouse through the lower windward side area and inside PG it generates two flow movements, one of horizontal type and parallel to the airflow movement outside, this flow is shown in direction and exit through the lower leeward side, the other flow is of convective type appears over the upper central part and between this region and the cover of PG, this is generated by the movement of warm, lower density air that occurs floating from the floor of the greenhouse to the roof region (López et al., 2011), Some of this air flow exits the structure through the roof ventilation areas (Figure 7). For PG scenarios S2, S3 and S4 the same flow pattern is observed with different speed due to the speed of the outside air in the simulated case, air movement for PG under this scenario shows inflow through the windward side and roof ventilation areas to the leeward side ventilation areas where they eventually exit to the outside (Figure 7).
For the simulated MPG scenarios, the effects of the air collection towers on the flow pattern inside the structure are observed (Figure 7). It identifies a pattern of flow that enters the greenhouse through the ventilation tower, this airflow is accelerated and directed towards the ground and then enter MPG and move horizontally over the area where the crops are established, it is also observed that as air flow inside MPG suffers a slowdown, these results are consistent with those reported in the study developed por Li & Mak (2007) and Pakari & Ghani. (2019b).
In the case of low velocity it is observed that the flow pattern leaves the greenhouse through the leeward side ventilation tower and the ventilation areas located in the roof area of the greenhouse, on the contrary for the remaining scenarios it is observed that the airflow entering  from the windward side collides with an airflow entering from the leeward tower causing two convective cells of different size and exiting the greenhouse through the roof ventilation areas (Figure 7). These airflow patterns differ from those reported by Pakari & Ghani. (2019b). This can be generated by the buoyancy effects and the presence of anti-insect nets considered in our study. The quantitative behavior of the air movement was made through the calculation of the mean speed (MV) and the normalized speed (NV) for two heights above ground level 1 and 2 meters, this parameter is the relation of the speed of the internal air flow with the speed of the external wind (Table 3). MV values for PG ranged from 0.27 ± 0.05 to 0.54 ± 0.22 m s -1 for the 1 m height, while for the 2 m height these values were 0.22 ± 0.06 and 0.66 ± 0.35 m s -1 characteristic velocities in sideventilated and roof-ventilated passive greenhouses (Teitel et al., 2006).
On the other hand, for MPG the values of MV were 0.34 ± 0.09 and 1.31 ± 0.39 m s -1 for the height of 1 meter, which represents an increase of 25.9 and 142.5% compared to PG, for the 2 meter height the MV values were 0.13 ± 0.09 and 0.78 ± 0.39 m s -1 which compared to PG are 40.4% lower for S1 and 18.1% higher for S4.
In general it can be mentioned that for MPG the indoor air flows have a higher velocity compared to PG and this is demonstrated by the VN value obtained in the cross section of the greenhouse at each of the two heights evaluated (Figure 8). Additionally, it is possible that the VN value in MPG presents a greater variation along the transversal axis and between heights analyzed, this is generated by the convective type movements that are presented in this MPG (Figure 7).
The VN values show a reduction in airflow velocity compared to outdoor velocity of 82.3 and 79.6% for the most critical case in PG and MPG respectively, while the lowest reduction was 53.9 and 44.2% for PG and MPG, this difference in speed between the exterior and interior flow of the greenhouses is caused by the presence of the insect protection meshes, which generate a restriction to the movement of air towards the interior of the greenhouse, the values obtained in PG are similar to those reported by Baeza et al. (2009) in a natural ventilation study developed in a greenhouse equipped with anti-insect screens.

Indexes of hourly renewal
The air exchange rate of a naturally ventilated structure is a parameter that allows measuring the efficiency of the ventilation system, this rate is defined as the relationship between the volumetric rate of air flow exiting the structure and the volume of the structure, this value was calculated by integration method mass flow rate (Senhaji et al., 2019;Villagrán et al., 2019).
The results found allow us to verify that the air renewal in the MPG greenhouse has a higher value than that of PG   ( Figure 9). The values for PG ranged from 20.1 to 74.9 vol h -1 , while those for MPG were 31.7 and 107.5 for S1 and S4 respectively, The renewal rate in MPG is 57.3 and 43.4% higher than PG for these two cases. In the case of MPG it can be said that for speeds e" at 1 ms-1 the renewal rates are above the minimum value recommended for naturally ventilated greenhouses which is 45 vol h -1 , value that allows to regulate the temperature and relative humidity conditions inside the greenhouse (ASHRAE, 2009;Villagrán et al., 2012).

Thermal behaviour
The spatial distribution patterns inside PG and MPG, allowed to observe the differentiated behavior that exists between greenhouses, which allowed to verify the benefits obtained from the implementation of the ventilation towers. In MPG, by obtaining higher renewal rates and better airflow movement, an additional thermal distribution of lower magnitude and more homogeneous compared to PG is obtained (Figure 10). This is a relevant result since sometimes a higher ventilation rate does not necessarily generate a reduction of the temperature value inside the structure and a greater homogeneity, this happens when the optimization strategies of natural ventilation do not generate adequate air flows inside the greenhouse, especially in the area where the crop will be established (Sase, 2006;Senhaji et al., 2019). The average values of temperature (MT) and thermal differential between outside and inside (∆T m ) were calculated for each of the simulated scenarios and for the two heights analyzed (Table 4). The MT value obtained for PG ranges from 36.3 ± 0.6 to 37.1 ± 1.1 °C, giving values of ∆T m of 1.5 ± 0.6 and 2.3 ± 1.1 °C respectively. On the other hand, for MPG the MV value varies between 35.0 ± 0.2 and 35.9 ± 0.3 which gives values of ∆T m 0.2 ± 0.2 and 1.0 ± 1.5 °C.
Likewise, the standard deviation values of the parameters evaluated show us that the spatial behavior of the temperature in this greenhouse is more homogeneous, a factor that is undoubtedly key in the management of the microclimate since this variable intervenes in the physiological and biological processes of the plants and any heterogeneous distribution will directly influence the uniformity in the growth of the crops, which affects the quality of the final product obtained (Ma et al., 2019 a,b).
Additionally, to analyze the spatial distribution of the thermal differential (∆T), it was plotted for all simulated scenarios, in a cross-sectional profile of the greenhouse (Figure 11). For MPG it can be seen that the distribution of ∆T is very uniform in the profiles analyzed, in the case of the scenario S1 it is observed as a large part of the transversal profile (0-15 m) presents a ∆T of approximately 1.0 °C, while for the scenarios S2, S3 and S4 the value of ∆T shows values lower than 0.8 °C in a large part of the transversal axis of the greenhouse with the exception of a region between 12.5 and 15 m, this region is the one influenced by the two convective cells analyzed in the air flow patterns (Figure 7).
On the contrary, for scenario S1 of PG the greatest heterogeneity found in this study is observed where the  Figure  11). The greater homogeneity of the temperature in MPG can be explained by the convective-type airflow patterns, which usually generate a better air mixture that if complemented by adequate ventilation areas generate homogeneous thermal conditions (Mesmoudi et al., 2017).

CONCLUSIONS
The implementation of ventilation towers proved to be an interesting alternative for the optimization of the natural ventilation of passive greenhouses. This microclimatic optimization alternative can be applied to the various models of passive greenhouses used in tropical countries.
Air flow patterns in Colombian greenhouses equipped with side air capture towers showed higher velocity and convective type movement which together generate ventilation rates of up to 57.3% compared to the standard greenhouse. The greenhouse with ventilation towers showed a highly homogeneous thermal behavior with thermal differentials between inside and outside below 1.1 °C.
Under the methodological approach presented in this work, future research can be proposed to evaluate the behavior of natural ventilation in this same greenhouse prototype with the presence of crops of commercial interest such as; tomatoes, cucumbers and paprika. To study the combined effect of natural ventilation together with evaporative cooling strategies. This methodology can be applied to other types of established greenhouses in various climatic conditions. Finally, it would also be interesting to reinforce the validation method by means of a sonic anemometry or hot-wire analysis, or in its absence by a study in the wind tunnel.