%% Modelado de no-linealidad de sector como incertidumbre, diseño de controladores lineales con prestaciones garantizadas ante no-linealidad
% © 2019, Antonio Sala Piqueras, Universitat Politècnica de València.
% Presentación en vídeo en http://personales.upv.es/asala/YT/V/nlfml.html
% Tabla de Contenidos
% Planteamiento del problema
% 1. Solución: modelado
% 2. Solución: Diseño de control (tma. pequeña ganancia escalado)
%
%% Planteamiento del problema
% Un sistema no-lineal tiene un modelo dado por la ecuación de estado no lineal:
% \frac{dx_1}{dt}= \phi(x_1)+x_2+d \\
% \frac{dx_2}{dt}= 2\phi(x_1)-x_1-5x_2+3u
% La no-linealidad viene dada por \phi(x_1)=sin(x_1).
% La variable manipulada (acción de control) es u, siendo d una perturbación (entrada exógena) cuyo efecto debe ser cancelado. En concreto, ante una perturbación con un cierto contenido máximo en frecuencia, se desea minimizar la norma ponderada de (y,u) viniendo y determinada por la ecuación de salida y=x_1.
%
%% 1. Solución: modelado
% Una no-linealidad phi=f(z) de sector (esto es, que verifique mingan\leq \frac{f(z)}{z}\leq maxgan) puede expresarse como una "ganancia incierta", que tras escalar/normalizar la "incertidumbre", puede escribirse como este diagrama de bloques, donde \|\Delta\|_\infty\leq 1 (entendiendo la expresión como \|m\|_2\leq \|q\|_2):
%

%zona de operación
z=(-pi/2):(.001+1e-6):(pi/2);
phi=sin(z);
mingan=min(phi./z), maxgan=max(phi./z)
medgan=(mingan+maxgan)/2, incgan=maxgan-medgan
%matriz de transf. entre entradas (m,z) y salidas (q,phi).
cosaarriba=[0 1;incgan medgan];

% Integrémoslo todo en una única representación interna (excepto "Delta")
%

A=[0 1;-1 -5];
Cperf=[1 0]; %salidas con objetivos de control
C2=[1 0];%variables medidas

%Efecto de la no-linealidad
%phi=f(z); 
Cz=[1 0]; 

Bgrande=[1 1 0;...
   2 0 3]; %entradas phi,d,u
Cgrande=[Cz;Cperf;0 0;C2];%salidas z,y,u,y
sys_sinnolin=ss(A,Bgrande,Cgrande,[0 0 0;0 0 0;0 0 1;0 0 0]);

sys3x3=lft(cosaarriba,sys_sinnolin,1,1);
sys3x3.OutputName={'q','y','u','ycopia'};
sys3x3.InputName={'m','d','u'};
eig(sys3x3) %con la ganancia "media" es inestable

%% 2. Solución: Diseño de control (tma. pequeña ganancia escalado)

s=tf('s');
Wy=100/(s+1);Wu=1;
Wd=1/(s+2);

rangomultiplicador=logspace(-2,2.5,30);
resultado=zeros(size(rangomultiplicador));
for i=1:length(rangomultiplicador)
    multiplicador=rangomultiplicador(i);
    Wout=blkdiag(multiplicador,Wy,Wu,1);
    Win=blkdiag(1/multiplicador,Wd,1);
    GPW=minreal(Wout*sys3x3*Win);
    [K,CL,GAM,INFO]=hinfsyn(GPW,1,1);
    GAM
    resultado(i)=GAM;
end
%simulemos el proceso con la ganancia "media" (delta=0)
syssim=lft(tf(0),lft(sys3x3,K));
sigma(syssim*Wd), grid on
step(syssim*Wd)