高頻撲動(dòng)微撲翼飛行器多目標(biāo)優(yōu)化設(shè)計(jì)

蔣進(jìn)，鄭祥明，馮卓群，沈歡

(南京航空航天大學(xué) 飛行器先進(jìn)設(shè)計(jì)技術(shù)國防重點(diǎn)學(xué)科實(shí)驗(yàn)室，南京 210016)

0 引 言

微撲翼飛行器因其體積小、重量輕、隱蔽性好、成本低等特點(diǎn)，在偵察、監(jiān)測、搜救等領(lǐng)域具有較好的應(yīng)用前景，受到越來越多關(guān)注[1]。近年來，針對(duì)微撲翼飛行器的設(shè)計(jì)與優(yōu)化，國內(nèi)外做了相關(guān)的研究[2-3]。發(fā)現(xiàn)微撲翼機(jī)翼與機(jī)身最大載荷隨著撲動(dòng)頻率的增加急劇增加，嚴(yán)重影響飛行控制效果和飛行器壽命[4-5]。為解決飛行載荷過大的問題，有學(xué)者提出在設(shè)計(jì)過程中使用高強(qiáng)度、疲勞特性較好的材料以減小對(duì)飛行器壽命的影響[6]；同時(shí)也有學(xué)者提出對(duì)飛行器加裝減震裝置和使用多傳感器修正姿態(tài)的方法，以減小載荷對(duì)飛行控制的影響[7-8]。分析可知，上述解決方案一定程度上都會(huì)增加整機(jī)重量，且不能同時(shí)滿足要求。

本文針對(duì)一種高頻微撲翼飛行器，提出以改善載荷在時(shí)間域上分布情況為目標(biāo)的多目標(biāo)優(yōu)化模型，基于NSGA-Ⅱ算法以飛行器所受升力峰值和對(duì)翼根合力矩峰值最小為目標(biāo)，得到一組懸停狀態(tài)下的Pareto最優(yōu)解，以期為總體設(shè)計(jì)提供數(shù)據(jù)選擇依據(jù)。

1 微撲翼曲柄搖桿機(jī)構(gòu)設(shè)計(jì)及力學(xué)分析

1.1 微撲翼設(shè)計(jì)

設(shè)計(jì)的微撲翼如圖1(a)所示，由電機(jī)驅(qū)動(dòng)，轉(zhuǎn)動(dòng)經(jīng)帶輪和曲柄搖桿機(jī)構(gòu)轉(zhuǎn)換為機(jī)翼的擺動(dòng)運(yùn)動(dòng)(如圖1(b)所示)，再經(jīng)過迎角限制裝置使機(jī)翼保持一定迎角(如圖1(c)所示)，從而產(chǎn)生升力和阻力。分析可知，機(jī)翼的運(yùn)動(dòng)主要由曲柄搖桿機(jī)構(gòu)決定，并且電機(jī)和帶輪可用曲柄轉(zhuǎn)速代替，故選取曲柄搖桿桿機(jī)構(gòu)參數(shù)作為主要研究對(duì)象。

(a) 微撲翼實(shí)物圖

(b) 微撲翼設(shè)計(jì)圖

(c) 機(jī)翼迎角限制裝置

1.2 曲柄搖桿機(jī)構(gòu)運(yùn)動(dòng)模型

曲柄搖桿傳動(dòng)機(jī)構(gòu)簡圖如圖2所示，曲柄逆時(shí)針轉(zhuǎn)動(dòng)，O1，O2為機(jī)架上的支點(diǎn)；角θ1是曲柄與水平線的夾角，θ2是連桿與水平線的夾角，θ3是搖桿與水平線的夾角；l1是曲柄長度，l2是連桿長度，l3是搖桿長度，l0是支點(diǎn)間的距離，T為擺動(dòng)幅角。

圖2 曲柄搖桿機(jī)構(gòu)示意圖

設(shè)曲柄角速度為ω1，則曲柄與水平線夾角為

θ1(t)=ω1t

(1)

根據(jù)曲柄搖桿機(jī)構(gòu)各個(gè)參數(shù)關(guān)系得到搖桿與水平線夾角數(shù)學(xué)模型為[9-11]：

(2)

傳動(dòng)角γ公式如下：

(3)

搖桿的角速度和角加速度分別為θ3(t)對(duì)時(shí)間的一階導(dǎo)和二階導(dǎo)。

1.3 微撲翼氣動(dòng)力慣性力模型

為了使機(jī)翼始終保持較大迎角，設(shè)計(jì)機(jī)翼具有較大剛度，忽略柔性變形。機(jī)翼迎角α(t)≥60°，超過失速迎角，在背流面出現(xiàn)脫體現(xiàn)象，環(huán)量存在條件被破壞。通過計(jì)算可以得到機(jī)翼在運(yùn)動(dòng)時(shí)雷諾數(shù)較小，此時(shí)可以采用平板大迎角繞流氣動(dòng)力近似計(jì)算公式分析氣動(dòng)力[12]。

則易得機(jī)翼升力L(t)、阻力D(t)和機(jī)翼對(duì)翼根的升力矩ML(t)，阻力矩MD(t)為

L(t)

(4)

(5)

ML(t)

(6)

(7)

式中：u(x)=xω3(t)為距離機(jī)翼根部x長度處的來流速度；c(x)為弦長；x1和x2分別為機(jī)翼翼面的起點(diǎn)和終點(diǎn)；Cπ/2為當(dāng)迎角為π/2；雷諾數(shù)為104～106范圍內(nèi)時(shí)，垂直于流動(dòng)方向的二維平板的壓差阻力系數(shù)；b為修正系數(shù)，取奇數(shù)。

經(jīng)計(jì)算，機(jī)翼撲動(dòng)過程中雷諾數(shù)約為30 000，滿足雷諾數(shù)要求，考慮到機(jī)翼展弦比較小，同時(shí)在撲動(dòng)過程中會(huì)發(fā)生柔性變形，與理想情況差距較大，需進(jìn)行修正，此處Cπ/2取1.68[13]，b取1[12]。

機(jī)翼翻轉(zhuǎn)的過程迎角α隨時(shí)間改變，為了準(zhǔn)確模擬撲動(dòng)時(shí)狀態(tài)，使用高頻相機(jī)進(jìn)行拍攝，拍攝裝置擺放如圖3～圖4所示。

圖3 高頻拍攝實(shí)驗(yàn)裝置擺放示意圖

得到攝影視頻后，將單幀拍攝照片上機(jī)翼特征點(diǎn)處(AB=BC)弦長(BD)像素長度與展長(EF)像素長度對(duì)比，再根據(jù)實(shí)際展長即可得到機(jī)翼弦長的投影長度，結(jié)合實(shí)際長度即可得到機(jī)翼迎角。

圖4 拍攝視頻單幀

根據(jù)實(shí)驗(yàn)數(shù)據(jù)，進(jìn)行數(shù)據(jù)擬合，提出機(jī)翼翻轉(zhuǎn)迎角近似公式。

(8)

式中：mod為取余數(shù)運(yùn)算；ai、bi、ci等系數(shù)的取值如表1所示。

表1 擬合公式參數(shù)表

實(shí)驗(yàn)結(jié)果與近似曲線對(duì)比如圖5所示。

圖5 機(jī)翼迎角實(shí)驗(yàn)測量結(jié)果與近似曲線對(duì)比

慣性力主要由機(jī)翼撲動(dòng)角加速度引起，只受機(jī)翼重量影響。假設(shè)機(jī)翼質(zhì)量分布沿展向分布均勻，則機(jī)翼慣性力Finer(t)和對(duì)翼根處慣性力矩Miner(t)可以表示為[14]

(9)

(10)

式中：m為機(jī)翼質(zhì)量；a3(t)為機(jī)翼擺動(dòng)角加速度。

根據(jù)以上模型，可得到合力矩∑M(t)表達(dá)式為

(11)

2 多目標(biāo)優(yōu)化模型

本文優(yōu)化目標(biāo)為：最小化飛行過程中的機(jī)翼機(jī)身所受載荷峰值。經(jīng)過計(jì)算可得到在兩個(gè)撲動(dòng)周期內(nèi)的機(jī)翼機(jī)身所受力和力矩變化情況，如圖6所示。其中，升力在高頻撲動(dòng)飛行過程中對(duì)傳感器測量產(chǎn)生較大的噪聲，從而影響飛行控制，而合力矩主要影響飛行器壽命，故優(yōu)化目標(biāo)可以轉(zhuǎn)化為使飛行過程中升力峰值和合力矩峰值最小。同時(shí)為了統(tǒng)一量級(jí)，對(duì)應(yīng)量的絕對(duì)值平均值進(jìn)行無量綱化。選取曲柄搖桿機(jī)構(gòu)的四桿長度l1、l2、l3、l0作為優(yōu)化變量。約束函數(shù)對(duì)平均升力Lave、四桿長度關(guān)系、傳動(dòng)角γ和擺動(dòng)角T進(jìn)行限制。最終多目標(biāo)優(yōu)化數(shù)學(xué)模型如下：

(12)

設(shè)計(jì)變量及變化范圍：6mm≤l1≤7mm，16.5mm≤l2≤20.5mm，6.5mm≤l3≤8.5mm，17mm≤l0≤21mm。

約束函數(shù)：g1: Lave≥13.32g，g2: l0+l1≤l2+l3，g3: γ≥40°，g4: T≥105°。

其中g(shù)1為飛行器最低升力要求，13.32g為設(shè)計(jì)微撲翼飛行器全機(jī)重量的1/4；g2為四連桿曲柄存在約束條件；g3是曲柄機(jī)構(gòu)傳動(dòng)角限制，以保證傳動(dòng)效率；g4是擺動(dòng)角限制，105°由綜合多種高頻撲動(dòng)昆蟲鳥類翅膀擺動(dòng)角度得到[15]。

(a) 機(jī)翼機(jī)身受力變化曲線

(b) 機(jī)翼根部受力矩變化曲線

該優(yōu)化問題為一個(gè)多變量、多約束的非線性多目標(biāo)優(yōu)化問題。典型的多目標(biāo)進(jìn)化算法有：MOGA(基于排序的適應(yīng)度賦值多目標(biāo)遺傳算法)、NPGA(小生境Pareto遺傳算法)、NSGA(非劣解排序遺傳算法)等。其中NSGA算法計(jì)算性能較好，采用適應(yīng)度共享策略，有利于Pareto前沿上的個(gè)體分布均勻，維持種群多樣性，防止過早收斂。但由于其本身局限性，算法復(fù)雜度高達(dá)o(MN3)[16-19]。本文采用NSGA-Ⅱ算法提出的快速最優(yōu)非劣解排序算法，降低算法復(fù)雜度至o(MN2)，引入精英保留策略，大大提高計(jì)算速度[20]。

3 優(yōu)化結(jié)果與分析

運(yùn)用Matlab平臺(tái)編寫NSGA-Ⅱ程序，代入上述優(yōu)化模型，在撲動(dòng)頻率為48Hz的情況下，對(duì)升力和慣性力在撲動(dòng)周期內(nèi)的分布情況進(jìn)行優(yōu)化。

設(shè)置種群規(guī)模為100，交叉概率為0.9，變異概率為0.25，進(jìn)行1 000代進(jìn)化后得到該優(yōu)化的Pareto最優(yōu)解集[21]，如圖7所示，可以看出：兩個(gè)目標(biāo)函數(shù)之間是互相沖突的，在減小一個(gè)目標(biāo)函數(shù)的同時(shí)必然會(huì)犧牲另一個(gè)目標(biāo)函數(shù)，因此該優(yōu)化結(jié)果需要進(jìn)行取舍。

圖7 NSGA-Ⅱ算法計(jì)算結(jié)果

根據(jù)微撲翼飛行器設(shè)計(jì)的具體結(jié)構(gòu)尺寸要求，從最優(yōu)解集中選取出最合適的幾組解如表2所示。

表2 優(yōu)化前后設(shè)計(jì)參數(shù)對(duì)比

表3 優(yōu)化前后結(jié)果對(duì)比

為驗(yàn)證NSGA-Ⅱ算法優(yōu)化效果，選取表中優(yōu)化后第1組尺寸，得到優(yōu)化前后機(jī)翼擺動(dòng)角、擺動(dòng)角速度、升力、合力矩在撲動(dòng)頻率為48Hz(此時(shí)飛行器處于懸停狀態(tài))情況下的變化曲線對(duì)比，如圖8所示。

(a) 優(yōu)化前后機(jī)翼擺動(dòng)角對(duì)比

(b) 優(yōu)化前后機(jī)翼擺動(dòng)角速度對(duì)比

(c) 優(yōu)化前后升力曲線對(duì)比

(d) 優(yōu)化前后慣性力矩曲線對(duì)比

從圖8可以看出：優(yōu)化前后機(jī)翼擺動(dòng)幅度增大了2.0°，一定程度上增大了平均升力；機(jī)翼翼面第一次翻轉(zhuǎn)變得更為平緩，而第二次翻轉(zhuǎn)速度加快，使翻轉(zhuǎn)過程中產(chǎn)生的載荷分布更為均勻，同時(shí)減小了翻轉(zhuǎn)產(chǎn)生的最大載荷；優(yōu)化使得升力峰值明顯降低，升力分布更為均勻，減小飛行過程中由于升力變化引起的飛行性能下降。

為驗(yàn)證優(yōu)化結(jié)果對(duì)其他飛行頻率也適用，計(jì)算得到0～60Hz情況下優(yōu)化前后平均升力變化曲線(如圖9(a)所示)，升力和慣性力變化云圖(如圖9(b)所示)是優(yōu)化前升力云圖，優(yōu)化后升力云圖如圖9(c)所示，優(yōu)化前力矩云圖如圖9(d)所示，優(yōu)化后力矩云圖如圖9(e)所示。

(a) 優(yōu)化前后平均升力變化曲線

(b) 優(yōu)化前升力云圖

(c) 優(yōu)化后升力云圖

(d) 優(yōu)化前慣性力矩云圖

(e) 優(yōu)化后慣性力矩云圖

從圖9可以看出：優(yōu)化后平均升力曲線幾乎未變，滿足飛行升力需求；同時(shí)在其他頻率下，升力峰值與力矩峰值也得到降低，故優(yōu)化結(jié)果適用于其他撲動(dòng)頻率。

從表1、表2及圖9可以看出：表明優(yōu)化后的微撲翼飛行器飛行過程中所受最大載荷明顯降低，飛行性能和飛行壽命在一定程度上得到改善。后續(xù)設(shè)計(jì)過程可繼續(xù)使用該優(yōu)化方法，直至達(dá)到設(shè)計(jì)目的。

4 結(jié) 論

(1) 建立優(yōu)化模型，經(jīng)優(yōu)化后，升力峰值與慣性力峰值得到明顯降低，飛行器飛行過程中所受載荷分布更為均勻，提高了飛行性能，延長了飛行壽命。

(2) 該優(yōu)化方法還可為設(shè)計(jì)者提供多組最優(yōu)數(shù)據(jù)以供選擇，減輕勞動(dòng)強(qiáng)度，提高設(shè)計(jì)效率，符合工程設(shè)計(jì)需求。