某全機跨聲速顫振模型顫振特性仿真與試驗驗證

錢 衛(wèi)，楊國偉，張桂江，鄭冠男

（1.沈陽飛機設(shè)計研究所，沈陽 110035；2.中國科學院力學研究所，北京 100190）

0 引 言

對于高機動飛機，經(jīng)典顫振邊界在跨聲速區(qū)往往存在一個凹坑，因此跨聲速顫振特性是飛機設(shè)計限制的關(guān)鍵設(shè)計參數(shù)。根據(jù)規(guī)范［1］要求，在設(shè)計初期就應進行顫振模型風洞試驗，以驗證規(guī)定的速度安全余量，并驗證變參數(shù)顫振的分析結(jié)果。

美國NASA蘭利研究中心TDT風洞從1960年至今，對各種飛機型號進行了大量的跨聲速顫振風洞試驗［2］，主要解決飛機型號的跨聲速顫振設(shè)計問題和跨聲速顫振仿真程序的評估與驗證。尤其是在F－16，B－777，F(xiàn)18E／F等飛機的顫振風洞試驗中取得的數(shù)據(jù)，對發(fā)展顫振仿真程序作出很大貢獻?？屏_拉多大學對F－16的飛行狀態(tài)進行了完整的CFD／CSM耦合的仿真，結(jié)構(gòu)FEM模型有168799自由度，CFD網(wǎng)格有403919空間點，63044表面點（小攻角歐拉方程仿真）［3－6］。對于CFD跨聲速顫振仿真，波音公司進一步發(fā)展了NASA基于歐拉／N－S方程的CFD程序CFL3D，進行跨聲速顫振仿真，并使用HSR項目的TDT風洞模型試驗結(jié)果對程序進行了修正，取得比較大的進展［7］，波音公司還使用CFL3D對B－1飛機的極限環(huán)振動（LCO）進行的全面的仿真研究［8－10］。美國在戰(zhàn)斗機的工程顫振設(shè)計中從F－18到F－22仍然采用線化理論偶極子網(wǎng)格法（DLM）進行顫振計算，同時通過跨聲速顫振風洞試驗進行驗證［11－12］。近年，ZONA公司發(fā)展了簡化的歐拉方程顫振仿真軟件ZEUS Code，可以進行飛機的跨聲速顫振仿真［13］。

國內(nèi)基于偶極子網(wǎng)格法（DLM）的顫振設(shè)計工作已形成比較完整的體系［14］，受跨聲速風洞試驗條件的限制，跨聲速顫振設(shè)計的準確性受到一定的限制。近年國內(nèi)在跨聲速CFD顫振仿真上取得很大的進展，對全機狀態(tài)，實現(xiàn)了跨聲速顫振仿真［15］，但一直缺乏對跨聲速顫振仿真結(jié)果的驗證。

飛機跨聲速顫振評估的難點在于強的氣動力非線性造成顫振耦合機理比較復雜，目前基于非定常CFD的仿真成熟度不高，基于縮比彈性相似模型的全機跨聲速顫振風洞試驗技術(shù)難度很大，尤其是全機彈性結(jié)構(gòu)相似模型，國內(nèi)在本文之前沒有先例。

針對某全機復雜耦合的顫振現(xiàn)象，本文成功設(shè)計了全機縮比結(jié)構(gòu)相似跨聲速顫振模型，在FL－26風洞開展了全機顫振風洞試驗，同時使用基于N－S方程的CFD仿真程序進行了跨聲速顫振仿真。通過仿真結(jié)果與跨聲速風洞試驗結(jié)果的對比，驗證了一種全機復雜耦合的跨聲速顫振形式。

1 全機跨聲速顫振風洞試驗

1.1 全機結(jié)構(gòu)相似跨聲速顫振模型

中國空氣動力發(fā)展中心FL－26風洞是一座試驗段尺寸為2.4m×2.4m的引射式、半回流、暫沖式跨聲速增壓風洞。該風洞具備一套全機顫振模型的懸浮式支持系統(tǒng)，是目前國內(nèi)唯一具備全機顫振風洞試驗條件的跨聲速風洞。

根據(jù)該風洞試驗段的條件，設(shè)計了全機縮比彈性相似顫振模型，模型主要使用玻璃纖維復合材料＋碳纖維復合材料＋泡沫組成主要結(jié)構(gòu)，局部連接使用鋁，鋼等金屬材料。具體見圖1。

圖1 全機縮比彈性相似顫振模型Fig.1 Full aircraft structuresimilar flutter model

該模型主要結(jié)構(gòu)傳力形式與飛機相似，利用玻璃鋼、碳纖維和泡沫塑料的混合結(jié)構(gòu)，實現(xiàn)了全機縮比彈性結(jié)構(gòu)相似。蒙皮主要使用陰模，采用濕法在常溫下加工分區(qū)變厚度的蒙皮構(gòu)件。機翼上的梁、肋和蒙皮可以組成多閉室薄壁結(jié)構(gòu)，承受機翼的載荷。梁、肋和蒙皮主要依靠粘接形成結(jié)構(gòu)，見圖2。

模型安裝在FL－26風洞的全機懸浮式顫振風洞試驗系統(tǒng)上，模型在風洞中的安裝見圖3。

圖2 全機縮比彈性相似顫振模型的制作Fig.2 Manufacture of the flutter model

1.2 模型的FEM仿真與地面GVT試驗

對模型進行了地面共振試驗（GVT），測量模型的頻率、模態(tài)，在試驗中把模型安裝在風洞懸浮式支持系統(tǒng)上，試驗采用多點隨機激勵的方法進行，見圖4。

對風洞模型進行了結(jié)構(gòu)動力學限元模型（FEM）仿真，仿真模型具有約1000000自由度，見圖5。模型主要頻率、模態(tài)計算結(jié)果見圖6。

圖3 全機顫振模型在風洞中的安裝Fig.3 Full aircraft flutter model in wind tunnel

圖4 模型地面GVT試驗Fig.4 GVT for the flutter model

圖5 結(jié)構(gòu)動力學仿真有限元模型（FEM）Fig.5 FEM of the flutter model

圖6 主要模態(tài)計算結(jié)果Fig.6 Modes of FEM simulation

模型的主要模態(tài)試驗結(jié)果與仿真計算結(jié)果振型基本一致，頻率對比見表1。

表1 仿真計算與GVT結(jié)果Table 1 Result of simulation and GVT

使用亞聲速偶極子網(wǎng)格法（DLM）進行了顫振計算，軟件采用 MSC／NASTRAN－Aeroelastic，非定常氣動力網(wǎng)格見圖7。

圖7 DLM法非定常氣動力網(wǎng)格Fig.7 Mesh of DLM

2 跨聲速顫振CFD仿真計算

采用CFD／CSD耦合的方法對此模型的風洞試驗狀態(tài)進行了仿真，非定常氣動力通過求解Navier－Stokes方程（混合網(wǎng)格）獲得，結(jié)構(gòu)變形由廣義結(jié)構(gòu)運動方程計算。將流體方程與結(jié)構(gòu)運動方程分別構(gòu)造為含子迭代的二階隱式時間推進格式，通過子迭代提高整個計算的時間精度。計算模型的氣動外形和風洞模型外形完全一致。

本文采用的方法里將所有階模態(tài)都進行參數(shù)統(tǒng)計，根據(jù)計算出的廣義位移增長率最大的那階模態(tài)的增長率等于0時判斷為顫振臨界點。而隨著時間步的推進，所有模態(tài)的頻率逐漸靠近為一個統(tǒng)一的頻率，就是顫振頻率。

本文采用含有附面層網(wǎng)格的混合網(wǎng)格進行計算，計算中采用的機翼表面和空間非結(jié)構(gòu)網(wǎng)格如圖8。其中模型表面的三角形網(wǎng)格數(shù)為555190個，附面層為9層網(wǎng)格，空間四面體網(wǎng)格數(shù)為10094920個。

為提高計算中繁重的流體／結(jié)構(gòu)之間的數(shù)據(jù)轉(zhuǎn)換效率，將FEM模型的結(jié)構(gòu)模態(tài)數(shù)據(jù)插值到飛機物面的計算網(wǎng)格點上，通過結(jié)構(gòu)運動方程計算物面的廣義位移也轉(zhuǎn)換到物面的計算網(wǎng)格點上。

對非結(jié)構(gòu)網(wǎng)格，本文采用一種基于彈簧模型的動網(wǎng)格生成技術(shù)，通過迭代得到內(nèi)網(wǎng)格每一網(wǎng)格點的變形量，將變形量疊加在原網(wǎng)格點的坐標值上，得到下一時間步的計算網(wǎng)格。

首先在馬赫數(shù)為0.2時的計算顫振邊界，主要與NASTRAN的DLM方法對比，主要采用變來流速度方法，計算結(jié)果見圖9。然后根據(jù)風洞試驗中的氣流參數(shù)特性完全仿真風洞試驗過程，計算馬赫數(shù)取為0.9，來流速度為281.04m／s，進行變來流密度顫振計算，通過改變來流密度計算得到了振動等幅時的來流密度，見圖10。

圖8 模型CFD網(wǎng)格Fig.8 CFD mesh for the flutter model

圖9 M＝0.2變來流速度的顫振計算結(jié)果Fig.9 Result of CFD simulation using variable airflow speed（M＝0.2）

圖10 M＝0.9變來流密度的顫振計算結(jié)果Fig.10 Result of CFD simulation using variable airflow density（M＝0.9）

3 跨聲速顫振風洞試驗及綜合分析

在FL－26風洞進行了全機模型顫振風洞試驗，風洞試驗采用固定M 數(shù)改變來流速壓的吹風方式進行，通過應變測量和加速度傳感器測量模型的動響應判斷模型是否發(fā)生顫振，模型在風洞中的安裝見圖3。在馬赫數(shù)0.9，速壓48.1kPa，發(fā)生了顫振現(xiàn)象，通過風洞緊急關(guān)車使模型退出了顫振。顫振發(fā)生時的動響應見圖11，所有傳感器都達到同一個振動頻率。

綜合以上的計算和試驗結(jié)果，CFD仿真在M＝0.2時與DLM法的結(jié)果顫振頻率基本一致，顫振速度有一定的誤差，原因可能因為DLM法對于小阻尼型顫振在阻尼為0的點誤差比較大，一般取1%～3%左右的阻尼點作為顫振邊界，具體對比結(jié)果見表2。

圖11 M＝0.9顫振模型動響應Fig.11 Frequency response of flutter model（M＝0.9）

表2 CFD計算與DLM計算結(jié)果的對比（M＝0.2）Table 2 Comparison of CFD and DLM simulation（M＝0.2）

CFD仿真的顫振速壓與風洞試驗結(jié)果進行對比，在M 數(shù)0.9時誤差約為9%左右，這個結(jié)果已經(jīng)能夠滿足工程設(shè)計的需求，具體對比結(jié)果見表3。

表3 CFD計算與風洞試驗結(jié)果的對比（M＝0.9）Table 3 Comparison of CFD simulation and wind tunnel test（M＝0.9）

4 結(jié) 論

通過以上的計算與試驗研究獲得以下結(jié)論：

（1）全機縮比彈性結(jié)構(gòu)相似顫振模型，滿足全機模態(tài)的模擬要求；

（2）基于Navier－Stokes方程（混合網(wǎng)格）的CFD仿真，在采用固定密度，M數(shù)，改變來流速度的條件下，與DLM方法對比，結(jié)果符合比較好；

（3）基于Navier－Stokes方程（混合網(wǎng)格）的CFD仿真，對應風洞試驗狀態(tài)，采用固定M數(shù)，改變來流密度的方法，與風洞試驗結(jié)果符合比較好；

（4）通過全機縮比彈性相似顫振模型的跨聲速顫振風洞試驗以及基于Navier－Stokes方程（混合網(wǎng)格）的CFD仿真，驗證了一種全機復雜耦合的顫振形式，在飛機非定常CFD仿真上取得了進展，可以對全機復雜耦合的顫振特性進行仿真，用于飛機的工程顫振設(shè)計。

［1］GJB 67A7－2008.軍用飛機結(jié)構(gòu)強度規(guī)范第7部分：氣動彈性［S］.

［2］RIVERA J A，F(xiàn)LORANCE J R.Contribution of transonic dynamics tunnel testing to airplane flutter clearance［R］.AIAA 2000－1768，2010.

［3］ELI Livne.Future of airplane aeroelasticity［J］.Journal of Aircraft，2003，40（6）：1065－1092.

［4］FARHAT C，GEUZAINE P，BROWN G.Application of a threefield nonlinear fluid－structure formulation to the prediction of the aeroelastic parameters of an F－16fighter［J］.Computers and Fluids，2003，32：3－29.

［5］FARHAT C，PIERSON K，DEGAND C.Multidisciplinary simulation of the maneuvering of an aircraft［J］.Engineering with Computers，2001，17：16－27.

［6］GEUZAINE P，BROWN G，HARRIS C，et al.Aeroelastic dynamic analysis of a full F－16 con.guration for various flight conditions［J］.AIAA Journal，2003，41（3）：363－371.

［7］MOELJO S Hong.Evaluation of CFL3D for unsteady pressure and flutter predictions［R］.AIAA－2003－1923，2003.

［8］HARTWICH P M，DOBBS S K，ARSLAN A E，et al.Navier－Stokes computations of limit－cycle oscillations for a B－1like configuration［J］.Journal of Aircraft，2001，38（2）：239－247.

［9］ERICSSON L E.Comment on Navier－Stokes computations of limit－cycle oscillations for a B－1like configuration［J］.Journal of Aircraft，2002，39（3）：519.

［10］HARTWICH P M，DOBBS S K，ARSLAN A E，et al.Reply to comment on Navier－Stokes computations of limit－cycle oscillations for a B－1like configuration［J］.Journal of Aircraft，2002，39（3）：519，520.

［11］BRIAN Hayes W，CHARLES E.Goodman.F／A－18E／F Super Hornet flutter clearance program［R］.AIAA 2003－1940，2003.

［12］ANDERSON W D，MORTARA S.F－22 aeroelastic design and test validation［R］.AIAA 2007－1764，2007.

［13］CHEN P C，ZHANG Z C，AYAN S.Aeroelastic applications of a panel－model－based overset Euler solver to a twin－engine transport flutter model［R］.AIAA 2010－2550，2010.

［14］管德，主編.《飛機氣動 彈性力學手冊》［M］.航空工業(yè)出版社，1994.

［15］YANG G W.Numerical analyses of transonic flutter on an aircraft［J］.Chinese Journal of Theoretical and Applied Mechanics，2005，37（6）：769－776.