24小时热门版块排行榜    

查看: 1295  |  回复: 7

xwndf250

银虫 (小有名气)

[求助] matlab处理常微分方程作图问题

流行病模型 sir
dI/dt=a*S*I-b*I
dS/dt=-a*S*I
dR/dt=b*I
S+I+R=1且0 想要做一个横坐标t=[0,50], 纵坐标dI/dt、dS/dt、dR/dt的3条曲线,求具体程序,本人写的程序
m文件:function y=SIR(t,x)
a=0.2;b=0.1;
y=[a*x(1)*x(2)-b*x(1);
-a*x(1)*x(2);
b*x(3)];
end
命令:应该怎样写,才是t与dI/dt、dS/dt、dR/dt的图像(就是S、I、R的变化率随时间的变化),注意,不是t与S、I、R的图像!最好用matlab语言。
回复此楼

» 猜你喜欢

» 本主题相关价值贴推荐,对您同样有帮助:

已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

csgt0

荣誉版主 (著名写手)

彩色挂图

【答案】应助回帖

★ ★
感谢参与,应助指数 +1
dbb627: 金币+2, 感谢应助 2013-02-26 14:07:58
CODE:
function xwn
options = odeset('RelTol',1e-4,'AbsTol',[1e-4 1e-4 1e-5]);
intvalue=[0.1 0.4 0.5];       %t=0时的初值
[T,Y] = ode45(@rigid,[0 50],intvalue,options);
plot(T,Y(:,1),'-',T,Y(:,2),'-.',T,Y(:,3),'.')
title('T-Y图')
for i=1:length(T)     
vdy(:,i)=rigid(T(i),Y(i,:));
end
figure
plot(T,vdy(1,:),'-',T,vdy(2,:),'-.',T,vdy(3,:),'.')
title('T-dY图')
end

function dy = rigid(t,y)
dy = zeros(3,1);    % a column vector
a=0.2;
b=0.1;
dy(1) = a*y(2)*y(1)-b * y(1);
dy(2) = -a*y(2) * y(1);
dy(3) = b * y(1);
end

showmethemoney
2楼2013-02-26 11:06:11
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

xwndf250

银虫 (小有名气)

引用回帖:
2楼: Originally posted by csgt0 at 2013-02-26 11:06:11
function xwn
options = odeset('RelTol',1e-4,'AbsTol',);
intvalue=;       %t=0时的初值
= ode45(@rigid,,intvalue,options);
plot(T,Y(:,1),'-',T,Y(:,2),'-.',T,Y(:,3),'.')
title('T-Y图')
for i=1: ...

function dy = rigid(t,y)后面的我看懂了,但是前面function xwn没有看懂,我记得求导貌似用diff函数的,能不能详细讲解下odeset()函数,intvalue作用,“for i=1:length(T) 和vdy(:,i)=rigid(T(i),Y(i,);”的意思。我是小白,求指教。好的话我多给金币。
3楼2013-02-26 21:06:43
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

csgt0

荣誉版主 (著名写手)

彩色挂图

【答案】应助回帖

★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★ ★
xwndf250: 金币+15, ★★★很有帮助, 非常感谢您的帮助,给您金币,顺便问下,t0初值怎么确定的?也就是0.1,0.4,0.5是怎么确定的? 2013-02-27 10:00:25
fegg7502: 金币+1, 应助指数+1, 鼓励交流 2013-04-02 09:30:14
odeset()函数设置数值解的计算精度,其实一般不用也可以
intvalue是初值啊,t=0时的初值
“for i=1:length(T) 和vdy(:,i)=rigid(T(i),Y(i,);”的为了循环计算每个t下的3个导数值
showmethemoney
4楼2013-02-27 09:54:46
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

csgt0

荣誉版主 (著名写手)

彩色挂图

【答案】应助回帖


fegg7502: 金币+1, 应助指数+1, 鼓励交流 2013-04-02 09:30:24
引用回帖:
3楼: Originally posted by xwndf250 at 2013-02-26 21:06:43
function dy = rigid(t,y)后面的我看懂了,但是前面function xwn没有看懂,我记得求导貌似用diff函数的,能不能详细讲解下odeset()函数,intvalue作用,“for i=1:length(T) 和vdy(:,i)=rigid(T(i),Y(i,);”的 ...

就是t=0是的ISR啊,具体多少得看你的实际情况,我只是随意写的个数。
showmethemoney
5楼2013-02-27 10:37:52
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

xwndf250

银虫 (小有名气)

引用回帖:
5楼: Originally posted by csgt0 at 2013-02-27 10:37:52
就是t=0是的ISR啊,具体多少得看你的实际情况,我只是随意写的个数。...

好的 谢谢你
6楼2013-02-27 11:09:25
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

657801288

铁虫 (初入文坛)

引用回帖:
2楼: Originally posted by csgt0 at 2013-02-26 11:06:11
function xwn
options = odeset('RelTol',1e-4,'AbsTol',);
intvalue=;       %t=0时的初值
= ode45(@rigid,,intvalue,options);
plot(T,Y(:,1),'-',T,Y(:,2),'-.',T,Y(:,3),'.')
title('T-Y图')
for i=1: ...

你好,两种方法做出的图怎么不一样呢?求指教。
7楼2013-03-30 08:49:39
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖

csgt0

荣誉版主 (著名写手)

彩色挂图

【答案】应助回帖


fegg7502: 金币+1, 应助指数+1, 鼓励交流 2013-04-02 09:30:38
不是两种方法啊,前面画的是T-Y图,后面是你要的T-dY图。
showmethemoney
8楼2013-04-01 10:48:30
已阅   回复此楼   关注TA 给TA发消息 送TA红花 TA的回帖
相关版块跳转 我要订阅楼主 xwndf250 的主题更新
最具人气热帖推荐 [查看全部] 作者 回/看 最后发表
[考研] 335求调剂 +3 yuyu宇 2026-03-23 4/200 2026-03-23 19:03 by macy2011
[考研] 336求调剂 +4 收到VS 2026-03-20 4/200 2026-03-23 19:02 by macy2011
[考研] 一志愿中国石油大学(华东) 本科齐鲁工业大学 +4 石能伟 2026-03-17 4/200 2026-03-23 17:51 by 17862566385
[考研] 一志愿重庆大学085700资源与环境,总分308求调剂 +6 墨墨漠 2026-03-23 7/350 2026-03-23 16:09 by 墨墨漠
[考研] 298求调剂 +8 上岸6666@ 2026-03-20 8/400 2026-03-23 11:02 by laoshidan
[考研] 石河子大学(211、双一流)硕博研究生长期招生公告 +3 李子目 2026-03-22 3/150 2026-03-22 21:01 by 怎么释怀
[考研] 287求调剂 +8 晨昏线与星海 2026-03-19 9/450 2026-03-22 17:01 by i_cooler
[考研] 求调剂 +6 十三加油 2026-03-21 6/300 2026-03-22 17:00 by i_cooler
[考研] 一志愿 西北大学 ,070300化学学硕,总分287,双非一本,求调剂。 +3 晨昏线与星海 2026-03-20 3/150 2026-03-22 16:00 by ColorlessPI
[考研] 求调剂 +5 Zhangbod 2026-03-21 7/350 2026-03-22 13:13 by Zhangbod
[考研] 354求调剂 +7 Tyoumou 2026-03-18 10/500 2026-03-22 11:11 by 人来盛
[考研] 求调剂 +3 13341 2026-03-20 3/150 2026-03-21 18:28 by 学员8dgXkO
[考研] 307求调剂 +3 余意卿 2026-03-18 3/150 2026-03-21 17:31 by ColorlessPI
[考研] 303求调剂 +5 睿08 2026-03-17 7/350 2026-03-21 03:11 by JourneyLucky
[考研] 材料 336 求调剂 +3 An@. 2026-03-18 4/200 2026-03-21 01:39 by JourneyLucky
[考研] 304求调剂 +6 曼殊2266 2026-03-18 6/300 2026-03-21 00:32 by JourneyLucky
[考研] 295求调剂 +4 一志愿京区211 2026-03-18 6/300 2026-03-20 23:41 by JourneyLucky
[考研] 广西大学家禽遗传育种课题组2026年硕士招生(接收计算机专业调剂) +3 123阿标 2026-03-17 3/150 2026-03-20 15:58 by 飞行琦
[考研] 考研求调剂 +3 橘颂. 2026-03-17 4/200 2026-03-17 21:43 by 有只狸奴
[考研] 考研调剂 +3 淇ya_~ 2026-03-17 5/250 2026-03-17 09:25 by Winj1e
信息提示
请填处理意见