面向一體化設計的高超聲速飛行器變參建模方法

湯佳駿，曹瑞，劉燕斌，陸宇平，朱鴻緒，衣春輪

(1. 南京航空航天大學 自動化學院，江蘇 南京210016；2. 南京航空航天大學 航天學院，江蘇 南京210016)

高超聲速飛行器復雜的飛行環境以及自身特殊的外形結構，導致其呈現出強非線性、強耦合、模型不確定以及快速時變等復雜的系統特征，因而給飛行器外形和控制系統設計帶來了諸多困難與挑戰[1-2]。相比于傳統先外形后控制的設計思路，本體-控制一體化在設計初期就考慮了控制系統的性能優化，反向指導了飛行器的外形設計，節約了參數反復迭代的成本，設計效率高，成為了飛行器外形、控制等研究領域的熱點[3-4]。文獻[5]提出了融合控制的多學科優化設計思想，并以X-29的初始設計為例說明了控制對飛行器設計的限制。文獻[6]考慮了僅有操縱面大小一個本體參數的一體化設計問題，對飛行器進行了本體和控制器的一體化優化。文獻[7]涉及了高超聲速飛行器一體化設計的實例，研究了多個外形參數變化對本體設計和控制設計的影響。眾多的研究表明本體-控制一體化思想在飛行器的初期設計階段，對結構外形的優化設計具有重要的意義[8]。

高超聲速飛行器由于具有更強的非線性特性及參數不確定性，從實際工程應用的角度來看，在飛行器一體化設計優化框架中融入參數化建模是必要的。同時，進一步考慮氣動、隱身、控制等多學科需求的參數化建模必然需要引入更多的參數，參數的增加會使優化設計任務增大并降低優化效率[9-11]。因此，隨著本征系統模型復雜度的增加，引入代理模型技術是本體-控制一體化優化設計的必經之路。文獻[12]針對高超聲速飛行器后體/尾噴管進行一體化設計，通過試驗數據計算、代理模型建立和智能算法優化3個過程，完成了對飛行器的優化建模。文獻[13]關注系統建模、不確定性和本體系統參數間的耦合問題，分析了其對系統建模和控制系統設計的影響。文獻[3]分別采用靈敏度分析的數據擬合方法和基于機理推導方法建立高超聲速飛行器的代理模型并進行了變參性能分析。目前，許多的系統優化設計問題中使用的代理模型，都存在對本征系統動力學的大量簡化或是完全忽略動態特性的靜態分析，例如在固定巡航點的配平線性化模型進行優化設計，這類簡單的線性模型無法反映飛行條件的變化對設計參數帶來時變影響，也不能適應優化設計中本征系統的變化[14]。

綜上，本文提出了一種高超聲速飛行器面向一體化的變參代理建模方法。該方法從一體化設計的需求出發，綜合考慮高超聲速飛行器大包線飛行的特點，建立包含飛行條件、外形結構等多參數的變參代理模型。所建立的模型在高超聲速飛行器設計前期可以為面向本體-控制一體化的優化設計服務，在設計中后期使得增益調度控制、魯棒控制等線性系統設計方法能夠應用于非線性系統設計[15-16]，為后續的控制性能分析、控制系統設計提供基礎。

1 面向一體化的變參代理建模方法

本文提出面向一體化的高超聲速飛行器變參代理建模方法，方法包括3個步驟，如圖1所示。

1.1 多參數的非線性動力學模型建立

對于高超聲速飛行器這一復雜的動力學系統，其系統狀態與本體參數、飛行條件之間存在復雜的耦合關系。這些參數與飛行器運動狀態之間的關系是通過作用在飛行器上的力和力矩建立的。故而要在高超聲速飛行器的力和力矩計算中包含所研究的參數，假設高超聲速飛行器有n個可變的調度參數?1,?2,…,?n，n≥2，則高超聲速飛行器的動力學模型可表示為：

fmodel(G,L,D,M,T)=f(?1,?2,…,?n)

(1)

式中：調度參數?1,?2,…,?n既可以是飛行器本體設計參數p1,p2,…，如飛行器質量、質心位置、操縱面大小、前體壓縮角等，也可以是研究對象的飛行條件，如高度h和馬赫數Ma。

圖1 面向一體化的變參代理建模方法流程Fig.1 Flow chart of parameter varying modeling method for integrated design

1.2 高超聲速飛行器線性變參代理模型建立

對高超聲速飛行器結合間隙度量理論建立線性變參代理模型的過程中，本文提出調度參數的分類策略，以解決模型參數數量多而帶來的參數高維和計算復雜的問題。

利用間隙度量理論建立線性變參代理模型在理論上對調度參數數量并沒有限制[18]，但隨著調度參數維數的增加，區域維度和區域內狀態點數量也在增加，進而導致使用該方法建立飛行器變參模型時的計算量大大增加。為了解決這一問題，通過對系統間間隙度量值隨單個調度參數的變化趨勢進行分析，探索飛行器對象特性的特異性，并通過機理推導驗證個別調度參數對模型系統影響呈現出的規律性，以此分析形成分類策略，將調度參數分成了規律性參數和無規律性參數2類。即在之后的建模過程中可以省去規律性參數的計算過程，達到參數降維，減小計算量的目的，進而擴大建模方法的適用范圍。

圖2展示了建立線性變參代理模型的流程，這一過程主要又可以分為建立全包線線性時不變(LTI)系統、調度參數分類決策和建立線性變參代理模型3個部分，將在下文具體闡述。

圖2 建立多參數的線性變參代理模型流程Fig.2 Flow chart of linear multi-parameter varying modeling

1.2.1 建立全包線的LTI系統

1) 確定各個調度參數的變化區間[?imin,?imax]，i=1,2,…,n。

2) 對各參數在對應參數變化區間內進行均勻采樣，得到一個由所有不同參數的狀態點組成的n維參數變化區域，區域內的狀態點數量為:

(2)

3) 對各采樣狀態點進行線性化得到一系列線性時不變(LTI)系統集合，記作{Se|e=1,2,…,Nt}。

1.2.2 調度參數分類決策

取單個調度參數?i，分析系統間間隙度量值隨該參數的變化趨勢是否呈現出規律性，對變化趨勢呈現明顯規律的調度參數，推導線性系統量綱導數關于該調度參數的偏導數，判斷調度參數對模型的影響能否分離，若能，即為規律性參數，反之或間隙度量值隨參數的變化趨勢無明顯規律，該調度參數劃為無規律性參數。規律性參數的數量記為x。

下面以高度參數為例進一步說明調度參數分類決策及確定規律性參數的過程。平衡狀態點處飛行器小擾動線性化模型為:

式中：X=[Vαhqθ]T；U=[δeφ]T;ΔX=X-Xeq;ΔU=U-Ueq;下標“eq”表示平衡點處的值。

(4)

高度的變化主要影響矩陣A中的Xh、Zh/Veq和Mh3項，計算式分別為:

(5)

若高超聲速飛行器模型的氣動導數CL、CD、Cm、CT均與高度參數無關，則高度參數變化對模型氣動力和力矩的改變僅是影響了動壓Q，則式變為:

(6)

式中：Sr為飛行器機翼的參考面積;cr為飛行器機翼的平均弦長;由于高超聲速飛行器的轉動慣量大導致Kheq數值很小；而高超聲飛行器模型的質量及飛行速度也均較大，故系數KXeq與KZeq數值也較小且隨配平狀態的變化不大。因此，對于不同高度平衡點處的線性化系統，它們矩陣A之間的差異及系統間的間隙度量值δgap受上述3個量綱的影響就呈現與?Q/?h相關的變化規律:

(7)

故而，對氣動導數與高度參數關聯較小的高超聲速飛行器模型，可以將高度參數劃歸為規律性參數。參考上述分析，對于調度參數與氣動力、推力之間關聯簡單的特殊模型，亦或是無動力再入等特殊飛行過程，都可以參考上述思路。針對具體模型分析，提出調度參數的分類策略，尋找規律性參數簡化建模過程。

1.2.3 建立線性變參模型

1) 再次遍歷參數集，若調度參數?i為無規律性參數：

①遍歷其余調度參數?j，j≠i計算相鄰狀態點l:(?i,?j)和l+1:(?i+Δ?i,?j)間線性系統的間隙度量值，記作δgap,?i(Sl,Sl+1)，其中l=1,2,…,N?i-1。

②獲得隨調度參數?i變化的狀態點系統間間隙度量值序列，計算每組序列中間隙度量值的累加和:

∑l=∑l-1+δgap,?i,?j(Sl,Sl+1)

(8)

③確定衡量指標τ，每次當序列累加和滿足式要求時，取對應的?i,j,k作為參考區間邊界。其中，取值次數k=1,2,…,Nd,?i。

(9)

⑤對無規律性參數構建新的n-x維參數變化區域，并根據各參數區間內的邊界劃分整個n-x維區域，結果可得參數變化子區域的數量為:

(10)

⑥為保證標稱點線性系統與子區域內其他狀態點處線性系統的差異最小，遍歷參數變化子區域φm，m=1,2,…,N，子區域內包含Nm個狀態點。計算子區域內每個狀態點與其他狀態點系統間的間隙度量值并求和，再對總和取均值作為衡量指標:

(11)

選取指標最小的狀態點的調度參數值即為標稱點對應該參數的值?i,P。

若調度參數?i為規律性參數：

①按1.2.2節中分析出的間隙度量值隨該參數的變化規律，將參數變化范圍[?imin,?imax]劃分為Nd,?i個子區間。

b)取各子區間的中點作為參數變化子區域中標稱點對應該調度參數的值?i,P。

②確定整個n維參數變化區域的區域劃分及區域內的標稱點，標稱點記作Pm=(?1,P,?2,P,…,?n,P)。繪制所確定的標稱點線性系統的伯德圖，比較對應參數變化子區域內各狀態點的動態特性是否相似，判斷子區間和標稱點選取的合理性；若不合理，則減小的采樣間隔或指標τ，重復之前的流程。

3) 運用最小二乘法對標稱點處的線性系統的狀態矩陣(Ai,Bi)進行擬合得到關于調度參數的多項式矩陣為:

(12)

最終建立高超聲速飛行器的線性變參代理模型:

式中δ?i∈[-1,1]為各調度參數歸一化表示為：

(13)

1.3 面向一體化設計的代理模型校驗

本體設計參數是一體化設計關注的重點，代理模型的準確度，尤其是關于本體設計參數變化規律的保真度，對面向一體化設計進行變參代理建模工作具有重要的意義。因此，提出以下步驟對建立的變參代理模型進行一系列校驗。

1)驗證建立的代理模型是否準確反映了原模型的開環非線性特性。繪制高超聲速飛行器非線性模型及變參代理模型在相同參數下的階躍響應曲線，比較各模型開環響應是否一致。

2)驗證建立的模型是否保留了原模型隨本體參數變化的系統特性。繪制高超聲速飛行器非線性模型與變參代理模型隨本體參數變化下的配平狀態量，比較各模型隨本體參數的變化趨勢是否一致。為保證校驗準確，需要改變本體參數的采樣間隔，以確認在建模過程中未采樣的本體參數狀態下，模型之間是否仍呈現相似的特性。

3)驗證建立的代理模型是否可以替代原模型進行控制設計。設計經典的閉環控制結構用于高超聲速飛行器非線性模型與變參代理模型，繪制閉環階躍響應曲線，比較各模型在相同控制器作用下的閉環響應是否一致。

2 變參代理建模方法應用

高超聲速飛行器的質心位置與飛行器的縱向靜穩定性有著密切的關系，將質心位置作為設計變量進行本體-控制一體化設計十分必要[19]。在此按照本文所提出的面向一體化的變參代理建模方法對參考文獻[20]中的高超聲速飛行器模型建立以馬赫數、高度和質心位置為調度參數的變參代理模型，并仿真校驗建立的代理模型，說明了此代理建模方法的準確性和適用性。

2.1 高超聲速飛行器代理模型建立

所采用的高超聲速飛行器非線性動力學模型中包含馬赫數Ma,高度h和質心位置Xcg3個調度參數。

首先，高度和馬赫數是2個飛行條件調度參數，確定了所研究的高超聲速飛行器一體化優化設計問題覆蓋的包線范圍，取馬赫數范圍為[5,10]，高度范圍為[25 km，30 km]，文獻給出的模型參考質心位置在-16.8 m附近，確定模型質心位置的變化范圍為[-15.8 m,-17.8 m]，取平均采樣間隔ΔMa=0.5，Δh=0.5 km，ΔXcg=-0.2 m。計算相鄰系統間的間隙度量值，分析其隨馬赫數、高度及質心位置的變化趨勢，如圖3～5所示

圖3 相鄰系統間間隙度量值隨馬赫數變化趨勢Fig.3 Trend of gap metric changing with Mach number between adjacent systems

通過對趨勢圖的分析可知，系統間間隙度量值主要受馬赫數變化的影響，隨高度和質心位置的變化較小。結合式，馬赫數參數和質心位置參數與氣動力、力矩關聯復雜，變化規律難以推導，劃歸為無規律性參數；而參考1.2.2節中的分析，系數KXeq與KZeq均為負值，動壓的變化率隨高度增加而增大，驗證了系統間間隙度量值隨高度增加而減小，故將高度參數劃歸為規律性參數。

因此，根據規律性參數參考間隙度量值的變化規律將飛行區域的高度劃分為[25 km，26.5 km]和[26.5 km，30 km]2個區間，標稱點高度參數的確定以高度區間的中點為原則選取。馬赫數和質心位置2個無規律性參數選取指標τ=0.2作為劃分準則進行參數變化區域的劃分。最終將整個3維的參數變化區域劃分為16個子區域，各個分區的區域范圍及標稱點如圖6所示。其中，實線代表了整個參數變化區域的邊界，圓點對應的是平均化采樣的狀態點，虛線劃分了子區域，星號標記點為子區域內的標稱點。

圖4 相鄰系統間間隙度量值隨高度變化趨勢Fig.4 Trend of gap metric changing with altitude between adjacent systems

圖5 相鄰系統間間隙度量值隨質心位置變化趨勢Fig.5 Trend of gap metric changing with the position of centroid between adjacent systems

接著對每個子區域驗證標稱點選取的合理性，以第14個子區域為例，標稱點各參數分別為Ma=9.5、h=26 000 m、Xcg=-17.6 m。圖7中的圓圈表示子區域內任一狀態點與其他所有狀態點線性系統間間隙度量值的加和平均值，質心位置不同的狀態點在圖中同一高度和馬赫數位置疊加顯示，可以發現圓點對應的均值最小，該狀態點即為對應子區域的標稱點。

圖8為第14個子區域標稱點與參數變化區域內其他狀態點線性系統單輸入單輸出通道的頻域特性對比。標稱點線性系統頻域特性為紅色曲線，其他狀態點線性系統頻域特性為藍色曲線，可以看出紅色曲線位于藍色曲線包絡的中間位置，說明其他狀態點的線性系統與標稱點的線性系統有相似的動態特性，故而驗證選取的標稱點可以替代整個參數變化子區域。

圖6 參數變化區域劃分及標稱點選取Fig.6 Division of parameter varying area and selection of nominal points

將馬赫數、高度及質心位置參數歸一化表示為δMa、δh和δXcg。圖6中的16個標稱點組成了整個參數變化區域的平衡點集合，其分別對應16個線性系統(Ai,Bi)，其中i=1,2,…,16。采用最小二乘法擬合得到狀態矩陣關于δMa、δh和δXcg的解析表達式，最終得到高超聲速飛行器縱向運動的線性變參代理模型

(5)

式中：A(δMa,δh,δXcg)、B(δMa,δh,δXcg)均是由歸一化的調度參數δMa、δh和δXcg確定的矩陣。

2.2 高超聲速飛行器代理模型檢驗

2.2.1 開環特性檢驗

比較非線性模型與線性變參代理模型的開環響應，選取部分狀態點，以Ma=7、h=26 000 m的飛行狀態為例。調度參數中的質心位置作為本體設計參數在Xcg=-16 m、Xcg=-16.5 m、Xcg=-17 m、Xcg=-17.5 m之間變化時，非線性模型與線性變參代理模型在指定參數下的單通道階躍響應如圖9所示。

由圖9可得所建立的線性變參代理模型與非線性模型雖然存在由線性化方法產生的特性差異，但開環響應的變化趨勢相同，這保證所建立的變參代理模型具有與原對象一致的開環非線性特性。

圖9 2種模型在指定參數下單通道階躍響應對比Fig.9 Comparison of single-channel step response of two models under specified parameters

2.2.2 隨本體參數變化的系統特性檢驗

針對研究的質心位置這一本體設計參數，比較非線性模型和變參代理模型的配平狀態隨質心位置的變化趨勢。不同于之前的建模過程，為保證變參代理模型能夠捕捉到質心位置在區間未采樣部分變化的變化規律，在質心位置變化區間[-15.8 m,-17.8 m]取更小的間隔ΔXcg=-0.1 m。配平舵面偏轉角隨質心位置的變化如圖10所示。

圖10 不同模型配平舵面偏轉角隨質心位置變化Fig.10 Trim deflection angle of elevator changing with the position of the centroid among different models

可以看出非線性模型與變參代理模型的配平舵面偏轉角隨質心位置的變化趨勢一致，反映變參代理模型能夠有效表征本體參數在原模型中的特性。

2.2.3 閉環特性檢驗

參考文獻[18]對建立的變參代理模型取短周期部分，采用經典PI控制結構，以極點配置法設計了縱向增穩控制器。將控制器分別用于變參代理模型和原非線性模型，效果對比如下:

由圖11控制器階躍響應曲線的比較，非線性模型和變參代理模型在相同的控制器作用下的響應具有一致性。

綜合上述校驗，表明所建立的變參代理模型具有良好的保真度，可以代替原模型很好地用于本體-控制一體化方面的設計研究。此外，本文的變參建模方法優先考慮模型保真度的要求對確定性系統進行建模，沒有對系統不確定性因素作詳細分析。然而，不確定性系統可以由確定系統加上不確定性因素后擴維形成的確定性系統表示[21]。針對具有強不確定性的高超聲速飛行器，考慮不確定性作為模型參數使系統的維度增加，但本文的變參代理建模方法依舊適用。

圖11 縱向增穩控制器不同模型控制效果對比Fig.11 Comparison of control effects of longitudinal stabilization controller among different models

3 結論

1)本文提出的高超聲速飛行器變參代理建模方法面向本體-控制一體化優化設計，兼顧了本征系統的特性變化和大包線飛行的特點，拓展了線性變參代理模型反映本征參數對系統動態特性影響的能力，使模型既能適應本體參數優化設計迭代又能用于控制系統的分析與設計，具有工程實用價值。

2)針對面向一體化的多參數建模帶來的參數高維問題，結合間隙度量理論的分析和機理推導，對調度參數進行決策分類，改進了多參數代理建模過程中計算量大的問題，擴展多參數建模的適用范圍。

3)建模方法包含完整的模型校驗流程，通過實例驗證，證實該方法能夠保證模型具有較好的保真度，尤其是能良好地反映出本體特征參數的特性。