| 网站首页 | 业界新闻 | 小组 | 威客 | 人才 | 下载频道 | 博客 | 代码贴 | 在线编程 | 编程论坛
欢迎加入我们,一同切磋技术
用户名:   
 
密 码:  
共有 1937 人关注过本帖
标题:这个脉冲微分方程的数值模拟程序有什么问题.
只看楼主 加入收藏
bertqq
Rank: 1
等 级:新手上路
帖 子:2
专家分:0
注 册:2008-10-25
收藏
 问题点数:0 回复次数:2 
这个脉冲微分方程的数值模拟程序有什么问题.
function DrawTraj6()
PeriodT=1;  p=0.55;
NumOfStep=50;NumOfPer=5;
Init_x1=0.3; Init_x2=0.5;  Init_t=0;
NStep1=NumOfStep;
x1=zeros((NumOfStep+2)*NumOfPer,1);
x2=zeros((NumOfStep+2)*NumOfPer,1);
t=zeros((NumOfStep+2)*NumOfPer,1);
ppp=odeset('RelTol',1e-8,'AbsTol',1e-8);
for i=1:NumOfPer
    if i==1
        x1(1)=Init_x1;x2(1)=Init_x2;
    else
        Init_x1=x1((i-1)*(NumOfStep+2));
        Init_x2=(1-p)*x2((i-1)*(NumOfStep+2));
        x1((i-1)*(NumOfStep+2)+1)=Init_x1;
        x2((i-1)*(NumOfStep+2)+1)=Init_x2;
      
    end
   
    t((i-1)*(NumOfStep+2)+1)=(i-1)*PeriodT+Init_t;
    sol=ode45(@OdePP,[(i-1)*PeriodT+Init_t,i*PeriodT+Init_t],[Init_x1,Init_x2],ppp);
    j=linspace(PeriodT*(i-1)+Init_t,i*PeriodT+Init_t,NStep1);
    x1(((i-1)*(NumOfStep+2)+2):((i-1)*(NumOfStep+2)+1+NStep1))=deval(sol,j,1);
     x2(((i-1)*(NumOfStep+2)+2):((i-1)*(NumOfStep+2)+1+NStep1))=deval(sol,j,2);
      t(((i-1)*(NumOfStep+2)+2):((i-1)*(NumOfStep+2)+1+NStep1))=j;
     
 end
 
 save TrajData
 figure(1);
 plot(x1,x2,'k');
 box on
xlabel('x_1');
 ylabel('x_2');
   hgsave('Traj');
  figure(2)
   plot(t,x2,'k');
   xlabel('t');
   ylabel('x_2');
   hgsave('tseries_x1');
 function dxdt=OdePP(t,x)
   r1=8; b1=0.1; a1=2.5; r2=2; a2=0.58; k1=5; k2=5;
   dxdt=[x(1)*(r1-b1*x(1)-a1*x(2)/(x(1)+k1))
    x(2)*(r2-a2*x(2)/(x(1)+k2))
     ];
搜索更多相关主题的帖子: 微分方程 脉冲 数值模拟 
2008-10-25 17:06
pengchengxin
Rank: 2
等 级:论坛游民
帖 子:3
专家分:40
注 册:2012-2-24
收藏
得分:0 
.....
2012-04-07 18:50
快速回复:这个脉冲微分方程的数值模拟程序有什么问题.
数据加载中...
 
   



关于我们 | 广告合作 | 编程中国 | 清除Cookies | TOP | 手机版

编程中国 版权所有,并保留所有权利。
Powered by Discuz, Processed in 0.025388 second(s), 7 queries.
Copyright©2004-2024, BCCN.NET, All Rights Reserved