数学建模大作业.doc

数学实验报告实验名称数学建模与MATLAB 学 院 材料学院 专业班级材料1014 姓名徐萌 孔德成 戴思雨 学号41071046 41030400 410303992012年6月141、 问题的提出。 传染病是当今世界最严重的疾病之一，2009年4月26日世界卫生组织以确认，美国和墨西哥发生了甲型H1N1流感，随后疫情迅速蔓延，截止8月中旬，全球感染人数约5万人。因此，运用传染病的数学模型来描述传染病甲型H1N1流感的传播过程，分析受感染人数的变化规律，探索制止甲型H1N1蔓延的手段是值得关注的。2、 模型的建立。 考查中国内地疫情变化，在疾病传播期间不考虑人口的出生率和死亡率，人口总数不变，为常量。中国的疫情研究发现易感染人数多为2050岁的青壮年，故保守估计在此传染病系统的人数N=50000人。甲型HINI流感的传播途径是与病源的直接接触，患者与健康者接触时，都使健康者感染病变故将人群分为3类：健康者(易感染者人群)、患者(已被感染人群)、治愈者(研究期间6月14日8月14日间中国内地感染病毒死亡人数为0，故此处不考虑死亡者)三者在总人数中的比例分别为 : s(t),i(t),r(t)且s(t)+i(t)+r(t)=1,io,So分别为患者人数，健康人数的比例初始值 设每个患者每日感染健康者的平均人数为日感染率，记为j，则j=j日新增病例数/(j-1)日（累计确诊人数-累计出院人数）； 每日被治愈的患者人数占其总数的比例为日治愈率，记为j，则j=j日被治愈的人数/j日累计确诊病人数; 定义整个传染期内每个患者有效接触的平均人数为接触数,由s(t)+i(t)+r(t)=1可知，对于病愈免疫的治愈者而言应有dr/dt=i,因此考虑SIR传染模型，该模型的方程为 di/dt=si-i; ds/dt=-si （1） i(0)=io, s(0)=So;三、模型的求解 1、数值运算 由于在方程(1)中无法求出s(t)和i(t)的解析解，故先做数值运算据来自中国卫生部网站公布的2009年6月14日8月14日的疫情数据(见表1)包括日累计确诊病例、日累计治愈病例等其中缺失的部分数据，将以通过给定的数据拟合得到. 表1疫情原始数据日期新增病例确诊病例累计治愈累计新增治愈数6月14日20185736月15日4122686136月16日1123797116月17日27264114176月18日33297135216月19日31328160256月20日28356185256月21日58414199146月22日27441227286月23日49490251246月24日38528275246月25日42570321466月26日48618338176月27日60678373356月28日51729401286月29日37766445446月30日44810496517月1日56866554587月2日49915612587月3日45960660487月4日401000704447月5日401040749457月6日571097793447月7日541151870777月8日361187927577月9日361223985587月10日4012631035507月11日3913021085507月12日2613281110257月13日2613541134247月14日4513991166327月15日4514441197317月16日4114851230337月17日5215371263337月18日4415811293307月19日4416251323307月20日4316681355327月21日5217201404497月22日5217721454507月23日3818101529757月24日4218521604757月25日2618781663597月26日2619041722597月27日2619301781597月28日3719671817367月29日3620031853367月30日4320461883307月31日4420901912298月1日2021101937258月2日2121311962258月3日2121521988268月4日2921812031438月5日2922102074438月6日2722372098248月7日2722642122248月8日2822922137158月9日2823202152158月10日2823482167158月11日3823862203368月12日3924252240378月13日5724822261218月14日552537228322注：2009年疫情效据见文献8以6月15日为基日，当日累计确诊病例226例，累计出院者86例，故 s(0)=(50000-226+86)/50000=0.9972; I(0)=(226-86)/50000=0.0028;在研究期间，平均日感染率和平均日治愈率由每天相应数据平均求得设计程序为： 新增病例A 确诊病例累计B 治愈累计C 新增治愈数DA=41 11 27 33 31 28 58 27 49 38 42 48 60 51 37 44 56 49 45 40 40 57 54 36 36 40 39 26 26 45 45 41 52 44 44 43 52 52 38 42 26 26 26 37 36 43 44 20 21 21 29 29 27 27 28 28 28 38 39 57 55B=226 237 264 297 328 356 414 441 490 528 570 618 678 729 766 810 866 915 960 1000 1040 1097 1151 1187 1223 1263 1302 1328 1354 1399 1444 1485 1537 1581 1625 1668 1720 1772 1810 1852 1878 1904 1930 1967 2003 2046 2090 2110 2131 2152 2181 2210 2237 2264 2292 2320 2348 2386 2425 2482 2537C=86 97 114 135 160 185 199 227 251 275 321 338 373 401 445 496 554 612 660 704 749 793 870 927 985 1035 1085 1110 1134 1166 1197 1230 1263 1293 1323 1355 1404 1454 1529 1604 1663 1722 1781 1817 1853 1883 1912 1937 1962 1988 2031 2074 2098 2122 2137 2152 2167 2203 2240 2261 2283D=13 11 17 21 25 25 14 28 24 24 46 17 35 28 44 51 58 58 48 44 45 44 77 57 58 50 50 25 24 32 31 33 33 30 30 32 49 50 75 75 59 59 59 36 36 30 29 25 25 26 43 43 24 24 15 15 15 36 37 21 22E=A./(B-C) %日感染率e=sum(E)/61 %平均日感染率F=D./(B-C) %日治愈率f=sum(F)/61 %平均日治愈率运行结果： A = Columns 1 through 16 41 11 27 33 31 28 58 27 49 38 42 48 60 51 37 44 Columns 17 through 32 56 49 45 40 40 57 54 36 36 40 39 26 26 45 45 41 Columns 33 through 48 52 44 44 43 52 52 38 42 26 26 26 37 36 43 44 20 Columns 49 through 6121 21 29 29 27 27 28 28 28 38 39 57 55B = Columns 1 through 8 226 237 264 297 328 356 414 441 Columns 9 through 16 490 528 570 618 678 729 766 810 Columns 17 through 24 866 915 960 1000 1040 1097 1151 1187 Columns 25 through 32 1223 1263 1302 1328 1354 1399 1444 1485 Columns 33 through 40 1537 1581 1625 1668 1720 1772 1810 1852 Columns 41 through 48 1878 1904 1930 1967 2003 2046 2090 2110 Columns 49 through 56 2131 2152 2181 2210 2237 2264 2292 2320 Columns 57 through 612348 2386 2425 2482 2537C = Columns 1 through 8 86 97 114 135 160 185 199 227 Columns 9 through 16 251 275 321 338 373 401 445 496 Columns 17 through 24 554 612 660 704 749 793 870 927 Columns 25 through 32 985 1035 1085 1110 1134 1166 1197 1230 Columns 33 through 40 1263 1293 1323 1355 1404 1454 1529 1604 Columns 41 through 48 1663 1722 1781 1817 1853 1883 1912 1937 Columns 49 through 56 1962 1988 2031 2074 2098 2122 2137 2152 Columns 57 through 61 2167 2203 2240 2261 2283D = Columns 1 through 16 13 11 17 21 25 25 14 28 24 24 46 17 35 28 44 51 Columns 17 through 32 58 58 48 44 45 44 77 57 58 50 50 25 24 32 31 33 Columns 33 through 48 33 30 30 32 49 50 75 75 59 59 59 36 36 30 29 25 Columns 49 through 61 25 26 43 43 24 24 15 15 15 36 37 21 22e = 0.173970451885603 %平均日感染率=0.173970451885603 f = 0.164030384929960 %平均日治愈率=0.164030384929960 接触数： =/=0.173970451885603/0.164030384929960=1.060598936958463可得模型方程为： di/dt=0.173970451885603si-0.164030384929960i; ds/dt=-0.173970451885603si ； i(0)=0.9972, s(0)=0.0028; 然后用Matlab软件编程，解常微分方程做出患者人数比例i(t)-时间t/d关图，健康者比例s(t)-时间t/d关系图 ，患者人数比例i-健康者比例s图 。 采用常微分方程数值求解的方法，预测之后的疫情，并绘制出曲线。设计程序为：M文件function wo=bean(t,x)wo=-0.17397045*x(2） 0;0.17397045*x(2) -0.16403038*x;命令：t,x=ode45(bean,0,460,0.9972,0.0028)t= Columns 1 through 6 0 5.3145 10.6290 15.9435 21.2580 32.7580 Columns 7 through 12 44.2580 55.7580 67.2580 78.7580 90.2580 101.7580 Columns 13 through 18 113.2580 124.7580 136.2580 147.7580 159.2580 170.7580 Columns 19 through 24 182.2580 193.7580 205.2580 216.7580 228.2580 239.7580 Columns 25 through 30 251.2580 262.7580 274.2580 285.7580 297.2580 308.7580 Columns 31 through 36 320.2580 331.7580 343.2580 354.7580 366.2580 377.7580 Columns 37 through 42 389.2580 400.7580 412.2580 423.7580 435.2580 441.4435 Columns 43 through 45 447.6290 453.8145 460.0000x = Columns 1 through 3 0.9972, 0.0028 0.9946, 0.0020 0.9918, 0.0031 Columns 4 through 6 0.9889, 0.0032 0.9859, 0.0034 0.9790, 0.0036 Columns 7 through 9 0.9717, 0.0039 0.9640, 0.0041 0.9560, 0.0042 Columns 10 through 12 0.9479, 0.0043 0.9398, 0.0043 0.9318, 0.0042 Columns 13 through 15 0.9240, 0.0041 0.9166, 0.0039 0.9096, 0.0037 Columns 16 through 18 0.9031, 0.0034 0.8972, 0.0032 0.8918, 0.0029 Columns 17 through 21 0.8869, 0.0026 0.8826, 0.0023 0.8788, 0.0020 Columns 22 through 24 0.8755, 0.0018 0.8726, 0.0015 0.8700, 0.0013 Columns 25 through 27 0.8679, 0.0012 0.8660, 0.0010 0.8644, 0.0009 Columns 28 through 30 0.8631, 0.0007 0.8619, 0.0006 0.8609, 0.0005 Columns 31 through 33 0.8601, 0.0004 0.8594, 0.0004 0.8588, 0.0003 Columns 34 through 36 0.8583, 0.0003 0.8578, 0.0002 0.8575, 0.0002 Columns 37 through 39 0.8572, 0.0002 0.8569, 0.0001 0.8567, 0.0001 Columns 40 through 42 0.8565, 0.0001 0.8564, 0.0001 0.8563,0.0001 Columns 43 through 45 0.8562, 0.0001 0.8562, 0.0001 0.8561, 0.0001 i=x(:,2); plot(t,i,g-) xlabel(时间t/d),ylabel(患者比例i) title(i-t关系图) s=x(:,1); plot(t,s,k-) xlabel(时间t/d),ylabel(健康者比例s) title(s-t关系图) plot(s,i,r-) xlabel(健康者比例s),ylabel(患者比例i) title(i-s关系图) y=dsolve(Dy=1/x/1.06058937-1,y(0.9972)=0.0028,x)y=(235717995174701779256942769471657*log(x)/250000000000000000000000000000000-(235717995174701779256942769471657*log(2493/2500)/250000000000000000000000000000000 - x + 1 即为i=1/*ln(s/s0)-s+s0+i0 =1/1.06059894*ln(s/0.9972)-s+1 图一 患者比例i-时间t关系图图二 健康者比例s-时间t关系图 图三 i-s关系图由以上分析可知：(1) 不论初始条件So 、io如何，患者人数终将消失，即i=0从相轨线图形上看，无论从哪一点出发，i(t)终将与s轴相交(t充分大).(2) 最终未感染者的比例是S，在式(3)中令i=0，得到S是方程 (So+io)一S+1/*lnS/So=0 (4)在(0,1/)内的根，在图形上S是相轨线与s轴在(0,1/)内交点的横坐标，此例中S约为0.85(3)若So1/，则i(t)先增加，当S=1/时，i(t)达到最大值 I(max)=So+io-1/(1+InSo) (5)然后i(t)减小且趋于0，s(t)则单调减小至S；若So1/时传染病就会蔓延而减小传染期接触数即提高定值1/，使得知So1/，从式(4)、式(5)可以看出，减小时，S增加，I(max)降低，也控制了蔓延的程度 另外，s=s/是传染期内的一个患