Trabajo teórico práctico con Matlab

Enviado por Leonardo Villegas

*El desarrollo y análisis realizado en el presente trabajo, a pesar de que solo abarca algunos conceptos importantes sobre el inmenso ámbito del control de procesos, engloba nociones básicas para el conocimiento y comprensión de la teoría de control y su aplicación.

1- Introducción – Control Moderno

Las aplicaciones de control automático son sinónimos de la tecnología moderna, se encuentran dentro del ámbito de la robótica hasta en un simple tostador. El control moderno aborda el problema de obtener el comportamiento deseado de un sistema que trabaja por sí solo.

El control de procesos se basa en los fundamentos de la teoría de retroalimentación, análisis de sistemas, modelización de sistemas, conversión de señales, control óptimo, teoría de seguimiento, etc.

Un sistema de control es una interconexión de componentes que forman una configuración del sistema, la cual proporcionará una respuesta deseada del mismo sistema.

Debido a la complejidad en aumento de los sistemas de control y al interés de obtener funcionamientos óptimos, se consideran cada vez un mayor número de variables del proceso. Así se tienen los sistemas multivariables en relación con indicadores para mejorar su eficiencia en el desempeño.

2- Objetivos y alcances del trabajo

El siguiente trabajo pretende introducir en conjunto, temas teóricos (vistos en clases) en relación con un ejemplo real donde se puedan desarrollar e ilustrar conceptos tales como: los de un sistemas y su clasificación; modelo matemático asociado y sus formas de representación; discretización, retroalimentación y control (vía ubicación de polos); observación de variables del sistema por medio de su salida y control óptimo cuadrático en referencia a índices de optimización.

El desarrollo de este trabajo también procura ser una guía para el lector; por lo cual, para entender los temas a desarrollar, se han hecho algunas reseñas de los contenidos expuestos en el mismo, así como también, se utilizó algunos comandos de la herramienta Matlab para el desarrollo y resolución de algunos puntos(Control System Toolbox).

3- Ejemplo a Desarrollar

Con este ejemplo se tratará de ilustrar una serie de conceptos:

Versión muestreada de un sistema continuo; estabilidad, controlabilidad y observabilidad del sistema; reubicación de polos y control óptimo.

Sistema eléctrico de doble recinto -sistema térmico- que supondremos a modo de ejemplo y como hipótesis, que éste debe mantener y ajustar las temperaturas de los recintos acotándolas en ciertos valores, teniendo en cuenta algunos parámetros físicos de los recintos. También supondremos que estos recintos contienen fluidos que no afectan al modelo propuesto y que sus propiedades físicas y químicas no varían y tampoco hay entrada o salida de fluido dentro de los recintos en el período de ensayo.

Se toman en cuenta los siguientes parámetros:

u(t)es la cantidad de calor de entrada al recinto Tj.(calentador eléctrico)
To, Tj, Ti son las temperaturas en el exterior del horno (ambiente), en el primer recinto y en el recinto interior respectivamente.
Ai, Ao: son las superficies interiores y exteriores de los recintos respectivamente.
Ci, Cj: son la capacidades caloríficas del espacio interior y del exterior

respectivamente.

hi, ho son los coeficientes de transmisión de calor para las superficie interior y exterior respectivamente.

Consideraciones

Teniendo en cuenta una serie de consideraciones físicas, se supondrá que: los valores reales de superficie, los valores para las capacidades caloríficas de los recintos y sus coeficientes de transmisión de calor, permanecerán constantes en todo momento, bajo cualquier condición.

Se tomarán como variables de estado las temperaturas de los recintos.

x1 = Ti

x2 = Tj

y u(t) la señal de control.

Las ecuaciones propuestas del sistema son:

x´1= -2×1 + 2×2

x´2= 0,5×1 -0,75×2 + 0,5u + 0,25T0(amb)

Los valores propuestos anteriores, tienen como finalidad mostrar como los recintos podrían estar interactuando con la entrada (resistencia eléctrica) que solo afecta al recinto exterior (según figura); como se ve en el modelo, los dos recintos pierden temperatura, pero parte de la energía es absorbida por el otro recinto. Como se mencionó, el recinto exterior recibe la energía de la entrada, pero del mismo modo éste recibe la energía calórica que aporta el ambiente exterior.

Considerando constante la temperatura ambiente, se puede eliminar el correspondiente término compensándolo con un término de valor -0,25Tamb.

Se tomará como salida la temperatura del recinto interior X1. Es la de mayor interés hipotéticamente.

En Matlab escribiremos las ecuaciones de estado y de salida que conforman al sistema como sigue:

A = [-2 2 ;0.5 -0.75];

B = [0 ;0.5];

C = [1 0];

D = 0; % no existe transmisión directa entre la entrada y la salida

HC = ss(A,B,C,D)

La sintaxis es sencilla, se definen cada una de las matrices con nombres; para indicar las distintas filas se coloca "punto y coma". El comando "ss" nos permite representar a nuestro sistema en el espacio de estados, quedando con el nombre propuesto "HC". Más adelante se visualiza como quedan las instrucciones en el editor de líneas de comando "el Prompt" de Matlab.

4- Clasificación del Sistema Ejemplo bajo estos supuestos

Es un sistema de tiempo continuo (analógico), ya que las variables que se toman en cuenta pueden modificar sus valores (reales -de temperatura-) en cualquier momento.

Se trata de un sistema estacionario (invariante en el tiempo), ya que su respuesta a una misma entrada es la misma independientemente del instante en que se aplica la entrada. Es decir, como se verá mas adelante al discretizar la representación del sistema, en nuestro caso, las matrices del Modelo en Espacio de Estados asociadas al sistema no están afectadas por el período de muestreo (k); esto es porque las temperaturas registradas en los recintos no depende del momento en que son evaluadas.

Además, es un sistema lineal, ya que se verifica el Principio de Superposición, es decir, que su respuesta a una combinación lineal de entradas es la misma combinación lineal de las respuestas a las entradas individuales. La condición de linealidad implica que a una combinación lineal de entradas le corresponde la combinación lineal de salidas.

Dicho de otro modo, en nuestro caso tenemos una sola entrada, la cual pudiendo variar, hará variar linealmente a la salida (por ejemplo: si midiendo las dos temperaturas correspondientes a los recintos en relación a la temperatura del calentador u(t)) siempre se tiene una misma relación-ej: la suma de las temperaturas de los recintos es la mitad de la del calentador/es-, si tuviésemos 2 calentadores (2 entradas) debería suceder igual(proporcionalmente) siguiendo el mismo principio; siempre y cuando se cumplan las consideraciones propuestas inicialmente.

Por otro lado, es un sistema causal ya que la salida en un instante depende de las entradas pasadas hasta ese instante (no depende de valores futuros). Para un sistema en tiempo discreto una condición necesaria y suficiente para causalidad es:

Causalidad ó H(n) =0 con n <0 , donde H(n) es la función transferencia de pulso(S/E) del sistema. Por eso es, en nuestro caso, trabajamos con Transformada Z unilateral. (Que corresponde a las sucesiones que tienen todos los términos de la serie de potencias positivas nulos, es decir la serie sólo está compuesta por los términos de potencias negativas y el término independiente).

La estabilidad del sistema es estudiada mas adelante. Un Sistema estable es aquel cuya salida es acotada, es decir no diverge. A una entrada acotada le corresponde una salida acotada; en nuestro caso hay que tener en cuenta que se trata de un sistema térmico, por lo que se podrían considerar factores como la velocidad y temperatura de calentamiento máximo y la respuesta de los recintos (que no excedan de cierta temperatura, etc.) como parámetros de estabilidad.

Es un sistema S.I.S.O. por tener solo una entrada y una salida, a diferencia de los M.I.M.O. que tienen múltiples entradas y múltiples salidas.

5- Discretización del sistema continuo

Por ejemplo, si discretizamos el sistema mediante aproximaciones numéricas, de modo que y si tomamos como tiempo de muestreo T igual a 1 seg. Tendríamos:

x1(k+1)-x1(k)= -2×1(k) + 2×2(k)

x2(k+1)-x2(k)= 0,5×1(k) -0,75×2(k) + 0,5u(k)

Como se ve, es fácil este método de discretización (muy inexacto por cierto) tomando periodos de muestreo igual a la unidad.

– Matlab –

Pero para nuestro caso, utilizaremos un muestro de 0,25 seg. con el método de muestreo que trabaja por defecto el Matlab 5.3, el ‘zoh’-Zero-order hold- que considera constantes a las entradas durante el periodo de muestreo; así obtenemos:

Sistema discretizado (designado "HD" en el Matlab).

x1(k+1)= 0.62774×1(k) + 0.359712×2(k)+ 0.025105 u(k)

x2(k+1)= 0,089928×1(k) +0,85256×2(k) + 0.11503 u(k)

Además, Matlab, ofrece la posibilidad de intercambiar la representación del modelo utilizando distintos comandos:

Para representar al sistema como la función de transferencia de pulso, entre la salida y la entrada, se puede utilizar el comando "tf" como sigue: H = tf(HD); en este caso, partiendo de la representación del sistema en el espacios de estados(HD) se obtiene H.

También podemos representar nuestro modelo en función de sus ceros, polos y la ganancia, con el comando zpksys = zpk(HD).

En el presente trabajo, seguiremos analizando distintos conceptos utilizando la representación en el espacio de estados discretizado.

* Aclaración: las matrices que en tiempo continuo tienen asignado los nombres A,B,C,D deberían reasignarse -bajo el análisis de tiempo discreto- como G, H, C, D respectivamente, pero siguiendo la sintaxis del Matlab se omitirá dicha reasignación.

Ahora bien ¿Por qué discretizamos el sistema?

Ventajas de la discretización de variables

La razón principal del empleo de la variable discreta es que permiten el proceso y almacenamiento de la información (datos) en computadoras digitales. Para ello finalmente se reduce la información a códigos binarios; utilizando en este caso la Transformada Z que nos permite reducir Ecuaciones en Diferencias con coeficientes constantes en Ecuaciones Algebraicas Lineales, ahorrando el manejo de eAt(exponencial) que aparece en el estudio de sistemas y variables dentro del espacio continuo.

Desventajas de la discretización de variables

En la discretización de variables [analógicas], en cual se toman registros a intervalos discretos regulares de señales (eléctricas o de otra índole); donde estos valores se llaman muestras conformando una sucesión; se observa que estas muestras introducen errores de aproximación debido a la conversión de la señal. También podrían considerarse otros errores como los producidos al utilizar código binario para representar estas muestras (este código binario depende, a su vez, de la cantidad de bits que se utilicen para la representación en una computadora).

Por otro lado, en la reconstrucción de la función original, que es el proceso inverso para transformar los valores discretos en señales (eléctricas o de otra índole) de variable continua en el tiempo (representables por números reales) ;se obtienen funciones escalonadas según el período de muestro original, pudiéndose así, reconstruir infinitas funciones con los valores obtenidos en el muestreo(ya que éstos solo representan información de la función en un momento dado, y no la función en sí).

¿Como determinar el tiempo entre cada muestreo?

El tiempo entre muestras estará determinado por el problema en sí, se deben hacer varias consideraciones, como la velocidad de procesamiento de las señales, si éstas deben ser almacenadas o no, capacidad de almacenamiento, tiempo total del muestreo (si se está haciendo un seguimiento del sistema), etc.

De todas las consideraciones, una de ellas tiene notable importancia, ésta tiene que ver con la frecuencia de los cambios o fluctuaciones que se producen en el sistema. Como regla general se debe hacer un muestreo del sistema con una frecuencia que doble (al menos) a la frecuencia máxima con la que pueden cambiar alguna de sus variables. Esto es así, para poder recolectar correctamente todos los cambios producidos en las variables, de otro modo pueden suceder inconvenientes o efectos indeseados en la lectura de estas variaciones. (Efecto Doppler, etc.)

Para nuestro caso, se ha supuesto que las temperaturas de los recintos pueden variar de un grado (centígrado) a otro, como máximo, en medio segundo; por ello es que se ha escogido como período de muestra 0.25 segundos.

6- Análisis de Estabilidad del Sistema

A continuación se realiza el estudio de estabilidad del sistema utilizando comandos del ToolBox de control de sistemas que ofrece Matlab.

Luego podemos visualizar las soluciones del sistema y observar si éstas están dentro de mi región de estabilidad (temperaturas máximas) que defina oportunamente según el caso real.

El comando "dlyap" nos permite resolver el segundo método propuesto por Lyapunov para sistemas discretos.

Siguiendo este método para el análisis de estabilidad asintótica global;

si se halla una matriz X (hermítica) de modo que sea definida positiva, habiendo definido Q (hermítica) positiva de modo que:

AT X A – X = – Q

Entonces, queda demostrado que el origen x = 0 es un estado de equilibrio.

Lo que pretende Lyapunov, es utilizar conceptos de derivadas de funciones escalares para determinar si estas son convergentes.

Al definir de forma positiva una función (x) y su derivada primera está definida negativa, estamos en presencia de una función decreciente, por lo que tiene una cota superior y si extendemos esa función a lo largo del tiempo ésta tiende a tomar el valor nulo.

En este caso la derivada es interpretada como la diferencia de valores obtenidos de cada una de las muestras.

ΔV(x(kT))= V(x(k+1)T) – V(x(kT))

Se toma como posible función de Lyapunov : V(x(k))= A*(k)XA(k) -definida positiva-.

Y se requiere que ΔV(x(k)) sea definida negativa para la estabilidad asintσtica.

– Matlab –

Para nuestro caso:

Al aplicar el criterio de Silvester para probar la definición positiva de X, se encuentra que X es definida positiva; por lo tanto, el estado de equilibrio, el origen x=0, es asintóticamente estable global.

Además puede observarse la estabilidad en este ejemplo de sistema discreto verificando los autovalores de la matriz A, estos deben ser negativos en su parte real y para sistemas discretos, estrictamente, deben ser menores a la unidad.

Se debe notar que si estamos en presencia de un sistema discretizado que representa a un sistema continuo con polos complejos; puede ocurrir inestabilidad oculta en función del período de muestreo.

La siguiente figura representa la respuesta del sistema (evaluando solo a X1) frente a un impulso con las condiciones iniciales siguientes; trabajando en Matlab:

» x0 = [20 ; 100] % suponiendo temperaturas iniciales de 20 y 100 grados en cada recinto.

» initial(HD,x0)

Se observa como el recinto X1 es afectado por la temperatura mayor del otro recinto, pero éste pierde temperatura gradualmente.

En la siguiente figura podemos ver que al aplicar en la entrada valores de magnitud constantes igual a la unidad (escalón unitario), se observa que la salida tiende a equilibrarse en el valor 2 luego de 30 seg.

» step(HD)%comando para graficar la respuesta del sistema ante una entrada escalón unitario

7- Análisis de Controlabilidad del Sistema Muestreado

-Controlabilidad de los estados

Podemos decir que el sistema es controlable si todas sus variables de estados pueden ser controladas en un período finito de pasos por una señal de control no restringida.

En otras palabras, si nos concentramos en las matrices que representan las ecuaciones en el espacio de estado, lo que se busca es que el vector de entrada o control u(k) "afecte" directamente o indirectamente a los valores de la matriz A (matriz de estado) por medio de la matriz B (matriz de entrada) en cada una de las muestras durante el periodo de simulación – para este caso-.

Para verificar esa transmisión de control se arma la matriz de controlabilidad.

Suponiendo que el estado final es el origen en el espacio de estados y que el tiempo inicial es cero, se puede considerar que es condición necesaria y suficiente para la controlabilidad completa de estado que el rango de la matriz Co sea igual al número de variables de estado que posee el sistema (o su grado).

Al analizar el rango se busca que no existan vectores columnas (H,GH,..etc) nulos y que sean linealmente independientes unos de otros, de esta manera nos aseguramos que existirá la influencia del control en todo el espacio de estados de dicho sistema.

Utilizando el comando ctrb de Matlab (no es necesario que el sistema esté en la forma canónica controlable):

» Co = ctrb(HD.A,HD.B);% armo la matriz de controlabilidad

» nocont=length(A)-rank(Co) %nocont es el numero de estados no controlables

nocont = 0

El rango de la matriz Co es 2 por lo que el sistema es controlable.

También puede verificarse que con los datos obtenidos del rango, podemos decir que el determinante de Co es distinto de cero, por lo que estamos frente a una matriz no singular (tiene inversa); entonces, si tenemos en cuenta que el control u(k) es un escalar no acotado, podemos encontrar la secuencia única de señales de control u(k) siendo k=0,1,2,3 ..n-1; tal que cualquier estado inicial x(0) se transfiera al estado deseado final en n períodos de muestreo. El vector que contiene la secuencia de control puede obtenerse como se muestra a continuación:

8- Análisis de Observabilidad del Sistema Muestreado

Un sistema se dice que es completamente observable, si para cualquier estado inicial X(0) puede determinarse a partir de la observación de Y(kT) sobre un numero finito de periodos de muestreo.

Análogamente, el sistema es observable si cualquier transición de los estados afecta a todos los elementos del vector de salida.

La observabilidad se utiliza para calcular las variables de estado observando la salida en caso de no tener acceso a ésta. En nuestro caso y como en la mayoría de los sistemas reales, es difícil medir las variables de estado del sistema durante el transcurso del tiempo, ya sea porque resulta inaccesible su medición (físicamente, tecnológicamente, etc.) o ésta no puede llevarse a cabo mientras actúa el sistema.

Para verificar lo anterior se construye la matriz de observabilidad:

Esta matriz nos sirve para determinar si un sistema de n-esimo orden es o no completamente observable. El sistema es completamente observable si y solo si la matriz Ob es de rango máximo, es decir, si Ob tiene n filas linealmente independientes.

Una vez que a partir de la observación de la salida se puede determinar x(0) , también podrá determinarse x(kT), ya que u(0), u(T), … , u((k-1)T) son conocidas.

Utilizando el comando obsv de Matlab

» Ob = obsv(HD.A,HD.C);% armo la matriz observabilidad

» noobs=length(A)-rank(Ob) % es el numero de estados no observables

noobs = 0

Por lo que podemos decir que el sistema es completamente observable.

9-Aplicando control sobre el sistema

9.1 Diseño vía ubicación de polos

Una vez hecho el análisis de controlabilidad y haber determinado que el sistema es de estado completamente controlable; además, si tenemos en cuenta que nuestro sistema ejemplo es estable, podemos hacerlo más estable aún(o desestabilizarlo según el caso); para ello podemos comenzar introduciendo un control simple a través de una ganancia de retroalimentación negativa.

Esta retroalimentación compensa la entropía del sistema, es decir la tendencia que tiende el sistema (a converger, a divergir u oscilar); introduciendo entropía negativa por medio de la retroalimentación se puede controlar al mismo siempre y cuando éste sea un sistema controlable.

Si u(k) es una señal de control no acotada y hacemos u(k) = -Kx(k) donde K es la matriz de ganancia; el sistema se convierte en un sistema de control de lazo cerrado.

La ecuación de estados puede escribirse como:

x(k+1)=(G-HK)x(k)

La selección de los polos en lazo cerrado o de la ecuación característica deseada, dependerá del efecto esperado; pudiendo éste ser basado en la velocidad de respuesta, amortiguación deseada, minimización de efectos adversos, ruidos etc.

La determinación de la matriz de ganancia buscada puede obtenerse por distintos métodos, de forma que los autovalores de |G-HK| sean los polos en lazo cerrado buscados (polos del sistema dentro del espacio de estados).

Utilizando Matlab- con comando pole- obtenemos los polos dentro de la representación del espacio de estados.

» pole(HD)

ans =

0.5281

0.9522

Gráficamente podemos ver la ubicación del los polos y ceros del sistema con el comando pzmap.

» pzmap(HD);zgrid

zgrid permite visualizar las líneas de los factores de amortiguación constantes y las líneas de frecuencia natural dentro del círculo unitario.
Los polos del sistema están graficados con cruces y los ceros con círculos.

El comando rlocus

El comando "rlocus" calcula el diagrama de Evans, que muestra las trayectorias de los polos del sistema S.I.S.O. con retroalimentación negativa en función de la ganancia introducida.

Este comando es útil para el estudio de los efectos de distintas ganancias en un sistema retroalimentado. Indirectamente esas trayectorias muestran información del tiempo de respuesta en frecuencia.

» rlocus(HD);zgrid

Podemos ver que esas trayectorias salen del circulo unitario, como veremos mas adelante, si ubicamos los polos en distintos puntos de la trayectorias se obtiene respuestas diferentes, si estos se ubican fuera del circulo, se logra un efecto indeseado, el de desestabilizar al sistema.

Pero antes veamos de qué manera afectan al sistema los valores de los polos reguladores:

*El grafico representa la salida de un sistema con amortiguación.

Mientras más grande es el valor real, más amortiguación (sobregiro) hay en la respuesta.
Mientras más grande es el valor imaginario, hay mayor tiempo de establecimiento.

A través de la prueba y error, podemos ir ubicando los polos teniendo en cuenta la sobreamortiguación o la subamortiguación, así como también, el tiempo de respuesta de establecimiento-rapidez- del sistema.

Además, se debe tener en cuenta que los polos y ceros mas cercanos al eje imaginario son los que mas influyen en la respuesta (polos y ceros dominantes).

– Matlab –

En principio, Matlab sugiere la utilización del comando place para la reubicación de polos, y así obtener la matriz de ganancia K; este comando engloba distintos métodos y comandos como el acker (Método Ackermann).

Pero podemos utilizar comandos más ventajosos para la ubicación de los polos, con el comando rlocusfind podemos ubicar interactivamente los polos y así obtener la ganancia asociada a esa ubicación de polos.

Una vez obtenida esa ganancia, podemos aplicársela al sistema con el comando feedback, y por ultimo visualizar como responde el sistema (ej.: con un impulso unitario); mas adelante veremos una herramienta grafica que permite hacer todo lo anterior de manera aún mas intuitiva-el rltool-.

Veamos:

-Paso 1

Combinamos los comandos rlocus, zgrid y rlocfind para poder visualizar y señalar en la ventana gráfica donde queremos ubicar los polos.

» rlocus(HD);zgrid;k=rlocfind(HD)

Select a point in the graphics window

-Paso 2

Matlab nos pide que indiquemos en la pantalla con el puntero del ratón donde queremos ubicar a los polos (dentro de las trayectorias marcadas).

Elegimos los polos dentro del círculo unitario, éstos quedan marcados con un par de cruces rojas como se ve en la siguiente figura:

El sistema devuelve las coordenadas y la ganancia asociada:

selected_point =

0.4338 + 0.8613i

k =23.1136

-Paso 3

Retroalimentamos al sistema con la ganancia obtenida utilizando el comando feedback (lo convertimos en un sistema de lazo cerrado).

» HDC= feedback(HD,k);

-Paso 4

Visualizamos la respuesta del sistema (renombrado HDC) frente al escalón unitario:

» step(HDC)% grafica la respuesta del sistema ante una entrada escalón unitario

Visualizamos que ha cambiado el comportamiento de la salida, ésta reacciona con más rapidez; con algo de amortiguación pero con demasiada retroalimentación negativa la cual ha hecho disminuir la amplitud de la salida (temperatura).

Repitiendo los pasos anteriores podemos ir reubicando los polos de forma que encontremos la respuesta esperada.

– Matlab –

Utilizando el rltool

Cuando invocamos este comando, se inicializa una interfaz gráfica interactiva para el diseño de sistemas S.I.S.O. de lazo cerrado utilizando técnicas de reubicación de polos. Esta herramienta también permite importar, modificar y guardar modelos previamente diseñados, así como también, trabajar con modelos de Simulink*(modelador y simulador de sistemas dinámicos).

Para trabajar sobre nuestro sistema podemos utilizar el comando como sigue:

» rltool(HD)

Podemos obtener distinta información del sistema:

Ventana principal asociada a nuestro sistema HD.

Información del sistema dentro del espacio de estados

Información sobre los ceros y los polos.

En los siguientes gráficos se muestra al sistema con polos reubicados en sus nuevas coordenadas y el grafico de (por ejemplo) la salida con una entrada escalón unitario al sistema.

Gráfico de la salida con amortiguación de acuerdo a los polos del sistema en las coordenadas anteriores:

Gráfico de la salida correspondiente a la ubicación de polos fuera del círculo unitario, como puede observarse el sistema se vuelve inestable.

9.2 Control Óptimo Cuadrático

Todos los objetivos de respuesta temporal (amortiguamiento, tiempo de asentamiento, etc.) se pierden para centrarse en un único objetivo: la minimización(o maximización, según el caso) de un índice de desempeño J.

El problema que viene a solucionar el control óptimo es el de saber donde han de reubicarse los autovalores del sistema para conseguir unos fines determinados, como puede ser la minimización de una función de costo

escalar J. Esta función de costo, o índice de desempeño, suele estar relacionada con la energía y por lo tanto con el costo económico del proceso.

El control óptimo consiste en obtener el vector de las señales de control u(k), que haga mínimo el índice de desempeño J. A este vector se le denomina trayectoria óptima de control.

Un ejemplo de índice de desempeño cuadrático para un proceso de tiempo finito es:

La razón de incluir las variables de control u(k) en el índice de comportamiento se debe a razones energéticas como se mencionó anteriormente, lo cual es muy importante puesto que minimiza la energía de las señales de entrada que son las que realizan los actuadotes, y por lo tanto, son las que cuestan capital (recursos).

Si el problema es de estado estacionario, es decir, si se desea hacer un control de proceso durante todo momento sin restricciones (N ≈ ∞), podemos omitir el primer término del lado derecho de la ecuación del índice ya que tiende a cero (corresponde al costo terminal del estado final); esto es así porque en un problema de minimización se desea que x(N) esté tan cerca del origen como sea posible. También se debe considerar si el estado fina x(N) no es fijo.

La ventaja de utilizar el esquema de control óptimo cuadrático es que el sistema diseñado será asintóticamente estable, excepto para algunos casos especiales (ver pag. 7).

Normalmente las matrices H, Q y R son matrices diagonales con ij = 0 si i j y además son matrices definidas positivas. Cuando hay alguna variable que no se desea minimizar entonces se hace el correspondiente término ij = 0.

En general, en todos los índices de comportamiento los coeficientes de las matrices Q y R inciden directamente sobre los valores de las señales del sistema x(k) y u(k); cuanto mayores sean los coeficientes, más se optimiza la variable que corresponda. Es decir, por ejemplo, cuanto más se desee limitar la excursión de una cierta variable de control uj, mayor se escogerá el elemento

rii de la matriz R.

Si la matriz R aumenta, las ganancias disminuyen y también disminuye el costo energético de control; pero con valores grandes de la matriz R puede empeorar el comportamiento del sistema.

Si R >> Q la contribución del control u(k) es muy grande comparándola con la contribución del estado x(k). Dicho de otra forma, disminuye el costo del control pero empeora el comportamiento del sistema.

Si Q >> R mejora la respuesta del sistema pero aumenta el costo energético de control. El vector de estado estará muy controlado, pero las señales de control u(k) pueden ser muy grandes.

Aspectos sobre la implementación de control óptimo

Controles Permitidos: Se debe tener en cuenta que el vector de control obtenido, sea consecuente con el sistema real estudiado; hay señales de control obtenidas analíticamente que son imposibles de implementar por restricciones físicas del sistema. Por ello no se debe dejar de tener en cuenta: la naturaleza del índice de desempeño, el estado inicial o la salida inicial, el estado deseado o la salida deseada y las restricciones.

Un índice de desempeño puede ser apropiado para un problema dado, pero muy inapropiado bajo las consideraciones de otro problema. Esto es, por ejemplo, si suponemos un problema donde se requiere que un índice preserve al máximo la energía del sistema (pero se requiere más tiempo de asentamiento); bajo otro problema, donde se busca una respuesta rápida (que seguramente requerirá utilizar la energía del sistema sin demasiadas restricciones) el primer índice resultaría inadecuado.

La implementación del control óptimo suele ser off-line, desarrollándose previamente de forma analítica por métodos numéricos. La solución más apropiada puede consistir en aplicar una ganancia de control constante, que corresponda a la solución en régimen permanente (estacionario), generalmente esta alternativa es la de implementación más simple pero requiere que el sistema cumpla con la condición necesaria de ser controlable (es suficiente que sea observable).

Para nuestro caso, debemos tener en cuenta que se trata de un sistema térmico y que por ejemplo, entre otras consideraciones, debemos analizar que la resistencia eléctrica puede calentar a los recintos con una velocidad de respuesta máxima, la cual no ha sido contemplada bajo el modelo matemático (con señal de control no restringida) además se observa que la entrada (en el sistema real) no es capaz de absorber energía (enfriar) a los recintos. Estas consideraciones deberían hacerse al aplicar todo tipo de control al sistema en un problema real.

– Matlab –

Comando dlqr.

El comando dlqr permite diseñar un regulador cuadrático lineal para sistemas discretos utilizando técnicas de "programación dinámica".

Este comando trabaja con un índice discretizado y utiliza la ecuación de Riccati para resolver problemas de horizonte finito o infinito de acuerdo al los parámetros que se le pasen al comando.

Sintaxis

[K,S,e] = dlqr(a,b,Q,R)

[K,S,e] = dlqr(a,b,Q,R,N)

Si se omite la matriz N(o H de acuerdo a otras nomenclaturas) el comando puede utilizarse para problemas de estado de etapas infinitas y puede obtenerse la matriz K asociada en estado estacionario (constante). Además, se obtiene la matriz S(o P) y los autovalores (e) para cada instante.

A continuación puede estudiarse como afectan las matrices Q y R a los polos del sistema y la ganancia asociada para un problema de estado estacionario. Al omitirse la matriz N, Matlab realiza las iteraciones correspondientes hasta encontrar la matriz constante K.

a) Q>>R

» Q = [10 0 ; 0 20];

R = [1];

[K,S,e] = dlqr(HD.A,HD.B,Q,R) % extraemos a A y B de HD

K =

0.7862 3.0155

S =

16.8376 7.6034

7.6034 43.8472

e =

0.5071

0.6066

b) R>>Q

» Q = [10 0 ; 0 20];

R = [80];

[K,S,e] = dlqr(HD.A,HD.B,Q,R)

K =

0.0596 0.2198

S =

25.3199 38.6332

38.6332 157.4081

e =

0.5279

0.9256

Puede analizarse, en el caso a, como cambia la respuesta del sistema (observando e) pero el costo de control es alto (R es bajo ó K es alto); al incrementar R en el punto b, se observa como disminuye la ganancia nuevamente y los polos se reubican cerca de su posición original – ver pág. 13.

En los siguientes casos se analizan otros tipos de respuesta.

c) Q>>R

» Q = [50 0 ; 0 20];

R = [1];

[K,S,e] = dlqr(HD.A,HD.B,Q,R)

K =

1.7592 3.8596

S =

77.4645 20.4881

20.4881 53.5935

e =

0.4961 + 0.1480i

0.4961 – 0.1480i

d) R>>Q

»Q = [50 0 ; 0 20];

R = [60];

[K,S,e] = dlqr(HD.A,HD.B,Q,R)

K =

0.1646 0.4885

S =

98.8871 84.9434

84.9434 254.9485

e =

0.5307

0.8892

e) R>>Q

» Q = [50 0 ; 0 20];

R = [150];

[K,S,e] = dlqr(HD.A,HD.B,Q,R)

K =

0.0820 0.2511

S =

104.4415 103.9542

103.9542 320.3764

e =

0.5291

0.9202

Puede examinarse que (en los casos d y e): a medida que aumentamos R, la señal de control está cada vez mas acotada por lo que la respuesta del sistema casi no cambia (demasiado control de energía a pesar de haber control sobre las variables de estado).

También pueden obtenerse resultados similares utilizando el comando dare (discrete-time algebraic Riccati equations) que resuelve ecuaciones de Riccati para tiempo discreto.

O puede resolverse iterativamente haciendo:

»Q = [10 0 ; 0 20];

R = [1];

G = HD.A; %extraemos la matriz A del SS

H = HD.B; %extraemos la matriz B del SS

S =[0 0 ; 0 0]; % suponiendo problema de estado estacionario

for i=1:18, % delimitado por cálculos anteriores

S = Q + G'*S*G -G'*S*H* inv(R+H'*S*H)*H'*S*G

K= inv(R+H'*S*H)*H'*S*G

end;

Se obtiene el la primera iteración:

S =

10 0

0 10

K =

0.2293 0.9406

En la iteración 16 S permanece constante, -ver caso a -:

S =

16.8376 7.6034

7.6034 43.8472

K =

0.7862 3.0155

Una vez que calculada la matriz S constante, obtenemos la matriz K asociada en estado estacionario.

Sabiendo que u(k) = -Kx(k) y teniendo x(k+1)=(G-HK)x(k),luego podemos calcular la trayectoria óptima y los valores del control óptimo.

Como se mencionó anteriormente, ya que se deben analizar cada uno de los valores de las variables de estado durante el proceso; es que este método, la mayoría de las veces, resulta difícil de implementar mientras actúa el sistema. De otro modo, como se realiza en algunos casos, se aplica la ganancia de control constante previamente calculada.

Existen otros métodos para el cálculo de control óptimo, como el de Euler-Lagrange, el método de valores propios y vectores propios y el enfoque de Lyapunov para resolver problemas en estado estacionario, entre otros.

10- Conclusiones Finales

No se debe dejar de considerar que en la implementación, la representación de los modelos matemáticos y las herramientas analíticas deben aproximarse lo más posible al ámbito real aunque no puedan ser aplicados verazmente.

Por otro lado, queda a la vista, la gran versatilidad de Matlab para desarrollar, entender, visualizar y ampliar conceptos representados analíticamente.

Referencias Bibliográficas

Ogata, K. (1998). Ingeniería de Control Moderno (3ª ed.). Prentice-Hall

Kuo, B. (1996). Sistemas de Control Automático (7ª ed.). Prentice-Hall

Sitios de Internet

-MATLAB Ver. 5.3 http://www.mathworks.com/

– TUTORIALES DE CONTROL PARA MATLAB – UNIVERSIDAD DE MICHIGAN E.E.U.U. http://www.engin.umich.edu/group/ctm/

– PAGINA WEB SOBRE SISTEMAS DE CONTROL UNIVERSIDAD EAFIT – COLOMBIA – http://www.control-systems.net/

-APUNTES -CONTROL MODERNO – FACULTAD DE INGENIERÍA – UBA http://www.fi.uba.ar/materias/61107/Apuntes/TZ00.pdf

-CONTROL MODERNO – UNIVERSIDAD DE GUADALAJARA – MEXICO – http://proton.ucting.udg.mx/materias/moderno/principal.htm

-TEORÍA DE CONTROL AUTOMÁTICO – UNIVERSIDAD DE SEVILLA – ESPAÑA –

http://www.esi2.us.es/~vargas/docencia/tca/tca.html

-SISTEMAS DE CONTROL AVANZADO – UNIVERSIDAD DE ALCALÁ – ESPAÑA http://www.depeca.uah.es/docencia/ING-ECA/ctr_avz/

– CONTROL ÓPTIMO LINEAL – UNIVERSIDAD AUTÓNOMA DE NUEVO LEÓN – MEXICO http://gama.fime.uanl.mx/~salinas/CO.html

– CONTROL MODERNO – FACULTAD DE INGENIERÍA UNIVERSIDAD NACIONAL DE LA PLATA – ARG. http://davinci.ing.unlp.edu.ar/controlm/electronica/archivos.htm

– AUTOMÁTICA – DEPARTAMENTO DE INGENIERÍA Y CIENCIA DE COMPUTADORES – UNIVERSIDAD JAUME I ESPAÑA http://www3.uji.es/~mortiz/automatica/asignatura%20automatica.htm

Apéndice

Resumen de los comandos más importantes del Control

System Toolbox de Matlab.

Autor:

Leonardo Villegas