脈衝時滯微分方程matlab方程

2021-07-13 22:01:46 字數 3617 閱讀 1319

function maichon_shizhi1()

clear;clc

periodt=10;q=0.5;

%numofstep 是每半個週期取多少點畫圖,numofper 是畫多少個週期的相圖

numofstep=50;numofper=10;

%以下是初值和初始時刻

inx1=1;inx2=1;inx3=1;inxt=0;

for j=1:20

x1=zeros((numofstep+1)*numofper,1);

x2=zeros((numofstep+1)*numofper,1);

x3=zeros((numofstep+1)*numofper,1);

t=zeros((numofstep+1)*numofper,1);

temp1=0;temp2=0;te***=0;

%延遲時間值

lags=[0.5];

%下面是控制誤差的引數,越小則執行速度越慢

ppp=odeset('reltol',1e-8,'abstol',1e-8);

for i=1:numofper  %迴圈一次畫乙個週期的軌線

if i==1

x111= inx1;

x211= inx2;

x311= inx3;       

else

%以下三個語句確定脈衝後的值         

x111=temp1+2.1;

x211=temp2+0.5;

x311=te***;

inx1=x111;

inx2=x211;

inx3=x311;

endsol1 = dde23(@myddefun,lags,[inx1,inx2,inx3],[(i-1)*periodt+inxt,(i-0.5)*periodt+inxt]);

%     

%     sol1.x;

%     sol1.y(1,:);

%     sol1.y(2,:);

%     sol1.y(3,:);   

temp1=sol1.y(1,end);

temp2=sol1.y(2,end);

te***=sol1.y(3,end);

%以下三個語句確定脈衝後的值    

inx1 = 0.5 * temp1;

inx2 = temp2 + 5*temp1;

inx3 = te***;

x1((i-1)*(numofstep+1)+1)=inx1;

x2((i-1)*(numofstep+1)+1)=inx2;

x3((i-1)*(numofstep+1)+1)=inx3;

t((i-1)*(numofstep+1)+1)=(i-0.5)*periodt+inxt;

sol=ode45(@odepp,[(i-0.5)*periodt+inxt,i*periodt+inxt],[inx1,inx2,inx3],ppp);

j=linspace(periodt*(i-0.5+1/numofstep)+inxt,i*periodt+inxt,numofstep);

x1(((i-1)*(numofstep+1)+2):(i*(numofstep+1)))=deval(sol,j,1);

x2(((i-1)*(numofstep+1)+2):(i*(numofstep+1)))=deval(sol,j,2);

x3(((i-1)*(numofstep+1)+2):(i*(numofstep+1)))=deval(sol,j,3);

t(((i-1)*(numofstep+1)+2):(i*(numofstep+1)))=j;

t1=t(((i-1)*(numofstep+1)+1):(i*(numofstep+1)));

x11=x1(((i-1)*(numofstep+1)+1):(i*(numofstep+1)));

x21=x2(((i-1)*(numofstep+1)+1):(i*(numofstep+1)));

x31=x3(((i-1)*(numofstep+1)+1):(i*(numofstep+1)));

temp1=x11(end,1);

temp2=x21(end,1);

te***=x31(end,1);

t1_zhuanzhi= (t1).';

x11_zhuanzhi=(x11).';

x21_zhuanzhi=(x21).';

x31_zhuanzhi=(x31).';

%     figure(1)      

% %     plot(sol1.x,sol1.y(1,:),(t1).',x11_zhuanzhi);     

%     ylabel('s1(t)');

%     hold on

%     axis([(i-1)*periodt,i*periodt,0,2]);

% %     figure(2)    

%     plot(sol1.x,sol1.y(2,:),(t1).',x21_zhuanzhi);

%     ylabel('s2(t)');

%     hold on

%     figure(3)    

%     plot(sol1.x,sol1.y(3,:),(t1).',x31_zhuanzhi);

%     ylabel('x(t)');

%     hold on

figure(4)

plot3(horzcat(sol1.x,t1_zhuanzhi),horzcat(sol1.y(1,:),x11_zhuanzhi),horzcat(sol1.x,t1_zhuanzhi),horzcat(sol1.y(3,:),x31_zhuanzhi),horzcat(sol1.x,t1_zhuanzhi),horzcat(sol1.y(2,:),x21_zhuanzhi));

xlabel('s1(t)');

ylabel('x(t)');

zlabel('s2(t)');

hold on

%     figure(5)

%     plot3(horzcat(sol1.y(1,:),x11_zhuanzhi),horzcat(sol1.y(3,:),x31_zhuanzhi),horzcat(sol1.y(2,:),x21_zhuanzhi));

% %     xlabel('s1(t)');

% %     ylabel('x(t)');

% %     zlabel('s2(t)');

%     hold on

endinx1=inx1+0.5;

inx2=inx2+0.3;

inx3=inx3+0.2;

endfunction dx= myddefun(t,x,z)

dx=[

2.5-x(1);

-x(2)-(x(2)*x(3))/(0.5*(1+x(2)));

-x(3)+exp(-0.5)*z(2,1)*z(3,1)/(1+z(2,1))];

function dxdt=odepp(t,x)

%以下引數是系統引數

dxdt=[ 2.5-x(1)

-x(2)

-0.5*x(3) ];

matlab求解時滯微分方程

matlab求解時滯微分方程,dde23呼叫格式 sol dde23 ddefun,lags,history,tspan ddefun函式控制代碼,求解微分方程y f t,y t y t 1 y t k 必須寫成下面形式 dydt ddefun t,y,z 其中t對應當前時間t,y為列向量,近似於y...

matlab解微分方程

1.dsolve函式 這是最簡單的一種求解微分方程的一種方法 符號解法。一般來說,在matlab中解常微分方程有兩種方法,一種是符號解法,另一種是數值解法。在本科階段的微分數學題,基本上可以通過符號解法解決。用matlab解決常微分問題的符號解法的關鍵命令是dslove命令。該命令中可以用d表示微分...

解方程及微分方程 MATLAB

x1,x2,x3,solve eq1 eq2 eq3 x1 x2 x3 s solve eq1 eq2 eq3 x1 x2 x3 第一種方式,對solve的括號中的x1,x2,x3的次序會有要求 與前面的一致 第二種無要求。eq都是等式,不是表示式。最好都打引號,不打引號也可以,但是一定要預先用sy...