旋轉的Rayleigh-Bénard問題Lorenz模型的動力學行為及數值仿真

王賀元, 陳相霆

(沈陽師范大學 數學與系統科學學院, 沈陽 110034)

1963年,Lorenz[1]將截譜方法作用于Rayleigh-Bénard,在流模型上推導出了三維Lorenz方程。后來,國內外的眾多學者研究了Lorenz模型的高維數推廣、Lorenz吸引子的存在性及性態[2-3]、Couette-Taylor流的力學機理等[4-10]。Bhattacharjee等[11]和張銀[12]先后對旋轉Rayleigh-Bénard問題的四維Lorenz模型進行了討論。本文對他們的模型進行了進一步討論,進行了全局穩定性分析并對系統的動力學行為進行了數值仿真。

1 數學模型

旋轉的Rayleigh-Bénard問題無量綱化后的擾動方程可通過如下的偏微分方程來描述[13]:

(1)

邊界條件為自由邊界?zu|z=0,1=?zv|z=0,1=w|z=0,1=θ|z=0,1=0。式(1)中:u=(u,v,w)表示速度的擾動;θ表示溫度的擾動;p表示壓強;R是Rayleigh數;T是Taylor數,物理意義是旋轉的大小;k=(0,0,1);Pr是Prandtl數。

對方程組(1)化簡后展開,選取模式后代入式(1)得到方程組(2)[13]:

(2)

方程組(2)中的X與對流的強度成比例,Y與上下層流體間的溫度差成比例,Z與垂直溫差的非線性強度成比例,G與流體渦旋強度成比例。

2 全局穩定性分析

2.1 全局吸引子的存在性

對于方程(2)作代換xx,yy,zz+r+Pr后作運算取u(t)=(X(t),G(t),Y(t),Z(t)),H=R4,由上述變量表示可知,參數r,Pr,b均為正數。

2.2 全局穩定性分析

下面利用李雅普諾夫函數法對系統(2)進行全局穩定性分析。構造李雅普諾夫函數

V(x,g,y,z)=13x2+13g2+5y2+5(z-63)2=K>0

很明顯,當K是常數時,上式表示一個四維橢球面,把這個橢球面所包圍的單連通區域記做E,K越大,橢球面越大。在方程組(2)的混沌區(Pr=5,b=8/3,r=50,τ=1)求V的導數:

顯然,下式可以視為一個四維橢球面,記為U:

由李雅普諾夫定理[14]的分析可知,E外的軌線都會進入E內。由此可知,E是這個Rayleigh-Bénard系統的捕捉區。雖然該系統的平衡點都不穩定,但是系統依然具有全局穩定性。

3 數值仿真

改變參量r的取值,借由數值仿真給出的7種指標分析系統的混沌現象。取Pr=5,b=8/3,τ=1,在參數r變化時,對系統的動力學行為進行仿真。

圖1給出了系統(2)隨參數r變化的最大李雅普諾夫指數圖像,可以看出系統的分岔過程對應李雅普諾夫指數的變化。圖2給出了系統(2)隨r變化的分岔圖,它展示了系統分岔和混沌演變的全過程。系統發生混沌后,出現了明顯的倒分岔現象。

圖1 最大李雅普諾夫指數圖像Fig.1 Maximum Lyapunov exponent image

圖2 分岔圖Fig.2 Bifurcationdiagram

如圖3和圖4所示,當 0≤r≤1 000時,平衡點逐漸不穩定,從全局吸引子開始生成2個不穩定的極限環,并且軌線條數隨r的增大而逐漸增多,最終產生了奇怪吸引子。之后吸引子在奇怪吸引子、擬周期吸引子和極限環之間變換,最終變為擬周期吸引子。

圖3 r=20時的吸引子圖像Fig.3 Attractor image at r=20

圖4 r=900時的吸引子圖像Fig.4 Attractor image at r=900

以r=50時為例,給出圖5～8的混沌指標。參考最大李雅普諾夫指數(圖1)和分岔(圖2),可知此時系統(2)處于混沌狀態。

圖5 r=50時的龐加萊截面Fig.5 Poincare section at r=50

圖6 r=50時的時間序列Fig.6 Time series at r=50

圖7 r=50時的功率譜Fig.7 Power spectrum at r=50

圖8 r=50時的返回映射Fig.8 Return mapping when r=50

4 結 語

通過對系統的全局穩定性分析可知,旋轉的Rayleigh-Bénard問題的四維Lorenz模型是全局穩定的;通過對系統模型(2)的數值仿真可知,Rayleigh-Bénard系統的動力學行為隨著參量r取值而變化。在[0,280],奇怪吸引子、擬周期軌道和極限環并存,系統存在倒分岔過程;在[280,500]吸引子穩定為極限環的形式。這些現象說明系統的穩定性是增加的,與文獻[13]中的結論一致。