dydtfraliu.m

用来进行分数阶Liu系统数值仿真的源代码——This is the program to simulate fractional-order Liu system

%function dy=DYDtFraLiu
q1=1;q2=1;q3=0.95;
h=0.01;N=3500;
a=10;b=4;c=2.5;d=1;
x0=0.1;y0=0;z0=0.2;
M1=0;M2=0;M3=0;
x(N+1)=[0];y(N+1)=[0];z(N+1)=[0];
x1(N+1)=[0];y1(N+1)=[0];z1(N+1)=[0];
x1(1)=x0+a*h^q1*(y0-x0)/(gamma(q1)*q1);
y1(1)=y0+h^q2*a*x0*(b-z0)/(gamma(q2)*q2);
z1(1)=z0+h^q3*d*(b*x0*x0-c*z0)/(gamma(q3)*q3);

x(1)=x0+a*h^q1*(y1(1)-x1(1)+q1*(y0-x0))/gamma(q1+2);
y(1)=y0+a*h^q2*(x1(1)*(b-z1(1))+q2*x0*(b-z0))/gamma(q2+2);
z(1)=z0+d*h^q3*(b*x1(1)*x1(1)-c*z1(1)+q3*(b*x0*x0-c*z0))/gamma(q3+2);
for n=1:N
    M1=(n^(q1+1)-(n-q1)*(n+1)^q1)*(y0-x0);
    M2=(n^(q2+1)-(n-q2)*(n+1)^q2)*x0*(b-z0);
    M3=(n^(q3+1)-(n-q3)*(n+1)^q3)*(b*x0*x0-c*z0);
    N1=((n+1)^q1-n^q1)*(y0-x0);
    N2=((n+1)^q2-n^q2)*x0*(b-z0);
    N3=((n+1)^q3-n^q3)*(b*x0*x0-c*z0);
for j=1:n
  M1=M1+((n-j+2)^(q1+1)+(n-j)^(q1+1)-2*(n-j+1)^(q1+1))*(y(j)-x(j));
  M2=M2+((n-j+2)^(q2+1)+(n-j)^(q2+1)-2*(n-j+1)^(q2+1))*x(j)*(b-z(j));
  M3=M3+((n-j+2)^(q3+1)+(n-j)^(q3+1)-2*(n-j+1)^(q3+1))*(b*x(j)*x(j)-c*z(j));
  N1=N1+((n-j+1)^q1-(n-j)^q1)*(y(j)-x(j));
  N2=N2+((n-j+1)^q2-(n-j)^q2)*x(j)*(b-z(j));
  N3=N3+((n-j+1)^q3-(n-j)^q3)*(b*x(j)*x(j)-c*z(j));
end
x1(n+1)=x0+a*h^q1*N1/(gamma(q1)*q1);
y1(n+1)=y0+a*h^q2*N2/(gamma(q2)*q2);
z1(n+1)=z0+d*h^q3*N3/(gamma(q3)*q3);
x(n+1)=x0+a*h^q1*(y1(n+1)-x1(n+1)+M1)/gamma(q1+2);
y(n+1)=y0+a*h^q2*(x1(n+1)*(b-z1(n+1))+M2)/gamma(q2+2);
z(n+1)=z0+d*h^q3*(b*x1(n+1)*x1(n+1)-c*z1(n+1)+M3)/gamma(q3+2);
end
%m5=xcorr(x,y,'coeff');
%size(m5);
%subplot(1,1,1),plot(x,y);Xlabel('X'),ylabel('Y');%,Zlabel('Z');
subplot(1,1,1),plot3(x,y,z);xlabel('X'),ylabel('Y'),zlabel('Z');
%subplot(2,2,1),plot(x,y);Xlabel('X'),ylabel('Y');
%subplot(2,2,2),plot(x,z);Xlabel('X'),Zlabel('Z');
%subplot(2,2,3),plot(y,z);ylabel('Y'),Zlabel('Z');
%subplot(2,2,4),plot3(x,y,z);Xlabel('X'),ylabel('Y'),Zlabel('Z');
%A=[max(x) max(y) max(z)];
%B=[min(x) min(y) min(z)];
%axis([0 4000 -1 1]);
%fn=fopen('liu.txt','wt');
%fwrite(fn,m,'integer*1');
%fclose(fn)
grid on