Considere el siguiente sistema:
10
G ( s )= 2
s +0.1 s+ 3
Donde se desea hacer el seguimiento de un objeto como el que se muestra a
continuación.
1.- Bajo el diseño de una ley de control tipo LQR obtener el mejor desempeño
2.- Bajo el diseño de una ley de control tipo MPC obtener el mejor desempeño
Nota 1: El examen debe entregarse el miércoles 15 de noviembre en la hora de clase
Nota 2: El examen debe contener todos los desarrollos correspondientes, incluyendo
gráficas y código de programas.
Nota 3: LA ley de control MPC debe ser en su forma completa.
Nota 4: En Teams se debe compartir los programas correspondientes.
Desarrollo:
Para empezar, se debe pasar la función de transferencia a modelo de estados y bueno
este depende como se realice, yo lo realice de esta manera.
10
G ( s )= 2
s +0.1 s+ 3
y (s ) 10
=
u(s) s2 +0.1 s +3
y (s)(s¿¿ 2+0.1 s+3)=10 u(s)¿
2
Ys +0.1 Ys+3 Y =10 U
1 1
Y = [ 10 U −3 Y ] −[ 0.1Y ]
s
2
s
Y ≔ x1 → y ≔ x1
1 1
S x 1 → ẋ1= [ 10 u−3 y ] − [ 0.1 y ]
s
2
s
∴ ẋ 1=x 2−0.1 x 1
S x 2 → ẋ2 =[ 10u−3 x 1 ]
ẋ (t )=
[−0.1
−3 0 ]
1
x ( t )+
[ 100 ] u(t )
y ( t ) =[ 1 0 ] x (t)
Sin embargo, al comprobar en Matlab si el resultado era correcto me arrogo lo
siguiente:
Por lo que decidí tomar los valores que arrogo Matlab
Por otra parte, de tenia que modificar la señal de referencia de nuestro bezier para que
se viera como la que tenia que seguir el sistema para lo que se hizo lo siguiente.
El Bezier que se compartió a lo largo de las clases es el que se ve en la imagen de
arriba donde sus parámetros son los que se observan en la imagen de abajo, sin
embargo, se pueden modificar para llevar la curva de respuesta a lo que se requiere.
La respuesta del Bezier sisn modificar sus parámetros nos da la curva que se observa
abajo por lo que se debe modificar los parámetros.
De igual forma se quitaron las señales de velocidad y aceleración del Beizer ya que no
interesaban para este caso, quedando como se muestra en la figura de arriba.
Propuse los parámetros que se ven en la imagen superior porque son los que más
hacían que se pareciera su curva a la que se tenia que llegar del sistema, así que deje
esos valores.
En la imagen superior se puede apreciar la respuesta de Bezier con los parámetros que
propuse y se aproximaban al seguimiento que debe seguir el sistema, por lo que esta
es la señal de referencia.
Desarrollo de la ley de control tipo LQR
Como ya se tenia el espacio de estados solo era de meter valores en el programa que
se vio durante clases y modificar Q y R para que la respuesta del control LQR se
aproximara a la señal de referencia.
Del mismo modo el diagrama a bloques que se utilizo en simulink se tuvo que modificar
quitando las señales de seno y escalón, como se ve en la imagen inferior así quedo el
diagrama a bloques.
Del mismo modo se modificó los parámetros del espacio de estados como se obtuvo
anteriormente.
Para obtener la salida deseada se toma el valor que se obtiene del LQR en este caso
como es un valor de una matriz y un vector se selecciona esa función y se escribe en
ganancia k(2) para el valor que necesitamos.
Primer propuesta Q=1 y R=1
Segunda propuesta Q=2 y R=1
Tercera propuesta Q=2 y R=0.5
Cuarta propuesta Q=2 y R=0.435
Quinta propuesta Q=4 y R=1
Sexta propuesta Q=4.5 y R=1
Séptima propuesta Q=4.69 y R=1
Puedo concluir después de las propuestas de Q y R que hice que fueron más no solo
las que se registraron en el documento, sin embargo, a lo largo de las simulaciones se
pudo observar que tener valores de Q=1 y R=1 la respuesta tenia un valor de
estabilidad dejos del de referencia, por lo que decidí a variar valores de Q sin embargo
se puedo ver que se aproximaba el valor pero llego un momento en que se comenzaba
a desestabilizar y se noto que la respuesta del error se volvía negativa lo que indicaba
que la señal de respuesta ya había sobrepasado la señal de respuesta, también se
notó que algunas propuestas de valore parecía que eran las definitivas para la
estabilidad del sistema o que llegaran a la señal de respuesta sin embargo al aumentar
el tiempo de simulación se podía apreciar que una señal estaba por arriba de la otra y
en la respuesta del error se podía apreciar esto también.
Del mismo modo se observó que para valores de Q iguales, pero cambiando R a
valores menor a uno se observo que el error tendía ser negativo y que las señales de
respuesta no lograban llegar a la señal de referencia.
Hablando de la acción de control se ve cuando se cambiaba el valor de R esta
disminuía sin embargo al cambiar Q tendía a consumir más energía de la acción de
control.
Al igual pude observar que fue un sistema sencillo de encontrar sus valores de Q y R a
diferencia del que se trabajo en clase que ningún valor se aproximaba a la señal de
referencia ni con valores pequeños ni grandes, si no hubiera sido sencillo se tenia que
ocupar un sistemas que aproximara los valores.
Desarrollo de la ley de control tipo MPC
Para el desarrollo de esta estrategia de control primeramente se establecieron algunos
parámetros como el horizonte de predicción, el horizonte de control, el tiempo de
muestreo y las iteraciones.
Había propuesto un tiempo de muestreo de 10 o 15 segundos, ya que se nota en la
grafica que en ese intervalo de tiempo se realiza un cambio significativo, para el
horizonte de predicción y de control, propuse valores de 12 y 10 respectivamente y
para las iteraciones establecí 10, sin embargo, al colocar esos valores en el programa
note que no hacían lo esperado y que al volver proponer valores algunas veces no
cambiaba el comportamiento o seguía guardando ese comportamiento.
El programa que utilice fue el siguiente:
clc
close all
A=[-0.1 -3; 1 0];
B=[1; 0];
C=[0 10];
D=[0];
Ts=15;
sys=ss(A,B,C,D);
sysd=c2d(sys,Ts);
Ad=sysd.A;
Bd=sysd.B;
Cd=sysd.C;
Nh=12;
Nc=6;
xk=[1.3 0]';
u=[1;2;0.5;-0.5;1;3;0];
ukmenos1=0;
uk=0;
%[s1,s2]=Solartrayectory ();
tray=Br;
%solucion analitica
iter=15000;
k=0.02;
[F,H]=PE(Ad,Bd,Cd,Nh,Nc);
tic
for i=1:iter
xk1=Ad*xk+Bd*uk;
yk=Cd*xk1;
fc=F*[xk;uk];
yref=tray(i:1:(Nh+i)-1)*pi/180;
U=inv(H'*H+k*eye(Nc,Nc))*H'*(yref-fc);
deltau=U(1);
xk=xk1;
ukmenos1=uk;
uk=ukmenos1+deltau;
sala(i)=yk;
refea(i)=yref(1);
controla(i)=uk;
end
toc
figure
set(gca,'fontsize',14)
plot(controla*100)
grid on
axis([0 60 0 20])
ax=gca;
figure
set(gca,'fontsize',14)
plot(sala*(180/pi)), hold, plot(refea*(180/pi),'r')
grid on
axis([0 60 0 20])
ax=gca;
figure
set(gca,'fontsize',14)
plot((refea-sala)*(180/pi))
grid on
axis([0 60 0 20])
ax=gca;
El programa d la función objetivo que utilice fue el siguiente:
function[F,H]=PE(Ad,Bd,Cd,Nh,Nc)
%Matrices originales
No=size(Cd);
Ni=size(Bd);
Ne=size(Ad);
Ne=Ne(1);
Ni=Ni(2);
No=No(1);
%Matrices aumentadas para modelo incremental
ceroM=zeros(Ni,Ne);
unoM=eye(Ni);
ceroQ=zeros(No,Ni);
M=[Ad Bd;ceroM unoM];
S=[Bd;unoM];
Q=[Cd ceroQ];
k=1;
for i=1:No:Nh*No
F(i:i+No-1,:)=Q*M^k;
k=k+1;
end
H=zeros(Nh*No,Nh*Ni);
k=0;
for i=1:No: Nh*No;
if i==1
H(i:i+No-1,i:i+Ni-1)=Q*M^k*S;
else
H(i:i+No-1,1:Ni*Nh)=[Q*M^k*S H(i-No:i-1,1:Nh*Ni-Ni)];
end
k=k+1;
end
H=H(:,1:Nc*Ni);
Al correr el programa se pueden observar tres figuras una es la acción de control, la
segunda es la respuesta y la señal de referencia y la tercera es el error.
Como mencione anteriormente no se podía cargar el Workspace del Bezier, sin limpiar
los datos de los workspace, por lo que los datos que se almacenan alteraban las
respuestas de salida.
Esta fue la gráfica del error de los valores que primero propuse sin embargo la
respuesta no se parecía en nada a la de referencia por lo que tuve que ir modificando
los valores para ver cómo se comportaba.
Para valores propuestos de: tiempo de muestreo 100, Nh=16, Nc=8, iter=1, k=500.
Acción de control:
Respuesta
Error
Valores propuestos tiempo de muestreo 150, Nh=12, Nc=6, iter=150, k=2000.
Acción de control
Respuesta del sistema
Error
Valores propuestos en otro programa, tiempo de muestreo 150, Nh=16, Nc=8, iter=1,
k=500, R=10
Respuesta del sistema.
Conclusión:
Durante todas las simulaciones que hice puedo concluir que es difícil identificar los
parámetros de un sistema y más debido a que estas van a interactuar con otros
parámetros y su respuesta no es lo que se esperaba, como en la primera simulación
que hice esperaba algo como la señal de referencia de respuesta sin embargo su
respuesta fue lo contrario. Del mismo modo puedo concluir que con el programa
completo añadiendo R, ayuda al sistema a mejorar su comportamiento se puede ver en
la última simulación que no fue la única que hice, pero solo puse esa que ya empezaba
a tomar mas la forma del sistema que se quería sin embargo se debía seguir ajustando
puesto que algunos valores no eran los esperados.
El implementar un MPC no es algo sencillo pues se requiere conocer muy bien los
parámetros de un sistema y algunas veces lo físico y lo simulado no será lo mismo, si
bien solo dará una aproximación del comportamiento.
Puedo decir que la diferencia del LQR y la implantación del MPC que es más fácil
implementar el LQR y que su respuesta se encuentra mas rápido a diferencia de MPC,
de igual forma puedo decir que la acción de control del MPC es menor diferencia del
LQR, dejando claro que la función de costos que se explico durante el parcial es cierta
y aquí se vio como es su comportamiento en un sistema.