A Volterra model of the greenhouse temperature using natural ventilation

A Volterra model of the greenhouse temperature using natural ventilation J.K. Gruber ∗,F.Rodr´ ıguez ∗∗,C. Bord´ ons ∗,J.L. Guzm´ an ∗∗, M. Berenguel ∗∗,E. F.Camacho ∗ ∗Dep. de Ingenier´ ıa de Sistemas yAutom´ atica, Universidad de Sevilla, C. de los Descubrimientos s/n 41092 Sevilla, Spain. E-mail: [jgruber,bordons,eduardo]@cartuja.us.es ∗∗Dep. de Lenguajes yComputaci´ on, Universidad de Almer´ ıa, Ctra. Sacramento s/n 04120 Almer´ ıa, Spain. E-mail:[frrodrig,joguzman,beren]@ual.es Abstract: Model-based control techniques are commonly applied to control the greenhouse climate. As well-known, these techniques require accurate models for adequate results. Several first-principle models havebeen developed for the greenhouse climate control problem considering all the physical and physiological processes. However,these models are too complexto beused for control purposes. On the other hand, empirical models based on input/output real data allowtoobtain better results and less complexmodel structures. In the Mediterranean areas the main problem is cooling the greenhouse and this leads to natural ventilation as astandard tool. This paper presents the development and the results of aVolterra model for the greenhouse temperature including the crop effect and using natural ventilation. ©Copyright IFAC2008. 1. INTRODUCTION Crop growth requires greenhouse climate conditions achievable by means of advanced control techniques based on reliable models of the greenhouse climate. In the southeast (SE) of Spain, one of the main problems is greenhouse cooling by natural ventilation. During recent years, advanced control techniques havebeen developed, such as adaptivecontrol (Sigrimis et al., 1999), feedforward control (Rodr´ ıguez et al., 2001), optimal control (van Straten et al., 2002), feedback linearizing control based on physical (Pasgianos et al., 2003) and empirical (Berenguel et al., 2006) models. This list of references is evidently limited and manyimportant papers on the temperature control (Tantau, 1993; Seginer and Sher, 1993; Martin-Clouaire et al., 1996; Fuchs et al., 1997), etc. are not mentioned due to space constraints. Twowidely differing strategies for modelling are those based on either first principles or input/output experimental data (Pearson, 1995). These techniques use a model of the greenhouse climate which links the output variable (temperature) to the control variable (vent position) and to the disturbances (mainly outside weather and crop). There exists a correlation between vent position and output temperature that can be ascribed to the ventilation rate. As acontributing factor, other researchers haveidentified this relationship by constructing semi-physical models, such as (Boulard and Baille, 1995) and (Kittas et al., 1997), and input/output empirical models (Rodr´ ıguez et al., 1999; Parra et al., 2004, 2006), and even by using an on-line estimator of the ventilation rate based on an unknown input observer with an output linearising feedback, as done by (Bontsema et al., 2006). The authors havepreviously used both physical and empirical models. In the first case, a greenhouse climate model based on energy balances wasdeveloped. Fundamentally,it isadynamic model represented by asystem of differential equations dX/dt =f(X,U,P,V,C,t), with X(ti)=Xi,where X(t),U(t),P(t),V(t)and Care respectively vectors of state variables, input variables, disturbances, system variables and system constants, tis the time, Xiis the known initial state at the initial time ti.Hence, f=f(t)is a nonlinear function based on mass and heat transfer balances. The number of equations describing the system characteristics depends on the greenhouse elements, the installed control actuators and the type of cultivation method. The physical processes included in the energy and mass balances are solar and thermal radiation absorption, heat convection and conduction, crop transpiration, condensation and evaporation (Rodr´ ıguez, 2002). In order to model the natural ventilation rate, φ ven,ofagreenhouse equipped with sidewall and roof continuous actuators, some proposals from other researchers such as (Boulard and Baille, 1995) and (Kittas et al., 1997) havebeen tested. These models suggest that the ventilation process is induced by the stack effect due to thermal buoyancy(produced by the difference between internal and external air temperature, Xt,aand P t,erespectively) and the wind effect due to the wind speed (Pws,e)action on the greenhouse structure. The best results were obtained with the following proposal of (Kittas et al., 1997) φ ven =cven,d2cgcven,h Xt,a−P t,e ¯ Xt · Vven,ar−latVven,ar−rf V2 ven,ar−lat+V2 ven,ar−rf   2 (1) +Vven,ar−lat+Vven,ar−rf 22 cven,wP2 ws,e+ φ losses where cven,dis the discharge coefficient of the ventilators, cven,wis the wind pressure coefficient, cgis the acceleration of gravity,¯ Xtis the mean absolute temperature, cven,his the vertical distance between the midpoints of the lateral and roof vents. Proceedings of the 17th World Congress The International Federation of Automatic Control Seoul, Korea, July 6-11, 2008 978-3-902661-00-5/08/$20.00 © 2008 IFAC 2925 10.3182/20080706-5-KR-1001.1866 Vven,ar−latand Vven,ar−rf are the areas of the sidewall and roof ventilation openings, givenby: Vven,ar−lat=cven,long−latcven,wid−latUven (2) Vven,ar−rf =2cven,long−rfcven,wid−rfsin(Uven/2) where cven,long−latand cven,long−rf are the lengths of the lateral and roof ventilation, cven,wid−latand cven,wid−rf are the widths of the lateral and roof ventilation, and Uven is the control signal of the ventilation openings in [%]. The leakage term, φ losses, represents loss through leaks when the vents are closed, approximated by an imposed threshold between high and lowwind speed. It wasclear that the direct application of this kind of model to model-based control strategies entailed severe complications for techniques that required model inversion in order to obtain vent aperture (Uven)from ventilation rate estimates (Rodr´ ıguez et al., 2001). Thus, empirical ventilation rate models were developed, both artificial neural-network based models (Rodr´ ıguez et al., 1999) and static empirical ones (Parra et al., 2006) to try to avoid the mentioned problems. In this last case, some studies were conducted in the greenhouse used in this work, where the well established logarithmic transformations of the mass balance equation for the decay-rate of atracer gas(N2O) were employed by (Baptista et al., 1999), although that method cannot be used during crop production. One of the main results obtained by applying regression techniques to actual measurements on the effect of combined roof and sidewall ventilator actions wasthe fact that the expected ventilation rate per unit ventilator area and wind speed, followed an exponential expression of the combined aperture set point, φ (Uven)=Pws,e α U β ven + φ losses,where Uven is the percentage aperture of the ventilators. The quantities α and β are tuning parameters which, according to actual measurements, show subtle variations between leeward and windward ventilation. φ losses is the calculated inherent ventilation leakage due to greenhouse aerodynamics and layout. Further studies were carried out in order to achievean approximate, butwell-tuned practical model, suitable for automatic control. Anyway,as pointed out by (Bontsema et al., 2006), since in practice the greenhouse climate is neverinsteady state, due to the external disturbances, the error is large due to neglecting the dynamic storage term. In all the previous models, it can be seen howthe models of the natural ventilation rate are quite difficult to obtain, so that, an approach based on Volterra models has been developed in order to find atrade-offbetween the difficulty for obtaining the model and the potential use for control purposes, applied for diurnal temperature control and applied to the typical parral greenhouse structure of the SE of Spain. Volterra models represent the logical extension of convolution models which havebeen applied successfully in linear model-based predictive control (MPC) (Doyle et al., 2002; Maner et al., 1996). These models exhibit generically agood behaviour and their structure can be exploited in the design of controllers, especially in the case of second order models. Asecond order Volterra model, with the truncation of terms (truncation orders N1and N2), can be defined as: y(k)=h0+ N1 ∑ i=1 a(i)u(k−i) + N2 ∑ i=1 N2 ∑ j=1 b(i,j)u(k−i)u(k−j)(3) which corresponds to alinear convolution model with anonlinearity as additional and additiveterm, as described in section 2. In that model y(k)and u(k)represent the last measured output and input to the system, respectively (kis the actual sampling instant). Notice that Volterra models havebeen frequently used to model bilinear systems in such away that it seems to be agood idea to use this formulation for modelling greenhouse temperature dynamics, including the disturbances in the nominal formulation of second order Volterra models. In aprevious work, apreliminary Volterra model wasdeveloped to model the inside temperature of an empty greenhouse (without crop) in order to evaluate the behavior of this modelling technique for this kind of systems (Gruber et al., 2007). Several PRMS (Pseudo Random MultilevelSequence) tests were performed using natural ventilation to obtain adequate data for identification purposes, because typical PRBS (Pseudo Random Binary Sequence) and RBS tests do not sufficiently excite nonlinear systems (Braun et al., 2001). The resulting model adequately fitted the real data butthe number of parameters was excessively high. Furthermore, the crop has an important effect on the greenhouse temperature and thus it is akeyfactor to be included in the system model. Therefore, the present work faces the previous drawbacks reducing the number of model parameters and considering the crop effects. The paper is organized as follows. It begins with materials and methods used where both equipment and the identification and modelling procedures are explained. Next, representative results are outlined. Finally,the conclusions section will summarize the achievements and propose intended applications. 2. MATERIALS AND METHODS 2.1 Greenhouse The data were acquired in aParral greenhouse located in Las Palmerillas Experimental Station (El Ejido, Almer´ ıa, Spain, Figure 1). It is asymmetric curved flat roof greenhouse with five North-South oriented spans of 7.56×23.2 m(877 m2total soil surface) and height from 2.8 mto4.4 m. The covering material is a200-micron thick PE film, laid on agalvanized steel structure 1.The control actuators installed in the greenhouse are automated flap roof and rolling lateral continuous ventilators, the former with amaximum aperture angle of 45oand the latter with 37 mlength and 1.2 maperture. Both haveinsect screens installed in the ventilation openings. Soil temperature measurements were made using semiconductor sensors at different depths, i.e., on both side of the mulch, immediately under the soil surface layer and at adepth of 50 mm. The greenhouse air temperature thermoresistance sensor and the air relative humidity capacitivesensor were placed at the top of the crop. A meteorological station wasinstalled outside at aheight of 6m for measurements of temperature, relativehumidity,solar and photosynthetically activeradiations, rain, and wind speed and direction. Some additional crop data were logged, i.e., substrate water content, water supply,amount of drainage, electrical conductivity of the substrate, plant transpiration measured with an electronic weight, leaf temperature and CO2concentration. A uniform sampling time of one minute wasestablished. 1http://aer.ual.es/CJPROS/engindex.php?Opcion=5 17th IFAC World Congress (IFAC'08) Seoul, Korea, July 6-11, 2008 2926 Fig. 1. Greenhouse used for the experiences. 2.2 Models, data and parameter identification As the main objectiveof the work presented in this paper was to evaluate the potential use of Volterra models to account for the main dynamics describing changes in inside air temperature to outside weather using only natural ventilation, the data used for identification purposes wascollected with crop (tomatoes) inside the greenhouse (April/May,2006). The influence of the crop wastaken into account as adisturbance to the greenhouse temperature. The mentioned disturbance wasmodeled by means of the the Leaf Area Index(LAI). Thus, the main variables considered for modelling purposes are: •Output: Xt,ainside air temperature. •Input: alat,rf aperture of the roof and lateral ventilations. •Disturbances: P t,eoutside temperature, Pws,eoutside wind speed, P t,ss soil surface temperature, Psol,eoutside global solar radiation, PLAI Leaf Area Index. During the identification process, additional inputs were tested (not included in this document), as the product alat,rfPws,e,but no large improvements in the parameter identification stage was detected while the number of parameters did increase. Moreover,the selection of these modulated variables will prevent the use of model predictivecontrollers. Asecond order truncated Volterra models that represented in equation (3). In the case of autoregressive(AR) models, the structure is: y(k)=h0+ Ny ∑ i=1 h1(i)y(k−i)+ N1 ∑ i=1 a(i)u(k−i) + N2 ∑ i=1 N2 ∑ j=i b(i,j)u(k−i)u(k−j)(4) where h0is the offset and h1(·),a(·)and b(·),the AR parameters, the linear parameters and the nonlinear second order parameters, respectively.Second order AR terms havenot been included, as the selection allows obtaining adequate results. In the theoretical case the upper limit of the convolution sum will be ∞,but in astable system as the greenhouse, the number of terms can be truncated (sum limits Ny,N1and N2)toreduce the number of parameters to identify,insuch away that h1(i)= 0,∀i>Ny,a(i)=0,∀i>N1,b(i,j)=0,∀i>N2,j>N2.Diagonal models are often used to reduce the number of parameters to estimate (b(i,j)=0,∀i=j), but inthe case treated in this paper only non-diagonal second order Volterra models havebeen used. Notice that, although the disturbances could be distinguished from the manipulable variables (for instance, when using these models for control purposes), there is no difference between them for identification purposes, and then, both manipulable variables and disturbances will be denoted by u(·)in this work without loss of generality. As the system under consideration is MISO (Multiple Input - Single Output), the model in the AR case (in the non-AR case h1(i)=0,∀i)isgiven by: y(k)=h0+ Ny ∑ i=1 h1(i)y(k−i)+ n ∑ l=1 N1,l ∑ i=1 al(i)ul(k−i) + n ∑ l=1 N2,l ∑ i=1 N2,l ∑ j=i bl(i,j)ul(k−i)ul(k−j)(5) In order to identify the parameters of the Volterra models, data in the period from 19 April to 02 May,2006 havebeen used (Figure 2, only aperiod of 5days is shown, time scale in days), while data from 03 May to 07 May,2006 havebeen used for validation purposes (Figure 3). Obviously,more data should be necessary to obtain more reliable results, butthe selected data are enough for assessing the validity of the technique for greenhouse temperature modelling. The sample time was 2minutes and the data were adequately filtered through first order filters before using them for calibration purposes. Before performing the identification and validation, the filtered data are also normalized, 17th IFAC World Congress (IFAC'08) Seoul, Korea, July 6-11, 2008 2927 truncation orders   Name Method Order alat,rf P t,ePws,ePsol,eP t,ss PLAI Xt,a∑parameters ε id ε val lb,ub A,b 0lsq 1st 30 1 10 5 10 1 0 234 0.6024 0.6994 =[],=[]=[],=[] 2nd 15 1 10 0 0 0 0 Afmin 1st 30 1 10 5 10 1 0 234 0.7710 0.9450 =[],=[]=[],=[] 2nd 15 1 10 0 0 0 0 Bfmin 1st 20 1 2 5 10 1 0 99 0.8297 1.1358 =[],=[]=[],=[] 2nd 10 1 2 0 0 0 0 Efmin 1st 20 1 2 5 10 1 0 99 0.8909 1.1680 =[],=[]=[],=[] 2nd 10 1 2 0 0 0 0 Car 1st 1111111 13 0.5976 0.7357 =[],=[]=[],=[] 2nd 1111100 Dar 1st 2222222 30 0.5961 0.7264 =[],=[]=[],=[] 2nd 2222220 Table 1. Main numerical results 0 1 2 3 4 5 0 50 100 t [d] alat,rf [%] 0 1 2 3 4 5 15 20 25 30 t [d] Pt,e [ºC] 0 1 2 3 4 5 0 5 10 15 t [d] Pws,e [m/s] 0 1 2 3 4 5 0 500 1000 t [d] Psol,e [W/m2] 0 1 2 3 4 5 20 22 24 26 28 t [d] Pt,ss [ºC] 0 1 2 3 4 5 1.6 1.8 2 2.2 t [d] PLAI [ºC] 00.5 11.5 22.5 33.5 44.5 5 15 20 25 30 35 t [d] Xt,a [ºC] Fig. 2. Rawdata for identification. 0 1 2 3 4 5 0 50 100 t [d] alat,rf [%] 0 1 2 3 4 5 15 20 25 30 t [d] Pt,e [ºC] 0 1 2 3 4 5 0 5 10 15 t [d] Pws,e [m/s] 0 1 2 3 4 5 0 500 1000 t [d] Psol,e [W/m2] 0 1 2 3 4 5 20 22 24 26 28 t [d] Pt,ss [ºC] 0 1 2 3 4 5 2.4 2.6 2.8 t [d] PLAI [ºC] 00.5 11.5 22.5 33.5 44.5 5 15 20 25 30 35 t [d] Xt,a [ºC] Fig. 3. Rawdata for validation. Vnor id =Vid −Vmin id Vmax id −Vmin id ;Vnor val =Vval −Vmin id Vmax id −Vmin id where Vis each one of the used input/disturbance/output signals both in the identification (id)and validation (val)data sets, supra-indexmin means the minimum value in that set and max the maximum one. In this way,the identification data are normalized to the interval [0,1] and the validation data to an interval close to that. With these normalized values, the output of the Volterra model is also normalized (ynor), in such away that the real output is given by y=ynor(Xmax t,a−id −Xmin t,a−id)+ Xmin t,a−id.Toidentify the parameters of the Volterra models, both least squares methods and constrained nonlinear optimization using sequential quadratic programming havebeen used. 3. REPRESENTATIVE RESULTS In this section several representativeresults are shown, both in tabular and graphical ways. Both non-AR and AR models havebeen identified and validated. The meaning of the abbreviations used in Table 1are: lsq (least squares method), fmin (constrained nonlinear optimization), ar (AR model), lb, ub (lower and upper bounds of the parameters) A,b(matrices defining linear constraints of the form Ax ≤B), ε id, ε val (mean square error for identification and validation: ε =∑n i=1(Xt,a(i)−ˆy(i))2 n where ˆyis the output of the model. Table 1shows the main results of the identification and validation stages. 3.1 Non-AR models Identification test 0:the first identification wasperformed using the least squares method. One of the initial problems was to select the truncation orders for the six variables influencing the greenhouse inside temperature. Due to this reason, in a first step identifications were performed in such away that for each variable along truncation order wasselected, while the truncation order for the rest waszero. In these identifications the truncation orders of the second order terms were also zero. In asecond step, six newidentifications were carried out to determine the second order truncation orders, using the truncation orders for the linear terms selected in the first step. Finally,once the truncation orders were determined, the identification using the least-squares method wasdone. The used truncation orders are shown in Table 1, as well as the number of parameters (234, including the offset) and the mean square errors. Identification test A:Asthe plots of the parameters obtained in the previous case did not showaconsistent convergence (many parameters havenoisy profiles), constraints havebeen included in the identification process, then using nonlinear optimization techniques. The number of parameters is the same than in the 17th IFAC World Congress (IFAC'08) Seoul, Korea, July 6-11, 2008 2928 Order Const. alat,rf P t,ePws,ePsol,eP t,ss PLAI 1st lb ––000– ub 0––––– 2nd lb 0––––– ub ––0––– Table 2. Constraints for identification A previous case (see Table 1) butwith the constraints shown in Table 2. The results of the identification and validation are shown in Table 1. With the constraints, the linear parameters havemore reasonable values, and the parameters corresponding to the second order terms that are not zero are close to the main diagonal (not shown for the sakeofspace). Identification test B:Inthe previous experience, it wasdetected that several of the linear terms parameters (specially those related to Pws,eand alat,rf)were equal to zero, in such away that the corresponding truncation orders could be decreased to 20 and 2, respectively.Inthe same way,the truncation order of the second order terms of Pws,ewere reduced to 2, providing atotal number of 99 parameters. The same constraints were used and the mean square errors were slightly higher than in test A. Identification test E:Inthe three previous cases, the linear parameters associated to variable alat,rf did not achievecoherent values, in such away that additional conditions were added. The same truncation orders and constraints than in the B case were used, but it was imposed that aalat,rf (k)≥aalat,rf (k− 1)fork=2,...,20, this conditions being quite easily implementable in the nonlinear constrained optimization problem. In an analogous waythe conditions aPsol,e(k)≤aPsol,e(k−1)fork= 2,...,5and aP t,ss (k)≤aP t,ss (k−1)fork=2,...,10 for the outside global solar radiation and the soil surface temperature were used. Aplot of the results of the identification and validation is shown in Figure 4. The resulting linear parameters are shown in Figure 5, while the second order term parameters are shown in Figure 6. It can be seen that the linear terms corresponding to the variables Psol,e,P t,ss and alat,rf are smooth curves. A further reduction of the truncation order for the corresponding variables is not possible anymore. Furthermore, with the imposition of the additional conditions for the linear parameters, the nonlinear parameters of the model input alat,rf showacoherent and reasonable shape. Diagonal models were also tested, but the mean square errors increased. 3.2 AR models Twoidentification tests are shown using AR models with nonlinear open loop identification. Identification test C:Inafirst step, all the truncation orders were set to 1for the linear parameters of the six variables influencing the greenhouse temperature. Forthe disturbance PLAI no second order terms were used, for the remaining five variables anonlinear truncation order of 1was chosen. With regards to the AR terms, atruncation order equal to 1was used for the linear part and 0for the nonlinear one. The offset has been also added, with atotal of 13 parameters (Table 1). Identification test D:Finally,asecond identification is shown where the truncation orders for the linear and nonlinear parameters related to the input and disturbance variables are set to 2, except for the second order terms of PLAI which were not used as in the previous identification. The AR part of the model has 00.5 11.5 22.5 33.5 44.5 5 15 20 25 30 35 t [d] Xt,a, y [ºC] real model 00.5 11.5 22.5 33.5 44.5 5 15 20 25 30 35 t [d] Xt,a, y [ºC] real model Fig. 4. Test E:identification (up) and validation (down) results 0 1 2 3 4 5 0 0.05 0.1 0.15 0.2 0.25 k [−] aPsol,e (k) 0 2 4 6 8 1 0 0 1 2 3 4x 10−3 k [−] aPt,ss (k) 0 2 4 6 8 10 12 14 16 18 2 0 −0.015 −0.01 −0.005 0 k [−] aalat,rf (k) Fig. 5. Test E:some linear parameters 0 2 4 6 8 10 1 2 3 4 5 6 7 8 9 10 0 1 2 3 4 5 6 7 8 x 10−3 k [−] j [−] balat,rf (j,k) Fig. 6. Test E:some second-order parameters atruncation order with 2parameters in the linear part and 0 parameters in the second-order terms (nonlinear AR terms have been neglected). Again, the offset is included, providing atotal 17th IFAC World Congress (IFAC'08) Seoul, Korea, July 6-11, 2008 2929 of 30 parameters. The main results are shown in Table 1and in Figure 7(plot of the results of the identification and calibration with the AR model). It can be seen that the mean square error is reduced to ε id =0.5961 and ε ver =0.7264. 00.5 11.5 22.5 33.5 44.5 5 15 20 25 30 35 t [d] Xt,a, y [ºC] real model 00.5 11.5 22.5 33.5 44.5 5 15 20 25 30 35 t [d] Xt,a, y [ºC] real model Fig. 7. Test D:identification (up) and validation (down) results 4. CONCLUSIONS This paper has presented results on the application of Volterra models for the modelling of greenhouse temperature evolution using natural ventilation and including the crop effect. The main interest wastosee howthese models cope with the nonlinear behaviour inherent in the relationship between temperature and vents aperture, through the ventilation rate, which is one of the most difficult dynamics to be modelled in the greenhouse as discussed in the introduction section. Future works will include the application of model predictivecontrol schemes for greenhouse temperature control using the models obtained in this work. ACKNOWLEDGEMENTS This work wassupported by the Spanish CICYT-FEDER grants DPI2004-07444-C04-01/04 and DPI2007-66718-C04- 01/04, by Fundaci´ on Cajamar and by the Consejer´ ıa de Innovaci´ on, Ciencia yEmpresa de la Junta de Andaluc´ ıa. REFERENCES B.J. Baptista, B.J. Bailey,J.M. Randall, and J.F.Meneses. Greenhouse ventilation rate: Theory and measurements with tracer gastechniques. J.Agri. Eng.Res. 72, pp. 363-374, 1999. M. Berenguel, F.Rodr´ ıguez, J.L. Guzm´ an, D. Lacasa, and J. P´ erez Parra. Greenhouse diurnal temperature control with natural ventilation based on empirical models. Acta Horticulturae 179, pp. 205-211,2006. J. Bontsema, E.J. vanHenten, J. Hemming, J. Budding, and T.Rieswijk. On-line estimation of the ventilation rate of greenhouses. ISHS Acta Horticulturae 718,2006. T.Boulard and A. Baille. Modelling of air exchange in a greenhouse equipped with continuous roof vents. J.Agri. Eng.Res. 61, pp. 37-48,1995. M.W.Braun, D.E. Rivera, and A. Stenman. Amodel-on- demand identification methodology for nonlinear process systems. Int. J.of Control, 74(18), pp. 1708-1717,2001. F.J. Doyle, R.K. Pearson, and B.A. Ogunnaike. Identification and Control Using VolterraModels.Springer,London, 2002. M. Fuchs, E. Dayan, D. Shmuel, and Y.Zipori. Effects of ventilation on the energy balance of agreen house with bare soil. Agricultural and Forest Meteorology,86, pp. 273-282, 1997. J. K. Gruber,F.Rodr´ ıguez, C. Bord´ ons, J. L. Guzm´ an, and M. Berenguel. Asurveyof computer-based approaches for greenhouse climate management. IFACAgricontrol,pages 41–46, 2007. C. Kittas, T.Boulard, and G. Papadakis. Natural ventilation of agreenhouse with ridge and side openings: Sensitivity to temperature and wind effects. Trans. of the ASAE. 40(2), pp. 415-425,1997. B.R. Maner,F.J. Doyle, B.A. Ogunnaike, and R.K. Pearson. Nonlinear model predictivecontrol of asimulated multivariable polymerization reactor using second-order volterra models. Automatica, 32(9), pp. 1285-1301,1996. R. Martin-Clouaire, P.J. Schotman, and M. Tchamitchian. A surveyof computer-based approaches for greenhouse climate management. Acta Horticulturae,406:409–423, 1996. J. P´ erez Parra, E. Baeza, J.I. Montero, and B.J. Bailey.Natural ventilation of parral greenhouses. Biosystems Engineering 87(3), pp. 355-366,2004. J. P´ erez Parra, M. Berenguel, F.Rodr´ ıguez, and A. Ram´ ırez- Arias. Ventilation rate models of mediterranean greenhouse for control purposes. Acta Horticulturae 179, pp. 197-204, 2006. G.D. Pasgianos, K.G. Arvanitis, P.Polycarpou, and N. Sigrimis. Anonlinear feedback technique for greenhouse environmental control. Computersand Electronics in Agriculture,40, pp. 153-177,2003. R.K. Pearson. Nonlinear input/output modelling. J.Proc. Control, 5(4), pp. 197-211,1995. F.Rodr´ ıguez, M.R. Arahal, and M. Berenguel. Application of artificial nn for greenhouse climate modeling. European Control Conference,1999. Karlsruhe, Germany. F.Rodr´ ıguez, M. Berenguel, and M.R. Arahal. Feedforward controllers for greenhouse climate control based on physical models. European Control Conference,2001. Porto, Portugal. F.Rodr´ ıguez. Modelling and hierarchical control of greenhouse crop production. PhD thesis, University of Almera, Spain, 2002. I. Seginer and A. Sher.Optimal greenhouse temperature trajectories for amulti-state-variable tomato model,chapter in: The computerized greenhouse, pages 153–172. Academic Press, 1993. N. Sigrimis, P.N. Paraskevopoulos, K.G. Arvanitis, and N. Rerras. Adaptivetemperature control in greenhouses based on multirate-output controllers. Proc. 14th IFACWorld Congress. pp. 571-576,1999. Beijing, China. H. J. Tantau. Optimal Control for plant production in greenhouse,chapter in: The computerized greenhouse, pages 139– 152. Academic Press, 1993. G. vanStraten, L.G. vanWillengenburg, and F.Tap. The significance of crop co-states for receding horizon optimal control of greenhouse climate. Control Eng.Practice,10, pp.625-632,2002. 17th IFAC World Congress (IFAC'08) Seoul, Korea, July 6-11, 2008 2930