[文学]实验2差分方程和数值微分

实验2 差分方程和数值微分题目3某湖泊每天有的河水流入，河水中污物浓度为，经渠道排水后湖泊容积保持不变，现测定湖泊中污物浓度为，建立差分方程计算湖泊中1年内逐月（按30天计算）下降的污物浓度，问要多长时间才能达到环保要求的浓度？为了把这个时间缩短为一年，应将河水中污物浓度降低到多少？【模型及其求解】记湖泊总水量，每天流入湖泊的河水量，湖泊最初污物浓度，目前河流的污物浓度，假定外部环境没有任何其他的污染物进入湖泊，并且湖泊本身不会制造污染物，流入湖泊的河水的污物浓度不变，不考虑每天任意时刻河水与湖水的混合过程，另外，假定湖很大，河流入水口与湖泊出水口相距很远（认为当天从湖泊流出的水的浓度是不变的），同时，湖水的总量任意时刻都保持不变（不会出现先蓄积再排出的现象）。基于以上假定，则第k天湖水的污物浓度只和第k-1天湖水的污物浓度及流过的河水污物浓度有关，可以得到如下关系式：记，则：代入题目中给定的数据，在MATLAB中计算。源程序为：%第二章第3题源程序clear all V0 = 200*104; %参数初始值设定Vin = 104; K = 0; pr0 = 0.02; x0 = 0.2;a = (V0-Vin)/V0; P = 0 0.01 0.02 0.04 0.06 0.08 0.1 0.14 0.16 0.2 0.24 %河水不同的污物浓度for j = 1:11 pr = P(j); b = Vin*pr/V0; x(1) = x0; x(2) = a*x(1) + b; for k = 2 : 720 x(k) = a*x(k-1) + b; %差分关系式 C(j,k) = x(k); %湖水污物浓度变化矩阵 end C(j,1) = x0; if(j = 3) %题目给出的河水污染物浓度下的情况 m = 1:12; B(m+1) = x(30*m); B(1) = x0; endendfigure, plot(0:12,B,b*,0:12,B,r-);%一年内按月份污物浓度变化(pr=0.02) title(图1.一年内湖泊污染物浓度(按月) xlabel(月份) ylabel(污物浓度(g/m3) gtext(pr=0.02);figure, plot(1 : 720,C(1,:),m,1 : 720,C(2,:),b,.%pr取不同值污染浓度曲线比较1 : 720,C(3,:),y,1 : 720,C(4,:),r,1 : 720,C(5,:),g,. 1 : 720,C(6,:),c,1 : 720,C(7,:),k,1 : 720,C(8,:),m,.1 : 720,C(9,:),g,1 : 720,C(10,:),r,1 : 720,C(11,:),k),.axis(0 720 0 0.25); title(图2.河水污物浓度不同时湖泊污染物浓度变化情况(按天); xlabel(天) ylabel(污物浓度(g/m3),gtext(pr=0); gtext(pr=0.01);gtext(pr=0.02);gtext(pr=0.04); gtext(pr=0.06); gtext(pr=0.08); gtext(pr=0.1); gtext(pr=0.14); gtext(pr=0.16); gtext(pr=0.2); gtext(pr=0.21); a = 0.9950; k = 360; xk = 0.04; x1 = 0.2; %如果一年内污染物浓度达要求求解prX b = (1-a)*(xk-x1*a(k-1)/(1-a(k-1); prX = b*V0/Vin;问题解决u Step1：首先，得到一年内112月份的湖水污物浓度数据，如表1。表1. 一年内湖水污物浓度（月份） （浓度单位：g/m3）MATLAB绘制出曲线如图1：由此曲线可以清楚地看到湖水的污染物浓度不断降低，但是一年时间仍然不能使得湖水污物浓度达到环保要求。利用Matlab进一步计算可以得到如表2所示结果：表2. 河水污物浓度0.02g/m3时，湖水污物浓度达标时间 （浓度单位：g/m3）湖泊污物浓度将在第439天达到u Step2：由差分方程逐步递推得到关系式：从而得到b的表达式：按照题目要求，如果要求把湖水污染物浓度降低到环保标准（即）所用时间减少到一年，则可以通过上式求出。参数为=0.9950，k=360（1年），代入式中，求得：再由关系式：计算得到： ，也就是说当河水的污染物浓度为时，按照模型计算，可以在一年时间使得湖水达到环保要求。【结果分析与讨论】在建立本题目的模型过程中，首先作了一些假定，归结起来主要有：1 外部环境没有任何其他的污染物进入湖泊，并且湖泊本身不会制造污染物；2 流入湖泊的河水的污物浓度不变；3 不考虑每天任意时刻河水与湖水的混合过程；4 假定湖很大，河流入水口与湖泊出水口相距很远（认为当天从湖泊流出的水的浓度是不变的）；5 湖水的总量任意时刻都保持不变（不会出现先蓄积再排出的现象）。在这五条假定的前提下，建立了形如：的差分方程数学模型，并对问题相应的求解。利用MATLAB画出了不同的河水污染物浓度pr的情况下，湖水污染物浓度随时间的变化曲线，如图2：从图上曲线的变化趋势可以很清楚的看到，河水的污染物浓度是影响湖水质量的重要因素。在低于湖水的污染物浓度的情况下（），河水的污染物浓度越低，湖水的污染物浓度下降得越快，湖水质量提高的越快；当达到0.2时，湖水的浓度会保持不变；超过0.2时，湖水的污染物浓度会越来越高。以上这三点实验结论都是符合客观事实的。在现实中，河水是不断流动的，而湖水是不流动的，所以湖泊的污染情况直接取决于河水的污染情况，因此该模型能够明确的给出一种湖泊类不流动水系污染治理的思路，即应将治理工作的重点集中在治理同湖泊相连的河流的治理。在实际应用中，该模型还应该考虑到湖泊原本积累的污染物会导致新的污染，从而加剧湖水的污染状况，因此要在差分方程中引入一个阻滞项，这样上面的曲线斜率就会降低，即湖泊污染程度改善速度会降低。题目5据报道，某种山猫在较好、中等及较差的自然环境下，年平均增长率分别为1.68%，0.55%和-4.50%，假定开始时有100只山猫，按以下情况讨论山猫数量逐年变化过程及趋势：（1）3种自然环境下25年的变化过程（作图）；（2）如果每年捕获3只，会发生什么情况？山猫会灭绝吗？如果每年捕获一只呢？（3）在较差的自然环境下，如果使山猫数量稳定在60只左右，每年要人工繁殖多少只？【第(1)问模型建立及求解】记第k年该山猫的数量为xk，自然环境下的年平均增长率为r，且a=1+r，则第k+1年山猫的数量为xk+1=axk+b， a=1+r， k=0,1,2, (1)在较好、中等及较差的自然环境下，以r=0.0168，0.0055，-0.0450以及x0=100，b=0代入，用MATLAB计算并作图，程序如下：%-作业题2_5函数M文件源程序-function x=Exf2_5(x0,n,r,b) % 建立名为Exf2_5的函数M文件， x0,n,r,b可调节a=1+r;x=x0; % 赋初值for k=1:n x(k+1)=a*x(k)+b; % 按照（1）式进行迭代计算end%-作业题2_5脚本 M文件源程序-clear all x0=100;n=25;b=0; % x0:山猫初始数量； n：年数； b:人工繁殖数（若捕获则为负数）k=(0:n);y1=Exf2_5(x0,n,0.0168,b); % 给定x0,n,r,b,调用Exf2_5函数计算y2=Exf2_5(x0,n,0.0055,b);y3=Exf2_5(x0,n,-0.0450,b);round(k,y1,y2,y3), % 对结果四舍五入取整plot(k,y1,k,y2,:,k,y3,-), % 将3条线画在一个图上xlabel(年数)ylabel(山猫数量)title(山猫数量在3种自然环境下25年的变化过程（b=0）) % 加入X轴标记，Y轴标记和标题gtext(r=0.0168), gtext(r=0.0055), gtext(r=-0.0450), % 在图上标记各图线的增长率得到的结果如下面的表1和图1所示：表1：山猫数量在3种自然环境下25年的变化过程（b=0）kxk（r=1.68%，b=0）xk（r=0.55%，b=0）xk（r=-4.50%，b=0）01001001001102101962103101913105102874107102835109103796111103767112104728114104699116105661011810663111201066012122107581312410755141261085215128109501613110948171331104618135110441913711142201401124021142112382214411336231471133524149114332515211532图1：山猫数量在3种自然环境下25年的变化过程（b=0）【第(1)问结果分析及理论计算】由Matlab计算出来的数据表格和变化趋势图可见，无人工干预时，在较好（r=1.68%）和中等（r=0.55%）自然环境下，该山猫数量逐年增多（但在较好环境下增长速度更快）；而在较差（r=-4.50%）自然环境下，山猫数量逐年减少，趋于灭绝。理论分析如下：在自然环境下（1）式的解的形式为xk=akx0， a=1+r， k=0,1,2显然当a1（即r0时），xk，山猫数量逐年增多；而a1（即r0且a1时，xk，表明山猫数量不断增长；当c=x0-b1-a1时，xk-，表明山猫数量逐年减小，濒于灭绝；当a1时，xkb1-a0，表明山猫的数量趋于一个恒定值，但该恒定值在b0的前提下（即每年捕获山猫而不是人工繁殖）是小于零的，所以山猫在该情况下仍然濒于灭绝。分别将具体数值r=0.0168，0.0055，-0.0450以及b=-3，b=-1代入以上条件，可得如下结果：b rr=0.0168r=0.0055r=-0.0450b=-3灭绝灭绝灭绝b=-1灭绝灭绝逐年增长可见，理论分析结果与Matlab计算结果是一致的。【第(3)问理论计算及Matlab求解】由第2问的分析可得，山猫数量满足xk=akx0+b1+a+ak-1=akx0+b1-ak1-a, k=0,1,2=cak+b1-a, k=0,1,2 其中c=x0-b1-a由该式可见，在每年人工繁殖b只山猫的条件下，要使山猫数量稳定在60只左右，需满足a1b1-a=60 ,将r=-0.0450代入，可得b=2.7。取最接近的整数得b=3。用Matlab计算较差环境中b=3条件下的山猫数量变化趋势，并使年数充分大（这里取200），可得上述条件下的山猫变化趋势如图：图6：每年人工繁殖3只山猫时，数量在较差环境下的变化趋势可见，理论计算结果与Matlab计算结果是一致的。取k=200代入，可得x200=67，也与变化趋势图一致。故可得，为使较差环境下山猫数量稳定在60只左右，平均每年需人工繁殖2.7（取整后为3）只山猫。【本题分析、总结】(1)在无人工干预的前提下，该种山猫在较好的和中等的自然环境条件下数量逐年增长，并且在较好自然条件下的增长速度较快。这与直观感觉和理论分析均是一致的。然而在较差的自热环境下，山猫将趋向灭绝，必须加以人工干预或者提高环境质量等措施防止其灭绝。山猫在自然环境中是否濒于灭绝的具体临界条件已在第（1）问的理论计算中加以分析。(2)如果每年捕获三只山猫，无论处于何种自然环境中，山猫都将灭绝，并且环境条件越差，灭绝速度越快；若每年捕获一只山猫，则较好自然条件下，山猫不会灭绝（但增长速度变慢），而中等及较差自然条件下山猫仍将濒于灭绝。山猫面对人类捕获时在各种自然环境中是否灭绝的具体条件已在第（2）问的理论计算中加以分析。(3) 为使较差环境下山猫数量稳定在60只左右，平均每年需人工繁殖2.7（取整后为3）只山猫。再结合前两问的结果可见，在不同的自然环境及人工干预下，山猫的数量变化趋势将大不相同。我们可以根据不同的自然环境通过计算分析，选择合适的人工干预的程度，从而使生态稳定。题目6下表是美国加州和其他各州1955年至1960年人口出生、死亡和迁移（只考虑加州和其他各州之间的净迁移）的情况（单位：千人），根据这些数据建立差分方程模型，计算并预测直至2010年的人口状况。（1）以加州和其他各州的人口数量为对象；（2）以加州和其他各州在美国总人口中的比例为对象，观察趋势。表美国加州和其他各州人口状况【模型假定】建立模型之前，现作如下假定：（以5年为单位）1加州人口增长率和死亡率都是常数；2其他州的人口增长率和死亡率都是常数；3人口的迁移率和加州人口和其他州人口数成正比；【模型建立和求解】记加州人口的增长率为，死亡率为，其他州人口增长率，死亡率；以5年为单位，第k1时段加州人口数量为，其他州人口数量为；第k时段加州人口数量为，其他州人口数量为，迁移率系数q。（1）若以加州人口和其他州人口作为研究对象：可以得到如下差分方程：，（）（2）若以加州和其他州和其他州在美国人口中的比例为研究对象（由于比例相加之和为1，因此只需要研究加州的比例即可）记为加州第k年人口比例，a为该比例的增长系数，则可列差分方程：，（）参数的估计：，使用MATLAB编程解决：（1）若以加州人口和其他州人口作为研究对象：源程序如下：clear allr1 = 0.1315; d1 = 0.0473;r2 = 0.1282; d2 = 0.0488;q = 5.6904*10(-7); a = 1.0842;x0 = 12988; y0 = ; x(1) = x0; y(1) = y0;for k = 1 : 20 x(k+1) = x(k) + (r1 - d1)*x(k) + q*x(k)*y(k); y(k+1) = y(k) + (r2 - d2)*y(k) - q*x(k)*y(k);endB = x;y;n = 1955:5:2055;plot(n,B(1,:),b*,n,B(2,:),r,n,B(1,:),-,n,B(2,:),-), axis(1955 2050 0 6*105);title(美国加州和其他各州人口数量变化)xlabel(时间/年)ylabel(人口数量/人)gtext(加州)gtext(其他各州)运行程序得到曲线如下：从曲线变化趋势可以看到，按照所建立的差分模型，加州的人口将不断增加，而其他州人口将在2020年左右达到峰值，随后将呈现衰减现象。利用MATLAB计算该模型，预测2010年加州人口将达到千人。（2）若以加州和其他州和其他州在美国人口中的比例为研究对象MATLAB源程序如下：clear alla = 1.0842;z0 = 0.0787; % z0 = 12988/z(1) = z0;for k = 1:20 z(k+1) = a*z(k);endn=1955:5:2055;plot(n,z,ro,n,z,-),axis(1955 2050 0 0.4);title(加州人口占美国人口比例变化)xlabel(时间/年)ylabel(人口比例100%)从图上可以清楚看到，按照所建立的模型分析，加州人口的比例将逐渐升高。【结果分析与讨论】在建立模型的过程中，忽略了很多因素，例如人口增长率一般不会为常数，会一定程度上遵循Logistic阻滞增长模型；人口的死亡率也不可能为常数；人口的迁移率于加州和其他各州人口的数量关系定义过于简单等。从相关参数的估计上，只选取了五年的数据，因此带有很大的偶然性，必然会导致预测结果的偏差。题目7与蛛网模型稍有差别，设第k+1与k时段商品上市数量之差是第k时段价格的线性增函数，系数为a；第k+1时段与k时段商品价格之差是第k时段数量的线性减函数，系数为b。又已知当商品数量为500、价格为200时，处于平衡状态。建立差分方程模型描述商品数量和价格的变化规律，对以下情况作图讨论其变化趋势。（1）设a=0.2，b=0.1，开始时商品数量和价格分别在500和100附近。（2）利用特征根讨论变化趋势。【第(1)问模型建立及求解】记商品第k时段的上市数量为xk，价格为yk。因为第k+1时段与k时段商品价格之差是第k时段数量的线性减函数，则需求函数f可表示为yk+1-yk=-bxk+c又因为第k+1与k时段商品上市数量之差是第k时段价格的线性增函数，则供应函数h可表示为xk+1-xk=ayk+d令k，将平衡状态下x=500，y=200代入，可得c=500b，d=-200a。故需求函数和供应函数变为yk+1-yk=-b(xk-500)xk+1-xk=a(yk-200) （3）(1)以a=0.2，b=0.1代入，可得yk+1-yk=-0.1(xk-500)xk+1-xk=0.2(yk-200) （4）并取x0=500，y0=100，n=100，用MATLAB计算并作图，程序如下：%-作业题2_7函数M文件源程序-function x,y=Exf2_7(x0,y0,a,b,n) % 建立名为Exf2_7的函数M文件，x0,y0,a,b,n可调节x=x0; % 赋初值y=y0;for k=1:n y(k+1)=y(k)-b*(x(k)-500); x(k+1)=x(k)+a*(y(k)-200); % 按照（4）式进行迭代计算end%-作业题2_7脚本 M文件源程序-clear all x0=500;y0=100; % x0:商品初始数量； y0:商品初始价格a=0.2;b=0.1;n=100;k=(0:n);x,y=Exf2_7(x0,y0,a,b,n); % 给定x0,y0,a,b,n,调用Exf2_7函数计算round(k,x,y), % 对结果四舍五入取整plot(x,y), % 以商品数量为横坐标，商品价格为纵坐标作图xlabel(商品数量x)ylabel(商品价格y)title(商品数量与价格变化趋势（x0=500,y0=100,a=0.2;b=0.1）)% 加入X轴标记，Y轴标记和标题gtext(x0,y0)=(500,100) % 在图上标记初始点得到的结果如下面的表3和图7所示。表3：商品数量与价格变化趋势（x0=500,y0=100,a=0.2,b=0.1,n=100）kxkykkxkykykxkyk050010033696211674973941480100346981916853639424601023569617169574391344010636690152706133834422112376801337164937254041203866711572684357638812939650987371533973741414062983747433178362153416067075766293935316742580597678526610346182435525177798238113421974452246788062081234221345491447980817713344229464604580803147143502444742949817921161535925948399568277687163712734937066837536017386286503437984725341840329851319958569212194223075229811386654-7204443155328113387613-23214673215426715588568-34224913245525817889521-41235163255625420290473-43245413235725422791425-40255663195826025292377-33265893135927027693330-20276123046028529994286-328633293613053209524518296512796232934096209433066726463357357971777331680248643883719815210532690230654223829913314066459390100121176图7：商品数量与价格变化趋势（x0=500,y0=100,a=0.2;b=0.1,n=100）还可以分别单独作出商品数量x和商品价格y随时段k变化趋势的图形如下：图8：商品数量随时段变化趋势（x0=500,y0=100,a=0.2;b=0.1,n=100）图9：商品价格随时段变化趋势（x0=500,y0=100,a=0.2;b=0.1,n=100）【第(1)问结果分析及理论计算】由Matlab计算出来的数据表格和变化趋势图可见，在x0=500,y0=100,a=0.2，b=0.1的条件下，商品数量和价格的振荡的振幅越来越大，若没有外界如政府的干预等措施，经济将趋于崩溃。理论分析如下：由之前的分析，有yk+1-yk=-0.1xk-500xk+1-xk=0.2(yk-200)消去yk，可得xk+2-2xk+1+1.02xk=10其对应的特征方程为2-2+1.02=0其特征根为1，2=10.i在差分方程中令xk+2=xk+1=xk=x，可知平衡点为x=x0=500,故方程的解可表示为xk=c11k+c22k+x0, k=1,2,c1, c2由初始值x1,x2确定。同理可以求得，商品价格yk的解的形式也具有相同的形式。由解的形式可知，xkx0(k)(经济趋向稳定)的条件是1，21 ，故经济趋向于不稳定。可见，理论分析结果与Matlab计算结果是一致的。【第(2)问分析、求解】由第一问的分析，有yk+1-yk=-0.1xk-500xk+1-xk=0.2(yk-200)从上述两式中消去yk，可得xk+2-2xk+1+1+abxk=500ab （5）二阶差分方程（5）的的特征方程为2-2+1+ab=0 （6）算出其特征根为1，2=1abi （7）在差分方程中令xk+2=xk+1=xk=x，可知平衡点仍为x=x0=500,方程（5）的解可表示为xk=c11k+c22k+x0, k=1,2, (8) c1, c2由初始值x1,x2确定。同理可以求得，商品价格yk的解的形式也具有相同的形式。由（8）式可知，xkx0(k)(经济趋向稳定)的条件是1，20,b0,所以1，2恒大于一，不满足xkx0(k)(经济趋向稳定)的条件，故经济趋向于不稳定。必须采取外界措施如政府干预来防止经济崩溃。【本题总结、分析】对于题中所给的经济模型，无论a，b取何值（但必须a0,b0，即满足第k+1时段与k时段商品价格之差是第k时段数量的线性减函数，而第k+1与k时段商品上市数量之差是第k时段价格的线性增函数），与模型对应的二阶差分方程的两个特征根的模均大于1，故方程的解趋向于不稳定，商品数量和价格的振荡的振幅越来越大，若没有外界如政府的干预等措施，经济将趋于崩溃。题目8在某种环境下猫头鹰的主要食物来源是田鼠，设田鼠的年平均增长率为r1，猫头鹰的存在引起的田鼠增长率的减少与猫头鹰的数量成正比，比例系数为a1；猫头鹰的年平均减少率为r2；田鼠的存在引起的猫头鹰减少率的增加与田鼠的数量成正比，比例系数为a2。建立差分方程模型描述田鼠和猫头鹰共处时的数量变化规律，对以下情况作图给出50年的变化过程。（1）设r1=0.2，r2=0.3，a1=0.001，a2=0.002，开始时有100只田鼠和50只猫头鹰。（2）设r1，r2，a1，a2同上，开始时有100只田鼠和200只猫头鹰。（3）适当改变参数a1，a2（初始值同上）。（4）求差分方程的平衡点，它们是稳定吗？问题分析 猫头鹰与田鼠组成了食物链，它们的数量是相互制约的。由题意知，随着猫头鹰数量的增加，被猫头鹰捕获的田鼠数量逐渐增多，那么田鼠数量就逐渐减少了，又由于缺少了食物的来源，随之改变的是猫头鹰的数量减少。因此，可以根据题目条件建立差分方程来确定猫头鹰与田鼠的数量关系，进而求出50年田鼠与猫头鹰共处时的数量变化规律。模型假设 根据题意，设第k年田鼠的数量为xk，猫头鹰的数量为yk。田鼠的年平均增长率为r1，猫头鹰的年平均减少率为r2。鼠增长率的减少与猫头鹰的数量成正比，比例系数为a1，猫头鹰减少率的增加与田鼠的数量成正比，比例系数为a2。随时可以改变a1，a2来观察猫头鹰和田鼠的变化规律。本题中为了与实际相结合，任何一种生物繁殖都种群需要满足最小的数量，这一数量即为阈值，故认为猫头鹰或田鼠的数量为个位数是，即可认为该种群灭绝。模型建立 建立猫头鹰与田鼠数量关系的差分方程，田鼠的数量变化方程为，猫头鹰的数量变化方程为。联立两个差分方程来求解每年猫头鹰和田鼠的数量。这样通过这两个相互制约的递推关系，由MATLAB计算可以得到田鼠和猫头鹰共处时50年的数量变化，并作出图形。算法设计 由于本题的相关参数有多次变化，故可以编写函数文件和脚本文件，通过改变脚本文件中的不同参数就能做出不同的图形。函数文件中主要关于差分方程的迭代运算，特别是要求x和y两个行向量（代表田鼠和猫头鹰的数量）要大于零且为整数。程序 用MATLAB编程计算，源程序如下所示：%以下为计算田鼠和猫头鹰50年数量的函数文件function x y=quantity(r1,r2,a1,a2,a,b)x(1)=a;y(1)=b;for k=1:49 x(k+1)=round(x(k)+(r1-a1*y(k)*x(k); y(k+1)=round(y(k)+(a2*x(k)-r2)*y(k); if x(k+1)0 x(k+1)=0; end if y(k+1)0 y(k+1)=0; endend%脚本文件，作图观察田鼠和猫头鹰的变化趋势k=1:50;x1 y1=quantity(0.2,0.3,0.001,0.002,100,50);%第一小题x2 y2=quantity(0.2,0.3,0.001,0.002,100,200);%第二小题x3 y3=quantity(0.2,0.3,0.1,0.1,100,200);%以下为改变a1,a2x4 y4=quantity(0.2,0.3,0.0001,0.0001,100,200);x5 y5=quantity(0.2,0.3,0.001,0.003,100,200);x6 y6=quantity(0.2,0.3,0.1,0.0001,100,200);subplot(2,1,1),plot(k,x1,r-*,k,y1,b-),gridxlabel(年份),ylabel(数量)title(田鼠和猫头鹰50年的变化规律（田鼠100只，猫头鹰50只）)legend(田鼠的数量,猫头鹰的数量)subplot(2,1,2),plot(k,x2,r-*,k,y2,b-),gridxlabel(年份),ylabel(数量)title(田鼠和猫头鹰50年的变化规律（田鼠100只，猫头鹰200只）)legend(田鼠的数量,猫头鹰的数量)pausesubplot(2,1,1),plot(k,x3,r-*,k,y3,b-),gridxlabel(年份),ylabel(数量)title(田鼠和猫头鹰50年的变化规律（a1=0.1,a2=0.1）)legend(田鼠的数量,猫头鹰的数量)subplot(2,1,2),plot(k,x4,r-*,k,y4,b-),gridxlabel(年份),ylabel(数量)title(田鼠和猫头鹰50年的变化规律（a1=0.0001,a2=0.0001）)legend(田鼠的数量,猫头鹰的数量)pausesubplot(2,1,1),plot(k,x5,r-*,k,y5,b-),gridxlabel(年份),ylabel(数量),axis(0 50 0 300)title(田鼠和猫头鹰50年的变化规律（a1=0.001,a2=0.003）)legend(田鼠的数量,猫头鹰的数量)subplot(2,1,2),plot(k,x6,r-*,k,y6,b-),gridxlabel(年份),ylabel(数量)title(田鼠和猫头鹰50年的变化规律（a1=0.1,a2=0.0001）)legend(田鼠的数量,猫头鹰的数量)计算结果 （1）当r1=0.2，r2=0.3，a1=0.001，a2=0.002，开始有100只田鼠和50只猫头鹰。（2）当r1=0.2，r2=0.3，a1=0.001，a2=0.002，开始有100只田鼠和200只猫头鹰。它们的50年来的数量关系如下所示：对于问题（1），田鼠和猫头鹰的数量开始都相继，后来当猫头鹰的数量很高时田鼠的数量就趋近于零。这样，猫头鹰和田鼠就会灭绝。对于问题（2），由于猫头鹰初始数量的不同，老鼠和猫头鹰的数量变化比问题（1）中要平缓，但是在50末附近也会出现猫头鹰的数量激增至800多只，这样就造成了老鼠数量的急剧减少，逐渐灭绝，猫头鹰无食物可吃也会灭绝。（3）改变a1=0.1，a2=0.1，开始有100只田鼠和200只猫头鹰；改变a1=0.0001，a2=0.0001，开始有100只田鼠和200只猫头鹰。它们的50年来的数量关系如下所示：若先改变参数a1=0.1，a2=0.1，同时增大了a1与a2，其实际意义是田鼠和猫头鹰之间的实际制约关系加强了，但在实际中是不符合现实状况的。在这样的情况下，由于猫头鹰的数量一开始是田鼠的两倍，所以猫头鹰的数量大幅度增长，老鼠的数量减小的很快，进而很快就灭绝了。猫头鹰也相继灭绝。再观察下同时减小a1与a2的情况，令a1=0.0001，a2=0.0001，则与上面的那种情况相反，田鼠和猫头鹰之间的实际制约关系减弱了，在两个物种数量不大的情况下从图上可以看出开始时田鼠的数量增多，猫头鹰的数量减小，由于猫头鹰对田鼠的数量影响不是很大，田鼠的数量会出现大幅度的增长。但是在此又出现了不符合实际的情况，因为田鼠数量的剧增，本来数量只剩下个位数的猫头鹰竟然不合理的增加。就如上面提到过的阈值，猫头鹰的数量太少了，它们是没办法进行大量繁殖的，所以从图上看来，只有25年以前的情况是合理的，后面25年的情况因为忽略了生物繁殖方面的因素，只单纯地考虑数值上的变化，这样就产生了与实际不相符的错误估计。改变a1=0.001，a2=0.003，开始有100只田鼠和200只猫头鹰；改变a1=0.1，a2=0.0001，开始有100只田鼠和200只猫头鹰。它们的50年来的数量关系如下所示：再改变参数a1=0.001，a2=0.003，现在这两个参数的取值很精妙，它们的取值是根据田鼠和猫头鹰的初始数量来决定的。根据两个反映田鼠和猫头鹰数量变化的差分方程，可以知道x=100，y=200是其中的一个平衡点。故从图表上反映出的现象是猫头鹰和田鼠的数量一直分别维持在200只和100只。在这样的情况下，虽然从差分方程上认为它们之间是相互影响的，但是二者的初值却决定了它们数量一直不变。自然界中有很多不可预计的因素，显然这样的情况在自然界中也是不可能的。最后调整参数a1=0.1，a2=0.0001，这样取值的意义是加强猫头鹰对田鼠的影响并减缓猫头鹰减少率的变化。从图上可以看出猫头鹰和田鼠的数量都在减少，直至为零，两个种群都相继灭绝。结果分析 根据题目所给参数和两个差分方程，可以令，这样求出差分方程的平衡点为（0，0）和（150，200）。对于这种非线性差分方程组，可以通过作图和一些数据来进行讨论。对于（0，0），虽然在很多情况下可以到达平衡点（0，0）,但是若猫头鹰先于田鼠灭绝，田鼠就能大量繁殖且不会灭绝，则易知（0，0）不是稳定平衡点。另外（150，200）只有在初值为x=150，y=200时才能达到平衡，则（150，200）也不是稳定平衡点。实际意义 该题在数学建模过程中首先作了一些假定，归结起来主要有以下几点：1.只有一条简单的食物链田鼠猫头鹰；2.田鼠和猫头鹰的数量只是互相影响，不与其他任何自然因素有关；3.田鼠和猫头鹰的繁殖中有一个最小的阈值，若它们其中某一种群的数量太少，就可以认为是灭绝。根据题意，通过以上几点假定作图分析的结果不太符合实际，条件若能再限制多些，比如考虑进大自然的影响和增加食物链，就能使得该数学建模更加符合实际的结果。结论 虽然该题条件不太符合实际，但是有一定的启示。若有两个相互影响的物种，由于初始数量的不同，这也是能够影响它们数量未来的发展趋势的；另外，若一种群的天敌灭亡，则该种群是能够得到大量繁殖的。通过以上MATLAB的作图分析即可看出以上结论。题目9研究将鹿群