实验传染病模型-微分方程模型

tttt实

验

六实验项：染病模型——微分方程模型实验实验目：进一步巩固、加强微分方程模型的建模、求解能力；2.学习掌握用数学软件包求解微分方程数值解的相关命令。实验内：建模实例，传染病问题等；2.编程求解。一模实-----传病型•

问题:有一种传染病（如SARS甲型H1N1）正在流行。现在希望建适当的数学模型用已经掌握的一些数据资料对该传染病进行有效地研究以期对其传播蔓延进行必要的控制减少人民生命财产的损失考虑如下的几个问题建立适当的数学模型并进行一定的比较分析和评价展望。•1、不虑环境的限制，设单位时间内感染人数的增长率是常数，建立模型求时刻的感染人数。•2假设环境条件下所允许的最大可感染人数为m单位时间内感染人数的增长率是感染人数的线性函数最大感染时的增长率为零建立模型求时刻的感染人数。实方与骤1、题析、这是一个涉及传染病传播情况的实际问题，其中涉及传染病感染人数随时间的变化情况及一些初始资料，可通过建立相应的微分方程模型加以解决。b、问题表中已给出了各子问题的一些相应的假设。c、在实际中，感染人数是离散变量，不具有连续可微性，不利于建立微分方程模型但由于短时间内改变的是少数人口这种变化与整体人口相比是微小的因此为了利用数学工具建立微分方程模型我们还需要一个基本假设：感染人是时的连续微函2、题解2.1问题一的解答—模型A、型假、感染人数是时间的连续可微函数；2)、单位时间内感染人数的增长率是常数，或单位时间内感染人数的增长量与/

当时的感染人数成正比。B、模型构成设t时刻的感染人数为x(t)始时刻

t0

的感染者人数为

，感染者的增长率为

r

，根据单位时间内感染人数的增长率是常数的假设，到t

时间内感染人数的增量为：t)(t)因此，

xt)

满足如下的微分方程：

,

x(0)xC、型求MATLAB算求解（介绍完MATLAB求解微分方程数值解的相关命令后再运行）x=x0*exp(r*t)即

(t)rt

tD、型分由上述解的形式，可以看出，感染人数将随着时间的增长按指数规律无限增长。特别地当时间趋向于无穷时感染人数也将趋向于无穷大这显然是不符合现实的说明该模型不可能用于传染病的长期预报同时也说明迫切需要对该模型进行必要的修正。E、改进方向单位时间内感染人数的增长率不是常数而是逐渐下降的原因感染人数增长到一定数量后，环境条件、人口总数等因素将对感染者数量的增长起阻滞作用，且阻滞作用随感染者数量增加而变大长率是感染人数的减函数染者越多，增长率越低。2．问题二的答——模型A、型假•1)、感染人数是时间的连续可微函；

•2)、感染人数受环境条件的限制，一个最大的可感染人数。•3)、单位时间内感染人数的增长率感染人数有关，是其线性函数，最大感染时对应增长率为零。B、模型构成/

rr0rr0•

t仍然设时刻的感染人数为

x(t

，初始时(

t0

)的感染者人数••

为，感染者人数为0时，感染人数的增长率为。根据单位时间内感染人数的增长率和感染人数有关是其线性函数的假设可得增长率关于感染者人数的线性函数关系式：r()r进一步，由最大感染时对应的增长率为零可确定参数的值为：0mxt)因此，在该模型的假设下，感染人数应满足如下的微分方程：

r(x)x(1)x,dtmx(0)0C、型求解：MATLAB计算求解（介绍完MATLAB求解微分方程数值解的相关命令后再运行）,'x(0)=x0''t')即D、型分

(tm1me

t•根据前述微分方程作出

/dt~

的曲线图，见图1-1，这是一条抛物线由该图可看出感染人数增长率随感染人数的变化规律增长率随着感染人数的增加而先增后减，在

/2

时达到最大。这预示着传染病高潮的到来是医疗卫生部门关注和需要密切注意的时刻因为感染人数增长率在一定程度上代表了医疗卫生水平增长率越小卫生水平越高所以改善保健设施、提高卫生水平可以推迟传染病高潮的到来。/

tt.根据模型求解得到的结果作出~曲线见图这是一条S曲线。由该图可看出感染人数随时间的变化规律以看出间趋于无穷时(t)

趋于

，且对一切

t

，

()

。此性质说明感染者数量不可能达到最大容量，但可无限趋近于最大容量。二利MATLAB求微方数解相命1指令数及调格式1.1指令函数：注：此指令函数用于求解微分方程(组)的符号(解析)解。1.2单变量常微分方程的调用格式：（,‘cond’‘v’)注：此调用格式用于求符号微分方程的通解或特解，其中eq代表微分方程，cond代表微分方程的初始条件（若缺少，则求微分方程的通解为指定自变量（如未指定，系统默认t自变量1.3常微分方程组的调用格式：/

’,‘eq2’,…‘’‘’,‘’…‘condn,‘v1’,‘v2’,…,vn’)注：此调用格式用于求解符号常微分方程组。其eq1…，代表n微分方程构成的微分方程组；cond1，…，表微分方程组的初始条件（若缺少，则求微分方程组的通解v1，，指定自变量（如未指定，系统默认t为自变量1.4记述规定：MATLAB中，用D(注意：一定是大写记述微分方程中函数的导数。当y是因变量时用‘表示‘y的n阶导数’如Dy表示y一阶导数y'，Dny表示yn导数。Dy(0)=5表示y(0)=5。表示微分方程y'''+y''+y'-x+5=0。2实例示例1、求微分方程

'xy

的通解命令输入：y=dsolve('Dy+2*x*y=2*x*exp(-x^2)','x')得结果为：y=若输入命令：>>y=dsolve('Dy+2*x*y=2*x*exp(-x^2)')则系统默认t为自变量而把真正的自变量x当作常数处理把y当作的函数，得到错误的结果：y=exp(-2*x*t-x*(x-2*t))+exp(-2*x*t)*C1例2、求微分方程命令输入：

2x2dt

的通解得结果为：x=C1*exp(5/2*t)+C2*exp(5/2*t)*t%系统默认为自变量例3、求微分方程下的特解。

y''y'

在条件

y

x

6,y

x

命令输入：>>y=dsolve('D2y+5*Dy-4*y+10=0','y(0)=6','Dy(0)=4','x')得结果为：y=41^(1/2))*x)*(7/4-51/164*41^(1/2))+5/2/

例4、求下述微分方程组的解dy

,x.命令输入：[xy]=dsolve('Dx=-3*y',得结果为：x=1/2*6^(1/2)*(C1*cos(6^(1/2)*t)-C2*sin(6^(1/2)*t))y=C1*sin(6^(1/2)*t)+C2*cos(6^(1/2)*t)3上机练习（可选择其中一部分习题）1、将讲授的例子中的命令输入MATLAB命令窗口，执行命令，观察输出结果并体会MATLAB在求解微分方程符号解方面的功能。2、解下列分方程（组（1du