麻岳敏，曹樹謙,3，郭虎倫

(1. 天津大學(xué) 機(jī)力學(xué)系，天津 300354； 2. 天津市非線性動(dòng)力學(xué)與控制重點(diǎn)實(shí)驗(yàn)室，天津 300354；3. 天津大學(xué) 力學(xué)國家級(jí)實(shí)驗(yàn)教學(xué)示范中心，天津 300354)

葉片作為航空發(fā)動(dòng)機(jī)、燃?xì)廨啓C(jī)等旋轉(zhuǎn)機(jī)械的重要部件之一，它的穩(wěn)定性直接關(guān)系到發(fā)動(dòng)機(jī)能否正常工作。它在高的離心荷載、氣動(dòng)荷載及振動(dòng)交變荷載等荷載的惡劣工作環(huán)境中，是極易發(fā)生故障的。已有研究表明，葉片故障占發(fā)動(dòng)機(jī)故障的62%以上，而葉片的損壞絕大多數(shù)是由于葉片振動(dòng)產(chǎn)生的動(dòng)應(yīng)力過大所致。因此，研究葉片的振動(dòng)問題具有重要的意義。

從已有文獻(xiàn)看，對(duì)旋轉(zhuǎn)葉片的研究主要集中在將不帶葉冠的葉片簡(jiǎn)化為旋轉(zhuǎn)懸臂梁來分析。Hodges等[1-2]建立了三自由度懸臂梁模型，研究了振動(dòng)頻率，分析了顫振速度與極限環(huán)響應(yīng)，并通過實(shí)驗(yàn)與理論結(jié)果進(jìn)行對(duì)比。Bercin等[3-4]考慮剪切和翹曲等因素對(duì)旋轉(zhuǎn)梁彎扭耦合振動(dòng)進(jìn)行了分析。張偉等[5-6]通過多尺度等方法研究了受活塞氣動(dòng)力預(yù)扭轉(zhuǎn)薄壁旋轉(zhuǎn)葉片的非線性動(dòng)力學(xué)問題，指出隨著氣流流速的增加，系統(tǒng)呈現(xiàn)周期運(yùn)動(dòng)和混沌運(yùn)動(dòng)等復(fù)雜動(dòng)力學(xué)行為。Eskandary等[7-8]用數(shù)值方法研究了Timoshenko旋轉(zhuǎn)梁的非線性振動(dòng)。Saxton等[9]建立了復(fù)合材料葉片動(dòng)力學(xué)方程，分析系統(tǒng)的動(dòng)力學(xué)行為。姜靜等[10]建立了考慮離心力且由于偏心發(fā)生彎扭耦合的非線性動(dòng)力學(xué)方程。Kashani等[11]針對(duì)Euler-Bernoulli梁分析了梁自由端受軸向載荷和力矩對(duì)系統(tǒng)彎曲扭轉(zhuǎn)振動(dòng)特性和穩(wěn)定性的影響。任勇生等[12]在不考慮葉片幾何和氣動(dòng)力非線性情況下，對(duì)風(fēng)力機(jī)葉片的彎扭耦合線性氣彈穩(wěn)定性進(jìn)行研究，并用特征值法求穩(wěn)定性。Hansen[13]對(duì)風(fēng)機(jī)葉片的兩大類氣動(dòng)彈性問題進(jìn)行了綜述，列出了在對(duì)葉片氣動(dòng)彈性穩(wěn)定性分析時(shí)所需要關(guān)注的重要問題。王丹等[14]建立彎扭耦合的葉片截面模型，指出葉片顫振過程中振動(dòng)幅值隨著來流速度的增加而增大。Hosseini[15]研究了旋轉(zhuǎn)葉片在空氣熱彈性載荷下的動(dòng)力學(xué)行為。Ladge等[16]通過實(shí)驗(yàn)研究了渦輪機(jī)葉片在動(dòng)態(tài)流作用下的穩(wěn)定性。Pourazarm等[17]考慮氣體作用下的變截面旋轉(zhuǎn)葉片的耦合振動(dòng)。隨著高速旋轉(zhuǎn)機(jī)械的發(fā)展，葉冠對(duì)葉片的振動(dòng)特性以及穩(wěn)定性的影響是我們不能忽視的一個(gè)問題。

李劍釗等[18]從結(jié)構(gòu)、特性分析和相關(guān)試驗(yàn)介紹了汽輪機(jī)帶冠葉片。周傳月等[19-20]通過有限元方法，研究了帶冠葉片振動(dòng)特性。Ansari等[21]研究了加質(zhì)量塊的旋轉(zhuǎn)梁的彎扭振動(dòng)頻率，得到質(zhì)量塊質(zhì)量的增加引起頻率的增加，而末端帶質(zhì)量塊梁長(zhǎng)度的增加降低了頻率的結(jié)論。目前，國內(nèi)外學(xué)者針對(duì)旋轉(zhuǎn)帶冠葉片受到氣動(dòng)載荷和離心載荷等多種載荷共同作用時(shí)的非線性氣動(dòng)彈性問題相對(duì)較少，有待進(jìn)一步研究。

本文主要研究旋轉(zhuǎn)帶冠葉片的非線性氣動(dòng)彈性問題?？紤]發(fā)生彎扭耦合變形，受離心力和氣動(dòng)力的影響，利用Hamilton原理建立柔性葉片的非線性偏微分方程，對(duì)其無量綱化，運(yùn)用Galerkin離散得到兩自由度耦合的常微分方程。通過諧波平衡法以及數(shù)值模擬研究帶冠旋轉(zhuǎn)葉片的動(dòng)力學(xué)行為。

1 帶冠葉片動(dòng)力學(xué)模型

將旋轉(zhuǎn)帶冠葉片簡(jiǎn)化為如圖1所示的無預(yù)扭轉(zhuǎn)角、軸向不可伸長(zhǎng)、末端加質(zhì)量塊的懸臂梁?？紤]葉片安裝在半徑為R，角速度Ω旋轉(zhuǎn)的軸上，受到Y(jié)方向的氣動(dòng)荷載。其中，梁的質(zhì)心軸與彈性軸距離為e，質(zhì)量塊質(zhì)心與梁末端的剛心重合。

為了簡(jiǎn)化建模過程，主要假設(shè)如下：

(1) 葉冠簡(jiǎn)化為集中質(zhì)量點(diǎn)，并忽略葉冠對(duì)系統(tǒng)結(jié)構(gòu)剛度的影響;

(2) 帶冠葉片為各向同性材料，本構(gòu)關(guān)系滿足胡克定律。

(1)

該點(diǎn)的速度可表示為

(2)

自由端質(zhì)量塊的速度表示為

(3)

系統(tǒng)的動(dòng)能為

(4)

式中：L為梁的長(zhǎng)度；ρ為梁密度；A為梁的橫截面面積；M為質(zhì)量塊質(zhì)量。

系統(tǒng)勢(shì)能由葉片變形和旋轉(zhuǎn)產(chǎn)生，表示為

(5)

式中：EI為抗彎剛度；GJ為抗扭剛度；κ(x)為曲率；Fc為葉片所受離心力；Fm為質(zhì)量塊所受離心力，其表達(dá)式分別為

(6)

Fm=M(R+L)Ω2

(7)

由材料力學(xué)可知，當(dāng)懸臂梁的變形較小時(shí)，曲率可近似表達(dá)為

(8)

由于葉片看作柔性梁，有相對(duì)較大的變形，此時(shí)曲率的選取就不能近似等于式(8)，而表達(dá)為

(9)

設(shè)由彎曲變形引起的橫截面的轉(zhuǎn)角為φ(x,t)，由材料力學(xué)的知識(shí)可知，其值是一個(gè)相對(duì)較小的量，故可近似認(rèn)為

(10)

(11)

由于φ比較小，將式(9)和式(11)代入式(5)，簡(jiǎn)化為

(12)

活塞理論是一種非定常氣動(dòng)理論，這種理論被廣泛應(yīng)用于超聲速和高超聲速翼面氣動(dòng)彈性分析中。根據(jù)一階活塞理論氣動(dòng)力模型[22]，得到二自由度系統(tǒng)受到的氣動(dòng)力La和氣動(dòng)力矩Ma為

(13)

(14)

系統(tǒng)虛功為

(15)

運(yùn)用Hamilton原理

(16)

將式(4)、式(12)、式(15)代入式(16)得動(dòng)力學(xué)方程

(17)

(18)

邊界條件

當(dāng)x=0時(shí)，

w(0,t)=θ(0,t)=0,w′(0,t)=0

(19)

當(dāng)x=L時(shí),

(20)

引入無量綱量

式(17)、式(18)變?yōu)?/p>

(21)

(22)

其中，

(23)

(24)

運(yùn)用Galerkin方法對(duì)無量綱方程做一階截?cái)啵x取一階振型離散

(25)

(26)

其中，

(27)

(28)

β=1.875

(29)

將式(25)～式(29)代入式(21)和式(22)中，然后進(jìn)行Galerkin正交化，整理得到兩自由度的非線性動(dòng)力學(xué)方程

(30)

系數(shù)Ai和Bj(i=1,2,,7；j=1,2,,6)為系統(tǒng)參數(shù)和轉(zhuǎn)速ω(無量綱化的角速度)的函數(shù)，其表達(dá)式詳見附錄。

2 數(shù)值響應(yīng)分析

2.1 臨界顫振速度

分析彎曲位移和扭轉(zhuǎn)位移在不同轉(zhuǎn)速下的時(shí)間響應(yīng)曲線，探討系統(tǒng)穩(wěn)定性。圖2為轉(zhuǎn)速ω=5時(shí)，彎曲位移和扭轉(zhuǎn)位移時(shí)間歷程圖。由圖可以看出葉片振動(dòng)收斂到零，是穩(wěn)定的。圖3和圖4分別為轉(zhuǎn)速ω=7時(shí)，彎曲位移和扭轉(zhuǎn)位移的時(shí)間歷程圖和相圖。由它們可以看出此時(shí)葉片收斂到穩(wěn)定極限環(huán)(Limit Cycle Oscillation, LCO)，發(fā)生極限環(huán)顫振。這說明在轉(zhuǎn)速5～7存在一個(gè)臨界轉(zhuǎn)速，系統(tǒng)發(fā)生Hopf分岔，由穩(wěn)定的零解分岔出穩(wěn)定的極限環(huán)。

圖2 ω=5時(shí)，無量綱彎曲位移和扭轉(zhuǎn)位移隨時(shí)間 變化的曲線Fig.2 The time history of displacement bending and torsion (dimensionless) at the speed ω=5

圖3 ω=7時(shí)，無量綱彎曲位移和扭轉(zhuǎn)位移隨時(shí)間 變化的曲線Fig.3 The time history of displacement bending and torsion (dimensionless) at the speed ω=7

圖4 ω=7時(shí)，無量綱彎曲和扭轉(zhuǎn)的相圖Fig.4 The phase diagram dimensionless bending and torsion at the speed ω=7

為了研究旋轉(zhuǎn)葉片的顫振臨界解，引入向量

(31)

則式(30)可表示為矩陣形式

(32)

式中：β為系數(shù)矩陣；C5和D5為系數(shù)，具體表達(dá)式見附錄。

這樣可以通過式(32)的特征值問題來求系統(tǒng)顫振臨界解。一般特征值都是復(fù)數(shù)，當(dāng)特征值的實(shí)部都是負(fù)數(shù)時(shí)，旋轉(zhuǎn)葉片的振動(dòng)是衰減振動(dòng)，即系統(tǒng)是穩(wěn)定的。當(dāng)特征值的實(shí)部有一個(gè)等于零時(shí)，葉片作簡(jiǎn)諧運(yùn)動(dòng)，處于臨界狀態(tài)。當(dāng)特征值的實(shí)部有一個(gè)是正數(shù)，則旋轉(zhuǎn)葉片振動(dòng)發(fā)散。通常情況下，旋轉(zhuǎn)葉片系統(tǒng)的所有特征值實(shí)部都是負(fù)數(shù)，隨著轉(zhuǎn)速的增加，系統(tǒng)特征值的實(shí)部有一個(gè)會(huì)提前變?yōu)榱悖@時(shí)，旋轉(zhuǎn)葉片失穩(wěn)。而使特征值最大實(shí)部達(dá)到零的轉(zhuǎn)速就是旋轉(zhuǎn)葉片發(fā)生顫振的臨界轉(zhuǎn)速。圖5表示特征值最大的實(shí)部隨轉(zhuǎn)速的關(guān)系圖，可以看出轉(zhuǎn)速ω=6.149時(shí)，特征值最大實(shí)部等于零，這時(shí)解有一對(duì)純虛根，說明達(dá)到臨界條件，發(fā)生顫振。

圖5 特征值最大實(shí)部值隨轉(zhuǎn)速ω的變化關(guān)系Fig.5 The relationship between the maximum real value of the eigenvalue and the speed ω

2.2 極限環(huán)響應(yīng)

非線性因素的存在，使得氣動(dòng)彈性系統(tǒng)在顫振點(diǎn)處發(fā)生了Hopf分岔。轉(zhuǎn)速大于顫振轉(zhuǎn)速時(shí)，非線性系統(tǒng)將會(huì)出現(xiàn)極限環(huán)，而不是線性系統(tǒng)的氣動(dòng)彈性失穩(wěn)。根據(jù)Hopf分岔定理，分析極限環(huán)響應(yīng)。圖6為彎曲運(yùn)動(dòng)和扭轉(zhuǎn)運(yùn)動(dòng)的LCO幅值隨轉(zhuǎn)速的變化圖，系統(tǒng)在K點(diǎn)(ω=6.149)發(fā)生Hopf分岔后，出現(xiàn)了極限環(huán)，并且極限環(huán)的幅值隨轉(zhuǎn)速ω的增大而增大。

圖6 彎曲運(yùn)動(dòng)和扭轉(zhuǎn)運(yùn)動(dòng)LCO幅值(Aq和 Ap)隨轉(zhuǎn)速ω的變化關(guān)系Fig.6 The relationship between bending and torsional amplitude of LCO (Aq and Ap)and the rotational speed ω

諧波平衡法是求解非線性振動(dòng)問題常用的一種近似解析法，它的求解過程歸結(jié)為代數(shù)方程組的求解。本文運(yùn)用諧波平衡法研究葉片的極限環(huán)響應(yīng)?？紤]方程是自由振動(dòng)方程，無外激勵(lì)，所以可以設(shè)解的模式為

q=a1cos?t+b1sin?tp=a2cos?t

(33)

將式(33)代入式(30)得

(34)

(-B1a1?2-B2a1?2+B4a2+B5b1?+B6a1)cos?t+ (-B2b1?2-B3a2?-B5a1?+B6b1)sin?t=0

(35)

略去3次諧波，由cos?t及sin?t項(xiàng)的系數(shù)等于零可得

(36)

(37)

-B1a2?2-B2a1?2+B4a2+B5b1?+B6a1=0

(38)

-B2b1?2-B3a2?-B5a1?+B6b1=0

(39)

由式(36)～式(39)得到a1,b1,a2,?的近似解。

對(duì)解進(jìn)行分析，與數(shù)值解進(jìn)行對(duì)比。如圖7所示，諧波平衡法的結(jié)果與數(shù)值結(jié)果吻合較好。兩種方法都說明在轉(zhuǎn)速ω=6.149時(shí)葉片發(fā)生顫振，出現(xiàn)極限環(huán)。而且隨著轉(zhuǎn)速的增加，極限環(huán)幅值A(chǔ)q的增長(zhǎng)趨勢(shì)也相似。

圖7 解析解與數(shù)值解的結(jié)果對(duì)比Fig.7 Comparison of analytical solution and numerical solution

3 結(jié)構(gòu)參數(shù)對(duì)臨界顫振轉(zhuǎn)速及極限環(huán)響應(yīng)的影響

3.1 對(duì)臨界顫振轉(zhuǎn)速的影響

取l=1.5,n=1時(shí)，分析質(zhì)量塊質(zhì)量對(duì)顫振轉(zhuǎn)速的影響。從圖8可知，由于質(zhì)量塊的存在，改變了原有系統(tǒng)的穩(wěn)定性。0.1≤S≤0.2時(shí)，顫振轉(zhuǎn)速ωf急劇下降，說明當(dāng)質(zhì)量塊質(zhì)量占梁質(zhì)量達(dá)到某一值時(shí)，開始影響系統(tǒng)的穩(wěn)定性。隨著S的增大導(dǎo)致顫振轉(zhuǎn)速ωf減小，系統(tǒng)容易發(fā)生顫振。

圖9為顫振轉(zhuǎn)速隨質(zhì)量塊位置的變化曲線。取l=1.5,S=0.16分析。當(dāng)0.5≤n≤0.6時(shí)，發(fā)生顫振的轉(zhuǎn)速ωf有一個(gè)大的變化，質(zhì)量塊開始影響系統(tǒng)的穩(wěn)定性； 當(dāng)n≥0.7時(shí)，臨界顫振轉(zhuǎn)速隨n的增大而減小，近似呈現(xiàn)線性關(guān)系。不難看出隨著質(zhì)量塊越靠近梁自由端，越容易發(fā)生顫振。

圖8 顫振轉(zhuǎn)速ωf隨質(zhì)量塊質(zhì)量S(無量綱)的變化關(guān)系Fig.8 The relationship between the flutter speed ωf and the mass S(dimensionless) of the mass

圖9 顫振轉(zhuǎn)速ωf隨質(zhì)量塊位置n(無量綱)的變化關(guān)系Fig.9 The relationship between the flutter speed ωf and the mass position n

圖10表示無質(zhì)量塊時(shí)，改變梁的長(zhǎng)度對(duì)顫振速度的影響。一開始長(zhǎng)度較短時(shí)，不會(huì)發(fā)生顫振。達(dá)到一臨界值后，開始發(fā)生顫振。之后隨著梁長(zhǎng)度的增大，發(fā)生顫振的幾率就隨之增大。Resor等[23]也指出隨著葉片長(zhǎng)度增長(zhǎng)，顫振速度會(huì)逐漸減小。分別取S=0.16和S=0.40(相對(duì)于梁L=0.15 m時(shí))，n=1，分析有質(zhì)量塊時(shí)，梁長(zhǎng)度對(duì)顫振速度的影響，如圖11所示。由于質(zhì)量塊的存在，在葉片較短時(shí)會(huì)發(fā)生顫振，并且隨著梁長(zhǎng)度的增加，顫振速度也隨之增加。但是，當(dāng)葉片達(dá)到一定長(zhǎng)度，質(zhì)量塊的影響就基本不存在了。這時(shí)和無質(zhì)量塊時(shí)的情況一樣，會(huì)隨著葉片長(zhǎng)度的增大，而導(dǎo)致容易發(fā)生顫振。說明由于質(zhì)量塊存在，短梁會(huì)受到一定的影響。

圖10 無質(zhì)量塊時(shí)，顫振轉(zhuǎn)速ωf隨梁長(zhǎng)度l的變化關(guān)系Fig.10 No mass, the relationship between the flutter speed ωf and the beam length l

圖11 有質(zhì)量塊時(shí)，顫振轉(zhuǎn)速ωf隨梁長(zhǎng)度l的變化關(guān)系Fig.11 There is mass, the relationship between the flutter speed ωf and the beam length l

圖12 顫振轉(zhuǎn)速ωf隨和的變化關(guān)系Fig.12 The relationship between the flutter speed ωf

3.2 對(duì)極限環(huán)響應(yīng)的影響

圖13和圖14是分析質(zhì)量塊對(duì)極限環(huán)響應(yīng)的幅值影響。固定l和n，分別取S為0.16, 0.25, 0.50, 0.75，得到圖13所示的四條曲線，可以看出質(zhì)量塊重量越大，越容易發(fā)生Hopf分岔，極限環(huán)幅值越大。同理，固定l和S，分別取n為0.85, 0.90, 0.95, 1.00，圖14說明質(zhì)量塊位置越靠近末端，越容易發(fā)生Hopf分岔，同時(shí)極限環(huán)幅值越大。

圖13 質(zhì)量塊不同質(zhì)量下的LCO幅值A(chǔ)q 隨轉(zhuǎn)速ω的變化關(guān)系Fig.13 The relationship between the LCO amplitude Aq and the speed ω at different quality of the mass

圖14 質(zhì)量塊不同位置下的LCO幅值隨轉(zhuǎn)速ω的變化關(guān)系Fig.14 The relationship between the LCO amplitude Aq and the speed ω at different positions of the mass

4 結(jié) 論

本文建立了帶冠葉片的兩自由度耦合非線性動(dòng)力學(xué)方程。由特征值理論得到了顫振臨界轉(zhuǎn)速，用諧波平衡法和數(shù)值方法對(duì)比分析了極限環(huán)響應(yīng)，結(jié)果發(fā)現(xiàn)兩者吻合的很好。研究了旋轉(zhuǎn)狀態(tài)下，質(zhì)量塊的質(zhì)量和位置、梁的長(zhǎng)度和寬度等參數(shù)對(duì)振動(dòng)響應(yīng)的影響規(guī)律。結(jié)果表明，質(zhì)量塊對(duì)顫振速度和極限環(huán)幅值有明顯的影響，質(zhì)量塊的質(zhì)量越大，位置越靠近末端，越容易發(fā)生顫振，極限環(huán)幅值越大。同時(shí)質(zhì)量塊的存在會(huì)對(duì)短梁有一定影響。本文主要針對(duì)不考慮葉片之間耦合的單個(gè)帶冠葉片的超音速、高超音速非線性氣動(dòng)彈性的穩(wěn)定性分析，對(duì)多個(gè)葉片之間的耦合問題將在后續(xù)研究中給與考慮。