一種剛性桿-彈簧擺模型的混沌動力學(xué)行為研究

趙　聰，于洪潔

(上海交通大學(xué)工程力學(xué)系，上海 200240)

?

一種剛性桿-彈簧擺模型的混沌動力學(xué)行為研究

趙聰，于洪潔

(上海交通大學(xué)工程力學(xué)系，上海 200240)

針對具有大范圍運(yùn)動慢變量和小幅度振蕩快變量的強(qiáng)非線性剛-柔耦合多體系統(tǒng)，建立一種剛性桿-彈簧擺模型。給出了該雙時間尺度變量系統(tǒng)的無量綱動力學(xué)方程，以頻率比、擺長比作為控制參數(shù)，對系統(tǒng)在不同初始條件下的非線性動力學(xué)行為進(jìn)行了數(shù)值模擬和分析。由于快、慢變量之間的相互耦合，動力學(xué)方程表現(xiàn)出強(qiáng)非線性的特點(diǎn)，對數(shù)值方法提出了更高要求。采用一種高精度的三次Lagrange插值精細(xì)積分法進(jìn)行數(shù)值求解，并給出了系統(tǒng)不同的運(yùn)動狀態(tài)對應(yīng)的參數(shù)范圍。數(shù)值分析結(jié)果表明，系統(tǒng)變量在不同的控制參數(shù)和初始條件下，呈現(xiàn)出了復(fù)雜的混沌動力學(xué)行為，快變量顯示了經(jīng)由準(zhǔn)周期環(huán)面破裂分岔通往混沌的途徑。

剛-柔耦合系統(tǒng)；雙時間尺度；彈簧擺；插值精細(xì)積分法；混沌

剛-柔耦合多體系統(tǒng)是目前常見的力學(xué)模型，在工程中有著非常廣泛的應(yīng)用，例如直升機(jī)旋翼、柔性機(jī)械臂、渦輪機(jī)葉片等等。該類系統(tǒng)的大范圍剛體運(yùn)動與彈性小變形運(yùn)動的耦合使得系統(tǒng)的動力學(xué)性態(tài)十分復(fù)雜，其動力學(xué)方程具有強(qiáng)非線性，隨之帶來的動力剛化現(xiàn)象、數(shù)值求解的困難等是建立該類模型時需要重點(diǎn)考慮和解決的問題。

平面雙擺剛體系統(tǒng)是混沌動力學(xué)中的經(jīng)典模型，Stachowiak T等[1]和Nikitina N V[2]對平面雙擺系統(tǒng)的混沌運(yùn)動等動力學(xué)行為作了相關(guān)研究。彈簧擺系統(tǒng)也是比較典型的非線性運(yùn)動模型之一，在不考慮彈簧的扭振和扭轉(zhuǎn)變形，只考慮彈簧縱向伸長時，彈簧擺是一個簡單的二自由度系統(tǒng)。近年來，國內(nèi)外學(xué)者對彈簧擺作了許多方面的研究，例如彈簧擺的內(nèi)共振現(xiàn)象[3-5]、混沌運(yùn)動[6-7]等。文獻(xiàn)[8]研究了彈簧擺系統(tǒng)在頻率比以及初值條件變化時的混沌、擬周期運(yùn)動等復(fù)雜動力學(xué)行為。對于將平面雙擺剛體系統(tǒng)與彈簧擺系統(tǒng)結(jié)合起來構(gòu)成的更復(fù)雜的剛性桿—彈簧擺系統(tǒng)的研究目前還鮮有發(fā)表。

鐘萬勰[9-10]針對線性定常結(jié)構(gòu)動力系統(tǒng)提出了一種精細(xì)時程積分法，這種方法具有高精確度、高效、穩(wěn)定的特點(diǎn)，它的解在積分點(diǎn)處數(shù)值上逼近于精確解。張洵安[11]等指出非線性精細(xì)積分算法可有效地計算多自由度系統(tǒng)的混沌運(yùn)動。精細(xì)積分法經(jīng)眾學(xué)者的改進(jìn)，已經(jīng)得到了長足發(fā)展和廣泛應(yīng)用，例如更新精細(xì)積分方法[12]、增維精細(xì)積分法[13-15]，三次Lagrange插值精細(xì)積分法[16]等等。其中，本文作者提出的三次Lagrange插值精細(xì)積分法[16]是一種單步顯式、預(yù)估-校正的具有四階精度的算法，具有自起步、精度高及計算量小的特點(diǎn)，適用于求解強(qiáng)非線性動力學(xué)方程。

本文結(jié)合了平面雙擺系統(tǒng)和彈簧擺系統(tǒng)的特點(diǎn)，構(gòu)造了更為復(fù)雜的剛性桿-彈簧擺模型，該模型是三自由度雙時間尺度變量系統(tǒng)。其中，兩個慢變量分別為剛性桿和彈簧的擺角，對應(yīng)剛-柔耦合系統(tǒng)的大范圍運(yùn)動；快變量為彈簧伸長率，對應(yīng)剛-柔耦合系統(tǒng)的彈性小變形。以參數(shù)頻率比描述快、慢變量之間的時間尺度差異。由于變量之間的耦合，系統(tǒng)動力學(xué)方程具有強(qiáng)非線性的特點(diǎn)，采用本文作者提出的一種三次Lagrange插值精細(xì)積分法[16]對無量綱化的動力學(xué)方程進(jìn)行數(shù)值求解。結(jié)合Poincaré映射和最大Lyapunov指數(shù)等數(shù)值手段，對系統(tǒng)的非線性動力學(xué)行為進(jìn)行分析。

1　模型的建立及動力學(xué)方程推導(dǎo)

如圖1所示，建立剛性桿-彈簧擺平面運(yùn)動模型。質(zhì)量分布均勻的剛性細(xì)桿在x-y平面內(nèi)繞支點(diǎn)O擺動，彈簧擺與剛性桿的下端鉸接。假設(shè)彈簧不發(fā)生扭轉(zhuǎn)變形和扭振，并且不計質(zhì)量，小球可視為質(zhì)點(diǎn)。系統(tǒng)受y方向的重力作用作平面運(yùn)動，不計一切摩擦力。其中：剛性桿的質(zhì)量為m1，長度為l1，相對于豎直方向的擺動角度為θ1；彈簧的剛度為k，原長為l0，動態(tài)擺長為r，相對于豎直方向的擺角為θ2，小球質(zhì)量為m2。

設(shè)系統(tǒng)靜平衡時彈簧的長度為l，由受力分析得到：l=l0+m2g/k。根據(jù)第二類Lagrange方程建立系統(tǒng)的動力學(xué)方程：

(1)

(2)

方程(2)描述了一個三自由度雙時間尺度變量剛—柔耦合系統(tǒng)，其中剛性桿的擺角θ1和彈簧擺的擺角θ2是慢變量，彈簧伸長率x是快變量，它們之間相互耦合；該方程為無量綱方程，具有強(qiáng)非線性，含有兩個控制參數(shù)：頻率比γ1、擺長比γ2。

2 三次Lagrange插值精細(xì)積分法

動力學(xué)系統(tǒng)微分方程的一般形式為

(3)

其中,x是n維向量，n是自由度數(shù)；r(t)是時變的載荷向量；α和β是已知的初始條件。

由于精細(xì)積分法易于處理一階方程，因而將式(3)化成下列一階常微分方程

(4)

其中，v是2n維向量，H是2n×2n階的滿秩常矩陣，f(t)是由r(t)構(gòu)成的已知列向量。

vk=v(kτ)=exp(Hkτ)v0=exp(Hτ)vk-1=Tvk-1

(5)

其中，矩陣T可由文獻(xiàn)[1]給出的方法精確計算出來，進(jìn)而計算出解向量，得到各離散時刻的解vk。

(6)

P3(t,v(t))=r0+r1(t-tk)+r2(t-tk)2+r3(t-tk)3

(7)

插值基點(diǎn)取

r0=f(tk,v(tk))≡fk

(8)

用Euler公式取vk+1的預(yù)測值[16]，將式(8)代入式(7)求出P3(t,v)，將其代入式(6)，令t=tk+1，積分得

(9)

(10)

式(10)是具有強(qiáng)非線性的無量綱一階微分方程組。由式(4)，構(gòu)造出非線性函數(shù)f(t,v)、H0分別為

3　剛性桿-彈簧擺系統(tǒng)非線性動力學(xué)行為分析

本節(jié)主要研究剛性桿-彈簧擺系統(tǒng)在不同的初始條件下隨著控制參數(shù)的變化而表現(xiàn)出的不同動力學(xué)行為。剛-柔耦合多體系統(tǒng)具有快、慢變量及多參數(shù)的特點(diǎn)?？紤]到工程實(shí)際，參數(shù)頻率比γ1的取值不宜過大或過小。根據(jù)頻率比的定義，γ1值的大小一定程度上反映了系統(tǒng)快、慢變量之間時間尺度相差程度的大小。若頻率比值取得過大，系統(tǒng)表現(xiàn)出很強(qiáng)的剛性，難以和多剛體系統(tǒng)區(qū)分開來；反之，則表現(xiàn)出很強(qiáng)的柔性。因此，首先設(shè)定一組固定參數(shù)：頻率比參數(shù)為定值γ1=50，擺長比參數(shù)γ2=1，剛性桿初始擺角θ10=45°，彈簧擺初始時刻為豎直狀態(tài)，即θ20=0°，彈簧初始伸長率x0=0.001。通過數(shù)值計算，觀察系統(tǒng)變量的時程曲線、相軌跡圖，并結(jié)合Poincaré映射和最大Lyapunov指數(shù)等數(shù)值手段對系統(tǒng)變量的運(yùn)動狀態(tài)進(jìn)行分析。

圖2a、b、c分別為系統(tǒng)慢變量剛性桿擺角θ1、彈簧擺擺角θ2和快變量彈簧伸長率x在θ10=45°時的Poincaré映射圖。圖2a顯示剛性桿擺角θ1的Poincaré映射點(diǎn)組成帶狀的環(huán)面；圖2b、2c顯示彈簧擺擺角θ2和彈簧伸長率x的Poincaré映射點(diǎn)為模糊一片的密集點(diǎn)集，表現(xiàn)出混沌運(yùn)動的特征；從圖中還可看出，彈簧擺擺角θ2的最大值小于0.15 rad (8.6°)，盡管剛性桿的初始擺角達(dá)到45°，彈簧擺的擺動仍處于小角度范圍內(nèi)。經(jīng)計算得出，此時系統(tǒng)最大Lyapunov指數(shù)λ1=0.909 014 7，可判定出現(xiàn)混沌運(yùn)動。

下面考慮剛性桿初始擺角θ10的值在5°～50°大范圍內(nèi)變化，其他參數(shù)條件不變。圖3a～d分別顯示當(dāng)θ10=5°,20°,40°,50°時，系統(tǒng)快變量彈簧伸長率x的Poincaré映射圖。從圖中可以看出，當(dāng)θ10=5°時，快變量x的Poincaré映射點(diǎn)近似于準(zhǔn)周期運(yùn)動對應(yīng)的一條封閉曲線，隨著θ10值的增大，Poincaré映射環(huán)面逐漸破裂，最終形成一大片模糊點(diǎn)集。圖3a～d四種情況下對應(yīng)的系統(tǒng)最大Lyapunov指數(shù)λ1分別為1.828 782×10-5、2.036 442×10-5、0.087 738 55、0.946 333 2。由此可推測快變量彈簧伸長率x可能經(jīng)由準(zhǔn)周期環(huán)面破裂分岔這一途徑通往混沌。

隨后，設(shè)定頻率比參數(shù)為定值γ1=20，考慮擺長比參數(shù)γ2的變化對系統(tǒng)狀態(tài)變量動力學(xué)行為的影響。仍然取彈簧初始伸長率x0=0.001，彈簧擺初始時刻豎直，即θ20=0°，令γ2在0.25～4范圍內(nèi)變化，剛性桿初始擺角θ10在5°～50°大范圍內(nèi)變化，觀察參數(shù)組合條件下的系統(tǒng)變量Poincaré映射圖，并計算相應(yīng)的系統(tǒng)最大Lyapunov指數(shù)。表1給出在部分γ2值下系統(tǒng)最大Lyapunov指數(shù)的變化情況，圖4是相應(yīng)的變化曲線。

由表1可以看出，在所研究的參數(shù)范圍內(nèi)，計算得到的系統(tǒng)最大Lyapunov指數(shù)均為正數(shù)。隨著剛性桿初始擺角θ10的增大，系統(tǒng)最大Lyapunov指數(shù)有增大的趨勢。擺長比γ2越小，系統(tǒng)最大Lyapunov指數(shù)增大到10-1量級時對應(yīng)的θ10臨界值越小，結(jié)合圖4，即擺長比γ2越小，系統(tǒng)最大Lyapunov指數(shù)隨剛性桿初始擺角θ10增長的越快。

4　結(jié)論

構(gòu)建了一種三自由度雙時間尺度變量剛—柔耦合系統(tǒng)——剛性桿—彈簧擺模型，體現(xiàn)了剛—柔耦合多體系統(tǒng)中大范圍運(yùn)動慢變量與彈性小變形快變量耦合的特點(diǎn)。對于此類容易產(chǎn)生剛性問題的強(qiáng)非線性動力學(xué)方程，采用一種三次Lagrange插值精細(xì)積分法進(jìn)行數(shù)值求解，該方法具有單步顯式、預(yù)估—校正、自起步、四階精度等特點(diǎn)。將頻率比、擺長比以及剛性桿初始擺角作為控制參數(shù)，研究了慢變量、快變量不同的動力學(xué)行為。發(fā)現(xiàn)了慢變量的混沌行為和快變量經(jīng)由準(zhǔn)周期環(huán)面破裂通往混沌的途徑。工程實(shí)際中包含很多類似的剛—柔耦合多體系統(tǒng)，例如柔性機(jī)器人手臂等，研究本文模型大范圍運(yùn)動和自身形變在合理參數(shù)范圍內(nèi)的耦合非線性動力學(xué)行為，對工程中針對此類模型進(jìn)行設(shè)計和控制具有相當(dāng)?shù)膮⒖純r值。

表1　系統(tǒng)最大Lyapunov指數(shù)λ1隨控制參數(shù)變化的情況

[1]Stachowiak T, Okada T. A numerical analysis of chaos in the double pendulum[J]. Chaos, Solitons & Fractals, 2006, 29(2): 417-422.

[2]Nikitina N V. Ultimate energy of a double pendulum undergoing quasiperiodic oscillations[J]. International Applied Mechanics, 2007, 43(9): 1035-1042.

[3]司麗榮, 張競夫. 彈簧擺內(nèi)共振現(xiàn)象的實(shí)驗(yàn)研究[J]. 物理實(shí)驗(yàn), 2002, 22(3): 9-12.

Si Lirong, Zhang Jingfu. Study on the experiment of the autoparametric resonance of a spring fendulum[J]. Physics Experimentation, 2002, 22(3): 9-12.

[4]Lee W K, Park H D. Chaotic dynamics of a harmonically excited springpendulum system with internal resonance[J]. Nonlinear Dynamics, 1997, 14(3): 211-229.

[5]Zheng J L, Yu X W. Analysis of the autoparametric resonance of a spring pendulum[J]. Physics and Engineering, 2010, 20(2): 13-16.

[6]Van der Weele J P, De Kleine E. The Order-chaos-order sequence in the spring pendulum[J]. Physica A: Statistical Mechanics and its Applications, 1996, 228(1-4): 245-272.

[7]Li Y S, Shu X F. Internal resonance and chaotic motion of a spring pendulum[J]. Journal of Taiyuan University of Technology, 1998, 29(6): 555-559.

[8]于洪潔, 張靖姝, 洪嘉振. 剛?cè)狁詈蠌椈蓴[的復(fù)雜動力學(xué)行為[J]. 科技導(dǎo)報, 2013, 31(18): 32-38.

Yu Hongjie, Zhang Jingshu, Hong Jiazhen. Complex dynamical behavior of rigid-flexible coupling spring pendulum[J]. Science and Technology Review, 2013, 31(18): 32-38.

[9]鐘萬勰. 結(jié)構(gòu)動力方程的精細(xì)時程積分法[J]. 大連理工大學(xué)學(xué)報, 1994, 34(2): 131-136.Zhong Wanxie. On precise time-integration method for structural dynamics[J]. Journal of Dalian University of Technology, 1994, 34(2): 131-136.

[10] Zhong W X, Williams F W. A precise time step integration method[J]. Journal of Mechanical Engineering Science, 1994, 208: 427-430.

[11] 張洵安, 楊人風(fēng). 結(jié)構(gòu)混沌運(yùn)動的非線性精細(xì)積分計算方法[J]. 西安公路交通大學(xué)學(xué)報, 2000, 20(2): 99-101.

Zhang Xunan, Yang Renfeng. The method of nonlinear precise integration for chaotic motion of structure[J]. Journal of Xi’an Highway University, 2000, 20(2): 99-101.

[12] Wang M F, Zhou X Y. Renewal frecise time step integration method of structural dynamic analysis[J]. Acta Mechanica Sinica, 2004, 36(2): 191-195.

[13] Gu Y X, Chen B S, Zhang H W. Precise time-integration with dimension expanding method[J]. Acta Mechanica Sinica, 2000, 32(4): 447-456.

[14] Zhang S Y, Deng Z C. Increment-dimensional precise integration method for nonlinear dynamic equation[J]. Chinese Journal of Computational Mechanics, 2003, 20(4): 423-426.

[15] Ge G, Wang H L, Tan G. Improved increment-dimensional precise integration method for the nonlinear dynamic equation with multi-degree-of-freedom[J]. Journal of Tianjin University, 2009, 42(2): 114-117.

[16] 呂和祥, 于洪潔, 裘春航. 精細(xì)積分的非線性動力學(xué)積分方程及其解法[J]. 固體力學(xué)學(xué)報, 2001, 22(3): 303-308.

Lv Hexiang, Yu Hongjie, Qiu Chunhang. A integral equation of non-linear dynamics and Its solution method[J]. Acta Mechanica Solida Sinaca, 2001, 22(3): 303-308.

(責(zé)任編輯耿金花)

On the Chaotic Dynamic Behaviour of the Rigid Rod-Spring Pendulum Model

ZHAO Cong，YU Hongjie

(Department of Engineering Mechanics, Shanghai Jiaotong University, Shanghai 200240, China)

A planar rigid rod-spring pendulum model was constructed and dimensionless dynamic equation was given. We numerically simulated and analyzed the dynamical behavior of the two-time-scale system while the frequency ratio and length ratio and initial conditions vary. The dynamic equation is strongly nonlinear as the fast and slow variables couple each other. A cubic interpolation precise integration method was applied to solve nonlinear dynamic equation. We also employ the Poincaré maps and the maximum Lyapunov exponent methods. Numerical simulation results demonstrate that the system presents complex chaotic motion in different parameters conditions. It is also found that the fast variable may transform to chaos via the quasi-periodic torus breakdown.

rigid-flexible coupling system; two-time-scale; spring pendulum; interpolation precise integration method; chaos

1672-3813(2016)03-0097-06；DOI:10.13306/j.1672-3813.2016.03.013

2014-09-16；

2015-01-04

國家自然科學(xué)基金重點(diǎn)項(xiàng)目(11132007)

趙聰(1990-)，男，遼寧鞍山人，碩士研究生，主要研究方向?yàn)榉蔷€性動力學(xué)。

于洪潔(1968-)，女，遼寧大連人，副教授，博士研究生，主要研究方向?yàn)閺?fù)雜系統(tǒng)非線性動力學(xué)，混沌同步與控制。

O313， O322

A