基于時間有限元的非線性系統周期響應求解及穩定性分析

鄭永進, 汪 利, 劉祚秋

(中山大學 航空航天學院,深圳 518107)

非線性現象廣泛存在于實際結構的動力分析之中,主要有幾何非線性、材料非線性與接觸非線性三種來源[1]。在早期的非線性研究中,一般通過線性逼近的方式將非線性問題轉化為線性問題進行求解,對于簡單問題所得近似解能夠滿足工程界的精度要求,但是隨著工程問題的復雜程度增大、分析時間增加,系統內部的非線性因素會使得分析和計算中出現非線性系統特有的現象,如多值響應、極限環、弛緩、分岔以及混沌等[2]。比如,多值響應現象主要出現在受強迫振動的非線性系統中,強非線性系統在諧波激勵下可能有多個解,利用改變激勵頻率或幅值的方式,可以發現系統振幅出現了跳躍現象,也意味著多解的存在[3];Bunonmo[4]在求解Van der Pol系統時,對極限環以及弛緩現象進行了分析研究,在不同的非線性系數下系統所產生的振動現象也不相同;Bernardo等[5]詳細描述了由分岔引起的復雜動力學行為;Venkatesan等[6]對Duffing-van der Pol系統中的分叉與混沌現象進行了分析。

為了解析這些特有的非線性現象,通常需要求得非線性系統的穩態響應(包括周期解),進而分析解的分岔和穩定特性。根據頻域或時域上的分析特點,主要出現了兩類非線性系統穩態響應求解方法,即頻域的諧波平衡法[7]和時域的打靶法[8]。諧波平衡法主要根據穩態響應的頻域特性,使用Fourier級數近似表示周期解或準周期解。該方法能夠直接獲取解的頻響特征,且對于光滑或弱非線性系統,使用少量的級數項就可得到高精度的解,具有極高的計算效率[9];但是,在求解非光滑以及強非線性問題時,諧波平衡法需要選用大量的諧波項,非線性項的諧波系數計算也十分困難,導致計算效率急劇降低。打靶法則是在通用數值積分方法上,引入周期邊界條件求得穩態響應的初值。打靶法能夠適用于各種問題,包括接觸、沖擊和摩擦等非光滑問題[10],且能夠直接得到系統的傳遞矩陣,進而判定周期解的穩定性;但是在面對大規模弱非線性及光滑問題時,計算過程相對來說十分繁瑣、代價較為高昂,且無法直接得到頻域特性。

與以上兩類方法不同,本文旨在提出一種基于時間有限元的非線性系統周期響應求解新方法。時間有限元是一種基于伽遼金變分公式和分段多項式插值的方法,與一般的基于有限差分的數值積分方法相比,時間有限元求解精度高,計算量適中[11],并且在面對激勵快速變化問題時,方法尤為精確[12]。在研究過程中,許多學者基于不同的變分公式、不同問題提出了相應的時間有限元方法[13-15]。Hulbert建立了結構動力學方程的時間不連續伽遼金有限元法,并證明了算法的穩定性和收斂性[16];Wang等[17]建立了一種Petrov-伽遼金型的時間有限元法,證明了在相同的時間步長下,針對線性問題時間有限元法可以提供比Newmark-β法更精確的響應解,周宇等[18]研究表明算法的周期延長率幾乎為0,意味著時間有限元的計算結果不會發生周期漂移,對實際工程的應用十分重要;Qin等[19]提出了一種新型的時間求積單元,顯著提高了計算效率。以上的各類時間有限元法主要面向線性動力系統的瞬態響應分析,本文則是提出一種新型的時間有限元,以求解非線性系統的穩態響應,顯著拓展了時間有限元的應用范疇。該方法在穩態響應伽遼金變分公式基礎上,代入位移響應的三次樣條插值,并結合牛頓迭代法進行求解。由于時間有限元具有精度高、能處理快速變化的激勵等特點,所提方法能高效、高精度地求解非光滑周期荷載下的穩態響應。此外,考慮到時間有限元也能求解瞬態響應,可從系統矩陣中得到周期系統的傳遞矩陣,進而直接根據傳遞矩陣的特征值判斷解的穩定性。

1 非線性系統周期響應的時間有限元方法

1.1 基本原理

考慮周期荷載作用下的結構非線性動力系統,其控制方程的一般形式為

(1)

(2)

為便于后續分析,令外荷載f∈L2(IT)在時間區間IT=[0,T]上平方可積,此時系統的速度和位移也是連續的,考慮周期邊界條件,系統位移u從屬于如下空間

(3)

式中,H2(IT)為二階導數平方可積函數的空間。

(4)

(5)

(6)

(7)

求解線性問題(6)即可從已有位移解uk(t)得到迭代更新解uk+1(t),反復迭代直至滿足收斂條件,最終得到滿意的結果。

(8)

(9)

將插值函數代入線性方程,可建立單元剛度矩陣Di和單元荷載向量Fi

(10)

通過單元組裝,可得到整體剛度矩陣D和荷載向量F。為了便于后續分析,對剛度矩陣和荷載向量作分塊處理,即

(11)

(12)

其中：

(13)

綜上所述,時間有限元的求解過程為：

步驟1給定初值{u}0,設置收斂精度ε=10-6,以及迭代步數k=1;

步驟2基于上一步的解{u}k-1計算整體剛度矩陣D和荷載向量F,引入周期邊界條件得到如式(12)的剛度方程。解剛度方程得到迭代更新解{u}k;

步驟3若不滿足收斂條件,即

(14)

步驟4若滿足收斂條件,則終止迭代,得到問題的周期解。

1.2 周期解的穩定性

(15)

令傳遞矩陣A的特征值中具有最大模的特征值為λ,該特征值也稱Floquet乘子,可直接用于判定系統的穩定性：|λ|<1,周期解是穩定的;|λ|>1,周期解是不穩定的;|λ|=1則通常處于臨界分岔狀態,根據特征值λ實部與虛部的取值,可能出現倍周期分岔、Hopf分岔、鞍結分岔等。

從上可知,周期解穩定性分析的關鍵在于得到傳遞矩陣A。為了獲取傳遞矩陣A,需要計算初值擾動下的系統響應,求解下述線性擾動方程

(16)

(17)

(18)

2 數值算例

本章主要以Duffing系統的受迫振動為例,驗證所提時間有限元在非線性系統周期響應求解和穩定性分析上的有效性和正確性,尤其要突出所提方法在非光滑周期荷載作用下的可行性。

2.1 受余弦荷載作用的單自由度系統

考慮一個含阻尼的單自由度Duffing系統,對其施加簡諧外激勵,其控制方程如下

f(t)=0.2cos(0.8t)

(19)

式中：m=1 kg;c=0.05 N·s/m;k=1 N/m;α=0.1 N/m3;T=2.5πs。

為了進行時間有限元分析,將時間區間[0,T]劃分為N=24個均勻分布的單元。同時,為了驗證解的正確性,以打靶法的解作為參考。圖1給出了時間有限元所求的周期響應,圖2給出了響應的相圖,兩圖的結果均與打靶法的結果進行對比。可以看出,時間有限元的結果與打靶法的結果幾乎重合,驗證了時間有限元在非線性系統周期響應求解中的正確性。此外,還可根據時間有限元計算得到系統的Floquet乘子|λ|=0.140 1,同時,為了對比,還可以計算打靶法的Floquet乘子|λ|=0.140 3。可以看出,時間有限元給出了正確的Floquet乘子計算結果,該結果表明周期解是穩定的。為了進一步探討非線性系數α對系統穩定性的影響,通過弧長法可得到基頻振幅與參數α的關系,如圖3所示,其中實線表示穩定解,虛線表示不穩定解。選取α=-0.1,此時系統存在3個周期解,其中2個為穩定解,1個不穩定。著重考慮不穩定解,此時時間有限元和打靶法的計算結果如表1所示。以打靶法為參考,時間有限元顯然給出了正確的不穩定周期解及穩定性判斷結果。時間有限元顯然可以得到不穩定周期解,并準確判斷出其穩定性。

圖1 單自由度非線性系統位移和速度響應

圖2 單自由度非線性系統相圖

圖3 基頻振幅-非線性項系數關系圖

(20)

表1 打靶法與時間有限元法計算結果對比

圖4給出了初始狀態誤差隨單元尺寸τ的變化趨勢。可以看出,隨著網格加密,時間有限元的初始狀態誤差越來越小,且該誤差達到了2階收斂速度,即error∝τ2。混合荷載示意圖如圖5所示。

圖4 初始狀態誤差隨單元尺寸τ的變化趨勢

圖5 混合荷載示意圖

2.2 受階躍和沖擊周期荷載作用的多自由度系統

考慮如下兩自由度Duffing方程

(21)

圖6 時間有限元法和Newmark法求多自由度非線性系統位移和速度響應的對比圖

使用時間有限元,給定單元數目N=24,可以得到兩個自由度x,y的位移速度響應,如圖6所示。同時,以時間有限元所得到的初始狀態作為初值,使用參數為(γ,β)=(1/2,1/6)的Newmark-β法計算系統的振動響應,經過50個周期后,解達到穩定,取此穩定后的周期解作為參考解。對比時間有限元與參考解(見圖6),可以看出,時間有限元給出了準確的周期響應計算結果,能夠很好地處理沖擊荷載引起的速度突變。同時,通過計算,還可得到Floquet乘子|λ|=0.986 0<1,而參考解是從瞬態響應經過長時間穩定得來,它是穩定的,有限元的解與參考解吻合,也表明了穩定性判斷的正確性。

為了進一步探討外激勵頻率ω對系統穩定性的影響,通過弧長法可得到基頻振幅與外激勵周期頻率ω的關系,最終可得x自由度基頻幅值曲線,如圖7所示,其中實線表示穩定解,虛線表示不穩定解。取ω=2.5 rad/s,此時系統存在3個解,分別如圖7中所標A,B,C三點。A,B,C三點對應的x坐標相位曲線及相應的Floquet乘子如圖8～10所示。可以看出,A點,C點對應的兩個穩定解與參考解完美相符,而B點對應的不穩定解則明顯偏離參考解,這正好驗證了周期響應求解和穩定性判斷的正確性。

3 結 論

本文在伽遼金變分公式的基礎上,提出了一種非線性系統周期響應求解的時間有限元法,并基于時間有限元的系統矩陣,直接求出了系統的傳遞矩陣,進而快速判斷解的穩定性。結合弧長法,還可得到非線性系統的頻響曲線。數值算例表明,所提方法精度高,能有效處理非光滑周期激勵(如階躍或沖擊周期荷載),

能給出正確的穩定性判斷結果,為非光滑問題的周期響應分析提供了一種新的思路。