MATLAB数学建模2乒乓球的弹跳和罗基斯帝模型

1乒乓球的弹跳罗基斯第模型[问题]罗基斯第模型一个乒乓球离球拍的高度为h0，落在球拍上反弹，设恢复系数为e，不计空气阻力。(1)如果e为常数，讨论球的高度变化的规律。如果e2与高度hn成线性关系e2=μ(1–hn/H0)(2.1)其中H0是最大高度，μ是参数。对于不同的参数讨论小球高度的变化规律。(2)当参数连续变化时，分析最后分布的高度。(3)计算前几个分岔点。(4)用李雅普洛夫指数判断混沌的发生。[解析](1)当球从高度hn下落到球拍上之前速度为2nnvgh(2.2)球与球拍碰撞后反弹的速度为v'n=evn(2.3)球反弹的高度为hn+1=e2hn(2.4)如果e1，则球的反弹高度随次数不断减小；如果e=1，则球反弹后始终保持初始高度；如果e1，例如球拍每次加一个向上的冲击力，则球的高度随次数不断增加。e2与高度的线性关系说明：如果球的高度较大，则恢复系数较小，反之较大。设相对高度为xn=hn/H0，则下一次上升的相对高度为xn+1=μ(1–xn)xn，(n=0,1,2,…)(2.5)这是著名的罗基斯第模型。由于相对高度0≤xn≤1，而(1–xn)xn的最大值为1/4，所以参数的值在0到4之间。球的高度强烈依赖参数。[算法](1)先取一个参数，再取一个相对高度，通过迭代算法计算下一次碰撞后的高度，画出高度点，依此类推。再取另一高度参数，重新通过迭代算法计算高度，画出高度点，依此类推。[程序]MATH2_1.m如下。%乒乓球与球拍的碰撞高度clear%清除变量u=input('请输参数(参考值:0.5,2,3.25,3.5,3.56,3.8):');%键盘输入初始相对高度(1)xn=0.9;%第1个的初始相对高度(2)figure%开创图形窗口plot(0,xn,'.')%画高度点text(0,xn,num2str(xn),'FontSize',16)%标记第1个的初始高度gridminor%加细网格title(['乒乓球与球拍的碰撞高度(\it\mu\rm=',num2str(u),')'],'FontSize',16)%标题n=50;%迭代次数axis([0,n,0,1])%坐标范围holdon%保持图像forj=1:n%按次数循环xn=u*(1-xn)*xn;%计算下一次的相对高度(3)plot(j,xn,'.')%画高度点2end%结束循环xn=0.1;%取初始相对高度(4)plot(0,xn,'ro')%画高度点text(0,xn,num2str(xn),'FontSize',16)%初始高度forj=1:n%按次数循环xn=u*(1-xn)*xn;%计算下一次的相对高度(5)plot(j,xn,'ro')%画高度点end%结束循环[说明](1)程序执行时要用户用键盘输入参数，提供6个参数选择。(2)取第1个较大的初始高度。(3)迭代计算下一个高度。(4)取第2个较小的初始高度。在说明混沌时，将此句改写如下，使第2个高度比第1个高度大一点。e=1e-8;%小量xn=0.9+e;%取初始相对高度text(0,0,num2str(e),'FontSize',16)%初始高度(5)同样迭代计算下一个高度。M2.1a图M2.1b图[图示](1)如M2.1a图所示，当参数μ为0.5时，如果初始相对高度取0.9，球与球拍碰撞之后高度不断降低，最终的高度为零。即使初始相对高度取0.2，球与球拍碰撞之后高度也不断降低，最终的高度为零。M2.1c图M2.1d图(2)如M2.1b图所示，当μ为2时，如果初始相对高度取0.9，球第一碰撞之后高度降低，以后碰撞则高度升高，最后碰撞保持一定的高度。如果初始相对高度取0.2，则碰撞高度不3断增加，最后稳定在一定的高度。高度稳定前的过程称为过渡过程或暂态过程，过渡过程与初始高度有关，但是最后高度的稳定与初始高度无关。高度值称为不动点，即重复自身轨迹的点。(3)如M2.1c图所示，当μ为3.25时，不论初始高度如何，经过过程期后，球最后在2个高度之中交替变化。不动点的个数随参数的增加而增加。(4)如M2.1d图所示，当μ为3.5时，不论初始高度如何，经过过程期后，球最后在4个高度之中交替变化，只是过渡过程稍微长一点。不动点的个数随参数的增加而进一步增加。(5)如M2.1e图所示，当μ为3.57时，不论初始高度如何，经过过渡期后，球最后在8个高度之中交替变化。(6)如M2.1f图所示，当μ为3.8时，第1个初始相对高度取0.9，第2个初始相对高度取0.2，球的高度杂乱无章地变化。(7)如M2.1g图所示，μ不变，第1个初始相对高度不变，第2个初始相对高度只增加10-5，以后的高度变化也迥然不同。(8)如M2.1h图所示，μ不变，第1个初始相对高度不变，第2个初始相对高度稍大一点(10-8)，以后的高度变化也迥然不同。这种对于初始条件十分敏感的运动称为混沌运动。M2.1e图M2.1f图M2.1g图M2.1h图[解析](2)当参数连续变化时，同样利用(2.5)式计算高度。当迭代次数n足够多的时候，对于周期性的不动点，xn就代表稳定值x∞。取μ为自变量，取x=x∞为函数，可作μ-x曲线。[算法](2)μ从0到4连续取值，先通过迭代算法筛去过渡值，继而用迭代算法获取迭代的结果，画出迭代图。[程序]M2_2.m如下。%罗基斯第模型的倍周期分岔和混沌图clear%清除变量4x=0.2;%初始值(可任取)u=0:0.0001:4;%参数向量(1)n=1000;%迭代次数(2)fori=0:n%按迭代次数循环x=u.*(x-x.^2);%迭代计算消除暂态过程(3)end%结束循环figure%开创图形窗口gridon%加网格xlabel('\it\mu','FontSize',16)%横坐标ylabel('\itx','FontSize',16)%纵坐标title('罗基斯第模型的倍周期分岔和混沌图','FontSize',16)%标题holdon%保持图像cc='bgrk';%颜色代码(4)fori=1:2^6%再按迭代次数循环x=u.*(x-x.^2);%迭代一次plot(u,x,['.'cc(mod(i,4)+1)],'MarkerSize',1)%画点图(5)end%结束循环[说明](1)参数向量的间隔很小，可当作连续分布的。(2)进行1000次迭代，对于周期性运动，可筛去不稳定的点(相对高度)。(3)变量x的初值是一个数值，第一次迭代之后就变成与u同样大小的向量。x的每一个元素都代表u的对应元素的点(相对高度)，这就是用向量的好处。(4)取4种颜色符号。M2.2图(5)在循环中画点时，对于周期运动，画出的同一高度；对于混沌运动，则画不同高度。5每循环4轮用同一颜色画点。[图示]如M2.2图所示，当参数μ从0到1时，高度为零；当μ从0到3时，高度有一个不为零的值；当μ3时，高度首先有两个值，然后分岔为4个值，再分岔为8个值，…，这种情况称为倍周期分岔；当μ达到某一值时，系统进入混沌状态。混沌图还有复杂的结构。[解析](3)在倍周期分岔中，分岔点划分了周期的范围。设二元函数f(μ,x)=μ(1–x)x(2.6)对于周期1不动点，当n→∞时，有xn+1→x∞，xn→x∞，x∞是不动点，用x表示x∞，可得x=f(μ,x)=μ(1–x)x(2.7)由此解得x(1)=0，x(2)=1–1/μ(2.8)不动点x(1)与参数无关，称为平凡不动点。不动点x(2)与参数有关，称为本征不动点。由于0≤xn≤1，所以x(2)≥1。函数对自变量的导数为(,)(12)xffxxx(2.9)不动点的稳定条件是|fx|1(2.10)对于平凡不动点，由于fx(μ,x(1))=fx(μ,0)=μ(2.11)可知：当μ1时，x(1)是稳定不动点；当μ1时，x(1)是不稳定的不动点，或者说x(1)失稳。μ1=1是一个分岔点。对于本征不动点，由于(2)1(,)(,1)2xxfxf(2.12)只有满足-1fx1条件的点才是稳定的，所以当1μ3时，x(2)是稳定的不动点。当μ3时，fx-1，x(2)失稳，因此μ2=3是一个分岔点。分岔值为x(2)=2/3(2.13)对于周期2不动点，则有xn+2→xn→x∞，用函数表示为f[μ,f(μ,x)]=x，简记为f2(μ,x)=x(2.14)指数表示函数嵌套。根据(2.6)式可得μ[1–μ(1–x)x]μ(1–x)x=x分解因式得x[μx–(μ-1)][μ2x2–(μ+1)μx+μ+1]=0(2.15)方程除了x(1)和x(2)两个解之外还有两个解(3)(4)1[1(1)(3)]2x(2.16)可见：当μ3时，x才有不相等的实数解，就是产生周期2的不动点。不动点稳定条件是|fx(μ,x(3))fx(μ,x(4))|1(2.17)设g(μ)=fx(μ,x(3))fx(μ,x(4))，利用(2.9)式可得g(μ)=μ2(1–2x(1))(1–2x(2))=μ2[1–2(x(1)+x(2))+4x(1)x(2)]再利用(2.16)式可得g(μ)=-μ2+2μ+46当g(μ)=-1时，解得16。由g(μ)-1可得[(16)][(16)]0由于610，必有316(2.18)当g(μ)=1时，解得μ=3和-1。由g(μ)1可得μ3。因此，在μ2μμ3范围内有周期2的稳定的不动点。当μμ3时，x(3)和x(4)失稳，因此μ3是分岔点，分岔值为x(3)=0.440，x(4)=0.8499(2.19)[算法](3)当迭代次数n足够多的时候，对于周期性的不动点，xn就代表稳定值x∞。μ从0到μ3连续取值，通过迭代算法最后获得稳定点x，可与解析解进行比较。μ从μ3到4连续取值，先通过迭代算法筛去过渡值，继续用迭代算法获取迭代的结果。[程序]M2_3.m如下。%罗基斯第模型的倍周期分岔图clear%清除变量x=0.2;%初始值(可任取)n=1000;%迭代次数fs=16;%字体大小u3=1+sqrt(6);%分岔参数(1)u=linspace(0,u3);%参数向量fori=0:n%按迭代次数循环x=u.*(x-x.^2);%迭代计算形成不动点end%结束循环figure%开创图形窗口plot([u,u],[x,u.*(x-x.^2)],'.')%画两支迭代曲线(点)holdon%保持图像plot([0,1],[0,0],'r')%画平凡不动点u=1:0.1:3;%分岔前的参数向量x1=1-1./u;%分岔前的不动点plot(u,x1,'r')%画曲线u=linspace(3,u3);%分岔后的参数向量d=(u+1).*(u-3);%根的判别式x21=(1+u+sqrt(d))./(2*u);%分岔后的上支不动点x22=(1+u-sqrt(d))./(2*u);%分岔后的下支不动点plot(u,x21,'r',u,x22,'r')%画2倍分岔后的曲线gridon%加网格xlabel('\it\mu','FontSize',fs)%横坐标ylabel('\itx','FontSize',fs)%纵坐标title('罗基斯第模型2倍周期分岔前后的不动点','FontSize',fs)%标题legend('迭代点','解析曲线',2)%图例axis([0,u3,0,1])%坐标范围x=0.2;%初始值(可任取)7u=3.401:0.0001:4;%参数向量(2)fori=0:n%按迭代次数循环x=u.*(x-x.^2);%迭代计算消除暂态过程end%结束循环figure%开创图形窗口gridon%加网格xlabel('\it\mu','FontSize',fs)%横坐标ylabel('\itx','FontSize',fs)%纵坐标title('罗基斯第模型的倍周期分岔和混沌图','FontSize',fs)%标题holdon%保持图像cc='bgrk';%颜色代码fori=1:2^6%再按迭代