ode45求解倒立摆模型一致控制问题

小小鲁班 · 发表于 2020-10-23 13:50

EDA365欢迎您登录！

您需要登录才可以下载或查看，没有帐号？注册

x

%%下面是画图的函数
x0=[pi/3;pi/4;-pi/3;-pi/5;
0;0;0;0;
0;0;0;0;
0;0;0;0;
0;0;0;0;
0;0]; %%初始条件
[t,x]=ode45(@example180928,[0:30],x0);%求解微分函数值
figure(1)
plot(t,x(:,1),'-.b' ,t,x(:,2),'-m' ,t,x(:,3),'--g' ,t,x(:,4),':k' ,t,x(:,21),'-.r','LineWidth',2.0);
%axis([0 20 -2 5]);
set(gca,'FontSize',10,'LineWidth',1.2); set(gcf,'Color','White');
xlabel('Time (seconds)','FontSize',16.0);
legend('x_{1}','x_{2}','x_{3}','x_{4}','x_{0}');
title('The trajectories of five agents','FontSize',16.0 );
%%%%下面是函数定义的程序
function [ dx] = example180928( t,x )
A=[0 1 0 1;%拓扑
   1 0 1 0;
   0 1 0 1;
   1 0 1 0];
H=diag([sum(A(1,

) sum(A(2,

) sum(A(3,

) sum(A(4,:))]);
L=H-A;
B=diag([1 0 0 0]);
%%%实际系统参数
f=zeros(4,1);
d=zeros(4,1);
g=zeros(4,1);
for i=1:4
li0=1+0.04*i;li1=0.5*li0;mi=0.4*li0;
XIi=(cos(x(i))+li0/li1)^2*li1^2/li0;
f(i)=li1*sin(x(i))*(1+li1/li0*cos(x(i)))*x(i+4)^2/XIi-9.8*sin(x(i))*(1+cos(x(i)))/XIi;
g(i)=0.4/(mi*li0*XIi);
d(i)=0.2*i*cos(t);
end
d=[d(1); d(2); d(3); d(4)];
f=[f(1); f(2); f(3); f(4)];
gg=diag([g(1) g(2) g(3) g(4)]);
%%%%%%%5%x1变量x,v,参考轨迹
q=[x(1) x(2) x(3) x(4)]';%所有智能体的位置
v=[x(5) x(6) x(7) x(8)]';%所有智能体的速度
Wf=[x(9) x(10) x(11) x(12)]';%未知动态
k=[x(13) x(14) x(15) x(16)]';%nussbaum
alpha=[x(17) x(18) x(19) x(20)]';%未知外部扰动误差
vd=x(22);
%dvd=x(22);
%%%%%%%%%%%%%%%%跟踪误差z
xi=zeros(4,1);
for i=1:4
S=0;
for j=1:4
      S=A(i,j)*(x(i)-x(j))+S;
end
xi(i)=S+B(i,i)*(x(i)-x(21));
end
xi=[xi(1); xi(2); xi(3); xi(4)];%\xi

dxi=zeros(4,1);
for i=1:4
T=0;
for j=1:4
      T=A(i,j)*(x(i+4)-x(j+4))+T;
end
dxi(i)=S+B(i,i)*(x(i+4)-x(22));
end
dxi=[dxi(1); dxi(2); dxi(3); dxi(4)];%\dxi
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%模糊逻辑系统%%%%%%%%%
s1=zeros(16,1);
for i=1:4
phi1=[exp(-(x(i)-3+1)^2/16);
      exp(-(x(i)-3+2)^2/16);
      exp(-(x(i)-3+3)^2/16);
      exp(-(x(i)-3+4)^2/16)];
ph1i=phi1/sum(phi1);
s1(4*(i-1)+1:4*i)=ph1i;
end
Sf1=s1(1:4);Sf2=s1(5:8);Sf3=s1(9:12);Sf4=s1(13:16); %%\varphi_{i}
Sf=[Sf1 Sf2 Sf3 Sf4];
%%%%%%参数,
Beta=[15 0 0 0 ;
      0 15 0 0 ;
      0 0 15 0 ;
      0 0 0 15];
Mu=[1 0 0 0 ;
   0 1 0 0 ;
   0 0 1 0 ;
   0 0 0 1];
Rho=[1 0 0 0 ;
   0 1 0 0 ;
   0 0 1 0 ;
   0 0 0 1];
  z=xi+dxi;
  z=[z(1);z(2);z(3);z(4)];
  dWf=Mu*Sf*z;
  bf=Sf*Wf;
  theta=exp(-0.05*t);
  pi=[z(1)*x(17)/(z(1)^2+theta^2)^0.5;
z(2)*x(18)/(z(2)^2+theta^2)^0.5;
z(3)*x(19)/(z(3)^2+theta^2)^0.5;
z(4)*x(20)/(z(4)^2+theta^2)^0.5];
  dalpha=Rho*[z(1)^2/((z(1)^2+theta^2)^0.5);
         z(2)^2/((z(2)^2+theta^2)^0.5);
         z(3)^2/((z(3)^2+theta^2)^0.5);
         z(4)^2/((z(4)^2+theta^2)^0.5)];
dk=diag([z(1) z(2) z(3) z(4)])*(Beta*z+bf+pi);
Nk1=-exp(0.5*x(13)^2)*(x(13)^2+2)*sin(x(13));Nk2=-exp(0.5*x(14)^2)*(x(14)^2+2)*sin(x(14));
Nk3=-exp(0.5*x(15)^2)*(x(15)^2+2)*sin(x(15));Nk4=-exp(0.5*x(16)^2)*(x(16)^2+2)*sin(x(16));
Nk=diag([Nk1 Nk2 Nk3 Nk4]);
u=Nk*(Beta*z+bf+pi);
%%%%%%%%%%%%%%%%%%动态方程的表达式
dq=v;
dv=f+gg*u+d;
dxd=vd;%%%leader
dvd=cos(t);
%%%%%%%%%%%%%%%%%%解微分方程
dx=[dq;dv;dWf;dk;dalpha;dxd;dvd];
x
end

图像只出坐标轴不显示图像，x的结果到一定时间之后是NAN，不知道要怎么调试参数，求助

shuddkk · 发表于 2020-10-23 14:37

帮你顶一下

大小的小 · 发表于 2020-10-23 16:23

来学习，蹲一个大神

zzz.dan · 发表于 2020-10-23 16:32

帐号		自动登录	找回密码
密码			注册

ode45求解倒立摆模型一致控制问题

EDA365欢迎您登录！

推荐内容 /1