倒立擺PID控制及MATLAB仿真

本文轉載自查看原文 2017-10-16 01:01 4082

說明：

全文共四個參考：

終極參考，程序+仿真圖+結果，但是缺少原理和文字建模過程。
參考0和1重在文字部分，模型建立怎么敘述，但無實際程序。
參考2重在程序部分，實際的matlab程序部分，直接上機調試。
每個大參考分別給了幾個小參考，以供選擇。

具體看個部分內容前先全部游覽下。

終極參考：

https://wenku.baidu.com/view/df64a8ece009581b6bd9eb2f.html

參考0 這是個題目，但是直接告訴你數學模型了，具體分析看參考1中的論文

參考一重在借鑒論文各部分敘述，但無實際程序，不用擔心，程序請看參考二。

1 倒立擺PID控制及MATLAB仿真

關鍵圖：

2 基於MATLAB的倒立擺控制及參數優化

關鍵圖截取：

參考二完整程序，自己修改數據和數學符號，配合參考一文字敘述模塊，整個問題就解決了。

下面給了兩個參考程序，第二個參考程序易懂，標注清楚。你都跑一下，看看效果。

程序1 http://ask.wenku1.com/question/%E5%80%92%E7%AB%8B%E6%91%86%E5%BB%BA%E6%A8%A1.htm

%主程序：

clear all
close all
global A B C D

g=9.8;%重力加速度
M=1.0;%小車質量
m=0.1;%桿的質量
L=0.5;%桿長度的一半
Fc=0.0005;%小車相對導軌的摩擦系數
Fp=0.000002;%桿相對於小車的摩擦系數

I=1/12*m*L^2;
l=1/2*L;
t1=m*(M+m)*g*l/[(M+m)*I+M*m*l^2];
t2=-m^2*g*l^2/[(m+M)*I+M*m*l^2];
t3=-m*l/[(M+m)*I+M*m*l^2];
t4=(I+m*l^2)/[(m+M)*I+M*m*l^2];

A=[0 1 0 0;
t1 0 0 0;
0 0 0 1;
t2 0 0 0];

B=[0;t3;0;t4];

C=[1 0 0 0;
0 0 1 0];

D=[0;0];

Q=[100 0 0 0;
0 10 0 0;
0 0 1 0;
0 0 0 1];

R=[0.1];
K=lqr(A,B,Q,R);

e1_1=0;
e2_1=0;
e3_1=0;
e3_1=0;

u_1=0;
xk=[-10/57.3,0,0.20,0];

ts=0.02;

for k=1:1:1000
time(k)=k*ts;
Tspan=[0 ts];

para=u_1;
[t,x]=ode45('daolibai_PID_Controller_SubFunction',Tspan,xk,[],para);
xk=x(length(x),:);

r1(k)=0.0;
r2(k)=0.0;
r3(k)=0.0;
r4(k)=0.0;

x1(k)=xk(1);
x2(k)=xk(2);
x3(k)=xk(3);
x4(k)=xk(4);

e1(k)=r1(k)-x1(k);
e2(k)=r2(k)-x2(k);
e3(k)=r3(k)-x3(k);
e4(k)=r4(k)-x4(k);

S=1;
if S==1
u(k)=K(1)*e1(k)+K(2)*e2(k)+K(3)*e3(k)+K(4)*e4(k);
elseif S==2
de1(k)=e1(k)-e1_1;
u1(k)=-50*e1(k)-10*de1(k);
de2(k)=e2(k)-e2_1;
u2(k)=-10*e2(k)-10*de2(k);
de3(k)=e3(k)-e3_1;
u3(k)=-10*e3(k)-10*de3(k);
de4(k)=e4(k)-e4_1;
u4(k)=-10*e4(k)-10*de4(k);
u(k)=u1(k)+u2(k)+u3(k)+u4(k);
end

if u(k)>=10
u(k)=10;
elseif u(k)<=-10
u(k)=-10;
end
e1_1=e1(k);
e2_1=e2(k);
e3_1=e3(k);
e4_1=e4(k);
u_1=u(k);
end

figure(1);
subplot(411);
plot(time,r1,'k',time,x1,'k');
xlabel('time(s)');ylabel('Angle');
subplot(412);
plot(time,r2,'k',time,x2,'k');
xlabel('time(s)');ylabel('Angle rate');
subplot(413);
plot(time,r3,'k',time,x3,'k');
xlabel('time(s)');ylabel('Car Position');
subplot(414);
plot(time,r4,'k',time,x4,'k');
xlabel('time(s)');ylabel('Car rate');

figure(5);
plot(time,u,'k');
xlabel('time(s)');ylabel('Force');

%子程序
function dx=dym(t,x,flag,para)
global A B C D
u=para;
dx=zeros(4,1);
dx=A*x+B*u;
%剛學，望多多交流，多多批評指正

0

程序2 http://www.verysource.com/code/6426077_1/untitled202.m.html

clear;
M = 1.096;
m = 0.109;
b = 0.1;
f=0;
I= 0.0034;
g = 9.8;
l = 0.25;
D = I*(M+m)+M*m*l^2; %denominator for the A and B matricies
A = [0      1              0           0;
    0  -(I+m*l^2)*b/D (m^2*g*l^2)/D    -m*l*f/D;
    0      0              0           1;
    0   -(m*l*b)/D    m*g*l*(M+m)/D  -(m+M)*f/D]
B = [     0; 
     -(I+m*l^2)/D;
          0;
        -m*l/D];
C = [1 0 0 0;
     0 0 1 0];
D = [0;
     0];
disp('The rank of controllability matrix')
rc=rank(ctrb(A,B))
% p=[-2-3*i,-2+3*i,-20,-21];
p=[-20,-21,-2-3*i,-2+3*i];

k=place(A,B,p)

Cn=[1 0 0 0 ];
Nbar=rscale(A,B,Cn,0,k)
t0=0;
tspan=0.006;
R=1;
y0=[0 0 0.1 0];
n=1000;
F=R*Nbar-k*y0';
for i=1:n
    tf=t0+tspan;
    [t,y]=ode45(@xprime10,[t0 tf],y0',[],F);
    t0=tf;
    y1=y(end,1);
    y2=y(end,2);
    y3=y(end,3);
    y4=y(end,4);
    y0=[y1 y2 y3 y4];
    F=R*Nbar-k*y0';
    h(i,1)=y1;
    h(i,2)=y2;
    h(i,3)=y3;
    h(i,4)=y4;
    h(i,5)=F;
end 
     T=0.006:0.006:6;
     figure(1)
     plot(T,h(:,1),':b')
     legend('小車位移')
     grid
     figure (2)
     plot(T,h(:,3),'-r')
     legend('擺桿擺角')
     grid
     figure(3)
     plot(T,h(:,5),'--g')
     legend('F')
     grid

　程序3　有網友說Matlab自帶一個例程，我電腦沒安裝matlab，你去找找。