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

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

1 引言

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

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

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

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

2 理論基礎(chǔ)

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

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

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

2.2 CST翼型參數(shù)化

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

3.2 旋翼流場CFD仿真與實驗驗證

基于ANSYS ICEM進(jìn)行網(wǎng)格繪制，應(yīng)用ANSYS CFX進(jìn)行CFD仿真計算。流域網(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 仿真-實驗驗證結(jié)果Fig.6 Comparison of Thrust between CFD Results and Experimental Data

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

圖5 實驗臺展示Fig.5 Experimental Device

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

依照2.1章節(jié)中的參數(shù)化可以得到有關(guān)旋翼三維幾何形狀的22個設(shè)計參數(shù)，包括：6個結(jié)構(gòu)參數(shù)（記為Pairfoil）和16個葉尖和最大弦長處翼型的CST翼型控制參數(shù)（記為Pplatform）。隨后進(jìn)行拉丁超立方實驗設(shè)計，并進(jìn)行章節(jié)2.2中的敘述的CFD仿真得到訓(xùn)練集和測試集數(shù)據(jù)建立懸疑扭矩TQ和推力TH的擴展自適應(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)化算法的搜索時間；另外，還考慮了CST二維翼型參數(shù)化以及翼型最小厚度tmax的約束，以保證能夠產(chǎn)生真實有效的旋翼外形。

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

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

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

優(yōu)化前后參數(shù)對比，如表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對比了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)化模型的對比Fig.8 Comparison between the Optimized and Baseline Model

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

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

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

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

5.2 流固耦合驗證

為了確保優(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)過程中受到的空氣作用力（此處忽略葉片重力），葉片的運動與受力分析圖，如圖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é)論

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