MATLAB数学建模2-乒乓球的弹跳与罗基斯帝模型.pdf

1 乒乓球的弹跳罗基斯第模型 问题罗基斯第模型 一个乒乓球离球拍的高度为 h0，落在球拍上反弹，设恢复系数为 e，不计空气阻力。 (1)如果 e 为常数，讨论球的高度变化的规律。如果 e2与高度 hn成线性关系 e2= (1 hn/H0)(2.1) 其中 H0是最大高度，是参数。对于不同的参数讨论小球高度的变化规律。 (2)当参数连续变化时，分析最后分布的高度。 (3)计算前几个分岔点。 (4)用李雅普洛夫指数判断混沌的发生。 解析(1)当球从高度 hn下落到球拍上之前速度为 2 nn vgh(2.2) 球与球拍碰撞后反弹的速度为 vn= evn(2.3) 球反弹的高度为 hn + 1= e2hn(2.4) 如果 e 1，例如球拍每次加一个向上的冲击力，则球的高度随次数不断增加。 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个的初始高度 grid minor%加细网格 title(乒乓球与球拍的碰撞高度(itmurm=,num2str(u),),FontSize,16)%标题 n=50;%迭代次数 axis(0,n,0,1)%坐标范围 hold on%保持图像 for j=1:n%按次数循环 xn=u*(1-xn)*xn;%计算下一次的相对高度(3) plot(j,xn,.)%画高度点 2 end%结束循环 xn=0.1;%取初始相对高度(4) plot(0,xn,ro)%画高度点 text(0,xn,num2str(xn),FontSize,16)%初始高度 for j=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%清除变量 4 x=0.2;%初始值(可任取) u=0:0.0001:4;%参数向量(1) n=1000;%迭代次数(2) for i=0:n%按迭代次数循环 x=u.*(x-x.2);%迭代计算消除暂态过程(3) end%结束循环 figure%开创图形窗口 grid on%加网格 xlabel(itmu,FontSize,16)%横坐标 ylabel(itx,FontSize,16)%纵坐标 title(罗基斯第模型的倍周期分岔和混沌图,FontSize,16)%标题 hold on%保持图像 cc=bgrk;%颜色代码(4) for i=1:26%再按迭代次数循环 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 + 1x，xnx，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。 函数对自变量的导数为 ( , )(1 2 ) x f fxx x (2.9) 不动点的稳定条件是 |fx| 1(2.10) 对于平凡不动点，由于 fx(,x(1) = fx(,0) = (2.11) 可知：当 1 时，x(1)是不稳定的不动点，或者说 x(1)失稳。 1= 1 是一个分岔点。 对于本征不动点，由于 (2) 1 ( ,)( ,1)2 xx fxf (2.12) 只有满足-1 fx 1 条件的点才是稳定的，所以当 1 3 时，fx 3 时，x 才有不相等的实数解，就是产生周期 2 的不动点。不动点稳定条件是 |fx(,x(3)fx(,x(4)| -1 可得 (16)(16)0 由于6 1 0 ，必有 3 16 (2.18) 当 g() = 1 时，解得 = 3 和-1。由 g() 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);%参数向量 for i=0:n%按迭代次数循环 x=u.*(x-x.2);%迭代计算形成不动点 end%结束循环 figure%开创图形窗口 plot(u,u,x,u.*(x-x.2),.)%画两支迭代曲线(点) hold on%保持图像 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倍分岔后的曲线 grid on%加网格 xlabel(itmu,FontSize,fs)%横坐标 ylabel(itx,FontSize,fs)%纵坐标 title(罗基斯第模型2倍周期分岔前后的不动点,FontSize,fs)%标题 legend(迭代点,解析曲线,2)%图例 axis(0,u3,0,1)%坐标范围 x=0.2;%初始值(可任取) 7 u=3.401:0.0001:4;%参数向量(2) for i=0:n%按迭代次数循环 x=u.*(x-x.2);%迭代计算消除暂态过程 end%结束循环 figure%开创图形窗口 grid on%加网格 xlabel(itmu,FontSize,fs)%横坐标 ylabel(itx,FontSize,fs)%纵坐标 title(罗基斯第模型的倍周期分岔和混沌图,FontSize,fs)%标题 hold on%保持图像 cc=bgrk;%颜色代码 for i=1:26%再按迭代次数循环 x=u.*(x-x.2);%迭代一次 plot(u,x,. cc(mod(i,4)+1),markersize,1)%画点图 end%结束循环 说明(1)平凡解和二倍周期分岔的范围比较宽，单独画出来比较理想。 (2)多倍周期分岔和混沌的范围比较小，单独画出来比较理想。 图示(1)如 M2.3a 图所示，用迭代法和用解析法计算的结果完全相同。在1= 1 处发生 一次分岔。在2= 3 处发生 2 倍周期分岔，这是树枝分岔。 (2)如 M2.3b 图所示，随着参数的增加，在 2 倍周期分岔后又发生 4 倍周期分岔，8 倍 周期分岔，这种连续的分岔又叫做费根包姆分岔。当周期无限增加时，就失去周期性， 倍周期分岔就走向混沌。混沌之中又有分岔，最明显的是周期 5 和周期 3 分岔，每个分支又 通过倍周期分岔重新走向混沌。在 3.678 处，两个大分支开始相遇。相遇之前的上支与下 支相似，两支又与总体相似。不仅如此，更小的局部分支都与整体相似。这种相似称为自相 似。 M2.3a 图M2.3b 图 解析(4)对于一维映射 xn + 1= f(xn)(2.20) 可用初值 x0和附近值 x0+ x0来计算分离。作一次迭代后，距离为 0 10000 d () ( )() d f x xf xxf xx x (2.21) 再作一次迭代后，距离为 8 1 21111100 d ( ) ( )( )( )() d f x xf xxf xxfx fxx x (2.22) 经过 m 次迭代后，距离为 1 0 0 () m mn n xxfx (2.23) 令 0 exp() mm xxm(2.24) 则得 1 =1 0 11 ln()ln| ()| m m mn n x fx mxm (2.25) 李雅普诺夫指数为 1 =1 1 limln| ()| m n m n fx m (2.26) 当 0 时，表示轨道分离，即对值 具