数学实验报告——利用MALTAB计算常微分方程数值解.docx

实验二 常微分方程数值解一、火箭飞行器问题描述小型火箭初始质量为1400kg，其中包括1080kg燃料，火箭竖直向上发射时燃料燃烧率为18kg/s，由此产生32000N的推力，火箭引擎在燃料用尽时关闭，设火箭上升时空气阻力正比于速度的平方，比例系数为0.4kg/m，求引擎关闭瞬间火箭的高度，速度，加速度及火箭达到最高点时的高度和加速度，并画出高度、速度、加速度随时间变化的图形。方法与公式1、简要分析本题的求解需要用到常微分方程，而整个过程又被分为两个阶段：火箭加速上升阶段和燃料燃尽后减速的阶段。由题目易知第一个阶段持续时间T1=108018=60s列出第一阶段的方程组：设M0为火箭本身质量，m为燃料质量，g为重力加速度 = 9.8m/s2,燃料燃烧率为a，空气阻力的比例系数为k，F为推进力。M0 = 1400-1080 = 320kg; v=F-kv2-M0+mgM0+m m=-a初值v=0，m=1080。由以上各式可以求出t=T1时火箭的速度。再求解第二阶段： v=-kv2-M0gM0 m=0可以求出火箭速度降为0的时刻。将整个过程中的时间向量以及速度向量联合起来，利用第三章所学插值与数值积分的方法可以求得任意时刻火箭的近似高度。2、方法求解常微分方程时，我分别采用了自己编写的欧拉公式、改进欧拉公式、4级4阶龙格-库塔公式，以及MATLAB自带的龙格-库塔方法，求解数值积分时采用辛普森公式。由于Matlab自带的Simpson公式是自适应的，因此需要使用自己在上一次实验时所编的Simpson公式。结果与分析1、各种公式的对比首先，我作出了各种不同公式计算得到的火箭速度随时间变化的图像，图如下：从图中可以看出，各种公式计算得到的结果基本一致，为确定其区别，将图像放大，放大约2000 倍后，得到下图：分析：从图中可以看出，自编欧拉公式距离MATLAB自带龙格-库塔公式最远，精度最差；自编的改进欧拉公式和自编的龙格-库塔公式结果基本一致，两者中自编龙格-库塔公式距MATLAB自带龙格-库塔公式的结果稍近。与之前的分析基本一致。然而产生自编龙格-库塔公式与MATLAB自带龙格-库塔公式之间的差距的原因还未知。由于MATLAB自带龙格-库塔公式精度较高，因此以下各项实验不再计算其它几项公式的结果。2、第一阶段火箭关闭瞬间的速度v=267.261240773261 m/s关闭前瞬间的加速度a=F-M0g-kv2M0=0.914286475421320 m/s2此时火箭的高度h=12189.7791507305 m3、第二阶段由初始速度以及常微分方程可以求得火箭达到最高点的时间约为71.3s；此时火箭的高度h=13115.7148739887 m加速度a=-9.80000406052111m/s24、整体过程下图为火箭加速度与时间的关系。由图可以看出，火箭一开始的加速度很大，随着时间的推移，火箭的燃料有所减少，与此同时速度有所上升。两者中前者使火箭的加速度增大，后者使其减小，综合作用，最终体现为加速度先略有上升，然后慢慢减小。当火箭中燃料燃尽时，火箭丧失推动力，因而加速度急剧减小为负值。此后火箭速度不断减小，导致火箭所受阻力逐渐减小，因而加速度的绝对值有所减小，直到最终火箭速度降为零，火箭不受阻力，仅受重力，加速度为重力加速度。 下图为火箭速度与时间的关系： 此图可以看作是由第一幅图对时间积分所得结果，本图的斜率对应第一幅图的值。第一阶段火箭加速度为正，因此速度不断增加，只是增加的速度不断减慢。燃料燃尽后，加速度变为负值，因此速度开始急剧下降，与此同时下降的速率不断减小。最终火箭速度降为0。 下图为火箭高度与实践的关系：此图可以看作是由第二幅图对时间积分所得结果，本图的斜率对应第一幅图的值。一开始或减速度较小，因此高度缓慢增加，之后增加的速度不断提升，直到火箭的燃料耗尽，此后速度不断减小，但仍为正值，因此火箭继续向上飞行，只是高度提升的速度逐渐变慢。知道最终火箭速度降为零。3、总结这是自己接触数学实验以来第一次解决相对复杂的实际问题，总体上还是有些不适应。虽然最终得以成功解决，但是编程的逻辑以及脚本之间的层次没有考虑清楚，需要在之后的学习中不断加强相关能力。当然也有不小收获：(1)了解了MATLAB解决常微分方程组的原理与步骤；(2)进一步熟悉了MATLAB的相关操作，包括画图、调试、编写脚本等；(3)加强了对之前插值与数值积分的理解，学习将其应用到实际问题的解决中。程序清单1、第一阶段的微分方程组function dx = fly(t,y)M0 = 320;a = 18; k = 0.4;g = 9.8;F = 32000;dx = (-k*y(1)2 + F - (M0 + y(2)*g)/(M0+y(2) ; - a*sign(y(2);2、第二阶段微分方程组function dx = fly1(t,y)M0 = 320; k = 0.4;g = 9.8;F = 0;a = 0;dx = (-k*y(1)2 + F - (M0 + y(2)*g)/(M0+y(2) ; -a;3、自编欧拉公式function y = Euler(fun,ts,x0)NumOfLine = numel(ts);NumOfRow = numel(x0);h = ts(2)-ts(1);y = zeros(NumOfLine,NumOfRow);y(1,:) = x0;for i = 1:NumOfLine-1 dy = feval(fun,ts(i),y(i,:); y(i+1,:)=y(i,:)+rot90(dy)*h;end4、自编改进欧拉公式function y = ImprovedEuler(fun,ts,x0)NumOfLine = numel(ts);NumOfRow = numel(x0);h = ts(2)-ts(1);y = zeros(NumOfLine,NumOfRow);y(1,:) = x0;for i = 1:NumOfLine-1 dy = feval(fun,ts(i),y(i,:); y1 = y(i,:)+rot90(dy)*h; dy1 = feval(fun,ts(i)+h,y1); y(i+1,:)=y(i,:)+h/2*(rot90(dy+dy1);end5、自编龙格-库塔公式function y = myRungeKutta(fun,ts,x0)NumOfLine = numel(ts);NumOfRow = numel(x0);h = ts(2)-ts(1);y = zeros(NumOfLine,NumOfRow);y(1,:) = x0;for i = 1:NumOfLine-1 k1 = feval(fun,ts(i),y(i,:); dy1 = y(i,:)+h/2*rot90(k1); k2 = feval(fun,ts(i)+h/2,dy1); dy2 = y(i,:)+h/2*rot90(k2); k3 = feval(fun,ts(i)+h/2,dy2); dy3 = y(i,:)+h*rot90(k3); k4 = feval(fun,ts(i)+h,dy3); y(i+1,:)=y(i,:)+h/6*rot90(k1+2*k2+2*k3+k4);end6、“自启动”脚本1e =0.1;ts1 = 0 : e : 60;y00 = 0,1080;y01 = Euler(fly,ts1,y00);y11 = ImprovedEuler(fly,ts1,y00);y21 = myRungeKutta(fly,ts1,y00);t1,y1 = ode45(fly,ts1,y00);ts2 = 60 : e :80;y00 = y1(60/e+1),1),0;y000 = y01(60/e+1),1),0;y001 = y11(60/e+1),1),0;y002 = y21(60/e+1),1),0;y02 = Euler(fly1,ts2,y000);y12 = ImprovedEuler(fly1,ts1,y001);y22 = myRungeKutta(fly1,ts1,y002);t2,y2 = ode45(fly1,ts2,y00);ts = ts1,ts2; t = t1;t2;i=1;while(y2(i,1)0) i = i+1;endy = y1;y2(1:i-1,:);y0 = y01;y02(1:i-1,:);y1 = y11;y12(1:i-1,:);y2 = y21;y22(1:i-1,:);j = numel(y)/numel(y00);plot(t(1:j),y(:,1),k,ts(1:j),y0(:,1),b,ts(1:j),y1(:,1),y,ts(1:j),y2(:,1),r);grid on;legend(MATLAB自带龙格-库塔公式,自编的欧拉公式,自编的改进欧拉公式,自编的龙格-库塔公式);7、自启动脚本2e =0.01;ts1 = 0 : e : 60;y0 = 0,1080;t1,y1 = ode45(fly,ts1,y0);ts2 = 60 : e :80;y0 = y1(60/e+1),1),0;t2,y2 = ode45(fly1,ts2,y0);ts = ts1,ts2;t = t1;t2;i=1; while(y2(i,1)0) i = i+1;endy = y1;y2(1:i-1,:);y(1,2) = (32000-(320+1080)*9.8)/(320+1080);for j = 2:60/e+i h(j) = mySimpson(t(1:j),y(1:j,1); if(j60/e+1) y(j,2) = (-0.4*y(j,1)2+32000-(320+y(j,2)*9.8)/(320+y(j,2); else y(j,2) = (-0.4*y(j,1)2-320*9.8)/320; endendh = h(1:j);t = t(1:j);plot(t,h,k);grid on;二、小船过河问题描述一只小船渡过宽度为d的河流，目标是起点A正对着的另一岸B点。已知水流速v1与船在静水中的速度v2之比为k。建立小船航线的数学模型，求解析解以及数值解。方法与公式1、数学模型与解析解由题图，以B为原点，BA为y轴负半轴建立直角坐标系，设任意时刻小船坐标为（x,y），则易得下列公式 x=v1+-x*v2x2+y2y=-y*v2x2+y2此微分方程组可用于作数值计算。为求解析解，采用极坐标系分析，设小船到原点的距离为r，与x轴正半轴的夹角为 r=-v2+v1cosr=-v1sin()两式相除，分离变量，积分，得lnr=v2v1lntan2-lnsin-C当r=d时，=90，代入得r=dsin(tan2)v2v1此式即为解析解，经过坐标变换可以得到x,y的关系。2、数值解 利用上一小节得到的微分方程组，通过龙革-库塔方法即可得到方程组的数值解，此外还可以与解析解所得数值进行比较。3、各种其他情况当流速为0,0.5,1.5,2时，可以现根据解析解进行判断，然后再进行数值计算。结果与分析1、d = 100m, v1=1m/s, v2=2m/s(1)使用(x(1)+0.00001计算函数M文件如下：function dz = Cross(t,x)v1=1;v2=2;r = sqrt(x(1)+0.00001)2+x(2)2);dz = v1-x(1)*v2/r;-x(2)*v2/r; 其中写为x(1)+0.00001是为了确定更充分地测试小船到达的时间，因为当(x,y)趋近于(0,0)时，Matlab的Ode45函数会因为r趋于0而无法顺利算出结果。为了说明这样做并没有带来多少误差，下面会用x(1)+0.00002再次计算，观察两者结果的差别。主方法如下：ts = 0:0.1:70;x0 = 0,-100;t,x = ode45(Cross,ts,x0);从测试结果得到，t=66.8s时，x= 9.91158339584917e-06，y= 3.10706640759385e-07，近似于原点。故可认为所需时间为66.8s。任意时刻位置如下： t（s） x（m） y（m） 0 0 -1.000000000000000e+002 1.000000000000000e+000 9.899329169404653e-001 -9.800003384356300e+001 2.000000000000000e+000 1.959459981808873e+000 -9.600027487106881e+001 3.000000000000000e+000 2.908165884197600e+000 -9.400094195010428e+001 4.000000000000000e+000 3.835625011732558e+000 -9.200226757752603e+001 5.000000000000000e+000 4.741402426855176e+000 -9.000449892891601e+001 6.000000000000000e+000 5.625050746668156e+000 -8.800789678100108e+001 7.000000000000000e+000 6.486110777329467e+000 -8.601273864220910e+001 8.000000000000000e+000 7.324113001651424e+000 -8.401932041779020e+001 9.000000000000000e+000 8.138577579100687e+000 -8.202795640981674e+001 1.000000000000000e+001 8.929014345798258e+000 -8.003897931718340e+001 1.100000000000000e+001 9.694922823250010e+000 -7.805274025331977e+001 1.200000000000000e+001 1.043579221204906e+001 -7.606961290713413e+001 1.300000000000000e+001 1.115109539532280e+001 -7.408999387621988e+001 1.400000000000000e+001 1.184029353402413e+001 -7.211430428718745e+001 1.500000000000000e+001 1.250283877098684e+001 -7.014299264573725e+001 1.600000000000000e+001 1.313817423092567e+001 -6.817653483665951e+001 1.700000000000000e+001 1.374573402043624e+001 -6.621543412383441e+001 1.800000000000000e+001 1.432494325546332e+001 -6.426022116961551e+001 1.900000000000000e+001 1.487522003112953e+001 -6.231145960144670e+001 2.000000000000000e+001 1.539596223181069e+001 -6.036974973213587e+001 2.100000000000000e+001 1.588655786949937e+001 -5.843572897362697e+001 2.200000000000000e+001 1.634639227881750e+001 -5.651007437986918e+001 2.300000000000000e+001 1.677484811701635e+001 -5.459350264681689e+001 2.400000000000000e+001 1.717130536397652e+001 -5.268677011242966e+001 2.500000000000000e+001 1.753514139504863e+001 -5.079067275926418e+001 2.600000000000000e+001 1.786573788019881e+001 -4.890605168463580e+001 2.700000000000000e+001 1.816245032926312e+001 -4.703380942270495e+001 2.800000000000000e+001 1.842463162884271e+001 -4.517490393331951e+001 2.900000000000000e+001 1.865164996561456e+001 -4.333034628944002e+001 3.000000000000000e+001 1.884288882633143e+001 -4.150120067713974e+001 3.100000000000000e+001 1.899774699782184e+001 -3.968858439560458e+001 3.200000000000000e+001 1.911563875342150e+001 -3.789366775974959e+001 3.300000000000000e+001 1.919601251292106e+001 -3.611767068597514e+001 3.400000000000000e+001 1.923827436902908e+001 -3.436192717024910e+001 3.500000000000000e+001 1.924184405457128e+001 -3.262785320447557e+001 3.600000000000000e+001 1.920619988005019e+001 -3.091692040884441e+001 3.700000000000000e+001 1.913087873364515e+001 -2.923065603183126e+001 3.800000000000000e+001 1.901547608121231e+001 -2.757064295019759e+001 3.900000000000000e+001 1.885964648940898e+001 -2.593851906196922e+001 4.000000000000000e+001 1.866315381901801e+001 -2.433591835407299e+001 4.100000000000000e+001 1.842564378334017e+001 -2.276475253244378e+001 4.200000000000000e+001 1.814679487434363e+001 -2.122706856439528e+001 4.300000000000000e+001 1.782644456008798e+001 -1.972490660448889e+001 4.400000000000000e+001 1.746458928472416e+001 -1.826029999453374e+001 4.500000000000000e+001 1.706138446849452e+001 -1.683527526358664e+001 4.600000000000000e+001 1.661714657296644e+001 -1.545184790711154e+001 4.700000000000000e+001 1.613254393856277e+001 -1.411153523845545e+001 4.800000000000000e+001 1.560764738725192e+001 -1.281711081948180e+001 4.900000000000000e+001 1.504249262731523e+001 -1.157180092254236e+001 5.000000000000000e+001 1.443757721814349e+001 -1.037835103360374e+001 5.100000000000000e+001 1.379386057023698e+001 -9.239025852247460e+000 5.200000000000000e+001 1.311276394520541e+001 -8.155609291669929e+000 5.300000000000000e+001 1.239616788478941e+001 -7.129406425249822e+000 5.400000000000000e+001 1.164550647085605e+001 -6.162400305079939e+000 5.500000000000000e+001 1.086152883779034e+001 -5.258361410743773e+000 5.600000000000000e+001 1.004575408440289e+001 -4.420111721734861e+000 5.700000000000000e+001 9.200389899411650e+000 -3.649193151059897e+000 5.800000000000000e+001 8.328011629192378e+000 -2.946397463538269e+000 5.900000000000000e+001 7.430349226982941e+000 -2.315074927076017e+000 6.000000000000000e+001 6.509732292564531e+000 -1.757524897325378e+000 6.100000000000000e+001 5.569261527276989e+000 -1.273637564234637e+000 6.200000000000000e+001 4.611316377221014e+000 -8.662361322252641e-001 6.300000000000000e+001 3.638932779284402e+000 -5.357652477978345e-001 6.400000000000000e+001 2.654994008077677e+000 -2.838279528432279e-001 6.500000000000000e+001 1.662619419689339e+000 -1.109391788034044e-001 6.600000000000000e+001 6.649434118513309e-001 -1.772322378297511e-002 6.700000000000000e+001 9.968875530977945e-006 -6.404494258255037e-008 6.800000000000000e+001 9.955642978032652e-006 -2.335004420014544e-007 6.900000000000000e+001 9.914062890159047e-006 3.157597756744254e-0077.000000000000000e+001 9.992844098020244e-006 -5.478066431501524e-007航行曲线如下：(2) 用x(1)+0.00002再次计算得出t=66.8s时，x=1.99973197904343e-05，y=-2.77389619481921e-07，趋近于原点，于是可以看出在公式中增加微小的一项并不会对本题中的目标造成多少误差。(3)与解析解进行对比由于r、与t的关系非常复杂，因此这里不进行数值上的对比，仅将两条曲线画在同一幅图中加以比较。增加代码：fi = 0 :pi/2/100:pi/2;r = (100./sin(fi).*(tan(fi/2).2);x1 = r.*cos(fi);y1 = -r.*sin(fi);图如下：由图可以看出，数值解与解析解极其相似，两者几乎化为一条曲线。可见数值计算的精度很高。2、其它情况(1)v1=0此时相当于船在静水中前进，所以航线为一条直线，到达对岸所用时间为50s。(2)v1=0.5此时水流速度较小，因此船航行曲线不像原图那样弯曲。计算得t=53.5s时，x=3.27480360132058e-06，y=-3.39925703953842e-07，趋近于原点。（注：横坐标已改为0,12）(3)v1=1.5m/s 此时水流速度较大，小船航线较以往弯曲。t=114.4s时，x=2.99731260021995e-05，y=-3.84983039618783e-07，趋近于原点。(4)v1=2m/s此时水流速度相当大，可以预测小船无法顺利到达目的地。由实验结果可以看出小船将近似停留在(50,0)的位置，这与解析解得出的结果相吻合趋于/2时，r=50。航线如下图：3、总结a.小船航线的弯曲程度、到达对岸的时间与水流速度成正相关，也就是说，水流速度越大，航线越弯曲、到达对岸的时间越长；b.当水流速度达到2m/s时，小船将永远无法到达目的地，最终停留在(50,0)的位置；c.采用龙革-库塔进行的数值计算的结果与解析解的差别很小，基本无法看出两者的差别。程序清单1、微分方程组function dz = Cross(t,x)v1=2;v2=2;r = sqrt(x(1)+0.00001)2+x(2)2);dz = v1-x(1)*v2/r;-x(2)*v2/r;2、“自启动”脚本ts = 0:0.1:200;x0 = 0,-100;t,x = ode45(Cross,ts,x0); fi = 0 :pi/2/1000:pi/2;r = (100./sin(fi).*(tan(fi/2);x1 = r.*cos(fi);y1 = -r.*sin(fi); plot(x(:,1),x(:,2),b,x1,y1,r);axis(0,60,-100,20);grid on;xlabel(x);ylabel(y);legend(,);三、种群竞争问题描述给出两个种群相互竞争的模型： xt=r1x1-xn1-s1yn2yt=r2y1-s2xn1-yn2(1)给定一组初值以及各个参数值，说明t充分大以后x(t),y(t)的变化趋势；(2)改变各个参数的值，保证s1、s2不变，计算并分析所得结果；再改变s1、s2；(3)试验当s1=0.8、s2=0.7时会有什么结果；s1=1.5、s2=1.7时会有什么结果；方法与公式 使用龙革-库塔方法，对所给微分方程组进行数值分析。结果与分析1、r1=r2=1，n1=n2=100，s1=0.5，s2=2，x0=y0=10.函数：function dz = Compete(t,x)r1=1;r2=1;n1=100;n2=100;s1=0.5;s2=2;dz = r1*x(1)*(1-x(1)/n1-s1*x(2)/n2);r2*x(2)*(1-s2*x(1)/n1-x(2)/n2);得到图形：下图为x,y随t的变化图下图为y随x的变化图由图可以看出，时间足够长后，x与y数量悬殊变大，最终一方繁盛，一方灭绝。2、改变系数(1)保持s1=0.5，s2=2.改变r1、r2a.同时增大，r1=r2=2由图可以看出，x、y构成的系统比原来更快地趋于稳定。b.同时减小，r1=r2=0.5由图可以看出，x、y构成的系统趋于稳定的速度变低了。c. x增大，y减小，r1=2，r2=0.5由图可以看出，x、y构成的系统趋于稳定的速度与r1、r2同时增大效果相近。d. y增大，x减小，r1=0.5，r2=2由图可以看出，x、y构成的系统趋于稳定的速度变低了，而且比r1、r2同时减小的情况趋于稳定的速度还低，在一开始时出现了y大于x的时间段。e.分析 r1、r2很大程度上决定了系统变化的速率物种的增长率越高，生长得越快，系统演变的速度就越快。改变n1、n2a.同时增大，n1=n2=200由图可以看出，x、y构成的系统趋于稳定的速度基本没有改变，但是稳定时的状态改变了x的数量由100变为了200；与此同时，y的最大值比以前有所提升。b.同时变小，n1=n2=50x、y构成的系统趋于稳定的速度仍旧基本没有改变，但是稳定时的状态改变了x的数量由100变为了50；与此同时，y的最大值也比之前有所降低。c. x增大，y减小 n1=200，n2=50x、y构成的系统稳定时的状态与n1、n2同时增大时大体相同。d. y增大，x减小 n1=50，n2=200x、y构成的系统稳定时的状态与n1、n2同时减小时大体相同。e.分析 n1、n2很大程度上决定了系统的最终状态，相当于是系统的“特解”。当然它并不确定解的形式。改变x0、y0a.同时增大x0=y0=50. 系统趋于稳定的速度变化不大，稳定时的状态也没有什么变化，不过这次y没有经过波动便直接趋向灭绝，原因应为：开始时x的数量较大，使得供养乙的资源很快被耗尽，导致y很快趋于灭绝。b.同时减小，x0=y0=5. 系统趋于稳定的速度变化不大，稳定时的状态也没有什么变化，乙的数量出现了波动，且比原本的状态波动大。原因应为：x、y的初始数量都较小，因此供养乙的资源没有被很快地耗尽，因此y出现了一个生长的高峰。c. x增大，y减小，x0=50，y0=5.效果和原因与，x0、y0同时增大时相似。d. x减小，y增大，x0=5，y0=50.效果和原因与，x0、y0同时增大时相近。e.分析初始数量决定了系统刚开始一段时间大致状态，但是对系统最终的影响几乎为零。从系统的角度上讲，它并不影响系统的特解。(2)改变s1、s2，使s1=0.8，s2=0.7其余不变 由图可以看出，x、y的地位发生了翻转，即x很快趋于灭绝，y则长盛不衰。在此基础上改变r1、r2、n1、n2、x0、y0实验所得曲线与之前所得曲线相近，但x、y所处地位翻转。分析 s的含义是：对于供养一方的资源而言，单位数量另一物种的消耗为单位数量本物种消耗的s倍。3、改变s1、s2(1)当s1=0.8，s2=0.7时其余不变 可以看出，x、y最终达到了平衡，但是平衡时两者都没有达到各自的最大容量。原因应为：两者对于对方的供养资源的占用率都没有彼此对自己资源的占用率高，这也就保证了每个品种都有自己赖以生存的一部分资源，从而也得以维持恒定。在此基础上改变r1、r2、n1、n2、x0、y0各项对于系统的作用仍旧如之前的分析：r对应系统变化的速率、n对应系统最终稳定时各物种的绝对数量、初始物种数量对应系统在开始一段时间的状态。分析这一次没有出现之前一直出现的一方昌盛一方灭绝的情况，而是两个物种并存，数量有所差别，这是由于每个物种对自己赖以生存的资源都没有比另一方更有竞争力，这就保证了每个物种都不会因为资源短缺而趋于灭绝。(2)当s1=1.5，s2=1.7时其余不变 可以看出，系统再次出现了一方灭亡一方昌盛的局面，但是此时系统是否和初始时的状态相同呢？改变na. n1=100, n2=200而从实验增大n2的结果可以看出，这一项并没有引起什么变化。b. n1=100, n2=20从图中可以看出，x、y的地位颠倒了。这可以这样解释：尽管初始时两者数量相当，但由于n2较小，也就是说y所占的比重较大，这在很大程度上影响了s的作用，即，n2的减小增大了s2的作用，这就使y占领了优势地位，从而获得竞争的“胜利”。改变ra.r1=1， r2=2由于原本局势对y不利，因此可以预见增大r1将仅使系统更快地达到之前的效果；而从实验增大r2的结果可以看出，x与y的地位再次颠倒。这种情况可以这样解释：由于s大于1，所以x、y都有置对方于灭绝之境的“潜力”，此时，拥有更高的增长率的一方将在角逐中获胜，占领竞争的优势地位，最终“获胜”。b .r1=1， r2=1.2由图可以看出，y并没有因为拥有较高的增长率而最终占据优势地位，相反却在短暂领先后趋于灭绝。可见s仍旧在很大程度上决定了双方竞争的优劣地位。改变初值