一種仿真數(shù)據(jù)驅(qū)動(dòng)的小型無人機(jī)旋翼優(yōu)化設(shè)計(jì)框架

郝佳瑞,呂利葉,孫 偉,宋學(xué)官

(大連理工大學(xué)機(jī)械工程學(xué)院,遼寧 大連116024)

1 引言

無人機(jī)應(yīng)用日益廣泛，故提升其的續(xù)航時(shí)間和負(fù)載能力具有重要工程意義。目前，計(jì)算流體力學(xué)（Computational Fluid Dynamics，CFD）被廣泛應(yīng)用于無人機(jī)旋翼優(yōu)化問題的仿真計(jì)算中。文獻(xiàn)［1］基于CFD對(duì)不同外形的無人機(jī)旋翼在懸停狀態(tài)下的推力進(jìn)行了對(duì)比，結(jié)論是弧型的旋翼投影形狀可以得到更大的推力。文獻(xiàn)［2］結(jié)合CFD仿真和遺傳算法（Genetic Algorithm，GA）確定了無人機(jī)葉片的弦長和安裝角的最佳分布。文獻(xiàn)［3］提出了CFD仿真驅(qū)動(dòng)的序列二次規(guī)劃（Sequential Quadratic Programming）算法來優(yōu)化亞音速機(jī)翼的結(jié)構(gòu)參數(shù)，取得了較高效的翼型。但是由于葉片參數(shù)化模型的限制，上述工作沒有將翼型參數(shù)和結(jié)構(gòu)參數(shù)的結(jié)合，因此設(shè)計(jì)自由度不夠大。

為了解決上述問題，具有較強(qiáng)翼型描述性能的Class/shape Function Transformation（CST）［4］被應(yīng)用于描述沿翼展方向的二維翼型變化，并結(jié)合結(jié)構(gòu)參數(shù)實(shí)現(xiàn)對(duì)三維無人機(jī)旋翼葉片的外形控制。同時(shí)，構(gòu)建基于CFD仿真數(shù)據(jù)的代理模型作為無人機(jī)旋翼的真實(shí)物理模型的替代能夠大大減少計(jì)算量，被廣泛應(yīng)用［5］。文獻(xiàn)［6］提出的一種高保真擴(kuò)展自適應(yīng)混合代理模型（Extended Adaptive Hybrid Functions，E-AHF），對(duì)參數(shù)量大的工程黑盒問題具有更高的擬合精度，適用于無人機(jī)旋翼的復(fù)雜幾何外形優(yōu)化。

綜上，建立了如圖1的一種新的快速優(yōu)化框架：首先，建立了基于翼型E387原無人機(jī)旋翼的高精度CFD模型，在所有工作轉(zhuǎn)速下誤差在4.8%以內(nèi)；然后，基于CST函數(shù)對(duì)葉片三維外形進(jìn)行了參數(shù)化，參數(shù)共22個(gè)，使用拉丁超立方采樣［7］進(jìn)行實(shí)驗(yàn)設(shè)計(jì)并建立了上述［6］中的擴(kuò)展自適應(yīng)混合代理模型；最終，在懸停狀態(tài)下，應(yīng)用GA遺傳算法進(jìn)行了以旋翼推力為約束，以最小化旋翼扭矩為目標(biāo)的優(yōu)化，并對(duì)優(yōu)化結(jié)果進(jìn)行了流固耦合分析，保證了其可行性。

圖1 基于E-AHF混和代理模型的葉片優(yōu)化設(shè)計(jì)框架Fig.1 E-AHF Model-Based Optimization Framework of Rotor Blades

2 理論基礎(chǔ)

2.1 計(jì)算流體力學(xué)（CFD）

旋翼的空氣動(dòng)力學(xué)模型假設(shè)為忽略重力三維可壓縮粘性流動(dòng)，其控制方程組（Reynolds Average Navier Stokes Equations）寫作：

式中：i,j=1,2,3—x,y,z三個(gè)方向；ν—速度；ρf—流體密度；p—流體微元體上的壓力；E—流體的能量；τij—因分子粘性作用而產(chǎn)生的作用在微元體表面上的粘性應(yīng)力的分量；δij—克羅內(nèi)爾符號(hào)（Kronecker delta）；qi—流體熱通量的分量。

2.2 CST翼型參數(shù)化

CST方法［8］是由Kulfan提出的一種描述翼型上下表面坐標(biāo)點(diǎn)的參數(shù)化方法，對(duì)翼型的擬合精度較高，且能產(chǎn)生連續(xù)光滑的幾何外形。該參數(shù)化方法分為類函數(shù)（Class function）和形狀函數(shù)（Shape function）兩個(gè)部分，其擬合公式為：

式中：c—弦長，x、z—翼型沿弦長方向的坐標(biāo)及其垂向坐標(biāo)，ψ、ζ—x、z坐標(biāo)的無因次表達(dá)；ζU、ζL—翼型上下表面的無因次坐標(biāo)；公式最后一項(xiàng)ψ·ΔζU、ψ·ΔζL定義了翼型上下表面的尾緣厚度—描述一行基本輪廓的類函數(shù)。

2.3 擴(kuò)展自適應(yīng)混和代理模型（E-AHF）

通過對(duì)單一代理模型進(jìn)行誤差分析，擴(kuò)展自適應(yīng)混和代理模型將不同的代理模型賦予不同的權(quán)重進(jìn)行融合，權(quán)重系數(shù)和各模型的內(nèi)置參數(shù)由數(shù)據(jù)集訓(xùn)練得到，最終獲得高于任何子代理模型的精度。首先，采用高斯隨機(jī)過程進(jìn)行誤差估計(jì)：

設(shè)均值函數(shù)為0，協(xié)方差由貝葉斯線性回歸求得：

其次進(jìn)行模型評(píng)估，引入核函數(shù)作為輸出響應(yīng)的函數(shù)：

式中：ybest(xj)—組分代理函數(shù)在輸入變量xj處的最好預(yù)測響應(yīng)，通過基準(zhǔn)函數(shù)獲得。

最終計(jì)算局部權(quán)重：

模型評(píng)估采用決定系數(shù)R2以及均方根誤差RMSE：

式中：yi—測試點(diǎn)xi處的真值，測試點(diǎn)xi處的代理模型所預(yù)測響應(yīng)為—所有測試點(diǎn)的真值的平均，Ntest—測試點(diǎn)個(gè)數(shù)。

3 無人機(jī)旋翼參數(shù)化與仿真

原旋翼采用E387翼型工作參數(shù)如表1，實(shí)物以及3D模型，如圖2所示。

圖2 原模型實(shí)物（左）及其三維模型示意（右）Fig.2 Baseline Rotor and 3-D Model

表1 無人機(jī)旋翼工作參數(shù)Tab.1 Parameters of UAV Operating

3.1 旋翼幾何參數(shù)化

由于旋翼葉片幾何外形的復(fù)雜性，其參數(shù)化過程分為兩個(gè)部分：首先是對(duì)旋翼葉片結(jié)構(gòu)（Platform）參數(shù)化，即確定弦長、安裝角的分布，前掠或后掠，尖梢比等參數(shù)；其次是對(duì)二維翼型（Airfoil）的參數(shù)化描述，方法諸如B樣條插值，PARSEC參數(shù)化，Hick-Henne參數(shù)化，CST參數(shù)化等等。如圖3結(jié)合旋翼結(jié)構(gòu)參數(shù)化，三階CST翼型參數(shù)化方法能夠用較少的變量去描述盡可能大的設(shè)計(jì)空間，實(shí)現(xiàn)對(duì)旋翼的快速優(yōu)化設(shè)計(jì)。

圖3 旋翼幾何參數(shù)化示意Fig.3 Platform Parameters

旋翼半徑R=45mm，葉尖、最大弦長處的位置固定，運(yùn)用6個(gè)參數(shù)描述旋翼葉片的結(jié)構(gòu)。φtip、ctip為葉尖處的安裝角和弦長；φmid、cmid為葉片最大弦長處的安裝角和弦長；dmt表征了最大弦長處的位置；ο表示二維翼型幾何扭轉(zhuǎn)中心，參數(shù)ε描述：當(dāng)ε=0時(shí)，ο位于翼型前緣，前掠為0。原模型結(jié)構(gòu)參數(shù)如下表2。

表2 原旋翼結(jié)構(gòu)參數(shù)Tab.2 Values of the Baseline Platform Parameters

對(duì)原模型二維翼型E387的描述采用三階CST參數(shù)化方法，通過Matlab R2018b優(yōu)化工具箱中的fmincon函數(shù)最小化幾何描述誤差，確定了上下翼型表面的CST參數(shù)。最終，控制懸疑葉片結(jié)構(gòu)的參數(shù)6個(gè)，葉尖翼型和最大弦長處翼型參數(shù)各8個(gè)，參數(shù)總計(jì)22個(gè)（詳見章節(jié)4.1與優(yōu)化結(jié)果一并列出）。

3.2 旋翼流場CFD仿真與實(shí)驗(yàn)驗(yàn)證

基于ANSYS ICEM進(jìn)行網(wǎng)格繪制，應(yīng)用ANSYS CFX進(jìn)行CFD仿真計(jì)算。流域網(wǎng)格為混合網(wǎng)格，外流域采用結(jié)構(gòu)網(wǎng)格尺寸如圖4a所示；內(nèi)流域采用非結(jié)構(gòu)網(wǎng)格。邊界條件設(shè)置如圖4b，外流域上下表面以及圓柱面為Open；側(cè)面為交界面（Interface），根據(jù)轉(zhuǎn)速設(shè)置周期旋轉(zhuǎn)邊界條件（Rotational Periodicity）；槳轂部分、葉片表面設(shè)置為Wall。湍流模型采用SST（Shear Stress Transport）模型［9］。

圖4 流域網(wǎng)格劃分以及邊界條件示意Fig.4 CFD Meshing and Boundary Conditions

圖6 仿真-實(shí)驗(yàn)驗(yàn)證結(jié)果Fig.6 Comparison of Thrust between CFD Results and Experimental Data

為了驗(yàn)證CFD仿真結(jié)果的準(zhǔn)確性，建立旋翼運(yùn)動(dòng)的高保真模型，搭建了實(shí)驗(yàn)平臺(tái)，如圖5所示，對(duì)原旋翼的工作轉(zhuǎn)速范圍內(nèi)的6個(gè)轉(zhuǎn)速下的推力進(jìn)行測量，并與仿真結(jié)果進(jìn)行對(duì)比如圖4b，其最大誤差不超過5%。

圖5 實(shí)驗(yàn)臺(tái)展示Fig.5 Experimental Device

4 基于代理模型的優(yōu)化設(shè)計(jì)

依照2.1章節(jié)中的參數(shù)化可以得到有關(guān)旋翼三維幾何形狀的22個(gè)設(shè)計(jì)參數(shù)，包括：6個(gè)結(jié)構(gòu)參數(shù)（記為Pairfoil）和16個(gè)葉尖和最大弦長處翼型的CST翼型控制參數(shù)（記為Pplatform）。隨后進(jìn)行拉丁超立方實(shí)驗(yàn)設(shè)計(jì)，并進(jìn)行章節(jié)2.2中的敘述的CFD仿真得到訓(xùn)練集和測試集數(shù)據(jù)建立懸疑扭矩TQ和推力TH的擴(kuò)展自適應(yīng)混合代理模型，得到精度如表3。并根據(jù)所建立的代理模型進(jìn)行如方程組（10）下定義的單目標(biāo)優(yōu)化以提升旋翼的懸停效率。

表3 E-AHF代理模型精度Tab.3 R2 and RMSE of E-AHF Model

除了旋翼推力TH的約束之外，還將旋翼葉片安裝角的負(fù)扭作為優(yōu)化問題的約束之一，這樣做有利于提高懸停效率，節(jié)省優(yōu)化算法的搜索時(shí)間；另外，還考慮了CST二維翼型參數(shù)化以及翼型最小厚度tmax的約束，以保證能夠產(chǎn)生真實(shí)有效的旋翼外形。

優(yōu)化算法采用Matlab R2018b內(nèi)置GA遺傳算法［10］。

5 優(yōu)化結(jié)果討論與流固耦合驗(yàn)證

5.1 優(yōu)化結(jié)果與討論

優(yōu)化前后參數(shù)對(duì)比，如表4所示。

表4 優(yōu)化前后的參數(shù)比較Tab.4 Parameter Values of Optimized Model and Baseline Model

經(jīng)優(yōu)化后的葉片正投影形狀以及壓力分布如圖7所示，其中CT、CQ和FM分別為旋翼的推力系數(shù)、功率（扭矩）系數(shù)和懸停效率。如圖8a對(duì)比了Y=0.3R、0.7R、0.9R處的截面翼型輪廓并展示了元模型和優(yōu)化后模型的3D效果圖，圖8b則比較了二者沿翼展方向的弦長和安裝角的分布。

圖7 旋翼葉片壓力云圖Fig.7 Pressure Contour of Baseline Blade and Optimization Blade

圖8 原模型與優(yōu)化模型的對(duì)比Fig.8 Comparison between the Optimized and Baseline Model

可見，優(yōu)化后的旋翼葉片低壓面壓力分布更加均勻且低壓分布面積較廣，葉片的有效面積更小，由葉根到葉尖的安裝角負(fù)扭更大；且相對(duì)于原旋翼翼型E387，優(yōu)化后的旋翼翼型具有較大的彎度。優(yōu)化后的旋翼葉片的扭矩減少12.5%，懸停效率提高了8.32%，總重相對(duì)減少24.23%。

為了檢測優(yōu)化結(jié)果是否在不同的飛行狀態(tài)也表現(xiàn)出良好性能，對(duì)優(yōu)化后旋翼進(jìn)行了不同轉(zhuǎn)速下的CFD仿真，得到了不同狀態(tài)下的優(yōu)化旋翼的工作性能參數(shù)，并與原模型比較結(jié)果，如圖8所示。

圖9 非優(yōu)化轉(zhuǎn)速下原模型與優(yōu)化模型性能對(duì)比Fig.9 Performance of Optimized Rotor in Off-design Conditions

顯然，優(yōu)化后的旋翼不僅在懸停狀態(tài)下表現(xiàn)出良好性能，在其他非優(yōu)化設(shè)計(jì)點(diǎn)也能減小無人機(jī)旋翼扭矩從而提高效率。

5.2 流固耦合驗(yàn)證

為了確保優(yōu)化旋翼在極限工作狀態(tài)（17000RPM）下的結(jié)構(gòu)安全性，基于ANSYS WORKBENCH執(zhí)行了其單向流固耦合分析，結(jié)合高保真CFD分析和FEM（Finite Element Method）有限元靜力學(xué)分析確定其最大應(yīng)力以及最大變形，判斷其是否失效。旋翼葉片主要受到由于自身質(zhì)量產(chǎn)生的離心力以及在旋轉(zhuǎn)過程中受到的空氣作用力（此處忽略葉片重力），葉片的運(yùn)動(dòng)與受力分析圖，如圖10所示。單向流固耦合的網(wǎng)格信息以及結(jié)果如表5，可見，優(yōu)化后的旋翼仍然能保證極限工作狀態(tài)下的結(jié)構(gòu)安全性。

圖10 旋翼葉片受力分析Fig.10 Loads on UAV Rotor Blade

表5 E-AHF代理模型精度Tab.5 Data for FSI of Baseline Model and Optimized Model

6 結(jié)論

針對(duì)小型無人機(jī)旋翼外形提出了一種基于擴(kuò)展自適應(yīng)混合代理模型（E-AHF）的快速高保真優(yōu)化設(shè)計(jì)框架，并且驗(yàn)證了其可行性?；贑ST方法提出了一種簡潔的旋翼參數(shù)化方法，設(shè)計(jì)變量共計(jì)22個(gè)，進(jìn)行了經(jīng)過實(shí)驗(yàn)驗(yàn)證的高精度CFD仿真計(jì)算，得到其懸停狀態(tài)下的旋翼拉力和扭矩，分別建立E-AHF代理模型，并執(zhí)行以旋翼拉力、葉片負(fù)扭、翼型厚度等為約束，以旋翼扭矩為目標(biāo)的GA遺傳算法，得到的優(yōu)化后的旋翼葉片的扭矩減少12.5%，懸停效率提高了8.32%，總重相對(duì)減少24.23%。最終，我們驗(yàn)證了在不同轉(zhuǎn)速下優(yōu)化后旋翼的良好性能，并進(jìn)行單向流固耦合應(yīng)力分析，驗(yàn)證了其在極限工作轉(zhuǎn)速下的結(jié)構(gòu)安全性。通過比較原旋翼和優(yōu)化后旋翼的幾何外形可見，在低雷諾數(shù)下，適當(dāng)?shù)脑黾右硇蛷澏群腿~片安裝角的負(fù)扭能夠提高小型無人機(jī)旋翼的效率。