基于LPV模型的顫振試飛響應預測及激勵優化

寇寶智， 雷 鳴， 盧曉東

(中國飛行試驗研究院 飛機所，西安 710089)

顫振試飛是一種高風險試飛課目，保障試飛安全和提升試飛效率是工程研究和實施的主要方向。試飛安全監控上普遍采用響應的時域監控和準實時數據處理結合的方法，時域監控主要關注響應的收斂情況，準實時數據處理一般在激勵動作結束后，通過遙測數據進行模態參數辨識，獲取阻尼變化情況。為了獲取良好的響應數據，顫振試飛需要采用有效的激勵手段，將線性掃頻激勵信號通過電傳系統施加到操縱面上是目前大部分電傳飛機顫振試飛的主要激勵手段之一[1]。

雷鳴等[2-3]提出了全機動力學有限元模型為基礎的激勵響應預測方法，在民機氣動伺服彈性和顫振試飛安全監控上得到應用。該方法屬于從設計角度出發的正向方法，依賴有限元模型的準確性。頻域子空間的方法由Mckelvey等[4-6]提出并進行了各種改進和擴展，這種面向控制的多入多出辨識模型，以其簡潔的狀態空間模型優勢，得到了各種系統建模中的廣泛應用，也是基于LPV(linear parameter-varying)模型插值建模的基礎，該辨識算法在顫振試飛數據處理中得到了較為成熟的應用[7]。Caigny等[8]提出了局部估計的狀態空間模型插值算法(state-space model interpolation of local estimates，SMILE)，適用于基于狀態空間的變參數系統插值建模。唐煒等[9]從試飛數據出發，基于LPV模型實現了氣動彈性魯棒顫振邊界預測，并在風洞試驗中進行了驗證。Baldelli等[10]基于飛行試驗建立了基于LPV模型的顫振預測模型，并通過ATW(aerostructure test wing)進行了飛行試驗數據的驗證，驗證了該建模方法的有效性。Desforges等[11]研究了顫振試飛中掃頻信號的優化問題，給出了分段線性掃頻相較于傳統簡單線性掃頻的優勢。

為了提升實時監控的及時性及試飛數據的信噪比，開發基于局部氣動彈性系統模型的響應預測和激勵優化方法是一種有效的途徑。在顫振試飛速度擴展過程中，通過已飛狀態點數據建模，實現關鍵部位響應的預測，在試飛激勵動作實施過程中，將預測響應作為時域數據安全監控的基準，可同時直觀考慮響應收斂情況及量值大小。若實際試飛響應偏離預測值太多，則考慮阻尼大的突變，及時實施安全措施，不必等待準實時數據處理結果，提高安全監控的及時性。通過建模仿真，預測模態頻帶范圍，在結構最大響應限制情況下，通過優化激勵信號，提升試飛響應數據質量，同時盡量縮減每個激勵動作的時間，提高試飛效率。

1 基于試飛數據的局部氣彈系統建模

1.1 顫振試飛響應預測及激勵優化的流程

建立全局LPV模型來描述某一飛機氣動彈性系統，面臨依賴參數較多和模型規模太大帶來的不確定性等問題。顫振試飛一般采用等高度變速法進行，對于顫振試飛包線擴展，針對在局部試飛狀態，選取少量局部關鍵結構響應，簡化為隨速壓參數改變的局部線性變參數模型，通過不斷迭代更新數據，提高準確度及工程可操作性。

顫振試飛響應預測與激勵優化隨試飛狀態推進的流程，如圖1所示。從已飛試飛狀態點數據出發，通過頻域子空間結合最小二乘復頻域的兩步辨識算法建立各狀態點有限階線性時不變(linear time-invariant，LTI)。然后以速壓為變參數進行插值，通過局部估計的狀態空間模型插值算法[12]，建立隨速壓變化的LPV模型。最后通過所建模型的仿真進行激勵信號優化和響應預測，在試飛實施過程中以預測響應為基準進行安全監控。飛行結束后將新飛狀態點數據迭代到模型進行模型更新后接著進行下一點的預測，直到速度擴展結束。

圖1 顫振試飛響應預測與激勵優化流程圖Fig.1 Flowchart of response prediction and excitation optimization for flight flutter test

顫振試飛掃頻激勵信號一般面臨激勵時間限制和最大結構響應限制這兩個限制，所以激勵優化包含兩方面：一種是對掃頻頻率變化率優化；一種是依據響應最大限制值進行激勵幅值優化。頻段優化一方面為了避免模態遺漏，應覆蓋一定的頻帶范圍，同時可根據預測下一狀態點的模態頻率，采用分段線性掃頻進行激勵信號優化，在選取的頻點范圍內增加掃頻時間，提高關心頻帶激勵能量和響應的信噪比。顫振試飛中關鍵部位振動限制值的設置是為了保障安全，避免太大的激勵值帶來飛機靜強度與操縱的問題，所以對激勵幅值需要進行一定的限制。

1.2 已飛狀態點LTI狀態空間模型辨識

(1)

可等價表示為

H～l=A～lB～lC～lD～l

(2)

假設連續傳遞函數為Gc(s)，雙線性變換給出的離散時間傳遞函數在頻域上為

(3)

(4)

以頻響函數G(ejωkT)為系統的輸出，單位陣為系統的輸入，p個輸出、m個輸入的一個離散時間LTI系統的狀態方程表示為

(5)

經過各頻點ωk(k=1,2,…,N)處的重復迭代及整合后得到如下形式

G=OrX+ΓrW

(6)

其中，

Γr為下三角Toeplitz塊矩陣

式(6)的實數矩陣形式為

Gre=OrXre+ΓrWre

(7)

式中：Gre=[Re(G) Im(G)]；Wre=[Re(W) Im(W)]； Re為取實部； Im為取虛部。

對矩陣[(Wre)T(Gre)T]進行QR分解

(8)

(9)

(10)

其中,

J1=[I(r-1)p0(r-1)p×p],
J2=[0(r-1)p×pI(r-1)p],
J3=[Ip0p×(r-1)p]

(11)

(12)

上述子空間方法中，通過奇異值大小確定系統階次，在氣彈系統辨識中，由于有噪聲干擾，往往不容易給出清晰的階次分界。辨識系統的阻尼取值對顫振試飛有重要意義，所以在響應估計中需要兼顧極點的精確性，文獻[15]提出結合頻域最小二乘法的兩步辨識方法，通過LSCF法辨識出極點[16-18]，結合上述子空間的零點辨識，提高系統模態參數的精度，有利于進行插值建模時兼顧響應幅值和模態參數預測精度。在待辨識氣彈系統結合地面試驗和理論計算得出辨識階次較為明確的情況下，也可以采用其他辨識方法給出系統的極點。

(13)

(14)

(15)

系統在兩步法下的連續時間狀態空間估計為

H～l=A～lB～lC～lD～l

(16)

1.3 氣彈系統局部LPV模型估計

H(q,q·,…)=A(q,q·,…)B(q,q·,…)C(q,q·,…)D(q,q·,…)

(17)

一般顫振試飛是以一系列離散速壓為狀態點進行，為簡化分析，認為保持一定高度和速度范圍下的激勵響應測試是一個等速壓狀態，可以只保留式(17)LPV模型的靜態響應部分來描述局部變速壓氣彈系統，即式(17)可簡化為

H(q)=A(q)B(q)C(q)D(q)

(18)

(19)

為第l個LTI的觀測陣，通過各相干性轉換矩陣建立的相干性局部LTI等價模型為

Hl=TlA～lT-1lTlB～lC～lT-1lD～l

(20)

H^(q)=∑Ni=1fi(q)H^i=∑Ni=1 fi(q)A^ifi(q)B^ifi(q)C^ifi(q)D^i

(21)

考慮如下代價函數

(22)

(23)

式中，F(h)為關于h的線性函數，插值建模過程轉變為一個線性最小二乘優化問題。

(24)

2 仿真示例

采用Texas A & M University帶副翼二元翼段俯仰-沉浮二自由度風洞試驗件的線性數值模型進行仿真驗證，模型簡圖如圖2所示，具體參數見文獻[19-21]。該模型輸入為副翼偏度，符合一般顫振試飛舵面激勵情況，顫振試飛中響應測量一般為加速度信號，在該翼段前緣處設置加速度計，將原模型的俯仰角及沉浮位移輸出轉變為前緣法向加速度信號的響應輸出。

圖2 二元翼段模型Fig.2 2DOF wing section model

該二元翼段氣動彈性方程為

(25)

式中：M∈2×2為廣義質量矩陣；C,K∈2×2為廣義阻尼矩陣和廣義剛度矩陣，為速壓的q的函數；F為輸入矩陣；u=β為舵面偏度輸入；x=[hα]T為俯仰位移與俯仰角組成的向量；E=[1db]為輸出矩陣；y為加速度傳感器輸出。

將氣動彈性方程寫成狀態空間形式如式(26)，理論數值計算顫振速壓為89.9 Pa。

x·x··y=0I0-M-1K-M-1CF-EM-1K-EM-1CEFxx·u=H(q)xx·u

(26)

由狀態空間方程進行仿真，獲取在速壓20 Pa，25 Pa，30 Pa，35 Pa下該系統的舵面激勵及加速度響應。激勵為0～5 Hz，幅值0.1 rad，時間60 s的線性掃頻信號，響應信號添加信噪比10 dB的白噪聲。

估計以上4個狀態點的傳遞函數，通過頻響函數，經兩步法辨識獲得系統在4個速壓點的LTI模型估計。以在速壓20 Pa下的情況為例，仿真得到的激勵輸入與響應輸出經濾波后，如圖3所示，使用最小二乘法獲得系統穩態圖，如圖4所示，頻域子空間和兩步法獲得的局部LTI估計的頻響與原始頻響對比，如圖5所示。

圖3 激勵輸入與響應輸出(20 Pa)Fig.3 Excitation and response (20 Pa)

圖4 LSCF法獲得的系統極點穩態圖(20 Pa)Fig.4 Poles stable chart obtained by LSCF method (20 Pa)

圖5 辨識LTI模型與原模型頻響對比(20 Pa)Fig.5 Comparison of frequency response between identified LTI model and original model (20 Pa)

(27)

圖6 預測響應與原數值模型響應對比(60 s，40 Pa)Fig.6 Comparison of predicted response and original numerical model response (60 s，40 Pa)

圖7 預測響應與原數值模型響應對比(10 s，40 Pa)Fig.7 Comparison of predicted response and original numerical model response (10 s，40 Pa)

頻帶上優化為0～0.8 Hz、掃頻5 s，0.8～3.0 Hz、掃頻45 s，3.0～5.0 Hz、掃頻10 s，的三段線性掃頻信號；由于是線性系統，依據上述的仿真響應最大值與激勵信號幅值比例，可確定在理想響應值下，激勵信號優化后的幅值為0.137 rad。優化后激勵信號,如圖8所示，使用該信號激勵下的系統仿真輸出，如圖9所示，可見響應最大值為1.5 m/s2，對比圖6響應，兩個模態頻率附近激勵時間增多，更有利于事后數據分析。

圖8 優化后激勵信號(40 Pa)Fig.8 Optimized excitation signal (40 Pa)

圖9 使用優化后激勵獲得的系統響應(40 Pa)Fig.9 Response using optimized excitation (40 Pa)

根據圖1所示流程圖，40 Pa后進行速壓點擴展，每增加5 Pa為一個速壓點，在目標響應值保持1.5 m/s2的情況下基于LPV預測模型優化后激勵選取幅值與理論模型計算所需激勵幅值對比，如圖10所示。

圖10 數值模型計算與LPV模型預測所需激勵幅值Fig.10 Required excitation amplitudes from calculation of numerical model and prediction of LPV model

速度點擴展過程中，可依據所建LPV模型進行顫振邊界預測，輔助判斷下一速度點穩定性及余量，隨著速度點增加，即用來插值的LTI模型增加，預測精度逐步提高，如圖11所示。

圖11 顫振速度預測值隨插值點數的變化Fig.11 Predicted flutter velocity varying with modeling interpolation points

3 試飛應用

3.1 顫振試飛安全監控

傳統的顫振試飛實時安全監控主要關注飛機的飛行狀態參數和結構的時域振動響應，在飛機達到要求的飛行狀態后，施加幾十秒的激勵信號，試飛工程師通過遙測數據監控結構響應信號，觀察掃頻下時域響應的峰值變化情況，結合自身監控經驗及該型機前期飛行監控情況做出定性判斷，待激勵結束后，方可對記錄在本地計算機中的遙測數據進行準實時分析，獲得關鍵模態的頻率及阻尼。這種實時監控對監控人員的經驗要求較高，且在整個試飛周期內會盡可能安排同一工程人員對同一組關鍵響應進行持續監控，提高對響應變化的敏感性，同時準實時分析對遙測信號的質量要求較高，斷點和跳點等均對結果影響較大。使用已飛數據建立的變速壓模型進行響應預測，則能夠為激勵期間的實時安全監控提供一個定量的參考值，通過觀察實測各峰值與預測峰值之間的偏差判斷目前系統穩定性是否在預算范圍內，及時發現阻尼較為迅速的變化趨勢，縮短應急中斷試驗的反應時間。

采用某型機顫振試飛的9個試飛狀態點數據，驗證本文響應預測模型的工程適用性，狀態點飛行高度一致，采用操縱面掃頻進行激勵。在通過前3個狀態點試飛后，開始構建基于變速壓模型預測下一狀態點激勵響應。依據圖1所示進行實時安全監控設計，每飛一個狀態點，進行模型更新，把下一點的響應預測作為飛行安全監控參考值。

以最后一個試飛狀態點為例，速度點VD試飛前，進行機翼尖和垂尾尖部響應預測。首先辨識獲得前8個狀態點的局部LTI模型，如機翼尖響應在副翼反對稱激勵下，前8個狀態點辨識所得LTI模型的幅頻和實測幅頻,如圖12所示，其中,第5個狀態點的辨識LTI模型頻響與原頻響對比,如圖13所示，可見兩步法比單子空間法精確度高，有利于后續響應預測準確度。接著通過插值算法獲得系統變速壓LPV模型，預測第9狀態點處的機翼尖部及垂尾尖部響應。最后一個狀態點，副翼反對稱激勵下機翼尖部預測響應與實際試飛響應對比,如圖14所示，副翼對稱激勵下機翼尖部預測響應與實際試飛響應對比,如圖15所示，方向舵激勵下垂尾尖部預測響應與實際試飛響應對比,如圖16所示。對比量值及變化趨勢，關鍵的幾個模態在時域響應上的變化均得到較好預測，設計合理的安全系數后，預測響應可滿足顫振試飛安全監控要求。

圖12 前8個狀態點辨識與實測幅頻Fig.12 Amplitude frequency responses of identified LTI model and measured system for preceding 8 state points

圖13 辨識LTI模型頻響與原頻響對比圖(13.2 kPa)Fig.13 Frequency response of identified LTI model compared with that measured (13.2 kPa)

圖14 機翼尖在反對稱激勵下的預測與實際響應值Fig.14 Predicted and actual acceleration responses at wing tip under anti-symmetric aileron excitation

圖15 機翼尖在對稱激勵下的預測與實際響應值Fig.15 Predicted and actual acceleration responses at wing tip under symmetric aileron excitation

圖16 垂尾尖預測與實際響應值Fig.16 Predicted and actual acceleration responses at vertical tail tip

以垂尾響應為例，垂尾尖部響應最大值與預測響應最大值在各狀態點下對比,如圖17所示。可見預測值準確度隨著插值狀態點的增加而增加。

圖17 真實及預測垂尾響應最大值隨速壓擴展的變化Fig.17 The real and predicted maximum response values at vertical tail tip varying with the extension of dynamic pressure

基于響應預測可以直觀預測下一試飛狀態點是否穩定，同時也可以基于已知狀態點所建LPV模型進行外推，得到下一狀態點的穩定性預估及是否具有足夠余量，確保下一狀態點安全。該型飛機危險模態為垂尾彎扭耦合，基于變速壓模型進行該機垂尾部結構的顫振邊界預測。如采用8個已飛狀態點數據構建垂尾激勵下垂尾響應的局部LPV模型，模型預測的速度-阻尼-頻率曲線與各狀態點實測值關系，如圖18所示，預測速壓達到21.3 kPa時，發生垂尾彎扭顫振，扭轉阻尼趨于0。第9個試飛狀態點速壓為18.53 kPa,具有適當余量，可進行該試飛點的試飛。實際試飛中結合邊界預測結果通過控制速度增加步長來確保試飛安全。若下一狀態點速壓比預測邊界速壓大，為保障試飛安全，則縮短試飛速度步長，使下一試飛狀態點速壓減小到預測邊界速壓以下，待試飛后，依據試飛數據更新模型后進行進一步預測分析。試飛中重點關注已知數據建模后預測的下一試飛狀態點穩定性及余量，這種基于外插法的邊界預測雖然數據點少時不太準確，但隨著試飛速度擴展及接近邊界，精度逐步提高。

圖18 垂尾局部LPV模型頻率-阻尼-速壓曲線Fig.18 Frequency-damping-dynamic pressure diagram of local LPV model for vertical tail

3.2 激勵優化

某型飛機進行顫振試飛，由于飛機油量、性能等限制，無法在高速度狀態點下保持足夠多激勵時間，隨著顫振試飛逐步進行，依據已知狀態點基于LPV模型建模后，逐步進行激勵優化，縮減激勵時間，提高響應信噪比。

以副翼激勵為例，原始副翼激勵設計為幅值0.8°，2～40 Hz，38 s的線性掃頻激勵，隨著飛行速度增加，狀態保持時間減少，需要縮減激勵時間。依據已飛的幾個試飛狀態點，預測該響應中包含的3個主要模態下一狀態點位于11～24 Hz。在響應不超過4g,保持掃頻范圍的情況下，依據所建LPV模型仿真后，優化激勵信號為幅值0.8°，2～11 Hz、掃頻2 s，11～24 Hz，掃頻13 s，24～40 Hz、掃頻4 s，共計19 s的激勵信號。為了驗證效果，在下一狀態點分別實施原激勵信號與優化后激勵信號，原激勵信號掃頻頻率隨時間變化及機翼翼尖響應，如圖19所示，優化后激勵信號的分段掃頻頻率隨時間變化及機翼尖響應，如圖20所示。可知優化后激勵信號在關注的頻帶范圍內仍然保持了1 Hz/s的頻率變化率，縮減激勵時間則有利于完成更多的激勵動作。估算原始激勵下機翼尖響應的信噪比為6.13 dB，優化后激勵下機翼尖信噪比估算為13.57 dB。圖21為本次副翼激勵優化前后自譜對比，可見優化后激勵提高了關心頻帶范圍內的激勵能量。

圖19 優化前線性掃頻頻率變化及機翼尖響應(38 s)Fig.19 Linear sweep frequency and wing tip response before optimization(38 s)

圖20 優化后分段線性掃頻頻率變化及機翼尖響應(19 s)Fig.20 Multi-section linear sweep frequency and wing tip response after optimization(19 s)

圖21 第一次激勵信號優化前后激勵自譜對比Fig.21 The comparison of PSD of excitation signals before and after the first optimization

繼續隨著包線擴展接近右邊界，速壓進一步增大，噪聲也進一步增加，且隨著數據積累及模型迭代，幾個主要關注模態的頻率范圍進一步穩定，進一步需要優化激勵信號為幅值0.8°，5.0～13.5 Hz、掃頻4.25 s，13.5～23.5 Hz、掃頻15 s，23.5～25.0 Hz、掃頻0.75 s，共計20 s的激勵信號。整體頻帶范圍進一步縮小，且關心頻帶范圍頻率變化率提高到1.5 Hz/s。試飛實施中作為對比組的幅值0.8°，頻率范圍5.0～25.0 Hz，20 s的簡單線性掃頻激勵下掃頻頻率變化及機翼尖響應，如圖22所示，優化后掃頻頻率及機翼尖響應，如圖23所示，估算的信噪比由8.29 dB提高到10.83 dB。該次副翼激勵優化前后自譜對比，如圖24所示，關心頻帶激勵能量進一步提高。

圖22 優化前線性掃頻頻率變化及機翼尖響應(20 s)Fig.22 Linear sweep frequency and wing tip response before optimization (20 s)

圖23 優化后分段線性掃頻頻率變化及機翼尖響應(20 s)Fig.23 Multi-section linear sweep frequency and wing tip response after optimization (20 s)

圖24 第二次激勵信號優化前后自譜對比Fig.24 The comparison of PSD of excitation signals before and after the second optimization

4 結 論

基于LPV模型使用顫振試飛數據進行飛機氣彈系統的局部建模，形成了模型預測支撐試飛，試飛數據迭代修正模型的顫振試飛流程，在部分型號飛機的顫振飛行試驗中得到較好的應用。隨著插值建模的試飛狀態點數據增加，預測模型多次迭代后的精度逐步提高，在臨近右邊界可以提供較為準確的預測值，可有效解決顫振試飛時域監控參考值問題，滿足工程安全監控及時性的需求。在激勵時間、頻帶需求和關鍵部位響應限制情況下，基于預測模型進行仿真計算所優化的激勵信號，可有效提升響應數據信噪比，避免激勵過程中響應超限，有利于提高試飛效率。

基于狀態空間構建動力學模型的方法由于有較為明確的激勵和響應數據要求，適用于能夠采用操縱面掃頻激勵得到主要模態的顫振試飛。文中建模方法在亞音速、高亞音速飛機的顫振試飛中得到較好驗證，對于馬赫數效應明顯的超音速飛機的顫振試飛，在實際建模時不能只考慮速壓一個變參數，同時需要考慮馬赫數和其他可能的影響因素。后續進一步研究適用于超音速顫振試飛的多變參數LPV預測模型，同時研究該方法在氣動伺服彈性試飛中的應用。