基于matlab的Lorenz系統(tǒng)的仿真研究

1、MATLAB課程期末作業(yè) 以下報(bào)告完成的是大作業(yè)第七題：7. Simulink仿真在高等數(shù)學(xué)課程中的應(yīng)用 21130223 宋沛儒基于MATLAB/Simulink對(duì)Lorenz系統(tǒng)仿真研究21130223 宋沛儒1.引言1963年Lorenz通過(guò)觀察大量大氣現(xiàn)象并進(jìn)行數(shù)值實(shí)驗(yàn)和理論思考，得到了一系列混沌運(yùn)動(dòng)的基本特征，提出了第一個(gè)奇異吸引子Lorenz吸引子1 ，Lorenz通過(guò)計(jì)算機(jī)模擬一個(gè)由三階微分方程描述的天氣模型時(shí)發(fā)現(xiàn)，在某些條件下同一個(gè)系統(tǒng)可以表現(xiàn)出非周期的無(wú)規(guī)則行為。Lorenz揭示了一系列混沌運(yùn)動(dòng)的基本特征，成為后人研究混沌理論的基石和起點(diǎn)，具有非常重要的意義。Lorenz系統(tǒng)

2、方程如下： （1）其中，a，b，c為正的實(shí)常數(shù)。本人利用了數(shù)學(xué)工具matlab，對(duì)Lorenz系統(tǒng)進(jìn)行了仿真研究，加深了對(duì)其的認(rèn)知。2.matlab求解Lorenz系統(tǒng)首先創(chuàng)建文件“Lorenz.m”定義Lorenz方程，假設(shè)固定a=10，b=2.6667，c=30,程序如下：function dx=Lorenz(t,x)dx=-10*(x(1)-x(2);30*x(1)-x(2)-x(1)*x(3);x(1)*x(2)-2.6667*x(3);end然后利用ode45（Runge-Kutta算法）命令求解Lorenz方程并繪制圖形，初值取x=y=z=0.1,程序如下：>> clf

3、>> x0=0.1,0.1,0.1;>> t,x=ode45('Lorenz',0,100,x0);>> subplot(2,2,1)>> plot(x(:,1),x(:,3)>> title('(a)')>> subplot(2,2,2)>> plot(x(:,2),x(:,3)>> title('(b)')>> subplot(2,2,3)>> plot(x(:,1),x(:,2)>> title('(c

4、)')>> subplot(2,2,4)>> plot3(x(:,1),x(:,2),x(:,3)>> title('(d)') 運(yùn)行后，得如下波形：圖中,（a）為L(zhǎng)orenz混沌吸引子在x-z平面上的投影，（b）為其在y-z平面上的投影，（c）為其在x-y平面上的投影，（d）為L(zhǎng)orenz混沌吸引子的三維圖。四張圖都類似于“8”字形。3. Lorenz系統(tǒng)對(duì)初值的敏感性此時(shí)因?yàn)楣潭▍?shù)a=10，b=2.6667，c=30時(shí)，為混沌系統(tǒng)，對(duì)初值具有敏感性，初值很小的差異會(huì)引起系統(tǒng)的大變化。例如在上例中取初值x=z=0.1，y=0.11

5、，繪制此時(shí)混沌吸引子在x-z上的投影，并與x=y=z=0.1在同一張圖比較。（初值為x=y=z=0.1時(shí)投影用藍(lán)色，初值為x=z=0.1，y=0.11時(shí)投影用紅色）程序如下：>> clf>> x0=0.1,0.1,0.1;>> t,x=ode45('ex_lorenz',0,100,x0);>> plot(x(:,1),x(:,3)>> hold on>> x0=0.1,0.11,0.1;>> t,x=ode45('ex_lorenz',0,100,x0);>> pl

6、ot(x(:,1),x(:,3),'r*')得到圖形如下：可以看到初值y僅變化0.01，圖中紅色與藍(lán)色不重合出明顯。證明了Lorenz系統(tǒng)的敏感性。4.matlab對(duì)Lorenz系統(tǒng)的仿真 由文獻(xiàn)1可知在上述方程組（1）中，令，當(dāng)c>1時(shí)，系統(tǒng)有三個(gè)平衡點(diǎn)：， ，。當(dāng)c=1時(shí)，系統(tǒng)在原點(diǎn)失去穩(wěn)定。當(dāng)c<1時(shí)，原點(diǎn)是唯一的平衡點(diǎn)并且是匯點(diǎn)。利用matlab的Simulink功能，搭建Lorenz系統(tǒng)模型，并探討參數(shù)對(duì)Lorenz系統(tǒng)的影響。仿真模型如圖：在仿真模型中，取參數(shù)a=10，b=8/3，觀察參數(shù)c取不同值時(shí)系統(tǒng)的運(yùn)行狀態(tài)。根據(jù)文獻(xiàn)1的分析，當(dāng)參數(shù)0<c&

7、lt;1時(shí)，只有一個(gè)穩(wěn)定平衡點(diǎn)O（0,0,0）。取初值為x=y=z=2,參數(shù)c=0.5，仿真停止時(shí)間取為50，運(yùn)行仿真。得到x、y、z的相圖以及x-z，y-z，x-y的圖形依次如下所示：可見(jiàn)，系統(tǒng)很快地趨向并穩(wěn)定在O（0,0,0），驗(yàn)證了前面所述。當(dāng)c>1時(shí)，系統(tǒng)有三個(gè)平衡點(diǎn)：原點(diǎn)O(0,0,0)和S+，S-。此時(shí)原點(diǎn)的特征值中有正值，因此原點(diǎn)為鞍點(diǎn)，是不穩(wěn)定平衡點(diǎn)。當(dāng)1<c<13.926時(shí)，不穩(wěn)定流形最終螺旋地趨于與之同側(cè)的平衡點(diǎn)S+或S-；當(dāng)c=13.926時(shí)，不穩(wěn)定流形剛好無(wú)限趨于原點(diǎn)O，即出現(xiàn)同宿軌；當(dāng)c>13.926時(shí)，不穩(wěn)定流形將繞到另一側(cè)，最終趨于與之異側(cè)

8、的S+或S-?？梢?jiàn)，c是一個(gè)同宿分岔點(diǎn)。因此，取初值x=y=z=2，c=8，仿真停止時(shí)間為50，運(yùn)行仿真，得到x、y、z的相圖以及x-z，y-z，x-y的圖形依次如下所示：可以看到，系統(tǒng)趨于與之同側(cè)的平衡點(diǎn)S+或S-。取初值x=y=z=2，c=18，仿真停止時(shí)間為50，運(yùn)行仿真，得到x、y、z的相圖以及x-z，y-z，x-y的圖形依次如下所示： 可以看到，系統(tǒng)趨于與之同側(cè)的平衡點(diǎn)S+或S-。為了觀察c=13.926的同宿分岔點(diǎn)現(xiàn)象，在c=13.926附近不斷嘗試，最終在c= 15.39682328時(shí)觀察到比較明顯的過(guò)渡跡象。取初值x=y=z=2，c=15.39682328，仿真停止時(shí)間為50，

9、運(yùn)行仿真，得到x、y、z的相圖以及x-z，y-z，x-y的圖形依次如下所示：可以看到，雖然最終軌線趨向于與之同側(cè)的平衡點(diǎn)S+或S-，但有著明顯的過(guò)渡跡象?？梢酝茰y(cè)，當(dāng)c取15.39682328到15.39682330間的某一個(gè)數(shù)值時(shí)，會(huì)出現(xiàn)同宿軌現(xiàn)象。根據(jù)文獻(xiàn)1，當(dāng)c>24.74時(shí)，S+與S-變?yōu)椴环€(wěn)定的，也就是說(shuō)系統(tǒng)進(jìn)入“混沌區(qū)”。此時(shí)三個(gè)平衡點(diǎn)O、S+、S-都不穩(wěn)定。取初值x=y=z=2，c=30，仿真停止時(shí)間為100，運(yùn)行仿真，得到x、y、z的相圖以及x-z，y-z，x-y的圖形依次如下所示：可以看到，上述圖形中，軌線繞著S+若干圈后，又繞著S-若干圈，如此循環(huán)，符合文獻(xiàn)1的描述。

10、為了觀察由系統(tǒng)趨向于與之異側(cè)的平衡點(diǎn)向系統(tǒng)的混沌狀態(tài)的過(guò)渡現(xiàn)象，在c=24.74附近反復(fù)不斷嘗試，最終發(fā)現(xiàn)當(dāng)c=23.299時(shí)，可以觀察到明顯的過(guò)渡跡象。因此，取初值x=y=z=2，c=23.299，仿真停止時(shí)間為100，運(yùn)行仿真，得到x、y、z的相圖以及x-z，y-z，x-y的圖形依次如下所示：可以看到，在上圖中，軌線看起來(lái)穩(wěn)定在一條圍繞與之異側(cè)的平衡點(diǎn)的軌道上。僅從仿真運(yùn)行的這段時(shí)間，無(wú)法判斷系統(tǒng)是處于混沌狀態(tài)還是會(huì)趨向于與之異側(cè)的平衡點(diǎn)，可以看出明顯的過(guò)渡跡象。5.結(jié)論 本文初步了解了Lorenz系統(tǒng)，并簡(jiǎn)單觀察了Lorenz混沌系統(tǒng)對(duì)初值的敏感性，比較分析在不同參數(shù)下的Lorenz系統(tǒng)仿真結(jié)果，通過(guò)使用matlab的simulink對(duì)Lorenz系統(tǒng)仿真，直觀地觀察到了L