基于伴隨理論的無人機氣動彈性優化

黃胤錚，王曉喆，王柳青

（1.北京航空航天大學 飛行學院；2.北京航空航天大學 無人系統研究院，北京 100083）

0 引言

無人機因具有體積小、成本低、對惡劣環境適應性好等優點，在民用和軍用領域均有巨大的應用潛力［1］。但無人機的作戰能力會受到飛機結構、氣動布局、任務載荷、巡航時間等因素影響，無人機任務完成度以及生存概率也與其飛行性能和作戰性能密切相關。隨著航空技術的發展，未來空域的環境復雜度都會有不同程度的上升，對無人機的整體性能提出了更高要求。

隨著有效載荷不斷增加，無人機的結構重量系數不斷下降，結構柔度不斷增加，氣動、結構等學科的耦合關系進一步增強，若在優化過程中采用串列優化方法會造成學科間的反復迭代，在降低設計效率的同時，也影響了整體性能的提升。氣動彈性優化充分考慮了氣動與結構之間的耦合效應，并以此挖掘設計潛力［2］，實現飛行器在復雜環境下經濟性和使用性能的提升［3］。在大型客機［4-5］和超音速客機［6］的設計中，氣動彈性優化已經得到了較為廣泛的應用。基于敏度信息的優化算法在求解大規模的多學科優化問題時，收斂較快［7-8］，而以遺傳算法為代表的進化式算法需要大量樣本，不適合直接應用于多設計變量的氣動/結構耦合優化問題中［9］。1988 年，Jameson［10］首次將伴隨方程法應用于氣動目標函數梯度的求解過程中，并對機翼進行優化設計。由于伴隨理論的計算消耗與設計變量的數量無關，適合處理具有大規模設計變量的多學科優化問題，因此近年來得到了廣泛關注［11-14］。Anderson 等［15］采用多學科建模優化方法和伴隨理論，實現了風力渦輪機葉片的高保真結構優化和負載應力最小化；Wang 等［16］采用基于伴隨方法的優化設計框架對機翼受動載荷進行了分析優化；Batay 等［17］采用多學科設計優化工具，結合伴隨求解器，對風力渦輪機進行了并行氣動設計優化；Wilke［18］采用自主開發的多目標優化框架實現了對直升機旋翼槳葉的優化設計；李潤澤等［19］在多目標優化過程中結合伴隨方法，提高了超臨界機翼氣動優化設計效率；劉曉東等［20］采用基于伴隨理論的氣動優化設計方法，實現了飛翼布局飛行器的氣動優化。

隨著現代無人機對高機動以及長航時性能的需求不斷提高，采用基于伴隨理論的氣動彈性優化方法深度挖掘多學科耦合所帶來的設計潛力是必要且迫切的。本文結合渦格法和有限元分析法，通過松耦合的方式實現氣動彈性分析，并基于序列二次規劃和伴隨理論開展非線性敏度優化。分別針對長航時和高機動兩類典型無人機，以最小燃油消耗為目標，對機翼的氣動外形和空間梁結構等參數進行優化設計，并對優化設計的結果進行對比分析。

1 氣動彈性仿真

1.1 氣動渦格法

對于飛行器的氣動分析本質上是對空間N-S 方程的求解。渦格法在位勢流方程的求解過程中較為常見，與空間N-S 方程的求解相比，利用基于求解位勢流的氣動計算方法能夠以更高的計算效率對啟動導數進行比較準確的求解［21］。

根據邊界不可穿透定理，可以在控制點處添加馬蹄渦與氣動表面相切的條件，如式（1）所示。其中，CC為控制點處誘導速度的系數矩陣，Γ為每個渦線渦強組成的列向量，V∞為控制點來流速度矩陣，N為控制點處單位法向量。

通過求解式（1）可得到每個渦的渦流強度和氣動力，如式（2）所示。其中，Fn為渦段上的升力，ρ∞為來流密度，V∞為來流速度，Vn為當前渦段中點誘導速度，Γn為當前渦段上的渦強。

渦格法渦線布置方法如圖1所示。

Fig.1 Vortex lattice method vortex line layout method圖1 渦格法渦線布置方法

1.2 結構有限元法

對于組成機翼的各部分而言，機翼上下表面的蒙皮主要承載機翼上的彎曲載荷，機翼的翼梁和副梁主要承載剪切載荷，蒙皮和梁共同承載扭轉載荷。在分析過程中，將機翼翼梁簡化為空間梁單元的有限元模型進行分析。空間梁單元的有限元模型如圖2所示［22］。

Fig.2 Finite element model of space beam element圖2 空間梁單元有限元模型

通過對空間梁模型的分析，可以獲得結構分析部分的控制方程如下：

其中，f為作用在節點上的力和力矩，d為節點位移列陣。通過求解式（3），可以得到各節點上的靜力學載荷分布。

1.3 仿真模型

氣動彈性的耦合求解采用松耦合方法，氣動與結構之間的數據傳遞如圖3、式（4）和式（5）所示［23］。

Fig.3 Aerodynamic load transfer圖3 氣動載荷傳遞

其中，i=1，2，分別代表氣動網格對應有限元模型中的左右兩個結構節點，T為氣動網格上通過渦格法計算得到的合力，rsai為結構節點指向氣動力作用點的方向向量，Fi、Mi分別為通過該牽引方案得到的節點力和節點扭矩。

在完成氣動載荷傳遞后，通過有限元分析可以得到對應的結構位移，并傳遞回氣動網格。

其中，i=1，2，分別代表氣動網格對應有限元模型中的左右兩個結構節點，di為平動位移，θi為轉動位移，rsai同式（5），為結構節點指向氣動力作用點的方向向量。在此基礎上，結合開源框架OpenMDAO 建立無人機的氣動、結構以及氣動彈性耦合分析模型［24］。

2 優化算法及框架

2.1 伴隨理論

傳統的導數計算方法以有限差分法和復變量差分法為主，在梯度優化中，都需要對控制方程反復求解，計算量大。而伴隨理論算法可以在保持精度的前提下，一次性求解出目標函數對所有設計變量的導數［3］。其計算量僅與目標函數和約束條件的數量有關，與設計變量的數量無關，適合多設計變量的優化設計問題。

本文的無人機氣動彈性優化以氣動結構耦合系統作為分析系統，以氣動參數和結構參數作為設計變量。優化目標函數I和控制方程殘差R可設定為：

其中，x表示設計變量；w、u表示氣動/結構耦合系統狀態變量，其中w具體表現為渦格法中網格單元的氣動負載，u具體表現為有限元分析中結構有限元分析得到的結構位移；Ra對應氣動分析的殘差，Rs對應結構分析的殘差。

令優化目標函數I和控制方程殘差R分別對設計變量x進行求導，并利用恒等變換和伴隨算子ψ得到伴隨方程。通過求解伴隨方程，可以一次性求解出目標函數針對所有設計變量的導數，如式（8）所示：

2.2 序列二次規劃

常用的梯度優化方法包括BFGS 擬牛頓算法、共軛梯度法、序列二次規劃算法（Sequential Quadratic Programming，SQP）等。其中，序列二次規劃法對大規模數據的計算效率較高，在工程中得到了廣泛應用。具體形式如式（9）所示：

其中，f（x）為目標函數，h（x）為等式約束函數，g（x）為不等式約束函數。

序列二次規劃法是一種用于求解非線性最優化問題的方法，其基本思路是將原問題轉化為一系列二次規劃（Quadratic Programming，QP）問題進行迭代求解。在每次迭代中，通過求解二次規劃問題得到下一次迭代的值，直到序列｛xk｝收斂于極值點。具體來說，對于某一次迭代值xk，NLP 會被近似為該點處的二次規劃問題。通過求解該二次規劃問題，可以得到下一次迭代的值xk+1。隨著迭代次數k的增加，序列｛xk｝將收斂于極值點［25］。本文采用SLSQP 算法（Sequential Least SQuares Programming），該算法是由Kraft［26］于1988 年提出的一種非線性約束優化算法，通過將優化問題轉化為一系列線性或二次規劃子問題，用于求解無約束或約束的非線性優化問題。具體來說，SLSQP 算法采用牛頓法求解每個子問題，并利用輔助函數法處理約束條件。

2.3 氣動彈性優化框架

本文所采用的氣動彈性優化框架如圖4 所示。具體流程如下：①設定優化目標并設置目標參數，明確目標函數、約束函數、設計變量與狀態變量等相關變量及數據；②根據初始數據建立氣動和結構的物理模型；③對目標函數和約束函數關于所有設計變量的總梯度進行計算；④利用序列二次規劃算法進行尋優；⑤通過不斷更新氣動和結構的物理模型，直至得到收斂解，完成優化設計。

Fig.4 Basic process of optimization design based on gradient optimization algorithm圖4 基于梯度優化算法的優化設計基本流程

3 無人機氣動彈性優化

3.1 大過載無人機

在有人機上，8g是受過專業訓練的航天員一般可以承受的最大過載［27］。而無人機沒有相關限制，可以做出更多大過載的機動動作來規避潛在危險。

3.1.1 工況計算

大過載無人機具體的飛行工況如表1 所示。選取NACA0015 翼型，機翼初始參數如表2 所示，并選取7075 鋁合金作為翼梁結構材料。

Table 1 High overload UAV flight state parameters表1 大過載無人機飛行狀態參數

Table 2 Initial shape parameters of high overload UAV wing表2 大過載無人機機翼初始外形參數

目標函數為巡航狀態下的燃油消耗量，對于飛機燃油消耗，采用Breguet 距離方程進行計算。具體方程形式如下：

其中，Wf為燃油重量，W0為飛機空重，Ws為結構重量，R為最大航程，V為巡航速度，SFC為燃油消耗率。設計變量為機翼梢根比、機翼后掠角、機翼扭轉角、機翼翼梁厚度、機動狀態下的機翼迎角。設計變量具體數量及數值范圍如表3 所示。約束條件為：機動狀態下機翼結構應力接近但不超過應力極限、機動狀態下升力等于重力。

Table 3 Parameter range and number of design variables for high overload UAV表3 大過載無人機設計變量參數范圍及數量

3.1.2 結果優化

經過優化設計后，可得到如表4 所示的3 種優化設計結果。無人機機翼的外形經過優化設計后的變化如圖5所示。

Table 4 Optimal design results of high overload UAV表4 大過載無人機優化設計結果

Fig.5 Optimal design results of UAV with different wing root chord lengths圖5 不同機翼翼根弦長無人機優化設計結果

通過對不同機翼翼根弦長情況下機翼優化設計結果的綜合分析，可以看出，在機翼經過優化設計后，大過載無人機的燃油消耗量相較于初始狀態均降低了20%～30%。以機翼翼根弦長5.5 m 為例，如圖6 所示，在經過優化設計后，可以看到馮米塞斯應力在巡航狀態下遠低于應力極限，符合基本邏輯和設計要求。在扭轉角方面，機翼翼尖和翼根的相對扭轉角為2.89°，符合機翼外形設計邏輯，改善了升力分布，有利于提升飛機的巡航性能［28］。設計方案總體上較為合理，符合飛行器機翼外形的設計邏輯。

Fig.6 Optimal design results of wing with root chord length of 5.5m圖6 5.5m機翼翼根弦長優化設計結果

3.2 長航時無人機

3.2.1 工況計算

對于長航時無人機，具體的飛行工況如表5 所示。選取NACA0015 翼型，采用如圖7 所示的矩形機翼，初始參數如表6所示，并選取7075鋁合金作為翼梁結構材料。

Table 5 Flight state parameters of long-endurance UAV表5 長航時無人機飛行狀態參數

Table 6 Initial parameters of long-endurance UAV wing表6 長航時無人機機翼初始參數

目標函數為燃油消耗量，設計變量為機翼梢根比、機翼后掠角、機翼扭轉角、機翼翼梁厚度、巡航狀態下的機翼迎角。設計變量的具體數量及數值范圍如表7 所示。約束條件為：機翼結構應力不超過應力極限、巡航狀態下升力等于重力。

Table 7 Parameter range and number of high overload UAV flight status表7 高過載無人機飛行狀態參數范圍及數量

3.2.2 結果優化

長航時無人機優化后的機翼幾何外形如圖8所示。

Fig.8 Optimized wing geometry of long-endurance UAV圖8 長航時無人機優化后的機翼幾何外形

初始狀態下的燃油消耗量為3 898 kg，長航時優化后的燃油消耗量為3 004 kg，優化后的結構重量為1 761 kg。與初始狀態相比，燃油消耗量降低了22.93%。優化設計后的具體參數如表8 所示。長航時情況下的優化設計結果如圖9 所示。在扭轉角方面，機翼翼尖和翼根的相對扭轉角為5.15°，處于合理范圍內。在馮米塞斯應力方面，由于需要獲得最佳的優化設計結果，因此需要在設計過程中盡可能逼近結構極限。本文的優化設計方案未超出應力極限，符合飛行器外形設計邏輯。

Table 8 Optimized design results of long endurance UAV表8 長航時無人機優化設計結果

Fig.9 Optimized design results of long endurance UAV圖9 長航時情況下的優化設計結果

4 結語

本文提出一種現代無人機的高效氣動彈性優化方法，并通過大過載和長航時兩類典型無人機進行了驗證，具體如下：

（1）提出基于伴隨理論的氣動彈性優化方法，通過氣動渦格法與有限元分析法進行氣動彈性松耦合分析，并基于伴隨理論和序列二次規劃法進行非線性敏度優化。

（2）針對大過載無人機進行結構優化設計，以燃油消耗量為目標函數，以機翼外形參數和機翼翼梁結構參數為設計變量，結合機動狀態下氣動、結構兩方面的約束函數，使無人機的燃油消耗量分別減少了21.66%、27.93%、32.05%。

（3）針對長航時無人機進行結構優化設計，以燃油消耗量為目標函數，以機翼外形參數和機翼翼梁結構參數為設計變量，結合巡航狀態下氣動、結構兩方面的約束函數，使無人機的燃油消耗量減少了22.93%。

綜上所述，本文提出的方法高效、穩健、可靠，時間成本和學習成本低，可為無人機氣動彈性優化提供理論和方法參考。