MATLAB數學建模2-乒乓球的彈跳與羅基斯帝模型.pdf

1 乒乓球的彈跳羅基斯第模型 問題羅基斯第模型 一個乒乓球離球拍的高度為 h0，落在球拍上反彈，設恢復系數為 e，不計空氣阻力。 (1)如果 e 為常數，討論球的高度變化的規(guī)律。如果 e2與高度 hn成線性關系 e2= (1 hn/H0)(2.1) 其中 H0是最大高度，是參數。對于不同的參數討論小球高度的變化規(guī)律。 (2)當參數連續(xù)變化時，分析最后分布的高度。 (3)計算前幾個分岔點。 (4)用李雅普洛夫指數判斷混沌的發(fā)生。 解析(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%開創(chuàng)圖形窗口 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%按次數循環(huán) xn=u*(1-xn)*xn;%計算下一次的相對高度(3) plot(j,xn,.)%畫高度點 2 end%結束循環(huán) xn=0.1;%取初始相對高度(4) plot(0,xn,ro)%畫高度點 text(0,xn,num2str(xn),FontSize,16)%初始高度 for j=1:n%按次數循環(huán) xn=u*(1-xn)*xn;%計算下一次的相對高度(5) plot(j,xn,ro)%畫高度點 end%結束循環(huán) 說明(1)程序執(zhí)行時要用戶用鍵盤輸入參數，提供 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 斷增加，最后穩(wěn)定在一定的高度。高度穩(wěn)定前的過程稱為過渡過程或暫態(tài)過程，過渡過程與 初始高度有關，但是最后高度的穩(wěn)定與初始高度無關。高度值稱為不動點，即重復自身軌跡 的點。 (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)當參數連續(xù)變化時，同樣利用(2.5)式計算高度。當迭代次數 n 足夠多的時候， 對于周期性的不動點，xn就代表穩(wěn)定值 x。取為自變量，取 x = x為函數，可作-x 曲線。 算法(2)從 0 到 4 連續(xù)取值，先通過迭代算法篩去過渡值，繼而用迭代算法獲取迭代 的結果，畫出迭代圖。 程序M2_2.m 如下。 %羅基斯第模型的倍周期分岔和混沌圖 clear%清除變量 4 x=0.2;%初始值(可任取) u=0:0.0001:4;%參數向量(1) n=1000;%迭代次數(2) for i=0:n%按迭代次數循環(huán) x=u.*(x-x.2);%迭代計算消除暫態(tài)過程(3) end%結束循環(huán) figure%開創(chuàng)圖形窗口 grid on%加網格 xlabel(itmu,FontSize,16)%橫坐標 ylabel(itx,FontSize,16)%縱坐標 title(羅基斯第模型的倍周期分岔和混沌圖,FontSize,16)%標題 hold on%保持圖像 cc=bgrk;%顏色代碼(4) for i=1:26%再按迭代次數循環(huán) x=u.*(x-x.2);%迭代一次 plot(u,x,. cc(mod(i,4)+1),MarkerSize,1)%畫點圖(5) end%結束循環(huán) 說明(1)參數向量的間隔很小，可當作連續(xù)分布的。 (2)進行 1000 次迭代，對于周期性運動，可篩去不穩(wěn)定的點(相對高度)。 (3)變量 x 的初值是一個數值，第一次迭代之后就變成與 u 同樣大小的向量。x 的每一個 元素都代表 u 的對應元素的點(相對高度)，這就是用向量的好處。 (4)取 4 種顏色符號。 M2.2 圖 (5)在循環(huán)中畫點時，對于周期運動，畫出的同一高度；對于混沌運動，則畫不同高度。 5 每循環(huán) 4 輪用同一顏色畫點。 圖示如 M2.2 圖所示，當參數從 0 到 1 時，高度為零；當從 0 到 3 時，高度有一個 不為零的值；當 3 時，高度首先有兩個值，然后分岔為 4 個值，再分岔為 8 個值， 這種情況稱為倍周期分岔；當達到某一值時，系統進入混沌狀態(tài)?；煦鐖D還有復雜的結構。 解析(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) 不動點的穩(wěn)定條件是 |fx| 1(2.10) 對于平凡不動點，由于 fx(,x(1) = fx(,0) = (2.11) 可知：當 1 時，x(1)是不穩(wěn)定的不動點，或者說 x(1)失穩(wěn)。 1= 1 是一個分岔點。 對于本征不動點，由于 (2) 1 ( ,)( ,1)2 xx fxf (2.12) 只有滿足-1 fx 1 條件的點才是穩(wěn)定的，所以當 1 3 時，fx 3 時，x 才有不相等的實數解，就是產生周期 2 的不動點。不動點穩(wěn)定條件是 |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)失穩(wěn)，因此3是分岔點，分岔值為 x(3)= 0.440，x(4)= 0.8499(2.19) 算法(3)當迭代次數 n 足夠多的時候，對于周期性的不動點，xn就代表穩(wěn)定值 x。從 0 到3連續(xù)取值，通過迭代算法最后獲得穩(wěn)定點 x，可與解析解進行比較。從3到 4 連續(xù)取 值，先通過迭代算法篩去過渡值，繼續(xù)用迭代算法獲取迭代的結果。 程序M2_3.m 如下。 %羅基斯第模型的倍周期分岔圖 clear%清除變量 x=0.2;%初始值(可任取) n=1000;%迭代次數 fs=16;%字體大小 u3=1+sqrt(6);%分岔參數(1) u=linspace(0,u3);%參數向量 for i=0:n%按迭代次數循環(huán) x=u.*(x-x.2);%迭代計算形成不動點 end%結束循環(huán) figure%開創(chuàng)圖形窗口 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%按迭代次數循環(huán) x=u.*(x-x.2);%迭代計算消除暫態(tài)過程 end%結束循環(huán) figure%開創(chuàng)圖形窗口 grid on%加網格 xlabel(itmu,FontSize,fs)%橫坐標 ylabel(itx,FontSize,fs)%縱坐標 title(羅基斯第模型的倍周期分岔和混沌圖,FontSize,fs)%標題 hold on%保持圖像 cc=bgrk;%顏色代碼 for i=1:26%再按迭代次數循環(huán) x=u.*(x-x.2);%迭代一次 plot(u,x,. cc(mod(i,4)+1),markersize,1)%畫點圖 end%結束循環(huán) 說明(1)平凡解和二倍周期分岔的范圍比較寬，單獨畫出來比較理想。 (2)多倍周期分岔和混沌的范圍比較小，單獨畫出來比較理想。 圖示(1)如 M2.3a 圖所示，用迭代法和用解析法計算的結果完全相同。在1= 1 處發(fā)生 一次分岔。在2= 3 處發(fā)生 2 倍周期分岔，這是樹枝分岔。 (2)如 M2.3b 圖所示，隨著參數的增加，在 2 倍周期分岔后又發(fā)生 4 倍周期分岔，8 倍 周期分岔，這種連續(xù)的分岔又叫做費根包姆分岔。當周期無限增加時，就失去周期性， 倍周期分岔就走向混沌?；煦缰杏钟蟹植?，最明顯的是周期 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 時，表示軌道分離，即對值 具