周健斌，章俊杰，孟 光

（1.中國商用飛機有限責(zé)任公司 上海飛機設(shè)計研究院，上海 201210；2.上海交通大學(xué) 機械系統(tǒng)與振動國家重點實驗室，上海 200240）

0 引 言

具有陀螺效應(yīng)的彈性體在航空、航天、機械等領(lǐng)域有廣泛的應(yīng)用。D′Eleuterio［2－6］和 Hughes［7］針對自旋航天器研究了一類具有分布陀螺力矩的彈性系統(tǒng)的動力學(xué)特性并提出了陀螺彈性系統(tǒng)的概念。Peck［8－9］等在此基礎(chǔ)上提出了運用嵌入在結(jié)構(gòu)中的角動量來實現(xiàn)結(jié)構(gòu)自適應(yīng)調(diào)諧，如頻移、模態(tài)耦合和解耦、相位調(diào)整等。針對一類帶有旋轉(zhuǎn)部件的機械臂，Yamanaka［10－11］分析了端部帶有轉(zhuǎn)子的均質(zhì)梁的動力學(xué)和穩(wěn)定性問題，其轉(zhuǎn)子轉(zhuǎn)軸沿梁的軸向。Li［12］等研究了端部帶有任意方向旋轉(zhuǎn)軸的轉(zhuǎn)子的柔性連接的穩(wěn)定性問題。

國內(nèi)外關(guān)于陀螺效應(yīng)對機翼或飛機顫振特性的影響機理研究文獻(xiàn)較少。在實踐中常根據(jù)經(jīng)驗將陀螺效應(yīng)等效成一定的系統(tǒng)阻尼，在ω法或p－k法計算時采用阻尼矩陣計入陀螺效應(yīng)進(jìn)行顫振計算［1］，取得較好的效果。盡管陀螺力矩在形式上與系統(tǒng)速度有關(guān)，但其功率恒等于零，即其本身并不是阻尼力［13］，因此將陀螺效應(yīng)等效成系統(tǒng)阻尼的方法缺少一定的理論支撐。王彬文［14］等用有限元方法研究了轉(zhuǎn)子陀螺效應(yīng)對帶有翼吊發(fā)動機系統(tǒng)振動特性的影響。

本文從機翼結(jié)構(gòu)動力學(xué)和轉(zhuǎn)子動力學(xué)基本理論出發(fā)，建立了計及陀螺效應(yīng)的機翼彎扭顫振方程，忽略發(fā)動機吊掛柔性等影響，分析了陀螺效應(yīng)對系統(tǒng)顫振特性的影響規(guī)律，為進(jìn)一步分析機翼－發(fā)動機系統(tǒng)顫振特性提供一定的理論依據(jù)。

1 顫振方程的建立

如圖1（a）所示，將機翼－發(fā)動機系統(tǒng)簡化為一端固支的彎扭耦合Bernoulli－Euler懸臂梁，忽略發(fā)動機與機翼間的連接剛度。機翼翼展為l，半弦長為b，右手坐標(biāo)系o－xyz原點o位于機翼根部，ox軸與機翼剛軸重合，oy軸沿飛機航向。圖1（b）所示為機翼變形后x＝x0處弦向截面示意圖。o′為機翼剛軸與截面交點，截面位移為v、w，截面繞o′轉(zhuǎn)角為θ。Cw為截面重心，沿弦向與剛軸距離為e。發(fā)動機重心位置Cs，坐標(biāo)為（x0，y0，z0），發(fā)動機重量為 ms，轉(zhuǎn)子轉(zhuǎn)速為Ω，旋轉(zhuǎn)矢量方向與剛軸垂直，來流速度為U∞。

圖1 機翼－發(fā)動機系統(tǒng)示意圖Fig.1 Illustration of the Wing－Engine system

由Hamilton變分原理

可推導(dǎo)得系統(tǒng)動力學(xué)方程，其中δU、δTw、δTs和δW分別為機翼勢能、機翼動能、發(fā)動機動能和外力做功的變分。由彎扭耦合懸臂梁動力學(xué)分析，機翼勢能的變分為：

其中，EI2為梁的面外彎曲剛度，EI3為梁的面內(nèi)彎曲剛度，GJ為梁的扭轉(zhuǎn)剛度。

本文通過建立動車組單輛車整體稱重模型，提出了3個調(diào)平指導(dǎo)參數(shù)以及高度調(diào)整閥調(diào)平條件。3個調(diào)平指導(dǎo)參數(shù)分別反映了影響動車組單輛車整體稱重調(diào)平的3個獨立因素：前轉(zhuǎn)向架不平、后轉(zhuǎn)向架不平以及車輛重心偏心。結(jié)合現(xiàn)場數(shù)據(jù)驗證了高度調(diào)整閥調(diào)平條件的正確性。在此基礎(chǔ)上可以利用調(diào)平指導(dǎo)參數(shù)來進(jìn)行判斷，并明確指出應(yīng)在轉(zhuǎn)向架加設(shè)墊片以滿足調(diào)平條件。應(yīng)用此稱重調(diào)平規(guī)律和3個調(diào)平指導(dǎo)參數(shù)以及高度調(diào)整法調(diào)平條件來指導(dǎo)車輛現(xiàn)場稱重，可以有效避免工人的盲目操作，顯著提高車輛現(xiàn)場稱重效率，還為計算機編寫計算動車組單輛車整體稱重調(diào)平程序提供了必要的理論支持，同時也為車輛稱重調(diào)平的數(shù)字化、智能化提供指導(dǎo)方向。

機翼動能的變分可表示為

其中，mw為機翼單位長度質(zhì)量，e為機翼截面重心與剛軸的距離（重心在剛軸前為正），σ為機翼單位長度截面繞剛軸的慣性半徑。

對發(fā)動機動能分析可知，發(fā)動機的動能的變分由兩部分組成，δTs＝δTm＋δTr，其中δTm為發(fā)動機集總質(zhì)量的動能的變分，δTr為發(fā)動機轉(zhuǎn)子的轉(zhuǎn)動動能的變分。

δD＝Dirac（x0－x）。當(dāng)機翼變形時，向量Ω的偏轉(zhuǎn)角在oxy和oyz平面的投影分別為v′和θ，由轉(zhuǎn)子動力學(xué)理論可知，發(fā)動機轉(zhuǎn)子的轉(zhuǎn)動的動能可表示為：

其中，Jp為轉(zhuǎn)子繞其轉(zhuǎn)軸的極轉(zhuǎn)動慣量，Jd為轉(zhuǎn)子的赤道轉(zhuǎn)動慣量。對Tr求變分可得：

非定常氣動力做虛功為：

其中，La和Ma分別為氣動升力和力矩，根據(jù)Theodorsen非定常氣動力理論：

Lw、Lθ、Mw和Mθ為氣動系數(shù)［15］。

將式（2）～（4）和（6）、（7）代入式（1）可得如下計及發(fā)動機轉(zhuǎn)子效應(yīng)的動力學(xué)方程：

本文采用extended Galerkin［16］對上述系統(tǒng)動力學(xué)方程進(jìn)行離散化并求解顫振特性，該方法僅需基函數(shù)滿足幾何邊界條件的要求，而將自然邊界條件包含在控制方程中。因此，令

選取如下滿足幾何邊界條件的函數(shù)作為基函數(shù)：

代入方程（10）～（12），積分可得離散后的系統(tǒng)顫振方程：

其中q＝｛q11，…，q1N1，q21，…，q2N2，q31，…，q3N3｝T，質(zhì)量矩陣M為對稱正定矩陣，剛度矩陣K、氣動力矩陣Qhh為對稱矩陣，由于系統(tǒng)陀螺力項與和相關(guān)，陀螺矩陣G為非對稱矩陣。對系統(tǒng)動力學(xué)方程（19）運用p－k法可求得系統(tǒng)的顫振特性。

2 算例分析

采用Runyan機翼模型［17］結(jié)構(gòu)參數(shù)，對系統(tǒng)顫振特性進(jìn)行分析，各主要參數(shù)如下：l＝1.2192m，b＝0.1016m，EI2＝403.76N·m2，EI3＝403.76N·m2，GJ＝198.76N·m2，mw＝1.2942kg·m－1，ms＝1.578kg，ρ＝1.224kg·m－3，σ＝0.053m，Jp＝0.01kg·m2，Jd＝0.0185kg·m2。定義無量綱慣性矩用以表征轉(zhuǎn)子慣性矩與梁結(jié)構(gòu)剛度、質(zhì)量的相對大小，式中圖2所示為發(fā)動機轉(zhuǎn)速Ω＝0rpm時，發(fā)動機位于不同機翼展向位置x0／l時的系統(tǒng)彎扭顫振速度和顫振頻率，圖中同時給出了文獻(xiàn)［17］的試驗結(jié)果，從圖可以看出，采用本文簡化模型計算得到的結(jié)果與文獻(xiàn)中結(jié)果吻合得較好。

圖2 Ω＝0rpm時發(fā)動機不同展向位置的系統(tǒng)顫振速度和顫振頻率Fig.2 Flutter speeds and frequencies of the system with different spanwise locations andΩ＝0rpm

如圖3所示為發(fā)動機轉(zhuǎn)子轉(zhuǎn)速發(fā)動機位于機翼展向位置x0／l＝22.9%且Ω＝0rpm（即沒有陀螺效應(yīng)作用）時系統(tǒng)顫振V－g、V－ω圖。當(dāng)Ω＝0rpm，系統(tǒng)彎扭顫振速度為91.18m／s，顫振頻率為18.17Hz。當(dāng)計及發(fā)動機陀螺效應(yīng)時，在陀螺力矩的作用下系統(tǒng)固有頻率和固有振型發(fā)生變化，如圖4所示為系統(tǒng)固有頻率隨發(fā)動機轉(zhuǎn)速（無量綱慣性矩）變化的Campbell圖。根據(jù)陀螺力矩理論，結(jié)合系統(tǒng)結(jié)構(gòu)特點可知，轉(zhuǎn)子轉(zhuǎn)動方向隨著機翼結(jié)構(gòu)振動而變化，陀螺力矩作用在機翼面內(nèi)彎曲運動和扭轉(zhuǎn)運動方向，從而使得涉及上述兩運動的系統(tǒng)模態(tài)相互耦合。

圖3 x0／l＝22.9%、Ω＝0rpm時的V－g圖和V－ω 圖Fig.3 V－gand V－ω graph of the system with x0／l＝22.9%andΩ＝0rpm

圖4 系統(tǒng)固有頻率Campbell圖Fig.4 Campbell diagram of the system natrual frequencies

陀螺效應(yīng)對系統(tǒng)固有特性的耦合效應(yīng)將進(jìn)一步影響其顫振特性。如圖5所示為系統(tǒng)顫振特性隨發(fā)動機轉(zhuǎn)速變化曲線，由圖可知，隨著發(fā)動機轉(zhuǎn)速的提高，系統(tǒng)顫振速度和顫振頻率都逐漸增大。

圖5 發(fā)動機轉(zhuǎn)速對系統(tǒng)顫振速度和顫振頻率的影響曲線Fig.5 Influence of gyro effects on the flutter velocity and frequency

進(jìn)一步考察Ω＝5730rpm（ˉH＝0.257）時系統(tǒng)的顫振特性，其顫振V－g、V－ω圖如圖6所示，此時顫振速度為92.96m／s，顫振頻率為18.52Hz。對比圖3和圖6可知，當(dāng)不考慮系統(tǒng)吊掛剛度時，在陀螺效應(yīng)作用下，機翼－發(fā)動機系統(tǒng)的V－g、V－ω圖基本相似，即系統(tǒng)顫振發(fā)生機理未發(fā)生根本性變化，因此，發(fā)動機陀螺效應(yīng)可等效成一定的系統(tǒng)結(jié)構(gòu)阻尼geq。等效阻尼可由如下方法計算：設(shè)在一定的發(fā)動機轉(zhuǎn)速下系統(tǒng)顫振速度為Vf，在不計陀螺效應(yīng)時系統(tǒng)顫振V－g曲線的相應(yīng)顫振模態(tài)分支上空速為Vf時對應(yīng)的阻尼比即為該轉(zhuǎn)速所對應(yīng)的等效阻尼。

根據(jù)上述等效阻尼的定義及計算方法可得Ω＝5730rpm（ˉH＝0.257）時，轉(zhuǎn)子陀螺效應(yīng)等效阻尼為geq＝0.022，如圖3中水平虛線所示。圖7給出了系統(tǒng)等效阻尼隨著發(fā)動機轉(zhuǎn)速提高的變化曲線。由圖可知，等效阻尼亦隨著發(fā)動機轉(zhuǎn)速的增加呈拋物線型增大。

圖6 x0／l＝22.9%、＝5730rpm時的V－g圖和V－ω 圖Fig.6 V－gand V－ω graph of the system with x0／l＝22.9%andΩ＝5730rpm

圖7 轉(zhuǎn)子效應(yīng)等效阻尼曲Fig.7 Curve of the gyroscopic equivalent damping

3 結(jié) 論

本文通過建立計及發(fā)動機轉(zhuǎn)子陀螺效應(yīng)的機翼彎扭顫振分析方程，在不考慮發(fā)動機吊掛剛度的條件下，通過分析得到如下結(jié)論：

（1）在發(fā)動機轉(zhuǎn)子陀螺效應(yīng)的作用下，機翼彎扭型顫振速度和顫振頻率隨轉(zhuǎn)子轉(zhuǎn)速的增加而增大，且發(fā)動機轉(zhuǎn)子陀螺效應(yīng)不影響機翼顫振形式；

（2）盡管陀螺力矩本身并非阻尼力，但由于其在機翼顫振特性中具有上述特點，仍可以將其對機翼顫振特性的影響等效成一定的系統(tǒng)阻尼，且該等效阻尼隨著發(fā)動機轉(zhuǎn)速的增加而增大。

［1］管德.飛機氣動彈性力學(xué)手冊［M］.航空工業(yè)出版社，1994.

［2］D′ELEUTERIO G M T.On the theory of gyroelasticity［J］.Journal of Applied Mechanics，Transactions ASME，1988，55（2）：488－489.

［3］D′ELEUTERIO G M T，HUGHES P C.Dynamics of gyroelastic spacecraft［J］.Journal of Guidance，Control，and Dynamics，1987，10（4）：401－405.

［4］D′ELEUTERIO G M T，HUGHES P C.General motion of gyroelastic vehicles in terms of constrained modes［C］.AIAA／ASME／ASCE／AHS 26th Structures，Structural Dynamics and Materials Conference.1985，384－390.

［5］D′ELEUTERIO G M T，HUGHES P C.Dynamics of gyroelastic continua［J］.Journal of Applied Mechanics，Transactions ASME，1984，51（2）：415－422.

［6］D′ELEUTERIO G M T.Dynamics of gyroelastic vehicles［R］.UTIAS Report（University of Toronto，Institute for Aerospace Studies），1986.

［7］HUGHES P C，D′ELEUTERIO G M T.Modal parameter analysis of gyroelastic continua［J］.Journal of Applied Mechanics，1986，53（4）：918－924.

［8］PECK M A，CAVENDER A R.Structural tuning through embedded angular momentum［A］.44th AIAA／ASME／ASCE／AHS／ASC Structures，Structural Dynamics，and Materials Conference［C］，2003，2：1462－1469.

［9］PECK M A，CAVENDER A R.Practicable gyroelastic technology［A］.27th Annual AAS Rocky Mountain Guidance and Control Conference［C］，2004，118：239－253.

［10］YAMANAKA K，HEPPLER G R，HUSEYIN K.On the dynamics and stability of a beam with a tip rotor［A］.The 35th AIAA／ASME／ASCE／AHS／ASC Structures，Structural Dynamics，and Materials Conference［C］，1994，2 1031－1038.

［11］YAMANAKA K，HEPPLER G R，HUSEYIN K.Stability of gyroelastic beams［J］.AIAA Journal，1996，34（6）：1270－1278.

［12］LI L，HEPPLER G R，HUSEYIN K.Stability of a flexible link with an arbitrarily oriented tip rotor and a conservative tip load［C］.IEEE Interational Conference on Robotics and Automation.2000，2：1472－1477.

［13］陳予恕，曹登慶，黃文虎.近代機械非線性動力學(xué)與優(yōu)化設(shè)計技術(shù)的若干問題［J］.機械工程學(xué)報，2007，43（011）：17－26.

［14］王彬文，孫俠生，齊丕騫.轉(zhuǎn)子陀螺效應(yīng)對系統(tǒng)振動特性的影響［J］.振動工程學(xué)報，2010，23（增刊）：63－67.

［15］HODGES D H，PIERCE G A.Introduction to structural dynamics and aeroelasticity［M］.Cambridge：Cambridge University Press，2002.

［16］FLETCHER C.Computational Galerkin methods［M］.New York：Springer－Verlag，1984.

［17］RUNYAN H，WATKINS C.Flutter of a uniform wing with an arbitrarily placed mass according to a differential－equation analysis and a comparison with experiment［R］.NACA Rep.1848，1950.