喻世杰，周興華，黃銳

南京航空航天大學(xué) 機(jī)械結(jié)構(gòu)力學(xué)及控制國家重點實驗室，南京 210016

近年來，對航空航天飛行器執(zhí)行復(fù)雜任務(wù)、適應(yīng)不同飛行環(huán)境的能力要求越來越高，而變體飛行器因具有可根據(jù)任務(wù)類型和飛行環(huán)境的要求自主改變結(jié)構(gòu)外形和氣動布局的優(yōu)點，受到了廣泛關(guān)注。

變彎度機(jī)翼技術(shù)是一類典型的變體形式，可自適應(yīng)改變翼型彎度改善飛行器氣動性能。隨著智能材料與結(jié)構(gòu)技術(shù)的發(fā)展，變彎度機(jī)翼可望成為極具發(fā)展前景的變體技術(shù)。當(dāng)前，國內(nèi)外對變彎度機(jī)翼技術(shù)進(jìn)行了一系列研究。例如，美國國家航空航天局（National Aeronautics Space and Administration，NASA）和美國空軍研究實驗室（Air Force Research Laboratory，AFRL）聯(lián)合開展任務(wù)自適應(yīng)機(jī)翼（Mission Adaptive Wing，MAW）［1-6］項目，改裝了F-111 飛機(jī)的機(jī)翼后緣；在20 世紀(jì)末又啟動了智能機(jī)翼（Smart Wing）［7-8］項目，利用基板彎曲控制后緣偏轉(zhuǎn)；之后研究重心轉(zhuǎn)移到柔性材料與機(jī)械驅(qū)動的組合，開展了自適應(yīng)柔性后緣（Adaptive Compliant Trailing Edge，ACTE）［9］項目，在灣流III 飛機(jī)上應(yīng)用并成功試飛。歐盟FP7（7th Framework Programme）計劃包含了基于可變后緣的增升裝置的研究［10］，瑞士聯(lián)邦材料科學(xué)與技術(shù)研究所的Campanile 通過改變輻條與蒙皮的間隔實現(xiàn)了機(jī)翼的變彎度［11］，德國宇航院（Deutsches Zentrum für Luftund Raumfahrt，DLR）的Sinapius 等［12］提出手指方案對機(jī)翼變彎度實現(xiàn)了單自由度驅(qū)動。中國關(guān)于變彎度機(jī)翼的研究仍處于發(fā)展之中［13］，西北工業(yè)大學(xué)劉世麗等［14］研究了機(jī)翼分布式柔性結(jié)構(gòu)方案，南京航空航天大學(xué)李飛［15］設(shè)計了基于形狀記憶合金的差動驅(qū)動自適應(yīng)機(jī)翼結(jié)構(gòu)，哈爾濱工業(yè)大學(xué)黃建［16］把新型零泊松比蜂窩結(jié)構(gòu)應(yīng)用在變彎度機(jī)翼中。

上述研究主要集中于變彎度機(jī)翼的驅(qū)動和結(jié)構(gòu)設(shè)計等問題，少有關(guān)注變體過程的動力學(xué)問題，尤其是變體過程的氣動力變化及參變氣動彈性穩(wěn)定性問題少有研究［17］。西北工業(yè)大學(xué)倪迎鴿［18］、南京航空航天大學(xué)趙永輝等［19-20］均采用偶極子網(wǎng)格法（Doublet-Lattice Method，DLM）構(gòu)建了折疊機(jī)翼的參數(shù)化氣動彈性模型［21］，分析了折疊角對變體機(jī)翼顫振特性的影響規(guī)律；北京航空航天大學(xué)楊寧等［22］利用DLM 對考慮結(jié)構(gòu)非線性折疊翼進(jìn)行了顫振分析。但DLM 是一種基于平面假設(shè)的方法，該方法無法考慮機(jī)翼彎度影響，更無法研究具體變彎度過程的氣動特性變化規(guī)律，因此DLM 無法用于變彎度機(jī)翼的氣動彈性研究。目前較為通用的方法是基于CFD 的變彎度機(jī)翼非定常氣動力計算分析。例如，西北工業(yè)大學(xué)劉艷［23］研究了采用連續(xù)變彎度后緣舵面系統(tǒng)的彈性翼飛機(jī)氣動彈性計算和結(jié)構(gòu)優(yōu)化問題，北京航空航天大學(xué)吳優(yōu)等［24］基于Open-FOAM 模擬了連續(xù)變彎度翼型的動態(tài)氣動特性。雖然CFD 方法具有高保真度，但對于研究變體過程中引發(fā)的參變氣動彈性穩(wěn)定性方面，存在顫振邊界丟失等問題，且計算效率低，需耗費大量的時間［25-29］。在變彎度機(jī)翼的參數(shù)化建模，國內(nèi)外研究者采用自由型面變形（Free Form Deformation，F(xiàn)FD）技術(shù)實現(xiàn)了變彎度參數(shù)化建模［30-32］，但需耦合基于RANS（Reynolds Averaged Navier-Stokes）方程的CFD 求解器，依舊存在CFD 固有缺陷，也未明確在參變條件下機(jī)翼顫振邊界的變化趨勢及其中的顫振機(jī)理。因此，對于變彎度機(jī)翼的參變氣動彈性問題，迫切需要一種高效、準(zhǔn)確的非定常氣動力計算方法，并耦合參變結(jié)構(gòu)動力學(xué)模型，進(jìn)而高效、準(zhǔn)確預(yù)測變體過程的氣動彈性特性［33］。非定常渦格法（Unsteady Vortex Lattice Method，UVLM）是一種基于時間步長的非定常時域氣動力計算，可動態(tài)捕捉附著渦和尾渦對分布式氣動載荷的影響，能分析機(jī)翼在彎度變化下的參變氣動彈性特征，計算速度較傳統(tǒng)方法有很大提升，是一種兼顧計算效率和精度的折衷方案。對于柔性機(jī)翼的低速非定常氣動力問題，可采用UVLM 進(jìn)行氣動力建模和響應(yīng)求解［34］。

綜上所述，為解決如何高效計算、分析全參數(shù)域下，變彎度機(jī)翼的氣動彈性穩(wěn)定性問題，克服DLM 方法無法考慮變彎度效應(yīng)，而直接CFD計算方法成本過高等困難，本文提出一種基于流形切空間插值和非定常渦格相結(jié)合的連續(xù)變后緣機(jī)翼參數(shù)化氣動彈性建模方法，避免了傳統(tǒng)方法中針對不同參數(shù)需重復(fù)建模的問題，提高了計算效率，且能精確高效地捕捉機(jī)翼在全參數(shù)空間內(nèi)的氣動特性的變化，并分析、預(yù)測了機(jī)翼在翼型彎度變化下的氣動特征變化機(jī)理及趨勢。

1 變彎度機(jī)翼參數(shù)化氣動彈性建模

變彎度機(jī)翼的幾何外形及內(nèi)部結(jié)構(gòu)如圖1 所示，機(jī)翼的后緣變彎度角θ在0°～30（°最大變形狀態(tài)）變化。其中可變后緣部分占機(jī)翼弦長的56%，兩端的柔性引導(dǎo)面可向下延展，以保持變形過程中的翼面連續(xù)性。變形過程中只有可變后緣部分的上下蒙皮發(fā)生彎曲變形，維持了變彎度時的結(jié)構(gòu)和氣動的連續(xù)性，有效地避免了機(jī)翼主要承力結(jié)構(gòu)發(fā)生彈性變形，可在變形過程中保持良好的穩(wěn)定性。

圖1 機(jī)翼模型結(jié)構(gòu)設(shè)計Fig.1 Structural design of wing model

1.1 參數(shù)化結(jié)構(gòu)動力學(xué)建模

設(shè)定初始未變形狀態(tài)的機(jī)翼模型后緣變彎度角θ=0°，建立初始結(jié)構(gòu)模型，對其離散化獲得有限元模型，結(jié)果如圖2 所示。

圖2 機(jī)翼有限元模型（ θ=0°）Fig.2 Finite element model of the wing（ θ=0°）

有限元模型主要采用殼、桿和梁單元構(gòu)建，來傳遞翼面壓力，承受機(jī)翼氣動力和傳遞氣動力矩。網(wǎng)格模型共有1 632 個節(jié)點，考慮翼根固支邊界條件的影響，節(jié)點自由度N=9 588。對模型不同部位賦予不同的材料參數(shù)如表1 所示。

表1 模型參數(shù)Table 1 Parameters of the model

變彎度機(jī)翼運動方程依賴于機(jī)翼的參數(shù)變化，即受后緣變彎度角θ的影響。后緣彎度發(fā)生變化時，建立運動方程所需的質(zhì)量與剛度矩陣都會發(fā)生變化，因此需要對機(jī)翼結(jié)構(gòu)參數(shù)化建模。在傳統(tǒng)的建模方法中，在每個不同變彎度角θ下，都要重新構(gòu)建結(jié)構(gòu)模型，計算量巨大，增加了時間。因此本文利用流形切空間插值（Manifold Tangent Space Interpolation，MTSI）方法，通過有限個不同參數(shù)的初始模型擬合出全部參數(shù)變化條件下的結(jié)構(gòu)動力學(xué)模型［35-37］。

采用節(jié)點坐標(biāo)變換的方式獲取變彎度角分別為θ1=0°、θ2=15°和θ3=30°的可變后緣節(jié)點坐標(biāo)，進(jìn)而得到3 個初始位置的整體有限元模型，結(jié)果如圖3 所示。

圖3 節(jié)點坐標(biāo)變換Fig.3 Coordinate transformation of node

通過有限元軟件分析求解3 個初始位置有限元模型的質(zhì)量矩陣M(θ1)、M(θ2)、M(θ3)和剛度矩陣K(θ1)、K(θ2)、K(θ3)，從而獲得包含前6 階模態(tài)的結(jié)構(gòu)振型Φr(θ1)、Φr(θ2)、Φr(θ3)，再對矩陣降階。

質(zhì)量矩陣的參數(shù)化建模：上文中得到的Mr(θ)均為對稱正定（Symmetric Positive Definite，SPD）矩陣，張成的空間為SPD 流形，記作MSPD。取Mr(θ3)為插值原點M0，MSPD的對數(shù)映射為

式中：logm 為矩陣對數(shù)；Γi為Mr(θi)由MSPD映射到其切空間，記作TM0MSPD的投影點所對應(yīng)的矩陣。

切空間TM0MSPD是“平直”向量空間，在切空間上可以利用Lagrangian 插值法對Γi插值，得到任意變彎度角θk的模型在TM0MSPD中的插值矩陣Γk為

MSPD的指數(shù)映射為

式中：expm 為矩陣指數(shù)；Mr(θk)為Γk由TM0MSPD映射到MSPD的投影點所對應(yīng)的矩陣，即目標(biāo)參數(shù)降階質(zhì)量矩陣。同理可得到目標(biāo)參數(shù)降階剛度矩陣Kr(θk)。MTSI 插值方法的示意圖見圖4。

圖4 流形切空間插值示意圖Fig.4 Schematic diagram of manifold tangent space interpolation

振型矩陣的參數(shù)化建模：由于振型矩陣不處于SPD 流形，引入Φr(θi)T的Moore-Penrose 廣義逆Ψ(θi)作為約束矩陣，滿足

式中：I6為6 階單位矩陣。

Grassmann 流形是y維歐幾里得空間中所有通過原點的x維平面集合，記作G(x，y)，Φr(θ)的第j列（即第j階(j=1，2，…，6)模態(tài)對應(yīng)的振型）張成的特征子空間可以視作流形G(1，N)上一點。取Φr(θ3)的第j列作為插值原點Φ0j，G(1，N)的對數(shù)映射為

式中：Φij表示振型矩陣Φr(θi)對應(yīng)的第j列。

同理利用式（3）得到變彎度角θk下的振型矩陣單列在切空間中插值結(jié)果Γkj。G(1，N)指數(shù)映射為

采用Gram-Schmidt 正交化方法使得Φkj滿足正交約束條件

式中：Ψkn表示約束矩陣Ψ(θk)對應(yīng)的第n列。

利用同樣的方法可以得到滿足正交約束條件的Ψkj，且有=1。最后，通過縮放使Φkj和Ψkj的2 范數(shù)大小與插值原點保持一致，得到參數(shù)化振型矩陣Φr(θk)=[Φk1Φk2…Φk6]。后文所用的振型Φ均為截斷振型，省去下標(biāo)r。

1.2 非定常氣動力建模

如圖5 所示，基于非定常渦格法，劃分機(jī)翼升力面的氣動網(wǎng)格并標(biāo)注非定常參數(shù)。

圖5 機(jī)翼升力面非定常參數(shù)Fig.5 Unsteady parameters of wing lift surface

位于(i，j) 處的渦環(huán)Q(i，j，Γij) 在任意點P(x，y，z)引起的誘導(dǎo)速度為

式中：fiv為誘導(dǎo)速度計算函數(shù)。

利用2 個相鄰時間步的后緣渦格后節(jié)點的位置創(chuàng)建尾渦渦格，重復(fù)此過程以模擬尾渦脫落過程如圖6 所示。

則有

式中：Γwt為當(dāng)前時間步脫落的尾渦強度；Γtet-Δt為前一時間步后緣附著渦強度。

在第K個控制點（即第i，j個渦格）上指定無穿透邊界條件為

式中：n為渦格法向量；m為附著渦總數(shù)；U∞為來流速度；vη為升力面在該控制點的變形速度；(uw，vw，ww) 為尾渦誘導(dǎo)速度（通過式（11）和式（12）計算并求和得出）。則第K個附著渦的渦強為

由伯努利方程可得分布壓力差的計算公式：

式中：ρ為空氣密度；Δci，j和Δbi，j分別是對于第i，j個渦格的i和j方向的長度（見圖5）；τii，j和τji，j分別是對于第i，j個渦格的i和j方向的切向向量。

最終得到分布?xì)鈩虞d荷：

式中：ΔSij為第i，j個渦格的面積。

與傳統(tǒng)氣動網(wǎng)格平面相比，UVLM 方法不再局限于平面無彎度翼型，為將參數(shù)化結(jié)構(gòu)動力學(xué)建模和UVLM 計算方法相匹配，需建立一個可隨參數(shù)變化的氣動網(wǎng)格模型，本文選取參數(shù)變化后結(jié)構(gòu)的中弧面。利用排序與分類算法提取出有限元結(jié)構(gòu)的上下表面節(jié)點順序序列，并根據(jù)其順序重新構(gòu)造上下曲面網(wǎng)格，得到當(dāng)前變彎度角下的氣動網(wǎng)格模型。

獲得的各變彎度角下的氣動網(wǎng)格模型如圖7所示。

圖7 不同變彎度角下的氣動模型Fig.7 Aerodynamic model with different camber angles

1.3 流-固耦合插值建模

對機(jī)翼結(jié)構(gòu)進(jìn)行有限元分析可獲得結(jié)構(gòu)網(wǎng)格節(jié)點處的模態(tài)振型，通過非定常氣動力建模方法可獲得渦格控制點的分布?xì)鈩虞d荷。為進(jìn)一步得到氣動網(wǎng)格節(jié)點處的模態(tài)振型及其等效空氣動力，利用無限板樣條插值（Infinite Plate Spline，IPS）方法來實現(xiàn)對于模態(tài)振型和氣動力的插值［38］。即對于結(jié)構(gòu)模態(tài)振型向量us和氣動模態(tài)振型向量ua，由載荷向量的平衡關(guān)系可得樣條插值矩陣G，滿足的條件為

針對氣動模型模態(tài)振型求解：通過式（19）構(gòu)造出機(jī)翼結(jié)構(gòu)的中弧面網(wǎng)格，對于任意節(jié)點(xk，yk)都能在機(jī)翼結(jié)構(gòu)網(wǎng)格的上下表面分別找到與之唯一對應(yīng)節(jié)點(xi，yi)與(xi＇，yi＇)，有xk=xi=xi＇，yk=yi=yi＇。因此可得中弧面的模態(tài)振型：

式中：Φs為有限元結(jié)構(gòu)振型。

因為氣動面與中弧面構(gòu)造方法一致，所以G為單位矩陣，通過式（22）計算氣動網(wǎng)格的模態(tài)振型為

針對速度和氣動載荷的求解，由式（20）得到氣動網(wǎng)格節(jié)點向渦格控制點插值的樣條插值矩陣Ga。設(shè)氣動網(wǎng)格節(jié)點的速度矩陣為va，則有

式中：vc為控制點的速度矩陣；Fa為氣動網(wǎng)格節(jié)點的氣動力矩陣。

1.4 參數(shù)化氣動彈性建模

前面已經(jīng)獲得了任意變彎度角θ下機(jī)翼有限元結(jié)構(gòu)模型的降階質(zhì)量矩陣、降階剛度矩陣及截斷振型。忽略阻尼的影響，構(gòu)建參變氣動彈性方程的一般形式為

式中：Fg為廣義氣動力，其表達(dá)式為

在參變條件下，式（11）～（式18）給出了非定常氣動力建模方法、式（22）～式（24）給出了流-固耦合插值建模方法、式（25）～式（27）給出了氣動彈性建模法，綜合上述方法計算不同變彎度角的機(jī)翼模型在不同來流速度下時域響應(yīng)。通過每一時刻的模態(tài)位移響應(yīng)可構(gòu)造得隨時間變化的氣動網(wǎng)格及其渦格，通過每一時刻的模態(tài)速度響應(yīng)可計算隨時間變化的誘導(dǎo)速度，這兩項的改變造成了氣動力的變化。最后利用穩(wěn)定狀態(tài)時域響應(yīng)的斂散性來確定機(jī)翼的顫振邊界，分析不同變彎度角對于機(jī)翼顫振特性的影響及其規(guī)律，實現(xiàn)了在給定風(fēng)速條件下對機(jī)翼模型的時域響應(yīng)分析，整體流程圖如圖8 所示。

圖8 氣動彈性建模流程圖Fig.8 Flow diagram of aeroelastic modeling

2 數(shù)值仿真與驗證

2.1 變彎度機(jī)翼參變結(jié)構(gòu)動力學(xué)特性與驗證

為分析不同變彎度角下機(jī)翼的結(jié)構(gòu)動力學(xué)特性，利用流形切空間插值方法，在0°～30°間隔2.5°選取共13 個后緣彎度位置的參數(shù)化矩陣，以此獲得不同變彎度角下模型的前6 階固有頻率及固有振型，并利用節(jié)點坐標(biāo)變換方式直接計算得到相應(yīng)的頻率和振型作為對比。固有頻率對比結(jié)果如圖9 所示，變化曲線基本一致擬合效果較好，最大誤差出現(xiàn)在θ=5°的第2 階模態(tài)處為3.6%。

圖9 固有頻率對比Fig.9 Comparison of natural frequencies

θ=5°時的固有振型的對比結(jié)果如圖10 所示（黑色為未變形狀態(tài)，紅色為MTSI 插值模態(tài)振型，藍(lán)色為直接計算模態(tài)振型）。由圖10 可知：由于機(jī)翼展弦比較大，前6 階模態(tài)中包含面內(nèi)運動（一階面內(nèi)彎曲），不過面內(nèi)振動幾乎不引起非定常氣動力，對機(jī)翼的氣動彈性特性影響不大。

圖10 MTSI 方法計算的固有振型Fig.10 Natural mode shapes by MTSI method

模態(tài)置信準(zhǔn)則（Modal Assurance Criterion，MAC）值是振型向量之間的點積，用于檢驗?zāi)B(tài)振型的相關(guān)性。不同變彎度角下參數(shù)化模態(tài)振型與直接計算模態(tài)振型的MAC 值如表2 所示，作為依據(jù)驗證參數(shù)化模態(tài)振型與直接計算模態(tài)振型的誤差。

表2 固有振型MAC 值Table 2 MAC value of natural mode shapes

由圖10 和表2 可得：參數(shù)化模態(tài)振型與直接計算模態(tài)振型吻合較好，在高階模態(tài)（第6 階）及θ=5°的第2 階模態(tài)處存在一定誤差，與圖9 中的頻率誤差分布相近，但數(shù)值均不低于0.996，滿足工程指標(biāo)要求，驗證了參數(shù)化結(jié)構(gòu)動力學(xué)建模的準(zhǔn)確性。

2.2 流-固耦合插值驗證

為了驗證IPS 插值方法的準(zhǔn)確性，本文將θ=5°的IPS 插值氣動模態(tài)振型與圖10 中的MTSI 結(jié)構(gòu)模態(tài)振型進(jìn)行對比驗證，結(jié)果如圖11所示（紅色為IPS 插值氣動模態(tài)振型，藍(lán)色為MTSI 結(jié)構(gòu)模態(tài)振型）。由圖11 可知，插值得到的氣動網(wǎng)格模態(tài)振型與參數(shù)化結(jié)構(gòu)模態(tài)振型基本一致，驗證了插值方法的準(zhǔn)確性。

圖11 插值模態(tài)振型對比Fig.11 Comparison of interpolated mode shapes

2.3 非定常氣動力驗證與計算

為驗證非定常氣動力計算方法的準(zhǔn)確性，采用Theodorsen 非定常氣動力理論與本文非定常氣動力方法比較，計算二元翼段俯仰角非定常氣動力，簡諧運動為

翼段具體參數(shù)見表3。利用本文方法建立相同的氣動模型，增大展弦比以忽略機(jī)翼三維效應(yīng)的影響，其余參數(shù)同表3。對比結(jié)果如圖12所示。

表3 翼段參數(shù)Table 3 Parameters of wing segments

圖12 Theodorsen 與UVLM 非定常氣動力對比Fig.12 Comparison of unsteady aerodynamic forces calculated by Theodorsen theory and UVLM

由圖12 可知，UVLM 算法計算得到的非定常氣動力與Theodorsen 氣動力吻合良好。為進(jìn)一步研究UVLM 算法對氣動力的計算效果，本文在給定條件下（不考慮彈性變形，給定風(fēng)速10 m/s，初始迎角0°），計算各變彎度角下機(jī)翼的定常氣動力系數(shù)的展向分布及隨彎度角變化的總體氣動力系數(shù)的變化情況，結(jié)果如圖13 和圖14 所示。

圖13 升力系數(shù)展向分布Fig.13 Spanwise distribution of lift coefficient

圖14 總升力系數(shù)隨彎度角的變化Fig.14 Variation of total lift coefficient with camber angle

2.4 顫振計算方法驗證

為驗證本文所用顫振計算方法的準(zhǔn)確性，構(gòu)建與機(jī)翼模型投影相似的平面機(jī)翼模型（忽略翼型的影響以滿足驗證算法的條件），分別利用本文方法與偶極子網(wǎng)格法進(jìn)行顫振分析對比。

構(gòu)建圖15 所示的機(jī)翼模型。有限元單元僅使用殼單元，邊界條件為根部固支。將結(jié)果提交至軟件中利用偶極子網(wǎng)格法計算，根軌跡結(jié)果如圖16（a）所示（沿箭頭方向來流速度逐漸增大），顫振速度18.4 m/s，顫振頻率3.15 Hz。

圖15 平面機(jī)翼網(wǎng)格Fig.15 Mesh of wing in plane

圖16 DLM 與UVLM 根軌跡對比Fig.16 Comparison of root locus by DLM and UVLM

再利用本文方法獲得該模型在顫振點附近的時域響應(yīng)。為在頻域中直觀地分析顫振速度與顫振模態(tài)，通過動態(tài)模態(tài)分解（Dynamic Mode Decomposition，DMD）方法［39］，辨識不同風(fēng)速下系統(tǒng)的特征值并排序，得到機(jī)翼各階模態(tài)的阻尼比和頻率隨來流速度變化情況，同樣繪制出根軌跡，結(jié)果如圖16（b）所示，顫振速度18.8 m/s，顫振頻率3.15 Hz。結(jié)果表明，本文方法計算所得的機(jī)翼顫振特性與偶極子網(wǎng)格法計算所得的吻合良好。

2.5 參變顫振計算

對不同變彎度角下機(jī)翼模型的顫振邊界計算與分析。對于給定變彎度角，判斷時域模態(tài)響應(yīng)在不同風(fēng)速下的斂散性，找出此時模型的顫振區(qū)間上下邊界如圖17 所示（以θ=5°的時域響應(yīng)為例），以此為基準(zhǔn)可進(jìn)一步確定該角度下模型的精確顫振速度。

圖17 θ=5°機(jī)翼顫振區(qū)間Fig.17 Flutter interval at θ=5°

通過上文所述的DMD 方法，對參變條件下的時域響應(yīng)進(jìn)行頻域分析，求得風(fēng)速變化條件下不同變彎度角的根軌跡，以此分析參變條件模型顫振前后的模態(tài)特征變化，結(jié)果如圖18 所示。由圖18 可知，機(jī)翼模型在各彎度角下的顫振均由1 階扭轉(zhuǎn)模態(tài)失穩(wěn)引起。對于1 階面內(nèi)彎曲模態(tài)，由于其面內(nèi)振動的特性，在顫振被激發(fā)時并不參與振動，根軌跡均位于左半平面且無序分布，因此圖18 中并未畫出。對于其他階模態(tài)，隨著風(fēng)速的增加，阻尼比數(shù)值均在減小，且大部分模態(tài)頻率向顫振頻率靠近。

圖18 參變根軌跡Fig.18 Root locus with parameter variation

為進(jìn)一步分析各變彎度角下顫振發(fā)生時的模態(tài)成分，本文利用不同變彎度角下模型在各自臨界顫振狀態(tài)下的模態(tài)響應(yīng)做頻譜分析。取穩(wěn)定狀態(tài)的時域響應(yīng)數(shù)據(jù)進(jìn)行快速傅里葉變換處理，結(jié)果如圖19 所示。

圖19 不同變彎度角下各階模態(tài)位移的頻譜分析Fig.19 Spectral analysis of each modal displacement at different camber angles

圖19 中均出現(xiàn)了2 個明顯峰值，較小值與1 階面內(nèi)彎曲模態(tài)的固有頻率接近，較大值為該變彎度角下的顫振頻率。由圖19 分析可知，在后緣偏轉(zhuǎn)角度較小時，2 階彎曲和1 階扭轉(zhuǎn)模態(tài)為主要顫振模態(tài)，除1 階面內(nèi)彎曲外各階模態(tài)的主要頻率均為顫振頻率。隨著后緣彎度角的增大，1 階扭轉(zhuǎn)模態(tài)成分占比逐漸增大，1、2 階彎曲模態(tài)成分占比逐漸減小。隨著可變后緣的偏轉(zhuǎn)，機(jī)翼后緣段的定常氣動升力逐漸增大，機(jī)翼的顫振形式由彎曲向扭轉(zhuǎn)過度，所以1 階扭轉(zhuǎn)模態(tài)對顫振的影響逐漸增大，1、2 階彎曲的影響逐漸減小。機(jī)翼模型的顫振速度和顫振頻率隨彎度的變化如圖20 和圖21 所示。

圖20 顫振速度隨彎度角的變化Fig.20 Variation of flutter speed with camber angle

圖21 顫振頻率隨彎度角的變化Fig.21 Variation of flutter frequency with camber angle

仿真驗證表明顫振速度擬合效果較好，顫振頻率擬合存在一定差異。在低偏轉(zhuǎn)角時，隨著后緣彎度的增大，機(jī)翼顫振速度逐漸減小，顫振頻率逐漸增大，機(jī)翼可變后緣部分的翼型彎度逐漸增大，機(jī)翼的定常氣動力變大，使其顫振邊界降低。同時由于后緣偏轉(zhuǎn)，增大了氣動力對機(jī)翼扭轉(zhuǎn)的影響，使機(jī)翼1 階扭轉(zhuǎn)模態(tài)對顫振的貢獻(xiàn)度增大，進(jìn)而降低了顫振速度，顫振頻率也向1 階扭轉(zhuǎn)模態(tài)頻率靠近，提高了顫振頻率。在高偏轉(zhuǎn)角之后，1 階扭轉(zhuǎn)模態(tài)完全主導(dǎo)顫振，顫振速度也不再因彎度角的增大而產(chǎn)生較大變化，解釋了此時顫振頻率變化趨于平緩的原因。

3 結(jié) 論

針對變彎度機(jī)翼的參數(shù)化建模以及顫振計算的問題，提出了基于流形切空間插值和非定常渦格法相結(jié)合的變彎度機(jī)翼參數(shù)化氣動彈性建模方法，高效準(zhǔn)確地分析了變彎度機(jī)翼的參變結(jié)構(gòu)特性和氣動彈性力學(xué)特性，并對機(jī)翼模型的氣動性能、非定常氣動力計算、插值方法、顫振計算均進(jìn)行了數(shù)值仿真與驗證，主要結(jié)論如下：

1）實現(xiàn)了連續(xù)變彎度機(jī)翼在全參數(shù)空間內(nèi)的精確化建模，包含結(jié)構(gòu)參數(shù)化建模、非定常氣動力建模、流-固耦合插值和氣動彈性模型建模?？梢詼p少不同參數(shù)模型重復(fù)構(gòu)建的次數(shù)，降低結(jié)構(gòu)模態(tài)計算的時間消耗。同時準(zhǔn)確地模擬機(jī)翼在不同外形參數(shù)下的結(jié)構(gòu)動力學(xué)特性，可在時域分析中綜合考慮機(jī)翼彎度、機(jī)翼彈性變形帶來的影響，對參變氣動性能進(jìn)行有效分析。

2）利用該方法對變彎度機(jī)理分析，分別從時域和頻域兩個角度深入研究了不同彎度角下機(jī)翼的參變特性。探討了機(jī)翼隨彎度角變化引起的參變顫振特性，計算并分析了不同變彎度角模型的顫振邊界和顫振模態(tài)成分，研究了機(jī)翼隨變彎度角變化引起的模態(tài)特征的變化趨勢及其原因：隨著變彎度角的增大，機(jī)翼的顫振速度增大而顫振頻率減小，但在一定角度后，由于1 階扭轉(zhuǎn)模態(tài)的主導(dǎo)性，機(jī)翼顫振受其他各階模態(tài)的影響較小，導(dǎo)致顫振速度的變化出現(xiàn)閾值，速度和頻率的改變趨于平穩(wěn)。

3）本文的參數(shù)化建模方法可擴(kuò)展應(yīng)用到其他結(jié)構(gòu)可變機(jī)翼構(gòu)型（例如可折疊機(jī)翼、變后掠機(jī)翼等），在參變條件下對飛行器進(jìn)行氣動彈性建模與分析，具有貼合工程實際的普適性意義。但本文提出的參數(shù)化建模方法仍有一些值得深入研究的問題。例如，仍受到氣動力模型和計算效率的限制，仍無法考慮翼型厚度的影響等。此外，需要深入探索該方法對于時變問題的應(yīng)用，進(jìn)而為變體飛行器的時變氣動彈性力學(xué)研究提供高效、準(zhǔn)確的建模方法。