大展弦比機翼氣動顫振的有限元分析

楊瓊梁，史曉鳴，2，唐國安

(1.復旦大學 力學與工程科學系，上海 200433;2.上海機電工程研究所)

大展弦比機翼氣動顫振的有限元分析

楊瓊梁1，史曉鳴1，2，唐國安1

(1.復旦大學 力學與工程科學系，上海 200433;2.上海機電工程研究所)

在Theodorsen二元氣動力的基礎上，建立非定常氣動力時域內積分形式的表達式或者等價的頻域表達式，利用粘彈性結構振動分析中對積分方程的等價變換將其寫成與結構動力學方程一致的二階常微分方程，將氣動力的影響作為對結構有限元模型質量陣、剛度陣和阻尼陣的補充，保留了結構原有的所有動力學特性，并且能夠直接用計算結構動力學的通用有限元軟件進行空氣－結構耦合的整體動力學分析，適合應用于具有復雜結構的氣彈問題。氣動力模型的建立可以利用各種試驗及數值方法得到的氣動力數據，適用性強。算例給出了大展弦比機翼的顫振邊界計算結果。

氣動彈性;非定常氣動力;顫振;有限元

氣動彈性是研究彈性體與氣動力之間相互作用的一門學科，是現代飛行器設計中的一個重要環節。其中，顫振是最引人關注的氣動彈性動力學問題，顫振通常能導致飛行器結構的毀滅性破壞［1］。由于氣彈分析中空氣與彈性體相互耦合，而目前還沒有一種分析工具能很好的對流體和彈性體兩者同時進行分析計算，這也使得氣彈分析變得較為復雜。

目前常見的氣彈分析方法有:CFD/CSD耦合計算方法［2］和對氣動力進行簡化后的快速工程分析方法，例如基于當地流活塞理論的分析法［3－5］。其中 CFD/CSD耦合法計算精度較高，但是采用CFD計算非定常氣動力時計算量大，且要在每個時間步內交換CFD和CSD兩者的計算數據，造成整體的計算效率低;為了克服用CFD計算所帶來的計算困難，實際應用中通常采用一些快速工程算法對非定常氣動力計算做簡化，如利用當地流活塞理論計算氣動力，但活塞理論只在來流馬赫數M∞?1的情況下有較為滿意的結果，計算范圍有所限制。

在上世紀90年代，Silva提出將Volterra級數理論用于非線性非定常氣動力的建模，建立起具有積分形式的氣動力降階模型，氣動力計算精度高，并在一定程度上實現了氣動力與結構的解耦［6－7］，但是積分形式的氣動力不屬于傳統的結構激勵載荷，從而造成計算復雜。本文的目的是在以積分形式表示的非定常氣動力的基礎上，在時域內將非定常氣動力表示為與結構動力學方程一致的二階常微分方程，從而將非定常氣動力對結構的影響表示為對結構的附加質量、剛度和阻尼，并以此為依據在物理坐標下直接對結構進行修改，實現氣動力與結構的解耦，在保留結構完整的動力學特性的條件下，應用傳統的通用有限元軟件進行空氣－結構耦合的整體動力學分析。

1 非定常氣動力建模及氣動力與結構耦合

有限元離散模型表示的彈性機翼在非定常氣動力作用下的振動方程為:

其中u為結構振動的節點位移向量;Ms，Cs和Ks為結構的質量、阻尼和剛度矩陣;－fα為作用在結構上的非定常氣動力，通常可以表示為結構位移、速度和加速度的函數。

Fung［8］給出了機翼在不可壓縮流中積分形式的升力系數和力矩系數，將氣動力表示為:

其中Ma，Ca，Ka，Vp和Ap為相關的氣動力參數。利用降階法［9］可以通過各種氣動力計算方法或者試驗得到頻域內的氣動力數據點，并擬合成有理多項式，經過Laplace變化就能得到如(2)所示的非定常氣動力表達式。

將非定常氣動力式(2)代入式(1)，得到結構與氣動力的耦合方程:

高淑華［10］在分析粘彈性結構振動時，證明了方程(3)表示的微分－積分型動力學方程可以等價地變換成純微分型方程:

其中U1i、∑i等是矩陣Vp(p=1，…，k)的奇異值分解Vp=U1p∑pU2p。

周文博［11］利用這種等價變換的方法，分析了帶有控制面的三元機翼顫振問題。證實了通過等價變換、可以利用通用有限元結構分析程序實現氣動彈性耦合的顫振計算。但文獻［11］僅僅考慮了簡單的機翼模型，機翼是用彈簧和扭簧連接的二自由度剛體，控制面為由扭簧連接到機翼的單自由度剛體。這種模型簡化較大，與實際彈性機翼存在很大差異。

2 大展弦比三維彈性機翼與氣動力耦合模型

用 彈性平板有限元模型表示大展弦比三維彈性機翼，如圖1所示。 氣 彈計算將采用MSC.Nastran有限元程序，平板表示為80個CQUAD單元、節點數為99。

圖1 大展弦比彈性機翼模型Fig.1 Model of the high-aspect ratio elastic wing

考慮到大展弦比機翼，機翼片條ab的弦向彈性彎曲可以忽略不計、弦向變形視為剛性，即直線ab變形后仍近似保持直線，其變形可以由某個節點的z向位移w(c)和繞x軸轉角θ(c)z來表示。作用在整個機翼上的氣動力可以近似地認為是由一系列作用在機翼片條上的氣動力的累加，而片條上連續分布的氣動力又可以簡化為作用在節點c上的一對集中力(升力)L(c)和集中力矩M(c)。

其中ΔL(c)為片條長度，系數見附錄。

3 大展弦比均勻直機翼顫振實例計算

根據片條理論將直機翼沿展向分為10塊，利用上述方法通過改變飛行速度V，得到每一塊的氣動力附加矩陣。由于=0，即氣動力和力矩的作用點取在翼面中點，利用上述得到的氣動力附加矩陣修改該點的質量、剛度、阻尼，最后用有限元分析程序(MSC.Nastran)的復特征值分析功能，計算空氣－結構耦合情況下的復特征值，以特征值實部穿越零點作為依據來計算得到機翼的顫振邊界。文獻［12］中，將機翼近似為梁。計算結果表明第三階振型發生顫振，相應的復特征值頻率的實部變化如圖2所示。本文計算所得到的顫振結果基本與文獻結果如表2所示，誤差小于6%，由于對氣動力擬合時會引起一定誤差，增加擬合階數將有助于提高計算精度。

表1 大展弦比機翼參數Tab.1 Parameter of high －aspect ratio wing

表2 臨界顫振速度與頻率Tab.2 Critical flutter speed and frequency

圖2 第三階頻率實部Fig.2 Real part of the third mode

4 結論分析

利用粘彈性結構振動分析中對積分方程的等價變換方法，推導了積分形式的非定常氣動力與結構動力學方程耦合的氣動彈性分析方法，并以大展弦比直機翼有限元模行為算例，計算得到了其顫振頻率和顫振速度，驗證了方法的可行性。

將氣動力作為對有限元中結構質量陣、剛度陣和阻尼陣的補充，完整保留了結構原有的動力學特性，避免了模態降階帶來的誤差。對復雜結構具有較高的適應性，不僅可以用于計算平板單元的顫振，對其它各類有限元單元，如桿系、殼體等單元同樣適用。直接利用有限元軟件進行顫振分析，計算效率高。

文中的氣動力模型來自于Theodorsen理論，但方法同樣適用于以一階Volterra級數表示的氣動力降階模型。只要通過CFD全流場計算或者試驗流場分析等手段得到氣動力的卷積表達式，就能對復雜結構進行氣動－彈性耦合計算，具有很好的工程適用性。

［1］Dowell E H，Clark R，Cox D，et al.A Mordern Course in Aeroelasticity［M］.New York:Kluwer Academic Publishers，2004:1－5.

［2］ Chen X Y，Zhan G C，Hu Z J.Numerical Simulation of Flow Induced Vibration Based on Fully Couple Fluid-Structual Interactions［R］.AIAA-2004 －2240.

［3］楊炳淵，宋偉力.用當地流活塞理論計算大功角翼面超音速顫振［J］.振動與沖擊，1995，14(2):60－63.

［4］楊炳淵，宋偉力.前緣激波脫體的大迎角翼面顫振工程計算方法［J］.振動與沖擊，2002，21(4):100－103.

［5］Zhang W W，Ye Z Y，Zhang C A.Supersonic Flutter Analysis Based on a Local Piston Theory［J］.AIAA，2009，47(10)，2321－2328.

［6］張偉偉，葉正寅.基于非定常氣動力辨識技術的氣動彈性數值模擬［J］.航空學報，2006，27(4):579－584.

［7］張偉偉，葉正寅.基于CFD的氣動力建模及其在氣動彈性中的應用［J］.力學進展，2008，38(1):77－86.

［8］Fung Y C.An Introduction to the Theory of Aeroelasticity［M］.New York，John Wiley＆ Sons，Inc，1955.

［9］Piergiovanni M，Liviu L，Walter A.S.Volterra Series Approach for Nonlinear Aeroelastic Response of 2-D Lifting Surfaces［R］.AIAA －2001－1459.

［10］高淑華，趙 陽，等.粘彈性結構動力學分析的等效粘性阻尼算法［J］.振動與沖擊，2005，24(1):18－21.

［11］周文博，陳力奮，楊瓊梁，等.基于結構動力學方法的氣動彈性分析［J］.振動與沖擊，2011，30(4):135－138.

［12］趙永輝.氣動彈性力學與控制［M］.北京:科學出版社，2007，134－135.

附錄

式(6)中的系數

Flutter analysis for a high-aspect ratio wing with finite element method

YANG Qiong-liang1，SHI Xiao-ming1，2，TANG Guo-an1

(1.Department of Mechanics and Engineering Science，Fudan University，Shanghai 200433，China;2.Shanghai Electro-Mechanical Engineering Institute)

Expressions for unsteady aerodynamics in time domain or frequency domain were established based on Theodorsen theory.The expressions with integration were converted into second order ordinary differential equations.The effects of unsteady aerodynamic force on a structure were regarded as comlements of mass matrix，stiffness one and damping one of finite element model of a structure.The whole dynamic characteristics were retained.General finite element software could be used to calculate the flutter boundary of an air-structure coupled model.The method was suitable for solving aeroelastic problems of complex structures.The unsteady aerodynamic data obtained with any numerical method or test could be used to establish an unsteady aerodynamic model.The proposed method was applicable to different cases.The flutter boundary of a high-aspect ratio wing was calculated.

aeroelasticity;unsteady aerodynamics;flutter;finite element method

V211.47

A

國家自然科學基金資助項目(90716001)

2010－05－21 修改稿收到日期:2010－08－20

楊瓊梁 男，博士生，1984年7月生

唐國安 男，教授，1962年生