T-S模糊模型与状态反馈控制及Matlab仿真（附代码）

出门两双鞋

已于 2022-03-16 18:55:02 修改

阅读量6.9k

点赞数 19

分类专栏：笔记文章标签：控制器

于 2021-12-09 16:11:40 首次发布

本文链接：https://blog.csdn.net/qq_29710939/article/details/121828964

版权

笔记专栏收录该内容

6 篇文章 4 订阅

订阅专栏

(一) 仿射非线性系统建模

以overhead crane system为例，

忽略摩擦，我们可以得到如下的状态空间方程，
$\dot x=\begin{pmatrix}x_2\\\frac{mlx_4^2\;\sin\;x_3+mg\;\sin\;x_3\cos\;x_3\;}{M+m\;\sin^2\;x_3}\\x_4\\-\frac{mlx_4^2\;\sin\;x_3\;\cos\;x_3+(M+m)g\;\sin\;x_3\;}{l(M+m\;\sin^2\;x_3)}\end{pmatrix}+\frac u{M+m\;\sin^2\;x_3}\begin{pmatrix}0\\1\\0\\-\frac{\cos\;x_3}l\end{pmatrix}$

(二) 计算T-S模糊模型子系统

取工作点。一般地，我们会选取平衡点。
$x_0=\begin{pmatrix}0\\0\\0\\0\end{pmatrix}$ , $x_1=\begin{pmatrix}5\\0\\\pi/6\\\pi/20\end{pmatrix}$ , $x_2=\begin{pmatrix}5\\-1\\-\pi/6\\0\end{pmatrix}$
$A_0=J(x_0)$ , $B_0=g(x_0)$
$A_i=J(x_i)+\frac {1}{{x_i}^Tx_i}[f(x_i)-J(x_i)x_i]{x_i}^T$ , $B_i=g(x_i)$ , $i = 1, 2$

(三) 建立推理，验证开环特性

if x is about $x_0$ , then $\dot x=A_0x+B_0u$
if x is about $x_1$ , then $\dot x=A_1x+B_1u$
if x is about $x_2$ , then $\dot x=A_2x+B_2u$
取 $h_0(x_3)=\left\{\begin{array}{l}\frac{\sin\;x_3}{x_3},x_3\neq0\\1,x_3=0\end{array}\right.$ ， $h_1(x_3)=h_2(x_3)=\frac{1-h_0(x_3)}{2}$

(四) 极点配置，验证闭环特性

选取极点 $p_0=\begin{pmatrix}-1.1 -1.5 -1.8 -2.1\end{pmatrix}$ ， $p_1=\begin{pmatrix}-1.3 -1.6 -1.9 -2.2\end{pmatrix}$ ， $p_2=\begin{pmatrix}-1.5 -1.6 -1.8 -1.9\end{pmatrix}$

(五) 使用LMI验证稳定性

$\exists P=\begin{pmatrix}0.4016&0.1427&0.5521&0.5098\\0.1427&0.1903&-0.0489&0.7069\\0.5521&-0.0489&5.7520&-0.1442\\0.5098&0.7069&-0.1442&2.8079\end{pmatrix}$
$\left\{\begin{array}{l}{(A_i-B_iK_i)}^TP+P(A_i-B_iK_i)<0,\;(i=1,2,3)\\{(A_i-B_iK_j+A_j-B_jK_i)}^TP+P(A_i-B_iK_j+A_j-B_jK_i)<0,(i<j\leq3)\end{array}\right.$

clear all;close all;
M=100;m=100;l=5;g=9.82;u=100;
syms x1 x2 x3 x4;
x00=[0 0 0 0]';x01=[5 0 pi/6 pi/20]';x02=[5 -1 -pi/6 0]';
x=[x1;x2;x3;x4];
f=[x2;(m*l*x4^2*sin(x3)+m*g*sin(x3)*cos(x3))/(M+m*sin(x3)*sin(x3));x4;-(m*l*x4^2*sin(x3)*cos(x3)+(M+m)*g*sin(x3))/l/(M+m*sin(x3)*sin(x3))];
gx=1/(M+m*sin(x3)*sin(x3))*[0;1;0;-cos(x3)/l];
J=jacobian(f,x);
Ai=J+(f-J*x)*x'/(x'*x);
A0=double(subs(J,x,x00));B0=double(subs(gx,x,x00));
A1=double(subs(Ai,x,x01));B1=double(subs(gx,x,x01));
A2=double(subs(Ai,x,x02));B2=double(subs(gx,x,x02));

a0=[-1.1 -1.5 -1.8 -2.1];
a1=[-1.3 -1.6 -1.9 -2.2];
a2=[-1.5 -1.6 -1.8 -1.9];

K0=place(A0,B0,a0);K1=place(A1,B1,a1);K2=place(A2,B2,a2);

setlmis([]);
P=lmivar(1,[4,1]);
lmiterm([1 1 1 P],1,A0-B0*K0,'s');
lmiterm([2 1 1 P],1,A1-B1*K1,'s');
lmiterm([3 1 1 P],1,A2-B2*K2,'s');
lmiterm([4 1 1 P],1,A0-B0*K1+A1-B1*K0,'s');
lmiterm([5 1 1 P],1,A1-B1*K2+A2-B2*K1,'s');
lmiterm([6 1 1 P],1,A2-B2*K0+A0-B0*K2,'s');
lmiterm([-7 1 1 P],1,1);
lmisys=getlmis;
[tmin,xfeas]=feasp(lmisys);
tmin
pmat=dec2mat(lmisys,xfeas,P)

$t_{min}=-0.0163$
$P=\begin{pmatrix}4.2687&3.1372&6.7828&14.3290\\3.1372&6.5107&-10.3561&29.7178\\6.7828&-10.3561&151.4136&-37.1458\\14.3290&29.7178&-37.1458&145.7880\end{pmatrix}$

下载链接
https://download.csdn.net/download/qq_29710939/84993530

trolly_sim.m

clear all;close all;
t0=0;tf=30;step=0.1;
x0(1:4)=[10 1 pi/10 0.1]';
global M m l u g;
global A0 A1 A2 B0 B1 B2 x3 K0 K1 K2;
M=100;m=50;l=4;g=9.81;u=100;
syms x1 x2 x3 x4;
x00=[0 0 0 0]';x01=[5 0 pi/6 pi/20]';x02=[5 -1 -pi/6 0]';
x=[x1;x2;x3;x4];
f=[x2;(m*l*x4^2*sin(x3)+m*g*sin(x3)*cos(x3))/(M+m*sin(x3)*sin(x3));x4;-(m*l*x4^2*sin(x3)*cos(x3)+(M+m)*g*sin(x3))/l/(M+m*sin(x3)*sin(x3))];
gx=1/(M+m*sin(x3)*sin(x3))*[0;1;0;-cos(x3)/l];
J=jacobian(f,x);
Ai=J+(f-J*x)*x'/(x'*x);
A0=double(subs(J,x,x00));B0=double(subs(gx,x,x00));
A1=double(subs(Ai,x,x01));B1=double(subs(gx,x,x01));
A2=double(subs(Ai,x,x02));B2=double(subs(gx,x,x02));

[t,x_fuzzy]=ode45('trolly_fuzzy',[t0:step:tf], x0);
[t,x_affine]=ode45('trolly_affine',[t0:step:tf], x0);
figure(1);
subplot(2,2,1);plot(t, x_fuzzy(:,1),'b-.',t,x_affine(:,1),'b:');xlabel('时间(秒)');ylabel('x1');
subplot(2,2,2);plot(t, x_fuzzy(:,2),'b-.',t,x_affine(:,2),'b:');xlabel('时间(秒)');ylabel('x2');
subplot(2,2,3);plot(t, x_fuzzy(:,3),'b-.',t,x_affine(:,3),'b:');xlabel('时间(秒)');ylabel('x3');
subplot(2,2,4);plot(t, x_fuzzy(:,4),'b-.',t,x_affine(:,4),'b:');xlabel('时间(秒)');ylabel('x4');

a0=[-1.1 -1.5 -1.8 -2.1];
a1=[-1.3 -1.6 -1.9 -2.2];
a2=[-1.5 -1.6 -1.8 -1.9];
K0=place(A0,B0,a0);K1=place(A1,B1,a1);K2=place(A2,B2,a2);

[t,x_fuzzy_ctrl]=ode45('trolly_fuzzy_ctrl',[t0:step:tf], x0);
[t,x_affine_ctrl]=ode45('trolly_affine_ctrl',[t0:step:tf], x0);
figure(2);
subplot(2,2,1);plot(t, x_fuzzy_ctrl(:,1),'b-.',t,x_affine_ctrl(:,1),'b:');xlabel('时间(秒)');ylabel('x1');
subplot(2,2,2);plot(t, x_fuzzy_ctrl(:,2),'b-.',t,x_affine_ctrl(:,2),'b:');xlabel('时间(秒)');ylabel('x2');
subplot(2,2,3);plot(t, x_fuzzy_ctrl(:,3),'b-.',t,x_affine_ctrl(:,3),'b:');xlabel('时间(秒)');ylabel('x3');
subplot(2,2,4);plot(t, x_fuzzy_ctrl(:,4),'b-.',t,x_affine_ctrl(:,4),'b:');xlabel('时间(秒)');ylabel('x4');

trolly_affine.m

function xdot = trolly_affine(t, x);
global M m l u g;
f=[x(2);(m*l*x(4)^2*sin(x(3))+m*g*sin(x(3))*cos(x(3)))/(M+m*sin(x(3))*sin(x(3)));x(4);-(m*l*x(4)^2*sin(x(3))*cos(x(3))+(M+m)*g*sin(x(3)))/l/(M+m*sin(x(3))*sin(x(3)))];
gx=1/(M+m*sin(x(3))*sin(x(3)))*[0;1;0;-cos(x(3))/l];
xdot=f+gx*u;

trolly_fuzzy.m

function xdot = trolly_fuzzy(t, x);
global A0 A1 A2 B0 B1 B2 u M m;
if(x(3)==0)
    h0=1;
else
    h0=sin(x(3))/x(3);
end
h1=(1-h0)/2;h2=h1;
xdot=h0*(A0*x+B0*u)+h1*(A1*x+B1*u)+h2*(A2*x+B2*u);

trolly_affine_ctrl.m

function xdot = trolly_affine_ctrl(t, x);
global M m l u g K0 K1 K2;
f=[x(2);(m*l*x(4)^2*sin(x(3))+m*g*sin(x(3))*cos(x(3)))/(M+m*sin(x(3))*sin(x(3)));x(4);-(m*l*x(4)^2*sin(x(3))*cos(x(3))+(M+m)*g*sin(x(3)))/l/(M+m*sin(x(3))*sin(x(3)))];
gx=1/(M+m*sin(x(3))*sin(x(3)))*[0;1;0;-cos(x(3))/l];
h0=M*(sin(x(3)))^2/(M+m*sin(x(3))*sin(x(3)));
h1=M*(cos(x(3)))^2/(M+m*sin(x(3))*sin(x(3)));
h2=m*(sin(x(3)))^2/(M+m*sin(x(3))*sin(x(3)));
xdot=f+gx*(u-h0*K0*x-h1*K1*x-h2*K2*x);

trolly_fuzzy_ctrl.m

function xdot = trolly_fuzzy_ctrl(t, x);
global A0 A1 A2 B0 B1 B2 u M m K0 K1 K2;
if(x(3)==0)
    h0=1;
else
    h0=sin(x(3))/x(3);
end
h1=(1-h0)/2;h2=h1;
u0=u-h0*K0*x-h1*K1*x-h2*K2*x;
xdot=h0*(A0*x+B0*u0)+h1*(A1*x+B1*u0)+h2*(A2*x+B2*u0);