数学模型与数学建模4.2-常微分方程组模型

数学模型安徽大学数学科学学院4.2常微分方程组模型由几个微分方程联立而成的方程组称为微分方程组。本节建立传染病模型、男生追女生模型和种群增长模型，并判断方程组的稳定性。我国目前的法定传染病有甲、乙、丙三类，共39种。传染病的特点是有病原体，有传染性和流行性，感染后常有免疫性。有些传染病还有季节性或地方性。传染病的分类尚未统一，有人按病原体分类，有人按传播途径分类。传染病的预防应采取以切断主要传播环节为主导的综合措施。传染病的传播和4.2.1传染病模型流行必须具备3个环节，即传染源（能排出病原体的人或动物）、传播途径（病原体传染他人的途径）及易感者（对该种传染病无免疫力者）。若能完全切断其中的一个环节，即可防止该种传染病的发生和流行。各种传染病的薄弱环节各不相同。在预防中应充分利用。除主导环节外对其他环节也应采取措施，只有这样才能更好地预防各种传染病。不同类型的传染病，其传播过程有着各自不同的特点，了解这些具体的传染病的传播过程需要了解其病理知识，这里不可能从医学角度一一进行分析，而主要按照一般的传播机理建立几类一般的传染病模型分析受感染人数的变化规律，讨论终止传染病蔓延的方法和手段。模型1.用表示时刻的病人数量，假设病人一旦与健康人群接触就会使健康人群患病，且单位时间内每个病人能够使健康人患病的人数为。初始时刻的病人数为。利用微元法，考虑内病人数的变化，则有，上式两边同时除以，并令，得到（4.2.1）()xtt0(0)xx[,]ttt()()()xttxtxttt0t0(0)dxdtxxx；式（4.2.1）的求解程序及结果如下：dsolve('Dx-lambda*x=0','x(0)=x0','t')ans=x0*exp(lambda*t)即式（4.2.1）的解为（4.2.2）式（4.2.2）为指数形式，故称模型（4.2.1）为指数增长模型。根据式（4.2.2），当时，，即所有人都会患病。很显然，这与事实不符。0()txtxetx模型2.（SI模型）考虑以下假设：（1）不考虑人口的出生、死亡和迁移等种群动力因素，在疾病传染期所考察地区内总人数始终保持常数不变；（2）人群分健康人群（易感染者，Susceptible）和病人（已感染者，Infective），时刻这两类人在总人数中所占比例分别为和，即有，并设初始时刻患病人数比例为；t()st()it()()1stit0(0)ii（3）病人一旦与健康人群接触，就必然具有一定的传染性。假设时刻单位时间内每个病人能传染的易感染者数量与该地区易感染者比例成正比，比例系数为，即时刻单位时间内每个病人可使个易感染者患病，而病人总数为，因而时刻单位时间内共有个易感染者患病。利用微元法，内患病人数的变化量为：t()stt()st()itNt()()stitN[,]ttt()()()()NittitstitNt上式两段同时除以，并令，有。结合假设（2），有（4.2.3）先观察的图像，建M文件myfun1.m。命令如下：functiony=myfun1(i);y=0.01*i*(1-i);再输入命令：fplot('myfun1',[0,1]);fplot函数的命令格式为fplot(‘fun’,lims)，表t0tdidtsi0(1);(0)didtiiii/didt示绘制字符串fun制定的函数在lims=[xmin,xmax]或lims=[xmin,xmax,ymin,ymax]上的图形。但是fun必须是M文件的函数名或是独立变量为x的字符串。图形如下：图4.2.1图像由图4.2.1可以看出，当时，达到最大值，此时病人数增长得最快，意味着传染病高潮的到来。模型（4.2.3）的求解程序如下：symsialphat;dsolve('Di-alpha*i*(1-i)','i(0)=i0','t')结果为：ans=1/(1-exp(-alpha*t)*(-1+i0)/i0)0.5i/didt这里syms是Matlab中定义多个变量的。根据运行结果，有式（4.2.3）的解：（4.2.4）由式（4.2.4），当时，，这是传染病最高峰时刻，此时与成反比。因为反映了单位时间每个病人的传染数量，所以也被称为单位时间接触率，它直接反映了当地的卫生000()11(1)/1111ttiteiiie0.5i0ln11tit水平可以推迟传染病高潮的到来。函数的作图命令及图形如下（这里取，）：fplot('1/(1-exp(-0.5*t)*(-1+0.01)/0.01)',[0,30]);xlabel('t');ylabel('i');0.500.01i图4.2.2患病人数比例的图像由图4.2.1可以看出，单调递增，且在时下凸，而在时上凸，且有渐近线，即病人数量一直递增，在病人比例不到时增长速度非常快，而在病人比例超过时增长速度放缓，最终病人的比例将接近于100%。显然这与事实不符，究其原因，是因为模型（4.2.3）中没有考虑病人被治愈的情况。()it0.5i0.5i1i模型3.（SIS模型）增加假设：假设被治愈的人不具有免疫能力，将成为易感染者，且时刻单位时间内治愈者变化率与病人数量成正比，比例系数为，利用微元法，内患病人数的变化量为，故有(4.2.5)t[,]ttt()()()()()NittitstitNtitNt0(1);(0)didtisiiiiii显然，当，即没有病人治愈时，模型(4.2.5)变为模型(4.2.3)。令，则有（4.2.6）因为反映了单位时间每个病人的传染数量，而反映了单位时间内每个病人的治愈率，则可表示为平均每个病人的治愈时间，为传染期内平均每个病人的有效传染数量。0/(11/)didtii1//图4.2.3（a）时的图像图4.2.3（b）时的图像由图4.2.3（a）和（b）可以看出，当时，随着的增大呈现先增后减的趋势，且在处为0，即此时达到最高峰；而当时，小于0且随着的增大一直递减，即单调递减且递减的速度越来越快。模型（4.2.6）的程序为：dsolve('Di=-alpha*i*(i-(1-1/sigma))','i(0)=i0','t');1/didti11/()it1/didti()it运行结果为：ans=(sigma-1)/(sigma-exp(-alpha*(sigma-1)*t/sigma)*(-sigma+1+i0*sigma)/i0/(sigma-1)*sigma+exp(-alpha*(sigma-1)*t/sigma)*(-sigma+1+i0*sigma)/i0/(sigma-1))图4.2.4（a）时患病人数比例的图像图4.2.4（b）时患病人数比例的图像图4.2.4（c），，时患病人数比例的图形。110.5200.6iiii图4.2.4（a）为，，时的图形，图4.2.4（b）为，，时的图形，图4.2.4（c）为，，时的图形。可以看出，当时，从点出发，单调递减趋于0，即此时病人数会越来越小最终趋于0；当时，的单调性取决于的大小，但是不论初始值还是，均有:0.50.200.01i0.5200.0111/i0.5200.611/i1()it0i1()it0i011/i011/ilim()11/tit该极限值随着的增大而增大，也就是说传染期内平均每个病人的有效传染数量越多，最后患病人数也会越多。模型3中假设治愈者还可以被传染，实际上，大多数传染病如流感、肝炎、天花等治愈后均有一定的免疫力，治愈者既不属于健康人群（易感染者），也不属于患病人群（已感染者）。因此需要单独考虑这类人群。模型4.（SIR模型）假设（1）总人数不变，人群分为健康人群、患病人群和病愈具有免疫能力的移出者（Removed）。时刻这三类人在总人数中所占比例分别为、和，即有，并设初始时刻患病人数比例为，，；Nt()st()it()rt()()()1stitrt0(0)0ii0(0)0ss0(0)0rr（2）单位时间内每个病人的接触率为，单位时间内每个病人的治愈率为，为传染期内平均每个病人的有效传染数量；利用微元法，内移出者的变化量为,等式两边同除以，并令，则有（4.2.7）/[,]ttt()()()rttrtNitNtNt0tdrdti故有SIR模型：（4.2.9）且因模型（4.2.9）无法得到解析解，因此考虑作数值分析，先建立M文件：functiony=ill(t,x);alpha=1,beta=0.3;y=[alpha*x(1)*x(2)-beta*x(1),-alpha*x(1)*x(2)]';000,(0),(0);,(0)0.dsdtsissdidtisiiidrdtirr；　()()()1.stitrt这里假设，，然后输入命令：ts=0:50;x0=[0.02,0.98];[t,x]=ode45('ill',ts,x0);plot(t,x(:,1),t,x(:,2)),grid,pause;plot(x(:,2),x(:,1)),grid得到、以及图形如下：10.3()st()it()is图4.2.5（a）、的图像图4.2.5（b）的图像()st()it()is下面分析模型（4.2.9）中的、和的变化情况：（1）由知恒有，即单调递减，且对任意时刻有；（2）可表示为，若，由知恒有，即患病人数始终不增加，故传染病不再传染。也就是说只()st()it()rt/dsdtsi/0dsdt()st0()sts/didtisi/(1/)didtis/01/s0()sts/0didt有在时传染病才传播，因此反映了传染病能否传播的一个临界值。若存在某时刻，使得，因为为单调递减函数，所以当时，，当时，，即有（4.2.10）01/s1/0t0()1/st()st0tt0()1/st0tt0()1/st0000,(1/)0,0,.ttdidtistttt由式（4.2.10）知，在之前，患病人数在增加，而在之后患病人数在减少，即传染病传播基本得到控制。（3）由（4.2.9）知，若存在使得，而时总有，故有，进而有，即传染病传播完全终止，此时所有的患病者全部被移出，故传染病传播的最终结果0t0t1t1()0it1tt()0it110ttttdsdtdidt10ttdrdt是，.（4）将（4.2.9）中的第一个方程和第三个方程相除得到（4.2.11）求解得到.令，得到，该式表明传染病终止后并非所有的人都要患病，然后治愈成为移出者，总有未曾患病的人存在。1()0it11()()1strt00/;(0);(0)0.dsdrsssrr0rsse1tt1()000rtssese（5）由得到（4.2.12）式（4.2.12）反映了传染病传播终止后，被移出人数占群体总人数的比率，可作为衡量传染病传染强度的一种指标，该比值越大传染强度也越大。11()()1strt11()1()rtst（6）将式（4.2.9）中的第二个方程和第一个方程相除得到（4.2.13）对式（4.2.13）进行求解，得到（4.2.14）即（4.2.15）00/1/()1;();/didssisi000lniissss000()()lnitstisss令，且有，则满足.即：（4.2.16）而由式（4.2.13）知，结合（4.2.10）和（4.2.15），当时，患病人数达到最大值（4.2.17）1tt1()0it1()s