含間隙非線性機翼跨聲速顫振時滯反饋控制

聶雪媛，鄭冠男，楊國偉

(中國科學院力學研究所 流固耦合系統力學重點實驗室，北京 100190）

氣動力與柔性結構相互作用會產生氣動彈性，顫振是氣動彈性領域中最危險的一類動不穩定現象，極易引發災難性事故。顫振主動控制技術是目前研究最多的顫振抑制方法，一般是通過在機翼上布置多個控制面，控制其偏轉改變作用在機翼上的氣動力，以達到抑制顫振的目的。

主動控制氣動彈性系統包含作動器、傳感器、控制器和數字濾波器等元器件，其動態特性會導致最終作用于結構的控制力產生時滯［1-2］。時滯會破壞控制器的控制性能，甚至導致被控系統失穩。在氣動彈性領域，已有學者對此現象展開了研究。國內，Zhao［3］對不可壓流場反饋通道存在延時的翼型氣彈穩定性進行了研究，指出在飛行器氣動伺服彈性設計時，不可忽略時滯影響。Huang等［4］分析了輸入延時對不可壓流場的飛行器穩定性的影響，并提出了一種最優控制方法抑制顫振。Cai等［5-6］針對前向通道存在時延，采用控制和滑模變結構控制方法對二元翼型顫振進行控制，指出不考慮時滯所設計的控制器無法有效抑制延時系統的顫振。Xu等［7］研究了超聲速下時滯對氣動彈性系統顫振邊界穩定性的影響。國外，Ramesh和Narayanan［8］對超聲速下反饋通道存在時滯的翼型，設計了狀態反饋控制方法。Marzocca等［9］研究了考慮時滯存在時，對線性和非線性控制器作用下的二維翼型氣動彈性系統穩定性的影響。Araujo和Santos［10］采用Smith預估方法研究了對控制通道存在不同時滯的顫振抑制效果。

就目前已有文獻來看，針對考慮時滯的顫振主動控制研究幾乎都集中在亞、超聲速域，氣動力計算基于線化氣動力模型，結構視為線性結構。然而，現代民用／軍用飛機大都以跨聲速巡航，激波位置對結構響應影響敏感，較為精確的模擬氣動力非線性的方法是計算流體力學(CFD），但該方法計算量太大，且難以用于控制器設計。此外，機翼控制面的鉸鏈處普遍具有間隙，會引起剛度非線性現象［2］，間隙非線性對氣動彈性有顯著影響，在分析帶控制面的機翼顫振問題時，必須加以考慮。

針對以上問題，本文基于氣動力降階模型(Reduced Order Model，ROM）技術，以含間隙非線性的二維翼型為對象，對輸入信號存在延時的控制系統展開了跨聲速機翼顫振時滯反饋主動控制方法的研究。通過狀態變換，將時滯被控系統轉換為不顯含時滯的狀態方程，進行最優時滯反饋控制。在數值仿真中，采用自適應時間步長以準確捕捉結構間隙切換點。

1 基于CFD／ROM 的氣動彈性模型

在跨聲速域，CFD方法能提供精確的非定常氣動力。本文通過CFD／CSD耦合方法，以白噪聲信號作為結構位移輸入，通過系統辨識技術建立非定常氣動力的自回歸滑動平均(Auto Regressive Moving Average，ARMA）降階模型。

1.1 控制方程

氣動力計算采用基于RANS的三維Navier-Stokes控制方程，守恒型的流動方程可表示為

式中：Q為守恒向量；Gc和Gv分別為對流通量和黏性通量；S為控制體V的邊界面積；n為面的法向量；t為物理時間。

方程(1）對無黏項離散采用Roe格式，黏性項離散采用二階中心差分格式，時間推進采用雙時間步。

1.2 結構運動方程

本文以含間隙非線性的帶控制面二維翼型(見圖1）作為研究對象，俯仰方向結構剛度含有間隙非線性，非線性廣義“位移-力”關系如圖2所示。

圖1 帶控制面的二元翼型Fig.1 Two-dimensional airfoil with control surface

圖2 間隙非線性Fig.2 Free-play nonlinearity

式中：M代表回復力矩，ˉM(α）為分段線性函數，α為俯仰角；kα為俯仰方向的剛度系數；δ為間隙角閾值。

二元翼型的非線性控制方程可寫為

式中：h為沉浮位移；b為機翼半弦長；kh為沉浮方向剛度系數；m為機翼質量；Iα為機翼慣性矩；xα為彈性軸到機翼重心的無量綱距離；L和M 分別為氣動力和力矩。

其中：μ＝m／(ρπb2）；Ma∞＝U∞／a∞；Cl和Cm分別為升力和俯仰力矩系數；ρ為來流密度。

定義xs＝［ξ α］T，則式(4）可轉換為狀態空間形式：

1.3 基于ARM A的開環氣動彈性模型

對于模型辨識而言，最重要的是激勵信號的選取。一般來說，要求激勵信號能夠激起被辨識系統所關心的所有頻段的信號。本文選取隨機信號作為結構位移，與CFD耦合，計算得到在該激勵下的氣動力輸出。對輸入白噪聲，輸出非定常氣動力的信號采用系統辨識方法，選用ARMA差分模型對氣動力建模。該模型具有時間域離散形式，如下：

式中：fa代表氣動力向量；η為結構位移輸入向量；A、B為要辨識的系統模型參數；na和nb分別為輸出和輸入的延遲階數。

方程(6）通過變量代換，可以化為狀態方程離散形式(見(7）），具體推導過程可詳見文獻［11］。

由于結構方程為連續系統，需要對氣動力離散結構轉換為連續系統，通過雙邊變換，可得氣動力的連續系統狀態空間為

將方程(5）與方程(8）耦合，可得開環氣動彈性模型為

值得注意的是，從式(9）可以看出，該方程中系數矩陣Asa包含分段線性函數(α），因此該方程為非線性方程。由于間隙非線性屬于不光滑非線性，本文所研究的結構具有2個切換點±δ(見圖2），若采用時間步長固定的龍格庫塔方法求解，將無法準確捕捉到間隙非線性切換點的變化。為此，本文采用二分法［12］搜索切換點確定自適應時間步長的龍格庫塔法對式(9）進行求解。

2 時滯反饋主動控制設計

2.1 時滯控制系統

對式(9）所描述的氣動彈性模型，可通過控制面偏轉，其目的是引入附加氣動力來抑制顫振。引入控制作用后的閉環氣動彈性結構框圖如圖3所示(開關撥到1點，不考慮控制信號的延時）。

圖3 氣動彈性系統閉環控制框圖Fig.3 Block diagram of closed-loop control of aeroelastic system

圖3中：β為控制面偏轉角，代表控制器輸出量，fca為控制面偏轉角偏轉引起的氣動力。λ為時滯。該氣動力模型采用1.3節方法建立了控制面偏轉氣動力降階模型：

式中：xc為控制面偏轉狀態變量；Aca、Bca、Cca、Dca為該模型的系數矩陣。

將式(9）和式(10）合并，即可得廣義被控對象模型為

當控制器存在時滯λ時(即在如圖3所示的控制系統中開關撥到點2），廣義被控對象模型可描述為

2.2 考慮時滯的控制方法設計

定義變量z(t）如下：

將該變量代入式(12a），由文獻［13］可知，式(12a）可變為不顯含時滯的系統：

可以證明，只要系統(12a）是完全可控的，則系統(14）也是完全可控的。

當時滯λ已知時，式(14）所描述的系統為分段線性系統(考慮結構間隙非線性，系數Asac為分段線性矩陣），本文對分段子系統采用最優反饋控制，通過極小化式(15）所示的目標函數，獲得控制器最優控制量β。

式中：Γ和R分別為狀態變量和控制變量的加權矩陣。

控制器的輸出最優控制量為

式中：P為Ricatti方程的解。

對式(14）所示的分段線性系統采用式(16）所示的最優控制，即可實現顫振主動控制。

由式(13）可知，實際系統狀態變量xsac的求解中包含了積分項s，對于任意時滯λ，總可以寫成λ＝lT－n，T為連續系統采樣周期，l為大于0的整數，0≤n＜T。在實際控制中，控制器輸出的控制量通過保持器將數字電路中的離散信號轉換為連續信號，以零階保持器為例：

基于上述思路，通過積分變量替換，式(13）中的積分項可以用不同歷史采樣時刻的控制量來表示［14］。

由式(18）可以看出，系統當前時刻的狀態與控制量的歷史數據相關，而歷史數據的長度則取決于時滯λ的大小。

3 算例與分析

以圖1所示的含間隙非線性二元機翼為例［15-16］，考慮控制信號存在時滯，對跨聲速機翼顫振進行主動控制。該模型無量綱參數設置選自文獻［15］，如表1所示。計算工況為馬赫數Ma∞＝0.8，迎角為0°。

表1 模型無量綱參數Table 1 Non-dimensional parameters of model

以濾波白噪聲作為結構位移，辨識得到非定常氣動力離散模型(7），通過雙邊變換，最終得到連續狀態方程(8）。

為驗證辨識模型的準確性，先以常用的3211信號作為位移η的給定輸入，進行CFD非定常氣動力計算，并與模型(8）的計算輸出結果進行比較，如圖4所示。

圖4 3211信號激勵下ROM和CFD計算的氣動力結果比較Fig.4 Comparison of aerodynamic forces with ROM and CFD excited by 3211 signals

進一步，按照文獻［15］設置初始擾動，即選取結構位移初值xs0＝［0 0.002 8 0 0.001］T，在Ma∞＝0.8U*＝1.48工況下，計算氣動力降階模型與結構耦合的開環氣彈模型(9）的響應，并與CFD／CSD直接耦合得到的結果進行比較，如圖5所示。可以看出，基于氣動力降階模型得到的開環模型(9）具有較高計算精度，能夠替代CFD／CSD的直接耦合計算。

圖5 ROM和CFD結算結構響應結果比較(Ma∞＝0.8U* ＝1.48）Fig.5 Comparison of structure response resultwith ROM and CFD(Ma∞＝0.8U* ＝1.48）

圖6給出了該工況下俯仰角的相軌跡。可以看出，俯仰角從初始位置(0.002 8，0.01）出發，在氣動力的作用下逐漸向原點位置收斂，結構是穩定的。繼續增大U*，當U*＝1.72時(文獻［15］在該工況初始條件下U*＝1.7），系統出現極限環振蕩，如圖7所示。

圖6 俯仰角相軌跡(Ma∞＝0.8U* ＝1.48）Fig.6 Phrase portrait of pitch angle at Ma∞＝0.8U* ＝1.48

針對圖7所發生的極限環振動，采用控制面偏轉的方式加以控制。圖8為給定控制面板偏轉3211信號、CFD計算結果與白噪聲信號辨識得到的氣動力模型(10）的結果比較。

圖7 俯仰方向極限環振蕩(Ma∞＝0.8U* ＝1.48）Fig.7 Limit-cycle oscillation of pitch at Ma∞＝0.8U* ＝1.48

圖8 3211信號控制面偏轉輸入控制面ROM和CFD計算的氣動力結果比較Fig.8 Comparison of control surface aerodynamic forces with ROM and CFD excited by 3211 deflection input signals

若控制通道不存在時滯，即圖3所示開關接通觸點1，此時采用常規的最優控制方法即可控制住顫振，如圖9所示。控制器輸出控制量的變化如圖10所示?？紤]到系統(11）狀態的可觀測性，本算例對式(15）中的狀態加權矩陣Q選為對角陣，對角線元素為100，控制加權矩陣R設為50。

圖9 無時滯系統顫振控制效果Fig.9 Flutter control result for system without delay

圖10 閉環控制系統控制面偏轉角Fig.10 Control surface deflection of closed-loop control system

若圖3所示系統控制信號存在時滯，即開關接通觸點2，繼續采用不考慮時滯的控制方法，在時滯λ較小時，控制方法是有效的。但當增大到一定值時(本算例中為無量綱時間0.37，T＝0.05），其控制后的系統響應相軌跡如圖11所示?？梢钥闯?，即使控制面偏轉到最大可允許偏轉角，但結構響應的相軌跡在整個相平面發散，系統出現失穩，控制方法失效。

圖11 時滯系統的無時滯控制相平面圖及控制面偏轉角Fig.11 Phase portrait of no delay control for system with time delay and corresponding control surface deflection

在這種情況下，采用本文提出的考慮延時設計的控制方法(16），對該系統進行顫振控制，此時l＝8，m ＝0.03，控制后的結構位移相軌跡如圖12所示。

從圖12可以看出，在控制作用下，翼型位移從初始狀態，出發，在相平面內其運動軌跡形成一條封閉曲線，收斂于各自的平衡態(0，0）。

圖12 考慮時滯的控制器對時滯系統顫振控制的相平面圖及控制面偏轉角(λ＝0.37）Fig.12 Phase portrait of flutter control for system with time delay exerted by controller designed with time delay and corresponding control surface deflection atλ＝0.37

繼續增大時滯，圖13為無量綱時滯λ＝1.6時，顫振控制的翼型位移相軌跡及控制面偏轉角?？梢钥闯?，考慮時滯設計的控制方法在大時滯下仍然有效。

圖13 考慮時滯的控制器對時滯系統顫振控制的相平面圖及控制面偏轉角(λ＝1.6）Fig.13 Phase portrait of flutter control for system with time delay exerted by controller designed with time delay and corresponding control surface deflection atλ＝1.6

4 結論

1）跨聲速下，以含間隙非線性的翼型為研究對象，對輸入信號存在時滯的閉環氣動彈性系統顫振主動控制方法進行了研究?；诜嵌ǔ鈩恿惦A模型結合間隙非線性具有分段線性的特點，建立了非線性被控氣動彈性模型，通過引入含積分項的狀態變換，將時滯被控系統轉換為無時滯系統，并在此基礎上進行最優反饋控制設計，所設計的控制量考慮了時滯的影響。

2）使用針對無時滯系統所設計的控制方法，對存在時滯的系統進行顫振主動控制，其控制效果隨時滯增大而減弱，直到增大到某一臨界值(本算例為無量綱時間0.37）時，控制方法失效；本文提出的時滯反饋控制方法能有效地處理控制通道的時滯，實現顫振主動抑制，可行性不受時滯大小的影響。

下一步工作將在此基礎上，考慮不確定時滯對氣動彈性系統穩定性的影響，時滯不再是事先確定的常數，應研究如何設計控制方法提高系統的穩定性。