Desarrollo e Implementaci´on de Algoritmos de Control para...

Desarrollo e Implementacion de Algoritmos de

Control para Turbinas de Gas Combustible

Br. Gamboa Arguello Genderson Enrique

Profesor Guıa: Prof. Addison Rıos Bolıvar

Proyecto de Grado presentado ante la ilustre Universidad de Los Andes

como requisito final para optar al Tıtulo de Ingeniero de Sistemas

Merida, Venezuela

Mayo 2007

Universidad de Los Andes c© Derechos Reservados 2007

CALIFICACION DE

PROYECTO DE GRADO

Los suscritos, Miembros del Jurado designado para evaluar el Proyecto

de Grado titulado “Desarrollo e Implementacion de Algoritmos

de Control para Turbinas de Gas Combustible” realizado por

Br. Gamboa Arguello Genderson Enrique, C.I. N◦ 13.356.847, presentado

como requisito parcial para la obtencion del grado de Ingeniero de Sistemas, han

resuelto, de acuerdo a lo pautado en el Artıculo 13 del Reglamento de Proyecto de

Grado, colocar la calificacion de

CALIFICACION:

Fecha: Mayo 2007

OBSERVACIONES:

Profesor Guıa:Prof. Addison Rıos Bolıvar

Jurado:Prof. Claudia S. Gomez Quintero

Prof. Juan Cardillo Albarran

Dedicado a Carmen Arguello de Gamboa,

Adelso Gamboa y todas aquellas personas

que hicieron posible el terminar este

largo camino

Indice

Indice de Figuras viii

Agradecimientos xii

Resumen xiii

Introduccion xv

1 Un Modelo de Turbina a Gas Combustible 1

1.1 Introduccion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1

1.2 Definicion de Turbina a Gas Combustible . . . . . . . . . . . . . . . . . 2

1.2.1 Antecedentes historicos . . . . . . . . . . . . . . . . . . . . . . . 3

1.2.2 Componentes Principales de una Turbina a Gas . . . . . . . . . 4

1.2.3 Funcionamiento de una Turbina a Gas Combustible . . . . . . . 4

1.2.4 Campos de Aplicacion . . . . . . . . . . . . . . . . . . . . . . . 5

1.3 Modelo de Turbina Prototipo HTDU . . . . . . . . . . . . . . . . . . . 5

1.3.1 Apreciacion General del Modelo . . . . . . . . . . . . . . . . . . 6

1.3.2 Diagramas de bloque en MatlabTM − Simulink . . . . . . . . . 6

1.3.3 Turbina . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7

1.3.4 Variables del Modelo y Unidades . . . . . . . . . . . . . . . . . 12

1.3.5 Ecuaciones Matematicas del Modelo . . . . . . . . . . . . . . . 14

2 Enfoques de Control Aplicados al Problema 20

2.1 Introduccion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

2.2 Definicion de un Sistema de Control . . . . . . . . . . . . . . . . . . . . 23

iv

2.2.1 Conceptos de Control . . . . . . . . . . . . . . . . . . . . . . . . 23

2.2.2 Clasificacion de sistema de control . . . . . . . . . . . . . . . . . 24

2.3 Fases para el Diseno de Sistemas de Control . . . . . . . . . . . . . . . 24

2.3.1 Modelado fısico . . . . . . . . . . . . . . . . . . . . . . . . . . . 24

2.3.2 Identificacion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25

2.3.3 Diseno de controladores . . . . . . . . . . . . . . . . . . . . . . 27

2.3.4 Sintonizacion y Validacion del diseno . . . . . . . . . . . . . . . 27

2.3.5 Implantacion . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27

2.3.6 Pruebas y Depuracion . . . . . . . . . . . . . . . . . . . . . . . 28

2.4 Criterios para el Analisis de Desempeno de un Sistema . . . . . . . . . 28

2.4.1 Respuesta transitoria . . . . . . . . . . . . . . . . . . . . . . . . 28

2.4.2 Respuesta Estacionaria . . . . . . . . . . . . . . . . . . . . . . . 29

2.5 Enfoque de Espacio Estado . . . . . . . . . . . . . . . . . . . . . . . . 29

2.5.1 Construccion del Vector de Estados . . . . . . . . . . . . . . . . 31

2.5.2 Observador de Estados . . . . . . . . . . . . . . . . . . . . . . . 32

2.5.3 Diseno de Sistemas por Seguimiento . . . . . . . . . . . . . . . . 35

2.5.4 Diseno de Control Optimo Cuadratico . . . . . . . . . . . . . . 38

2.6 Enfoque en el Dominio de Frecuencia . . . . . . . . . . . . . . . . . . . 41

2.6.1 Estabilidad Relativa . . . . . . . . . . . . . . . . . . . . . . . . 41

2.6.2 Sensibilidad de los Parametros . . . . . . . . . . . . . . . . . . . 41

2.6.3 Rechazo de Perturbaciones . . . . . . . . . . . . . . . . . . . . . 42

2.7 Tecnicas Clasicas de Control . . . . . . . . . . . . . . . . . . . . . . . . 42

2.7.1 Controlador por Adelanto de Fase . . . . . . . . . . . . . . . . . 42

2.7.2 Controlador por Atraso de Fase . . . . . . . . . . . . . . . . . . 46

2.7.3 Controlador PD (Proporcional Derivativo) . . . . . . . . . . . . 48

2.7.4 Controlador PI (Proporcional Integral) . . . . . . . . . . . . . . 50

2.7.5 Diseno del controlador PID . . . . . . . . . . . . . . . . . . . . 52

2.8 Esquema Anti-Rebote (Anti-Windup) . . . . . . . . . . . . . . . . . . . 54

2.8.1 Metodologıa General para el Diseno del esquema Antirebote

(Anti-windup) . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54

3 Identificacion de un Modelo de Turbina a Gas Combustible 56

3.1 Introduccion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56

3.2 Identificacion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57

3.2.1 Obtencion del Modelo Parametrico . . . . . . . . . . . . . . . . 58

3.2.2 Validacion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61

4 Diseno y Analisis de Resultados de las Estrategias de Control

Aplicadas 65

4.1 Introduccion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 65

4.2 Diseno de un Observador de Orden Completo . . . . . . . . . . . . . . 66

4.3 Diseno de Control por Seguimiento . . . . . . . . . . . . . . . . . . . . 69

4.4 Diseno del Control Optimo Cuadratico . . . . . . . . . . . . . . . . . . 73

4.5 Diseno de Controladores Clasicos . . . . . . . . . . . . . . . . . . . . . 77

4.5.1 Diseno de Controlador por Atraso de Fase . . . . . . . . . . . . 77

4.5.2 Diseno de un Controlador Proporcional Integral (PI) . . . . . . 81

4.5.3 Diseno de Controlador Proporcional Integral (PID) . . . . . . . 84

4.6 Diseno de Esquema Anti-Rebote (Anti-Windup) . . . . . . . . . . . . . 89

4.7 Analisis y Comparacion entre los Enfoques previamente mencionados . 93

5 Implantacion de Algoritmos de Control Bajo PLC. 97

5.1 Introduccion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 97

5.2 Transmision a traves de MatlabTM usando comunicacion serial RS-232

y/o TCP/IP. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98

5.2.1 Xpc . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99

5.2.2 Conexion de la Red para un Xpc . . . . . . . . . . . . . . . . . 99

5.2.3 Propiedades del Ambiente Xpc. . . . . . . . . . . . . . . . . . . 100

5.3 Plataformas Logix 5000 . . . . . . . . . . . . . . . . . . . . . . . . . . . 102

5.3.1 Caracterısticas de la Plataforma Logix 5000 . . . . . . . . . . . 102

5.3.2 Lenguajes de Programacion . . . . . . . . . . . . . . . . . . . . 103

5.4 Implatacion Digital de Algoritmos de Control . . . . . . . . . . . . . . 104

6 Conclusiones y Recomendaciones 106

Bibliografıa 108

A Codigos y Pantallas Usadas en el Editor de MatlabTM 110

A.1 Codigo Fuente de la Turbina de Gas Combustible HTDU . . . . . . . . 110

A.2 Codigo para generar las graficas de las diversas variables que maneja el

modelo de la Turbina de Gas Combustible HTDU . . . . . . . . . . . . 126

A.3 Codigo para generar el informe de los puntos de operacion de las diversas

variables que maneja el modelo de la Turbina de Gas Combustible HTDU 131

A.4 Codigo para la validacion mediante el error del modelo identificado de

la Turbina de Gas Combustible HTDU . . . . . . . . . . . . . . . . . . 136

A.5 Turbina HTDU Control por Seguimiento . . . . . . . . . . . . . . . . . 139

A.6 Turbina HTDU Control Optimo . . . . . . . . . . . . . . . . . . . . . . 140

A.7 Turbina HTDU Control por Atraso de Fase . . . . . . . . . . . . . . . 140

A.8 Turbina HTDU Control Proporcional Integral (PI) . . . . . . . . . . . 141

A.9 Turbina HTDU Control Proporcional Integral Derivativo (PID) . . . . 141

Indice de Figuras

1.1 Modelo de una turbina . . . . . . . . . . . . . . . . . . . . . . . . . . . 3

1.2 Modelo en Simulink de una turbina, (proporcionado por el Dr. Silvio

Simani) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7

1.3 Modelo en Simulink del Sub-sistema Turbina HTDU . . . . . . . . . . 11

1.4 Esquema del Modelo de la Turbina HTDU . . . . . . . . . . . . . . . . 13

1.5 Esquema del Modelo en Simulink HTDU . . . . . . . . . . . . . . . . . 14

2.1 Curva de respuesta de sistema de orden 2, ante una entrada escalon

unitario . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29

2.2 Diagrama de bloques del esquema de control realimentado con un

observador . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34

2.3 (a) Graficas de las curvas e(t) y u(t) que muestran una senal de control

diferente de cero cuando la senal de error es cero(control integral); (b)

graficas de las curvas e(t) y u(t) que muestran una senal de control de

cero cuando la senal de error es cero (control proporcional). . . . . . . . 35

2.4 Sistema de seguimiento de tipo 1 con un observador de estado . . . . . 38

2.5 (a) Aporte de fase positiva: Adelanto de fase. (b) Aporte de fase

negativa: Atraso de fase . . . . . . . . . . . . . . . . . . . . . . . . . . 43

2.6 Esquema de control Anti-windup . . . . . . . . . . . . . . . . . . . . . 55

3.1 Un modelo de turbina HTDU en laso abierto sin compensar . . . . . . 57

3.2 Senales entrada/salida del modelo sin compensar . . . . . . . . . . . . 58

3.3 Herramienta System Identification Toolbox de MatlabTM utilizada . . . 59

3.4 Metodos Parametricos usando el System Identification Toolbox MatlabTM 61

viii

3.5 Resultados de los Metodos Parametricos del System Identification

Toolbox MatlabTM . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 63

3.6 Error Absoluto entre el modelo real y el identificado . . . . . . . . . . . 63

3.7 Error Porcentual entre el modelo real y el identificado . . . . . . . . . . 64

4.1 Diagrama analogico del espacio estado identificado en MatlabTM −Simulink . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 67

4.2 Diagrama analogico del observador de orden completo en MatlabTM −Simulink . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69

4.3 Digrama del sistema de seguimiento en MatlabTM − Simulink . . . . . 71

4.4 Entrada/Salida del sistema compensado (control por seguimiento),

referencia: 1466,1 [rad/s] . . . . . . . . . . . . . . . . . . . . . . . . . . 72


referencia: 1660,5 [rad/s] . . . . . . . . . . . . . . . . . . . . . . . . . . 72


referencia: 1780,2 [rad/s] . . . . . . . . . . . . . . . . . . . . . . . . . . 73

4.7 Entrada/Salida del sistema compensado (control por optimo),

referencia: 1466,1 [rad/s] . . . . . . . . . . . . . . . . . . . . . . . . . . 75


referencia: 1660,9 [rad/s] . . . . . . . . . . . . . . . . . . . . . . . . . . 76


referencia: 1780,2 [rad/s] . . . . . . . . . . . . . . . . . . . . . . . . . . 76

4.10 Diagrama de Bode: a la izquierda la planta sin compensar mientras que

a la derecha esta compensada con Kc1 . . . . . . . . . . . . . . . . . . . 78

4.11 Entrada/Salida del sistema compensado (control por atraso de fase),

referencia: 1047,2 [rad/g] . . . . . . . . . . . . . . . . . . . . . . . . . . 80


referencia: 1256,6 [rad/s] . . . . . . . . . . . . . . . . . . . . . . . . . . 80


referencia: 1455,6 [rad/s] . . . . . . . . . . . . . . . . . . . . . . . . . . 81

4.14 Diagrama de bode de la planta sin compensar . . . . . . . . . . . . . . 82

4.15 Entrada/Salida del sistema compensado (control proporcional integral

(PI)), referencia: 1466,1 [rad/s] . . . . . . . . . . . . . . . . . . . . . . 83


(PI)), referencia: 1660,5 [rad/s] . . . . . . . . . . . . . . . . . . . . . . 83


(PI)), referencia: 1780,2 [rad/seg] . . . . . . . . . . . . . . . . . . . . . 84

4.18 Diagrama de Bode: a la izquierda, la planta sin compensar; a la derecha,

esta la planta compensada con Kp . . . . . . . . . . . . . . . . . . . . . 85

4.19 Diagrama de Bode compensado con Kp . . . . . . . . . . . . . . . . . . 86


derivativo(PID)), referencia: 1466,1 [rad/s] . . . . . . . . . . . . . . . . 87





4.23 Entrada/Salida del sistema compensado (Gobernador sin antirebote) y

la plata sin compensador . . . . . . . . . . . . . . . . . . . . . . . . . . 89

4.24 Entrada/Salida del sistema compensado (Gobernador con saturacion

[0.1,-0.1] sin Kcr) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90


[0.2,-0.2] sin Kcr) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 90


[0.1,-0.1] con Kcr=2000) . . . . . . . . . . . . . . . . . . . . . . . . . . 91


[0.1855,-0.1855] con Kcr=8) . . . . . . . . . . . . . . . . . . . . . . . . 91


[0.19,-0.19] con Kcr=8) . . . . . . . . . . . . . . . . . . . . . . . . . . . 92

4.29 Respuesta del sistema, comparacion entre controladores con referencia

= 1413,7 rads

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 93


= 1445,1 rads

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94


= 1660,9 rads

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 94


= 1780,2 rads

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95


= 1885 rads

. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 95

5.1 Xpc TargetBoxTM . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 99

5.2 Conexion entre PC, a traves de comunicion serial RS232 y TCPIP . . . 100

5.3 Ventana con propiedades de configuracion . . . . . . . . . . . . . . . . 101

5.4 Control Logix 5555 de Allen-Bradley . . . . . . . . . . . . . . . . . . . 102

5.5 Conexion con el Control Logix . . . . . . . . . . . . . . . . . . . . . . . 105

A.1 Pantalla secundaria de la turbina HTDU, estructura del Control por

Seguimiento . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 139

A.2 Pantalla secundaria de la turbina HTDU, estructura del Control Optimo 140

A.3 Pantalla secundaria de la turbina HTDU, estructura del Control Atraso

de Fase . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 140

A.4 Pantalla secundaria de la turbina HTDU, estructura del Control

Proporcional Integral (PI) . . . . . . . . . . . . . . . . . . . . . . . . . 141

A.5 Pantalla secundaria de la turbina HTDU, estructura del Control

Proporcional Integral Derivativo (PID) . . . . . . . . . . . . . . . . . . 141

Agradecimientos

A Dios Todopoderoso.

A la Ilustre Universidad de los Andes, por abrir sus puertas y ofrecerme la

oportunidad para mi formacion profesional y el financiamiento economico a traves del

CDCHT con el proyecto No I-1001-06-02-f.

Al Profesor Addison Rıos y la Profesora Claudia Gomez, por su incondicional

apoyo y guıa que hicieron posible la culminacion de este Proyecto.

A mi hermano Frank Gamboa y el Sr. Jesus Omar Leon por su continua ayuda

desde el inicio hasta el final de la carrera.

A la Sra. Juana Hernandez de Di’nucci que en vida me recibio como un Hijo,

deseo expresarle mi gratitud y enorme aprecio.

A mi familia y todas aquellas personas que de alguna manera brindaron su

colaboracion y ayudaron a hacer realidad esta meta.

xii

Resumen

Los sistemas de control se encuentran en la mayorıa de los sectores industriales,

tales como control de la calidad de los productos manufacturados, lıneas de ensamble

automatico, control de maquinas-herramientas, robotica y muchos otros.

El proposito del control automatico, a traves de sus metodos, consiste en

mantener en un determinado valor de operacion las variables del proceso, tales como:

temperaturas, presiones, flujos y velocidad. Los procesos son de naturaleza dinamica

y no lineales, en ellos siempre ocurren cambios y si no emprenden las acciones

pertinentes, las variables importantes del proceso, no cumpliran con las condiciones

de diseno, basadas en un estudio completo.

El estudio para el desarrollo de un modelo de control esta fundamentado en

esquemas existentes, para la construccion de controladores. La finalidad de estas

tecnicas es ser implementada en un sistema de turbinas a gas combustible, con el

proposito de alcanzar los objetivos especıficos del proyecto.

Para ello se disena una ley de control convencional de acuerdo a las

especificaciones de funcionamiento propuestas y las restricciones de los actuadores

presentes en el proceso. Su ejecucion es programada a traves de la herramienta

matematica Simulink−MatlabTM destinado a generar el codigo en texto estructurado

correspondiente a los enfoques de control.

Finalizando se explicara como podrıan ser implantados los controladores en el

Computador Logico Programable (PLC) bajo la Plataforma-Logix, cuyos elementos

xiii

y variables pueden ser modificados automaticamente a traves de un sistema de

supervision basado en esta tecnologıa.

Descriptores: Turbina a gas combustible, Enfoques de Control, Plataformas

Tecnologicas, optimizacion, modelos matematicos.

Cota: no tiene

Introduccion

Los enfoques de control, han experimentado diversos cambios a lo largo de la historia.

Desde tiempos antiguos el hombre ha necesitado crear mecanismos que le permitan

resolver problemas especıficos. En la actualidad se han desarrollado numerosas tecnicas

y principios a traves de los cuales se logra un desempeno optimo de determinados

procesos.

La principal inquietud del ser humano ha sido el diseno de sistemas capaces

de realizar sus tareas, de forma independiente, bajo un ambiente seguro, inclusive en

condiciones no aptas para el desenvolvimiento humano, donde la vida de las personas

no corra ningun tipo de riesgo, razon por la cual se deben implementar dichos metodos

en los procesos industriales.

En este proyecto, el proceso industrial de estudio tiene que ver con las turbinas

a gas combustible cuyo funcionamiento consiste en un motor de flujo continuo que se

caracteriza por presentar una baja relacion de peso-potencia y velocidad de giro muy

elevada. Los orıgenes de las turbinas a gas combustible se remontan a mucho antes

de que el desarrollo tecnologico, tanto de materiales como de procesos industriales,

permitiera su correcta implantacion.

Las turbinas de gas orientadas a la propulsion a reaccion se usan en la gran

mayorıa de aviones comerciales y militares, mientras las turbinas de gas orientadas

a la generacion de trabajo en un eje se usan tambien en buques, trenes, autobuses,

camiones y vehıculos, pero la funcion prioritaria es la produccion de energıa electrica.

xv

En general, los elementos basicos de una turbina estan conformados por:

compresor, camara de combustion, turbina de potencia y carga. La variable a ser

controlada es la velocidad de la turbina, manejada en unidades de radianes por

segundo o revoluciones por minuto; la variable manipulada es el flujo de combustible

y finalmente las variables restantes son las medidas.

Desde el punto de vista del diseno de sistemas de control para procesos

industriales, es necesario recurrir a enfoques, los cuales, mediante pruebas y

comparaciones, permiten alcanzar funcionamientos optimos. En particular se busca

una solucion de facil manejo y flexible a modificaciones.

Para implementar estos enfoques en la actualidad existe una variedad de

tecnologıas disponibles. Como es el caso de la Rockwell Automation que mediante un

PLC, disponible en el laboratorio de control de la Facultad de Ingenierıa, Escuela-

Sistemas, que permite la implantacion automatizada de leyes de control.

En este sentido, se debe alcanzar un objetivo principal como es la identificacion,

evaluacion, analisis e implementacion de polıticas de control para un modelo de

turbina prototipo High Temperature Demonstration Unit (HTDU) a gas combustible.

Especıficamente se pretende la Identificacion de las tecnologıas y metodologıas de

desarrollo disponibles para el control de turbinas a gas combustible que permitan

comparaciones con cualquier aplicacion de control existente en la industria y explicar

como pudiese ser la implantacion de dichos enfoques, en un PLC.

De esta manera aprovechar las potencialidades matematicas para la construccion

de algoritmos de control y de esta forma contribuir con una generacion de mecanismos

para la regulacion de turbinas a gas combustible bajo plataformas de automatizacion.

Capıtulo 1

Un Modelo de Turbina a Gas

Combustible

Este capıtulo concierne la teorıa de las turbinas de gas combustible, explicando todo

lo relacionado con la misma, desde ¿que es una turbina? hasta su funcionamiento,

descripcion, y campos de aplicacion. Tambien se explica en plenitud el modelo en

MatlabTM − Simulink de la turbina prototipo High Temperature Demonstration Unit

(HTDU) utilizado durante este estudio.

1.1 Introduccion

Las turbinas a gas, son las mas recientes. Si bien hubo intentos de fabricarlas desde

inicios del siglo pasado, el primer ensayo exitoso data de 1937. Los mismos difieren

de otras en el sentido de que la combustion se realiza dentro de la maquinas y, por lo

tanto, el fluido de trabajo son gases de combustion (de allı su nombre).

En los proximos parrafos se explora algunos aspectos fundamentales acerca de

la Turbina a Gas y su modelo matematico en MatlabTM − Simulink. Estos son:

• Definicion de la Turbina a Gas: descripcion general sobre la produccion de

1.2 Definicion de Turbina a Gas Combustible 2

gases para generar la energıa mecanica para luego transformarla en electrica.

• Aspectos historicos: algunos detalles sobre la evolucion historica de la turbina

a gas. Implementacion en la aviacion y para la produccion de la energıa electrica.

• Componentes principales: Compresor, turbinas y sistema de combustion.

• Funcionamiento de la Turbina a Gas: ciclo empleado y como se realiza el

ciclo de operacion.

• Campos de aplicacion: Usos orientados a la generacion de trabajo.

1.2 Definicion de Turbina a Gas Combustible

Una turbina a gas combustible es un motor rotativo de flujo continuo que se caracteriza

por presentar una baja relacion de peso potencia y velocidad de giro muy elevada.

La elevada velocidad de giro, que en funcion del tamano puede llegar a alcanzar

una velocidad de 40000 revoluciones por minuto, orienta su uso a una unidad de

generacion elevada de entalpıa que puede utilizarse para propulsion a reaccion o puede

ser encargada de accionar una turbina de potencia acoplada a un eje, en la que puede

acoplarse cualquier tipo de carga (Prado y Lemus 1987).

La turbina de gas esta conformada por dos elementos principales:

• El generador de gases.

• La unidad generadora de potencia.

El generador de gases esta formado a su vez por uno o varios compresores, la camara

de combustion, donde se mezclara el combustible con el aire y donde tendra lugar la

combustion y finalmente la o las turbinas de expansion de gases, que en este caso solo

obtendran la potencia necesaria para mover los compresores.


La unidad generadora de potencia es donde se obtendra la potencia util de la

maquina y dependiendo de la aplicacion sera otra turbina de expansion de gases o

bien una tobera de propulsion.

Figura 1.1: Modelo de una turbina

1.2.1 Antecedentes historicos

Los orıgenes de las turbinas a gas combustible se remontan a muchos antes del

desarrollo tecnologico. Las limitaciones esenciales provienen de las altas temperaturas

de trabajo de los materiales, el correcto equilibrio y la articulacion del rotor, por el

elevado regimen de giro del mismo (Ediciones UPC 2007).

La primera patente de una turbina de gas se obtiene en 1871 por John Barber,

pero no es hasta 1900 que se construye la primera turbina a gas que funciono realmente

y que fue disenada en Francia por Stolze, aunque los resultados obtenidos fueron

decepcionantes. El aporte de diversos cientıficos como Frank Whittle, permitio que en

1939 se construyese el primer avion del mundo propulsado por una turbina de gas, el

avion aleman He 178, ano en el que tambien se implemento la primera turbina a gas

para la produccion de energıa electrica.


1.2.2 Componentes Principales de una Turbina a Gas

Los principales accesorios de la turbina a gas son los componentes del rotor, formados

por el compresor de flujo axial y las ruedas de turbina, y los componentes estacionarios,

formados por la carcasa de la turbina y el combustor (Prado y Lemus 1987).

1. Compresor: Su funcion es la de proveer aire de alta presion a las camaras de

combustion para producir gases calientes, necesarios para hacer girar la turbina.

Tambien suministra aire de enfriamiento para las toberas, las ruedas de la

turbinas, las piezas de transmision y otras porciones de espacio recorrido por

los gases calientes.

2. Turbina: los alabes de turbina convierten la energıa cinetica de los gases de escape

en energıa mecanica y rotacional util, la cual produce la fuerza necesaria para

accionar el compresor de flujo axial y cumplir con los requisitos de carga.

3. Sistema de combustion: Provee energıa mediante el calor necesario para el ciclo

de turbina. Esto se logra quemando combustible en el flujo de aire del compresor

y diluyendo los productos de combustion con el aire excedente para lograr la

temperatura deseada.

1.2.3 Funcionamiento de una Turbina a Gas Combustible

El principio de funcionamiento de una turbina a gas combustible consiste en un

compresor dinamico que suministra aire a una camara de combustion en donde

se quema combustible, con exceso de aire, a presion constante. A este proceso lo

llamamos ciclo simple, solo significa que los productos de la combustion se mezclan

con un exceso de aire para producir gas con energıa a una temperatura lo bastante

baja para el tipo de materiales usados.

El gas energizado se expande en una turbina que impulsa el compresor de aire y

produce potencia adicional como salida mecanica. Como ultimo paso, los productos

1.3 Modelo de Turbina Prototipo HTDU 5

de la combustion se descargan en la atmosfera.

Una turbina de gas, en general, no funciona a plena potencia, por lo que la

variacion del rendimiento es un factor de gran importancia economica; para adaptar

una turbina de gas a una potencia dada se puede actuar sobre la temperatura de

admision de los gases en la turbina φ, o sobre el numero de r.p.m. del compresor n, lo

que implica una modificacion del gasto G y de la relacion α de compresion, que solo

es posible cuando el receptor pueda admitir variaciones de velocidad de gran amplitud

(Prado y Lemus 1987).

El punto de funcionamiento viene definido por la interseccion de las curvas

caracterısticas de cada maquina.

1.2.4 Campos de Aplicacion

Las turbinas de gas orientadas a la propulsion a reaccion se implementan en la gran

mayorıa de aviones comerciales y militares, mientras que las turbinas de gas orientadas

a la generacion de trabajo en un eje tambien se han utilizado en buques, trenes,

autobuses, camiones y vehıculos, pero la funcion prioritaria es la produccion de energıa

electrica, para cubrir las demandas de la dicha fuente energetica.

1.3 Modelo de Turbina Prototipo HTDU

Este es un modelo desarrollado por Alstrom Corporation el cual fue extraıdo de una

turbina prototipo HTDU, con todos los procesos internos modelados, mediante el uso

de tecnologıa de punta describiendo a plenitud su funcionamiento (Silvio Simani y

Cesare Fantuzzi b. 2006).


1.3.1 Apreciacion General del Modelo

El flujo de aire entra por el conducto hacia el compresor, este es sometido a grandes

presiones de aire, donde el compresor lo calienta y envıa hacia la camaras de

combustion y estas se propagan por la turbina de fase del compresor, donde una

valvula tipo Mariposa proporciona un medio de generar una presion posterior en la

turbina del compresor (no hay potencia presente en el modelo de la turbina).

El aire refrigerante manado por el paso de la valvula del compresor enfrıa

el estator y el rotor de la turbina. Bajo condiciones de estado estable, la potencia

generada por la turbina del compresor esta en equilibrio, de manera que el compresor

la asimila, junto a sus perdidas. Los elementos recesivos representan perdidas de

presion en los componentes, los elementos acumuladores representan volumenes en los

componentes.

1.3.2 Diagramas de bloque en MatlabTM − Simulink

Los componentes entre el conducto de la entrada y la valvula de paso de la turbina del

compresor, incluso las perdidas mecanicas, se representan por la etiqueta ’turbina de

HTDU ’.

El bloque de entrada se representa por la temperatura y presion ambiental

los cuales son ′t1′ y ′p1′, el flujo de combustible es ′ff ′ y la presion posterior de la

turbina es ′p7′. Las salidas de dicho bloque son el flujo de ventilacion refrigerante de

la turbina y la temperatura ′m10′ y ′t4′, el flujo de salida de la masa de la turbina y

temperatura ′m8′ y ′t7n′, y velocidad de la turbina ′Wt′. La valvula de paso de la

turbina representa el volumen del conducto y tiene por etiqueta ′c5′; los flujos ′m10′ y

′m8′ son elementos mixtos.

Siendo estas las entrada para el bloque ′c5′, las salidas de la turbina son ′m8′,

′m10′, ′t4′, ′t7n′ descritas anteriormente, junto con la presion del flujo de masa de la

valvula ′m9′. Los bloques sombreados mostrados en la parte superior izquierda inician


el modelo y procesan las diversas opciones tal como se muestra en la figura 1.2, un clic

del raton sobre el bloque apropiado ejecuta el archivo correspondiente en MatlabTM .

El bloque ′Init1′ inicializa los datos del modelo donde estan las variables que

conforman la dinamica de desempeno. Los parametros son la estructuracion de la

secuencia de la demanda de velocidad y el angulo de apertura de la valvula. El bloque

′Plot′ invoca rutinas que trazan la historia de tiempo de muchas de las variables de

los modelos (Silvio Simani y Ron J. Patton 2006).

Figura 1.2: Modelo en Simulink de una turbina, (proporcionado por el Dr. Silvio Simani)

1.3.3 Turbina

En esta seccion se describira cada componente, bloque y etiqueta que forma parte del

modelo que pueden apreciarce en la figura 1.3. Estos son:


1. Conducto de entrada: Las perdidas del conducto de entrada se representan

por el bloque de elemento recesivo ′R1′, y el volumen del conducto por el

elemento acumulador ′c1′.

2. Compresor: Se evalua el flujo de la masa y la validez isentropica como funcion

de presion de la entrada, la presion del paso de la valvula, temperatura de la

entrada y velocidad.

La temperatura de salida es calculada de la proporcion de presion, energıa,

temperatura de entrada. El torque del diagrama del compresor, calcula el flujo

de masa, escapes de temperatura, con media ′cp′ y el diagrama de velocidad.

3. Sistema de combustible: Todos los componentes desde la valvula de escape

del compresor hasta la entrada de la turbina son incluidos en este bloque.

4. Difusor: El bloque ′c2′ representa el volumen de gas del compresor, y promueve

el paso de la valvula del compresor. El bloques ′R2′ representa las perdidas del

difusor.

5. Camara de Combustion: Todas las combustiones se combinan y modelan

como un solo componente. El bloque ′c3′ representa el volumen asociado con los

casos externos de las camaras de la combustion.

El bloque ′R3′ representa las perdidas asociadas con el flujo de aire en las

camaras de la combustion. El bloque ’Combustion’ representa el proceso de la

combustion, para el cual ′t5c′ es la variable de salida de la temperatura, ′t5′ es

la entrada de la temperatura y dt el diferencial de temperatura. La variable dt

es el diferencial de temperatura. La salida de ’Combustion’ es la suma de los

flujos de las masa entrantes, el flujo de la masa de ventilacion es ′m5′ y flujo de


combustible es ′ff ′.

El bloque ′c4′ representa el volumen de gas de la camara de combustion y la ruta

del gas entre la camara de combustion y la entrada de la turbina.

6. Flujos de Aire Refrigerante que Escapan: El flujo de aire se desvıa de la

salida del compresor hacia el estator y el rotor para refrigerar la potencia de la

turbina. Los bloques ’Chk R10’ y ’Chk R6’ representan las trayectorias del flujo

del escape de ventilacion para la refrigeracion del compresor de la turbina.

El flujo a traves de estas trayectorias se considera bloqueado y una de las

funciones sobre la corriente son la temperatura y presion. Diversas trayectorias

en la ventilacion de enfriamiento son suministradas en el motor, estos se han

combinado y se han reducido a dos en el modelo.

La vıa de la trayectoria del bloque ’Chk R6’ se asume antes de la entrada de la

turbina y contribuye a la potencia del eje de la turbina. La vıa de la trayectoria

’Chk 10’ asume la entrada de la trayectoria de gas anterior de la turbina de

poder sin la contribucion del eje de potencia (Silvio Simani y Ron J. Patton 2006).

7. Compresor de la Turbina: El compresor de la turbina es modelado usando

las curvas genericas para calcular la eficacia de la turbina y el flujo de masa.

El rendimiento de la turbina es calculado de una serie de curvas parabolicas

modificadas que relacionan el rendimiento isentropico y la proporcion de

velocidad usando una curva de Hodge que permite calcular el flujo de masa

como una funcion proporcion a la presion (Silvio Simani y Ron J. Patton 2006).

Ambas curvas se encuentran a una escala para ajustar las condiciones del punto

de diseno. Las entradas al bloque son: la temperatura, presion, presion del paso


por la valvula, velocidad del eje. Las salidas del bloque son las siguientes: eje

del torque, potencia en el eje, rendimiento isentropico, temperatura de paso por

la valvula y flujo de masa.

8. Carga mecanica: Este bloque representa la relacion entre el torque de carga

y la velocidad de su eje rotacional. Es el esfuerzo de torsion de la carga y la

diferencia entre el compresor y el torque de la turbina, es la entrada al bloque

′c5′, y la velocidad del eje es la salida del bloque.

Dentro de las perdidas del bloque del torque, la velocidad proporcional asumida

se obtiene de la carga del torque dividida por la carga de inercia integrando

los resultados para proporcionar la velocidad del eje. Tambien es incluido

linealmente a la rampa la velocidad durante la parte inicial de la simulacion al

punto donde los modelos del componente son validos, de una manera similar al

juez de salida en el artefacto real.

Una caracterıstica que tambien se le incluye a la rampa es la velocidad lineal

durante la parte inicial de la simulacion, al punto donde los elementos de los

componentes del modelo esten validados, de manera similar al arranque en el

motor real.


Figura 1.3: Modelo en Simulink del Sub-sistema Turbina HTDU

9. Valvula presurizada: La valvula presurizada es de tipo Mariposa. El

segmento de obstaculo en la valvula estandar se une al paso de la valvula en la

turbina siguiendo la ruta del volumen del conducto, bloque ′c5′.

La valvula la constituye un orificio con dos estados, uno con paso obstruido y

otro de paso abierto. Este orificio se considera como una area eficiente siendo

funcion del angulo de la valvula y estando compensada con la presion.

Los bloques de entrada son: presion del caudal de subida, presion del caudal de

bajada (ambiente), temperatura del caudal de subida y angulo de apertura de

la valvula. La salida del bloque es el flujo de la valvula. El angulo de apertura

de la valvula es parte de una funcion escalon que puede ser programada para

introducir un escalon deseado para el angulo de la valvula durante el tiempo de

simulacion.


1.3.4 Variables del Modelo y Unidades

A continuacion se definen, en las siguientes tablas todas la variables pertenecientes al

modelo (Silvio Simani y Ron J. Patton 2006).

Variables Descripcion

Av Angulo de apertura

ff Flujo de combustible

t1 Temperatura del entorno

p1 Presion del entorno

Tabla 1.1: Definicion de las variables de entrada al modelo.


m1 Flujo masico hacia ′c1’

m3 Flujo masico hacia la camara de combustion

m4 Flujo masico hacia ′c3′

m5 Flujo masico hacia el proceso de combustion



m9 Salida de flujo masico presurizado

p2 Presion de entrada al compresor

p3 Presion de salida del compresor

p4 Presion de entrada ′R3′

p5 Presion de salida hacia la camara de combustion

p7 Presion de retorno a la turbina

Tabla 1.2: Definicion de las variables medidas en el modelo Parte I.



pt Presion de retorno a la turbina

qa Torque de la carga

qc Torque del compresor

qt Torque de la turbina

t3 Temperatura de salida del compresor

t4 Temperatura refrigeracion entrante en ′R3′



t7 Temperatura de entrada a la valvula presurizada

qc Torque del compresor

qt Torque de la turbina

Tabla 1.3: Definicion de las variables medidas en el modelo Parte II.

Esquema de variables segun las tablas 1.1, 1.2 y 1.3 que conforman el modelo

de la turbina de gas combustible.

Figura 1.4: Esquema del Modelo de la Turbina HTDU

en donde:

T1 = Condicion inicial fija de la temperatura del ambiente.

P1 = Condicion inicial fija de la presion del ambiente.


Av(t) = Presion de apertura de la valvula.

ff(t) = El flujo de combustible es la variable manipulada.

Wt(t) = La velocidad de la turbina es la variable a ser controlada.

Yi(t) = Variables medidas las cuales son descritas en las tablas 1.2 y 1.3.

Las variables del modelo en computadora, que estan en unidades SI, son las siguientes:

p - presion [ Nm2 ]

q - torque [N ·m]

m - flujo de masa [kgs]

ff - flujo de combustible [kgs]

t - temperatura [K]

Nt - velocidad de la turbina [rad/s]

av - angulo de apertura [deg]

pct - Potencia de la turbina [W ]

1.3.5 Ecuaciones Matematicas del Modelo

Las ecuaciones que se muestran a continuacion fueron extraıdas del modelo de la turbina

prototipo HTDU plasmadas en MatlabTM −Simulink. El siguiente esquema muestra

una apreciacion general del modelo.

Figura 1.5: Esquema del Modelo en Simulink HTDU


• Turbina HTDU : Este subsistema del modelo contiene la mayor parte de los

modulos que describen el desempeno de la turbina, cada modulo cumple una

labor diferente e importante para el proceso. A continuacion se mostraran las

ecuaciones de los modulos.

• Conducto(entrada, salida): La entrada y salida no tiene ecuaciones en el

modelo, al igual que la camara combustion, pero este modulo ocurre un proceso

en el cual la turbina de gas no sufre transformaciones termodinamicas, solo de

flujo de masa. Para realizar simulaciones dinamicas, su modelo es especialmente

importante porque el conducto es el lugar donde suceden fenomenos inerciales

debidos a la masa del fluido que contiene. Estos fenomenos son descritos en

termino de masa y momento por las ecuaciones (Silvio Simani 2002). En

particular, para la i-seccion del modulo, las siguientes son las ecuaciones que

lo describen:

k =cp

cv(1.1)

δpi

δt=

kRTi

V i(mi −m(i+1)) (1.2)

δmi

δt=

Ai

Li(p(i−1) − pi)−

λikRT(i−1)m2

AiDi(p(i−1) − pi)(1.3)

en donde:

A = area del conducto [m2].

L = longitud de conducto [m].

cp = calor especıfico a presion constante.

cv = calor especıfico a volumen constante.

k= constante de calor especıfico.

D = diametro hidraulico [m].

m = tasa de flujo masico [kg/s].


p = presion [ Nm2 ].

R = constante del gasto masico.

T= temperatura [K].

λ = coeficiente de friccion.

• Compresor: Las siguientes ecuaciones forman parte de esta porcion:

m3 = mp2√t1

(1.4)

t3 = (t2(p3

p2)nisc(α−1

α))(

1

Tcs + 1) (1.5)

qc =m3(t3− t2)

Ncp (1.6)

en donde:

m = flujo de masa [kg/s].

nisc = constante de eficiencia Isentropica del compresor.

cp = calor especıfico para una presion constante [cd].

Tcs= periodo en el compresor [s].

N= velocidad en el compresor [rad/s].

α= promedio gamma del calor especıfico.

Las variables restantes estan definidas en las tablas 1.1, 1.2 y 1.3

• Camara de Combustion: Este subsistema no contiene ecuaciones. En su

mayorıa usa tablas, procesos de interpolacion y extrapolacion de datos. El

Dr. Silvio Simani explica de manera general, lo siguiente en sus notas sobre la

procedencia matematica del proceso.

Combustor (Proceso de combustion interna): Cuando el modelo dinamico

del fluido incluye al combustor, las ecuaciones de balance de masa y momento


son integradas para calcular la presion y la tasa de fluido de masa de gas en la

salida del combustor para condiciones y geometrıas dadas.

La temperatura del gas Ti en la salida del combustor es calculada utilizando la

ecuacion de balance (1.7), bajo la hipotesis de que la combustion y salida de calor

son instantaneos ya que la inercia termica ha sido graduada respecto a la inercia

mecanica (Silvio Simani y Ron J. Patton 2002).

δρ

δt+

δv

δt= 0 (1.7)

ρδρ

δt+ ρv

δv

δt= −(

δp

δt) + F + ρg(

δz

δt) (1.8)

Ti ≈ Ti−1 +(LHV · ηcc)mf

mi · cp(1.9)

en donde:

ρ=coeficiente de la densidad.

v=velocidad [rad/s].

p=presion [ Nm2 ].

F=coeficiente de la fuerza de friccion.

g=gravedad [m/s2].

z=altura de la camara de combustion [m].

mf= flujo masico del combustible [kg/s].

LHV =Valor mınimo de calentamiento [cd]].

ηcc=eficiencia del combustor(camara de combustion).

Sub-Modulo Combustion: Este proceso forma parte del sistema de la

camara de combustion y representa el proceso de combustion de la turbina

cuya resultado es posible debido a la combinacion de combustible mas aire, las

siguientes ecuaciones forman parte del proceso.


mp5 = ff + m5 (1.10)

en donde:

mp5=Salida de flujo de masa [kg/s].

ff=flujo de combustible [kg/s].

m5=flujo de aire [kg/s].

• Turbina: Este modulo presenta las siguientes ecuaciones:

qt =(t6− t7n)cp ·m8

wt(1.11)

pt = (t6− t7n)cp ·m8 (1.12)

tisd = t6(1− (p7

p5)

(α−1)α ) (1.13)

en donde:

tisd = temperatura isentropica [K].


• Modulo Carga: Este modulo presenta las siguiente ecuacion:

qa = (qt− qc)− qloss (1.14)

qloss = tss2TI (1.15)


en donde: qloss=perdida en el torque [N ·m].

tss=tiempo de establecimiento en la turbina [s].

TI=2.018e-5 coeficiente de perdida en el torque [Nm]/[rad/s]2


• Valvula Presurizada: Este modulo presenta las siguientes ecuaciones:

Rg =cp(α− 1)

α(1.16)

mg = (cdv · p7 · Ae√

t7)(

√(α(α− 1)

2Rg)((

p8

p7)2/α − (

p8

p7)(1+α)/α)) (1.17)

en donde:

Ae = Area eficiente [m2].

cdv = 0.97 constante de valvula [cd].

Rg = constante de gas para una tempatura [K].

mg = flujo de masa del gas [kg/s].

Capıtulo 2

Enfoques de Control Aplicados al

Problema

En este capıtulo se dan algunos conceptos de control y se explican las fases para el

diseno de los controladores que se han implantado en este estudio. Se trataran diversos

enfoques tales como: espacio estado, dominio en frecuencia, tecnicas clasicas de control,

sistemas de seguimiento, control optimo cuadratico y anti-rebote.

2.1 Introduccion

La necesidad del hombre de entender y aprender todo lo que sucede en su medio, a

traves de la observacion, lo ha conducido a la necesidad de manipular el mismo, con

la finalidad de resolver los problemas que implica vivir en el.

Estos ambientes estan bajo la influencia de leyes fısicas, quımicas, mecanicas

o de alguna otra ındole, que permiten al hombre concebir mejor y en una forma mas

sistematica su entorno. Todos estos medios podemos llamarlos sistemas y, por tanto,

a traves de las propiedades que lo rigen, extraer modelos que representen de forma

idonea sus comportamientos.

2.1 Introduccion 21

Estos modelos constituyen una representacion abstracta realizada en terminos

de lenguaje y simbologıa. Una vez logrado este objetivo, es posible estudiar un

sistema cualquiera, proporcionando las bases necesarias para tener una idea de como

influenciar (controlar) el comportamiento del sistema real (Hebert Sira-Ramırez y

Francklin Rivas-Echeverrıa 1997).

En ultimo termino, este es el interes practico del modelo en sı. En ocasiones

suele ser complejo, difıcil y engorroso obtener el modelo matematico, por tal razon

se debe utilizar una tecnica llamada identificacion que consiste en la toma de

datos entrada y salida del proceso para luego manipularlos mediante el uso de una

herramienta matematica. De esta manera, se puede obtener un modelo que permitira

conocer e interpretar que sucede con la respuesta, siendo de gran utilidad al momento

de optimizar el desempeno.

Una vez conocidos y estudiados los requerimientos del sistema, se indaga la

tecnica conveniente para hacer control y posteriormente llevar a cabo el diseno. El

modelo de la planta y el controlador son descritos a traves de las ecuaciones que

representan su dinamica. En la actualidad existen diversas tecnicas de control de las

cuales concentraremos nuestra atencion a las siguientes:

1. Enfoque espacio estado: utiliza informacion directa relativa a sus estados;

cuando no esta disponible se debe estimar. La estimacion de semejantes datos o

variables por lo general se denomina observacion.

Una vez conocidos los valores estimados, podemos manipular la respuesta del

sistema a nuestra conveniencia, mediante la compensacion, que consiste en

la modificacion de la dinamica de un sistema, realizada para satisfacer las

especificaciones determinadas. Dicha compensacion puede hacerse bajo diferentes

esquemas, entre los cuales tenemos:

• Construccion del vector de estados: Se considera un sistema lineal o una

aproximacion al mismo para luego obtener la senal de control mediante un

estado instantaneo.

2.1 Introduccion 22

• Control por seguimiento: La respuesta debe seguir una senal de referencia

convirtiendose en un sistema regulador asintoticamente estable tal que e(t)

tienda a cero, dada cualquier condicion inicial e(t=0)(Ogata Katsuiko 1998).

• Control optimo cuadratico: Consiste en sistemas reguladores optimos

basados en ındices de desempeno usando la ecuacion de Riccati.

2. Tecnicas clasicas de control: Son suficientes para resolver problemas de

control en muchas aplicaciones industriales, particularmente cuando la dinamica

del proceso lo permite y los requerimientos de desempeno son modestos. Estas

tecnicas son:

• Familia de controladores Proporcional Derivativo Integral(PID) : Son

controladores cuyo proposito es hacer que el error en estado estacionario,

entre la senal de referencia y la senal de salida, sea cero de manera asintotica

en el tiempo, lo que se logra mediante el uso de la accion integral.

Ademas el controlador tiene la capacidad de anticipar el futuro a traves de

la accion derivativa que tiene un efecto predictivo sobre la salida del proceso.

• Controladores Adelanto-Atraso: Son una buena herramienta para ajustar

las ganancias de un sistema de control a fin de poder cumplir con las

especificaciones dadas. Existen dos maneras de calcular dichas tecnicas:

el diseno de sistemas de control mediante el lugar geometrico de las raıces

y el diseno de sistemas de control mediante la respuesta en frecuencia.

Un compensador en adelanto puede incrementar la estabilidad o velocidad

de respuesta de un sistema; un compensador en atraso puede reducir (pero

no eliminar) el error de estado estacionario. Dependiendo del efecto deseado,

uno o mas compensadores en adelanto y en atraso pueden usarse en varias

combinaciones.

2.2 Definicion de un Sistema de Control 23

3. Control anti-rebote (anti-windup): Este control es usado para evitar efectos

indeseables de saturacion en la actuacion, debido a lımites de operacion en

equipos, producto de la constante integracion del error.

2.2 Definicion de un Sistema de Control

Un sistema de control describe la conducta de un proceso con el fin de regular de

manera optima su funcionamiento.

2.2.1 Conceptos de Control

Unifiquemos terminos a fin de utilizar un lenguaje comun en este aspecto de la

tecnologıa, siempre teniendo en cuenta que nuestro objeto de estudio seran los sistemas

de control lineales.

• Proceso: Se define como una operacion continua, marcado por cambios

graduales que suceden uno al otro en forma fija, con el fin de obtener un resultado,

y en este caso particular se refiere a cualquier operacion sometida a control.

• Planta: Es el elemento fısico que se desea controlar. La planta puede ser un

motor, un circuito electrico, una vehıculo, un tanque, sistema de nivel de lıquidos

entre otros.

• Sistema: Consiste en un conjunto de elementos que actuan coordinadamente

para lograr un objetivo determinado.

• Sistema dinamico: Es aquel en el cual los efectos actuales (salidas) son el

resultado de causas actuales y previas (entradas).

• Variable controlada: Es la cantidad o condicion que se controla.

• Variable manipulada: Es la cantidad o condicion que el controlador modifica

para afectar el valor de la variable controlada.

2.3 Fases para el Diseno de Sistemas de Control 24

• Variable medida: Es la cantidad o condicion que se mide.

• Perturbacion: Es una senal que tiende a afectar adversamente el valor de la

salida de un sistema. Si la misma se genera dentro del sistema se la denomina

interna, de lo contrario se la denomina externa.

2.2.2 Clasificacion de sistema de control

Los sistemas de control se pueden describir de la siguiente manera:

• Sistemas de control lineales - no lineales: la mayorıa de los sistemas fısicos no

son lineales. Sin embargo, algunas veces pueden linealizarse dentro de un rango

para representar un comportamiento aproximado.

• Sistemas de control invariables en el tiempo - variables en el tiempo: los primeros

no son muy comunes en procesos reales debido a que presentan parametros que

no son funcion del tiempo, por el contrario los otros son muy comunes y poseen

parametros que varıan con el tiempo.

• Sistemas de control de tiempo continuo - tiempo discreto: en los primeros todas

las variables son funcion de un tiempo continuo ′t′, en cambio en los otros existen

una o mas variables que son conocidas solo en ciertos instantes de tiempo.

2.3 Fases para el Diseno de Sistemas de Control

En lo que sigue se explicara como se encaran los problemas reales planteando una serie

de pasos a partir de la concepcion del proceso hasta el diseno de los controladores.

2.3.1 Modelado fısico

El objetivo principal de la especificacion de un sistema y el modelado fısico consiste en

proporcionar una descripcion del sistema en estudio lo mas precisa posible, aunque lo

bastante sencilla como para permitir el analisis y diseno subsecuentes. Este modelo se


reduce a un conjunto de ecuaciones que representan su dinamica, tomando en cuenta

que no es unico.

La dinamica de muchos sistemas, tales como mecanicos, electricos, termicos,

economicos, biologicos, etc., se describen en terminos de ecuaciones diferenciales.

Dichas ecuaciones se obtienen a partir de leyes fısicas que gobiernan un sistema

determinado, como las leyes de Newton para sistemas mecanicos, leyes de Kirchoff

para sistemas electricos.

Debemos tener en cuenta que esta etapa es la parte mas importante de todo el

analisis. Para el modelado consideraremos los siguientes topicos:

• Construccion del modelo: Se debe descomponer. Una vez logrado, se

estudia la conducta de cada uno de los componentes y sus interacciones pueden

aproximarse utilizando teorıas conocidas.

• Identificacion de los componentes: Se debe trazar el diagrama de los

componentes que indique las entradas y salidas, ası como las interconexiones de

los componentes internos y los que estan al limite del sistema, se debe identificar

todas las variables, parametros necesarios, ası como sus sentidos y convencion de

signos.

• Relaciones entrada y salida: Describir el comportamiento de estos

componentes.

• Solucion del modelo: La solucion del modelo es fundamental para estudiar

su dinamica. Se requiere de ecuaciones diferenciales correctamente calculadas

que permitan determinar el comportamiento ante varias entradas y con distintas

configuraciones o valores de los parametros.

2.3.2 Identificacion

Habitualmente en la vida real no es sencillo obtener un modelo; esta labor puede

tornarse tediosa y hasta extremadamente complicada, por esta razon debemos utilizar


este metodo alterno para obtener un modelo.

Implementando la tecnica de mınimos cuadrados para generar estos modelos,

se logra una aproximacion al proceso real con alta precision. Esta tecnica consiste en

determinar la recta que mejor se aproxime cuando el error involucrado es la suma de

los cuadrados de las diferencias entre los valores de y en la recta de aproximacion y los

valores de y dados.

El problema general de ajustar la recta de mınimos cuadrados a un conjunto

de datos minimizando la ecuacion (2.1) con respecto a los parametros a y b se muestra

a continuacion (Burden Fairs 1996).

Ecuacion general:

m∑i=1

[yi − (axi + b)]2 (2.1)

Para que se presente un mınimo es necesario:

δ

δa

m∑i=1

[yi − (axi + b)]2 = 0 (2.2)

δ

δb

m∑i=1

[yi − (axi + b)]2 = 0 (2.3)

Por lo tanto (2.2) y (2.3) quedan:

2m∑

i=1

(yi − axi − b)(−xi) = 0 (2.4)

2m∑

i=1

(yi − axi − b)(−1) = 0 (2.5)

Las ecuaciones se simplifican a las ecuaciones normales:

am∑

i=1

x2i + b

m∑i=1

xi =m∑

i=1

xiyi (2.6)

a

m∑i=1

xi + b ·m =m∑

i=1

yi (2.7)


En la actualidad existe una variedad de herramientas cientıficas y tecnologicas

que contribuyen a facilitar este trabajo.

Un programa utilizado por los ingenieros para realizar este tipo de tareas

tiene por nombre MatlabTM (The MathWorks 2007), herramienta poderosa para las

matematicas que puede interpretar datos para luego generar un modelo matematico

del proceso.

2.3.3 Diseno de controladores

En el diseno de control automatico, se debe seguir una serie de tecnicas y algoritmos

matematicos cuya meta en estudio suele medirse por un conjunto de especificaciones de

rendimiento, el proceso de diseno es interactivo. Todos los comportamientos pertinentes

del sistema controlado, determinados durante la solucion del modelo, normalmente son

usados en el proceso de diseno.

2.3.4 Sintonizacion y Validacion del diseno

Existen diversos metodos para la sintonıa de controladores, algunos analıticos y otros

de tanteo. La eleccion correcta de los parametros permite el desempeno optimo y

deseado del sistema, para luego validarlo y asegurar una respuesta satisfactoria ante

diversas entradas o perturbaciones.

2.3.5 Implantacion

Consiste en la puesta en funcionamiento del sistema en un lugar donde desempenara su

labor. Para ello deben afrontarse tareas tan variadas como la formacion de los usuarios

o la carga inicial de los datos con los que trabajara.

2.4 Criterios para el Analisis de Desempeno de un Sistema 28

2.3.6 Pruebas y Depuracion

Son necesarias para comprobar la validez del modelo en un entorno no lineal, ofreciendo

las condiciones para estimar la robustez del sistema propuesto. Estas pruebas deben

practicarse en un lapso de tiempo considerable, para observar el desempeno del mismo,

y ası realizar la depuracion o afinacion de las ganancias.

2.4 Criterios para el Analisis de Desempeno de un

Sistema

El analisis del desempeno permite interpretar que sucede con los resultados obtenidos,

a partir de la respuesta de un sistema, de allı la premisa de optimizar el diseno.

La respuesta de un sistema de control en el tiempo se divide en dos etapas que se

explicaran a continuacion.

2.4.1 Respuesta transitoria

La respuesta temporal en lazo cerrado, normalmente es implementada ante una

entrada del tipo escalon, puesto que contiene un espectro lo suficientemente amplio de

frecuencias. Por lo general los requerimientos de la respuesta transitoria se describen

en terminos de polos y ceros, que se denominan polos dominantes del sistema que se

encuentran mas cercanos al origen en el plano complejo.

Las principales caracterısticas que modifican los polos dominantes son tiempo

de crecimiento (td), sobredisparo maximo (Mp), tiempo de establecimiento (ts),

tiempo pico (tp) y tiempo de retardo (tr) para sistemas de segundo orden, ya que las

especificaciones transitorias se expresan en funcion de estas caracterısticas.

2.5 Enfoque de Espacio Estado 29

La respuesta transitoria de un sistema de control practico exhibe con frecuencia

oscilaciones amortiguadas antes de alcanzar el estado estable (Ogata Katsuiko 1998).

Figura 2.1: Curva de respuesta de sistema de orden 2, ante una entrada escalon unitario

2.4.2 Respuesta Estacionaria

Por respuesta en estado estable, se refere a la manera en la cual se comporta la

salida del sistema conforme el tiempo tiende a infinito. Esta respuesta tiene como

caracterıstica principal seguir la referencia hasta un punto de equilibrio (set point) o

valor deseado de referencia y estabilizarse.

Entre mas precision exista entre la referencia y la salida, el error en estado

estable sera mas pequeno.

2.5 Enfoque de Espacio Estado

Este metodo se considera particularmente atractivo porque, proporciona la oportunidad

suficiente de aplicar la intuicion fısica sobre las causas y efectos en los componentes


de sistemas, siendo transparente para el tipo de sistema que se modela: mecanico,

quımico, electrico, etc. Este enfoque lo conforma:

• Variables de estado: Son el conjunto linealmente independiente de variables

que se utilizan para especificar el estado de un sistema.

El numero de variables de estado determina el orden del sistema, como se observa

en la siguiente ecuacion:

c1x1 + c2x2 + ... + cnxn (2.8)

• Ecuacion de estado: Representacion matematica que satisface particularmente

el requerimiento del valor inicial mediante una ecuacion diferencial ordinaria de

primer orden (para sistemas de tiempo continuo) o una ecuacion en diferencia

(para sistemas en tiempo discreto). Este conjunto de ecuaciones diferenciales

ordinarias de primer orden son representadas en (2.9),(2.10) y (2.11) (Ogata

Katsuiko 1998).

x1 = f1(x1, x2, ..., xn, u1, u2...ur, t) (2.9)

x2 = f1(x1, x2, ..., xn, u1, u2...ur, t) (2.10)

xn = f1(x1, x2, ..., xn, u1, u2...ur, t) (2.11)

• Espacio estado: Se describe como el espacio matematico (de dimension n) cuyas

coordenadas son variables de estado. Por tanto, en cualquier instante, el estado

del sistema esta representado por un punto en el espacio estado.

Esta descripcion explica el uso del termino enfoque de espacio estado, para

describir un sistema de estado como se observa en las ecuaciones (2.12) y (2.13).

x = Ax + Bu (2.12)

y = Cx (2.13)

en donde:

A=Matriz de estados.


B=Matriz de entrada.

C=Matriz de salida.

x=Variable de estado.

u=Senal de control.

y=Senal de salida.

Caracterısticas del enfoque espacio estado

Para la controlabilidad de los sistemas lineales, conocer los conceptos de

controlabilidad y observabilidad juega un papel importante debido a que

gobiernan la existencia de una solucion para un problema de control.

Controlabilidad: Se dice que el proceso es controlable si cada variable de

estado del proceso se puede controlar para llegar a un cierto objetivo en un

tiempo finito, a traves de algun control no restringido (Benjamin C. Kuo 1996).

Observabilidad: Un sistema es completamente observable si cada variable de

estado del sistema afecta alguna de las salidas (Benjamin C. Kuo 1996).

2.5.1 Construccion del Vector de Estados

El problema de reconstruccion del vector de estado, o, equivalentemente, de

construccion (diseno) de un observador para un sistema lineal fue resuelto a principios

de los anos 60 por un profesor americano de nombre David Luenberger (Salvador

Saucedo 2006).

En su honor, el observador determinıstico de estado recibe en la literatura el

nombre de ”Observador de Luenberger”. Para este diseno consideraremos un sistema

lineal en modelo de espacio estado y seleccionaremos la senal de control como:

u = −Kx (2.14)


en donde:

K=Ganancia del vector de estados.

Esto significa que la senal de control se determina mediante un estado

transitorio. Tal esquema se denomina realimentacion del estado. La matriz K

se denomina matriz de ganancias de la realimentacion del estado.

2.5.2 Observador de Estados

Un observador de estado estima las variables con base en las mediciones de salida y

de control. Las variables de estado que no estan disponibles son las que se necesitan

medir mediante este proceso llamado observacion.

Existen metodos para estimar las variables de estado que no se miden sin un

proceso de diferenciacion. Cuando se necesita captar todas las variables del sistema,

sin importar si algunas estan disponibles para una medicion directa, se denomina

observador de orden completo. Se debe tener en cuenta que el sistema en estudio

requiere ser controlable y observable.

Se dice que un sistema es controlable en el tiempo to si se puede llevar de

cualquier estado inicial x(to) a cualquier otro estado, mediante un vector de control

sin restricciones, en un intervalo de tiempo finito (Salvador Saucedo 2006). Se dice

que un sistema es observable en el tiempo to si, con el sistema en el estado x(to),

es posible determinar este estado a partir de la observacion de la salida durante un

intervalo de tiempo finito.

En muchos casos practicos el vector de estado observado x se usa en la

realimentacion del estado para generar el vector de control deseado.


Algoritmo de diseno

Suponga que se aproxima el estado mediante x. Considerando el sistema

definido mediante (2.12) y (2.13), entonces se tiene:

˙x = Ax + Bu + Ke(y − Cx) (2.15)

en donde:

x=Representa el vector de estado aproximado

Ke=Ganancia del observador de estado para manipular el error y asegurar una

dinamica asintoticamente estable.

Sustituyendo(2.13) en (2.15) y luego restando con (2.12) obtenemos:

(x− ˙x) = Ax− Ax + Bu−Bu + Ke(Cx− Cx) (2.16)

(x− ˙x) = (A−KeC)(x− x) (2.17)

el error de observacion es:

Ae = |A−KeC| (2.18)

L(s) = |sI − Ae| (2.19)

en donde:

Ae = Matriz que representa la dinamica del error.

I = Matriz identidad.

L(s)= Matriz de ganancias del observador.

Para calcular los valores Ke del observador, se debe igualar la matriz L(s) a un

polinomio deseado Pd(s).

• Polinomio con raıces reales negativas.

Pd(s) =n∏

i=1

(s + ai) (2.20)


• Polinomio con raıces complejas conjugadas.

Pd(s) =n∏

i=1

(s + ai ± bij) (2.21)

en donde:

a = Polo real en el plano complejo.

b = Polo imaginario en el plano complejo.

Finalmente, sustituyendo (2.18) en (2.19) e igualando al polinomio deseado se tiene.

Pd(s) = |sI − A + KeC| (2.22)

De tal manera, se pueden calcular los valores de la matriz Ke del observador. La

siguiente figura (3.2) ilustra el diagrama de bloques de la realimentacion de un sistema.

Figura 2.2: Diagrama de bloques del esquema de control realimentado con un observador


2.5.3 Diseno de Sistemas por Seguimiento

En el control proporcional de una planta, cuya funcion de transferencia no posee un

integrador, hay un error en estado estable, o desplazamiento, en la respuesta para

una entrada escalon. Tal error se elimina si se incluye la accion de control integral

en el controlador. En el control integral de una planta, la senal de salida a partir del

controlador, es en todo momento el area bajo la curva del error hasta tal momento.

La senal de control u(t) tiene un valor diferente de cero cuando la senal de

error e(t) es cero, como se aprecia en la figura 2.3(a). Esto es imposible en el caso del

controlador proporcional, dado que una senal de control diferente de cero requiere de

una senal de error diferente de cero. La figura 2.3(b) muestra la curva e(t) contra t y

la curva u(t) correspondiente contra t cuando el controlador es de tipo proporcional.

Figura 2.3: (a) Graficas de las curvas e(t) y u(t) que muestran una senal de control diferente de cero

cuando la senal de error es cero(control integral); (b) graficas de las curvas e(t) y u(t) que muestran una

senal de control de cero cuando la senal de error es cero (control proporcional).

Si la planta no tiene integrador (planta de tipo 0), el principio basico del diseno

de un sistema de seguimiento de tipo 1 es insertar un integrador en la trayectoria

directa entre el comparador de error y la planta, como se observa en la figura 2.4. Las

ecuaciones son descritas a continuacion (Ogata Katsuiko 1998).

u = −Kx + KIζ (2.23)

ζ = r − y = r − Cx (2.24)


en donde:

ζ=Salida del integrador.

r=Senal de referencia.

KI=Ganancia de integracion.

Algoritmo de diseno

1. Para evitar la posibilidad de que el integrador insertado se cancele por un cero

de la planta en el origen la funcion de transferencia de la planta, Gp(s) no debe

tener ceros en el origen.

2. Suponga que la entrada de referencia (funcion escalon) se aplica en t = 0. En

este caso, para t > 0, la dinamica del sistema se describe mediante una ecuacion

que es una combinacion de las ecuaciones (2.12),(2.13) y (2.23):

x(t) = Ax(t) + Bu(t)

ζ(t) = r(t)− Cζ(t) (2.25)

3. Se disenara un sistema asintoticamente estable, tal que x(t = ∞), x(ζ = ∞)

y x(u = ∞) tiendan a valores constantes, respectivamente. Ası, en un estado

estable, ζ(t) = 0 y obtenemos, y(∞) = r.

4. Considerando que r(t) es una entrada escalon, tenemos que r(∞) = r(t) = r

constante para todo t > 0. Por ello debemos restar las ecuacion (3.25) cuando

t = 0 y t = ∞. La idea basica es disenar un sistema regulador estable de (n

+ 1 )-esimo orden que lleve a cero el nuevo vector de error e(t), dada cualquier

condicion inicial e(t=0), de tal manera que:

˙xe(t) = Axe(t) + Bue(t) (2.26)

˙ζe(t) = r(t)− Cζe(t) (2.27)

ue(t) = −Kxe(t) + Kiζe(t) (2.28)


Luego se define la matriz dinamica del error y los siguientes vectores:

Ae =

[A 0

−C 0

](2.29)

e(t) =

xe(t)

ζe(t)

(2.30)

B =

B

0

(2.31)

K = [K −KI ]. (2.32)

Agrupando terminos y sustituyendo, resulta:

e(t) = (Ae− BK)e(t) (2.33)

que es, finalmente, la ecuacion de estado del error.

5. Para calcular los valores de K, se aplica la siguiente ecuacion:

Pd(s) = |sI − Ae + BK|. (2.34)

En el diseno actual, es necesario considerar varias matrices K diferentes

(que correspondan a varios conjuntos distintos de valores caracterısticos deseados) y

realizar simulaciones en computadora para encontrar aquella que resulten con el mejor

desempeno general del sistema. Como ocurre normalmente, no todas las variables de

estado se pueden medir en forma directa. En ese caso, necesitamos usar un observador

de estado. La figura 3.4 muestra un diagrama de bloques de un sistema de seguimiento

de tipo 1 con un observador de estado.


Figura 2.4: Sistema de seguimiento de tipo 1 con un observador de estado

2.5.4 Diseno de Control Optimo Cuadratico

El diseno de los sistemas de control optimo y los sistemas reguladores optimos basados

en ındices de desempeno cuadraticos se reduce a la determinacion de los elementos

de la matriz K. Una ventaja de usar el esquema de control optimo cuadratico es que

el sistema disenado sera estable, excepto en donde el sistema no es controlable. Al

disenar sistemas de control con base en la minimizacion de los ındices de desempeno

cuadratico necesitamos resolver las ecuaciones de Riccati, a fin de minimizar el ındice

de desempeno:

J =

∫ ∞

0

(xT Qx + uT Ru)dt (2.35)

en donde:

Q = una matriz simetrica real semidefinida positiva.

R = una matriz hermitiana o simetrica real definida positiva.

Las matrices Q y R determinan la importancia relativa del error y del gasto de este

sistema.

Ahora se mostrara que una funcion de Liapunov se usa efectivamente en la

solucion de este problema. Supongamos que:

xT Qx = − d

dt(xT Px) (2.36)


en donde P es una matriz hermitiana o simetrica real definida positiva.

En este caso, obtenemos:

xT Qx = −xT Px− xT Px = −xT AT Px− xT PAx = −xT (AT P + PA)x (2.37)

Mediante el segundo metodo de Liapunov, sabemos que, para una Q determinada,

existe P, si A es estable, tal que:

AT P + PA = −Q. (2.38)

Siguiendo el analisis obtenido al resolver el problema de optimizacion de parametros,

se establece:

x = Ax + Bu = Ax−BKx = (A−BK)x (2.39)

Suponiendo que la matriz en (2.39) es estable, o que los valores caracterısticos de esta

matriz tienen partes reales negativas:

xT (Q + KT RK)x = −xT Px− xT Px = −x[(A−BK)T P + P (A−BK)]. (2.40)

Comparando ambos miembros de esta ultima ecuacion y considerando que la misma

debe ser valida para cualquier x, requerimos que:

(A−BK)T P + P (A−BK) = −(Q + KT RK) (2.41)

en donde P es una matriz hermitiana o simetrica real definida positiva.

Mediante el segundo metodo de Liapunov, si en la ecuacion (2.39) la matriz es

estable, existe una matriz P definida positiva que satisface la ecuacion (2.41).

Dado que se ha supuesto que R es una matriz hermitiana o simetrica real

definida positiva, escribimos.

R = T T T (2.42)


en donde T es una matriz no singular.

Ası, la ecuacion (2.41) se puede escribe como:

(A−KB)T P + P (A−BK) + Q + KT RK = 0. (2.43)

Rescribiendo esta ultima, se obtiene:

AT P +PA+[TK− (T T )−1BT P ]T [TK− (T T )−1BT P ]−PBR−1BT P +Q = 0. (2.44)

La minimizacion de J con respecto a K requiere de la minimizacion, como:

xT [TK − (T T )−1BT P ]T [TK − (T T )−1BT P ]x. (2.45)

Dado que esta ultima expresion es no negativa, el mınimo ocurre cuando es cero, o

cuando:

TK = (T T )−1BT P. (2.46)

Por tanto:

K = T−1(T T )−1BT P = R−1BT P (2.47)

De la ecuacion (2.47) resulta la matriz de ganancia K. Ası, el problema de control

optimo cuadratico es lineal cuando el ındice de desempeno se obtiene mediante la

siguiente ecuacion:

u(t) = −Kx(t) = −R−1BT Px(t) (2.48)

La matriz P de la ecuacion anterior debe satisfacer la ecuacion (2.35) o la ecuacion

reducida siguiente:

AT P + PA− P (BR−1B)T P + Q = 0, (2.49)

la cual se denomina ecuacion matricial reducida de Riccati (Ogata Katsuiko 1998).

2.6 Enfoque en el Dominio de Frecuencia 41

2.6 Enfoque en el Dominio de Frecuencia

En el enfoque de la respuesta en frecuencia, se analiza el estado estacionario del

sistema ante perturbaciones (senoidales) en distintas frecuencias, con la expectativa de

deducir la dinamica del sistema. El sistema que se somete a una prueba de respuesta

en frecuencia tiene que ser asintoticamente estable, es decir, que todos los valores

propios se hallan en el plano complejo de la mitad izquierda, excluido el eje imaginario.

Para un sistema marginalmente estable (sistema estable con valores especıficos

en el eje imaginario), cualquier perturbacion persistira a menos que se imponga una

condicion inicial precisa que pueda cancelar con exactitud la perturbacion de inicio.

2.6.1 Estabilidad Relativa

El proposito del diseno mediante este enfoque es obtener una determinada velocidad

de respuesta del sistema en lazo cerrado. Los criterios mas usados para medir la

estabilidad relativa son el margen de fase y el margen de ganancia. Los valores para

un buen desempeno son un margen de ganancia mayor a 6 db (decibeles) y un margen

de fase entre 30◦y 60◦, por lo general, para sistemas de segundo orden.

2.6.2 Sensibilidad de los Parametros

La finalidad del diseno es lograr que el sistema en lazo cerrado se comporte de forma

deseada a pesar de errores (incertidumbres) o cambios en el modelo del proceso. La

sensibilidad es el grado con el cual los cambios en los parametros afectan el desempeno

global sistema. Los valores de ganancias de lazo altas contribuyen a la disminucion de

la sensibilidad, en lazo cerrado, debido a cambios en los parametros. Se debe buscar

un compromiso entre estabilidad y sensibilidad.

2.7 Tecnicas Clasicas de Control 42

2.6.3 Rechazo de Perturbaciones

La idea del diseno, en este caso, es lograr que la senal controlada sea igual a la senal de

referencia a pesar de la influencia permanente de las perturbaciones. La realimentacion

reduce el efecto de las perturbaciones y el ruido en el comportamiento del sistema.

2.7 Tecnicas Clasicas de Control

El objetivo principal es presentar los procedimientos clasicos de control que son de

gran uso en el campo de control de procesos. Actualmente son ampliamente utilizados

en la industria moderna, controlando un elevado porcentaje de procesos industriales

en lazo cerrado.

2.7.1 Controlador por Adelanto de Fase

En la presente seccion se explicara el procedimiento matematico para el diseno de

controladores adelanto-atraso . Consideremos la siguiente funcion de transferencia:

G(s) =s + zc

s + pc

(2.50)

en donde:

zc=Cero del controlador.

pc=Polo del controlador.

Gc(s)=Funcion de transferencia del controlador.


Figura 2.5: (a) Aporte de fase positiva: Adelanto de fase. (b) Aporte de fase negativa: Atraso de fase

la figura (2.5) nos muestra los aportes de fase del controlador respecto de un

punto en el plano complejo.

en donde:

θc=Fase aportada por el controlador.

θzc=Fase de zero del controlador.

θpc=Fase del polo del controlador.

Diseno del Controlador por Adelanto de Fase

Para una funcion de transferencia del controlador descrita por:

Gc(s) = Kc(s + 1

T)

(s + 1αT

)(2.51)

o tambien por :

Gc(s) = Kcα(Ts + 1)

(αTs + 1), (2.52)

Condicion:

0 < α < 1, (2.53)

en donde:

Kc(s) = ganancia del controlador, α = factor de atenuacion, T = periodo, la fase


aportada por el compensador es:

θc = θzc − θpc = tan−1(ωT )− tan−1(ωαT ). (2.54)

Diferenciando con respecto de ω, para encontrar la frecuencia donde ocurre el maximo

avance de fase, se tiene que:

δθ

δω=

T

1 + (ωT )2=

αT

1 + (αωT )2. (2.55)

Igualando a cero y resolviendo para ω, tenemos:

ωmax =1

T√

α. (2.56)

Sabiendo que:

θc = θzc − θpc ⇒ tan(θzc − θpc) =tan(θzc)− tan(θpc)

1 + tan(θzc) tan(θpc)(2.57)

Y tambien que:

θzc = arctan(ωT ) (2.58)

θpc = arctan(αωT ) (2.59)

se obtiene que el avance de fase maximo puede calcularse a traves de la siguiente

expresion:

tan(θm) =

1T√

αT − 1

T√

ααT

1 + 1T√

αT − 1

T√

ααT

=1− α

2√

α(2.60)

o, de manera directa, como:

θm = tan−1(1− α

2√

α) = sin−1(

1− α

1 + α) (2.61)


Algoritmo de Diseno

1. Determinar la ganancia requerida para eliminar el error en estado estable. Es

decir: es:

Ke = lim(s→0)

Gp(s) = lim(s→0)

K

∏mi=1(s + zi)∏nj=1(s + pj)

(2.62)

Ked = Kc1Ke (2.63)

en donde:

Ked=Constante deseada de error.

Kc1=Constante del controlador.

Ke =Constante de error de la planta.

K=Ganancia de la planta.

2. Dibujar el diagrama de Bode de la funcion de transferencia del sistema

compensado con la ganancia Kc1, esto es:

G(s) = Kc1K

∏ni=1(s + zi)∏mj=1(s + pj)

(2.64)

3. Calcular el valor del margen de fase actual y determinar el adelanto de fase que

sera aportado por el controlador. Esto es:

θc = θd − θp (2.65)

en donde:

θd=Margen de fase deseado. θp=Margen de fase de la planta.

4. Determinar el factor de atenuacion usando la ecuacion (2.56).

5. Determinar la frecuencia en la cual el sistema no compensado tiene un valor de

−20 lg 1√α. Seleccionar esta frecuencia como ωmax.

6. Encontrar el valor del cero del compensador, usando la ecuacion (2.61).

7. Determinar el valor del polo del compensador conocido el valor de 1/T y α.


8. Determinar el valor de Kc, sabiendo que Kc = Kc1/α.

Comentarios: Esta propuesta de diseno puede no alcanzar condiciones deseadas

sobre los tiempos de respuesta (ancho de banda, en el dominio frecuencial). En este

caso, podrıa ser necesario un rediseno del compensador.

Tomando en cuenta que el compensador ”empujara” la frecuencia de cruce de

ganancia actual hacia las altas frecuencias, sera, es necesario incorporar un factor

de correccion de error de (10◦ aprox.) sobre la fase que aportara el compensador.

Limitaciones del diseno: Si el sistema original es inestable o con un margen de

estabilidad baja (margen de fase pequeno), el adelanto de fase adicional que debe ser

aportado por el controlador puede ser excesivo, lo cual exige un aporte de fase mas

grande. Esto puede afectar la magnitud de la respuesta frecuencial alrededor de la

frecuencia que determina el ancho de banda (aumento de magnitud), incrementando

el ancho de banda y, en consecuencia , aumentando la sensibilidad a los ruidos que se

introducen por la entrada del sistema. En este sentido, si se requiere de un adelanto

de fase de mas de 90◦, se debe emplear un controlador de etapas multiples (Benjamin

C. Kuo 1996).

2.7.2 Controlador por Atraso de Fase

El principio diseno de este tipo de compensador requiere de un valor de ganancia que

satisfaga condiciones de error en estado estable, ubicando el cero de controlador para

lograr la atenuacion necesaria en la curva de magnitud en la frecuecia donde ocurre el

margen de fase deseado.

Diseno del Controlador por Atraso de Fase

la funcion de transferencia del controlador se describe como:

Gc(s) = Kc(s + 1

T)

(s + 1βT

)(2.66)


o, tambien como:

Gc(s) = Kcβ(Ts + 1)

(βTs + 1)(2.67)

Condicion:

β > 1 (2.68)

siendo beta el factor de atenuacion.

Algoritmo de Diseno

1. Determinar la ganancia de lazo que el requerimiento de error en estado estable.

Esto es:

Ked = Kc1Ke (2.69)

2. Dibujar el diagrama de Bode de la funcion de transferencia del sistema

compensado con la ganancia Kc1, es decir.

G(s) = Kc1K

∏ni=1(s + zi)∏mj=1(s + pj)

(2.70)

y ubicar la frecuencia deseada (ωgd) donde ocurre el margen de fase deseado (θd +

eφ). El margen de fase deseado esta asociado a una especificacion de respuesta

temporal. El termino (eφ) indica una correccion de error de fase entre 5◦ y 10◦

introducido por el compensador alrededor de ωgd.

3. Elegir un cero del compensador (frecuencia de cruce en 1/T ) una decada por

debajo de la frecuencia de cruce de ganancia deseada.

1

T=

ωgd

10(2.71)

4. Para lograr que (ωgd) sea la frecuencia de cruce de ganancia del sistema

compensado, se aprovecha la atenuacion introducida por la curva de magnitud

del compensador, la cual es igual a −20 lg 1√α. Ası debe lograrse que:

|G(jωgd)|dB = −20 lg(β). (2.72)

Usando la ecuacion anterior, se podra calcular el valor de β.


5. Obtener la ubicacion del polo del compensador (frecuencia de cruce en 1/βT ).

6. Determinar el valor de Kc, sabiendo que Kc = Kc1/β.

Diseno de controladores por etapas multiples: para este metodo, el uso del lugar

de las raıces puede resultar algo engorroso. Ası que es conveniente el diseno basado en

la respuesta frecuencial. En este caso, el diseno se traduce en compensar una parte de

la fase requerida con un compensador y a partir de este nuevo sistema compensado, se

terminan de satisfacer los requerimientos con el segundo compensador (Benjamin C.

Kuo 1996).

2.7.3 Controlador PD (Proporcional Derivativo)

Este tipo de controlador pertenece a la familia de controladores PID (Proporcional

Derivativo Integral), que aplican en forma casi general a la mayorıa de los sistemas

de control. Es un hecho bien conocido que los esquemas de control PID basicos

y modificados han demostrado su utilidad para aportar un control satisfactorio.

Interpretacion del controlador PD en el dominio frecuencial

El control PD se comporta como un filtro paso alto, agregando fase positiva (adelanto

de fase) al sistema no compensado. Esta caracterıstica permite:

• Mejorar al margen de fase de un sistema de control y, en consecuencia, mejorar

la estabilidad relativa del sistema.

• La curva de magnitud del controlador ”empuja” la frecuencia de cruce de ganancia

hacia el valor mas alto.

• Incrementar el ancho de banda lo que reduce los tiempos de respuesta pero

acentua los ruidos a alta frecuencia.


Diseno del Controlador PD

Consiste basicamente en localizar la frecuencia de corte del controlador ωgd = KpKd

para lograr una mejora en el margen de fase deseado en la nueva frecuencia de cruce del

sistema compensado. Para propositos de diseno en el dominio frecuencial, se representa

la funcion de transferencia del controlador de la siguiente manera:

Gc(s) = Kp + sKd = Kp(1 +Kd

Kp

) = Kp(1 +1

ωgd

) (2.73)

en donde Kp y Kd son las constantes del controlador.

u(t) = Kpe(t) + Kdδe(t)

δt(2.74)

U(s) = KpE(s) + sKdE(s) (2.75)

Algoritmo de Diseno

1. Fijar el valor de Kp para los satisfacer requerimientos de error en estado estable.

2. Dibujar el diagrama de bode del sistema compensado con Kp.

3. Encontrar el margen de fase del sistema, θp, y determinar la contribucion angular

del controlador θc, para alcanzar el margen de fase deseado θd. Recordar agregar

un poco mas de fase deseada en virtud del desplazamiento de la frecuencia de

corte del sistema no compensado por efecto de la incorporacion del cero del

controlador.

θc = θd − θp + eφ (2.76)

4. Elegir ωgd a la izquierda de la frecuencia de cruce de ganancia si θc > 45◦ o la

derecha si θc < 45◦.

5. Dibujar el nuevo diagrama de bode compensado y ajustar el controlador si es

necesario.


Con este diseno en el dominio frecuencial, tambien es posible satisfacer otra

especificacion del dominio temporal, tomando en cuenta que existe la siguiente

ecuacion que relaciona el margen de fase y la frecuencia donde este ocurre

(frecuencia de cruce de ganancia), en la respuesta frecuencial en lazo abierto,

con el tiempo de asentamiento del sistema en lazo cerrado:

tan θp =8

ωcgTs

, (2.77)

en donde:

ωcg=Frecuencia de cruce de ganancia.

Ts=Tiempo de estableciemiento.

2.7.4 Controlador PI (Proporcional Integral)

Aun cuando el controlador PI esta concebido para mejorar el error en estado

estacionario, este puede disenarse para lograr algunas mejoras en estado transitorio.

Sin embargo, teniendo en cuenta que el controlador PD puede ser usado para lograr

las mejoras en el regimen transitorio de manera adecuada, el PI sera disenado para

minimizar, en lo posible, el deterioro del desempeno transitorio, bien sea no compensado

o previamente compensado por un control del tipo PD (Benjamin C. Kuo 1996). Dicho

efecto es asegurado localizando el cero del controlador muy cerca del polo en origen,

con valores Ki y Kp relativamente pequenos.

Interpretacion del controlador PI en el dominio frecuencial

El control PI se comporta como un filtro pasa baja, agregando fase negativa (retraso

de fase) al sistema no compensado. Esta caracterıstica permite:

• Mejorar el margen de ganancia, pero puede perjudicar la estabilidad, segun la

posicion de la frecuencia de corte.

• Si Kp < 1, la curva de magnitud del controlador tiene un efecto atenuante,

”empujando” la frecuencia de cruce de ganancia hacia un valor mas bajo,


disminuyendo el ancho de banda. Esto filtra los ruidos a alta frecuencia pero

aumenta los tiempos de respuesta.

Diseno del Controlador PI

Para propositos de diseno en el dominio frecuencial, representamos a la funcion

de transferencia del controlador de la siguiente manera:

Gc(s) = Kp +1

sKi = (

Ki + Kp

s) =

Ki(1 + 1ωgd

)

s(2.78)

u(t) = Kpe(t) + Ki

∫ t

0

dt (2.79)

U(s) = KpE(s) +1

sKiE(s) (2.80)

Teniendo en cuenta que:

ωgd =Ki

Kp

(2.81)

y que Ki y Kp son las contantes del controlador, tal y como se describen en (2.79) y

(2.80).

Algoritmo de Diseno

1. Elegir el valor de Kp para proporcionar un factor de atenuacion adecuado sobre

la curva de magnitud del sistema no compensado con el fin de lograr el margen

de fase deseado.

2. Dibujar el diagrama de bode del sistema no compensado.

3. Encontrar el valor de frecuencia para el cual se logra el margen de fase deseado.

Denotar a esta frecuencia como la frecuencia de cruce deseada.


4. Elegir el valor de Kp para lograr una atenuacion de magnitud, tal que el sistema

compensado logre la frecuencia de cruce deseada. En efecto, el valor de Kp viene

dado por:

Kp = 10−|Gp(ωgd)|

20 (2.82)

Otra forma es:

Kp = 1/|Gp(ωgdj)| (2.83)

ωgd =Ki

Kp

=ωgd

10. (2.84)

Despejando se tiene:

Ki =(ωgdj)

10Kp. (2.85)

Recordemos que alrededor de ωgd, la curva de ganancia es afectada solo por

Kp. Por otro lado, en virtud de la eleccion de la frecuencia de corte del controlador,

en principio, a una decada por debajo de ωgd, generalmente se incorpora un termino

correctivo de error de fase de 5◦. Ası, ωgd se elige para una diferencia de fase de θd +5.

2.7.5 Diseno del controlador PID

El control PID se comporta como un filtro pasa banda, agregando fase negativa (retraso

de fase) en las bajas frecuencias y agregando fase positiva (adelanto de fase) en las altas

frecuencias (Benjamin C. Kuo 1996). Estas caracterısticas permiten:

• Mejorar el margen de ganancia.

• Mejorar el margen de fase.


• Minimizar los efectos no deseados del adelanto de fase (empujar la frecuencia de

cruce de ganancia hacia un valor mas alto), esto significa que pueden atenuarse

con el efecto del atraso de fase (empujar la frecuencia de cruce de ganancia hacia

un valor mas bajo), a traves de una ubicacion adecuada de las frecuencias de

corte del controlador PID.

Para propositos de diseno en el dominio frecuencial, representamos a la funcion

de transferencia del controlador de la siguiente manera:

Gc(s) = Kp +Ki

s+ Kds (2.86)

en donde Kp, Ki y Kd son las constantes del controlador.

u(t) = Kpe(t) + Kdδe(t)

δt+ Ki

∫ t

0

dt (2.87)

U(s) = KpE(s) + sKdE(s) +1

sKiE(s) (2.88)

Algoritmo de Diseno

Este enfoque de diseno es sencillo debido a que requiere de la construccion previa del

PD y el PI. Utilizado las siguientes ecuaciones obtenemos los parametros de diseno:

Kp = Kp2 + Kd1Ki2 (2.89)

Kd = Kd1Kp2 (2.90)

Ki = Ki2 (2.91)

en donde:

Kp2=Ganancia proporcional del controlador PI.

Kd1=Ganancia derivativa del controlador PD.

Ki2=Ganancia integral del controlador PI.

2.8 Esquema Anti-Rebote (Anti-Windup) 54

2.8 Esquema Anti-Rebote (Anti-Windup)

En la mayorıa de los problemas de control existen limitaciones en los actuadores. El

principal efecto indeseable se debe a la presencia de un integrador en el controlador

(suponiendo que haya alguno, como sucede con los controladores de la familia PID)

que continuara integrando aun mientras la entrada se encuentra saturada.

Ası, el estado del integrador en cuestion puede alcanzar valores excesivos, que

deterioraran la respuesta transitoria del sistema, generalmente produciendo grandes

sobrevalores (Virginia Mazzone 2002).

Si estas limitaciones se ignoran en la etapa de diseno, el desempeno real del

sistema de control puede sufrir una severa degradacion respecto al esperado si la senal

de control alcanza sus lımites.

2.8.1 Metodologıa General para el Diseno del esquema

Antirebote (Anti-windup)

1. Disenar un controlador PID o PI ideal.

2. Definir valores de saturacion (Maximo y Mınimo) previos a un estudio sobre los

lımites de desempeno del actuador.

3. Elegir un valor de ganancia Kcr, necesaria para compensar el error producto de

la accion integral, tal como lo muestra la figura 2.6. Este valor se calcula de la

siguiente manera:

• Kcr = 1Ki

, en el caso de un PI.

• Kcr =√

Ki ·Kd, en el caso de un PID.

4. Ajustar los valores de no lograr las especificaciones deseadas.

2.8 Esquema Anti-Rebote (Anti-Windup) 55

Cuando existe saturacion, se modifica el valor del efecto integral para que sea

igual a la salida lımite del actuador. Si el lazo no presenta saturacion el lazo no actua,

siendo Kcr la ganancia que indica la rapidez con la que se hace la correccion.

En caso que la salida no se pueda medir se toma un modelo matematico del

actuador (Jorge Vento 2005).

Figura 2.6: Esquema de control Anti-windup

Capıtulo 3

Identificacion de un Modelo de

Turbina a Gas Combustible

El siguiente capıtulo contiene el procedimiento para identificar el modelo de la turbina

facilitado por el Dr. Silvio Simani.

3.1 Introduccion

Los sistemas convencionales, basados en compresores, procesos de combustion,

valvulas y turbinas son procesos industriales que resultan tediosos y engorrosos al

momento de ser modelados fısicamente. Por tal razon, es usado el metodo de mınimos

cuadrados basado en datos entrada-salida de un proceso industrial.

Bajo la tecnica previamente mencionada la companıa Alstrom extrajo un

modelo de la turbina prototipo HDTU en Simulink−MatlabTM el cual fue facilitado

por el Dr. Silvio Simani, aportando informacion importante sobre el funcionamiento

del sistema. Para luego ser manipulado debido a su complejidad, desconectando el

controlador que tiene incluido. Posteriormente se inserta una senal del tipo escalon

como la entrada de flujo de combustible y lograr de esta manera la identificacion.

3.2 Identificacion 57

Para finalmente, usandando el Toolbox−MatlabTM obtener un modelo reducido

de segundo orden que se mostrara en esta seccion.

3.2 Identificacion

El primer paso para disenar cualquier esquema de control es la obtencion de un

modelo de la turbina a gas combustible en lazo abierto que aporte la informacion

correspondiente al comportamiento dinamico de la planta, para proceder al ajuste

de los parametros de los diversos controladores, de tal manera que satisfagan los

requerimientos de diseno establecidos.

Figura 3.1: Un modelo de turbina HTDU en laso abierto sin compensar

Para el caso particular de este proyecto fue un tanto laborioso lograr el

modelo en lazo abierto debido a que el modelo de referencia programado en

MatlabTM − Simulink fue facilitado con un gobernador incorporado cuya estructura

interna contempla un controlador tipo PI.


Antes de desconectar el bloque del gobernador se observo el valor estacionario de

la variable del flujo de combustible para luego insertarlo, a traves de una senal del tipo

escalon en cadena directa con el bloque de la turbina logrando los siguientes resultados.

Figura 3.2: Senales entrada/salida del modelo sin compensar

Observando la senal de entrada/salida en la figura (3.2) es claro notar una

respuesta libre de la accion de control y apta para ser identificada. Debido a la no

linealidad del proceso, la identificacion se realiza entorno a un punto de operacion en

el cual se estabiliza el sistema.

3.2.1 Obtencion del Modelo Parametrico

Este metodo consiste en usar la tecnica de mınimos cuadrados la cual, a traves de

cualquier herramienta matematica, puede ser implementada para la obtencion de un

modelo.

El principio de mınimos cuadrados, segun Gauss dice que: los parametros


desconocidos de un modelo se deben elegir de modo que ”la suma de los cuadrados

de las diferencias entre los valores observados realmente y los observados estimados,

multiplicada por numeros que midan el grado de precision, sea un mınimo”.

Este metodo puede aplicarse a una gran variedad de sistemas con excelentes

resultados.

Figura 3.3: Herramienta System Identification Toolbox de MatlabTM utilizada

El procedimiento aplicado, para identificar el sistema con la ayuda de MatlabTM

es:

1. Importar los datos.

2. Preprocesamiento: quitar la medias y separa los datos de identificacion y control.

3. Estimacion del modelo: seleccion del orden y especificacion de la estructura del

modelo.


Identificacion del Sistema

Por medio de la herramienta de identificacion de MatlabTM , utilizando un

periodo de muestreo de 0.001[s], 1100 muestras del modelo original, insertarlas al

System Identification Toolbox MatlabTM como vectores de entrada-salida dadas las

referencias 0.2186 [kg/seg] para el flujo de combutible y de 1660.9 [rad/s] para la

velocidad de la turbina como se ilustra en la figura (3.3) se obtuvieron los siguientes

resultados:

1. Metodo ARX: Para este metodo el modelo que se ajusto a la salida muestreada

con mayor porcentaje fue de grado 4, posee 4 ceros, 1 polos y ningun retardo,

este se ajusto con un 96.36 porciento.

2. Metodo OE: En este metodo el modelo que se ajusto con mayor porcentaje a la

salida muestreada fue de grado 2, posee 2 ceros, 1 polos y ningun retardo, este

se ajusto con un 97.86 porciento.

3. Metodo ARMAX: Para este metodo el porcentaje a la salida medida fue de grado

2, este siguio la senal en un 93.49 porciento.

4. Metodo Box-Jenkins(BJ): Para el metodo BJ el modelo con mas porcentaje de

ajuste a la salida medida fue de grado 2, el cual se aproximo a la respuesta real

del sistema con un 99.12 porciento.

La comparacion entre las respuestas adquiridas por los diversos modelos

parametricos resulto ser mejor el metodo Box-Jenkins (BJ) y la funcion de

transferencia obtnida es la ecuacion (5.1) y se demuestra en la siguiente figura (3.4):


Figura 3.4: Metodos Parametricos usando el System Identification Toolbox MatlabTM

Funcion de transferencia obtenida por el metodo (BJ):

G(s) =10025s + 41410

s2 + 9.0138s + 5.4382(3.1)

3.2.2 Validacion

De los 4 modelos identificados, es necesario seleccionar el que mejor se ajuste para

implementar la tecnica de control correspondiente; por tal razon se evalua el error

entre la respuesta del modelo identificado y los datos reales.

Se tomaron dos de los cuatro modelos por ajustarse con mayor porcentaje de

exactitud y poseer funciones de transferencia de orden 2, estos son los modelos (OE) y

el(BJ), para posteriormente realizar comparaciones las cuales arrojaron los siguientes

resultados:


Metodo Punto operacion Error Porcentual

BJ 1664.5456 [rad/seg] 0.00084108

OE 1664.5456 [rad/seg] 0.0033818

Tabla 3.1: Tabla de error porcentual de los metodos en el punto de operacion

Para validar el modelo se han variado los valores de la entrada escalon

incorporada para saber en que medida los modelos soportan los cambios.

Variacion de la referencia por debajo del punto de operacion.

Metodo 1ra Var Porc. 1er Error Porc. 2do Var. Porc. 2do Error Porc.

BJ 5 4.1418 10 8.6561

OE 5 8.6561 10 8.66

Tabla 3.2: Variacion porcentual de la referencia por debajo del punto de operacion

Variacion de la referencia sobre el punto de operacion.

Metodo 1ra Var. Porc. 1er Error Porc. 2do Var. Porc. 2do Error Porc.

BJ 5 3.3737 10 7.8316

OE 5 3.3737 10 7.8270

Tabla 3.3: Variacion porcentual de la referencia por arriba del punto de operacion


Las siguientes graficas ilustran los resultados.

Figura 3.5: Resultados de los Metodos Parametricos del System Identification Toolbox MatlabTM

Figura 3.6: Error Absoluto entre el modelo real y el identificado


Figura 3.7: Error Porcentual entre el modelo real y el identificado

A traves, de las figuras 3.5, 3.6 y 3.7, se observan los resultados de los metodos

parametricos, error absoluto y el error porcentual del sistema, donde se aprecia

claramente una gran diferencia en el regimen transitorio de las curvas, debido a que en

la identificacion se tomaron datos solo de un pequeno trozo de la curva transitoria y el

resto fue en regimen permanente lo que implica un error muy grande para el transitorio

y mınimo para el estado permanente.

Capıtulo 4

Diseno y Analisis de Resultados de

las Estrategias de Control Aplicadas

En el siguiente capıtulo se mostrara el proceso matematico para el diseno de observador,

diseno de controladores, resultados y analisis correspondientes a los algoritmos

mencionados en el capıtulo 2.

4.1 Introduccion

Este capıtulo busca probar la eficiencia y eficacia de los metodos matematicos

previamente mencionados, que permitan poner en funcionamiento un sistema de

turbinas a gas combustible que, sin arriesgar vidas humanas, logre realizar diversas

tareas.

Para ello se usan diversas tecnicas de compensacion, siendo una buena

herramienta el ajuste de ganancias de un sistema de control, y ası poder cumplir con

las especificaciones dadas.

Tambien se expondran los diversos calculos, respuestas y analisis de los sistemas

compensados con los diferentes metodos, a fin de apreciar las ventajas que ofrecen y de

elaborar una comparacion general para determinar cual es el mejor enfoque de diseno

4.2 Diseno de un Observador de Orden Completo 66

en cuanto a desempeno ante las diversas variaciones de los valores de referencia.

4.2 Diseno de un Observador de Orden Completo

Considerese el sistema de turbina de gas combustible que fue expuesto en el capıtulo

1. Este sistema esta descrito por el siguiente espacio estado proveniente de la ecuacion

(4-1):

[x1

x2

]=

[−9.0138 −5.4382

1 0

] [x1

x2

]+

[1

0

]· u (4.1)

y = 10025x1 + 41410x2 (4.2)

Es conocido que todo polo cuya parte real negativa es ai = ζωn, permite influenciar

el tiempo de establecimiento del sistema pues ts = 4ζωn

(criterio de 2 porciento de

presicion), de tal modo que el observador puede ser mas rapido o mas lento que el

sistema en estudio, segun que tan cerca esten los polos del polinomio deseado al origen.

Para tal sistema es necesario que el observador sea mas rapido que cualquier

accion de control, debido a que se debe estimar primero los estados antes de pasar

por el controlador; la aplicacion inversa de esta acotacion no funciona. Para el caso

senalado en este trabajo se escogieron valores de polos no muy cercanos al origen para

evitar esta clase de problemas en el futuro.

A partir del espacio estado, obtenemos el diagrama analogico necesario para la

simulacion en MatlabTM − Simulink. de las variables del sistema (figura 4.1)


Figura 4.1: Diagrama analogico del espacio estado identificado en MatlabTM − Simulink

Apreciando la salida del sistema se puede concluir que la estimacion debe

hacerse de los dos estados debido a que la salida del modelo identificado es

y = 10025x1 +41410x2, por lo tanto no se tiene ningun estado, lo que implica elaborar

un observador de orden completo.

Una vez aclarado este topico, obtenido el diagrama analogico y las ecuaciones,

se inicia con los calculos para el diseno del observador.

Usando la ecuaciones del capıtulo 2, tenemos:

x = Ax + Bu + Ke(y − Cx) (4.3)

en donde:

x=Representa el vector de estado aproximado

Por lo tanto el error de observacion es:

Ae = |A−KeC| (4.4)

L(s) = |sI − Ae| (4.5)


|sI−A+Ke| =

(s 0

0 s

)−

(−9.0138 −5.4382

1 0

)+

(Ke1

Ke2

) (10025 41410

)

|sI − A + Ke| = s2 + s(41410Ke2 + 9.0138 + 10025Ke1)+

9.0138Ke241410 + 5.4382− 5.4385Ke210025 + 41410Ke1 (4.6)

El polinomio caracterıstico deseado se describecomo:

Pd(s) = s2 + (a1 + a2)s + a1a2 (4.7)

donde:

a1 = 3 (4.8)

a2 = 8 (4.9)

Por lo tanto, queda:

Pd(s) = s2 + 11s + 24 (4.10)

Igualando terminos entre el polinomio caracterıstico y el polinomio deseado para

el calculo de ganancias del observador, se obtiene:

41410Ke2 + 9.0138 + 10025Ke1 = 11 (4.11)

9.0138Ke241410 + 5.4382− 5.4382Ke210025 + 41410Ke1 = 24 (4.12)

Resolviendo, se tiene que:

Ke1 = 1.1078 · 10−4 (4.13)

Ke2 = 2.1147 · 10−5 (4.14)

4.3 Diseno de Control por Seguimiento 69

Figura 4.2: Diagrama analogico del observador de orden completo en MatlabTM − Simulink

Los valores de ke1 y ke2 calculados e insertados en el diagrama analogico que se

ilustra con la figura 4.2, proporcionan una dinamica rapida y adecuada para el sistema,

debido a la posicion de sus polos descritos en las ecuaciones (4.8) y (4.9), de tal manera

que si se modifican el sistema podrıa ser mas lento, rapido o no cubrirıa adecuadamente

los requerimientos deseados.

4.3 Diseno de Control por Seguimiento

El objetivo de este control es poder seguir una senal de referencia, incrementar la

velocidad de respuesta del proceso y lograr una curva de suave, tambien cabe mesionar

que al implementar este enfoque el sistema se extiende a tercer orden y por lo tanto

se deben fijar 3 polos.

Para la seleccionar el polo a1 se tomo en cuenta la dinamica del observador,

para ello se ubica el polo deseado cerca del origen, de esta manera influenciar la


estabilidad del sistema y la velocidad del mismo mientras que el polo a2 esta alejado

por lo menos 8 veces para no afectar de manera significativa el primer polo y a su vez

el polo a3 esta 5 veces lejos, lo que permite al primer polo ser el dominante.

Usando la ecuaciones (2.33) y (2.34) del capıtulo 2, tenemos:

e(t) = (Ae− BK)e(t)

Con la ecuacion del error se calcula los valores de K , mediante:

Pd(s) = |sI − Ae + BK|s 0 0

0 s 0

0 0 s

−−9.0138 −5.4382 0

1 0 0

10025 41410 0

+

1

0

0

(K1 K2 −Ki

)

Pd(s) = s3 + s2(K1 + 9.0138) + s(10025Ki + K2 + 5.4382K2) + 41410Ki (4.15)

El polinomio caracterıstico deseado es:

Pd(s) = s3 + (a1 + a2 + a3)s2 + (a1a2 + a1a3 + a2a3)s + a1a2a3 (4.16)

donde:

a1 = 0.8086 (4.17)

a2 = 6.7085 (4.18)

a3 = 4.043 (4.19)

Por lo tanto, queda:

Pd(s) = s3 + 11.56s2 + 35.816s + 21.931 (4.20)

Igualando terminos entre el polinomio caracterıstico y el polinomio deseado para

el calculo de ganancias del control por seguimiento se obtiene:

K1 + 9.0138 = 11.56 (4.21)


41410Ki = 21.931 (4.22)

10025Ki + K2 + 5.4382K2 = 35.816 (4.23)

a partir de las cuales se obtiene:

Ki = 5.2961 · 10−4 (4.24)

K1 = 2.5471 (4.25)

K2 = 6.3877 (4.26)

La siguiente figura ilustra las conexion de un sistema de seguimiento

implementado con un integrador en MatlabTM − Simulink.

Figura 4.3: Digrama del sistema de seguimiento en MatlabTM − Simulink

Los resultados obtenidos para este esquema de control se ilustran en las

figuras 4.4, 4.5 y 4.6. Para poder apreciar los beneficios obtenidos por el sistema de

seguimiento se realizo una comparacion entre las graficas del sistema compensado y

sin compensar.


a. Valor de la referencia: 1466,1 [rad/s].

Figura 4.4: Entrada/Salida del sistema compensado (control por seguimiento), referencia: 1466,1 [rad/s]

b. Valor de la nueva referencia: 1660,5 [rad/s].


4.4 Diseno del Control Optimo Cuadratico 73

c. Valor de la nueva referencia: 1780,2 [rad/s].


Al comparar la respuesta simulada de la planta y la obtenida al implementar

el controlador mediante el sistema de seguimiento, se observa un desempeno mas

eficiente en la curva compensada y que manipula perfectamente la variable de flujo de

combustible (ff) para lograr cualquier valor de referencia.

Tambien se puede observar un pequeno sobresalto en el regimen transitorio

debido a la accion integral, pero este percance puede ser corregido a posteriori por el

control anti-rebote.

4.4 Diseno del Control Optimo Cuadratico

Este enfoque persigue un objetivo similar al esquema anterior excepto por las ventajas

que implica la aplicacion del mismo, las cuales garantizan la estabilidad basado en

ındices de desempeno cuadraticos.


Para este algoritmo, los criterios de diseno, deben cumplir las siguientes

condiciones:

Parametros de diseno:

Q=matriz hermitiana o simetrica real definida positiva (o semidefinida positiva).

R=matriz hermitiana o simetrica real definida positiva (o semidefinida positiva).

se tamara en cuenta la ecuacion de riccati (2.46) indicada en el capıtulo 2:

AT P + PA− P (BR−1B)T P + Q = 0

Haciendo uso de la herramienta MatlabTM se obtendra la matriz K, de valores

de ganancias correspondientes, a traves de la expresion (2.47).

Procedimiento general, para el diseno:

1. Declarar en la consola de MatlabTM las matrices Q y R correspodientes a los

parametros diseno.

2. Utilizar el comando Care que se encarga calcular una solucion unica X para la

ecuacion algebraica de Riccati en tiempo continuo.

[X, L,G] = CARE(A, B, Q, R, S, E)

Cuando los parametros R, S y E son omitidos, se asignan valores prefijados

R = I, S = 0 y E = I.

en donde:

X=Matriz de estados.

L=Matriz de autovalores.

G=Matriz de ganancias.

3. Modificar los valores de Q y R hasta lograr las especificaciones deseadas.

Para la diseno de control optimo en este trabajo se elaboraron una serie de

pruebas con diversas matrices Q y R hasta que se logro conseguir los valores indicados


para lograr resultados que permiten una respuesta sin sobresaltos del tipo amortiguada

y con un error que tiende a cero. Para lo cual se utilizaron las siguientes matrices:

Q =

(1 0

0 1

)(4.27)

R =

(900 0

0 900

)(4.28)

Para tales matrices se obtuvo una solucion de MatlabTM de K = [0.0720 · 10−3,

0.0720 · 10−3, 0.1022 · 10−3], donde:

Ki = 0.0720 · 10−3

K1 = 0.0720 · 10−3

K2 = 0.1022 · 10−3

los resultados se ilustran en las figuras 4.7, 4.8 y 4.9.


Figura 4.7: Entrada/Salida del sistema compensado (control por optimo), referencia: 1466,1 [rad/s]

4.5 Diseno de Controladores Clasicos 77

4.5 Diseno de Controladores Clasicos

En esta seccion se implementan redes adelanto-atraso de fase siendo una buena tecnica

para ajustar ganancias de un sistema de control, y de tal manera que cumpla con las

especificaciones que establecidas por el disenador. Tambien se mostrara el diseno de

una estructura de control que es casi universalmente utilizada en la industria. Se trata

de la familia de controladores de estructura fija llamada la familia de controladores

PID.

4.5.1 Diseno de Controlador por Atraso de Fase

La funcion primaria de un compensador en atraso consiste en atenuar la curva de fase

en el rango de alta frecuencia para dar al sistema suficiente margen de fase.

Usando las ecuacion (2.60) del capıtulo 2, tenemos:

Ked = Kc1Ke

Observando la constante del error Ke = 7614, 7 perteneciente a la planta sin compensar

se tomo el siguiente valor de Ked = 0.01 debido a que mejora el margen de fase y el

error en estado estable en un 95 porciento, por lo tanto:

Kc1 =Ked

Ke

=13.3958

7614, 7= 1.7592 · 10−3


La siguiente grafica ilustra ambos casos:

Figura 4.10: Diagrama de Bode: a la izquierda la planta sin compensar mientras que a la derecha esta

compensada con Kc1

Para desplazar la curva de magnitud hasta los cero decibeles y obtener una nueva

frecuencia de cruce de ganancia para mejorar el desempeno y la estabilidad del

sistema, es necesario utilizar la expresion: 20log 1β

= -20 decibeles, por lo tanto:

β = 10 (4.29)

Seleccionando el aporte de fase basado en la observacion y la experiencia de diseno, se

tiene:

θd + eφ = 50, 2◦ + 10◦ = 60, 2◦ (4.30)

A partir de este valor de margen de fase, se ubica en el diagrama de bode la frecuencia

de cruce deseada (ωgd):

ωgd = 2.22rad

s. (4.31)

Se elige el cero del compensador (frecuencia de cruce en 1/T ) una decada por debajo

de la frecuencia de gananacia de tal modo que resulta:


1

T=

ωgd

10= 0.222 (4.32)

Luego se calcula:

Kc =Kc1

β= 1.7592 · 10−4. (4.33)

Finalmente se calcula la funcion de transferencia del controlador:

Gc1(s) = Kc(s + 0.222)

(s + 0.0222). (4.34)

Debido a que el compensador no fue suficiente para llegar al desempeno deseado se

anade otro controlador por atraso.

Se asume el valor de Kc2 = 1, para elegir el siguiente aporte:

θd + eφ = 70◦ + 10◦ = 80◦ (4.35)

para una frecuencia deseada de:

ωgd = 2.5549rad

s. (4.36)

Se decrementa en una decada la frecuencia, obteniendo el cero del controlador por

atraso.

1

T=

ωgd

10= 0.2554 (4.37)

Finalmente la funcion de transferencia del segundo controlador es:

Gc2(s) =(s + 2.554)

(s + 0.2554)(4.38)

el Producto en cascada de los controladores resulta:

Gc(s) = 1.7592 · 10−4 (s + 0.222)

(s + 0.0222)· (s + 2.554)

(s + 0.2554)· (4.39)

Los resultados obtenidos al implantar este control se presentan en las figuras 4.10, 4.11

y 4.12.



Figura 4.11: Entrada/Salida del sistema compensado (control por atraso de fase), referencia: 1047,2

[rad/g]



[rad/s]




[rad/s]

Como se observa en las figuras anteriores, es claro notar una ganancia en tiempo

de la respuesta, la curva es mas suave y la manera eficiente como sigue la referencia,

pese a estas ventajas se paga con un sobresalto de magnitud aceptable.

4.5.2 Diseno de un Controlador Proporcional Integral (PI)

El principio de este metodo consiste en localizar la frecuencia de corte del controlador

en un valor que permita lograr ciertas especificaciones de error en estado estable.

En primer lugar se dibuja el diagrama de bode no compensado (figura 4.14).


Figura 4.14: Diagrama de bode de la planta sin compensar

Observando la figura anterior se puede notar claramente que para una fase de

θ = −60, 2 se tiene una frecuencia de cruce en el valor ωgd = 15, 6 · rads

por lo tanto

usando las ecuaciones (2.80) y (2.82) del capıtulo 2, tenemos:

Kp = 1/|Gp(ωgdj)| (4.40)

Kp = 1/(2.3377 · 103) = 6.6731 · 10−4 (4.41)

Por lo tanto, podemos obtener:

Ki =Gp(ωgdj)

10Kp

Ki =15, 6

106.6731 · 10−4 = 4.2777 · 10−4 (4.42)

La funcion de transferencia del controlador resultante es:.

Gc(s) =Ki + Kps

s=

(4.2777 · 10−4 + 6.6731 · 10−4s)

s· (4.43)


En las figuras 4.15, 4.16 y 4.17 y se muestran los resultados de la simulacion del sistema

con este tipo de control.


Figura 4.15: Entrada/Salida del sistema compensado (control proporcional integral (PI)), referencia:

1466,1 [rad/s]



1660,5 [rad/s]




1780,2 [rad/seg]

La respuesta para este control es satisfactoria debido a que, en estado estable y

para cualquier valor de referencia, el error tiende a cero. Tambien la respuesta es mas

rapida, aunque estan los sobresaltos por error de integracion no son un problema de

magnitud debido a que se pueden corregir con un sistema anti-rebote.

4.5.3 Diseno de Controlador Proporcional Integral (PID)

Para iniciar con el diseno de este controlador se debe tener en cuenta que se requiere

disenar un controlador del tipo PD, puesto que el PID es la suma de los enfoques: PD

para regimen transitorio y PI para regimen estacionario.

Diseno de Controlador PD

El principio de este metodo consiste en localizar la frecuencia de corte del

controlador en un valor que permita lograr mejoras en el regimen transitorio como se

explica en la seccion (2.4.2).


Con el diseno del PI se consiguio satisfacer mejoras de error en estado estable

de un 95 porciento. Es por ello que utilizaremos el mismo valor de Kp. En la figura

(4.18) se muestra el diagrama de bode compensado con dicho valor de Kp.

Figura 4.18: Diagrama de Bode: a la izquierda, la planta sin compensar; a la derecha, esta la planta

compensada con Kp

Luego, se determina la contribucion angular del controlador a partir del valor

de margen de fase mostrado en la grafica anterior.

θc = θd − θp + eφ = 130◦ − 109◦ + 9◦ = 30◦ (4.44)

Dado que θc < 45◦, se elige un valor de frecuencia deseada a la derecha de la frecuencia

de cruce del sistema.


Figura 4.19: Diagrama de Bode compensado con Kp

En la figura anterior se puede observar una fase de θ = −80, 1 tenemos una

frecuencia de cruce en el valor ωgd = 50, 6 · rads

Usando la ecuacion siguiente:

ωgd =Kp

Kd

(4.45)

tenemos:

Kd =Kp

ωgd

(4.46)

Kd =6.6731 · 10−4

50.6= 1.3188 · 10−5 (4.47)

Una vez obtenido los dos controladores, se pueden calcular los parametros del

PID con las ecuaciones (2.85), (2.86), (2.87) y (2.88) del capıtulo 2:

Kp = Kp2 + Kd1Ki2

Kp = 6.6732 · 10−4 (4.48)


Kd = Kd1Kp2

Kd = 8.8005 · 10−9 (4.49)

Ki = Ki2

Ki = 4.2777 · 10−4 (4.50)

Finalmente la funcion de transferencia del controlador es:

Gc(s) = Kp +Ki

s+ Kds

Gc(s) = 6.6732 · 10−4 +4.2777 · 10−4

s+ 8.8005s · 10−9· (4.51)

Debido que, al implementar el controlador en el modelo real de la turbina

no responde satisfactoriamente y genera errores de procesamiento, se ajustan los

parametros de diseno hasta lograr la siguiente funcion de transferencia con la cual

si alcanzar la respuesta deseada.

Gc(s) = 1.4620 · 10−5 +0.086

s+ 2.4854s · 10−10 (4.52)

Los resultados obtenidos al implantar este control se presentan en las figuras 4.20, 4.21

y 4.22.


Figura 4.20: Entrada/Salida del sistema compensado (control proporcional integral derivativo(PID)),

referencia: 1466,1 [rad/s]

4.6 Diseno de Esquema Anti-Rebote (Anti-Windup) 89

Como se aprecia en la figura anterior, la respuesta del sistema controlado en el

regimen transitorio muestras mejorıa, tales como: suavidad en la curva de respuesta, el

error estacionario tiende a cero y llega adecuadamente a la referencia con un pequeno

sobresalto pero no puede ocasionar danos fısicos relevante a la planta real por la sencillez

de su magnitud, ademas de poder corregirse con un sistema anti-rebote.

4.6 Diseno de Esquema Anti-Rebote (Anti-

Windup)

Para el siguiente diseno los parametros fueron obtenidos de manera empırica.

Implementando el esquema presentado en la figura 2.6 y siguiendo los casos:

1. Caso 1: Sin compensacion Kcr, ni saturacion.

Figura 4.23: Entrada/Salida del sistema compensado (Gobernador sin antirebote) y la plata sin

compensador

En donde:

Kcr= Ganancia de realimentacion del esquema antirebote.

Como se observa en la figura 4.23 hay un sobresalto producto del control PI al

integrar el error, en lo que sigue se aplicara el anti-rebote con diversos valores


hasta obtener los resultados deseados de la respuesta.

2. Sin compensacion Kcr con saturacion de [0.1,-0.1].

Figura 4.24: Entrada/Salida del sistema compensado (Gobernador con saturacion [0.1,-0.1] sin Kcr)

3. Sin compensacion Kcr con saturacion de [0.2,-0.2]

Figura 4.25: Entrada/Salida del sistema compensado (Gobernador con saturacion [0.2,-0.2] sin Kcr)

La respuesta de la velocidad saturada segun la figura 4.25 no muestra mejora,

siendo esta igual a la senal sin saturar.

4. Con compensacion de Kcr = 8, saturacion [0.1,-0.1]


Figura 4.26: Entrada/Salida del sistema compensado (Gobernador con saturacion [0.1,-0.1] con Kcr=2000)

En esta experiencia se calculo el valor de Kcr = 1Ki

= 2000, que al aplicarse no

se obtuvieron resultados satisfactorios (figura 4.27).

5. Con compensacion de Kcr = 8, saturacion [0.1855,-0.1855]

Figura 4.27: Entrada/Salida del sistema compensado (Gobernador con saturacion [0.1855,-0.1855] con

Kcr=8)

Observando la figura anterior la compensacion con los parametros establecidos

logra anular el error de integracion, limita la accion de potencia acoplada al eje

de la turbina y ofrece una respuesta mas rapida de establecimiento.


6. Con compensacion de Kcr= 8, saturacion [0.19,-0.19]

Figura 4.28: Entrada/Salida del sistema compensado (Gobernador con saturacion [0.19,-0.19] con Kcr=8)

La figura 4.28 indica que no hay cambios de ningun tipo en la respuesta del

sistema, lo que implica que los valores implementados no funcionaron.

Estudiando las graficas se demuestra que los parametros eficientes para aplicar

el algoritmo antirebote son los utilizados en el punto 5, debido a que elimina el error

en estado estable y mejora la velocidad de respuesta de manera muy aceptable.

4.7 Analisis y Comparacion entre los Enfoques previamente mencionados 93

4.7 Analisis y Comparacion entre los Enfoques

previamente mencionados

Esta seccion pretende ilustrar mediante graficas los resultados correspondientes a cada

controlador donde se observan las diversas respuestas de la velocidad de la turbina ωt

en rads

.

a. Para un valor de referencia de 13500 r.p.m o 1413,7 rads

.

Figura 4.29: Respuesta del sistema, comparacion entre controladores con referencia = 1413,7 rads


b. Para un valor de referencia de 13800 r.p.m o 1445,1 rads

.


c. Para un valor de referencia de 15860 r.p.m o 1660,9 rads

.



d. Para un valor de referencia de 17000 r.p.m o 1780,2 rads

.


e. Para un valor de referencia de 18000 r.p.m o 1885 rads

.

Figura 4.33: Respuesta del sistema, comparacion entre controladores con referencia = 1885 rads


Como se puede apreciar en las graficas anteriores, es sencillo clasificar al control

optimo como el unico que reune las siguientes carcaterısticas:

1. El error tiende a cero en estado estable para las variaciones de la referencia.

2. La curva no tiene sobresaltos siendo del tipo sobre-amortiguado e igual para las

diversas variaciones de la referencia.

3. Mejora el tiempo de respuesta para todos los cambios de la referencia.

Capıtulo 5

Implantacion de Algoritmos de

Control Bajo PLC.

Este capıtulo muestra como debe ser la implantacion de los algoritmos de control

usando la plataforma logix 5000 y estableciendo la comunicacion a traves de MatlabTM

con el Xpc y los puertos seriales RS − 232 o TCP/IP .

5.1 Introduccion

Para la implantacion de los algoritmos desarrollados podrıa usarse una plataforma

logix 5000 o un dispositivo para la comunicacion serial RS232 - TCPIP a traves del

Xpc para la transmision de datos, estos dispositivos son de suma importancia debido

a que permiten supervisar un proceso en tiempo real, esta seccion pretende servir de

guıa, explicando como seria la aplicacion real de los enfoques disenados en el capıtulo

4 e implementarlos bajo un PLC y ası controlar el desempeno de una turbina de gas

combustible.

Para elaborar tal implantacion se debe seleccionar mediante simulaciones,

pruebas y ajustes los algoritmos de control que mejor desempeno logren, para

luego discretizarlos y traducirlos a los lenguajes de programacion manejados por el

5.2 Transmision a traves de MatlabTM usando comunicacion serial RS-232 y/o TCP/IP.98

PLC debido a que este es un dispositivo digital que ademas posee 4 lenguajes de

programacion diferentes, los cuales son:

1. Funciones de Diagramas de Bloque: Aplicaciones donde se produce flujo de

informacion entre componentes de control.

2. Lenguaje de alto nivel: Similar a Pascal, es utilizado para desarrollar expresiones

complejas con un amplio rango en el manejo de datos.

3. Listado de instrucciones: Lenguaje de bajo nivel, basado en funciones de algebra

logica en donde solo es permitida una operacion por lınea.

4. Software de Programacion: La plataforma Logix de Rockwell Automation trabaja

con un software de programacion comun denominado RSLogix 5000, el cual

brinda herramientas para la programacion de aplicaciones a traves de lenguajes

estandares.

Con la implantacion de los algoritmos se pretende controlar la velocidad de

respuesta de la turbina y mejorar el desempeno de la misma.

5.2 Transmision a traves de MatlabTM usando

comunicacion serial RS-232 y/o TCP/IP.

El Xpc TargetBoxTM es distribuido por MathWorks. Es una PC Industrial con una

configuracion personalizada que provee una solucion de hardware de puesta en marcha

inmediata para las mas populares aplicaciones de prototipaje rapido.


5.2.1 Xpc

Este dispositivo provee a un sistema de tiempo real integrado, basado en un PC para

la ejecucion de codigo generado automaticamente desde Real-Time Workshop (The

MathWorks 2007).

Figura 5.1: Xpc TargetBoxTM

El Xpc soporta dos protocolos de conexion y comunicacion entre un computador

origen y destino: conexion serial y por internet. En la conexion serial la computadora

origen y destino estan conectadas directamente con un cable serial utilizando los puertos

RS232. Este cable es conectado a un enlace de MODEM nulo que puede ser de mas

de 5 metros y con una tasa de transferencia entre 1200 y 115200 baudios.

5.2.2 Conexion de la Red para un Xpc

Las computadoras origen y destino son conectadas a traves de una red. La red puede

ser del tipo LAN , INTERNET o una conexion directa utilizando un cable cruzado

para ETHERNET .


Los computadores origen y destino son conectados a la red con tarjetas para

ETHERNET utilizando el protocolo de comunicacion TCP/IP . La transferencia de

datos puede ser de 10 megabits por segundo o 100 megabits por segundo.

Figura 5.2: Conexion entre PC, a traves de comunicion serial RS232 y TCPIP

Debe instalarse un software destino Xpc en el computador de origen, pero no

es necesario. El software se puede bajar desde la web. Antes de instalar el Xpc, se

debe obtener un archivo de licencia valido o una clave personal.

El archivo de licencia o la clave de la licencia personal identifica los productos

comprados a MathWorks y los paquetes permitidos para uso e instalacion.

5.2.3 Propiedades del Ambiente Xpc.

Estas propiedades dan informacion al Xpc sobre los productos de software y hardware

con los que trabaja. Se pueden cambiar algunas propiedades, mientras que otras

raramente se pueden cambiar.

Procedimiento para la aplicacion del ambiente Xpc: en la ventana de comandos


de Matlab siga los siguientes pasos:

1. Escriba xpcsetup para que se abra el cuadro de dialogo.

2. La ventana del Xpc.

Figura 5.3: Ventana con propiedades de configuracion

3. Para la opcion Ccompiler: Se puede seleccionar Visual C o Watcom.

4. La opcion CompilerPath: Debe escribirse la ruta donde fueron instalados sus

compiladores de C y C + +.

5. HostTargetComm: Debe seleccionar RS − 232.

6. RS232HostPort: Permite seleccionar el COM1 o COM2 para la conexion al

computador destino. El Xpc determina el puerto COM que debe utilizar el

computador automaticamente.

5.3 Plataformas Logix 5000 102

7. RS232BaudRate: Permite seleccionar el rango de baudios de la conexion serial.

8. Una vez finalizados los cambios a las propiedades, oprimir el boton con la opcion

update.

El Xpc actualiza el ambiente de trabajo con las nuevas propiedades, no es

necesario salir y reiniciar Matlab, aun si se cambie la comunicacion entre el origen y

destino de RS − 232 a TCP/IP .

5.3 Plataformas Logix 5000

Las plataformas Logix de Allen-Bradley proporcionan una unica arquitectura de

control integrada para el control discreto de procesos, movimiento y variadores. Las

plataformas Logix proporcionan una maquina de control, un entorno de software

de programacion y compatibilidad para comunicaciones comunes a traves de varias

plataformas de hardware.

Figura 5.4: Control Logix 5555 de Allen-Bradley

Todos los controladores Logix funcionan con un sistema operativo de multitarea,

multiprocesamiento y aceptan el mismo conjunto de instrucciones en varios lenguajes

de programacion (Rockwell Automation 2005).

5.3.1 Caracterısticas de la Plataforma Logix 5000

El software RSLogix 5000 Enterprise Series esta disenado para funcionar con las

plataformas Logix de Rockwell Automation.

5.3 Plataformas Logix 5000 103

Tambien admite configuracion y programacion de ejes para el control de

movimiento, desarrollando las siguientes actividades:

• Crear y Modificar el codigo ejecutable.

• Supervisar datos en lınea.

• Configurar las comunicaciones.

• Configurar modulos I/O y dispositivos.

El software RSLinx (serie 9355) es un paquete de servidor de comunicacion

que proporciona conectividad de dispositivos a nivel de la planta para una amplia

variedad de aplicaciones. RSLinx permite que multiples aplicaciones de software se

comuniquen simultaneamente con una serie de dispositivos en muchas redes diferentes

(Jorge Vento 2005).

5.3.2 Lenguajes de Programacion

La Comision Internacional Electrotecnica (IEC) define una normativa de

estandarizacion a mediados de los anos 90, en su estandar IEC 61131, donde define 4

lenguajes normalizados (Rockwell Automation 2005).

1. Funciones de Diagramas de Bloque: Utilizados en aplicaciones donde se produce

flujo de informacion entre componentes de control, debido al manejo de los

elementos como bloques que pueden ser facilmente conectados entre sı en forma

similar a un circuito.

2. Lenguaje de alto nivel: Similar a Pascal, es utilizado para desarrollar

expresiones complejas con un amplio rango en el manejo de datos, ademas

posee instrucciones para ejecutar lazos de repeticion, estructuras condicionales y

funciones aritmeticas.

3. Listado de instrucciones: Lenguaje de bajo nivel, basado en funciones de algebra

logica en donde solo es permitida una operacion por lınea y cuyo manejo puede ser

5.4 Implatacion Digital de Algoritmos de Control 104

tedioso para grandes aplicaciones por ser similar al lenguaje ensamblador, cada

instruccion consta de un componente armonico de operacion y una direccion que

identifica elementos particulares dentro del controlador.

4. Software de Programacion: La plataforma Logix de Rockwell Automation

trabaja con un software de programacion comun denominado RSLogix 5000,

el cual brinda herramientas para la programacion de aplicaciones a traves de

lenguajes estandares, estableciendo un ambiente de programacion que eleva la

productividad, ofrece facilidad de mantenimiento y brinda estabilidad al sistema.

5.4 Implatacion Digital de Algoritmos de Control

El software permite organizar y configurar toda la informacion necesaria para controlar

el proceso en tiempo real mediante la creacion de un proyecto. La tabla 5.1 muestra

como configurar los datos para crear un proyecto nuevo.

Tipo del controlador Tipo del controlador

Tipo del Chasis 1756 10 Slot Control-Logix Chasis

Numero de Slot del controlador 5

Descripcion opcional

Nombre del controlador Control Optimo

Revision 11

Tabla 5.1: Propiedades del Proyectos

El controlador digital implementado consta de una tarea periodica con un solo

programa. El programa agendado en la tarea se compone de la rutina principal y

las subrutinas. En las subrutinas se deben implantar los algoritmos de control para

aplicar al sistema de turbinas de gas combustible, que debe ser editado en lenguaje

escalera.

5.4 Implatacion Digital de Algoritmos de Control 105

Para editar los algoritmos de control en lenguaje escalera es necesario discretizar

los controladores obtenidos lo que da como resultado la senal de control a envıar hacia

la valvula donde pasa el flujo de combustible.

Figura 5.5: Conexion con el Control Logix

Capıtulo 6

Conclusiones y Recomendaciones

Al finalizar esta tesis se han podido obtener importantes conclusiones, como el

desarrollo de algoritmos para el control de turbinas a gas combustible, que permiten

hacer comparaciones con cualquier aplicacion de control existente en la industria

logrando una solucion de facil manejo y flexible a modificaciones.

Se puede comprender la complejidad y la gran cantidad de tecnicas matematicas

que implica el desarrollo de las estrategias de control presentadas en este proyecto, por

lo cual este trabajo tambien pretende servir de aporte al estudio de los algoritmos de

control para turbinas bajo Plataformas - automatizadas y transmision de informacion

en tiempo real.

Al cabo de la construccion de los modelos de control, estos demostraron

resultados conforme a los objetivos planteados. La razon principal que conllevo a

los mismos fue la identificacion del sistema mediante la tecnica de mınimos cuadros

alrededor de un punto de operacion. A pesar que, el proceso es no lineal y los datos

no fueron extraıdos directamente de un modelo real sino de un simulador previamente

elaborado por Alstrom Company.

Las simulaciones e implementaciones en MatlabTM−Simulink de los esquemas

de control indican que cada algoritmo tiene su proposito particular y se piensa que

6 Conclusiones y Recomendaciones 107

los resultados constituyen un aporte significativo al momento de establecer cualquier

necesidad planteada.

Para la implantacion de los controladores en el PLC se menciona como deberıa

hacerce usando la tecnologıa de Rockwell Automation disponible en la escuela de

sistemas de la Universidad de los Andes, siendo una herramienta muy util debido ha

que cuenta con paquetes de software que permiten programar cualquier algoritmo

de control que dependa de operaciones matematicas como sumas, restas, productos,

divisiones, integraciones, derivaciones, incluso mas complejas, a traves de editores de

lenguaje escalera, texto estructurado, diagrama de bloques funcionales.

Otra herramienta util es el Xpc de MatlabTM para transmitir datos en tiempo

real de PC a PC por medio del puerto serial RS-232 o vıa TCP/IP, la facilidad

del mismo para adaptarse a cualquier tecnologıa actual en el campo del control de

procesos industrial es un motivo de relevancia para mencionarlo en este trabajo.

como recomendacion principal, se sugiere la toma de datos directamente

de un modelo de turbina real, para considerar otras variables que podrıan ser

importantes en la dinamica del proceso, ya que este trabajo no utilizo informacion

de una fuente real sino de un modelo previamente elaborado en MatlabTM−Simulink.

De tal manera de lograr una identificacion mas compleja para la continuacion

de este proyecto implantando los algoritmos bajo un PLC en cualquiera de sus

lenguajes disponibles, para su futuro uso en la industria o la universidad. A su

vez se propone implementar otros esquemas o enfoques de control para observar su

desempeno y ofrecer mas alternativas utiles.

Tambien se puede hacer enfasis en la implantacion de controladores bajo

Plataformas de automatizacion que manejen entornos de programacion libre, de tal

manera de evitar inconvenientes con licencias, permisos o costosas asesorıas por parte

de companıas privadas.

Bibliografıa

Benjamin C. Kuo (1996), Sistemas de Control Automatico, segunda edn, Prentice Hall.

Burden Fairs (1996), Analisis Numerico, segunda edn, Editorial Iberoamerica.

Ediciones UPC (2007), Pagina principal. http://www.edicionsupc.es/ftppublic/

pdfmostra/EM04001M.pdf.

Hebert Sira-Ramırez, R. y Francklin Rivas-Echeverrıa (1997), Control de Sistemas No

Lineales, Universidad de Los Andes.

Jorge Vento (2005), Implementacion de Tecnicas de Control Clasicas bajo Control-

Logix.

Ogata Katsuiko (1998), Ingenieria de Control Moderna, Prentice Hall, Universidad de

Minnesota.

Prado, L. y Lemus (1987), Polıticas de Mantenimiento en turbinas a Gas, Cadafe,

Vanezuela.

Rockwell Automation (2005), Controllogix plataformas. http://www.ab.com/logix/

controllogix/index.html.

Salvador Saucedo (2006), UbicaciOn de polos. http://www.prodigyweb.net.mx/

saucedo8/controlIV/ubicacion.pdf.

Silvio Simani (2002), Model-Based Fault Diagnosis in Dynamic Systems Using

Identification Techniques, PhD thesis, Universita degli Studi di Moderna e Reggio

Emilia.

http://www.edicionsupc.es/ftppublic/pdfmostra/EM04001M.pdf

http://www.edicionsupc.es/ftppublic/pdfmostra/EM04001M.pdf

http://www.ab.com/logix/controllogix/index.html

http://www.ab.com/logix/controllogix/index.html

http://www.prodigyweb.net.mx/saucedo8/controlIV/ubicacion.pdf

http://www.prodigyweb.net.mx/saucedo8/controlIV/ubicacion.pdf

BIBLIOGRAFIA 109

Silvio Simani, C. y Ron J. Patton (2002), Model-Based Fault Diagnosis in Dynamic

Systems Using Identification Techniques, Springer-Verlag.

Silvio Simani y Cesare Fantuzzi b. (2006), ‘Dynamic system identification and model-

based fault diagnosis of an industrial gas turbine prototype’, Mechatronics xxx

(2006) xxx–xxx p. 23. http://www.sciencedirect.com.

Silvio Simani y Ron J. Patton (2006), Identification and fault diagnosis of a simulated

model of an industrial gas turbine.

The MathWorks (2007), Pagina principal de matlab. http://www.mathworks.com.

Virginia Mazzone (2002), Estructuras practicas de control siso. http://iaci.unq.

edu.ar/caut1.

http://www.sciencedirect.com

http://www.mathworks.com

http://iaci.unq.edu.ar/caut1

http://iaci.unq.edu.ar/caut1

Anexo A

Codigos y Pantallas Usadas en el

Editor de MatlabTM

A.1 Codigo Fuente de la Turbina de Gas

Combustible HTDU

%file in1_htdu.m ’m’ file to initialise HTDU

%model: model htdu8.m and

later

R B Page: 26/9/96

Conditions to model dynamic run

% -------------- p l o t d a t a ---------------

Num_pts = 1e4; Dec_int = 100;

A.1 Codigo Fuente de la Turbina de Gas Combustible HTDU 111

% --------- I S O d a y c o n d i t i o n s

Tiso = 288.16; % ISO day temperature [K]

Piso = 101325; % ISO day pressure [N/m^2]

---------- table of Cp and Gamma for air -------------

table organised as [ T[K] Cp[kJ/kg] gamma ]

CpGmAir=[ 100 1002 1.401

200 1002 1.401

298.15 1005 1.400

400 1013 1.395

500 1030 1.387

600 1051 1.376

700 1075 1.364

800 1099 1.354

900 1121 1.344

1000 1141 1.336

1100 1159 1.329

1200 1174 1.323

1300 1188 1.318

1400 1200 1.314

1500 1211 1.311

1600 1221 1.308

1700 1229 1.305

1800 1237 1.302

2000 1250 1.298 ];


%---------- table of Cp and Gamma for combustion gas ---

% table organised as [ T[K] Cp[kJ/kg] gamma ]

% table below shows Cp as function of t from gate cycle

%data: note Cp

%also function of air fuel ratio

CpGmGas=[ 734 1121

861 1163

902 1162

940 1174

1080 1233

1143 1225

1190 1238

1289 1289

1332 1279

1399 1297

1561 1329

1619 1348 ];

% to smooth data fit polynomial

[p_ s_]=polyfit(CpGmGas(:,1),CpGmGas(:,2),2);

% generate new table from poly using same temperature

%range as for air

temp_=polyval(p_,CpGmAir(:,1));

CpGmGas=[CpGmAir(:,1) temp_];


% ----------- table of gamma for combustion gas -------

GmGas= [ 902 1.331

940 1.331

1143 1.314

1190 1.310

1332 1.299

1399 1.294

1561 1.286

1619 1.281];

% to smooth data fit polynomial

[p_ s_]=polyfit(GmGas(:,1),GmGas(:,2),2);

% generate new table from poly using same temperature

%range as for air

temp_=polyval(p_,CpGmAir(:,1)); CpGmGas=[CpGmGas temp_];

% cleanup

%clear p_ s_ temp_ GmGas;

% --------- gas constant for air --------------------

%Rgas = 287; % gas constant for air [J/kg/K]

% --------- C o m p r e s s o r d a t a -----------

% load date CompP1 CompN0 CompNmax CompM1 CompE1


load comp_dat;

% CompM1= matrix of compressor mass flows [kg/s]

% each column represents a particular compressor

%speed

%each row a particular compressor pressure ratio

%Pout/Pin

% CompN0= vector of compressor speeds [RPM]

%each element in n0 defines the speed associated with

%each column of m1

% CompP1= vector of pressure ratios

%each element in p1 defines the pressure ratio

%associated with each row of m1

% CompNmax= scaler 100% compresser speed [RPM]

% CompE1= polytropic efficiency [per unit]

CompM1= CompM1 *sqrt(Tiso)/Piso; % corrected mass flow

CompN0= CompN0 * (2*pi/60) / sqrt(Tiso);% corrected

%speed

%GcM = 1.4; % gamma for air at mean compressor

%temperature

TcS = .01; % compressor heat


%soak time constant [s]

%CpCm = 1000; % mean compresser Cp [J/kg/K]

% ------------- T u r b i n e D a t a --------------

% Hodge tables for turbine impedance;

%single stage turbine

Hodge1Pr=[0 1 1.003 1.006 1.02 1.041 1.07

1.09 1.113 1.1431.178 1.217 1.270

1.342 1.446 1.603 1.725 2.01e33];

Hodge1Wr=[0 0 .1 .2 .3 .4 .5 .55

.6 .65 .7 .75 .8 .85 .9 .95 .975 1

1];

Ke0 = 2.6; % Hodge efficiency curve ’shape’ coefficient

Ke1 = sqrt(Ke0)+1; % Hodge efficiency curve coefficient

EtD = .8; % Turbine efficiency at the design point

KsrD= 1.989; % Turbine design point pressure ratio

KfrD= 5.61e-4; % Turbine design point flow ratio

%CtCm = 1000; % Turbin mean Cp [J/kg/K]

%GtM = 1.4; % Turbine mean gamma

TtS =.01; % Turbine heatsoak time constant [s]

TtSic= 300; % Initial condition on heat soak

%time constant [K]

% ------ P r e s s u r i s i n g V a l v e ---------


% valve effective area tables:-

PrV = [.1 .274 .4]; % pressure ratio P2/P1

AnV = [0 5 10 15 20 25 30 35 40 45 50 55]; % valve

%angle [Deg]

AeV = [.0300 .0303 .0311 .0323 .0340 .0363

.0390 .0421 .0456 .0495

.0538 .0591 .0296 .0299

.0305 .0318 .0334 .0357

.0384 .0415 .0449

.0487 .0530 .0581 .0289

.0294 .0299 .0311 .0328

.0350 .0376 .0406

.0438 .0477 .0520 .0569]; % effct area [m^2]

%GmV = 1.3; % ratio of specific heats for air in valve

%Rv = 287; % gas constant for air in valve [J/kg/K]

CdV = 0.97; % Valve Cd

AvI = 40; % Initial valve angle, prior to time AvT [Deg]

AvF = 63.15; % Final valve angle, after time AvT [Deg]

AvT = 7; % Time at which valve angle

%changes from AvI to AvF [s]

% -------- C o m b u s t i o n C h a m b e r ------------

% temperature rise tables

TcI = [300 400 500 600 700 800 900];% Inlet air

%temperature [K]


AcfR = [.004 .016 .025]; % Air fuel ratio

TcR = [ 172 595 935

170 580 910

165 565 890

160 550 870

155 540 855

150 525 835

145 515 920];

% temperature rise liquid fuel [K]

TcR = TcR*49.0/42.6;

% correct tempreature rise: Liquid >> Gas [K]

%----- R e s i s t i v e p i p e s e c t i o n s -----

LcR1 = 1167; % loss coeff R1; compressor inlet duct loss

LcR2 = 2.28e5; % loss coeff R2; compressor o/p diffuser loss

LcR3 = 4.78e5; % loss coeff R3; combustion chamber loss

FcR4 = 4.85e-5; % flow coeff R4; cooling air bleed flow

FcR6 = 3.59e-5; % flow coeff R5; cooling air bleed flow

% ------------ inlet conditions ---------------------------

Pinlet = 1.009e5; % inlet air pressure [N/m^2]

Tinlet = 273.6; % inlet air temperature [K]

% ---- P l e n u m v o l u m e s ------------------------

% V1: i n l e t d u c t

%GmV1 = 1.4; % gamma for air in volume 1


PiV1 = 1e5; % initial pressure of air in volume V1 [N/m^2]

TiV1 = 300; % Initial temperature of air in volume V1 [K]

V1 = 9.75; % volume V1 [m^3]

% V2: c o m p r e s s o r + o u t l e t d i f f u s e r




V2 = 0.273; % volume V2 [m^3]

% V3: o u t e r c a s e o f c o m b u s t c h a m b e r




V3 = 0.032; % volume V3 [m^3]

% V4: c o m b u s t i o n c h a m b e r

%GmV4 = 1.32;% gamma for air in volume v4



V4 = 0.02; volume V4 [m^3]

V5: t u r b i n e o u t l e t d i f f u s e r

%GmV5= 1.35; % gamma for air in volume 5

PiV5 = 1e5; % initial pressure of air in volume V5 N/m^2]



V5 = 0.305; % volume V5 [m^3]

% -------------- CP of gas streams into gas junction -------

%Cp4 = 1000; % CP of bleed air flow [J/kg/K]

%Cp5 = 1300; % CP of gas in combustion chamber [J/kg/K]

%Cp6 = 1000; % CP of gas at turbine outlet [J/kg/K]

% -------------------- mechanical load ----------------------

Jl = 4.8; % total inertia [kg m^2]

Tl = 2.018e-5; % torque loss coeff [Nm]/[rad/s]^2

Ws1 = 7000*2*pi/60; % initial turbine speed at time t=0 [rad/s]

Ws2 = 10000*2*pi/60; % turbine speed at end of starter

%ramp at time t=Ts1[rad/s]

Ts1 = .5; % time at end of starter ramp [s]

Ts2 = .75; % time at end of starter hold time [s]

% --------------- gobernador ------------------------

Gi = .5e-3; % ganancia integral del gobernador

Gp = 1e-3; % ganancia proporcional del gobernador

Ts3 = 1.5; % rompiento del 2do tiempo del flujo de

%combustible[seg]

Ws4 = 11500*2*pi/60; % speed at which governor

%activated [rad/s]

Ws5 = 12286*2*pi/60;

% governor 1st speed set point [rad/s]

Ws6 = 17000*2*pi/60; %


% 15860 governor 2nd speed set point

Ts6 = 10;

% time at which 2nd speed set point activated

Fr1 = .05;

% initial fuel flow ramp rate [kg/s/s]

Fr2 = .05;

% final fuel flow ramp rate [kg/s/s]

Sdr = 1000*2*pi/60;

% speed demand rate limit [Rad/s/s]

%--------- Control por seguimiento ---------------------

TWf=10;

% tiempo de activacion del set point final

Winicio=12286*2*pi/60;

% 1ra Velocidad deseada de operacion [rad/seg]

Wfinal=17000*2*pi/60;

% 2da Velocidad deseada de operacion,

%set point final [rad/seg]

K1s=2.5471; K2s=6.3877; Kir=5.2961e-004;

%--------------Adelanto-Atraso ------------------------

Kc=1.7592e-004; beta=10; T1=4.5045; T2=0.3922;

%-------------- PID -----------------------------------


Kp=1.4620e-005; Ki=0.086; Kd=2.4854e-010;

%-------------- PI ------------------------------------

Kii=4.2777e-004 Kip=6.6731e-004;

%-------------Control Optimo ---------------------------

A=[-9.138 -5.4382;1 0];

B=[1 0;0 0]; C=[10025 41410]; D=[0];

S=eye(2)*0; QB=eye(2)*1; RB=eye(2)*900;

%QA=eye(2)*1;

%RA=eye(2)*2;

%Q1=eye(2)*1;

%R1=eye(2)*1;

%Q2=eye(2)*100;

%R2=eye(2)*100;

%Q3=eye(2)*1e-002;

%R3=eye(2)*1e-002;

%Q4=eye(2)*1e-003;

%R4=eye(2)*1e-003;

%Q5=eye(2)*1e-004;

%R5=eye(2)*1e-004;

%*************** Polinomio de Riccati


%[P,L,G] = care(A,B,QB,RB,S)

%*** Variacion de parametros Q y R de dise~no ********

%*** Variacion con R1 y Q1

%P = [0.0645 0.0912; 0.0912 1.1896];

%L = [-8.5562;-0.6462];

%G = [0.0645 0.0912;0 0];

%K1r = 0.0645;

%K1 = 0.0645;

%K2 = 0.0912;

%*************** Variacion con R1 y Q2

%P = [4.8403 5.9449; 5.9449 109.4215];

%L = [-13.1100;-0.8683];

%G = [4.8403 5.9449;0 0];

%K1r = 4.8403;

%K1 = 4.8403;

%K2 = 5.9449;


%P = [0.0006 0.0009;0.0009 0.0119];

%L = [-8.4986;-0.6400];

%G = 1.0e-003*[0.6477 0.9193;0 0];

%K1r = 1.0e-003*0.6477;


%K1 = 1.0e-003*0.6477;

%K2 = 1.0e-003*0.9193;


%P = [0.0006 0.0009;0.0009 0.0119];

%L = [-8.4986;-0.6400];

%G = 1.0e-003*[0.6477 0.9193;0 0];

%************** Variacion con R2 y Q1

%P = [0.0648 0.0919; 0.0919 1.1924];

%L = [-8.4986;-0.6400];

%G = 1.0e-003*[0.6477 0.9193; 0 0];


%P = [0.0648 0.0919;0.0919 1.1924]

%L = [-8.4986; -0.6400]

%G =1.0e-003 *[0.6477 0.9193;0 0]


%P = [0.0006 0.0009;0.0009 0.0119];

%L = [-8.4981;-0.6399]

%G = 1.0e-005*[0.6478 0.9194;0 0];


%K1r = 1.0e-005*0.6478;

%K1 = 1.0e-005*0.6478;

%K2 = 1.0e-005*0.9194;


%P = [0.0001 0.0001;0.0001 0.0012];

%L = [-8.4981;-0.6399]

%G =1.0e-006*[0.6478 0.9194;0 0];

%K1r = 1.0e-006*0.6478;

%K1 = 1.0e-006*0.6478;

%K2 = 1.0e-006*0.9194;


%X = 1.0e-003*[0.0065 0.0092;0.0092 0.1192]

%L = [-8.4981;-0.6399]

%G = 1.0e-007*[0.6478 0.9194;0 0]

%K1r = 1.0e-0014*0.6478;

%K1 = 1.0e-0014*0.6478;

%K2 = 1.0e-0014*0.9194;

%************** Variacion con RA y QA

%P = [0.0646 0.0916;0.0916 1.1910];

%L = [-8.5272;-0.6431];

%G = [ 0.0323 0.0458;0 0];

%K1r = 0.0323;


%K1 = 0.0323;

%K2 = 0.0458;


%P = [0.0484 0.0594;0.0594 1.0942]

%L = [-13.1100;-0.8683]

%G = [4.8403 5.9449;0 0]

%*************** Variacion con Q1 y R2

%P = [0.0648 0.0919;0.0919 1.1924];

%L = [-8.4986;-0.6400];

%G = 1.0e-003*[0.6477 0.9193;0 0];


%P = [0.0484 0.0594;0.0594 1.0942];

%L = [-13.1100;-0.8683];

%G = [4.8403 5.9449;0 0];


%P = [0.0246 0.0266;0.0266 1.0322];

%L = [-32.7363;-0.9802];

%G = [24.5785 26.6488;0 0];

%K1 = 14.5785;

%K2 = 16.6488;


modelo de la Turbina de Gas Combustible HTDU 126

%************************************************

%************************************************

%****************** Con QB y RB

%P = [0.0648 0.0919;0.0919 1.1924];

%L = [-8.4981;-0.6399];

%G = 1.0e-003*[0.0720 0.1022;0 0];

K1r=1.0e-003*0.0720; K1op=1.0e-003*0.0720;

K2op=1.0e-003*0.1022;

A.2 Codigo para generar las graficas de las diversas

variables que maneja el modelo de la Turbina

de Gas Combustible HTDU

file plt_t2.m plot commands for tests on HTDU_1

figure(1); clf; subplot(5,1,1),plot(t,p2),grid,

ylabel(’p2’)

title([’HTDU simulation’,’ : ’,date])

subplot(5,1,2),plot(t,p3),grid,ylabel(’p3’)






xlabel(’--------- time ------->’)

figure(2); clf; subplot(2,1,1),plot(t,p3-p5),grid,

ylabel(’p3-p4’)


subplot(2,1,2),plot(t,p4-p5),grid,ylabel(’p4-p5’)

xlabel(’--------- time ------->’)

figure(3); clf; subplot(4,2,1),plot(t,m1),grid,

ylabel(’m1’)


subplot(4,2,2),plot(t,m3),grid,ylabel(’m2’)







xlabel(’--------- time------->’)

figure(5); clf; subplot(4,2,1),plot(t,t3),grid,

ylabel(’t3’)




subplot(4,2,2),plot(t,t4n),grid,ylabel(’t4n’)

subplot(4,2,3),plot(t,t4),grid,ylabel(’t4’)






xlabel(’--------- time ------->’)

figure(6); clf; subplot(4,1,1),plot(t,qc),grid

,ylabel(’qc’)


subplot(4,1,2),plot(t,qt),grid,ylabel(’qt’)

subplot(4,1,3),plot(t,qa),grid,ylabel(’qa’)

subplot(4,1,4),plot(t,wt*60/2/pi),grid,ylabel(’Nt’)

xlabel(’--------- time ------->’)

figure(7); clf

compressor surge line

csl=[ 3 1.8

5 3

6.5 4.3

8 5.4



9.3 7

13 10.4

14 11.4

16 14

17.8 16.4

19 17 ];

plot(m3,p2.\p3),grid,ylabel(’p3/p2’),xlabel(’m3’)

hold on;

plot(csl(:,1),csl(:,2),’m’); title([’HTDU simulation:

compressor

operating point’,’ : ’,date]) clear csl;

figure(8); clf subplot(2,1,1),plot(t,ffa),grid, ;hold

on,plot(t,ffg,’r’),plot(t,ff,’k’),ylabel

(’Flujo de Combustible

[Kg/seg]’); title([’HTDU Simulacion’,’ : ’,date]);

subplot(2,1,2),,plot(t,wtg*60/2/pi,’k’),grid,;hold

on;plot(t,wta*60/2/pi,’r’)

,ylabel(’Velocidad [rad/seg]’)%,plot(t,av),

grid,ylabel(’av’);

plot(t,wts*60/2/pi,’g’),;hold on;grid,

plot(t,wtpid*60/2/pi,’k’),

plot(t,wtpi*60/2/pi,’r’)

,plot(t,wtco*60/2/pi,’c’),

ylabel(’Velocidad [rad/seg]’)%



plot(t,qc.*wt,’m’),ylabel(’pt’) ;

,plot(t,wtaa*60/2/pi,’y’)

legend(’Control Por Seguimiento’,’Control PID’,

’Control PI’,’Control

Optimo’,0);

xlabel(’--------- tiempo ------->’)

;hold on;

figure(9) title([’HTDU simulation’,’ : ’,date]);

plot(t,wt*60/2/pi),grid,ylabel(’wt’);

figure(10) title([’HTDU simulation’,’ : ’,date]);

plot(t,ff),grid,ylabel(’Flujo de Combustible’);

xlabel(’--------- time ------->’)



A.3 Codigo para generar el informe de los puntos

de operacion de las diversas variables que

maneja el modelo de la Turbina de Gas

Combustible HTDU

file printss.m

function to print HTDU steady state values;

disp(’ ’)

disp(’ ’)

disp(’ ’)

disp(’ ’)

%del htdu_ss.dry;

diary htdu_ss.dry; c_=clock;

tstr_=[’’,num2str(c_(4)),’:’,num2str(c_(5))];

disp([’HTDU simulation: steady-state values’ ,

date,tstr_]) formatshort e;

disp(’ ’) disp([’t= ’,num2str(t(length(t))),’ [s]’])

disp([’ Piso=’,num2str(Piso),’ [N/m^2]’])



disp([’ Tiso= ’,num2str(Tiso),’ [K]’])

disp(’ ’)

disp([’ p1 = ’,num2str(p1(length(p1))), ’ [N/m^2]’])







disp(’ ’)

disp([’ m1 = ’,num2str(m1(length(m1))), ’ [kg/s]’])


disp([’ m4 =’,num2str(m4(length(m4))), ’ [kg/s]’])

disp([’ m5p = ’,num2str(m5p(length(m5p))), ’ [kg/s]’])







disp([’ m10 = ’,num2str(m10(length(m10))),’ [kg/s]’])

disp([’ ff = ’,num2str(ff(length(ff))), ’ [kg/s]’])

disp(’ ’)

disp([’ t1 = ’,num2str(t1(length(t1))), ’ [K]’])


disp([’ t3n = ’,num2str(t3n(length(t3n))), ’[K]’])


disp([’ t4n = ’,num2str(t4n(length(t4n))),’ [K]’])


disp([’ t5n = ’,num2str(t5n(length(t5n))),’ [K]’])



disp([’ t5c = ’,num2str(t5c(length(t5c))),’ [K]’])


disp([’ t6 = ’,num2str(t6(length(t6))),’ [K]’])


disp([’ t7n = ’,num2str(t7n(length(t7n))), ’[K]’])

disp(’ ’)

disp([’ qc = ’,num2str(qc(length(qc))), ’ [Nm]’])

disp([’ qt = ’,num2str(qt(length(qt))), ’ [Nm]’])

disp([’ qa = ’,num2str(qa(length(qa))), ’ [Nm]’])

disp([’ qloss = ’,num2str(qloss(length(qloss))),’[Nm]’])

disp(’ ’)

disp([’ CpCm = ’,num2str(CpCm(length(CpCm))), ’

[J/kg/K]’])

disp([’ GmCm = ’,num2str(GmCm(length(GmCm))), ’

[]’])

disp([’ CpTm = ’,num2str(CpTm(length(CpTm))), ’

[Nm]’])



disp([’ GmTm = ’,num2str(GmTm(length(GmTm))), ’

[Nm]’])

disp(’ ’) pcta = CpCm.*m3.*(t3-t2);

ptta = CpTm.*m8.*(t6-t7n);

disp([’pcm = ’,num2str(qc(length(qt))*wt(length(wt))),’[W]’])

disp([’pct = ’,num2str(pcta(length(pcta))), ’ [W]’])

disp(’ ’)

disp([’ qtm = ’,num2str(qt(length(qt))*wt(length(wt))),

’ [W]’])

disp([’ ptt = ’,num2str(ptta(length(ptta))), ’

[W]’]) clear pcta ptta;

disp(’ ’)

disp([’qloss = ’,num2str(qloss(length(qloss))), ’ [Nm]’])

disp([’qloss = ’,num2str(qloss(length(qloss))*wt(length(wt))),

’[W]’])

disp(’ ’)

disp([’ wd = ’,num2str(wd(length(wd))), ’

[rad/s]’])

disp([’ wt = ’,num2str(wt(length(wt))), ’

[rad/s]’])

A.4 Codigo para la validacion mediante el error del modelo identificado de la

Turbina de Gas Combustible HTDU 136

disp([’ Nt = ’,num2str(wt(length(wt))*60/2/pi), ’

[RPM]’])

disp(’ ’)

disp([’ av = ’,num2str(av(length(av))), ’

[Deg]’])

disp(’ ’)

disp([’ CompE= ’,num2str(CompE(length(CompE))), ’

[]’])

disp([’ turbE= ’,num2str(turbE(length(turbE))), ’

[]’])

diary off;

A.4 Codigo para la validacion mediante el error

del modelo identificado de la Turbina de Gas

Combustible HTDU

clc close all

hold on plot(t,VelocidadTR(:,2))

plot(t,VelocidadTBJ(:,2),’r’)



plot(t,VelocidadTOE(:,2),’k’)

legend(’Velocidad Wt(rad/seg) Modelo

REAL’,’Velociadad Wt(rad/seg) Modelo IDENTICADO BJ’,

’Velociadad Wt(rad/seg) Modelo IDENTICADO OE’,0)

ylabel(’Velocidad (Rad/Seg)’)

title(’Respuesta de la Velocidad Wt, Turbina HTDU’)

xlabel(’Tiempo’);grid;hold off

figure(2) plot(t,(abs(wt-WtIdent)))

title(’Error Absoluto, Modelo

Real & Identificado bj, Turbina HTDU’)

xlabel(’Tiempo’) grid

figure(3) plot(t,(abs(wt-WtIdent)*100/1664.5456))

title(’Error Porcentual, Modelo Real & Identificado bj,

Turbina HTDU’)

ylabel(’Pocentaje del Error %’)

xlabel(’Tiempo’)

%grid



disp(’ ’)

disp(’ ’)

disp(’ ’)

disp(’ ’)

diary htdu_ss.dry; c_=clock;

tstr_=[’ ’,num2str(c_(4)),’:’,num2str(c_(5))];

disp([’Turbina HTDU : Valores ’,date,tstr_])

format short e;

disp(’ ’)

disp([’Valor del Set Point Wt de la Planta Real=

’,num2str(wt(length(wt))), ’ [rad/seg]’])

disp(’ ’)

disp([’Valor del Set Point Wt de la Planta Identificada=

’,num2str(WtIdent(length(WtIdent))),’ [rad/seg]’])

disp(’ ’)

disp([’Valor Maximo de Flujo de Combustible=

’,num2str(ff(length(ff))),’ [Kg/seg]’])

disp(’ ’)

disp([’Error Porcentual bj=

’,num2str(abs(wt(length(wt))-WtIdent(length(WtIdent)))

A.5 Turbina HTDU Control por Seguimiento 139

/wt(length(wt))*100 ),’ %’])

disp(’ ’)

disp([’Error Porcentual oe= ’,

num2str( abs(wt(length(wt))-

WtIdent1(length(WtIdent1)))

/wt(length(wt))*100 ),’ %’])

diary off;

A.5 Turbina HTDU Control por Seguimiento

Figura A.1: Pantalla secundaria de la turbina HTDU, estructura del Control por Seguimiento

A.6 Turbina HTDU Control Optimo 140

A.6 Turbina HTDU Control Optimo

Figura A.2: Pantalla secundaria de la turbina HTDU, estructura del Control Optimo

A.7 Turbina HTDU Control por Atraso de Fase

Figura A.3: Pantalla secundaria de la turbina HTDU, estructura del Control Atraso de Fase

A.8 Turbina HTDU Control Proporcional Integral (PI) 141

A.8 Turbina HTDU Control Proporcional Integral

(PI)

Figura A.4: Pantalla secundaria de la turbina HTDU, estructura del Control Proporcional Integral (PI)

A.9 Turbina HTDU Control Proporcional Integral

Derivativo (PID)

Figura A.5: Pantalla secundaria de la turbina HTDU, estructura del Control Proporcional Integral Derivativo

(PID)

Desarrollo e Implementaci´on de Algoritmos de Control para...

Documents

Transcript of Desarrollo e Implementaci´on de Algoritmos de Control para...