【文章內(nèi)容簡介】 beita=pi*d*sinsita./lamda。 %計(jì)算 β I=(sin(arfa).^2./(arfa.^2)).*... (sin(N.*beita).^2./(sin(beita).^2))。 %相對(duì)光強(qiáng)分布 figure。 %開辟圖形窗口 plot(x,I)。 %畫光強(qiáng)與觀察點(diǎn)位置關(guān)系圖 figure。view(0,90)。 %新開圖形窗口并在 xy 平 面內(nèi)觀察 hold on colormap(gray)。 %選擇灰度色圖 mesh(X,Y,I)。 %繪制衍射圖樣 運(yùn)行結(jié)果如圖 16(a)和圖 16(b)所示： 11 圖 16(a) 光強(qiáng)與位置的關(guān)系 圖 16(b) 光柵衍射的模擬圖樣 質(zhì)點(diǎn)在萬有引力作用下的運(yùn)動(dòng) 以萬有引力的固定不動(dòng)的施力質(zhì)點(diǎn) 0m 所在位置為坐標(biāo)原點(diǎn) O , 建立直角坐標(biāo)系 Oxy ，質(zhì)點(diǎn)的運(yùn)動(dòng)微分方程為 03Gmmmr rr??，分量方程為： 002 2 2 22 2 2 2( ) ( )G m G mxyxyx y x yx y x y? ? ? ?????， (14) 這兩個(gè)方程都是二階常微分方程，定義解矢量為 y ，令 (1 ) ( 2 ) ( 3 ) ( 4 )y x y x y y y y? ? ? ?， ， ， (15) 可將方程組 (14)化為： ? ?? ?0322203222( 1 )( 1 ) ( 2)( 2)( 1 ) ( 3 )( 3 )( 3 ) ( 4)( 4)( 1 ) ( 3 )G m ydy dyydt dt yyG m ydy dyydt dt yy?? ? ???? ? ??， (16) (1)編寫微分方程組函數(shù)文件 ： function ydot=yxlcfun(t,y,flag,p) %函數(shù)首行 ,p 為參量 Gm0 ydot=[y(2)。 p*y(1)/sqrt(y(1).^2+y(3).^2).^3。 y(4)。 p*y(3)/sqrt(y(1).^2+y(3).^2).^3]。 %建立微分方程組 (2)解微分方程的主程序 ： p=1。 %取 Gm0=1 y0=[10 6 。25 5 0。25 6 0]。 %三組不同初始條件 plot(0,0, 39。*r39。) %畫出 O 點(diǎn) for i=1:3 %分別以不同初始條件解 3 次方程 12 [t,y]=ode45(39。yxlcfun39。,[0::300],y0(i,:),[ ],p)。 hold on axis([25 25 20 20])。 %指定坐標(biāo)范圍 et(y(:,1),y(:,3)) %繪出質(zhì)點(diǎn)運(yùn)動(dòng)軌跡 (x,y) end %結(jié)束循環(huán) 解出的結(jié)果如圖 17 所示： 圖 17 萬有引力場中質(zhì)點(diǎn)運(yùn)動(dòng)軌跡 由上面例子，我們初步了解了 Matlab 解常微分方程的一般過程，首先是建立微分方程函數(shù)文件，文件的格式如下： fuction ydot=filename(t,y,p1,p2) %t,y 是積分區(qū)間和解矩陣 p1,p2 是參數(shù) ydot=[關(guān)于 t,y 的表達(dá)式 ]。 %ydot 表示 dy/dt 下面介紹 ode45 命令的用法， ode45 的一般調(diào)用格式為： [T,Y]=ode45(39。fun39。,tspan,y0,options,p1,p2,…) 其中含義如下表： 表 16 ode45 命令含義 Fun 求解的微分方程函數(shù)名 Tspan 單調(diào)遞增（減）的積分區(qū)間 [t0:tstep:tfinal] y0 初始條件矢量 Options 用 odeset 建立的優(yōu)化選項(xiàng)，一般用默認(rèn)值，為空矢量 “[ ]” 13 p1,p2 傳遞給 fun 函數(shù)的參數(shù) T,Y T 是輸出的時(shí)間列矢量，矩陣 Y 的每一個(gè)列矢量是解的一個(gè)分量 各個(gè)項(xiàng)在命令中的位置和順序不能顛倒，否則程序就會(huì)出錯(cuò) 。 167。 本章小結(jié) 本章首先對(duì) Matlab 進(jìn)行了簡單介紹，介紹了 Matlab 的發(fā)展及應(yīng)用前景， 然后詳細(xì)介紹了基本運(yùn)算功能，基本繪圖功能，數(shù)值分析功能，并且簡單介紹了 M 文件的編寫及 Matlab 的程序調(diào)試。 由于 Matlab 是集數(shù)值運(yùn)算、符號(hào)運(yùn)算、數(shù)據(jù)可視化、數(shù)據(jù)圖文字統(tǒng)一處理、系 統(tǒng)動(dòng)態(tài)仿真等功能于一體的數(shù)學(xué)軟件，所以為了加深對(duì) Matlab 的基本功能的理解，在本章第三節(jié)我們列舉了幾個(gè)簡單的應(yīng)用。 例 1 等量異號(hào)點(diǎn)電荷的電勢分布，應(yīng)用了 Matlab 三維網(wǎng)格作圖命令mesh(x,y,z)和基本函數(shù) 數(shù)值運(yùn)算功能。 例 2 光柵衍射， 應(yīng)用了 Matlab 的 基本運(yùn)算功能，基本繪圖功能，數(shù)值分析功能。 例 3 質(zhì)點(diǎn)在萬有引力作用下的運(yùn)動(dòng) ， 應(yīng)用了 Matlab 的 基本運(yùn)算功能，基本繪圖功能，在一個(gè)窗口下繪制多條圖形，并且利用 Matlab 的 ode45 命令求解常微分方程，最后 總結(jié)了 Matlab 的 ode45 命令 解常微 分方程的一般過程。 14 第二章 混沌行為與特性 167。 混沌理論 在現(xiàn)代物理的研究中，混沌理論的建立可能稱得上是最重要的成就之一?；煦绲陌l(fā)現(xiàn)被譽(yù)為繼相對(duì)論和量子力學(xué)后的第三次物理學(xué)革命，混沌的研究一直備受學(xué)術(shù)界的關(guān)注。經(jīng)典物理的確定論和近代量子物理的隨機(jī)論，雖然都非常成功地解決了許多自然現(xiàn)象，但是這兩種理論之間似乎存在著對(duì)立的矛盾，只適宜于不同的領(lǐng)域。在宏觀領(lǐng)域似乎只應(yīng)用經(jīng)典物理理論，而在微觀領(lǐng)域中更多是使用量子力學(xué)、統(tǒng)計(jì)物理的隨機(jī)性理論。 按照確定性理論，物體運(yùn)動(dòng)以后在任何時(shí)間狀態(tài) 只是一個(gè)點(diǎn)，而按照隨機(jī)性理論，物體的狀態(tài)不是點(diǎn)，而是由點(diǎn)組成的點(diǎn)云，點(diǎn)云的密度就表示物體出現(xiàn)這種狀態(tài)的概率大小。隨著物理理論的深人研究，人們發(fā)現(xiàn)在傳統(tǒng)的確定性理論領(lǐng)域中，一定的條件下也可以出現(xiàn)一定的隨機(jī)性現(xiàn)象。這種隨機(jī)性又不同于傳統(tǒng)的隨機(jī)性理論研究的隨機(jī)性問題，因?yàn)樗钟幸欢ǖ拇_定性。這種現(xiàn)象的發(fā)現(xiàn)和解決從而誕生了混沌理論。 在經(jīng)典物理中，物體的狀態(tài)變化規(guī)律都是用非線性方程來描述的，只要給出一定初始條件，就可以解出這個(gè)非線性方程的解，事實(shí)上，如果對(duì)這些非線性方程進(jìn)一步研究，就會(huì)發(fā)現(xiàn)，初始條件的變化，可以使這 些本來是確定性的解出現(xiàn)隨機(jī)性。 167。 簡單的數(shù)學(xué)游戲 (1)令 1 2nnXX? ? ，如果是利用二進(jìn)制，按 Modell（去掉整數(shù)，只留小數(shù)的操作）進(jìn)行迭代 [21]，我們來看看其結(jié)果如何。 取 0X = 0. 1000100100111001 1X = 2X = 3X = ?? 16X =0 也就是說，其結(jié)果最終趨向于零，或者說是趨向于一個(gè)點(diǎn)，一般我們把 15 這個(gè)點(diǎn)叫“匯”。而如果取 0X? = 001??如此循環(huán)下去，我們就會(huì)發(fā)現(xiàn)，其結(jié)果應(yīng)該有 16 種穩(wěn)定的點(diǎn)，我們叫“穩(wěn)定極限環(huán)”。 0X 和 0X? 雖然在數(shù)學(xué)上可以說是非常接近，只存在微小的不同，但是其計(jì)算的結(jié)果卻完全不同。這就說明，微小的初始條件的變化，可能引起不同的結(jié)果，這就是初始條件的敏感性。 (2)令 1 (1 )n n nX X X?? ??。 [0,1]nX ? 如果取 ? =2， 0X =，則 1X =， 2X =， 3X =， 4X = 999996，?? nX =，也就是其結(jié)果最終趨向于點(diǎn) X =0. 50。 如果取 ? =， 0X =，則 1X =,?? 8X =， 9X = 565。 如果繼續(xù)迭代下去，其結(jié)果仍然是在 8X 和 9X 上循環(huán)，這就是說，其迭代的結(jié)果可能是兩個(gè)，我們就叫“分岔”，像這樣有兩個(gè)結(jié)果的就叫“二周期”。 如果取 ? =， 0X =，則 1X =0. 5648，?? 13X =， 14X = 8230286， 15X = , 16X =。如果繼續(xù)迭代下去，就是進(jìn)行這四數(shù)字的循環(huán)，像這樣的結(jié)果就叫“四周期”。 如果取 ? =3. 9，那么我們就會(huì)發(fā)現(xiàn)其結(jié)果是無周期的，或者說是有無數(shù)的點(diǎn)，也就是進(jìn)人完全隨機(jī)狀態(tài)，我們就叫“混沌”（ chaos）。 167。 “蝴蝶效應(yīng)” 美國的著名氣象學(xué)家 Edward Lorenz 從旋轉(zhuǎn)的木桶實(shí)驗(yàn) [21]，總結(jié)出包括12 個(gè)方程的方程組，建立了一個(gè)仿真的氣象模型，他認(rèn)為盡管氣象變化萬千，但總是遵循經(jīng)典的物理定律，只要知道一定的初始條件，那么利用這些方 程總是可以把結(jié)果算出來的。這就是說按照傳統(tǒng)的確定性理論，他就可以確定將來的氣象變化的規(guī)律和任何時(shí)間的氣象狀態(tài)。這里需要說明的是，一般傳統(tǒng)的科學(xué)家都認(rèn)為，任何量的測量和獲得都不可能是完全精確的，都有一定的近似，所以在進(jìn)行計(jì)算的時(shí)候一般都采用一定的近似。因?yàn)樗麄冋J(rèn)為，極小的影響和變化、差別是可以忽略不計(jì)的，事物運(yùn)動(dòng)之中都具有一定的收斂性，極小的差別不會(huì)引起大的影響。 Lorenz 在利用計(jì)算機(jī)進(jìn)行計(jì)算的時(shí)候，一次，為了省時(shí)，他就把上次計(jì)算打印結(jié)果當(dāng)作初始值輸人了，然而，當(dāng)他一小時(shí)以后回來的時(shí)候突然發(fā)現(xiàn)其結(jié)果卻偏差 極大。開始他以為是計(jì)算機(jī)出了問題，后來經(jīng)過仔細(xì)的研究， 16 發(fā)現(xiàn)是由于初始值的微小差別導(dǎo)致其結(jié)果的極大偏差。因?yàn)槟菚r(shí)候的計(jì)算機(jī)還很簡單，存儲(chǔ)只是 6 位，但是打印出來的只是 3 位，例如輸入 0. 532020，只能打印出來 ，當(dāng)時(shí)他認(rèn)為這是極小差別，不會(huì)引起大的變化。但他的方程對(duì)這些微小的不同卻是極其敏感的，他把這種現(xiàn)象叫作“蝴蝶效應(yīng)”。意思就是：巴西的蝴蝶抖動(dòng)一下翅膀，就可能在德克薩斯引起一場風(fēng)暴。 蝴蝶效應(yīng)說明了初始條件的重要性，也說明了科學(xué)的嚴(yán)謹(jǐn)。任何隨意的忽略，都可能導(dǎo)致嚴(yán)重的后果。也正由此導(dǎo)致后來“混沌”理論的誕生。 167。 用 Matlab 演示混沌的基本性質(zhì) 在自然界中，絕大部分運(yùn)動(dòng)都是混沌運(yùn)動(dòng)，規(guī)則運(yùn)動(dòng)只在局部的范圍和較短的時(shí)間內(nèi)存在。從簡單的數(shù)學(xué)游戲入手，我們了解了混沌運(yùn)動(dòng)的產(chǎn)生，通過對(duì)“蝴蝶效應(yīng)”的介紹，我們了解了混沌運(yùn)動(dòng)的主要特征及性質(zhì)。在各種軟件中， Matlab 是非常適合混沌的演示和仿真實(shí)驗(yàn)的。本節(jié)將對(duì)如何使用Matlab 來演示混沌運(yùn)動(dòng)特征及性質(zhì)進(jìn)行研究。 167。 用 Matlab 產(chǎn)生標(biāo)準(zhǔn)的混沌信號(hào) 1963 年，美國氣象學(xué)家洛倫茲在《大氣科學(xué)雜志》上發(fā)表了著名的論文《確定性的非周期流》，文中指出：三階非線性自治系統(tǒng)中可能會(huì)出現(xiàn)混沌解[12]。洛倫茲提出了一個(gè)簡化的天氣預(yù)報(bào)模型，這就是著名的洛倫茲方程組： ()()x a y xy b z x yz xy cz????? ? ?????? (21) 這個(gè)簡化模型是一個(gè)完全確定的方程組。然而，當(dāng)方程組的三個(gè)參數(shù)取某些值時(shí) (比較常用的是 a =10,b =28,c =8/3)，方程組出現(xiàn)了混沌解。這是在耗散系統(tǒng)中，一個(gè)確定的方程能導(dǎo)出混沌解的第一個(gè)實(shí)例，它標(biāo)志著混沌學(xué)的涎生。 在 Matlab 中，可以用如下程序 產(chǎn)生洛淪茲信號(hào)，在對(duì)混沌信號(hào)的演示和處理中，洛淪茲信號(hào)是最常用到的標(biāo)準(zhǔn)混沌信號(hào)?；煦缦到y(tǒng)存在混沌吸引子，洛淪茲吸引子就是著名的蝶形圖。如圖 21 所示。 (1)洛倫茲函數(shù)程序： function dy=lorenz(t,y) 17 dy=zeros(3,1)。 dy(1)=10*(y(1)+y(2))。 dy(2)=28*y(1)y(2)y(1)*y(3)。 dy(3)=y(1)*y(2)8*y(3)/3。 (2)給定參數(shù)作圖程序： [t,y]=ode45(39。lorenz39。,[0 30],[12,2,9])。 plot3(y(:,1),y(:,2),y(:,3)) view([20,42])。 圖 21 洛淪茲信號(hào)的吸引子 167。 倍周期分岔 —— 通向混沌之路 倍周期分岔是許多非線性動(dòng)力學(xué)過程中常見的現(xiàn)象，也是進(jìn)人混沌的一種重要方式 [ 12]。可以用描述蟲日模型的 Logistic 方程來演示一個(gè)動(dòng)力學(xué)系統(tǒng)是如何通過倍周期分岔從規(guī)則運(yùn)動(dòng)進(jìn)人混沌運(yùn)動(dòng)的。 Logistic 差分方程為 1 (1 )n n nX X X?? ??，初值 0X 的取值范圍為 (0,1)， ? 的取值范圍為 [1,4]。由 Logistic 方程描述的系統(tǒng)的最終狀態(tài)取決于 ? 值， ? 在由1 變化到 4 的過程中，該系統(tǒng)通過不斷的倍周期分