|
楼主 |
发表于 2007-8-30 21:16
|
显示全部楼层
回复 #21 咕噜噜 的帖子
是贴错了.下面是对的.请问咕噜噜对我的模型有没有其他较合适的方法得到分支图?
function dy =epidemic(t,y)
global r a2 beta c c1 c2 k k1 k2
dy=zeros(3,1); % a column vector
dy(1)=r.*y(1).*(1-y(1)/k-y(2)/k)-beta.*y(1).*y(2);
dy(2)=beta.*y(1).*y(2)-c.*y(2)-c1.*y(3).*y(2)/(y(2)+k1);
dy(3)=(a2-c2.*y(3)/(y(2)+k2)).*y(3);
clc;clear all
global r a2 beta c c1 c2 k k1 k2
r=5;a2=1.2;beta=7;c=1.5;c1=2;c2=1;k=3;k1=0.6;k2=0.5;
options=odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-4 1e-4]);
[t,y]=ode45('epidemic',[0 100],[0.5 1 0.8],options);
figure(1);
plot(t,y(:,1),'b',t,y(:,2),'r');
figure(2);
plot(t,y(:,3),'g');
figure(3);
plot3(y(:,1),y(:,2),y(:,3),'r'); |
|