朱金文， 楊德慶

(上海交通大學(xué) 船舶海洋與建筑工程學(xué)院，上海 200240)

近年來，圓頻率函數(shù)依賴于坐標(biāo)的一次導(dǎo)數(shù)的非線性振動系統(tǒng)引起了學(xué)術(shù)界的廣泛關(guān)注。就單自由度系統(tǒng)而言，對應(yīng)的微分方程具有如下形式：

(1)

這類非線性振動微分方程奠基性的研究始于Mickens等[2]，他們考慮了式(1)兩個(gè)特例的一般性質(zhì)。之后，Beatty等[3]研究了式(1)中具有周期解的最簡單的情況，即:

(2)

針對式(2)Mickens等[4-7]，做了一些后續(xù)的擴(kuò)展性工作。應(yīng)當(dāng)注意到式(2)是描述具有靜態(tài)和慣性三次非線性的機(jī)械振子的一個(gè)特例。一般形式上，它們的微分方程是:

(3)

另一方面，當(dāng)式(3)中參數(shù)取為ν1=0，且初始條件與式(2)相同時(shí)，式(3)退化為:

(4)

其中：x是系統(tǒng)的位移，ω0是該系統(tǒng)線性情況下的固有頻率，和ν是非線性系統(tǒng)參數(shù)。該方程即在式(2)的基礎(chǔ)上考慮了三次恢復(fù)力項(xiàng)。

第一，在計(jì)算式(4)周期解的一階近似時(shí)，會遇到與式(2)相同的困難[4],即對于方程式(4)，傳統(tǒng)的攝動法或諧波平衡法得到頻率ω(A)=2π/T(A)在有限的幅值區(qū)間內(nèi)是奇異的或者說沒有定義。以諧波平衡法為例，為了求解式(4)的頻率-幅值關(guān)系，假設(shè)式(4)有如下形式的解：

x(t)=acosωt+bsinωt=Acos(ωt+θ)

(5)

將式(5)代入式(4)并令第一項(xiàng)諧波系數(shù)為0，得

(6)

第二，當(dāng)幅值A(chǔ)較大的時(shí)候，式(4)和式(2)的解均出現(xiàn)方波現(xiàn)象，文獻(xiàn)[6]描述了這種方波現(xiàn)象，但沒有探討原因。

然而，由于式(4)較之式(2)考慮了三次恢復(fù)力項(xiàng)，式(4)的復(fù)雜程度超過式(2)。例如，文獻(xiàn)[4]和[6]分別得到了式(2)的周期表達(dá)式，而式(4)的周期至今沒有文獻(xiàn)討論。本文將首先尋找式(4)的首次積分，進(jìn)而求得其周期表達(dá)式并研究其性質(zhì)，籍此探討幅值A(chǔ)較大時(shí)解方波現(xiàn)象及產(chǎn)生的原因。

本文安排如下：首先尋找式(4)的首次積分，得到振子的周期積分表達(dá)式，基于該積分表達(dá)式的可積性條件采用諧波平衡法求解式(4)一類初始條件下的精確解。然后研究式(4)的周期性質(zhì)，計(jì)算周期T(A)的近似解析表達(dá)式，解釋方波現(xiàn)象。最后，使用Hermite插值法求解式(4)，將解與數(shù)值解進(jìn)行對比以檢驗(yàn)該方法有效性。

1 方程的周期及一類初始條件下的精確解

(7)

易知式(7)的平衡點(diǎn)為(0,0)，該平衡點(diǎn)鄰域，有

(8)

式(8)中矩陣的特征值為λ=±ω0i，因此，式(7)的平衡點(diǎn)為中心[10]，故對于任意初始條件，式(4)的非平凡解均為周期解。

為了得到式(4)的首次積分，將其改寫為：

(9)

式(9)形式上等價(jià)于：

(10)

式(10)的解[11]為:

(11)

其中，c是任意積分常數(shù)。由于y2是x的偶函數(shù)，因此式(4)描述的振動必定是等幅值振動。另外，y代表該振子的運(yùn)動速度，因此周期為：

(12)

令式(11)中y2=0，并代入x=A，可求得c為

(13)

因此式(12)實(shí)際上是:

T(A)=

(14)

顯然，當(dāng)c=0時(shí)，可以得到T(A)的解析表達(dá)式

(15)

假設(shè)式(4)的精確解為:

x(t)=Acos(ωt+θ)

(16)

將式(16)代入式(4)，可以得到:

(17)

注意到，式(17)對任意時(shí)間t總是成立當(dāng)且僅當(dāng):

(18)

求解方程組(18)，得到:

(19)

與傳統(tǒng)諧波平衡取一項(xiàng)諧波不同，式(17)中取兩項(xiàng)諧波(即令兩項(xiàng)諧波系數(shù)均為零)，使得頻率與幅值均為定值，如式(19)所示。所以沒有出現(xiàn)前言中所述的奇異性現(xiàn)象。

因此，式(4)對一類初始條件，有如下精確解：

(20)

(21)

值得指出的是，式(20)可以用來檢驗(yàn)其他非線性近似方法在求解式(4)時(shí)的精度。

2 振子周期性質(zhì)及周期的近似解析表達(dá)式

A充分小時(shí)，有:

eνA2(1-u2)=1+νA2(1-u2)+O(A4)

(22)

將式(22)代入到式(14)并化簡，得到:

(23)

運(yùn)用以下不等式:

(24)

可以得到T(A)的下界:

(25)

為了考察T(A)是否具有上界，必須對其導(dǎo)數(shù)dT/dA進(jìn)行分析。由式(14)可以得到:

(26)

若式(26)中的積分變量換為w，那么有:

(27)

w隨著u在A2ν與0之間變化，即:

u=0?w=A2ν,u=1?w=0

(28)

將式(27)代入式(26)，得到:

(29)

為了公式的簡潔，式(26)與式(29)中積分函數(shù)分母中的一個(gè)u并未替換為w，這是因?yàn)閮H對積分函數(shù)中分子的符號感興趣。

應(yīng)用泰勒級數(shù)展開[12]，可以得到:

(30)

代入a,b，式(30)右端中的系數(shù)變?yōu)椋?/p>

(31)

由此

(32)

dT(A)/dA<0

(34)

因此T(A)無論如何都具有上界。需要指出的是，這里說A≥A*時(shí)，并不是說在0

由以上結(jié)果可知，T(A)為幅值A(chǔ)的有界函數(shù)，且最終隨著A的增大而衰減為0。

圖1 ω20，，和ν分別取0.1，5，和1時(shí)式(4)對應(yīng)的T(A)-A曲線

可以看到，T(A)隨著A的增長迅速的衰減到0。

可從另一個(gè)角度直觀的理解上述現(xiàn)象。由(11)及(13)可知，系統(tǒng)速度的平方為:

(35)

為了敘述的方便，令:

(36)

(37)

式中g(shù)(x)為開口向上的拋物線，f(x)為開口未定的鐘形線，如圖2所示。隨著A的增大，對應(yīng)的鐘形線由開口向上經(jīng)過水平位置變化到開口向下。其中水平位置，圖中的點(diǎn)劃線對應(yīng)著f(x)=0，亦即第1節(jié)中精確解的位置?？梢钥吹?，隨著A的增大，f(x)與g(x)之間的差值越大，與之對應(yīng)的平均速度迅速增大，而平均速度的增加量遠(yuǎn)遠(yuǎn)超過幅值A(chǔ)的增加量，故周期T(A)會隨著幅值A(chǔ)的增加而減小。

另外，由圖2可知，當(dāng)A取值較大時(shí)，由于幅值在A附近時(shí)振子的速度增長需要一定的時(shí)間，而0點(diǎn)周圍的振子的速度卻非常大，故振子表現(xiàn)出在±A處滯留，然后迅速穿過0點(diǎn)緊接著又在?A處滯留的特性，即由于ν的存在，振子的振動將表現(xiàn)出明顯的方波現(xiàn)象，如圖3所示，由于幅值A(chǔ)越大0點(diǎn)周圍的振子的速度就越大，則可以推斷當(dāng)A增大時(shí)，這種方波性質(zhì)將更明顯。此時(shí)，對于式(4)，數(shù)值解法可能會變得不可靠甚至無法求解，另一方面，要獲得式(4)的較高精度的近似解析解也將變的十分困難。

由式(23)易知，當(dāng)A充分小時(shí)，可以得到T(A)的一階近似表達(dá)式:

(38)

然而式(38)中并未包含ν，所以其精度必然較低，下面求解T(A)的二階近似表達(dá)式，由于：

(39)

將式(39)代入式(14)，可得:

(40)

積分式(40)可以用第一類完全橢圓積分K(q)表示出來，即:

T(A)=pK(q)

(41)

其中

(42)

圖，和ν分別取0.1，5，和1時(shí)方程(4)對應(yīng)的T(A)-A曲線

3 Hermite插值法求解方程

為了用Hermite插值法求解式(4)，將時(shí)間變量t轉(zhuǎn)變?yōu)棣?

τ=sinωt

(43)

該變換由Qaisi[13]提出冪級數(shù)法時(shí)引入，引入該變換，則當(dāng)t從t=0時(shí)開始增大到無窮大時(shí)，τ從τ=0開始以頻率ω在-1到+1之間振動，由此，無限的時(shí)間區(qū)域0≤t≤∞變?yōu)榱擞邢薜臅r(shí)間區(qū)域-1≤τ≤+1。相應(yīng)地，原來的方程及初始條件變?yōu)?

(44)

其中′表示對τ求導(dǎo)。則只需求得式(44)在τ∈[-1,1]之間的近似解，然后將τ替換為sinωt就可以得到式(4)的近似解。第一個(gè)步驟可以由Hermite插值法完成。

由于軌線上的初始狀態(tài)點(diǎn)由τ=0唯一確定，振子在“諧振時(shí)間τ”的半個(gè)周期內(nèi)必定已經(jīng)完成了一次振動，因此諧振時(shí)間τ的頻率ω必定是振子頻率Ω(A)的一半：

(45)

其中，T(A)由數(shù)值積分得到，且可以達(dá)到任意精度。

由于x(1)=-A，由式(44)可以得到:

(46)

為了得到式(44)在τ∈[-1,1]間的近似解，可以構(gòu)造如下的插值函數(shù):

(47)

其中ci(i=1,2)為待定常數(shù)。

進(jìn)一步，ci(i=1,2)由x(τ)在τ=1的函數(shù)值及導(dǎo)數(shù)值決定。令X(1)=x(1)，X′(1)=x′(1)，可以得到2個(gè)線性方程組成的線性方程組，求解該方程組得到ci(i=1,2):

(48)

只要將式(47)中的τ替換為sinωt，就得到了式(4)的近似解析解:

(49)

其中ci(i=1,2)由式(48)給出。

將(49)化簡得到:

X(t)=a1cosΩt+a3cos3Ωt

(50)

其中ai(i=1,3)是:

(51)

式(51)也可以寫成:

T=16.948 1,Ω=0.370 731

(53)

當(dāng)A=1.5時(shí)由數(shù)值積分得到振子的周期與頻率為:

T=1.525 83,Ω=4.117 88

(54)

圖5給出了A=0.1與A=1.5時(shí)近似解與數(shù)值解的對比，結(jié)果表明式(52)在A不很大時(shí)具有較高的精度。

當(dāng)振幅A較大時(shí)，如A=5時(shí)，此時(shí)

T=0.116 863,Ω=53.765 2

(55)

為了求解式(4)，可以取高階Hermite插值函數(shù)，如[10,5]階 Hermite插值函數(shù)解為:

X(t)=5-0.432 768sin2(Ωt/2)

-0.437 683sin4(Ωt/2)-0.595 948sin6(Ωt/2)

-0.910 523sin8(Ωt/2)-1.478 97sin10(Ωt/2)

-4 244.41sin12(Ωt/2)+18 236.8sin14(Ωt/2)

-31 945.5sin16(Ωt/2)+28 407.9sin18(Ωt/2)

-12 785.7sin20(Ωt/2)+2 324.77sin22(Ωt/2)

(56)

如圖6所示，可見精度令人滿意。

(a)A=0.1時(shí)的x(t)-t曲線 (b)A=1.5時(shí)的x(t)-t曲線

4 結(jié) 論

(57)

然后，基于周期的積分表達(dá)式研究了周期性質(zhì)，分析表明該振子的周期對于振幅的變化而言是有界的，且最終周期隨振幅的增大而衰減到0；得到了周期T(A)的二階近似解析表達(dá)式T(A)=pK(q)，其中p與q由式(42)給出，該近似表達(dá)式與數(shù)值積分的對比表明當(dāng)幅值A(chǔ)不很大時(shí)，具有較高的精度。最后提出一種Hermite插值法求解式(4)，與數(shù)值解的對比表明該方法是有效的；且當(dāng)幅值較大，即振子的振動表現(xiàn)出明顯的方波現(xiàn)象時(shí)，Hermite插值法仍可以使用，只要提高Hermite插值函數(shù)的階數(shù)，就可以達(dá)到所需的精度。

參 考 文 獻(xiàn)

[1]Hamdan M N， Shabaneh N H. On the large amplitude free vibrations of a restrained uniform beam carrying an intermediate lumped mass[J]. Journal of Sound and Vibration, 1997，199(5)：711-736.

[2]Mickens R E, Jackson D D. Oscillations in systems having velocity dependent frequencies[J]. Journal of Sound and Vibration, 1999, 223：329-333.

[8]Mickens R E. Oscillations in planar dynamic systems[M]. Singapore：World Scientific, 1996.

[9]Cveticanin L. Homotopy-perturbation method for pure nonlinear differential equation[J]. Chaos, Solitons & Fractals, 2006, 30(5)：1221-1230.

[10]Smale S, Hirsch M W, Devaney R L. Differential equations, dynamical systems, and an introduction to chaos[M]. Academic Press, 2003.

[11]Л.C.龐特里亞金.(林武忠、倪明康譯)常微分方程(第6版)[M]高等教育出版社, 2006.

[12]Courant R, John F. Introduction to calculus and analysis (Volume 1)[M]. Springer-Verlag, 1965.

[13]Qaisi M I. A power series approach for the study of periodic motion[J]. Journal of Sound and Vibration, 1996, 196(4):401-406.