当前位置:   article > 正文

传染病模型中作图与计算(matlab,数学模型)_1、在sis模型中设定不同的参数a,μ和初始值i(0),编程计算i(t)并作图,分析曲线

1、在sis模型中设定不同的参数a,μ和初始值i(0),编程计算i(t)并作图,分析曲线

一、实验要求

1.取k=0.1,画出dt/di~i的曲线图,求i为何值时dt/di达到最大值,并在曲线图上标注。试编写一个m文件来实现。

2.求出微分方程的解析解 i(t),画出如下所示的 i~t 曲线(i(0)=0.15, k=0.2, t=0~30)。试编写一个 m 文件来实现。(在图形窗口菜单选择 Edit/Copy Figure,复制 图形)

3. 取λ=0.1,σ=1.5,画出如下所示的dt/di~i曲线图。试编写一个m文件来实现。(在图形窗口菜单选择 Edit/Copy Figure,复制图形)

4.求出微分方程的解析解 i(t)。取λ=0.2,σ=3,t=0~40,画出如下所示的图形。试编写 一个m文件来实现。

5.设λ=1,μ=0.3,i(0)=0.02,s(0)=0.98。编写程序,求出(i(t),s(t)的数值计算结果;i(t),s(t)图形;i~s图形(相轨线))。

二、实验内容

1.

  1. clear
  2. clc
  3. fplot('0.1*x*(1-x)',[0 1.1 0 0.03]);
  4. x=fminbnd('-0.1*x*(1-x)',0,1);
  5. y=0.1*x*(1-x);
  6. hold on; %在上面的同一张图上画线(同坐标系)
  7. plot([0,x],[y,y],':',[x,x],[0,y],':');
  8. text(0,y,'(di/dt)m','VerticalAlignment','bottom');
  9. text(x,-0.001,num2str(x),'HorizontalAlignment','center');
  10. title('SI模型di/dt~i曲线');
  11. xlabel('i');
  12. ylabel('di/dt');

2.

  1. clear
  2. clc
  3. k=0.2;
  4. x0=0.15;
  5. x=dsolve('Dx=k*x*(1-x)','x(0)=x0');
  6. tx=0:0.15:30;
  7. for j=1:201
  8. t=tx(j);
  9. i(j)=eval(x);
  10. end
  11. plot(tx,i);
  12. title('图1SI模型的i~t曲线');
  13. xlabel('t(天)');
  14. ylabel('i(病人所占比例)');

3.

  1. clear
  2. clc
  3. fplot('-0.1*x*(x-(1-2/3))',[0 0.4 -0.0005 0.003]);
  4. hold on; %在上面的同一张图上画线(同坐标系)
  5. plot([0,0.4],[0,0]);
  6. title('SIS模型的 di/dt~i曲 线');
  7. xlabel('i');
  8. ylabel('di/dt');

4.

  1. clear
  2. clc
  3. k=0.2;
  4. x0=0.2;
  5. x=dsolve('Dx=-k*x*(x-(1-1/3))','x(0)=x0');
  6. tx=0:0.15:40;
  7. for j=1:267
  8. t=tx(j);
  9. i(j)=eval(x);
  10. end
  11. plot(tx,i);
  12. hold on; %在上面的同一张图上画线(同坐标系)
  13. plot([0,40],[2/3,2/3],'k:');
  14. x0=0.9;
  15. x=dsolve('Dx=-k*x*(x-(1-1/3))','x(0)=x0');
  16. tx=0:0.15:40;
  17. for j=1:267
  18. t=tx(j);
  19. i(j)=eval(x);
  20. end
  21. plot(tx,i,'r');
  22. title('SIS模型的 di/dt~i曲 线');
  23. xlabel('i');
  24. ylabel('di/dt');
  25. legend('i(0)=0.2','1-1/σ','i(0)=0.9');
  26. xlim([0 40]);
  27. ylim([0 1]);

5.

(此代码第二张图显示需要在命令栏输入任意字符)

  1. function y=f(t,x)
  2. y=[1*x(1)*x(2)-0.3*x(1),-1*x(1)*x(2)]';
  1. clear
  2. clc
  3. t=0:50;
  4. x=[0.02,0.98];
  5. [t,x]=ode45('f',t,x);
  6. [t,x]
  7. plot(t,x);
  8. pause
  9. plot(x(:,2),x(:,1));

三、实验结果及其分析

1.

2.

 

3.

4.

 

5.

 

         0    0.0200    0.9800

    1.0000    0.0390    0.9525

    2.0000    0.0732    0.9019

    3.0000    0.1285    0.8169

    4.0000    0.2033    0.6927

    5.0000    0.2795    0.5438

    6.0000    0.3312    0.3995

    7.0000    0.3444    0.2839

    8.0000    0.3247    0.2027

    9.0000    0.2863    0.1493

   10.0000    0.2418    0.1145

   11.0000    0.1986    0.0917

   12.0000    0.1599    0.0767

   13.0000    0.1272    0.0665

   14.0000    0.1004    0.0593

   15.0000    0.0787    0.0543

   16.0000    0.0614    0.0507

   17.0000    0.0478    0.0480

   18.0000    0.0371    0.0460

   19.0000    0.0287    0.0445

   20.0000    0.0223    0.0434

   21.0000    0.0172    0.0426

   22.0000    0.0133    0.0419

   23.0000    0.0103    0.0415

   24.0000    0.0079    0.0411

   25.0000    0.0061    0.0408

   26.0000    0.0047    0.0406

   27.0000    0.0036    0.0404

   28.0000    0.0028    0.0403

   29.0000    0.0022    0.0402

   30.0000    0.0017    0.0401

   31.0000    0.0013    0.0400

   32.0000    0.0010    0.0400

   33.0000    0.0008    0.0400

   34.0000    0.0006    0.0399

   35.0000    0.0005    0.0399

   36.0000    0.0004    0.0399

   37.0000    0.0003    0.0399

   38.0000    0.0002    0.0399

   39.0000    0.0002    0.0399

   40.0000    0.0001    0.0399

   41.0000    0.0001    0.0399

   42.0000    0.0001    0.0399

   43.0000    0.0001    0.0399

   44.0000    0.0000    0.0398

   45.0000    0.0000    0.0398

   46.0000    0.0000    0.0398

   47.0000    0.0000    0.0398

   48.0000    0.0000    0.0398

   49.0000    0.0000    0.0398

   50.0000    0.0000    0.0398

声明:本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有侵权的内容,请联系我们。转载请注明出处:【wpsshop博客】
推荐阅读
相关标签
  

闽ICP备14008679号