王杰方，安偉光，宋向華

（哈爾濱工程大學(xué)航天工程系，哈爾濱 150001）

一種改進(jìn)Mathieu方程動(dòng)力不穩(wěn)定邊界的方法

王杰方，安偉光，宋向華

（哈爾濱工程大學(xué)航天工程系，哈爾濱 150001）

基于Mathieu方程的臨界頻率方程式，提出了一種改進(jìn)Mathieu方程不穩(wěn)定邊界的方法，并獲得了比Bolotin近似邊界更精確的前三階收斂的不穩(wěn)定邊界。從改進(jìn)的不穩(wěn)定區(qū)域邊界表達(dá)式和Bolotin近似公式的對比中發(fā)現(xiàn)：兩種方法獲得的第一、二階不穩(wěn)定區(qū)域相差不大，但相較于Bolotin的第三階不穩(wěn)定區(qū)域，改進(jìn)的第三階不穩(wěn)定區(qū)域整體上移，且上移幅度隨著激發(fā)系數(shù)的增大而增大。當(dāng)激發(fā)系數(shù)μ取0．5時(shí)，上邊界上移幅度為8．61%，下邊界上移幅度為11．56%。對于受低頻載荷作用的動(dòng)力穩(wěn)定性問題，第三階不穩(wěn)定邊界公式的改進(jìn)具有重要的意義。

動(dòng)力穩(wěn)定性；Mathieu方程；臨界頻率；動(dòng)力不穩(wěn)定區(qū)域

Mathieu方程式常見于物理和工程的各個(gè)領(lǐng)域中：趙晶瑞［1］對一種新型的深海順應(yīng)式采油平臺（Spar平臺）的縱搖運(yùn)動(dòng)進(jìn)行了分析，得到了具有三次非線性項(xiàng)的有阻尼Mathieu方程。另外，在載流薄板的磁彈性動(dòng)力屈曲問題中也存在Mathieu方程［2］。Li等［3］利用Mathieu方程分析了諧振式慣性傳感器的動(dòng)力學(xué)特性。Mathieu方程有一個(gè)最重要的性質(zhì)：當(dāng)它們的系數(shù)間存在某種關(guān)系時(shí)，方程式具有無限增長的解，這些解在對應(yīng)的參數(shù)平面上完全布滿了許多區(qū)域，這些區(qū)域相當(dāng)于不穩(wěn)定區(qū)域［4］。因此，確定不穩(wěn)定區(qū)域是研究Mathieu方程問題的重點(diǎn)。

現(xiàn)有文獻(xiàn)中關(guān)于不穩(wěn)定區(qū)域的分析方法主要有攝動(dòng)法和傅里葉分析法［5］：攝動(dòng)法［6］適用于小參數(shù)的情況，傅里葉分析法［7－8］適用于大參數(shù)的情況，例如，在實(shí)際的海洋環(huán)境下，Mathieu方程中的參數(shù)不再是小參數(shù)［5，9］，攝動(dòng)法不再適用，需要采用傅里葉分析法。由于不穩(wěn)定區(qū)在任何工程或者物理問題中都屬于必須避開的禁區(qū)，因此其精確性是十分重要的。文中將對目前廣泛采用的傅里葉分析法進(jìn)行改進(jìn)，推導(dǎo)出比文獻(xiàn)［4］中的Bolotin近似公式更精確的不穩(wěn)定邊界公式。

本文思路是：將周期解按傅里葉級數(shù)展開并代入Mathieu方程中得到臨界頻率方程式，依次從臨界頻率方程式中提取n階前主子行列式等式，然后將其轉(zhuǎn)化為一元n次方程（n≤4），通過多次直接求解一元n次方程的根來獲得收斂的不穩(wěn)定邊界表達(dá)式。

1 結(jié)構(gòu)動(dòng)力穩(wěn)定性問題

1.1 結(jié)構(gòu)動(dòng)力穩(wěn)定性問題的Mathieu方程

以簡諧縱向力P＝P0（1+βcos（θt））作用下的鉸支直桿的動(dòng)力穩(wěn)定性問題為例［4］，其Mathieu方程為

從靜力屈曲角度：縱向力的最大值Pmax＝（1+ β）P0須小于臨界屈曲載荷p*1，即P0須小于P*1／（1+ β）。假設(shè)P0＝P*1／（1+β），代入μk中得，β取任意值時(shí)，激發(fā)系數(shù)μ恒等于0．5。由于μk是隨P0單調(diào)遞增的函數(shù)，所以激發(fā)系數(shù)μk的取值范圍為（0，0．5）。

Mathieu方程存在周期為T和2T的周期解，且周期相同的兩個(gè)解包圍著不穩(wěn)定區(qū)，周期不同的兩個(gè)解包圍著穩(wěn)定區(qū)。因此，找出使Mathieu方程具有周期T和2T的周期解的條件就可以確定不穩(wěn)定邊界［4］。

1.2 與周期為2T的周期解對應(yīng)的臨界頻率方程式

將Mathieu方程的解展開為

并代入Mathieu方程中，得到與周期2T的周期解對應(yīng)的臨界頻率方程式為

周期2T的周期解包圍的不穩(wěn)定區(qū)域由式（3）確定。

1.3 與周期為T的周期解對應(yīng)的臨界頻率方程式

同樣，將Mathieu方程的解展開為

得到與周期T的周期解對應(yīng)的臨界頻率方程式

周期T的周期解包圍的不穩(wěn)定區(qū)由式（5），式（6）確定。

2 改進(jìn)Mathieu方程不穩(wěn)定邊界的方法

2．1 與周期2T的周期解對應(yīng)的動(dòng)力不穩(wěn)定區(qū)域邊界

（1）若式（2）中的系數(shù)a1，b1存在非零解，其他系數(shù)為0，對應(yīng)于式（3）的左上角的對角線元素為0，從而獲得第一階動(dòng)力不穩(wěn)定區(qū)域的邊界表達(dá)式

式中：u表示上邊界，d表示下邊界，θ*iju表示第j次計(jì)算的第i階上邊界，j表示第j次計(jì)算的不穩(wěn)定邊界，下同。上式中的第一階動(dòng)力不穩(wěn)定區(qū)域邊界表達(dá)式是關(guān)于Mathieu方程存在周期為2T的周期解嚴(yán)格成立的條件，但它們并不是完整的第一階不穩(wěn)定區(qū)域的邊界。下面將會對第一階不穩(wěn)定區(qū)域的邊界表達(dá)式進(jìn)行進(jìn)一步的推導(dǎo)計(jì)算和修正，以獲得完整的、收斂的第一階不穩(wěn)定區(qū)域邊界。

（2）若式（2）中的系數(shù)a1，b1，a3，b3存在非零解，其他系數(shù)為0，對應(yīng)于式（3）中2階前主子行列式等于0，得到第二個(gè)第一階動(dòng)力不穩(wěn)定區(qū)域的邊界

圖1 第一次修正后的第一階不穩(wěn)定區(qū)域Fig．1 The first order unstable region after the first amend

從圖1可知，上邊界θ*12u在 θ*11u之上，下邊界θ*12d包含在θ*11d之內(nèi)，所以經(jīng)第一修正后的第一階不穩(wěn)定區(qū)域的邊界為 (θ*12u，θ*11d)。

（3）若式（2）中的系數(shù)a1，a3，a5和b1，b3，b5，存在非零解，其他系數(shù)為0，對應(yīng)于式（3）中3階前主子行列式等于0，將它們化為兩個(gè)標(biāo)準(zhǔn)的一元三次方程，

式中：x1＝x2＝θ2／4Ω2。

采用盛金公式求解上述兩個(gè)一元三次方程的根，即可獲得不穩(wěn)定區(qū)域的邊界公式。盛金判別式為

第三個(gè)第一階動(dòng)力不穩(wěn)定區(qū)域的邊界為

圖2給出了將不穩(wěn)定區(qū)（θ*32u，θ*32d）與（θ*31u，θ*31d）整合后得到的第一次修正的第三階動(dòng)力不穩(wěn)定區(qū)域。

圖2 第一修正后的第三階不穩(wěn)定區(qū)域Fig．2 The third order unstable region after the first amend

從圖2可知，第三階不穩(wěn)定區(qū)域的上邊界θ*32u在θ*31u之上，下邊界θ*32d包含在下邊界θ*31d之內(nèi)，所以，經(jīng)第一次修正后的第三階不穩(wěn)定區(qū)域的邊界為（θ*32u， θ*31d）。另外，經(jīng)計(jì)算發(fā)現(xiàn)θ*13u和θ*13d的值與上文中的θ*12u和θ*12d非常相近，這說明第一次修正后得到的邊界（θ*12u，θ*11d）是已經(jīng)收斂的邊界，因此并不需要對第一階不穩(wěn)定區(qū)進(jìn)行第二次修正。收斂的第一階不穩(wěn)定區(qū)域的邊界表達(dá)式為

（4）若式（2）中的系數(shù)a1，a3，a5，a7和b1，b3，b5，b7存在非零解，其他系數(shù)為0，即式（3）中4階前主子行列式等于0，將這兩個(gè)4階行列式等式化為標(biāo)準(zhǔn)的一元四次方程，

式中：x3＝x4＝θ2／4Ω2。

采用置換群解法求解上述兩個(gè)一元四次方程的根即可獲得不穩(wěn)定區(qū)邊界。置換群法的相關(guān)判別式為

第三個(gè)第三階動(dòng)力不穩(wěn)定區(qū)域的邊界

當(dāng)μ∈［0，0．5］時(shí)，θ*14u和θ*14d的值與上文中的θ*13u和θ*13d幾乎一致，這也再一次驗(yàn)證了式（14）是收斂的。當(dāng)μ∈［0，0．5］時(shí)，θ*33u的值略高于θ*32u（激發(fā)系數(shù)取0．5時(shí)，變化幅度僅為0．3%），因此第三階動(dòng)力不穩(wěn)定區(qū)域的邊界仍為（θ*32u，θ*31d）。最終得到收斂的第三階動(dòng)力不穩(wěn)定區(qū)域的上下邊界表達(dá)式為

2．2 與周期T的周期解對應(yīng)的動(dòng)力不穩(wěn)定區(qū)域邊界

（1）式（5）的左上角的對角線元素為0，式（6）中2階前主子行列式為0，得

（2）式（5）中的2階前主子行列式為0，式（6）中

3階前主子行列式為0，得到第二個(gè)第二階動(dòng)力不穩(wěn)定區(qū)域的邊界

將（θ*22u，θ*22d）與（θ*21u，θ*21d）整合起來得到第一次修正后的第二階動(dòng)力不穩(wěn)定區(qū)域見圖3。

圖3 第一次修正后的第二階不穩(wěn)定區(qū)域Fig．3 The second order unstable region after the first amend

從圖3可知，第二階不穩(wěn)定區(qū)域的上邊界θ*22u在θ*21u之上，下邊界θ*22d包含在下邊界θ*21d之內(nèi)，所以經(jīng)第一修正后的第二階不穩(wěn)定區(qū)域的邊界為（θ*22u，θ*21d）。

（3）式（5）中的3階前主子行列式為0，式（6）中4階前主子行列式為0，獲得兩個(gè)標(biāo)準(zhǔn)的一元三次方程，

式中：x5＝x6＝θ2／Ω2。

采用盛金公式求解上述方程的根，能夠獲得第三個(gè)第二階動(dòng)力不穩(wěn)定區(qū)域的邊界

式中：θi，Ai的含義與式（12）和式（13）相同。

計(jì)算表明，當(dāng)μ∈［0，0．5］時(shí)，θ*23u和θ*23d的所夾的不穩(wěn)定區(qū)域被（θ*22u，θ*21d）包裹，因此不需要對第二階動(dòng)力不穩(wěn)定區(qū)邊界進(jìn)行第二次修正。最終得到收斂的第二階動(dòng)力不穩(wěn)定區(qū)域的上下邊界表達(dá)式是

（4）式（5）中的4階前主子行列式為0，式（6）中5階前主子行列式為0，獲得兩個(gè)標(biāo)準(zhǔn)的一元四次方程，

式中：x7＝x8＝θ2／Ω2。

采用置換群法求解上述方程的根，即可獲得第四個(gè)第二階動(dòng)力不穩(wěn)定區(qū)域的邊界

式中：yi1，yi2，yi3，θi的含義與式（17）中相同。

計(jì)算表明，在μ∈[0，0．]5范圍內(nèi)，θ*24u和θ*24d的值與上文中的θ*23u和θ*23d幾乎一致，這也再一次驗(yàn)證了式（23）的收斂性。

3 改進(jìn)的不穩(wěn)定邊界與Bolotin近似邊界比較

為了獲得收斂的前三階不穩(wěn)定邊界，需要從臨界頻率方程式中提取的最高次方程為一元四次方程。經(jīng)驗(yàn)算，文中采用盛金公式求解一元三次方程的根和采用置換群法求解一元四次方程的根時(shí)，其誤差的數(shù)量級僅為10－15，完全可以認(rèn)為每一步的求解過程并不產(chǎn)生誤差。不穩(wěn)定邊界的誤差只與最終收斂情況有關(guān)，但邊界的收斂卻是可控的。文中的前三階不穩(wěn)定邊界（式（14）、式（23）和式（18））經(jīng)過多次修正已是收斂的邊界，因此采用本方法能獲得比Bolotin近似邊界更精確的表達(dá)式。

3.1 改進(jìn)的不穩(wěn)定邊界與Bolotin近似邊界的比較

將改進(jìn)的前三階不穩(wěn)定區(qū)域邊界表達(dá)式（14）、式（23）和式（18）與Bolotin前三階不穩(wěn)定區(qū)域近似邊界公式［4］進(jìn)行對比。圖4～圖6中分別給出了第一階、第二階、第三階改進(jìn)的不穩(wěn)定邊界和Bolotin近似邊界的曲線圖，實(shí)線表示改進(jìn)的不穩(wěn)定邊界，虛線表示Bolotin近似邊界。注意，由于本改進(jìn)的第一階和第二階不穩(wěn)定區(qū)域的下邊界表達(dá)式與Bolotin方法給出的邊界公式相同，因此，圖4和圖5中的下邊界只有一條實(shí)曲線。

圖4 第一階不穩(wěn)定區(qū)域的比較Fig．4 Comparison of the first order dynamic unstable region

圖5 第二階不穩(wěn)定區(qū)域的比較Fig．5 Comparison of the second order dynamic unstable region

圖6 第三階不穩(wěn)定區(qū)域的比較Fig．6 Comparison of the third order dynamic unstable region

從圖4和圖5中可知，對于第一階、第二階不穩(wěn)定區(qū)域的上邊界，兩種計(jì)算方法得出的結(jié)果相差都不大。計(jì)算表明，當(dāng)μ＝0．5時(shí)，兩種方法計(jì)算的第一階上邊界誤差僅為0．653%，第二階上邊界誤差也僅為－0．193%。這說明采用本方法獲得的第一階和第二階不穩(wěn)定區(qū)邊界表達(dá)式的精確度是合理的。

從圖6可知，相較于第三階Bolotin不穩(wěn)定區(qū)域，改進(jìn)的第三階不穩(wěn)定區(qū)域整體上移，且上移幅度隨著激發(fā)系數(shù)的增大而增大。當(dāng)μ分別取0．1、0．2、0．3、0．4、0．5時(shí)，改進(jìn)的第三階不穩(wěn)定區(qū)域的上邊界較Bolotin的第三階上邊界的上移幅度分別為1．003%、1．223%、2．492%、5．128%、8．609%；當(dāng)μ分別取0．1、0．2、0．3、0．4、0．5時(shí)，改進(jìn)第三階不穩(wěn)定區(qū)域的下邊界的上移幅度分別為0．0%、0．0%、0．326%、1．832%、11．558%，偏差較大，不能忽略。

3.2 第三階動(dòng)力不穩(wěn)定區(qū)域邊界改進(jìn)的意義

以文獻(xiàn)［10］中的水下超空泡高速航行器為例：超空泡尾部脈動(dòng)和尾渦脫落產(chǎn)生的空泡高頻特征頻率θ的取值范圍為80～160 Hz，航行器殼體的特征頻率ω＝149 Hz，分析航行器的動(dòng)力不穩(wěn)定區(qū)域與外載荷頻率之間的位置關(guān)系。

對于水下高速航行器的動(dòng)力屈曲問題而言，其第一階不穩(wěn)定區(qū)域位于2Ω附近，第二階不穩(wěn)定區(qū)域位于Ω附近，第三階不穩(wěn)定區(qū)域位于2Ω／3附近，用定值分量P0和臨界力P*表示為

其中，臨界力只與結(jié)構(gòu)參數(shù)有關(guān)，與周期力無關(guān)。α1，α2，α3依次為第一階、第二階、第三階不穩(wěn)定區(qū)域邊界的系數(shù)，且α1＝298，α2＝149，α3＝99．3。下面的第一階、第二階、第三階不穩(wěn)定區(qū)域上下邊界的系數(shù)依次用（α1u，α1d）、（α2u，α2d）、（α3u，α3d）來表示。由于μ＝0時(shí)，不穩(wěn)定區(qū)域的上下邊界相等，所以上式的系數(shù)省略了下角標(biāo)u和d。

表1給出了激發(fā)系數(shù)μ取不同值時(shí)，改進(jìn)的前三階不穩(wěn)定區(qū)域上下邊界系數(shù)（α1u，α1d），（α2u，α2d），（α3u，α3d）的值和第三階Bolotin不穩(wěn)定邊界系數(shù)（α3u′，α3d′）。

表1 前三階不穩(wěn)定區(qū)域邊界的系數(shù)Tab.1 The coefficient of the first three order unstable region

以μ＝0．5時(shí)的不穩(wěn)定區(qū)為例，結(jié)合表1中的數(shù)據(jù)，經(jīng)計(jì)算發(fā)現(xiàn)：

從以上四個(gè)結(jié)論可知：① 當(dāng)P0／P*的值較大時(shí)，第一二階不穩(wěn)定區(qū)與載荷頻率發(fā)生重合的可能性較大，是主要的不穩(wěn)定區(qū)。② 當(dāng)P0／P*的值較小時(shí)，第三階不穩(wěn)定區(qū)域與外載荷頻率區(qū)間發(fā)生重合的可能性較大，也可能是主要的不穩(wěn)定區(qū)域。且從“（3）”和“（4）”的對比中可知：采用改進(jìn)的不穩(wěn)定邊界表達(dá)式，P0／P*＞0．33時(shí)才能避開第三階不穩(wěn)定區(qū)域，而采用Bolotin的近似邊界表達(dá)式時(shí)，P0／P*＞0．21時(shí)就能避開第三階不穩(wěn)定區(qū)域。也就是說，采用Bolotin的近似邊界表達(dá)式認(rèn)為(0．21，0．33)是P0／P*的安全取值范圍，但采用改進(jìn)的第三階不穩(wěn)定表達(dá)式認(rèn)為(0．21，0．33)是P0／P*的危險(xiǎn)取值范圍。這說明采用精度更高的改進(jìn)不穩(wěn)定邊界使結(jié)構(gòu)避免了采用Bolotin近似公式時(shí)所隱藏的危險(xiǎn)。

4 超空泡運(yùn)動(dòng)體圓柱薄殼艙段的動(dòng)力穩(wěn)定性分析算例

超空泡運(yùn)動(dòng)體航行深度H＝10 m，流場密度ρw＝1 000 kg／m3，對于自然超空泡，空泡內(nèi)的飽和蒸汽壓pc＝2 350 Pa（20°C），標(biāo)準(zhǔn)大氣壓p＝101 325 Pa。圓柱薄殼艙段的幾何參數(shù)：半徑R＝0．2 m，長度L＝4 m，厚度h＝3 mm，空化器直徑dn＝0．2 m。材料物理參數(shù)：彈性模量E＝209 GPa，材料密度ρ＝7 850 kg／m3，泊松比ν＝0．3。

基于文獻(xiàn)［11］中圓柱薄殼動(dòng)力穩(wěn)定性微分方程，采用Bolotin方法和本方法計(jì)算超空泡運(yùn)動(dòng)體圓柱艙段的動(dòng)力不穩(wěn)定區(qū)。圖7～圖9給出了i＝1，k＝1（i，k分別為軸向和周向的半波數(shù)）時(shí)，依賴于航行速度V（m／s）和頻率f／Hz的前三階動(dòng)力不穩(wěn)定區(qū)域。

圖7 第一階不穩(wěn)定區(qū)域的比較Fig．7 Comparison of the first order dynamic unstable region

圖8 第二階不穩(wěn)定區(qū)域的比較Fig．8 Comparison of the second order dynamic unstable region

圖9 第三階不穩(wěn)定區(qū)域的比較Fig．9 Comparison of the third order dynamic unstable region

從圖7～圖9可知，隨著航行速度的增大，超空泡運(yùn)動(dòng)體圓柱薄殼艙段的不穩(wěn)定區(qū)間寬度也在不斷地增大。由圖7和圖8可知，對于第一階和第二階不穩(wěn)定區(qū)域，采用Bolotin方法和本方法獲得的結(jié)果相差甚小，可以忽略；但是，從圖9可知，采用本文方法獲得的不穩(wěn)定區(qū)較Bolotin方法獲得的不穩(wěn)定區(qū)整體上移，且航行速度越大上移幅度越大，這一結(jié)論與3．1節(jié)中的結(jié)論是一致的。

5 結(jié) 論

本文是采用多次直接求解一元n次方程的根來獲得不穩(wěn)定邊界表達(dá)式的。首先，利用盛金公式求解一元三次方程的根和置換群法求解一元四次方程的根時(shí)，其誤差的數(shù)量級僅為10－15，完全可以認(rèn)為每一步的求解并不產(chǎn)生誤差；其次，文中的前三階不穩(wěn)定邊界經(jīng)過多次修正已是收斂的邊界。以上兩點(diǎn)確保了本文改進(jìn)方法的精確度。并得出以下結(jié)論：

（1）從改進(jìn)的不穩(wěn)定區(qū)域邊界表達(dá)式和Bolotin近似公式的對比中發(fā)現(xiàn)：兩種方法獲得的第一階、第二階動(dòng)力不穩(wěn)定區(qū)域相差不大；相較于第三階Bolotin不穩(wěn)定區(qū)域，本文方法獲得的改進(jìn)的第三階不穩(wěn)定區(qū)域整體上移，且上移幅度隨著激發(fā)系數(shù)的增大而增大，當(dāng)μ＝0．5時(shí)，上邊界上移幅度為8．61%，下邊界上移幅度為11．56%。

（2）從水下超空泡高速航行器的動(dòng)力穩(wěn)定性分析算例中可知：當(dāng)P0／P*的值較小時(shí)，第三階不穩(wěn)定區(qū)域與外載荷頻率區(qū)間發(fā)生重合的可能性較大，也可能是主要的不穩(wěn)定區(qū)域。

（3）采用精度更高的改進(jìn)的不穩(wěn)定邊界使結(jié)構(gòu)避免了采用Bolotin不穩(wěn)定近似公式時(shí)所隱藏的危險(xiǎn)。

綜上所述，由改進(jìn)的傅里葉分析法確定的不穩(wěn)定區(qū)邊界的精度，為Mathieu方程的動(dòng)力不穩(wěn)定區(qū)提供了比Bolotin不穩(wěn)定區(qū)域更精確的邊界表達(dá)式。

［1］趙晶瑞，唐友剛，王文杰．傳統(tǒng)Spar平臺參數(shù)激勵(lì)Mathieu不穩(wěn)定性的研究［J］．工程力學(xué)，2010，27（3）：222－227．

ZHAO Jing-rui，TANG You-gang，WANGWen-jie．Study on the parametrically excited mathieu instability of a classic spar platform［J］．Engineering Mechanics，2010，27（3）：222－227．

［2］王平，王知人，白象忠．馬丟方程解的穩(wěn)定性在磁彈性屈曲中的應(yīng)用［J］．哈爾濱工業(yè)大學(xué)學(xué)報(bào)，2009，41（11）：222－224．

WANG Ping，WANG Zhi-ren，BAI Xiang-zhong．Applications of the stability of Mathieu equation's solution in magnetic-elastisity buckling problem［J］．Journal of Harbin Institute of Technology，2009，41（11）：222－224．

［3］Li Yan，F(xiàn)an Shang-chun，Guo Zhan-she．Mathieu equation with application to analysis of dynamic characteristics of resonant inertial sensors［J］．Commun Nonlinear Sci Numer Simulat，2013，18：401－410．

［4］Bolotin，V V．The Dynamic Stability of Elstic Systems［M］．Holen-Day San Francisco，1964．

［5］徐萬海，吳應(yīng)湘，鐘興福，等．海洋細(xì)長結(jié)構(gòu)參數(shù)激勵(lì)不穩(wěn)定區(qū)的確定方法［J］．振動(dòng)與沖擊，2011，30（9）：79－83．

XU Wan-hai，WU Ying-xiang，ZHONG Xing-fu，et al．Methods for parametric excitation instability analysis of slender flexible cylindrical structures in offshore engineering［J］．Journal of Vibration and Shock，2011，30（9）：79－83．

［6］Nayfeh A H．Introduction to Perturbation Technique［M］．New York：Wiley，1981．

［7］馮世寧，陳浩然．含分層損傷復(fù)合材料層合板非線性動(dòng)力穩(wěn)定性［J］．復(fù)合材料學(xué)報(bào)，2006，23（1）：154－160．

FENG Shi-ning，CHEN Hao-ran．Nonlinear dynamic instability behavior of delaminated composite laminates［J］．Acta Materiae Compositae Sinica，2006，23（1）：154－160．

［8］趙洪金，董寧娟，劉超，等．基于能量法的高溫（火災(zāi)）環(huán)境下軸心受壓格構(gòu)柱動(dòng)力穩(wěn)定性分析［J］．工 程 力 學(xué)，2011，28（12）：160－165．

ZHAO Hong-jin，DONG Ning-juan，LIU Chao，et al．Dynamic stability analysis on axial compression lattice column under high temperature condition using energy method［J］．Engineering Mechanics，2011，28（12）：160－165．

［9］王俊榮，謝彬．深水半潛式平臺Mathieu不穩(wěn)定問題研究［J］．工程力學(xué)，2012，29（10）：347－353．

WANG Jun-rong，XIE Bin．Mathieu instablity study of a deepwater semi-submersible platform［J］．Engineering Mechanics，2012，29（10）：347－353．

［10］黨建軍，張謀，郭芳，等．超空泡脈動(dòng)對水下高速航行器殼體強(qiáng)度的影響研究［J］．機(jī)械科學(xué)與技術(shù)，2011，30（11）：1930－1933．

DANG Jian-jun，ZHANG Mou，GUO Fang，et al．A study of the influence of the super-cavit pulsation on the intensity of the shell of the high-speed underwater vehicle［J］．Mechanical Science and Technology for Aerospace Engineering，2011，30（11）：1930－1933．

［11］王杰方，安偉光，宋向華．超空泡運(yùn)動(dòng)體圓柱薄殼動(dòng)力屈曲及可靠性分析［J］．振動(dòng)與沖擊，2014，33（8）：22－28．

WANG Jie-fang，AN Wei-guang，SONG Xiang-hua．Dynamic buckling and reliability analysis of a cylindrical thin shell for supercavitating vehicles［J］．Journal of Vibration and Shock，2014，33（8）：22－28．

Im proved method about dynam ic unstable boundary ofmathieu equation

WANG Jie-fang，ANWei-guang，SONG Xiang-hua
（Department of Aerospace Engineering，Harbin Engineering University，Harbin 150001，China）

An improvedmethod about unstable boundary of Mathieu equation was proposed according to the critical frequency equation，and the first three orders convergent unstable boundary was got，which ismore accurate than the Bolotin approximate boundary．Comparing the twomethods，it shows that their first two orders dynamic unstable region are almost the same，the third order unstable region of the improved method moves upward compared with the Bolotinmethod，and the range is amplified with the growth of excitation coefficient．When the excitation coefficientμis 0．5，the upper boundarymoves upward about8．61%and the lower boundary moves upward about 11．56%．With respect the dynamic stability problem caused by low frequency loading，the improvement on the third order unstable boundary expression has great significance．

dynamic stability；mathieu equation；critical frequency；dynamic unstable region

O342

A

10．13465／j．cnki．jvs．2015．12．031

國家科技部國際合作專項(xiàng)（2012DFR00070）

2014－04－28 修改稿收到日期：2014－05－27第一作者王杰方女，博士，1987年生

安偉光 男，教授，博士生導(dǎo)師，1943年生

郵箱：anweiguang＠hrbeu．edu．cn