馮進(jìn)鈐,劉亞妮，王迎宵，李玉婷

(西安工程大學(xué) 理學(xué)院，陜西 西安 710048)

Duffing振動(dòng)系統(tǒng)作為非線性系統(tǒng)的經(jīng)典模型，其全局動(dòng)力學(xué)的研究一直是學(xué)者關(guān)注的熱點(diǎn). 胞映射方法作為一種研究非線性系統(tǒng)全局動(dòng)力學(xué)的有效數(shù)值方法，最先由Hsu[1]在20 世紀(jì)80 年代初提出,隨后得到了廣泛發(fā)展，包括簡單胞映射[2]、插值胞映射[3-4]、廣義胞映射[5]、胞參照點(diǎn)映射[6-7]、圖胞映射等[8-12]. 胞映射方法的生成最常用的方法是采樣點(diǎn)方法. 胞映射方法不僅可以從空間范圍角度呈現(xiàn)系統(tǒng)的全局結(jié)構(gòu)，包括共存吸引子、吸引盆和吸引盆邊界等結(jié)構(gòu)，而且能捕捉到系統(tǒng)的非吸引不變集，如鞍型不動(dòng)點(diǎn)、周期鞍和混沌鞍等. 該方法的優(yōu)點(diǎn)是數(shù)值軌道都是根據(jù)原系統(tǒng)來生成的，但缺點(diǎn)是需要解決大量的常微分方程初值問題的模擬，導(dǎo)致計(jì)算效率較低.

鑒于此，論文基于攝動(dòng)法的思想，將每個(gè)小胞中的采用點(diǎn)看作是該胞中心點(diǎn)的初始擾動(dòng)，討論了每個(gè)胞中不同采用點(diǎn)作為初始點(diǎn)的常微分方程的快速求解，并將其應(yīng)用到圖胞映射方法中，大大提高了計(jì)算效率. 與文獻(xiàn)[13]中多項(xiàng)式逼近算法相比，論文方法公式推導(dǎo)簡單，便于程序?qū)崿F(xiàn). 對于隨機(jī)動(dòng)力系統(tǒng)而言，高斯短時(shí)截?cái)喾椒ū粡V泛采用[14-21]，大大提高了圖胞映射算法的效率. 但該方法受限于高斯白噪聲激勵(lì)下的隨機(jī)動(dòng)力系統(tǒng). 論文方法更具一般性，適用于更一般的Markov過程驅(qū)動(dòng)下隨機(jī)動(dòng)力系統(tǒng). 論文以典型Duffing振動(dòng)系統(tǒng)為應(yīng)用實(shí)例，分析了系統(tǒng)全局吸引子的空間結(jié)構(gòu), 與已有結(jié)果比較，驗(yàn)證了該方法的有效性.

1 軌道攝動(dòng)(orbit perturbation,簡稱OP)積分法

考慮如下常微分方程初值問題描述的非線性系統(tǒng)

dx/dt=f(t,x),x(t0)=x0，

(1)

其中：t∈R是獨(dú)立時(shí)間變量；向量x∈D?Rn表示系統(tǒng)的運(yùn)動(dòng)狀態(tài)，滿足Lipschitz條件；向量函數(shù)f(t,x)∈Rn充分光滑，描述系統(tǒng)的向量場.

對于系統(tǒng)(1)，當(dāng)有非線性項(xiàng)存在時(shí)，其精確解通常是無法得到的，可借助數(shù)值積分法來求解.常用的數(shù)值方法有泰勒級數(shù)法、Runge-Kutta(RK)單步法、Adams-Bashforth(AB)多步法.

首先將積分區(qū)間[t0,t0+T]等分為N等分，分點(diǎn)為tk=t0+kh,k=0,1,…,N,h=T/N為積分步長. 設(shè)φ(x0,t0;t)為向量場f(t,x)產(chǎn)生的流，x(tk)和xk分別為系統(tǒng)(1)在tk時(shí)刻的準(zhǔn)確值和近似值，記x(tk)=φ(x0,t0;tk).4階顯式積分公式對應(yīng)的局部截?cái)嗾`差可以表示為

εk+1=h5Φ(tk,x(tk))+O(h6)，

(2)

其中：函數(shù)Φ(t,x)=λd(5)x/dt(5)為局部誤差函數(shù)的主導(dǎo)項(xiàng)，λ為常數(shù).

由(2)式可以看到，縮小步長h和減小Φ都可以提高精度. 當(dāng)固定精度εk+1時(shí)，主導(dǎo)項(xiàng)h5Φ(tk,x(tk))中步長h與函數(shù)Φ成反比，即當(dāng)Φ越小時(shí)，可以采用較大的步長h，同樣可以獲得相同的精度；另一方面，采用較低精度的數(shù)值積分法也可得到較高精度的解.

基于上述討論，論文介紹軌道攝動(dòng)法. 假設(shè)對于給定的初值x0，系統(tǒng)(1)存在唯一解x(0)(t)=φ(x0,t0;t)，稱為參考軌道. 由方程(1)知

dx(0)/dt=f(t,x(0)),x(0)(t0)=x0.

(3)

對參考軌道的初值x0做微小的擾動(dòng)

得到的唯一解為

稱為目標(biāo)軌道，則

(4)

其中

稱為攝動(dòng)軌道.

由方程(3)和(4)可得其滿足如下微分方程

(5)

其中：g(t,δx)=f(t,x(0)+δx)-f(t,x(0)).

(6)

另外，當(dāng)需要計(jì)算大量的目標(biāo)軌道時(shí)，都可以共用一條參考軌道，利用大步長數(shù)值積分法求解方程(5)，然后借助(6)式即可快速得到目標(biāo)軌道的數(shù)值解. 僅損失存儲參考軌道的內(nèi)存，即可獲得運(yùn)算速度的大幅提升. 該算法可以推廣到一般的非光滑系統(tǒng)數(shù)值仿真中[17-19]，這里不再贅述.

2 Duffing振動(dòng)系統(tǒng)應(yīng)用實(shí)例

考慮如下的經(jīng)典非線性Duffing振動(dòng)系統(tǒng)

(7)

圖1 軌道攝動(dòng)法與RK直接法之間的誤差

從圖1可以看到，當(dāng)k值越大，兩種數(shù)值解的誤差逐漸變大，但在有限的積分時(shí)間內(nèi)，基本滿足一定的精度要求. 同時(shí)，隨著k值的增大，消耗的計(jì)算時(shí)間在逐漸減少，計(jì)算效率在不斷提高. 所得數(shù)值仿真的結(jié)果與上一部分的分析結(jié)果是一致的.

為了進(jìn)一步驗(yàn)證上述方法在解決多初值問題方面的有效性和高效性，取以(0,1)為中心長度為0.1的正方形單元格A，在該單元格內(nèi)隨機(jī)均勻選取10 000個(gè)初始點(diǎn). 對于每個(gè)初始點(diǎn)，分別用4階RK直接積分法和軌道攝動(dòng)法求解系統(tǒng)(7)，積分長度取為T=2π/ω,積分步數(shù)取為N=1 000.在實(shí)施軌道攝動(dòng)法時(shí)，其參考軌道的初始點(diǎn)取為單元格的中心(0,1).

(8)

此外，完成上述數(shù)值仿真RK直接法和軌道攝動(dòng)法所使用的CPU時(shí)間分別為26.79 s和1.24 s. 可見，對于大規(guī)模的胞映射算法的實(shí)現(xiàn)，軌道攝動(dòng)法不僅具有較高的精度，且大大提高了計(jì)算速度.

圖2 多初值情形下RK直接法和軌道攝動(dòng)法的比較圖

固定系統(tǒng)(1)的部分參數(shù)為a=1,b=0.15,c=1,ω=1，考慮外諧和激勵(lì)的強(qiáng)度對系統(tǒng)(7)全局動(dòng)力學(xué)的影響，選取Poincare截面為

(9)

將上述軌道攝動(dòng)積分算法運(yùn)用到廣義圖胞映射方法中，深入分析分叉和激變前后系統(tǒng)(7)的全局動(dòng)力學(xué)行為. 在廣義圖胞映射算法中，感興趣的狀態(tài)空間選取為D={(x1,x2):-2≤x1≤2,-2≤x2≤2}，將感興趣的區(qū)域D均勻劃分為200×200個(gè)胞，利用邊界內(nèi)部組合方法在每個(gè)胞內(nèi)取100個(gè)采樣點(diǎn). 圖3中，符號B(i1,i2,…,ik)表示第i1至第ik個(gè)吸引子的吸引盆邊界，S(i1,i2,…,ik)表示吸引盆邊界B(i1,i2,…,ik)上的鞍. 胞映射方法中采用RK直接法和軌道攝動(dòng)法所得結(jié)果幾乎一致，這里略去RK直接法所得的全局結(jié)構(gòu)圖.

圖3 全局結(jié)構(gòu)圖

3 結(jié)束語

多初值問題的數(shù)值積分是胞映射算法效率的瓶頸. 論文針對胞與胞之間的建立機(jī)制，借助攝動(dòng)法的思想，發(fā)展了軌道攝動(dòng)積分法，大大提高了胞映射算法的計(jì)算效率.

以典型的Duffing振動(dòng)系統(tǒng)為計(jì)算實(shí)例，討論了軌道攝動(dòng)積分法的有效性和計(jì)算效率. 同時(shí)，分析了系統(tǒng)的全局動(dòng)力學(xué). 研究表明：在一定的參數(shù)條件下，系統(tǒng)存在多吸引子共存現(xiàn)象，捕捉到了共存的穩(wěn)定周期吸引子和非吸引的混沌鞍；隨著外諧和力的變化，周期吸引子與邊界混沌鞍發(fā)生碰撞，導(dǎo)致該吸引子及其吸引盆一同消失.