MODELO DINÁMICO DEL CLIMA EN UN INVERNADEROPARA FLORES

FREDY ORLANDO RUIZ PALACIOS

PONTIFICIA UNIVERSIDAD JAVERIANAFACULTAD DE INGENIERÍAMAESTRÍA EN INGENIERÍA ELECTRÓNICABOGOTÁ D.C.2005

MODELO DINÁMICO DEL CLIMA EN UN INVERNADEROPARA FLORESFREDY ORLANDO RUIZ PALACIOSTrabajo de investiga ión de maestría

Dire torIng. CARLOS COTRINO M.S .PONTIFICIA UNIVERSIDAD JAVERIANAFACULTAD DE INGENIERÍAMAESTRÍA EN INGENIERÍA ELECTRÓNICABOGOTÁ D.C.2005

PONTIFICIA UNIVERSIDAD JAVERIANAFACULTAD DE INGENIERÍAMAESTRÍA EN INGENIERÍA ELECTRÓNICARECTOR MAGNÍFICO: R.P. GERARDO REMOLINA, S.J.DECANO ACADÉMICO: Ing. FRANCISCO REBOLLEDO M. M.S .DECANO DEL MEDIO UNIVERSITARIO: R.P. ANTONIO SARMIENTO, S.J.DIRECTOR DE MAESTRÍA: Ing. PEDRO VIZCAYA G. Ph.D.DIRECTOR DEL PROYECTO: Ing: CARLOS COTRINO B. M.S .

iii

ARTÍCULO 23 DE LA RESOLUCIÓN No. 13 DE JUNIO DE 1946�La universidad no se ha e responsable de los on eptos emitidos por sus alumnos ensus proye tos de grado.Solo velará porque no se publique nada ontrario al dogma y a la moral atóli a yporque los trabajos no ontengan ataques o polémi as puramente personales. Antesbien, que se vea en ellos el anhelo de bus ar la verdad y la justi ia.�

iv

Índi e general1. INTRODUCCIÓN 11.1. Revisión bibliográ� a . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11.2. Objetivos del proye to . . . . . . . . . . . . . . . . . . . . . . . . . . . 31.3. Organiza ión del do umento . . . . . . . . . . . . . . . . . . . . . . . . 42. MODELADO DEL INVERNADERO 52.1. Variables ambientales . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72.1.1. Radia ión . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72.1.2. Temperatura . . . . . . . . . . . . . . . . . . . . . . . . . . . . 82.1.3. Humedad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92.1.4. Dióxido de arbono . . . . . . . . . . . . . . . . . . . . . . . . . 112.1.5. Ventila ión . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 112.2. Modelo dinámi o . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122.2.1. Pro esos de temperatura . . . . . . . . . . . . . . . . . . . . . . 132.2.2. Pro esos de vapor de agua . . . . . . . . . . . . . . . . . . . . . 162.2.3. Modelo total . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 173. MÉTODOS DE IDENTIFICACIÓN DE SISTEMAS 183.1. Métodos no paramétri os . . . . . . . . . . . . . . . . . . . . . . . . . . 203.1.1. Análisis en tiempo . . . . . . . . . . . . . . . . . . . . . . . . . 203.1.2. Análisis en fre uen ia . . . . . . . . . . . . . . . . . . . . . . . . 213.2. Métodos paramétri os . . . . . . . . . . . . . . . . . . . . . . . . . . . 223.2.1. Métodos autoregresivos . . . . . . . . . . . . . . . . . . . . . . . 23v

3.2.2. Métodos por subespa ios . . . . . . . . . . . . . . . . . . . . . . 264. DISEÑO DE EXPERIMENTOS 344.1. Invernaderos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 344.2. Instrumenta ión . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 354.3. Experimentos preliminares . . . . . . . . . . . . . . . . . . . . . . . . . 364.3.1. Experimento 1. Gradientes espa iales . . . . . . . . . . . . . . . 374.3.2. Experimento 2. Comportamiento del follaje . . . . . . . . . . . 374.4. Experimentos para identi� a ión . . . . . . . . . . . . . . . . . . . . . 385. RESULTADOS 405.1. Experimentos preliminares . . . . . . . . . . . . . . . . . . . . . . . . . 405.1.1. Varia ión horizontal . . . . . . . . . . . . . . . . . . . . . . . . 415.1.2. Varia ión verti al . . . . . . . . . . . . . . . . . . . . . . . . . . 415.2. Pro esamiento previo . . . . . . . . . . . . . . . . . . . . . . . . . . . . 435.3. Modelo de aire . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 455.3.1. Estima ión del orden . . . . . . . . . . . . . . . . . . . . . . . . 455.3.2. Algoritmo de identi� a ión . . . . . . . . . . . . . . . . . . . . . 465.4. Modelo de aire y follaje . . . . . . . . . . . . . . . . . . . . . . . . . . . 495.4.1. Estima ión del orden . . . . . . . . . . . . . . . . . . . . . . . . 505.4.2. Algoritmo de identi� a ión . . . . . . . . . . . . . . . . . . . . . 515.5. Valida ión . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 535.6. Compara ión on otros resultados . . . . . . . . . . . . . . . . . . . . . 546. CONCLUSIONES 56

vi

Índi e de �guras2.1. Flujos de alor que se onsideran en el modelo del invernadero. Tomadode [16℄. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 52.2. Distribu ión general del espe tro solar dentro y fuera de la atmosfera.Tomado de [20℄. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72.3. Comportamiento de la radia ión dire ta en la ubierta. . . . . . . . . . 72.4. Modelo de ventila ión enital. Tomado de [1℄. . . . . . . . . . . . . . . 123.1. Diagrama en bloques de una estru tura ARX. . . . . . . . . . . . . . . 233.2. Diagrama en bloques de una estru tura ARMAX. . . . . . . . . . . . . 243.3. Diagrama en bloques de una estru tura ARARMAX. . . . . . . . . . . 243.4. Diagrama en bloques de una estru tura Box−Jenkins. . . . . . . . . . . 243.5. Diagrama en bloques de una estru tura en variables de estado. . . . . . 273.6. Interpreta ión de la proye ión obli ua en R

2. Tomado de [19℄ . . . . . 304.1. Plano de una nave on dimensiones. . . . . . . . . . . . . . . . . . . . . 354.2. Foto de la estru tura del invernadero. . . . . . . . . . . . . . . . . . . . 354.3. Vista frontal del te ho del invernadero. . . . . . . . . . . . . . . . . . . 354.4. Ubi a ión del invernadero usado en la primera parte del experimento . 374.5. Ubi a ión de los sensores en el primer experimento, parte a. . . . . . . 374.6. Ubi a ión del invernadero usado en la segunda parte del experimento . 374.7. Ubi a ión de los sensores en el primer experimento, parte b. . . . . . . 384.8. Distribu ión verti al de sensores. . . . . . . . . . . . . . . . . . . . . . 384.9. Ubi a ión de la instrumenta ión exterior. . . . . . . . . . . . . . . . . . 394.10. Distribu ión de los sensores a nivel de follaje. . . . . . . . . . . . . . . 39vii

5.1. Ilustra ión de la medida de desvia ión en la temperatura. . . . . . . . . 415.2. Comportamiento verti al de la temperatura. . . . . . . . . . . . . . . . 425.3. Comportamiento verti al de la humedad relativa. . . . . . . . . . . . . 425.4. Modelo re�nado del lima en el invernadero . . . . . . . . . . . . . . . 435.5. Valor máximo de humedad absoluta en fun ión de la temperatura . . . 435.6. Comportamiento del oe� iente de transmisión durante un día . . . . . 445.7. Autovalores de la proye ión obli ua on i = 5 . . . . . . . . . . . . . . 465.8. Autovalores de la proye ión obli ua on i = 10 . . . . . . . . . . . . . 465.9. Autovalores de la proye ión obli ua on i = 15 . . . . . . . . . . . . . 465.10. Error por entual de simula ión de la humedad para diferentes órdenes ytamaños de las matri es de Hankel. . . . . . . . . . . . . . . . . . . . . 475.11. Error por entual de simula ión de la temperatura para diferentes órdenesy tamaños de las matri es de Hankel. . . . . . . . . . . . . . . . . . . . 475.12. Humedad medida y simulada. . . . . . . . . . . . . . . . . . . . . . . . 495.13. Temperatura medida y simulada. . . . . . . . . . . . . . . . . . . . . . 495.14. Humedad relativa medida y simulada. . . . . . . . . . . . . . . . . . . . 495.15. Autovalores de la proye ión obli ua on i = 5 . . . . . . . . . . . . . . 505.16. Autovalores de la proye ión obli ua on i = 10 . . . . . . . . . . . . . 505.17. Autovalores de la proye ión obli ua on i = 15 . . . . . . . . . . . . . 515.18. Error por entual de simula ión de la humedad del aire para diferentesórdenes y tamaños de las matri es de Hankel. . . . . . . . . . . . . . . 515.19. Error por entual de simula ión de la temperatura del aire para diferentesórdenes y tamaños de las matri es de Hankel. . . . . . . . . . . . . . . 515.20. Error por entual de simula ión de la humedad del follaje para diferentesórdenes y tamaños de las matri es de Hankel. . . . . . . . . . . . . . . 525.21. Error por entual de simula ión de la temperatura del follaje para difer-entes órdenes y tamaños de las matri es de Hankel. . . . . . . . . . . . 525.22. Humedad del aire medida y simulada. . . . . . . . . . . . . . . . . . . . 525.23. Temperatura del aire medida y simulada. . . . . . . . . . . . . . . . . . 525.24. Humedad relativa del aire medida y simulada. . . . . . . . . . . . . . . 52viii

5.25. Humedad del follaje medida y simulada. . . . . . . . . . . . . . . . . . 525.26. Temperatura del follaje medida y simulada. . . . . . . . . . . . . . . . 525.27. Humedad relativa del follaje medida y simulada. . . . . . . . . . . . . . 525.28. Humedad del aire medida y simulada. . . . . . . . . . . . . . . . . . . . 535.29. Temperatura del aire medida y simulada. . . . . . . . . . . . . . . . . . 545.30. Humedad relativa del aire medida y simulada. . . . . . . . . . . . . . . 545.31. Humedad del follaje medida y simulada. . . . . . . . . . . . . . . . . . 545.32. Temperatura del follaje medida y simulada. . . . . . . . . . . . . . . . 555.33. Humedad relativa del follaje medida y simulada. . . . . . . . . . . . . . 55

ix

Índi e de uadros5.1. Desvia iones horizontales de temperatura y humedad en el aire, durantelos experimentos previos. . . . . . . . . . . . . . . . . . . . . . . . . . . 415.2. Errores obtenidos on el modelo de orden 2 para temperatura y humedad. 505.3. Errores obtenidos on el modelo de orden 3 para temperatura y humedadde aire y follaje. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 525.4. Errores obtenidos on el modelo de�nitivo para datos nuevos. . . . . . . 535.5. Errores obtenidos por otros autores en estudios similares. . . . . . . . . 55

x

Capítulo 1INTRODUCCIÓNAunque la agri ultura es una industria milenaria, sólo desde ha e algunos años se hanempezado a utilizar sistemas de informa ión y ontrol automáti o para optimizar losre ursos empleados en la produ ión agrí ola e in rementar la alidad y antidad de losprodu tos ultivados. En parti ular, desde la dé ada de los sesenta se han estudiado es-trategias de ontrol que permiten optimizar la produ ión de ultivos bajo invernadero.Para optimizar la produ ión en un invernadero se debe regular el lima al que essometido el ultivo, minimizando el osto de fun ionamiento del sistema. La solu ión alproblema planteado requiere ono er un modelo del omportamiento del lima dentrodel invernadero, des rito por la temperatura del aire, la humedad relativa, la ilumina ióny la on entra ión de CO2. Además se requieren modelos dinámi os de la evolu ión del ultivo, su re imiento, respira ión y fotosíntesis, entre otros.Este proye to se enfo a ha ia el primer paso en el desarrollo de un invernadero automa-tizado, la onstru ión de un modelo dinámi o que permita prede ir el omportamientodel lima dentro de éste. Esta herramienta es fundamental en las siguientes etapas deplanteamiento de estrategias de ontrol.1.1. Revisión bibliográ� aLos prin ipales estudios sobre modelamiento y ontrol de lima en invernaderos se hanrealizado en ondi iones limáti as y e onómi as diferentes a las olombianas. El grupo1

liderado por el profesor G. P. Bot, en la universidad de Wageningen (Holanda), hadesarrollado modelos dinámi os de invernaderos on base en pro esos biológi os, estosmodelos han sido validados en invernaderos ubi ados en Holanda [1℄. Adi ionalmentehan estudiado regímenes limáti os para ahorrar energía en los sistemas de ontrol[2, 6, 8℄.El grupo de P.F. Martínez, en la Universidad Polité ni a de Valen ia (España), hadesarrollado métodos de identi� a ión, basados en algoritmos genéti os, para en ontrarlos parámetros de un modelo físi o altamente no lineal del lima en invernaderos, y losha validado en invernaderos sometidos a lima mediterráneo [7℄.J.B. Cunha, del Centro de estudos e de te nologias do ambiente e vida de la Univer-sidade de Tras-Os-Montes e Alto Douro (Portugal), ha planteado té ni as re ursivasde identi� a ión, ha realizado estudios omparativos sobre diferentes modelos limáti- os de invernaderos y ha planteado estrategias de manejo óptimo del lima en inver-naderos [10�12℄.Mohammed Yassine El Ghoumari, en la Universidad Autónoma de Bar elona (España),desarrolló un sistema de ontrol predi tivo no lineal para invernaderos, manejandorestri iones. El modelo de lima fue desarrollado desde prin ipios físi os y alibradomediante Sequiential Quadrati Programming [4℄.El grupo de Peter Young, del Centre for resear h on environmental systems and sta-tisti s, en Lan aster University, ha desarrollado métodos estadísti os de identi� a ióny ontrol [14, 15℄, y los ha apli ado en el ontrol automáti o del lima en invernaderosubi ados en Bedford, Reino Unido.Los estudios anteriormente itados se han desarrollado en Europa, región en la que,debido a su ubi a ión geográ� a, las temperaturas extremas que se presentan en inviernoy verano ha en ne esario que los sistemas de ontrol optimi en el osto de opera iónde los sistemas de alefa ión y aire a ondi ionado, ne esarios para no someter a lasplantas a temperaturas extremas que ausen daños permanentes a los ultivos.2

1.2. Objetivos del proye toEste proye to tiene omo objetivo prin ipal Obtener y validar un modelo dinámi odel lima de un invernadero ubi ado en la sabana de Bogotá, y usarlo para generarinforma ión que sea útil en el manejo del invernaderoEste objetivo se ha subdividido en los siguientes objetivos espe í� os:Diseñar los experimentos y sele ionar los sensores que permitan adquirir infor-ma ión ade uada para la identi� a ión del modelo del lima.Obtener y validar un modelo dinámi o del lima en un invernadero ubi ado enla sabana de Bogotá, on base en los fenómenos físi os, quími os y biológi osexistentes.Plantear estrategias de manejo del invernadero, on base en el omportamientodel modelo obtenido, que eviten llevar el lima a ondi iones ríti as.A tualmente la gran mayoría de ultivos de �ores en Colombia no posee sistemas au-tomáti os de ontrol de lima, y su instala ión es inviable debido a limita iones e onómi- as. La onstru ión de un invernadero tradi ional en la sabana de Bogotá tiene un osto aproximado entre 3 y 5 dólares por metro uadrado, mientras que un invernaderoautomatizado puede ostar hasta 50 dólares por metro uadrado *.El prin ipal problema que en este momento poseen los ultivos de �ores en la sabana,rela ionado on el ontrol del lima, es el ontrol de los hongos Mildeo velloso y Mildeopolvoso. Estos hongos se reprodu en en ondi iones de alta humedad y temperatura,por lo que es ne esario evitar que el lima en el invernadero al an e estas ondi iones.Los modelos dinámi os de lima que se presentan en este proye to permiten prede irel omportamiento del lima en un invernadero en la sabana de Bogotá y usarlos paragenerar informa ión útil en el manejo del invernadero.*Informa ión suministrada por el Ingeniero William Torres del Departamento de Ingeniería deAmeri a�or Ltda3

1.3. Organiza ión del do umentoEn el apítulo dos se de�nen las variables que des riben el lima dentro del invernadero,los fenómenos que afe tan su evolu ión y su intera ión. Por último se presenta un mod-elo dinámi o del invernadero, basado en las leyes físi as que des riben los fenómenos.En el apítulo tres se presenta el problema de la identi� a ión de sistemas y algunasalternativas de solu ión tipo aja negra, usando modelos lineales y no lineales. Se de-s riben en detalle los algoritmos por sub-espa ios, ya que son la alternativa usada eneste proye to.En el apítulo uatro se des riben los invernaderos usados en el proye to, la instru-menta ión y el sistema de adquisi ión de datos. Además, se plantean los experimentosrealizados para validar las hipótesis del modelo y para obtener los datos requeridos porel algoritmo de identi� a ión.En el apítulo in o se presentan los resultados del trabajo, ini iando por los experimen-tos preliminares, que justi� an la aproxima ión de parámetros on entrados ha iendouna modi� a ión a la estru tura. El pro esamiento previo al que se someten los datosobtenidos en los experimentos de identi� a ión. La obten ión de un modelo dinámi oque des ribe la evolu ión de la temperatura y la humedad del aire en el invernadero.La obten ión de un modelo extendido que des ribe el omportamiento de temperaturay humedad del aire sobre el ultivo y del aire atrapado por el follaje. Por último, sepresenta la valida ión del modelo obtenido on un nuevo onjunto de datos y estosresultados se omparan on los reportados por otros autores.Finalmente en el apítulo seis se presentan las on lusiones del proye to y posiblesmejoras y trabajos futuros.

4

Capítulo 2MODELADO DEL INVERNADEROUn invernadero es una estru tura que en ierra una por ión de terreno, destinada a laagri ultura. Su fun ión es modi� ar las ondi iones ambientales a las que se somete el ultivo en bus a de un aumento en la alidad y rentabilidad de la a tividad agrí ola.El prin ipal bene� io que se obtiene on el uso de invernaderos es un mejor aprove hamien-to de una fuente gratuita de energía omo lo es el sol. Pero las ondi iones ambientalesdentro del invernadero también dependen del omportamiento del lima fuera del in-vernadero y las a tividades realizadas dentro de éste. La �gura 2.1 muestra los �ujosde alor presentes en el sistema.Otro fa tor que afe ta de manera importante la evolu ión del lima dentro del inver-nadero es el aire atrapado dentro de la estru tura, que presenta un inter ambio on elexterior mu ho menor al que permite un ambiente abierto.Para hablar de lima es ne esario de�nir ini ialmente las variables que lo des riben y ómo afe tan el re imiento del ultivo. La prin ipal variable limáti a usada en ultivosbajo invernadero es la temperatura, y la prin ipal fun ión de un invernadero es elevarla temperatura media dentro de la estru tura algunos grados por en ima de la mediaexterior. Este he ho permite ultivar espe ies en limas más fríos que aquellos en losque omúnmente re e la espe ie y además, permite in rementar la produ ión.Otra variable que des ribe el lima es la humedad. Su efe to sobre el ultivo no es menosimportante que el de la temperatura, ya que regula la transpira ión de las plantas,fa tor determinante de la alidad de las mismas. Existe un problema on el manejo de5

Figura 2.1: Flujos de alor que se onsideran en el modelo del invernadero. Tomadode [16℄.la humedad y es la existen ia de una amplia antidad de de�ni iones e interpreta ionesde esta variable; entre ellas la humedad relativa, usada omúnmente en meteorología;presión par ial de vapor, dé� it de satura ión, ontenido de humedad, entre otras.Todas estas de�ni iones se en uentran rela ionadas y son útiles en diferentes análisis.Más adelante se realizará una breve des rip ión de algunas de ellas y se indi ará ualesse usan en este proye to.La radia ión o �luz� es otro fa tor que determina el lima. Es la fuente de energía quepromueve la evolu ión de las variables ambientales. Aunque pueden instalarse fuenteslumíni as arti� iales, en este proye to sólo está presente la radia ión solar omo fuentede energía. Para el análisis de esta variable también existen multiples de�ni iones quese des ribirán más adelante.La presen ia de dióxido de arbono (CO2) en el ambiente es fundamental para el pro esode fotosíntesis, por ello es in luido por varios autores en sus modelos de lima. En esteproye to no se in luyó el CO2 omo variable ambiental debido a que el osto de la6

instrumenta ión es elevado. La ex lusión del CO2 del modelo no redu e la alidad delos resultados, dado que esta variable no afe ta la evolu ión de la temperatura o lahumedad, aunque depende de ellas.En la primera parte de este apítulo se des riben las variables ambientales que deter-minan el lima en el invernadero. En la segunda se desarrolla un modelo dinámi o del lima en el invernadero a partir de las e ua iones físi as que regulan el omportamientode ada variable.2.1. Variables ambientales2.1.1. Radia iónLa radia ión solar es la prin ipal fuente de energía en un invernadero. Es el motor dela fotosíntesis, fuente de energía térmi a y generador de inter ambio de aire ha ia elexterior.En la des rip ión del lima en invernaderos se divide el espe tro ele tromagnéti o endos bandas [1℄, [3℄. La radia ión de onda orta (luz visible) y la radia ión de onda larga(infrarrojo), ver �gura 2.2.La radia ión de onda orta orresponde a la banda entre los 400 y 700 nm, de interésdebido a que es la banda que a tiva la fotosíntesis. Esta radia ión es denominada PAR,por sus siglas en inglés Photosyntheti A tive Radiation y la unidad de medida esµmolm−2s−1.La radia ión de onda larga orresponde a la banda entre 700 y 3000 nm. Esta es laradia ión en infrarrojo, que sirve para alentar el invernadero, por lo tanto es la prin ipalfuente para el balan e de energía. La unidad de medida es Wm−2.Otra forma de lasi� ar la radia ión solar es por sus ara terísti as geométri as. Dea uerdo a su fuente y traye toria seguida, se lasi� a en radia ión dire ta y difusa [3℄ [4℄.La radia ión dire ta es aquella que llega dire tamente de la fuente lumíni a, en este aso, el sol. La transmisión de radia ión dire ta ha ia el invernadero se puede al ular, ono iendo la posi ión del sol y las ara terísti as ópti as de la ubierta, a través de laley de Snell, ver �gura 2.3. 7

Figura 2.2: Distribu ión general del espe tro solar dentro y fuera de la atmosfera. Toma-do de [20℄.

Figura 2.3: Comportamiento de la radia ión dire ta en la ubierta.8

La radia ión difusa es aquella que proviene de todo el ielo, el ejemplo más laro esla radia ión re ibida en un día nublado. Se ara teriza por no tener una dire iónparti ular, por lo que la ley de snell no des ribe su transmisión a través de la ubierta.Para modelar su transmisión en la ubierta se usa un oe� iente global de transmisión,independiente de la dire ión de arribo de la onda.En la sabana de Bogotá, la radia ión es prin ipalmente difusa, on sólo 5 horas diarias,en promedio, de radia ión dire ta [20℄. Por esta razón, en el modelo no es ríti o teneren uenta la varia ión temporal de oe� iente de transmisión, debido al ambio en laposi ión del sol.2.1.2. TemperaturaLa temperatura es un fa tor regulador de la a tividad metabóli a de las plantas. Existenregímenes de temperatura determinados para ada espe ie en los uales su desarrollo esóptimo. Además, existen límites biológi os, fuera de los uales, las plantas sufren dañospermanentes.La temperatura del aire dentro del invernadero se obtiene al realizar un balan e deenergía en el volumen limitado por la estru tura. Existen �ujos de alor por fenómenosde ondu ión, onve ión y radia ión. La úni a fuente de radia ión en este proye toes el sol. Por otro lado, se pierde alor a través del suelo, del aire que ir ula ha iael exterior y parte del alor es absorbido por el ultivo. Además, existe �ujo de alorlatente en la transpira ión de las plantas.En este modelo se asume que la temperatura del aire es uniforme en el volumen delinvernadero, hipótesis que ha sido validada en el experimento 1.En la literatura se en uentran diferentes aproxima iones en la onstru ión de modelosde temperatura de las plantas, que en prin ipio es diferente a la temperatura del aire.Bot [1℄ plantea un modelo en el que el follaje es una apa homogenea de masa que sólointera túa on el aire. Wang y Boulart [16℄ usan una apa uniforme de vegeta ión, peroesta intera túa on la ubierta, el aire y el suelo, generando un modelo más detallado.La mayoría de autores , omo Trigui [21℄,Ghoumari [4℄, Piñón [22℄, Herrero [7℄, entreotros, aproximan el omportamiento de la temperatura del follaje al de la temperatura9

del aire. Esta aproxima ión se evalúa en el experimento 2.2.1.3. HumedadLa humedad es la variable más onfusa en la literatura sobre ontrol de lima, debidoa la multipli idad de de�ni iones que existen de ella. En su de�ni ión más bási a, lahumedad es el ontenido de vapor de agua en el aire, medido en kg/m3.Humedad absoluta (χ): Masa de vapor de agua por unidad de volumen.χ :=

[kg

m3

] (2.1)Autores omo Herrero [7℄ expresan la humedad absoluta de forma adimensional, de�nién-dola omo la on entra ión de masa de vapor de agua por masa de aire se o.x =

kgH2O

kgaireseco

(2.2)Presión par ial de vapor (ea): La ley de Presiones par iales de Dalton expresa que lapresión total de una mez la de gases ideales es igual a la suma de las presiones par ialesde ada omponente. Esta ley gobierna todos los gases que omponen el aire, in luyendoel vapor de agua. Re ordando la ley de los gases ideales:PV = nRT (2.3)Donde P es la presión a la que es sometido el gas en Pa, V es el volumen que o upael gas en m3, n es el número de moles, R es la onstante de los gases ideales y T es latemperatura en K.Usando la e ua ión (2.3) es posible obtener la presión par ial de vapor de agua ea, ono iendo la humedad absoluta (2.1), la temperatura T y el peso mole ular del agua

Mv = 0,018kgmol−1, ver [3℄.ea =

χTR

Mv

=

[N

m2

]= [Pa] (2.4)Presión de satura ión(es): El vapor de agua presenta una ara terísti a que di�eredel omportamiento de los gases ideales; la existen ia de una antidad máxima de10

masa de agua, en forma de vapor, que el aire puede soportar; esta es una propiedadtermodinámi a que se aso ia a la antidad de energía interna (entalpía) de las molé ulasde vapor y es fun ión de la temperatura. La forma de expresarla es a través de la máximapresión par ial de vapor de agua que puede soportar el aire a una temperatura dada,se exprese en Pa.La rela ión entre la temperatura y la presión de satura ión no es simple, existen tablasque las rela ionan, y algunos autores han desarrollado e ua iones que aproximan lasrela iones obtenidas on datos empíri os, ver por ejemplo [3℄. La aproxima ión usadaen este proye to es la propuesta por Woodward y Sheehy en [23℄:ln es = A−

B

T− C lnT (2.5) on:

A = 57,96Pa

B = 6731PaK

C = 4,796PaDé� it de presión de vapor (V PD): Es la diferen ia entre la presión de satura ióny la presión par ial de vapor existente a una temperatura dada. Esta de�ni ión es muyimportante en el pro eso de transpira ión de las plantas, ya que es en fun ión de ésta quelas plantas liberan agua al ambiente en su sistema de transpira ión, fa tor determinatedel re imiento y alidad de las plantas.V PD = es − ea (2.6)Humedad relativa (RH): La humedad relativa es el por entaje o o iente, entre lapresión par ial de vapor existente y la presión de satura ión, a una temperatura dada.

RH = 100 ∗

[ea

es

]

T

= [ %] (2.7)La humedad relativa es la medida usada omúnmente ya que es la que per ibe la piel;aunque para la formula ión de modelos en in onveniente ya que es fun ión de dosvariables físi as, la temperatura y la humedad absoluta.11

2.1.4. Dióxido de arbonoEl dióxido de arbono o CO2 es la materia prima en el pro eso de fotosíntesis. Aunqueno es una variable ambiental, en los últimos años se ha venido onsiderando en lossistemas de optimiza ión y automatiza ión de invernaderos, debido a su in�uen ia enel desarrollo del ultivo, ver [3, 4℄.El prin ipal problema on el CO2 es que durante el día, en los periodos de mayorradia ión solar, la planta onsume grandes antidades de este gas, y si el invernaderono tiene un inter ambio de aire su� iente on el exterior, se redu e su on entra ión y ausa un de remento en el pro eso de fotosíntesis. Si la fotosíntesis se redu e, se detienetambién el re imiento de la planta, dañando así la produ ión.La antidad de CO2 libre en el aire se expresa usualmente en partes por millón o onsu presión par ial omo gas ideal.Existen diferentes métodos para determinar la on entra ión de CO2 en el aire, basadosen prin ipios quími os, térmi os, elé tri os y ópti os. Los sistemas más usados, debidoa su osto, son aquellos basados en prin ipios ópti os, que aprove ha las ara terísti asopa as del CO2 en el rango del infrarrojo (efe to invernadero).2.1.5. Ventila iónLa ventila ión es el �ujo de aire entre el exterior y el invernadero. Este �ujo es promovi-do por dos fenómenos: el gradiente de temperatura entre el interior y el exterior, y las orrientes de aire que ingresan al invernadero. Bot [1℄ presenta un modelo de ventila iónen el que uanti� a estos fenómenos. Ghoumari [4℄ presenta un modelo de ventila ión,basado en los trabajos de Kozai y Sase [5℄.Bot y Ghoumari modelan la ventila ión a través de una ventana enital de aperturavariable, omo la que se muestra en la �gura 2.4. A ontinua ión se reprodu en losresultados de Bot para los �ujos de ventila ión:Ventila ión debida al gradiente térmi o:12

Figura 2.4: Modelo de ventila ión enital. Tomado de [1℄.φT =

L0H320

3

(gβ

F0

) 12

(sinψ)32 (∆T )

12 (2.8)Ventila ión debida al viento:

φu =1

2L0H0

(Kf (ζ)

F0

) 12

u (2.9)donde:L0 es el an ho de la ventana [m].H0 es la longitud de la ventana [m].ζ es el ángulo de apertura de la ventana [rad].ψ es el ángulo de in lina ión del te ho [rad].β es el oe� iente de expansión térmi a [K−1].g es la a elera ión gravita ional [ms−2].F0 es el fa tor de forma de la ventana [adimensional].Kf (ζ) es el oe� iente de �u tua ión de presión de la ventana [adimensional].2.2. Modelo dinámi oEn este apartado se desarrollan las e ua iones dinámi as que des riben el ompor-tamiento del lima en el invernadero. Ini ialmente se modela el omportamiento de la13

temperatura del aire y del suelo, estas e ua iones se obtienen a partir de la ley de on-serva ión de la energía. A ontinua ión, se onstruye el modelo de vapor de agua en elaire, a partir de la ley de onserva ión de masa.En sentido estri to, la evolu ión del lima en ada punto del invernadero es diferente,siendo este un sistema de parámetros distribuidos, des rito por e ua iones diferen ialespar iales, on ondi iones de frontera dadas por ada interfaz on el medio exterior. Sinembargo, un modelo de este tipo produ e un sistema de in�nitas dimensiones, on unaaltísima omplejidad que no es útil para realizar análisis interpretativos ni en el diseñode sistemas de ontrol.En la revisión bibliográ� a realizada, se en ontró que la totalidad de los autores que hanrealizado estudios sobre modelos de lima han usado modelos de parámetros on en-trados omo aproxima ión al pro eso del lima, bajo la hipótesis de que las variablesambientales de interés son uniformes dentro del invernadero, tanto horizontal omoverti almente.Existe otra área de estudio, dentro del tema de lima, en la que se usan modelos deparámetros distribuidos, denominada dinámi a de �uidos omputa ional o CFD por sussiglas en inglés. La metodología usada en esta área se basa en la onstru ión de modelostotalmente teóri os y su simula ión por métodos numéri os. El objetivos de esta áreaes en análisis del omportamiento de la ventila ión en espa ios errados, usando unaaproxima ión netamente numéri a.2.2.1. Pro esos de temperaturaLa temperatura de un uerpo o volumen se obtiene al realizar un balan e de energíaentre el uerpo y su entorno. A ontinua ión se presentan los balan es de energía parael aire en el invernadero, el follaje y el suelo.Autores omo Wang [16℄ y Young [15℄ realizan adi ionalmente un balan e de energíapara la ubierta de la estru tura, tomando omo referen ia los modelos de Bailey yDeltour. La razón para adi ionar esta e ua ión es que las ubiertas usadas en la mayoríade invernaderos en el norte de Europa son onstruidas usando dos apas de vidrio, on unintermedio de aire entre ellas; aumentando así el aislamiento térmi o ha ia el ambiente14

y reando una apa itan ia térmi a importante.En el aso parti ular de los invernaderos modelados en este proye to esta e ua ión noes ne esaria, debido a que el material usado en las ubiertas es polipropileno, materialdelgado y de muy baja apa itan ia térmi a. Algunos grupos españoles y portugueses [7,11℄, que usan invernaderos on ubiertas plásti as, usan la misma aproxima ión seguidaen este proye to, modelando la ubierta omo una resisten ia térmi a, sin apa idad dealma enar energía.Balan e de energía en el aireEl volumen de ontrol, delimitado por el aire en errado en el invernadero, re ibe alorde la radia ión solar a través de la ubierta, pierde alor por ondu ión y onve ióna través del suelo, del aire que ir ula ha ia el exterior y del ultivo. La e ua ión quedes ribe este fenómeno es la siguiente:Caire

dTa

dt= Qrad −Qas −Qac −Qv (2.10)Donde:

Caire es la apa idad térmi a del aire del invernadero. Esta depende del volumen del in-vernadero, la densidad del aire y el alor espe í� o del aire [WK−1]. Una aproxima ión,usada por Herrero et. al. [7℄ es:Caire = V caireρaire (2.11)

Qrad Es el �ujo de energía que proviene de la radia ión solar, depende prin ipalmentede la radia ión re ibida en el exterior, es de ir el lima exterior, pero también dependede otros fa tores omo la in lina ión del sol, que varía on la hora y la épo a del año;la in lina ión relativa del te ho del invernadero y el alibre o olor del plásti o usadoen el re ubrimiento. [W ]

Qas Es el �ujo de alor que se presenta a ausa del gradiente térmi o existente entreel aire y el suelo. El �ujo se genera por ondu ión y onve ión natural. Un modelosu� iente para este fenómeno es:Qas = kas (Ta − Ts) (2.12)15

Qac Es el �ujo de alor a través de la ubierta, generado por el gradiente térmi o entreel aire exterior y el interior. El �ujo se genera por ondu ión y depende de la apa idadde ondu ión térmi a de la ubierta. Un modelo ade uado de este �ujo es:Qac = kac (Ta − Text) (2.13)

Qv es el �ujo de alor debido al inter ambio de aire entre el interior y el exterior.La ventila ión es ausada por el gradiente térmi o entre el interior y el exterior y las orrientes de aire en el exterior, omo se des ribió en la se ión 2.1.5. Cal ulando laenergía que transporta el aire ir ulante se tiene:Qv = φNρairecaire (Ta − Text) (2.14)Mientras la temperatura del aire dentro del invernadero sea mayor a la temperaturaexterna y a la temperatura del suelo, que es la ondi ión normal de opera ión, el úni o�ujo de alor que ingresa al volumen de ontrol es el de radia ión; los demás �ujos tienesigno negativo en (2.10), lo que indi a que son pérdidas de energía para el invernadero.Balan e de energía en el sueloEl volumen de ontrol para el modelo de temperatura del suelo es una apa de tierra queintera túa on el aire y on el subsuelo. Autores omo [1℄, usan una apa on un grosorde 15cm en sus modelos. La e ua ión que des ribe el omportamiento de la temperaturadel suelo es:

CsdTs

dt= Qas −Qss (2.15)Donde:

Cs es la apa itan ia térmi a del suelo, que depende de su volumen, densidad y alorespe í� o.Qas es el �ujo de alor que se presenta a ausa del gradiente térmi o existente entre elaire y el suelo, omo se des ribió en la se ión anterior.Qss es el �ujo de alor que se presenta a ausa del gradiente térmi o entre el suelo yel subsuelo. En este modelo se onsidera que el subsuelo es una variable de entrada,16

no manipulable y que evolu iona de manera independiente al lima del invernadero. El�ujo de alor se genera por ondu ión, siendo un modelo su� iente el siguiente:Qss = kss (Ts − Tss) (2.16)Algunos autores in luyen hasta 7 apas de subsuelo omo variables dependientes, adauna des rita on una e ua ión similar a (2.15), ver [16℄.De (2.15) se puede on luir que el suelo se alienta tomando alor del aire del inver-nadero, mientras éste presente una temperatura más alta, y pierde energía en su interfaz on el subsuelo, mientras que éste presente una temperatura inferior.2.2.2. Pro esos de vapor de aguaComo se expuso en la se ión 2.1.3, la humedad es la antidad de vapor de agua presenteen el aire; medida omo masa de vapor de agua por unidad de volumen (χ = [kg/m3]).Para uanti� ar la evolu ión de la humedad en un volumen limitado se realiza unbalan e de masa de vapor de agua en la frontera del volumen.En la literatura se en uentran múltiples formas de la e ua ión de onserva ión, formu-ladas a partir de diferentes variables. Por ejemplo Tap et. al [6℄ presenta un modelo enfun ión de la humedad absoluta, Herrero et. al [7℄ formula su modelo para la on en-tra ión de vapor de agua, esta di�ere de la formula ión de Tap en un fa tor de es ala.Bot [1℄ onstruye un modelo para la presión par ial de vapor que resulta diferente delos dos anteriores debido a su dependen ia de la temperatura.A ontinua ión se presenta la formula ión para la humedad absoluta propuesta por Tap.Balan e de masa de vapor de agua en el aireEl volumen de ontrol, delimitado por el aire en errado en el invernadero, re ibe vaporde las plantas produ to de la transpira ión, inter ambia vapor on el aire exteriorprodu to de la ventila ión, y pierde vapor uando se presenta ondensa ión en lasparedes del invernadero y en las hojas. Realizando un balan e de masa se obtiene lasiguiente e ua ión: 17

VH

dχ(t)

dt= φvent(χext − χ) + CsatE −Mc (2.17)Donde:

VH es el volumen efe tivo del invernadero.φvent(χext − χ) es el �ujo de vapor de agua ha ia o desde el exterior.E es la transpira ión del ultivo generada por fenómenos biológi os. Se en uentranmúltiples formula iones que uanti� an este �ujo. Cunha et. al. [13℄ presenta una aprox-ima ión simple, propor ional al de� it de presión de vapor, mientras que Herrero et.al. [7℄ usa un modelo altamente no lineal, dependiente, además del de� it de presión devapor, de la radia ión solar, el alor de vaporiza ión del agua, la ondu tan ia de lashojas, entre otros parámetros que varían on el lima.Csat Es la onstante de satura ión que modela la satura ión de agua en el aire. Csattoma un valor de 1 mientras el aire no esté saturado y onmuta a 0 uando el aire sesatura de vapor de agua, anulando el fenómeno de transpira ión.Mc Es la ondensa ión de vapor en las paredes del invernadero y en las hojas del ultivo.Este fenómeno solo se presenta uando la temperatura del aire des iende hasta al an-zar el punto de ro ío, prin ipalmente en la no he o uando hay lluvia. Es importantedesta ar que no ne esariamente todo el aire debe al anzar el punto de ro ío, solo aquelque se en uentra er a a una frontera on un gradiente térmi o importante.2.2.3. Modelo totalAl reunir las e ua iones que des riben el omportamiento de la temperatura y lahumedad, se obtiene un onjunto de e ua iones diferen iales a opladas que representanuna formula ión en variables de estados del lima dentro del invernadero.

dTa

dt= 1

Caire(Qrad −Qas −Qac +Qv)

dTs

dt= 1

Cs(Qas −Qss)

dχ(t)

dt= 1

VH(φvent(χext − χ) + CsatE −Mc)

(2.18)18

Capítulo 3MÉTODOS DE IDENTIFICACIÓNDE SISTEMASEn el apítulo anterior se presentó un modelo dinámi o del lima en un invernadero, onstruido a partir de prin ipios físi os y leyes fundamentales. El resultado de la formu-la ión es un onjunto de e ua iones diferen iales no lineales, on una gran antidad deparámetros que dependen de ara terísti as (geométri as, biológi as, et .) de los om-ponentes que integran el sistema. Este tipo de formula ión es ade uada para realizaranálisis interpretativos sobre los fenómenos que o urren dentro del sistema, aunque re-alizar la estima ión de los parámetros del modelo es una tarea ompleja y dispendiosaque requiere la realiza ión de mu hos experimentos espe í� os para determinar el valorde ada parámetro.Bot [1℄ presenta los resultados de su tesis do toral, en donde parte importante de susaportes es la estima de los parámetros de un invernadero ubi ado en Holanda, usandoprin ipios fundamentales y experimenta ión.Herrero et. al [7℄ usan algoritmos genéti os para en ontrar 15 parámetros del modelono lineal de la temperatura de aire y suelo, y humedad del aire de un invernaderoubi ado en Valen ia (España). La metodología usada en este trabajo tiene un alto osto omputa ional, requiriendo más de 80 horas para la eje u ión del algoritmo.Como se eviden ia en las referen ias anteriores, uando los sistemas que se quierenmodelar son de alta omplejidad, obtener modelos a partir de prin ipios físi os puede19

impli ar un esfuerzo no justi� ado o in luso ser inviable, omo es el aso de algunossistemas distribuidos. Además, los modelos ompletos y no lineales no son ade uadospara en análisis y diseño de sistemas de ontrol, ya que su alta omplejidad los ha e dedifí il manipula ión en el mar o de la teoría de ontrol.Otra metodología para la onstru ión de modelos dinámi os de fenómenos físi oses la identi� a ión de sistemas, de�nida por Ljung en [18℄ omo: �Construir mode-los matemáti os de sistemas dinámi os on base en datos experimentales obtenidos delsistema�. Esta dis iplina tiene una amplia gama de apli a iones, ya que los sistemasdinámi os se en uentran en múltiples ampos de la ien ia.La onstru ión de un modelo a partir de medi iones requiere tres elementos:Un onjunto de datos: Es el punto de partida para el pro eso de identi� a ión,la alidad de los datos determina en gran medida la oheren ia entre el modeloobtenido on el sistema real.Un onjunto de modelos: La mayoría de métodos de identi� a ión requierenque el usuario delimite el tipo de estru tura a la que debe pertene er el modeloque se bus a. Por ejemplo modelos autorregresivos, modelos lineales en variablesde estados, redes neuronales, et .Una norma o medida de alidad: Se debe espe i� ar una regla que permitamedir de forma uantitativa el desempeño de un modelo en su objetivo de repro-du ir el onjunto de datos. La medida más usada es el error medio uadráti o depredi ión, mostrado en 5.5.Jm =

1

N

N∑

k=1

([y(k) − y(k)]2) (3.1)Donde y(k) es la salida produ ida por el sistema y y(k) es la predi ión realizadapor el modelo, on base en las entradas y salidas de los tiempos anteriores.A ontinua ión se presenta una breve des rip ión de varias té ni as de identi� a iónbien estable idas, on el �n de ubi ar al le tor en el ontexto de los métodos usados eneste proye to. 20

3.1. Métodos no paramétri osSe denominan métodos no paramétri os a aquellas té ni as de identi� a ión de sistemasen las que no se de�ne un onjunto de modelos posibles, por lo que no existe un espa iode parámetros en el ual realizar la búsqueda del modelo que más se ajusta a los datos.Asumiendo que el onjunto de datos disponible para la identi� a ión fue generado porun sistema lineal, es posible obtener fun iones ara terísti as del sistema a partir de losdatos; en parti ular es posible obtener estimados de su respuesta impulso y su fun iónde transferen ia.3.1.1. Análisis en tiempoDado un sistema en tiempo dis reto:y(t) = G0(q)u(t) + v(t) (3.2)Donde q es el operador desplazamiento, u(t) es una señal de entrada manipulable y v(t)es una realiza ión de un pro eso esto ásti o des ono ido. Es posible obtener la respuestapulso del sistema 3.2 si se ono e la entrada apli ada. A ontinua ión se presentan losresultados para algunas entradas típi as. Para más detalles ver [18,24,25℄.Respuesta pulso: Si un sistema des rito por 3.2 es sometido a una entrada pulso deamplitud α, su� iente para que su efe to en la salida sea mayor al del ruido, es posibleestimar la respuesta impulso del sistema omo:

g(t) =y(t)

α(3.3)Respuesta paso: Si el sistema des rito por 3.2 es sometido a una entrada tipo pasode amplitud α, su respuesta impulso puede re uperarse a partir de la e ua ión:

g(t) =y(t) − y(t− 1)

α(3.4)Correla ión: Si al sistema des rito por 3.2 se le apli a ruido blan o on fun ión deauto orrela ión:

Ru(τ) = αδ(τ) (3.5)21

Es posible obtener la respuesta impulso del sistema a partir de:g(t) =

Ryu(τ)

α(3.6)Donde Ryu(τ) es un estimado de la orrela ión ruzada entre la entrada y la salida delsistema.3.1.2. Análisis en fre uen iaUna de las prin ipales ara terísti as de los sistemas lineales es su respuesta ante unaentrada sinusoidal, que resulta de la misma fre uen ia de la señal apli ada, variando suamplitud y su fase.Una primera aproxima ión a la estima ión de la fun ión de transferen ia de un sistemalineal es medir su respuesta ante un onjunto de entradas sinusoidales en la banda defre uen ias de interés.

G(jωi) =Y (ωi)

U(ωi)(3.7)Este es un pro eso tedioso y altamente sensible a ruido.Otro enfoque para determinar fun iones de transferen ia de sistemas lineales es el análi-sis espe tral. En este método se bus a obtener el o iente 3.7 en un solo experimento;para ello se apli a al sistema una señal on un amplio ontenido en fre uen ia y, uti-lizando la transformada de Fourier, se estima el ontendido espe tral de las señales deentrada y salida. Los métodos por análisis espe tral son menos sensibles al ruido quelos métodos en tiempo y la estima ión fre uen ia por fre uen ia.Una forma de mejorar los resultados del análisis espe tral es adi ionar a la hipótesisde linealidad del sistema, la hipótesis de que la fun ión de transferen ia es una señal ontinua y suave. En la prá ti a esto se logra onstruyendo un promedio de los estima-dos de la fun ión de transferen ia en fre uen ias er anas a través de una ventana depondera ión. Las ventanas más usadas son las de Barlett, Hamming y Parzen.Los métodos no paramétri os permiten obtener informa ión del omportamiento delsistema, pero no entregan un modelo matemáti o que lo aproxime. Estos métodos son22

útiles para realizar un primer a er amiento al omportamiento del pro eso que se quiereidenti� ar y ganar ono imiento sobre el mismo; se usa omo parte de una prueba previa.3.2. Métodos paramétri osCuando el objetivo del pro eso de identi� a ión es obtener un modelo matemáti o quedes riba de manera aproximada el omportamiento del sistema en estudio, el problemade identi� a ión se puede formular de la siguiente manera:Dada una estru tura de modelos M, parametrizada por un ve tor θ ∈ DM ⊂ Rd, onmodelos parti ulares M(θ); se de�ne el onjunto de modelos posibles omo:

M∗ = {M(θ)|θ ∈ DM} (3.8)Considérense los onjuntos de datos:Y N = [y(1) y(2) . . . y(N)] (3.9)UN = [u(1) u(2) . . . u(N)] (3.10)(3.11)Produ idos por un sistema dinámi o:y(t) = G(q, θ)u(t) +H(q, θ)e(t) (3.12)Con e(t) ruido blan o.El objetivo es de idir, on base en la informa ión ontenida en WN = [Y N , UN ], uál esel mejor ve tor de parámetros θN , y por lo tanto el mejor modelo M(θN), que expli alos datos.De manera formal, un algoritmo de identi� a ión es un operador que realiza un mapeo:

[Y N , UN ] → θN ∈ DM (3.13)23

Figura 3.1: Diagrama en bloques de una estru tura ARX.3.2.1. Métodos autoregresivosExiste un onjunto de algoritmos de identi� a ión, bien ara terizados y estudiados,que usan omo onjunto de modelos estru turas autoregresivas. Estos se diferen ianentre si por los grados de libertad que brindan para modelar el efe to del ruido.La estru tura más sen illa es la denominada autorregresiva on entrada externa o ARXpor su nombre en inglés (AutoRegressive with eXogeneous input). Esta estru tura estádes rita por la e ua ión de diferen ias:y(t)+a1y(t−1)+a2y(t−2)+. . .+ana

y(t−na) = b1u(t−1)+b2u(t−2)+. . .+bnbu(t−nb)+e(t)(3.14)El ve tor de parámetros para esta estru tura es:

θ = [a1 a2 . . . anab1 b2 . . . bnb

] (3.15)De�niendo los polinomios:A(q) = 1 + a1q

−1 + a2q−2 + . . .+ ana

q−na (3.16)B(q) = b1q

−1 + b2q−2 + . . .+ bnb

q−nb (3.17)Un diagrama en bloques del sistema representado por 3.14 se muestra en la �gura 3.1.La prin ipal desventaja de usar el modelo ARX es la falta de grados de libertad paramodelar las ara terísti as de los disturbios e(t).24

Figura 3.2: Diagrama en bloques de una estru tura ARMAX.Una estru tura que aumenta la �exibilidad en el tipo de disturbios que se pueden in luiren el modelo, es la auroregresiva on promedio móvil y entrada externa o ARMAX, porsus siglas en inglés. En este aso se asume que la entrada de disturbio es un promediomobil de un ruido blan o, omo lo muestra la e ua ión 3.18. Un diagrama de bloquesde esta estru tura se muestra en la �gura 3.2. El osto de agregar el promedio móvil alruido es el in remento de la dimensión del espa io de búsqueda.y(t) + a1y(t− 1) + · · · + ana

y(t− na) = b1u(t− 1) + b2u(t− 2) + · · · + bnbu(t− nb)+

c1e(t− 1) + c2e(t− 2) + · · · + cnce(t− nc) (3.18)El ve tor de parámetros para esta estru tura es:

θ = [a1 a2 · · · anab1 b2 · · · bnb

c1 c2 · · · cnc] (3.19)Existen otras estru turas, ada vez más omplejas, que permiten obtener dinámi as�exibles para el ruido. Entre ellas la ARARMAX, que adi iona un �ltro IIR para elruido, de�nido por los oe� ientes deD(q), antes de que ingrese al sistema, ver �gura 3.3.Y la estru tura Box−Jenkins que desa opla totalmente la dinámi a debida a la entradamanipulada, de�nida por [B(q) A(q)], de la debida al ruido, de�nida por [C(q) D(q)], omo se muestra en la �gura 3.4.Hasta el momento se han presentado las estru turas usadas para representar el on-junto de modelos posibles. Ahora se des riben los algoritmos de identi� a ión usados25

Figura 3.3: Diagrama en bloques de una estru tura ARARMAX.

Figura 3.4: Diagrama en bloques de una estru tura Box−Jenkins.26

omúnmente on estas estru turas.Mínimos uadrados: El algoritmo mas usado es el de mínimos uadrados, que bus aminimizar el error medio uadráti o de predi ión, de�nido en 5.5. De manera formal,este algoritmo se de�ne omo:θLS

N = arg mın

[1

N

N∑

k=1

1

2[y(k) − y(k)

2)

] (3.20)Para las estru turas autoregresivas, este riterio es una forma uadráti a en θ, quepuede ser minimizada analíti amente omo un problema sen illo de regresión lineal.En el mar o de las estru turas autoregresivas, orresponde al usuario de�nir el gradode ada uno de los polinomios que omponen el sistema. Por ejemplo, en el aso de unaestru tura ARX es el usuario quien debe de idir los valores de na y nb, on la úni arestri ión de que na ≥ nb. En el aso de estru turas más �exibles omo la Box−Jenkins,es ne esario de�nir uatro parámetros na, nb, nc y nd.Máxima verosimilitud: Otro enfoque en el planteamiento del problema de identi�- a ión es asumir que las medidas disponibles son realiza iones de un pro eso esto ásti o.Suponiendo que la fun ión de densidad de probabilidad de las medidas yN es:f(θ;x1, x2, . . . , xN) = fy(θ;xN) (3.21) on xi una variable aleatoria, que al realizar el experimento toma el valor yi∗, y θ el onjunto de parámetros determinísti os que des ribe la fun ión.Luego de realizado el experimento, la fun ión 3.21 toma la forma:

fy(θ, yN∗ )Esta es una fun ión determinísti a en θ y lo que se debe bus ar es un onjunto deparámetros θ, tal que el onjunto de medidas observadas sea el más probable posible.De manera formal, el estimado on máxima verosimilitud es:

θML(yN∗ ) = arg max

θfy(θ; y

N∗ ) (3.22)No es posible obtener una solu ión analíti a al problema planteado. Los algoritmosexistentes se aproximan al estimado óptimo de forma iterativa.27

Variables instrumentales: De manera ideal, el error de predi ión que genera unbuen modelo no debe estar orrela ionado on los datos pasados. De no ser así, es desuponer que existe más informa ión en los datos que no ha sido usada para onstruirel modelo. Es posible a�rmar que un buen modelo es aquel que produ e errores depredi ión independientes de los datos pasados.Aunque estri tamente, para veri� ar que el error es independiente del pasado se requierela totalidad del onjunto in�nito de datos, esto no es viable en la prá ti a. En lugar deesto, dada una se uen ia �nita de datos, es posible estimar la orrela ión del error depredi ión ε(t, θ), on una ierta transforma ión de los datos pasados ζ(t, θ). Es de ir:1

N

N∑

t=1

ζ(t, θ)α(ε(t, θ)) (3.23)Y el ve tor θ que es el mejor estimado on base en los datos disponibles es:θIV (yN

∗ ) = sol[ 1

N

N∑

t=1

ζ(t, θ)α(ε(t, θ)) = 0

]0 (3.24)3.2.2. Métodos por subespa iosOtro tipo de estru tura para formular sistemas lineales es la representa ión en el espa iode estados. La razón para usar este tipo de modelos es la fa ilidad para extender losmétodos de identi� a ión de manera natural a sistemas multivariable, y los algoritmosque solu ionan el problema de identi� a ión son numéri amente e� ientes y robustos.Matemáti amente, estos modelos están des ritos por el siguiente onjunto de e ua ionesde diferen ias de primer orden:

xk+1 = Axk +Buk + wk (3.25)yk = Cxk +Duk + vk (3.26) on:

E[

(wp

vp

) (wT

q vTq

)=

Q S

ST R

δpq (3.27)28

Las variables involu radas en el modelo son:uk ∈ R

m es el ve tor orrespondiente a las m entradas del pro eso en el instante k.yk ∈ R

l es el ve tor orrespondiente a las l salidas del pro eso en el instante k.xk ∈ R

n es el ve tor de estado del pro eso en el instante k.wk ∈ R

n y vk ∈ Rl son señales ve toriales esto ásti as, no medibles, que representan elruido de entrada y salida del sistema. Se asume que son se uen ias ve toriales de ruidoblan o, esta ionarias y de media ero.Las matri es se denominan:

A ∈ Rn×n es la matriz ara terísti a del sistema. Des ribe totalmente la dinámi a delpro eso, y está ara terizada por sus autovalores.

B ∈ Rn×m es la matriz de entrada que representa la transforma ión lineal entre elespa io de m entradas manipuladas y el espa io de n estados.

C ∈ Rl×n es la matriz de salida que des ribe omo los estados internos son transferidosa las l medidas yk.

D ∈ Rl×m es la matriz de transferen ia dire ta. Aunque en los sistemas ontinuos estetérmino usualmente es ero, no su ede en el aso de sistemas dis retos debido al efe todel muestreo.

Q ∈ Rn×n, R ∈ R

l×l y S ∈ Rn×l son las matri es de ovarianza de las se uen ias deruido wk y vk.La �gura 3.5 muestra un diagrama de bloques de una estru tura en variables de estado.Las requisitos que debe umplir el sistema a identi� ar son mínimos y su requerimientoresulta obvio. En primer lugar, el par {A, C} debe ser observable, lo que signi� a quetodos los modos del sistema pueden ser observados desde la salida yk. En segundolugar, el par {A, [B Q1/2]} debe ser ontrolable, lo que impli a que todos los modos delsistema pueden ser ex itados, ya sea por la entrada determinísti a uk o por la entradaesto ásti a wk.De las dos ondi iones expuestas se deriva el nombre del método, ya que el modeloidenti� ado solo ontendrá el subespa io observable y ontrolable del pro eso.La parametriza ión expuesta presenta varias ventajas sobre las parametriza iones au-toregresivas. Primero que todo, la formula ión es igual para sistemas SISO y para29

Figura 3.5: Diagrama en bloques de una estru tura en variables de estado.sistemas MIMO. No es ne esario re urrir a modi� a iones o extensiones para manejarestos sistemas.Además, en esta formula ión toda la dinámi a del sistema está ontenida en la matrizA. Los autovalores de esta matriz des riben todos los modos que los datos ontienen yque provienen del sistema real, sin importar que sean ex itados por las entradas deter-minísti as o por las esto ásti as. Esta es una ara terísti a que no poseen los modelosautorregresivos que dividen la dinámi a determinísti a, que expli an el sistema real, yla dinámi a del ruido, que expli an el efe to de las entradas esto ásti as. Por ejemplolos modelos de Box−Jenkins presentados en la se ión 3.14. El he ho de ondensar to-das las dinámi as en un solo modelo brinda omo resultado una realiza ión mínima delmodelo del sistema.A diferen ia de las estru turas autorregresivas, que requieren que el usuario �je var-ios parámetros, orrespondientes al orden de todos los polinomios de la estru tura; laparametriza ión del onjunto de modelos posibles es muy sen illa en los métodos porsubespa ios, y solo ontiene un parámetro ajustable que orresponde al orden del sis-tema. Además, este parámetro puede ser obtenido exa tamente en sistemas sin presen iade ruido o estimado fá ilmente en otro aso.El on epto en el que se basan los algoritmos por subespa ios es la posibilidad deen ontrar la se uen ia de estados xk a partir de los datos de entrada y salida medidos.30

Una vez que los estados son ono idos, el algoritmo de identi� a ión se redu e a unproblema de mínimos uadrados lineal para en ontrar las matri es del sistema. Se puedede ir que los algoritmos por subespa ios no identi� an modelos entrada−salida, sinomodelos entrada−estado−salida.A ontinua ión se presenta de manera general la formula ión del problema para el asodeterminísti o y su solu ión, tomada de [19℄.Problema de identi� a ión determinísti a:Dadas s medidas de la entrada uk ∈ Rm y la salida yk ∈ R

l generadas por un sistemalineal determinísti o des ono ido, de orden n:xk+1 = Axk +Buk (3.28)yk = Cxk +Duk (3.29)Determinar:El orden n del sistema des ono idoLas matri es del sistema A ∈ R

n×n, B ∈ Rn×l, C ∈ R

l×n y D ∈ Rl×m, bajo unatransforma ión similar.El he ho de que las matri es del sistema sólo se puedan re uperar bajo una transfor-ma ión similar se debe a que existen in�nitas bases que representan al mismo sistemay presentan la misma rela ión entrada−salida. Esto no genera problemas, ya que elobjetivo del pro eso de identi� a ión es pre isamente re uperar la rela ión entrada−salida.Proye ión obli ua: Antes de plantear la solu ión al problema de identi� a ión, esne esario de�nir la proye ión obli ua de ve tores, omo una extensión de la proye iónortogonal. Ya que es la herramienta fundamental en los algoritmos por subespa ios.El operador que permite proye tar el espa io �la de una matriz A ∈ R

p×j en el espa io�la de una matriz B ∈ Rq×j es:

ΠB , BT · (BBT )† ·B (3.30)31

Donde •† es la pseudoinversa de la matriz •. De manera resumida la proye ión ortogonalde A en B se nota:A/B , A · ΠB (3.31)Siendo B− el omplemento ortogonal del espa io �la de B, es posible des omponer lamatriz A omo una ombina ión lineal de las olumnas de B y de las olumnas de B−,asi:

A = A/B ·B + A/B− ·B− (3.32)En lugar de des omponer A omo una ombina ión de dos matri es ortogonales, esposible expresarla omo una ombina ión de dos matri es no ortogonales B y C, y del omplemento ortogonal de éstas. Esto es:A = LB ·B + LC · C + LB−,C− ·

B

C

− (3.33)La matriz Lc · C se de�ne omo la proye ión obli ua del espa io �la de A a lo largodel espa io �la de B en el espa io �la de C. La �gura 3.6 ilustra la proye ión obli uade un ve tor en R2.De manera formal, se de�ne la proye ión obli ua del espa io �la de A ∈ R

p×ja lo largodel espa io �la de B ∈ Rq×j en el espa io �la de C ∈ R

r×j omo:A/BC , A(CT BT ) ·

CCT CBT

BCT BBT

†primeras r olumnas · C (3.34)

A ontinua ión se realizan algunas de�ni iones ne esarias a la hora de presentar elalgoritmo de identi� a ión.Matri es de Hankel: Ini ialmente se de�nen algunas matri es de Hankel y otras ma-tri es rela ionadas on el sistema. 32

Figura 3.6: Interpreta ión de la proye ión obli ua en R2. Tomado de [19℄Matriz de Hankel de entrada:

U0|2i−1 ,

u0 u1 u2 . . . uj−1

u1 u2 u3 . . . uj

. . . . . . . . . . . . . . .

ui−1 ui ui+1 . . . ui+j−2

ui ui+1 ui+2 . . . ui+j−1

ui+1 ui+2 ui+3 . . . ui+j

. . . . . . . . . . . . . . .

u2i−1 u2i u2i+1 . . . u2i+j−2

(3.35)U0|2i−1 ,

U0|i−1

Ui|2i−1

,

Up

Uf

(3.36)

33

U0|2i−1 ,

u0 u1 u2 . . . uj−1

u1 u2 u3 . . . uj

. . . . . . . . . . . . . . .

ui−1 ui ui+1 . . . ui+j−2

ui ui+1 ui+2 . . . ui+j−1

ui+1 ui+2 ui+3 . . . ui+j

. . . . . . . . . . . . . . .

u2i−1 u2i u2i+1 . . . u2i+j−2

(3.37)U0|2i−1 ,

U0|i

Ui+1|2i−1

,

U+

p

U−f

(3.38)El número de �las (i) es de�nido por el usuario y debe al menos ser mayor que elmáximo orden del sistema. El número de olumnas (j) se elige igual a s− 2i+1, lo queimpli a que todos los datos disponibles son usados en la identi� a ión.Los subíndi es usados en las de�ni iones de las matri es ha en referen ia a la primeray la última �la de la matriz. Los subíndi es �p� y �f � ha en referen ia al �pasado� y�futuro�. Las matri es U+ y U− se de�nen desplazando la división entre el pasado y elfuturo un bloque de Hankel ha ia abajo.Las matri es de Hankel de salida se de�nen de manera similar a las matri es de entrada(3.35), (3.36), (3.37) y (3.38). Su nota ión es: Y0|2i−1, Yp, Yf , Y +

p y Y −f .Adi ionalmente, se de�nen otras matri es de Hankel formadas por las entradas y salidas:

W0|i−1 ,

U0|i−1

Y0|i−1

,

Up

Uf

, Wp

(3.39)De manera similar se de�ne W+

p omo:W+

p =

U+

p

Y +p

34

Otro elemento importante es la se uen ia de estados determinísti os, de�nida omo:Xd

i ,(xd

i xdi+i . . . x

di+j−2 x

di+j−1

) (3.40)Matri es del sistema:Se de�ne la matriz extendida de observabilidad Γi, on i > n, omo:Γi ,

C

CA

CA2

. . .

CAi−1

(3.41)La matriz reversa y extendida de ontrolabilidad ∆i se de�ne omo:

∆i =(Ai−1B Ai−2B . . . AB B

) (3.42)Por último se de�ne la matriz Toeplitz triangular inferior Hi omo:Hi ,

D 0 0 . . . 0

CB D 0 . . . 0

CAB CB D . . . 0

. . . . . . . . . . . . . . .

CAi−2B CAi−3B CAi−4B . . . D

(3.43)Algoritmo de identi� a ión:1. Cal ular las proye iones obli uas:

Oi = Yf/UfWp (3.44)Oi−1 = Y −

f /Uf−W+p (3.45)2. Realizar la des omposi ión en valores singulares, SVD, de la proye ión obli ua:

W1OiW2 = USV T (3.46)35

Con W1 y W2 matri es de pondera ión de�nidas por el usuario. W1 ∈ Rli×li es derango ompleto y W2 ∈ R

j×j es tal que rank(Wp) = rank(WpW2).3. Dividir el resultado de la des omposi ión en valores singulares de la siguientemanera:W1OiW2 = (U1 U2)

S1 0

0 0

V T

1

V T2

(3.47)Y determinar el orden del sistema inspe ionando el número de valores singularesdiferentes de ero en S.4. Determinar Γi y Γi−1 omo:

Γi = W−11 U1S

1/21 (3.48)

Γi−1 = Γi (3.49)5. Obtener las se uen ias de estados Xi y Xi−1 omo:Xi = Γ†

iOi (3.50)Xi−1 = Γ†

i−1Oi−1 (3.51)6. Resolver el sistema de e ua iones lineales para A, B, C y D: Xi+1

Yi|i

=

A B

C D

Xi

Ui|i

(3.52)El algoritmo des rito realmente abar a una familia de métodos de identi� a ión, or-respondientes a diferentes valores de las matri es de pondera ión W1 y W2. Ljung [18℄realiza un listado de diferentes valores de las matri es de pondera ión y el algoritmoque generan.

36

Capítulo 4DISEÑO DE EXPERIMENTOSEste proye to se desarrolló en la �n a Splendor Flowers El Rosal, ubi ada en el mu-ni ipio de El Rosal Cundinamar a. Esta �n a, de propiedad de la empresa Ameri a�orLtda., posee 122 invernaderos dedi ados al ultivo de rosas, distribuidos en 50 he táreas.La empresa Ameri a�or Ltda. posee un departamento de investiga ión y desarrollo queha venido estudiando los efe tos del lima en la evolu ión de las plantas y prin ipalmenteen la in iden ia de enfermedades. El grupo posee las instala iones y la instrumenta iónade uada para registrar la evolu ión del lima en un invernadero y alma enar las me-didas en un omputador.A ontinua ión se des riben los invernaderos usados en el proye to, la instrumenta iónutilizada y los experimentos realizados.4.1. InvernaderosLa unidad estru tural de un invernadero es la nave. Una nave es una estru tura re -tangular formada por 16 segmentos o uadros, on te ho a dos aguas, y uya alturaaumenta ha ia el entro del invernadero. La �gura 4.1 muestra el plano de una nave.El an ho de una nave es de 6.8 m y el alto en los extremos es de 2.9 m, al anzando 4.3m en el entro.Los invernaderos estándar se omponen de do e naves iguales, formando un bloque de77m de an ho por 82m de largo, ubriendo un área de 6314m2 y abar ando un volumen37

Figura 4.1: Plano de una nave on dimensiones.aproximado de 23,000m3.El invernadero está ubierto en su totalidad por plásti o semitransparente de polietileno.Es posible desplazar las paredes ha ia arriba mediante un sistema de adenas, ver �gura4.2.El te ho del invernadero está diseñado de tal manera que posee bajantes omo desagüesy permite la ventila ión enital, sin que ingrese agua al interior a ausa de la lluvia, ver�gura 4.3.El área de la apertura enital es un parámetro importante en el modelo de ventila iónnatural del invernadero, en el aso de los invernaderos utilizados, éstos poseen unaapertura de 370m2, que orresponden al 5.6% del área ubierta por el bloque.4.2. Instrumenta iónPara medir y alma enar las variables ambientales se utilizaron sistemas de registroWat hDog 450 *, produ idos por la empresa Spe trum Te hnologies y de propiedad del*http://spe mters. om 38

Figura 4.2: Foto de la estru tura del invernadero.departamento de investiga ión de Ameri a�or Ltda.Cada instrumento de medida está en apa idad de registrar temperatura y humedadrelativa on una pre isión de ± 0,6 ◦C y ± 3% HR, respe tivamente. Además, adaequipo uenta on dos anales de entrada adi ionales para instalar sensores de radia iónsolar, humedad en hoja, entre otros.La radia ión solar se midió on un sensor de radia ión de onda orta (PAR Quantumsensor), produ ido también por la empresa Spe trum Te hnologies. Como se des ribióen la se ión 2.1.1, este sensor mide la radia ión entre 400 y 700 nm. La medida deradia ión PAR se transformó en �ujo de energía en el espe tro ompleto(W/m2). Esta39

Figura 4.3: Vista frontal del te ho del invernadero.transforma ión es posible debido a que el espe tro ara terísti o de la luz solar es bien ono ido, ver [16,17℄.El intervalo de muestreo fue �jado en in o minutos en todos los experimentos. Lasele ión de este periodo se hizo on base en un estimado ini ial de la apa idadtérmi a del invernadero. Teniendo en uenta que el volumen de aire atrapado es de23,000m3, la densidad del aire se en uentra alrededor de 1,2kg/m3 y el alor espe í� odel aire es de 1000J/kg◦C; la apa idad térmi a del invernadero es de aproximadamente27,6 · 106J/◦C. El tiempo de muestreo sele ionado está dentro de las espe i� a ionesde los instrumentos utilizados y es omparable on los tiempos usados por otros autores,por ejemplo [9, 11,13,16℄.Para evitar desvia iones en las medidas de temperatura, ada equipo se instala dentro deuna aja prote tora, semejante a la ubierta de una esta ión meteorológi a. Su fun iónes proteger al dispositivo de la radia ión solar dire ta, de las orrientes de viento y de lalluvia, que ausan errores en la medida; el osto de introdu ir esta prote ión es obteneruna redu ión en la velo idad de respuesta del dispositivo, ya que se in rementan su apa itan ia y resisten ia térmi as.

40

Figura 4.4: Ubi a ión del invernadero usado en la primera parte del experimento4.3. Experimentos preliminaresPara obtener un modelo sen illo del lima en el invernadero, el primer paso es redu ir la omplejidad del modelo, aproximando su omportamiento al de un sistema de parámet-ros on entrados, bajo la hipótesis de que el lima es uniforme dentro del invernadero,tanto horizontal omo verti almente. Los primeros experimentos realizados bus an val-idar esta hipótesis, midiendo los gradientes limáti os dentro del invernadero.A ontinua ión se des riben los dos experimentos realizados.4.3.1. Experimento 1. Gradientes espa ialesEl objetivo del experimento es medir los gradientes en las variables limáti as a lolargo del invernadero; para ello se realizaron dos series de medi iones de temperaturay humedad, usando diferentes on�gura iones de sensores.El primer experimento se realizó entre los días 8 y 16 de agosto del 2005, en un inver-nadero tipo er ha, ubi ado omo se indi a en la �gura 4.4.Para la primera serie de medidas, los sensores se instalaron a una altura de 2.5 m sobre41

Figura 4.5: Ubi a ión de los sensores en el primer experimento, parte a.el suelo, por en ima del ultivo. A lo largo del invernadero, los sensores se ubi aronuniformemente espa iados, entre las naves 3-4, 6-7 y 9-10. La �gura 4.5 muestra unplano del invernadero on la ubi a ión detallada de los instrumentos.La segunda serie de medidas se realizó entre los días 6 y 12 de septiembre del 2005,en un invernadero tipo er ha, ubi ado omo se indi a en la �gura 4.6. Durante lasegunda serie de medidas, los sensores se instaron a una altura de 2.5 m sobre el suelo,por en ima del ultivo y se distribuyeron en uatro puntos dentro del invernadero, ver�gura 4.7.4.3.2. Experimento 2. Comportamiento del follajeEn un segundo experimento se midieron los gradientes de temperatura y humedad delaire dentro del follaje y el aire sobre el ultivo. Este experimento se realizó en el mismoinvernadero usado en la primera parte del experimento anterior, �gura 4.4.Se instalaron 6 sensores de temperatura y humedad, divididos en dos niveles y tres42

Figura 4.6: Ubi a ión del invernadero usado en la segunda parte del experimento

Figura 4.7: Ubi a ión de los sensores en el primer experimento, parte b.43

Figura 4.8: Distribu ión verti al de sensores.zonas. El primer nivel orresponde a una altura de 1.2 m sobre suelo, dentro del follaje.El segundo nivel a 2.5 m del suelo, por en ima del ultivo, ver �gura 4.8. A lo largo delinvernadero, los sensores se ubi aron uniformemente espa iados, entre las naves 3-4, 6-7y 9-10, de la misma forma que se instalaron en la primera parte del experimento 1, ver�gura 4.5.4.4. Experimentos para identi� a iónSe realizaron dos series de medidas omo experimento de identi� a ión. Los datosobtenidos en la primera parte del experimento, tomados entre el 6 y el 12 de septiem-bre, fueron usados para obtener los modelos dinámi os y los resultados de la segundaparte, tomados entre el 30 de septiembre y el 2 de o tubre, fueron usados omo datosde valida ión.Debido a la naturaleza de los invernaderos en estudio, no existen variables manipuladasdisponibles durante los experimentos. No existen sistemas de alefa ión o ventila iónforzada, y la ventila ión natural, aunque es variable, es a tivada por un sistema me áni- o frágil, tedioso y lento de manipular. Modi� ar la apertura de las ventadas laterales44

de un invernadero toma más de 30 minutos de trabajo de un operario y la ventila ión enital es �ja.Esto quiere de ir que el modelo a identi� ar orresponde a un sistema en malla abier-ta, sujeto a todos los disturbios ambientales y arente de entradas manipuladas. Sinembargo, dado que es posible registrarlas, las entradas del sistema son las ondi ionesambientales externas y la radia ión solar. Para medir estas variables se instaló un equipode registro, igual a los usados en los experimentos preliminares, en un punto de la �n aalejado de las instala iones, para evitar la in�uen ia en la medida de posibles fuentesde alor o humedad, que no re�ejen las ondi iones ambientales globales. La �gura 4.9muestra el plano de la �n a on la ubi a ión del registrador exterior.Dentro del invernadero, se ubi aron uatro instrumentos de medida a nivel de aire, omo se indi a en la �gura 4.7; y uatro instrumentos dentro del follaje, a lo largo deuna nave, ver �gura 4.10.Durante el experimento las ventanas laterales permane ieron abiertas un 50%, omose muestra en la �gura 4.2. La razón es que este es el punto normal de opera ión delinvernadero.

45

Figura 4.9: Ubi a ión de la instrumenta ión exterior.

46

Figura 4.10: Distribu ión de los sensores a nivel de follaje.

47

Capítulo 5RESULTADOSEl objetivo de este proye to es obtener un modelo dinámi o del lima dentro del inver-nadero, para lograrlo se realizaron los experimentos des ritos en el apítulo anterior.Los primeros experimentos permiten validar la hipótesis de que es posible modelar elsistema omo un pro eso de parámetros on entrados, aunque se deben agregar algunasvariables para re�nar el modelo, omo se expone más adelante.Los experimentos �nales permiten obtener la informa ión requerida por los algoritmosde identi� a ión para onstruir un modelo lineal en variables de estado de la evolu ióndel lima dentro del invernadero.Por último, los modelos obtenidos permiten realizar varios análisis sobre el pro eso enestudio y su omportamiento ante situa iones que no se presentaron durante el pro esode identi� a ión.5.1. Experimentos preliminaresComo se men ionó en el apítulo anterior, el objetivo de estos experimentos es validarla hipótesis de que el lima es uniforme dentro del invernadero y que este se puedemodelar omo un sistema de parámetros on entrados. Para veri� ar que el lima esuniforme, se miden los gradientes de temperatura y humedad en los ejes x, y y z, endos invernaderos tipo er ha.Para uanti� ar los gradientes, en ada instante de tiempo se evalúa la desvia ión en la48

medida, de�nida omo la diferen ia entre el punto más alto y el más bajo; y se obtienela media de la desvia ión durante ada experimento. La �gura 5.1 ilustra la medida.

Figura 5.1: Ilustra ión de la medida de desvia ión en la temperatura.5.1.1. Varia ión horizontalPara el nivel de aire, a 2.5 m del suelo, los gradientes a lo largo y an ho del invernaderoson pequeños, del mismo orden de magnitud que la pre isión de los instrumentos usados, omo se muestra en la tabla 5.1.En el peor aso, la temperatura presenta una desvia ión de 0,7 ◦C y la humedad presentauna desvia ión del 3.4%. Estos resultados veri� an la viabilidad de utilizar un modelo on entrado, usando el promedio de las variables en el invernadero omo salida.5.1.2. Varia ión verti alAunque la mayoría de autores onsideran que el ambiente dentro del invernadero esuniforme verti almente, ver [1,7,21,22℄; en este proye to se en ontró que para los inver-naderos en estudio se presentan gradientes importantes entre el aire que se en uentra er a a las plantas y el aire sobre el ultivo.49

Cuadro 5.1: Desvia iones horizontales de temperatura y humedad en el aire, durantelos experimentos previos.La �gura 5.2 muestra el omportamiento de la temperatura en un punto del inver-nadero. En la �gura se observa que el aire que se en uentra er a al follaje presentaun retardo, tanto en el alentamiento omo en el enfriamiento. Además los ambiosrápidos de temperatura que sufre el aire superior son �ltrados dentro del follaje. Estosdos elementos sugieren la hipótesis de que el follaje tiene un efe to de �ltro pasabajosque aumenta el orden del sistema. Al analizar el omportamiento de la humedad, elresultado es diferente. El follaje siempre presenta una humedad mayor a la del aire,debida a la evapotranspira ión de las plantas. Las varia iones temporales son similaresen los dos puntos, no se puede de ir que la humedad en el follaje sea una versión �ltradade la humedad del aire. Más bien, es la misma señal on un desplazamiento que varía on las ondi iones ambientales.La evapotranspira ión del ultivo es el �ujo de agua ha ia el ambiente en forma devapor. Existen varios modelos empíri os que intentan des ribir este �ujo. Todos ellosmodelan este �ujo en fun ión de la radia ión solar re ibida por la hoja, la temperaturade la hoja y el de� it de presión de vapor. Se debe re ordar que el de� it de presión devapor es la diferen ia entre la presión par ial de vapor existente y aquella para la ualse presenta satura ión.Cuando el de� it de presión de vapor se a er a a ero, le es más difí il a la hoja emitirvapor al aire. Cuando el de� it se ha e ero es imposible que exista �ujo de vapor ha ia50

Figura 5.2: Comportamiento verti al de la temperatura.el ambiente. Este omportamiento es oherente on el omportamiento observado enel experimento. En la �gura 5.3 se muestra el omportamiento de la humedad en unpunto del invernadero.Con base en el análisis anterior se on luye que una aproxima ión de parámetros on- entrados es viable, aunque existen gradientes importantes en dire ión verti al, queobligan a orregir el modelo presentado en el apítulo 2.Se propone un modelo re�nado, ompuesto por dos volúmenes de ontrol, el primervolumen es el aire atrapado en el invernadero fuera del follaje que orresponde al aireinterior usado en la mayoría de aproxima iones. Adi ionalmente, se in luye en el modelootra masa que orresponde al aire en errado por el follaje. En la interfaz entre los dosvolúmenes se presentan inter ambios de alor por ondu ión y onve ión; y �ujos demasa de vapor de agua a ausa de gradientes de presión par ial. La �gura 5.4 presentade manera esquemáti a el modelo re�nado que se propone.51

00:00 02:00 04:00 06:00 08:00 10:00 12:00 14:00 16:00 18:00 20:00 22:00 00:0040

50

60

70

80

90

100

Gradiente vertical de humedad

Humedad cerca del follaje

Humedad del aire

Figura 5.3: Comportamiento verti al de la humedad relativa.

Figura 5.4: Modelo re�nado del lima en el invernadero52

5.2. Pro esamiento previoAntes de realizar el pro edimiento de identi� a ión es ne esario a ondi ionar los datos,para que se adapten a los requerimientos del algoritmo y la rela ión entre las variablesde entrada y salida sea lo más lineal posible.El primer paso es transformar la humedad relativa a humedad absoluta, omo se de�nióen el apítulo 2. La razón es que la e ua ión de onserva ión de masa de agua se planteaen términos de la humedad absoluta, siendo la relativa una fun ión no lineal de ésta yde la temperatura. Este paso desa opla la medida de humedad de la temperatura. Lae ua ión (5.1) muestra la rela ión entre la humedad absoluta (χ) y la relativa (HR),y la �gura 5.5 muestra el omportamiento de la máxima humedad absoluta posible enfun ión de la temperatura. La e ua ión y la �gura re�ejan el a ople no lineal entrehumedad y temperatura.HR = 100 ∗

χ · TR

Mv exp (A−BT−1 − C lnT )(5.1)En donde T es la temperatura en grados kelvin, R es la onstante universal de losgases, Mv es el peso mole ular del agua, A, B y C son onstantes on las dimensionesade uadas.En segundo lugar, para realizar la aproxima ión de parámetros on entrados, se obtieneel promedio de las medidas de temperatura y humedad de los uatro puntos medidos anivel de aire, el resultado es el valor que se usa omo salida del sistema para el pro esode identi� a ión.Finalmente, dado que para medir la radia ión solar se usa un medidor de radia iónPAR (sensible solo en el espe tro visible entre 400 y 700 nm) es ne esario, a partir deesta medida, obtener la poten ia radiada en el espe tro ompleto. Para esto se parte dela hipótesis de que toda la radia ión proviene del sol, hipótesis válida debido a que noexiste ilumina ión arti� ial en el invernadero. La rela ión entre la radia ión PAR y laradia ión total es:

Qrad

[W

m2

]= 3,6 ∗QPAR

[µMol

m2s

] (5.2)53

0 5 10 15 20 25 300

0.5

1

1.5

2

2.5

3

3.5

Temperatura (°C)

Hum

edad

(kg

/m3 )

Valor máximo de humedad absoluta en función de la temperatura

Figura 5.5: Valor máximo de humedad absoluta en fun ión de la temperaturaSiendo 3.6 un fa tor de es ala on las unidades ade uadas, omo lo expone [17℄.Debido a la varia ión de la in lina ión del sol en el trans urso de un día, el ángulode in iden ia de la radia ión dire ta sobre la ubierta es variable, por lo que algunosautores onsideran el sistema omo variante on el tiempo. Para uanti� ar la varia ióndel oe� iente de transmisión de la ubierta se midió la radia ión in idente en el exteriory en el interior del invernadero, a 1,2m del suelo. La �gura 5.6 muestra el o iente entrela radia ión interior y la exterior, junto on su valor medio (trazo en rojo). Aunque larela ión varía, durante todo el día permane e alrededor de la media. Las varia ionesson debidas al paso de personas er a a los medidores y a la distan ia existente entreel sensor externo y el invernadero.Es omún en los pro esos de identi� a ión de sistemas, de�nir las variables in remen-tales eliminando el valor medio de las señales y las tenden ias lineales, on el objetivode eliminar la informa ión del punto de equilibrio y on entrarse en la informa ión deseñal pequeña. Para el sistema en estudio no es onveniente realizar este pro esamientoprevio debido a que los �ujos de alor y vapor de agua son ausados por los gradientestérmi os o de humedad, medidos omo la diferen ia de los valores absolutos de las vari-54

Figura 5.6: Comportamiento del oe� iente de transmisión durante un díaables. Si se elimina la media de ada variable, las diferen ias en pequeña señal se alterane in luso pueden ambiar de signo, re�ejando una dinámi a diferente a la que se quiereidenti� ar.5.3. Modelo de aireEn el primer aso se toman omo salidas del sistema la temperatura y la humedad, al uladas omo el valor medio de las medidas obtenidas por los uatro registradoresinstalados sobre el ultivo, omo se des ribió en el apartado 4.4. En total se tienen 1651medidas que se obtuvieron durante 138 horas ontinuas.5.3.1. Estima ión del ordenCon los datos prepro esados, el primer paso en el pro eso de identi� a ión es estimarel orden del sistema. Como se expuso en el apítulo 3, el orden del sistema a identi� ar,55

en los algoritmos por subespa ios, se obtiene a partir de un análisis de los valoressingulares de la proye ión obli ua del subespa io de las salidas futuras, a lo largo delsubespa io de las entradas futuras, sobre el espa io de las entradas y salidas pasadas.De manera intuitiva, el orden del subespa io resultado de la proye ión orresponde ala informa ión de las salidas futuras que no puede ser expli ado por las entradas futurasy que por lo tanto depende de las entradas y salidas pasadas.El úni o parámetro de ajuste que tiene este primer paso del algoritmo es el tamaño delos bloques de Hankel usados en la onstru ión de las matri es de entradas y salidaspasadas, denominado i en la se ión 3.2.2.Teniendo en uenta el modelo dinámi o que se presentó en el apítulo 2, el orden delsistema puede ser 2 o 3, de a uerdo on la in�uen ia de la apa de suelo en el lima. La�gura 5.7 muestra los valores singulares de la proye ión obli ua, obtenidos on i = 5, enla que se observa laramente que los primeros dos valores singulares son mu ho mayoresa los demás, por lo tanto el orden del sistema que des ribe los datos es 2.Al aumentar el tamaño de los bloques de Hankel, �jando i = 10, la diferen ia entre losprimeros dos valores singulares y el ter ero es menos mar ada, aunque aún es notorioel ambio de tenden ia luego del segundo valor, omo se ve en la �gura 5.8. Al �jari = 15 ya los valores singulares no presentan un ambio de pendiente luego del segundo,sino del ter er valor singular, ver �gura 5.9. Si se ontinua aumentando el orden de losbloques de Hankel, los valores singulares presentan un de aimiento monótono que nopermite obtener un estimado del orden.Dado que la estima ión del orden varía on el número de �las usadas en la proye iónobli ua, esta prueba no es on luyente, aunque los resultados indi an que los datosfueron generados por un sistema de orden dos o tres. Este resultado on uerda perfe -tamente on el modelo teóri o onstruido a partir de prin ipios físi os en el apítulo2.5.3.2. Algoritmo de identi� a iónPara obtener un modelo dinámi o del sistema se usa el algoritmo robusto ombinado,propuesto por De Moor y Overs hee [19℄. Este es un algoritmo de identi� a ión de56

Figura 5.7: Autovalores de la proye ión obli ua on i = 5

Figura 5.8: Autovalores de la proye ión obli ua on i = 10

57

Figura 5.9: Autovalores de la proye ión obli ua on i = 15sistemas on entradas determinísti as y esto ásti as, que de a uerdo on los autores ese� iente y presenta buenos resultados on datos que provienen de sistemas no linealesy ruidosos. Cuando los datos provienen de un sistema lineal e invariante on el tiempo,este algoritmo es un estimador no sesgado de los polos y eros del sistema.Para medir la alidad de los resultados se usaron 3 medidas de desempeño. La medidamás usada en el área de identi� a ión de sistemas es el error uadráti o medio, de�nido omo:J2 =

√√√√ 1

N

N∑

k=1

((y(k) − y(k))2) (5.3)Otra forma de medir el error de un modelo es el error absoluto medio:J1 =

1

N

N∑

k=1

|y(k) − y(k)| (5.4)De Moor y Overs hee de�nen en error por entual de simula ión omo:ǫ = 100 ∗

√√√√[∑N

k=1[y(k) − y(k)]2∑Nk=1 y(k)

2

] (5.5)Con el �n de omparar los resultados on los obtenidos por otros autores, en este58

3 4 5 6 7 8 9 103

3.5

4

4.5

5

5.5

6

i

Err

or p

orce

ntua

l de

hum

edad

Orden 1Orden 2Orden 3Orden 4

Figura 5.10: Error por entual de simula ión de la humedad para diferentes órdenes ytamaños de las matri es de Hankel.proye to se obtuvieron las tres medidas de error de�nidas atrás, para ada una de lasvariables de salida.Dado que no fue posible estable er un orden de�nido on los valores singulares de laproye ión obli ua, se obtuvieron modelos de varios órdenes, usando diferentes tamañosde las matri es de Hankel. Para ada uno de ellos se evaluó el error por entual desimula ión para la temperatura y la humedad. Las �guras 5.10 y 5.11 muestran el omportamiento del error para los asos evaluados. El eje x orresponde al tamañode las matri es de Hankel (i), el eje y al error por entual de la variable y ada trazo orresponde a un orden del modelo.En las �guras se eviden ia que el modelo que mejor expli a los datos es el de orden 2.Aunque para i ∈ [6, 7, 8] el error en la humedad es menor on un modelo de orden 3,este modelo no logra expli ar la temperatura ade uadamente. Para de�nir el tamañode la matriz se debe ha er un ompromiso entre las dos variables de salida, ya que elerror de temperatura disminuye on i, mientras que el error de humedad aumenta.Se de ide usar el resultado del algoritmo on i = 9; éste es el tamaño de matriz onel que se obtiene el menor error para la temperatura, que es la variable on mayor59

3 4 5 6 7 8 9 105

5.5

6

6.5

7

7.5

8

8.5

9

9.5

10

i

Err

or p

orce

ntua

l de

tem

pera

tura

Orden 1Orden 2Orden 3Orden 4

Figura 5.11: Error por entual de simula ión de la temperatura para diferentes órdenesy tamaños de las matri es de Hankel.error. El riterio usado en la sele ión del tamaño de matriz es la minimiza ión deuna norma in�nita, ara terizada por brindar robustez al resultado, minimizando elerror en el peor aso. Aunque este no es el mejor modelo para la humedad, el resul-tado es su� ientemente bueno. En total, el error obtenido es ǫT = 5,79 % y ǫH = 3,92 %.A ontinua ión se presentan los resultados obtenidos on el modelo de�nitivo, onstruido on i = 9 y orden 2, para los mismos datos usados en el pro eso de identi� a ión. Enla se ión 5.5 se evalúa el modelo on otra serie de datos.La �gura 5.12 muestra el omportamiento de la humedad absoluta, en esta se observaque el modelo es apaz de prede ir el omportamiento de esta variable en periodoslargos, aunque las os ila iones de alta fre uen ia las atenúa. Esta es la variable quemejor expli a el modelo, on un error menor al 4%, evaluado sobre los mismos datosusados en la identi� a ión.La �gura 5.13 muestra el omportamiento de la temperatura, en esta se observa omola evolu ión es muy similar, ex epto en los periodos de temperaturas altas, en los queel modelo predi e temperaturas inferiores a las que realmente se presentan. El error60

Cuadro 5.2: Errores obtenidos on el modelo de orden 2 para temperatura y humedad.por entual obtenido es inferior al 6%, evaluado sobre los mismos datos usados en laidenti� a ión.La �gura 5.14 muestra el omportamiento de la humedad relativa. Se observa omodurante el día (horas de baja humedad relativa) el modelo reprodu e ade uadamente el omportamiento de esta variable. En la no he, uando la humedad relativa es muy alta,al anzando in luso satura ión, la salida del modelo no orresponde a la evolu ión delsistema. La ausa de este fenómeno es la no linealidad existente en el pro eso uando seal anza satura ión, en ese momento la evapotranspira ión de las plantas esa y apare eun �ujo de vapor on signo negativo orrespondiente a la ondensa ión de agua en lasparedes y las hojas. El modelo lineal usado no es apaz de des ribir este fenómeno, portal razón la salida del modelo al anza valores de humedad de más del 100% durantealgunos intervalos.En la tabla 5.3.2 se resumen los resultados obtenidos on el modelo en ontrado. En éstase muestran los errores por entual, uadráti o medio y absoluto medio, para temper-atura y humedad.A ontinua ión se presenta el modelo obtenido para el sistema:24 x1(k + 1)

x2(k + 1)

35 =

24 0,9205 −0,1593

−0,016 0,7804

3524 x1(k)

x2(k)

35+

24 −0,1690 −174,6393 −0,0002

−0,2002 −245,2778 −0,0002

3526664 Text(k)

χext(k)

Qrad(k)

37775 (5.6)61

14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:000.006

0.007

0.008

0.009

0.01

0.011

0.012

0.013

0.014Resultado de la simulación de humedad

Tiempo (horas)

Hum

edad

(kg

/ m

3 )

Humedad simuladaHumedad medida

Figura 5.12: Humedad medida y simulada.

14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:000

5

10

15

20

25

30Resultado de la simulación de temperatura

Tiempo (horas)

Tem

pera

tura

(ºC

)

Temperatura simuladaTemperatura medida

Figura 5.13: Temperatura medida y simulada.62

14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:0020

30

40

50

60

70

80

90

100

110Resultado de la simulación de humedad relativa

Tiempo (horas)

Hum

edad

rel

ativ

a (%

)

Humedad simuladaHumedad medida

Figura 5.14: Humedad relativa medida y simulada.24 χ(k)

Ta(k)

35 =

24 −76,53µ −37,57µ

−0,6518 −0,3739

3524 x1(k)

x2(k)

35+

24 28,577µ 0,8852 0,1361µ

0,3238 −261,7601 71,38µ

3526664 Text(k)

χext(k)

Qrad(k)

37775 (5.7)5.4. Modelo de aire y follajeEn el segundo aso se in luyen omo salidas del sistema la temperatura y humedaddel aire a nivel de follaje, dado que se en ontraron gradientes importantes en el limade este nivel y el del aire a 2.5 m del suelo, omo se expuso en la se ión 5.1. Eneste aso no se toma omo salida el valor medio de las medidas obtenidas on los uatro sensores instalados, dado que el omportamiento del mi ro- lima dentro delfollaje depende onsiderablemente del tipo de siembra y la variedad de �or que seen uentre plantada. Por el ontrario, se toma un punto y se usa dire tamente en elpro eso de identi� a ión. La razón para usar este enfoque es que el modelo físi o indi aque la rela ión es prá ti amente unidire ional, el lima del aire sobre el ultivo afe ta63

el mi ro- lima a nivel de follaje, pero dado que el volumen de aire en errado es mu homenor al volumen total del invernadero, la in�uen ia en sentido inverso es muy pequeña,pudiendo ser despre iada.5.4.1. Estima ión del ordenAl obtener los valores singulares de la proye ión obli ua de las salidas futuras, a lo largode las entradas futuras, sobre las entradas y salidas pasadas, se obtiene un estimado delorden del sistema, que di�ere del resultado obtenido en la se ión anterior, en la que setomaron omo salidas solo la temperatura y humedad del aire.Al usar i = 5, omo tamaño de las matri es de Hankel, se obtiene un resultado queindi a que el orden del sistema es 3, ver �gura 5.15. El mismo resultado se obtiene alaumentar el tamaño de las matri es de Hankel on i = 10, ver �gura 5.16. Cuando se �jai = 15 el resultado ambia y se tiene uatro valores singulares de mayor valor, indi andoque el orden del sistema es 4, ver �gura 5.17. Si se ontinua aumentando el orden delos bloques de Hankel, los valores singulares presentan un de aimiento monótono queno permite obtener un estimado del orden.El modelo físi o y los experimentos preliminares indi an que el orden del sistema debeser 3. Dos e ua iones diferen iales des riben el omportamiento del aire sobre el ultivo, omo se expuso en la se ión anterior y se requiere una e ua ión adi ional para expli arel omportamiento de la temperatura del follaje, que se omporta omo una versión �l-trada de la temperatura del aire superior. La humedad, omo se indi ó en la se ión 5.1,presenta una dinámi a muy similar a la del aire superior, pero on un desplazamientovariable, que depende de la temperatura, el de� it de presión de vapor y la radia iónsolar de manera instantánea. Los resultados de la estima ión de orden indi an que elmodelo tiene un omportamiento muy similar al esperado on la formula ión a partirde prin ipios físi os.

64

0 5 10 15 20

0.5

1

1.5

2

2.5

3Valores singulares para j = 5

Orden

Val

ores

sin

gula

res

Figura 5.15: Autovalores de la proye ión obli ua on i = 5

0 5 10 15 20 25 30 35 40

0.5

1

1.5

2

2.5

3

3.5Valores singulares para j = 10

Orden

Val

ores

sin

gula

res

Figura 5.16: Autovalores de la proye ión obli ua on i = 10

65

0 10 20 30 40 50 60

0.5

1

1.5

2

2.5

3

3.5

Valores singulares para j = 15

Orden

Val

ores

sin

gula

res

Figura 5.17: Autovalores de la proye ión obli ua on i = 155.4.2. Algoritmo de identi� a iónAl igual que en el aso anterior, se usa el algoritmo robusto ombinado de De Moor yOvers hee [19℄, para diferentes tamaños de las matri es de Hankel y varios órdenes delsistema. A partir de lo ual, se sele iona el modelo que mejor des ribe los datos.Las �guras 5.18, 5.19, 5.20 y 5.21 muestran el omportamiento del error por entual delas uatro variables de salida que se tienen en este aso. Usando el riterio del peor aso(norma in�nita) se elige usar un sistema de orden 3 on i = 8, este modelo es el quemejor reprodu e el omportamiento de la temperatura, tanto del aire omo del follaje,redu iendo el error al 5.9%, manteniendo el error en la humedad menor al 4.5%.La tabla 5.3 resume los resultados obtenido on el modelo de orden 3 y i = 8. En estase apre ia que los errores de temperatura y humedad del aire aumentan levemente en ompara ión on el aso anterior en el que éstas eran las úni as salidas. Los resultadospara el lima a nivel de follaje son satisfa torios, ya que aumentando el orden del modeloen uno, es posible expli ar dos variables adi ionales on mayor pre isión, sin dañar losresultados anteriores.Las �guras 5.22, 5.23, 5.24, 5.25, 5.26 y 5.27 muestran la evolu ión de las variablessimuladas para los datos de identi� a ión, los resultados son similares a los obtenidos66

3 4 5 6 7 8 9 103

3.2

3.4

3.6

3.8

4

4.2

4.4

4.6

4.8

5

i

Err

or p

orce

ntua

l de

hum

edad

del

aire

Orden 1Orden 2Orden 3Orden 4

Figura 5.18: Error por entual de simula ión de la humedad del aire para diferentesórdenes y tamaños de las matri es de Hankel.

Cuadro 5.3: Errores obtenidos on el modelo de orden 3 para temperatura y humedadde aire y follaje. 67

3 4 5 6 7 8 9 105

5.5

6

6.5

7

7.5

8

8.5

9

i

Err

or p

orce

ntua

l de

tem

pera

tura

del

aire

Orden 1Orden 2Orden 3Orden 4

Figura 5.19: Error por entual de simula ión de la temperatura del aire para diferentesórdenes y tamaños de las matri es de Hankel.en el aso anterior.El siguiente es el modelo dinámi o obtenido para el sistema:26664 x1(k + 1)

x2(k + 1)

x3(k + 1)

37775 =

26664 0,9208 0,1108 −0,1312

−0,0105 0,8862 0,0240

−0,0143 0,0842 0,8535

3777526664 x1(k)

x2(k)

x3(k)

37775+

26664 −0,214 −189,9 −0,253m

0,102 87,5 0,0434m

−0,160 −174,9 −0,0989m

3777526664 Text(k)

χext(k)

Qrad(k)

37775(5.8)2666664 χaire(k)

Taire(k)

χfol(k)

Tfol(k)

3777775 =

2666664 −0,0495m −0,1520m −0,1740m

−0,5603 −0,0220 −0,3551

−0,0349 −0,0769 −0,3283

−0,5391 0,4153 −0,0995

377777526664 x1(k)

x2(k)

x3(k)

37775+

2666664 0,128m 0,825 0,078µ

0,372 −270,4 60,8µ

32,1µ 0,6274 0,152µ

71,5m −183,4 31µ

377777526664 Text(k)

χext(k)

Qrad(k)

37775(5.9)5.5. Valida iónPara veri� ar la validez de los modelos en ontrados, éstos se usan para simular el omportamiento del lima en el invernadero ante datos que no fueron usados en el68

3 4 5 6 7 8 9 103.8

3.9

4

4.1

4.2

4.3

4.4

4.5

i

Err

or p

orce

ntua

l de

hum

edad

del

folla

je

Orden 1Orden 2Orden 3Orden 4

Figura 5.20: Error por entual de simula ión de la humedad del follaje para diferentesórdenes y tamaños de las matri es de Hankel.

69

3 4 5 6 7 8 9 105

5.5

6

6.5

7

7.5

8

8.5

9

i

Err

or p

orce

ntua

l de

tem

pera

tura

del

folla

je

Orden 1Orden 2Orden 3Orden 4

Figura 5.21: Error por entual de simula ión de la temperatura del follaje para diferentesórdenes y tamaños de las matri es de Hankel.

14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:000.006

0.007

0.008

0.009

0.01

0.011

0.012

0.013

0.014Resultado de la simulación de humedad

Tiempo(horas)

Hum

edad

(kg/

m3 )

Humedad simuladaHumedad medida

Figura 5.22: Humedad del aire medida y simulada.70

14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:000

5

10

15

20

25

30Resultado de la simulación de temperatura

Tiempo (horas)

Tem

pera

tura

(ºC

)

Temperatura simuladaTemperatura medida

Figura 5.23: Temperatura del aire medida y simulada.

14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:0020

30

40

50

60

70

80

90

100

110Resultado de la simulación de humedad relativa

Tiempo(horas)

Hum

edad

rel

ativ

a(%

)

Humedad simuladaHumedad medida

Figura 5.24: Humedad relativa del aire medida y simulada.71

14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:000.006

0.007

0.008

0.009

0.01

0.011

0.012

0.013

0.014

0.015Resultado de la simulación de humedad

Tiempo(horas)

Hum

edad

(kg/

m3 )

Humedad simuladaHumedad medida

Figura 5.25: Humedad del follaje medida y simulada.

14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:000

5

10

15

20

25

30Resultado de la simulación de temperatura

Tiempo (horas)

Tem

pera

tura

(ºC

)

Temperatura simuladaTemperatura medida

Figura 5.26: Temperatura del follaje medida y simulada.72

14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:00 14:00 02:0020

30

40

50

60

70

80

90

100

110Resultado de la simulación de humedad relativa

Tiempo(horas)

Hum

edad

rel

ativ

a(%

)

Humedad simuladaHumedad medida

Figura 5.27: Humedad relativa del follaje medida y simulada.pro eso de identi� a ión. Para ello se registraron las variables ambientales del mismoinvernadero identi� ado durante tres días, posteriores al experimento de identi� a ión, omo se expli a en la se ión 4.3. Se uenta on 861 medidas que se obtuvieron durante72 horas.La tabla 5.4 muestra los errores del modelo ante datos nuevos. En esta se observa quela alidad de los resultados se redu e, la variable on mayor error es la temperatura delaire on un error por entual del 8.3%.Las �guras 5.28, 5.29,5.30, 5.31, 5.32 y 5.33 muestran el resultado de la simula ión paralos datos de valida ión. En ellas se apre ia omo la humedad obtenida en simula iónpresenta varia iones más lentas que las medidas y la temperatura estimada resulta sermás baja que la real. Además, se apre ian las desvia iones uando la humedad relativase a er a al 100%, ausada por la no linealidad de la satura ión del aire.73

Cuadro 5.4: Errores obtenidos on el modelo de�nitivo para datos nuevos.

08:00 20:00 08:00 20:00 08:00 20:000.006

0.007

0.008

0.009

0.01

0.011

0.012

0.013

0.014

0.015

0.016Humedad del aire, datos de validación.

Tiempo (horas)

Hum

edad

(kg

/m3 )

Humedad simuladaHumedad medida

Figura 5.28: Humedad del aire medida y simulada.74

08:00 20:00 08:00 20:00 08:00 20:005

10

15

20

25

30

35Temperatura del aire, datos de validación.

Tiempo(horas)

Tem

pera

tura

(ºC

)

Temperatura simuladaTemperatura medida

Figura 5.29: Temperatura del aire medida y simulada.

08:00 20:00 08:00 20:00 08:00 20:0020

30

40

50

60

70

80

90

100

110Humedad relativa del aire, datos de validación.

Tiempo (horas)

Hum

edad

rel

ativ

a (%

)

Humedad simuladaHumedad medida

Figura 5.30: Humedad relativa del aire medida y simulada.75

08:00 20:00 08:00 20:00 08:00 20:000.006

0.008

0.01

0.012

0.014

0.016

0.018

0.02Humedad del follaje, datos de validación.

Tiempo(horas)

Hum

edad

(kg

/m3 )

Humedad simuladaHumedad medida

Figura 5.31: Humedad del follaje medida y simulada.

08:00 20:00 08:00 20:00 08:00 20:005

10

15

20

25

30

35Temperatura del follaje, datos de validación.

Tiempo (horas)

Tem

pera

tura

(ºC

)

Temperatura simuladaTemperatura medida

Figura 5.32: Temperatura del follaje medida y simulada.76

08:00 20:00 08:00 20:00 08:00 20:0030

40

50

60

70

80

90

100

110Humedad relativa del follaje, datos de validación.

Tiempo(horas)

Hum

edad

rel

ativ

a (%

)

Humedad simuladaHumedad medida

Figura 5.33: Humedad relativa del follaje medida y simulada.5.6. Compara ión on otros resultadosComo último paso se realiza una ompara ión on los resultados obtenidos por otrosautores en estudios similares. Cabe anotar que durante el desarrollo de este proye tono se en ontraron reportes de identi� a ión de lima en invernaderos usando modeloslineales en variables de estados ni estru turas que onsideren de manera independienteel aire dentro del follaje y el aire sobre el ultivo.Una primera aproxima ión, expuesta por Wang et. al. [16℄ en 2000, presenta los resul-tados de un estudio en el que se estiman todos los parámetros de un modelo basadoen prin ipios físi os. En este usan un modelo de orden 6 y omo medida de desempeñoevalúan el error absoluto medio. Como resultados obtienen un error de 1,0 ◦C parala temperatura del aire y 5.3% para la humedad relativa. El error máximo absolutoobtenido es de 4,3 ◦C para la temperatura y 17.1% para la humedad relativa. Los re-sultados son omparables on los obtenidos en este estudio, on la diferen ia de que77

el modelo usado por Wang es onsiderablemente más omplejo y la estima ión de susparámetros requiere múltiples experimentos espe í� os y omplejos.Cunha presenta en 2003 una ompara ión entre varios métodos de identi� a ión apli a-bles a lima en invernaderos [10℄, entre ellos un modelo basado en prin ipios físi os, unmodelo lineal on estru tura ARX y un modelo no lineal basado en redes neuronales.Como medida de desempeño usa el error uadráti o medio y en uentra que el modelomás �el es el basado en prin ipios físi os on un error de 1,4 ◦C en temperatura y 4.94%en humedad relativa, En segundo lugar el modelo basado en redes neuronales presentaun error de 1,83 ◦C y 5.17%, por último el modelo lineal presenta un error de 2,0 ◦Cy 5.28%. Estos resultados son de menor alidad que los presentados en este trabajo yusan estru turas más omplejas en el aso de redes neuronales y prin ipios físi os.Herrero et. al. [7℄ presentan en 2003 un modelo basado en prin ipios físi os, uyos oe-� ientes son sintonizados por un algoritmo genéti o. Los resultados obtenidos, medidos omo el error absoluto medio de la temperatura y humedad relativa, son de 0,83 ◦C y5.28% respe tivamente. De nuevo estos resultados son muy er anos a los obtenidosen este trabajo y los algoritmos utilizados son onsiderablemente más omplejos om-puta ionalmente, ya que el algoritmo tarda aproximadamente 80 horas en sintonizar elmodelo.Salgado y Cunha [13℄ presentan en 2005 un modelo de aja gris en el que ombinanté ni as de lógi a difusa y prin ipios físi os para onstruir un modelo jerárqui o detemperatura y humedad del aire en un invernadero. En este trabajo obtienen un modelo on error uadráti o medio de 0,71 ◦C, 0,29g/m3 y 2.56%. Los errores obtenidos porSalgado y Cunha superan los resultados expuestos en este trabajo, usando un modeloaltamente no lineal omo lo son los sistemas difusos.La tabla 5.5 resume los resultados obtenidos por otros autores.78

Cuadro 5.5: Errores obtenidos por otros autores en estudios similares.79

Capítulo 6CONCLUSIONESEn este proye to se ha presentado el desarrollo de un modelo dinámi o del lima en uninvernadero ubi ado en la sabana de Bogotá. Los modelos obtenidos presentan erroresde simula ión en los mismos ordenes de magnitud que trabajos previos, on la diferen iade que la estru tura usada es lineal.Las medi iones realizadas y los modelos obtenidos permiten a�rmar que los invernaderosestudiados no aíslan el lima interior de las ondi iones externas. La resisten ia térmi aha ia el exterior es muy baja, por lo que el invernadero es in apaz de alma enar alorpara las no hes, bajando su temperatura hasta al anzar la exterior tan pronto esa laradia ión solar.La ondensa ión de agua sobre las hojas, ondi ión que se debe evitar para protegeral ultivo de la forma ión de hongos malignos, es una situa ión inevitable en los in-vernaderos en estudio. El bajo aislamiento on el exterior y al falta de sistemas de alefa ión no turna ha en que el aire de enfríe rápidamente al ano he er, redu iendosu apa idad de retener vapor de agua y al anzando satura ión, momento en el ual sepresenta ondensa ión en las hojas y las paredes.Para evitar la ondensa ión de agua en las hojas es ne esario modi� ar la estru turade los invernaderos , aumentando el aislamiento térmi o on el exterior sin sa ri� arventila ión natural, es de ir, si se ierran todas las ventanas laterales y enitales, seredu irá el inter ambio de aire on el exterior, lo ual in rementa la reten ión de alor,pero a la vez redu e el inter ambio de aire on el exterior, aumentando la humedad al80

interior debido a que las plantas están onstantemente transpirando. Esta alternativano pare e ser la ade uada.Otra op ión es alentar el aire interior, aumentando su apa idad de retener vaporde agua y evitando que se al an e la satura ión. Para que la energía invertida en el alentamiento sea baja es ne esario también in rementar el aislamiento térmi o on elexterior durante las horas de opera ión del sistema de alefa ión. A la vez, se debegarantizar una ade uada ventila ión durante el día, para evitar sobre argas de humedadpor transpira ión. El análisis lleva a la on lusión de que se requiere un invernadero on alefa ión y ventila ión arti� ial, la de isión de desarrollar e instalar un sistemade este tipo es más e onómi a que ientí� a, pues el osto de desarrollo y opera iónpuede ser prohibitivo en el ontexto olombiano.Aunque no se usaron dire tamente en el pro eso de identi� a ión, la onstru ión demodelos on base en prin ipios físi os es un paso fundamental en el pro eso de identi�- a ión, ya que ese paso permite entender la estru tura del sistema dinámi o en estudioy realizar interpreta iones de los resultados obtenidos en el pro eso de identi� a ión on mayor dominio del omportamiento físi o de la planta.La sele ión del onjunto de modelos posibles en un problema de identi� a ión es unpaso fundamental, y en este proye to la representa ión en variables de estados de sis-temas lineales resultó ser una parametriza ión ade uada. Además, los algoritmos porsubespa ios, son la herramienta ade uada para en ontrar los parámetros óptimos de unsistema dinámi o uando éste es multivariable.La modi� a ión propuesta de la estru tura del sistema, separando la dinámi a delaire atrapado por el follaje de la del aire superior, permitió in rementar la alidaddel sistema, sin sufrir un gran in remento en el osto omputa ional. La estru turautilizada no se ha en ontrado en la literatura utilizada y es un nuevo aporte que ofre eeste trabajo.Para mejorar la alidad del modelo es ne esario in luir más informa ión en el pro eso deidenti� a ión. En este proye to se dejaron de medir variables que in�uyen en el pro esodel lima omo la velo idad y dire ión del viento y la temperatura del suelo y subsuelo.La velo idad y dire ión del viento pueden llegar a ser parámetros ríti os en el modelo81

si son la prin ipal fuente de inter ambio de aire on el exterior. Esta situa ión se puededar si se aumenta el área de la ventila ión enital y se mejora el aislamiento térmi odel invernadero. Por ahora, el alor alma enado en el aire se pierde rápidamente ha iael exterior debido a la baja resisten ia térmi a de la ubierta.No fue posible medir la temperatura de suelo y subsuelo en este proye to debido ala aren ia de la instrumenta ión ade uada, sin embargo los resultados obtenidos sona eptables, lo ual indi a que estas variables in�uyen de manera marginal en la evolu ióndel lima, on base en los resultados obtenidos. Las razones para realizar esta a�rma iónson las mismas que se expusieron para el efe to del viento.Se debe realizar un estudio a fondo del mi ro lima dire tamente sobre las hojas, ya queeste es el lugar en el que re en los hongos que dañan el ultivo. Aunque en este proye tose tuvo un enfoque ma ro, en el que se estudió la evolu ión del lima promedio del inver-nadero, desde un enfoque �tosanitario no es relevante que el aire sobre el follaje al an esatura ión si se garantiza que no se presenta ondensa ión sobre las hojas. Este nuevoenfoque abre una nueva linea de investiga ión en la que se debe empezar por de�nirlas variables que des riben las ondi iones ambientales en la super� ie de una hoja,diseñar o sele ionar la instrumenta ión ade uada para medir las variables ne esarias,plantear modelos basados en prin ipios físi os y realizar experimentos de identi� a iónque permitan obtener modelos ade uados del pro eso. El �n último del estudio debeser el planteamiento de alternativas en el manejo del lima que sean e onómi amenteviables y se enfoque en evitar al ondensa ión de vapor de agua sobre las hojas.

82

Bibliografía[1℄ Bot, G. P. A., (1983). Greenhouse limate: from physi al pro ess to a dynami almodel. Tesis do toral, Agri ultural University of Wageningen, Holanda, 240 pági-nas.[2℄ Bot, G. P. A., (2001). Developments in indoor sustainable plant produ tion withemphasis on energy saving. Computers and Ele troni s in Agri ulture. Vol. 30,151-165.[3℄ Hanan, Joe J., (1998). Greenhouses: advan ed te hnology for proteted horti ulture.CRC Press. Bo a Raton, Florida.[4℄ El Ghoumari M. Y., (2003). Optimiza ión de la produ ión de un invernaderomediante ontrol predi tivo no lineal. Tesis do toral, Universidad Autónoma deBar elona, 181 páginas.[5℄ Kozai, T., Sase, S., (1978). A simulation of natural ventilation for multi-span green-house. A ta Horti ulturae. 87, 39-49.[6℄ Tap R. F., Van Willigenburg L.G., Van Straten G., (1996). Experimental resultsof re eding horizon optimal ontrol of greenhouse limate. Pro . Math. and Contr.Appli . in Agri . and Hort. 406, 229-238.[7℄ Herrero J. M., Martínez M., Blas o X., Ramos C. y Rodríguez S., (2003). Identi-� a ión paramétri a del modelo no lineal de un invernadero mediante AlgoritmosGenéti os. XXIV Jornadas de Automáti a, León, España, Sept. 2003.

83

[8℄ Körner O., (2003). Crop based limate regimes for energy saving in greenhouse ultivation. Tesis do toral, Agri ultural University of Wageningen, Holanda, 240páginas.[9℄ Körner O., Challa H., (2004). Temperature integration and pro ess-based humidity ontrol in hrysanthemum. Computers and ele troni s in agri ulture. Vol. 43, pag.1-21.[10℄ Cunha J. B., (2003). Greenhouse limate models: An overview, EFITA Conferen e,Budapest, Hungría, Jul. 2003.[11℄ Cunha J. B., Couto C. y Ruano A. E., (1997). Real-time parameter estimation ofdynami temperature models for greeenhouse environmetal ontrol. Control Engi-neering Pra ti e. Vol. 5, No. 10, pag. 1473-1481.[12℄ Cunha J. B. y De Moura J. P., (2003). Optimal management of greenhouse envi-ronmets. EFITA Conferen e, Budapest, Hungría, Jul. 2003.[13℄ Salgado P. y Cunha J. B., (2005). Greenhouse limate hieran hi al fuzzy modelling.Control Engineering Pra ti e. Vol. 13, pag. 613-628.[14℄ Young P., Chotai A. (2001). Data-Based me hanisti modelling, more asting and ontrol. Control systems magazine, IEEE, Vol. 21, No. 5, pag. 14-27.[15℄ Young P., Lees M., Chotai A., Ty h W. y Chalabi Z. (1994). Modelling and PIP ontrol of a glasshouse mi ro- limate. Control Engineering Pra ti e. Vol. 2, No. 4,pag. 1591-604.[16℄ Wang S., Boulard T.,(2000). Predi ting the mi ro limate in a naturally ventilatedplasti house in a mediterranean limate. J. Agri . engng Res. vol. 75, No.1, pag.27-38.[17℄ Kania S., Gia omelli G., (2002). Solar radiation avaibility for plant growth in Ari-zona ontrolled environment agri ulture systems, Pro . 30th National agri ulturalplasti s ongress, San Diego, USA, Feb. 2002.84

[18℄ Ljung L., (1999). System Identi� ation − Theory for the user, 2 ed., PTR Prenti eHall. Upper Saddle River.[19℄ Over hee V., De Mooor B., (1996). Subspa e identi� ation for linear systems, Kluw-er A ademi Publishers, Dordre ht, Holanda.[20℄ Leyva P., Editor, (2001). El medio ambiente en Colombia, Instituto de hidrología,meteorología y estudios ambientales (IDEAM), Bogotá, Colombia.[21℄ Trigui M., Barrington S., Gauthier L., (2001). A strategy for greenhouse limate ontrol, part I: Model development, J. Agri . engng Res. vol. 78, No.4, pag. 407-413.[22℄ Piñon S., Peña M., Ku hen B., (2002). Control predi tivo on restri iones para el lima de un inverná ulo, Dyna, vol. 135, pag. 11-25.[23℄ WoodwardF., Sheehy J., (1983). Prin iples and measurements in environmentalboilogy. Butterworths, Londres.[24℄ Franklin G., Powell J., Workman M., (1998). Digital ontrol of dynami systems,3 ed., Addison Wesley Longman, In ., Menlo Park, California.[25℄ Sinha N., Kuszta B., (1983). Modelling and identi� ation of dynami systems, VanNostrand Reinhold, New York.

85