平尾模型連續變速壓顫振試驗方法及數值計算研究

劉 南，易家寧，王 冬，侯良學，楊希明

(1.中國航空工業空氣動力研究院 氣動研究與試驗二部，沈陽 110034；2.高速高雷諾數氣動力航空科技重點實驗室，沈陽 110034)

顫振是一種典型的氣動彈性不穩定性現象，可能導致飛行器的結構遭到災難性破壞，在飛行包線內絕對不允許發生顫振現象且留有一定的裕度。風洞試驗是研究飛行器顫振問題的重要手段之一，例如我國國軍標GJB 67.7A—2008《軍用飛機結構強度規范第7部分》[1]中明確要求飛機必須通過顫振風洞試驗獲得其顫振特性，進而得出顫振設計結論和飛行包線。而跨聲速條件又由于流場非線性、顫振速壓低、試驗風險高，成為顫振風洞試驗研究的重點和難點[2]。

相比顫振計算研究[3-5]，目前國內跨聲速顫振試驗研究主要受到如下兩方面的約束：① 國內絕大部分高速風洞都是暫沖式，風洞試驗研究大多基于該類型風洞[5-7]，但是暫沖式顫振風洞試驗信號有效采樣時間短、風洞背景噪聲大、信噪比低，而且采用階梯變速壓方式進行顫振試驗耗氣量大，試驗效率低；② 試驗模型跨聲速顫振邊界評估手段不足，基于偶極子格網法的線化方法無法準確預測跨聲速顫振邊界，一般通過增加安全系數的方式進行處理，結果的精度和可靠性無法保證；而基于雷諾平均Navier-Stokes方程的高精度方法效率較低。

針對上述問題，本文基于自回歸滑動平均系統辨識方法建立了一套高效的模型顫振邊界預測方法，通過地面振動試驗(ground vibration test，GVT)結果對模態數據進行修正后，計算與風洞試驗吻合較好，顯著地提高了試驗的安全性和效率；基于補償解耦和專家PID建立了高速暫沖式風洞定馬赫數連續變速壓試驗方式，有效地提高了試驗效率。

1 試驗模型顫振邊界高效預測方法

1.1 ARMA辨識方法

首先簡單介紹本文采用的ARMA系統辨識方法[8]。該方法中廣義氣動力預測模型如下：

(1)

式中：fT為廣義氣動力向量；q為模態位移向量；na和nb分別為輸出、輸入延遲階數；k為離散時間步；Ai，Bi∈RN×N為ARMA模型參數，N為模態階數。

利用3211外部信號作為模態位移的激勵，在來流速壓QT條件下通過非定常計算得到廣義氣動力，進而獲得時域響應序列{q(k),ya(k)}(k=1,2,…,K)，其中K為迭代步數。通過最小二乘法可以得到矩陣{Ai}(i=1,2,…,na)和{Bi}(i=0,1,…,nb)，具體計算過程此處不再贅述，可以參考文獻[9-10]。

將式(1)改寫成矩陣形式：

(2)

式中各個量的具體定義可以參考文獻[4-5]。

將式(2)轉化為連續時域狀態空間方程如下：

(3)

結合結構模態運動方程：

(4)

式中：M,D,K分別是廣義質量、阻尼和剛度矩陣；q為模態位移，廣義力f與fT之間的關系如下：

(5)

式中，Q為來流動壓。結合式(3)～(5)可得：

(6)

其中

根據Lyapunov穩定性理論，根據矩陣L的特征值分布能夠判斷式(6)的穩定性，進而獲得試驗模型的顫振邊界。

1.2 驗證算例

利用AGARD445.6機翼對ARMA辨識方法進行驗證。該機翼于1961年在TDT風洞進行顫振試驗[11]，已成為跨聲速顫振研究和程序驗證的標準模型之一。

CFD/CSD計算結果(無黏和有黏)、ARMA辨識結果以及風洞試驗結果對比如圖1所示。由圖可見：CFD/CSD計算考慮黏性后，大幅提升與風洞試驗結果的吻合度，尤其是Ma=1.041和1.172工況；對于AGARD445.6機翼，ARMA預測結果與CFD/CSD計算結果較為吻合，顫振頻率和動壓的差別均不超過4%。效率方面，ARMA預測顫振邊界僅需一次非定常計算，無需在不同動壓下分析時域響應變化趨勢，能夠顯著提升顫振分析效率。

(a)顫振速壓

2 定馬赫數連續變速壓流場控制

FL-60暫沖式風洞流場輸入為主調壓閥開度、柵指位移和引射器調壓閥開度，輸出為總壓和馬赫數。其中引射器用于降低速壓運行下限，其運行方式在主調壓閥開啟前預先調節到特定引射壓力，在風洞主管路流場調節過程中，引射壓力保持恒定，可視作常量。因此，對于亞、跨聲速(Ma=0.4至1.2)流場控制系統可視作兩個存在耦合關系的閉環控制回路：① 總壓控制回路，以調壓閥為執行機構，被控對象是總壓；② 馬赫數控制回路，以柵指為執行機構，被控對象是馬赫數。其中，總壓調節過程中閥門開度對馬赫數產生耦合影響，馬赫數調節過程中柵指位移對總壓產生耦合影響。為了消除兩者之間的相互影響，首先通過補償解耦消除耦合影響，然后設計專家PID建立各個控制回路。

2.1 補償解耦

風洞總壓P0和馬赫數Ma的參數表達式為

(7)

式中：smv和scf分別是主調壓閥開度和柵指位移。由流場調試結果得出，柵指位移對總壓的影響近似視作為零，則式(7)可寫成如下矩陣形式：

(8)

引入補償因子并設計反饋解耦器，如下：

(9)

式中：λ12、λ21為補償因子。若實現完全解耦，則使系統矩陣對角線為0，即

(10)

聯立式(9)和(10)，有

(11)

但是風洞流場復雜且具有非線性特征，難以精確得到系統模型參數G11至G22，無法實現完全解耦。在實際工程應用中，往往只需要將耦合作用減弱到一定程度，便可以將耦合視作擾動，通過閉環控制消除耦合影響。

2.2 專家PID控制器

連續變速壓試驗需要在保證調節速度和調節精度符合需求的前提下盡量避免產生超調，但同一組PID參數的暫態性能和穩態性能相悖，也難以在不同的速壓范圍內均具備理想的控制效果。此外，風洞流場具有非線性、滯后和時變不確定性等特征，難以獲取精確的數學模型。只能將減弱之后的耦合視作擾動，并通過優化控制器性能來消除耦合帶來的不利影響。諸多因素使得常規的PID控制器難以滿足要求。

FL-60風洞所采用的控制器形式為反饋-前饋結合的符合控制器，其中反饋部分為增量式數字PI控制器，用于穩定被控參數；前饋部分為靜態前饋控制器，用于補償風洞氣源壓力變化，表達式為

u(k)=u(k-1)+KP[e(k)-e(k-1)]+

KIe(k)+KFFTS

(12)

式中：KFF為前饋補償系數；Ts為控制周期。

專家PID控制基于受控對象和控制規律的各種知識，根據專家規則、當前狀態和歷史狀態等智能調整PID控制器參數，使受控系統盡可能地優化和實用化，具備一定的自適應性和魯棒性。

在式(12)的基礎上引入各參數的加權因子，則控制算法如下：

u(k)=u(k-1)+αKP[e(k)-e(k-1)]+

βKIe(k)+ξKFFTS

(13)

式中：α、β、ξ分別為比例系數、積分系數和微分系數。

以總壓專家PID控制器為例，控制器結構見圖2。加權因子的大小由專家系統按照專家規則，根據被控參數的歷史輸入、輸出和其他相關條件進行在線修正。為避免過度修正引起失穩，加權因子取值范圍為α,β,ξ∈[0,2]。

圖2 風洞總壓專家PID控制器

2.3 試驗結果

在FL-60風洞開展定馬赫數連續變速壓流場調試試驗，分別在Ma=0.8和1.2工況下進行調試，設定總壓變化速率為3 kPa/s，對應速壓變化速率分別為0.88和1.2 kPa/s，得到來流總壓、馬赫數和速壓的變化趨勢如3所示。由圖可見：馬赫數控制誤差小于±0.003，取得較為理想的控制效果。

3 顫振風洞試驗

3.1 試驗模型

本文風洞試驗采用模型是依據某機平尾為基礎，按照FL-60風洞速壓運行范圍設計得到的縮比模型。幾何尺寸如下：展長650.6 mm，翼根弦長659.1 mm，翼根弦長220.3 mm，前緣后掠角34°。結構采用“承力蒙皮+梁”的結構形式，蒙皮、肋板和梁等主要結構形式采用玻璃纖維復合材料，翼根采用高強度鋁合金7075材料與轉軸連接，內部填充聚氨酯泡沫，翼梢加鉛塊配重并安裝加速度計測量振動響應。

(a)Ma=0.8

開展GVT試驗通過錘擊法獲得試驗模型安裝在風洞轉窗后的固有頻率和振型，如表1所示為計算與GVT試驗得到的固有頻率對比。由結果可見：前四階固有頻率誤差均在10%以內，有限元計算與GVT試驗符合性較好；但是前兩階固有頻率與GVT試驗結果存在一定差距，而且這兩個模態是顫振主導模態，因此需要在計算前對結構有限元模型進行修正，可以通過修改模型參數匹配GVT測得的固有頻率和結構振型。

表1 計算與GVT試驗得到固有頻率對比(以一階固有頻率計算值為1進行無量綱化處理)

3.2 信號采集和安全防護

試驗模型安裝在FL-60風洞試驗段如圖4所示，試驗過程中采用PXI動態系統進行振動信號采集，該系統是基于PXI總線的QTS2524動態測試系統，由CM4214放大器、CM4504低通濾波器、NI-PXI4472動態采集設備組成，系統共有24個通道，每個通道動態采樣率最高51.2 kHz，動態采集精度0.2%。

圖4 試驗模型安裝在FL-60風洞試驗段上

顫振試驗數據采集和安全防護軟件是基于PXI動態采集設備和Labview軟件環境開發的，如圖5所示，主要包括操縱面板、采集參數設置以及各個通道的時域響應。其中操縱面板提供人機交互功能，包括開始/停止數據采集、強行關車、顯示試驗狀態等功能，采集參數設置主要供使用人員在試驗前設置采樣率、各個采樣通道、安全防護監控通道和閾值等。

圖5 顫振試驗數據采集和安全防護軟件界面

試驗模型安全防護裝置如圖6所示。試驗模型根部安裝限位裝置，在風洞啟動、關車和模型即將發生顫振(連續多次加速度振幅超過預先設定閾值)等條件下，通過氣缸驅動限位裝置限制模型的扭轉位移，從而抑制顫振，保護模型。

圖6 試驗模型安全防護裝置

3.3 試驗方法

本文顫振試驗通過定馬赫數階梯變速壓和連續變速壓兩種試驗方法開展。前者在某個固定馬赫數條件下設置若干速壓階梯，每個階梯保持穩定5～10 s并記錄模型的振動信號；后者則是在某個固定馬赫數條件下速壓線性增大。試驗結束后，通過亞臨界響應分析方法分析振動信號，并外插得到顫振邊界。

本文亞臨界響應分析方法選擇Peak-Hold方法，該方法在國內外高速顫振風洞試驗中得到廣泛應用[12]。其原理如下：對于彈性系統置于風洞氣流中受到湍流激勵，系統的振幅可以表示為顫振微分方程的解，表示如下：

(14)

式中：Q為廣義力；A為廣義氣動力(是馬赫數和減縮頻率的函數)；V為來流速度；ρ為來流密度；下標f、s和t分別代表顫振、振動和湍流。

對上式左右兩邊求倒數，并將振動和湍流項分開，對于湍流激勵可得：

(15)

因此，當系統接近顫振臨界點時，1/|a|與速壓接近線性關系，因此可以利用譜峰值倒數進行顫振邊界外插。

3.4 試驗結果分析

針對顫振風洞試驗結果，首先將流場和翼梢加速度時域響應(由于加速度信號信噪比較低，本文利用低通濾波器已將200 Hz以上的高頻信號濾除)分段加窗，如圖7所示，然后計算各段的平均速壓及對應的加速度頻譜，畫出譜峰值隨速壓的變化趨勢及其擬合曲線如圖8所示，預測得到顫振邊界43.0 kPa，其中圖8(a)是3車階梯變速壓試驗結果的綜合展示，試驗過程中每車最多3個速壓階梯。

(a)馬赫數和速壓

圖8 譜峰值倒數隨速壓的變化趨勢及擬合得到的線性趨勢線(Ma=0.70)

綜上，通過本次試驗分別獲得兩種試驗方式得到的平尾模型顫振邊界如表2所示。結果表明：兩種試驗方式得到的顫振邊界非常接近，差別不超過2%。

表2 平尾模型兩種試驗方式得到的顫振速壓對比

3.5 顫振邊界對比

顫振速壓邊界計算與試驗結果對比如圖9所示，由圖可見：利用原始有限元模型進行顫振分析得到的顫振邊界較高，利用GVT試驗結果對有限元模型修正后得到顫振邊界與試驗值較為吻合。

圖9 平尾模型顫振速壓和頻率對比

由本文結果可見，對于一般的部件級顫振模型，通過數值計算在風洞試驗前可以得到較為可靠的顫振邊界。如果采用定馬赫數階梯變速壓試驗方式，根據計算結果在來流速壓遠低于顫振速壓時，可以適當增大速壓階梯，在靠近顫振邊界時適當減小速壓階梯。如果采用定馬赫數連續變速壓試驗方式，可以進一步改進設計分段線性速壓變化流場控制程序，在來流速壓遠低于顫振速壓時，可以適當增大速壓變化速率，在靠近顫振邊界時適當減小速壓變化速率。

4 結 論

針對跨聲速顫振裕度低、風險高、預測難度大的問題，本文建立了一套高效高精度顫振預測方法和連續變速壓試驗方法，得出結論如下：

(1)首先建立了基于CFD/CSD耦合的高精度顫振預測方法，但是該方法需要在時域內求解非定常流場和結構運動方程并且反復迭代，效率極低，在此基礎上，建立了基于ARMA系統辨識模型的高效顫振預測方法，每個馬赫數條件下僅需進行一次非定常計算即可，效率提升一個量級左右；

(2)通過AGARD445.6標模對CFD/CSD耦合計算方法和ARMA辨識方法進行驗證，計算結果與試驗吻合度很高；

(3)基于FL-60風洞建立了定馬赫數連續變速壓試驗能力，通過流場調試證明速壓變化線性度高、馬赫數偏差小(波動范圍小于±0.003)；

(4)通過平尾模型對上述成果進行綜合驗證，階梯變速壓和連續變速壓兩種試驗方式顫振速壓差別低于2%，計算結果與試驗吻合度較好，為試驗車次安排可以起到良好的指導作用。

以上試驗模型的設計結果固有頻率和結構振型與GVT試驗有一定差別，后續針對上述差異進一步研究，探索復合材料有限元分析、鋪層設計和加工工藝等影響；該試驗模型的跨聲速凹坑附近顫振速壓過低，無法通過暫沖式風洞試驗獲得，后續通過修改模型或改至連續式風洞(可負壓運行)進行試驗，獲得較為完整的顫振邊界，為數值計算驗證提供更加豐富的試驗數據。