薛 棟, 宋筆鋒, 宋文萍, 楊文青(西北工業(yè)大學(xué) 航空學(xué)院, 陜西 西安 710072)

0 引 言

自然界大約有100萬種昆蟲和9000種鳥類[1]，這些飛行生物向人類展示其卓越飛行能力的同時(shí)也啟發(fā)著人類研制類似的撲翼飛行器。

在最近20年里，人們已經(jīng)成功研制出一些仿昆蟲型和仿鳥型的撲翼飛行器(如圖1所示)，如美國(guó)哈佛大學(xué)的“果蠅”飛行器[2]、德國(guó)Festo 公司的“Smart bird”[3]、中國(guó)南京航空航天大學(xué)的“金鷹”[4]和西北工業(yè)大學(xué)的“信鴿”等[5]。

(a) 美國(guó)“果蠅” (b) 德國(guó)“Smart bird”

(c) 南航“金鷹” (d) 西工大“信鴿”

圖1國(guó)內(nèi)外具有代表性的撲翼飛行器

Fig.1RepresentativeFMAVs

未來的仿鳥型微型飛行器具有以下特點(diǎn)：氣動(dòng)性能好(能主動(dòng)適應(yīng)和利用不同氣流環(huán)境)、機(jī)動(dòng)性強(qiáng)(可穿越狹小復(fù)雜空間)、飛行效率高(低能量消耗獲得長(zhǎng)距離飛行)、易于隱身(采用仿生飛行的方式迷惑敵方)等優(yōu)勢(shì)，因而可克服微型固定翼及旋翼無人機(jī)成本高、能量利用率低、隱身能力弱等缺陷，具有巨大的發(fā)展?jié)摿蛻?yīng)用前景。

目前的仿鳥型微型飛行器主要從飛行方式進(jìn)行仿生，即單純靠撲動(dòng)翼?yè)鋭?dòng)同時(shí)產(chǎn)生升力和推力進(jìn)行飛行，但氣動(dòng)效率與鳥類相比仍然有很大差距，此外在傳感、驅(qū)動(dòng)和控制等方面鳥類飛行機(jī)理的研究仍然處于初級(jí)階段，所以亟需進(jìn)一步的深入研究才能使仿鳥型微型飛行器具有像鳥類一樣卓越的飛行能力。本文主要從氣動(dòng)、結(jié)構(gòu)與飛行穩(wěn)定性等方面進(jìn)行分析和總結(jié)。

仿鳥型撲翼飛行器和仿昆蟲型撲翼飛行器雖然都屬于撲翼飛行器，但它們的飛行機(jī)理卻有著很大區(qū)別。

如圖2[6]所示，仿鳥型撲翼飛行器像鳥類一樣以平飛為主要飛行方式，撲動(dòng)頻率一般在10～20 Hz，升力主要來源于前飛速度產(chǎn)生的升力。而仿昆蟲型撲翼飛行器則與昆蟲類似，其主要飛行狀態(tài)是懸停，撲動(dòng)頻率為20～600 Hz[7]，它們的升力主要靠撲動(dòng)翼的純撲動(dòng)產(chǎn)生。由于其具有很高的撲動(dòng)頻率，且撲動(dòng)翼質(zhì)量占總重量的比例較小(小于5%)，所以機(jī)體受撲動(dòng)翼的影響較小，飛行過程中不會(huì)產(chǎn)生明顯振動(dòng)，所以氣動(dòng)、結(jié)構(gòu)與飛行力學(xué)之間的耦合問題不是很突出。相反，仿鳥型撲翼飛行器的低撲動(dòng)頻率和高撲動(dòng)翼質(zhì)量占比，使得氣動(dòng)、結(jié)構(gòu)與飛行力學(xué)的耦合關(guān)系很強(qiáng)，是設(shè)計(jì)高效、高性能仿鳥型撲翼飛行器不能忽略的問題。

如圖3所示，仿鳥型撲翼飛行器的研制是一個(gè)多學(xué)科設(shè)計(jì)過程，目前仍面臨著以下問題：1) 鳥類和仿鳥型撲翼飛行的非定常氣動(dòng)機(jī)理尚不明確；2) 柔性撲動(dòng)翼結(jié)構(gòu)形式復(fù)雜、撲動(dòng)過程中會(huì)產(chǎn)生大幅變形等特點(diǎn)，針對(duì)該類氣動(dòng)結(jié)構(gòu)耦合問題目前缺少既高效又相對(duì)準(zhǔn)確的計(jì)算方法；3) 缺乏高精度、高效的氣動(dòng)、結(jié)構(gòu)與飛行動(dòng)力學(xué)耦合的計(jì)算方法，不能滿足控制律設(shè)計(jì)和優(yōu)化設(shè)計(jì)的輸入需求。只有解決好上述各個(gè)學(xué)科的難點(diǎn)問題和多學(xué)科之間的耦合問題，才能設(shè)計(jì)出像鳥類一樣具有卓越飛行能力的撲翼飛行器。

圖3 仿鳥型撲翼飛行器多學(xué)科設(shè)計(jì)過程Fig.3 Multidisciplinary design of bird inspired FMAV

本文的主要內(nèi)容包括：仿鳥型撲翼飛行器的氣動(dòng)機(jī)理研究現(xiàn)狀，包括非定常氣動(dòng)機(jī)理方面；氣動(dòng)結(jié)構(gòu)耦合研究現(xiàn)狀，及其尺度律的研究現(xiàn)狀；動(dòng)力學(xué)建模方法、飛行穩(wěn)定性分析方法的研究現(xiàn)狀；以及通過耦合氣動(dòng)、結(jié)構(gòu)與飛行力學(xué)來研究撲動(dòng)翼柔性對(duì)飛行穩(wěn)定性影響的研究進(jìn)展。

1 仿鳥型撲翼飛行器氣動(dòng)機(jī)理

1.1 無量綱氣動(dòng)參數(shù)

為方便論述，先介紹4個(gè)常用的描述撲翼飛行的無量綱氣動(dòng)參數(shù)。

雷諾數(shù)Re描述的是流體的慣性力和黏性力的比值，定義如下：

Re=ρfLrefUref/μ

(1)

式中：ρf為流體介質(zhì)密度，Lref為參考長(zhǎng)度，Uref為參考速度，μ為流體介質(zhì)的黏性系數(shù)。

自然界飛行生物的雷諾數(shù)范圍橫跨了6個(gè)數(shù)量級(jí)，從最小飛行昆蟲接近1的雷諾數(shù)到潛水獵鷹的100萬雷諾數(shù)[7]。鳥類和仿鳥型撲翼飛行器的雷諾數(shù)通常大于10000，此時(shí)邊界層容易經(jīng)歷層流/湍流轉(zhuǎn)捩，這些流動(dòng)狀態(tài)對(duì)該類飛行的氣動(dòng)特性起著關(guān)鍵作用。

前進(jìn)比J反映的是尾渦脫落對(duì)推進(jìn)效率的影響，定義如下：

J=U∞/(2fARcmφ)

(2)

式中：U∞為來流速度，f為撲動(dòng)頻率，AR為展弦比，cm為撲動(dòng)翼的平均氣動(dòng)弦長(zhǎng)，φ為撲動(dòng)幅度，當(dāng)J=0時(shí)表示撲翼飛行的懸停狀態(tài)。根據(jù)研究結(jié)果[8]，自然界中飛行生物前進(jìn)比的范圍約為0～5，具有較高推進(jìn)效率的前進(jìn)比范圍為1.25～2.5。

羅斯比數(shù)Ro代表的是慣性力和科氏力的比值，該參數(shù)影響著前緣渦的穩(wěn)定性，定義如下：

(3)

自然界中飛行生物的羅斯比數(shù)范圍是1～10, 前緣渦穩(wěn)定的條件是要求羅斯比數(shù)一般小于4[9]。

減縮頻率k表征的是撲動(dòng)運(yùn)動(dòng)的非定常程度，定義為流體擾動(dòng)的空間尺度與弦長(zhǎng)的比值：

k=πfcm/Uref

(4)

上述四個(gè)無量綱參數(shù)分別決定著仿鳥型撲翼飛行器在飛行過程中湍流的發(fā)生(雷諾數(shù)Re)、尾渦的形態(tài)(前進(jìn)比J)、前緣渦的穩(wěn)定性(羅斯比數(shù)Ro)，以及飛行過程中的非定常性(減縮頻率k)。

1.2 研究現(xiàn)狀

為能研制出像鳥一樣高飛行效率的撲翼飛行器，首先需要研究清楚鳥類飛行的空氣動(dòng)力學(xué)機(jī)理。

1909年Knoller和Betz發(fā)現(xiàn)二維翼型在撲動(dòng)過程中可以同時(shí)產(chǎn)生升力和推力，這一發(fā)現(xiàn)被稱為Knoller-Betz效應(yīng)[10-11]。Von Karman和 Burgers[12]首次提出撲動(dòng)翼型產(chǎn)生阻力或推力與尾渦的位置和方向有關(guān)，阻力和推力相應(yīng)的尾渦分別被稱為卡門渦街和反卡門渦街。1953年，Bratt[13]通過實(shí)驗(yàn)驗(yàn)證了 Von Karman 和 Burgers 的結(jié)論。

由于鳥類和仿鳥型撲翼飛行器的減縮頻率通常小于0.3，此時(shí)認(rèn)為非定常氣動(dòng)機(jī)理的影響不是很顯著[14]，可以使用準(zhǔn)定常的方法(Quasi-steady method)對(duì)這類撲翼飛行的氣動(dòng)力進(jìn)行計(jì)算。例如激勵(lì)盤理論[15]和Peters提出的Induced Inflow Theory(誘導(dǎo)入流理論)[16]。

雖然仿鳥型撲翼飛行器與仿昆蟲型撲翼飛行器在飛行機(jī)理上有很大區(qū)別，但針對(duì)兩者的數(shù)值計(jì)算方法和實(shí)驗(yàn)方法卻是相近相通的。由于以昆蟲為研究對(duì)象的數(shù)值計(jì)算和實(shí)驗(yàn)研究已經(jīng)很全面和深入，因此在分析仿鳥型撲翼飛行器氣動(dòng)特性時(shí)可以參照和借鑒，例如在數(shù)值計(jì)算方面。早在1998年，Liu等[17]就通過求解基于時(shí)間推進(jìn)的三維非定常Navier-Stokes方程，計(jì)算了昆蟲飛行的非定常流場(chǎng)結(jié)構(gòu)，并通過實(shí)驗(yàn)進(jìn)行了驗(yàn)證。其他針對(duì)昆蟲的氣動(dòng)特性的數(shù)值計(jì)算和實(shí)驗(yàn)研究還有很多，可以參考Shyy[6,18]等的綜述文章，這里不再贅述。

國(guó)內(nèi)對(duì)仿鳥型撲翼飛行過程中的空氣動(dòng)力學(xué)問題的研究起步較晚。在2003年, 南京航空航天大學(xué)的曾銳和昂海松[19]基于簡(jiǎn)化的綠頭鴨撲動(dòng)模型研究了撲動(dòng)幅度和撲動(dòng)頻率等參數(shù)對(duì)平均升力系數(shù)的影響。同年，西北工業(yè)大學(xué)的龔凱[20]通過求解Euler方程研究了不同的撲動(dòng)參數(shù)對(duì)氣動(dòng)特性的影響。2006年，楊淑利等基于改進(jìn)的片條理論方法研究了撲動(dòng)翼氣動(dòng)性能[21]。2008年，謝輝等[22]通過求解三維非定常N-S方程，開展了微型撲動(dòng)翼的數(shù)值模擬計(jì)算。2010年，楊文青等基于距離減縮法的嵌套網(wǎng)格技術(shù)，實(shí)現(xiàn)了對(duì)微型撲翼飛行器全機(jī)氣動(dòng)特性的研究[23]。

在國(guó)內(nèi)的撲翼氣動(dòng)實(shí)驗(yàn)方面，2007年，昂海松等針對(duì)撲翼微型飛行器開展了風(fēng)洞實(shí)驗(yàn)工作，研究了柔性撲動(dòng)翼產(chǎn)生升力和推力的流動(dòng)機(jī)理，為柔性撲翼飛行器的設(shè)計(jì)提供了思路[24]。同年，西北工業(yè)大學(xué)邵立民等開展了初步的微型撲翼飛行器風(fēng)洞實(shí)驗(yàn)[25]。2008年，王利光等通過對(duì)風(fēng)洞實(shí)驗(yàn)系統(tǒng)的改進(jìn)，研究了撲翼飛行器氣動(dòng)力和功耗的影響因素[26]。隨后在2014年，王利光等在風(fēng)洞實(shí)驗(yàn)系統(tǒng)的基礎(chǔ)上引進(jìn)了高速攝像系統(tǒng)，實(shí)現(xiàn)了柔性撲動(dòng)翼氣動(dòng)力和結(jié)構(gòu)變形量的同步測(cè)量[27]。

1.3 非定常氣動(dòng)機(jī)理

雖然昆蟲在飛行過程中有很多非定常氣動(dòng)機(jī)理(Clap and Fling機(jī)制，快速旋轉(zhuǎn)和尾跡捕捉等)[18]，但是在鳥類飛行過程中仍沒有發(fā)現(xiàn)類似的氣動(dòng)機(jī)理[28]。

隨著最近幾年關(guān)于前緣渦對(duì)鳥類和仿鳥型撲翼飛行影響研究的深入，發(fā)現(xiàn)前緣渦對(duì)該類飛行過程中氣動(dòng)力的產(chǎn)生不能被忽視，例如有研究發(fā)現(xiàn)有的鳥類在慢速前飛過程中前緣渦可以提供高達(dá)49%的升力[29]，如圖4所示。

圖4 一種叫Ficedula hypoleuca鳥在飛行速度為1 m/s時(shí)翅膀外段的前緣渦PIV云圖，渦量大小范圍從-2000/s (紅)至2000/s (白) [29]Fig.4 Air movements induced by the flapping wing of a flycatcher flying at V∞=1 m/s at the hand wing section. Color bars in PIV images refer to vorticity ranges, scaled from -2000/s (red) to 2000/s (white)[29]

那么影響自然界中飛行生物前緣渦形態(tài)有哪些因素呢？其中一個(gè)重要因素就是雷諾數(shù)。不同種類的飛行生物具有不同的雷諾數(shù)范圍，如圖5所示。在圖5的左上方分別顯示了在雷諾數(shù)為120和1400時(shí)的翼型剖面的速度場(chǎng)，同時(shí)在其上面疊加了展向速度的云圖，云圖范圍從0 m/s(白色)至0.5 m/s(黑色)，可以看出當(dāng)處于昆蟲的雷諾數(shù)范圍時(shí)(2000以下)，前緣渦的展向速度隨著雷諾數(shù)的增加而增加，而展向速度直接影響前緣渦的穩(wěn)定性[9]。圖5右上方則顯示當(dāng)雷諾數(shù)高于1400時(shí)，可以發(fā)現(xiàn)前緣渦開始出現(xiàn)分叉和破裂。

圖5 雷諾數(shù)對(duì)前緣渦的影響[9]Fig.5 Reynolds number (Re) effect on LEV formation [9]

很多研究者發(fā)現(xiàn)昆蟲能夠保持前緣渦的穩(wěn)定，從而使失速延遲[18]。不同于昆蟲，鳥類或仿鳥型撲翼飛行器的雷諾數(shù)通常大于10000，此時(shí)前緣渦會(huì)對(duì)氣動(dòng)特性有怎樣的影響呢？Hubel等[30]通過風(fēng)洞實(shí)驗(yàn)對(duì)一個(gè)仿鳥模型進(jìn)行了測(cè)試，發(fā)現(xiàn)在鳥類飛行的雷諾數(shù)和減縮頻率范圍內(nèi)普遍出現(xiàn)了前緣渦。并且前緣渦的出現(xiàn)、脫落和分離與雷諾數(shù)、減縮頻率有很大關(guān)系。由于展向上不同站位有著不同的流動(dòng)狀態(tài)，所以前緣渦對(duì)升力的有益影響和流動(dòng)分離造成失速的不利影響可能會(huì)相互抵消。同時(shí)前緣渦會(huì)增大阻力，即削弱撲動(dòng)產(chǎn)生的推力作用。此外，與在昆蟲飛行過程中發(fā)現(xiàn)穩(wěn)定的前緣渦現(xiàn)象不同，實(shí)驗(yàn)中仿鳥模型的前緣渦極易脫落。

雖然Hubel的實(shí)驗(yàn)沒有得到穩(wěn)定的前緣渦增升機(jī)制，但自然界中的鳥類有可能通過主動(dòng)控制翅膀展向的扭轉(zhuǎn)、彎度和有效攻角等參數(shù)實(shí)現(xiàn)對(duì)前緣渦的控制和利用。希望今后通過對(duì)這一方面的深入研究，實(shí)現(xiàn)對(duì)仿鳥型撲翼飛行器氣動(dòng)效率的提高。

2 氣動(dòng)結(jié)構(gòu)耦合研究

2.1 國(guó)內(nèi)外研究現(xiàn)狀

昆蟲和鳥類在飛行過程中，利用翅膀的有益變形可以在降低功耗的同時(shí)提高推力[31]。利用這一點(diǎn)，撲翼飛行器的撲動(dòng)翼通常采用柔性輕質(zhì)結(jié)構(gòu)。而輕質(zhì)結(jié)構(gòu)還可以降低撲動(dòng)翼的慣性力帶來的功耗。

國(guó)外研究者在柔性撲動(dòng)翼氣動(dòng)結(jié)構(gòu)耦合的數(shù)值計(jì)算和實(shí)驗(yàn)方面開展了很多工作。Shyy等[32]在2010年針對(duì)柔性撲動(dòng)翼氣動(dòng)結(jié)構(gòu)耦合的數(shù)值計(jì)算和實(shí)驗(yàn)方面進(jìn)行了比較詳細(xì)和全面的總結(jié)，這里僅介紹具有代表意義的研究工作。

在實(shí)驗(yàn)研究方面比較有代表意義的是Heathcote等人[33]和Wu等[34]的工作。Heathcote等通過水洞實(shí)驗(yàn)研究了柔性撲翼作沉浮運(yùn)動(dòng)時(shí)，弦向和展向的剛度對(duì)推力特性的影響。結(jié)果表明，柔性撲翼比剛性撲翼能產(chǎn)生更大的推力。該實(shí)驗(yàn)數(shù)據(jù)也經(jīng)常被用來驗(yàn)證數(shù)值計(jì)算的正確性。

Wu等通過更復(fù)雜的實(shí)驗(yàn)手段研究了撲動(dòng)翼柔性變形與氣動(dòng)特性的關(guān)系。具體實(shí)驗(yàn)方法是通過數(shù)字圖像相關(guān)技術(shù)測(cè)量結(jié)構(gòu)變形，利用六分量天平測(cè)量撲翼所受到的氣動(dòng)力和慣性力，流場(chǎng)信息則通過數(shù)字粒子圖像測(cè)速技術(shù)(DPIV)獲取。他們的研究結(jié)果表明：撲動(dòng)翼的結(jié)構(gòu)特性對(duì)氣動(dòng)特性的影響很大，過柔和過剛的撲動(dòng)翼結(jié)構(gòu)都不能得到好的氣動(dòng)特性。結(jié)構(gòu)特性與氣動(dòng)性能之間的具體關(guān)系將在下文中的尺度律的研究中進(jìn)行介紹。

在撲動(dòng)翼氣動(dòng)/結(jié)構(gòu)耦合的數(shù)值計(jì)算方面，Kim等[35]基于修正的片條理論模型開發(fā)了氣動(dòng)結(jié)構(gòu)耦合計(jì)算模型，并對(duì)柔性撲動(dòng)翼的氣動(dòng)特性進(jìn)行了數(shù)值模擬，計(jì)算結(jié)果與實(shí)驗(yàn)值吻合良好。Unger等[36]通過求解非定常雷諾平均N-S方程，并耦合基于有限元方法的非線性結(jié)構(gòu)求解器，研究了僅具有弦向柔性的撲動(dòng)翼的氣動(dòng)特性。

國(guó)內(nèi)對(duì)撲翼飛行器氣動(dòng)結(jié)構(gòu)耦合研究方面，西北工業(yè)大學(xué)的楊文青等[37]基于非定常雷諾平均N-S方程耦合靜態(tài)結(jié)構(gòu)方程，研究柔性撲動(dòng)翼的氣動(dòng)特性。陳麗利等[38-39]基于結(jié)構(gòu)動(dòng)力學(xué)方法開發(fā)了柔性撲動(dòng)翼的CFD/CSD耦合求解程序，并在計(jì)算過程中考慮了結(jié)構(gòu)慣性力的影響，實(shí)現(xiàn)了對(duì)柔性撲翼氣動(dòng)特性的更準(zhǔn)確模擬。

2.2 氣動(dòng)/結(jié)構(gòu)耦合的尺度律分析

由于柔性撲翼飛行的氣動(dòng)結(jié)構(gòu)耦合計(jì)算不僅涉及到復(fù)雜的非定常氣動(dòng)計(jì)算，還涉及各式各樣的撲動(dòng)翼氣動(dòng)外形、結(jié)構(gòu)布局以及運(yùn)動(dòng)參數(shù)，這些參數(shù)的繁雜增加了柔性撲動(dòng)翼氣動(dòng)性能分析的難度。而尺度律的研究可以較好地解決這個(gè)問題。

Kang等[31]對(duì)柔性撲動(dòng)翼的氣動(dòng)、結(jié)構(gòu)進(jìn)行了量綱分析，并引入三個(gè)無量綱的結(jié)構(gòu)參數(shù)用來分析撲動(dòng)翼的柔性對(duì)氣動(dòng)特性的影響。

其中參數(shù)ρ*表示結(jié)構(gòu)密度與流體密度的比值，定義如下:

ρ*=ρs/ρf

(5)

式中，ρs為各向同性材料組成的撲動(dòng)平板結(jié)構(gòu)的密度，ρf為流體介質(zhì)的密度。

h*表示撲動(dòng)翼厚度與參考弦長(zhǎng)的比值，定義如下：

h*=h/cm

(6)

式中，h為撲動(dòng)翼的厚度，cm為撲動(dòng)翼的平均氣動(dòng)弦長(zhǎng)。ρ*h*用來表征撲動(dòng)翼結(jié)構(gòu)的慣性力與流體慣性力的比例關(guān)系。

f*表示動(dòng)頻率與撲動(dòng)翼結(jié)構(gòu)的一階固有頻率的比值：

f*=f/f1

(7)

式中，f為撲動(dòng)頻率，f1為撲動(dòng)翼結(jié)構(gòu)的一階固有頻率。當(dāng)該比值很小時(shí)表示撲動(dòng)翼接近剛性，此時(shí)撲動(dòng)翼的變形很小，其位移主要由撲動(dòng)運(yùn)動(dòng)決定。

圖6 無量綱的變形量與平均升力的關(guān)系[40]Fig 6 Scaling law betweenγ and [40]

(8)

(9)

根據(jù)圖6可知，無論是昆蟲的翅膀還是人造的撲動(dòng)翼，它們的氣動(dòng)特性的表征參數(shù)與結(jié)構(gòu)特性的表征參數(shù)呈現(xiàn)很好的線性相關(guān)。因此，研究人員可以通過這種簡(jiǎn)單關(guān)系實(shí)現(xiàn)對(duì)柔性撲動(dòng)翼的快速性能分析。

由于該方法在計(jì)算過程只考慮了結(jié)構(gòu)的一階模態(tài)，氣動(dòng)計(jì)算方面只考慮了附加質(zhì)量力的影響，因此計(jì)算精度有待進(jìn)一步提高。此外該方法只考慮了各向同性材料構(gòu)成的平板結(jié)構(gòu)。由于撲動(dòng)翼的展向剛度通常需要比弦向的剛度大1～2個(gè)數(shù)量級(jí)[41-42]，因此人造撲動(dòng)翼一般通過碳桿-薄膜結(jié)構(gòu)來實(shí)現(xiàn)這一剛度分布要求，因此需要對(duì)這種各向異性材料構(gòu)成的復(fù)雜結(jié)構(gòu)的撲動(dòng)翼進(jìn)行尺度律分析。

3 撲翼飛行器的飛行動(dòng)力學(xué)研究

3.1 撲翼飛行器的運(yùn)動(dòng)方程

Gebert和Gallmeier[43]首次推導(dǎo)了撲翼飛行的運(yùn)動(dòng)方程，并考慮了撲動(dòng)翼慣性力、機(jī)體慣性力和慣性矩的變化影響。孫茂等[44]發(fā)現(xiàn)方程中的一些錯(cuò)誤，并且重新進(jìn)行了推導(dǎo)，得到：

(10)

(11)

3.2 仿鳥型撲翼飛行器動(dòng)力學(xué)建模過程中考慮撲動(dòng)翼慣性力的必要性

隨著飛行器的尺寸和質(zhì)量增加，仿鳥型撲翼飛行器的撲動(dòng)翼質(zhì)量占總重的比例通常會(huì)超過昆蟲和仿昆蟲型撲翼飛行器。那么質(zhì)量占比具體為多少時(shí)就需要考慮撲動(dòng)翼慣性力的影響呢？

Orlowski和Girard[45]對(duì)此問題展開了研究，發(fā)現(xiàn)如果直接忽略撲動(dòng)翼的慣性力可能在動(dòng)力學(xué)和控制問題中導(dǎo)致不一樣甚至完全錯(cuò)誤的結(jié)果。

Caetano等對(duì)此問題開展了定量研究[46]。研究對(duì)象選取了兩種布局的飛行器，一種是X布局四翼?yè)湟盹w行器“DelFly II”，另一種是雙翼?yè)湟盹w行器。同時(shí)采用了兩種動(dòng)力學(xué)模型進(jìn)行建模，一種是單剛體運(yùn)動(dòng)方程，另一種是考慮撲動(dòng)翼慣性力的多體動(dòng)力學(xué)模型?；谶@兩種模型，并通過飛行實(shí)驗(yàn)得到的位移和姿態(tài)信息來辨識(shí)飛行過程中的氣動(dòng)特性。

圖7為不同布局飛行器在垂直于機(jī)身軸線方向的升力的相關(guān)系數(shù)隨著質(zhì)量系數(shù)的變化歷程。兩組數(shù)據(jù)的相關(guān)系數(shù)越接近1，代表兩組數(shù)據(jù)之間誤差越?。环粗?，越接近0，代表兩組數(shù)據(jù)之間的誤差越大。

由圖7可知，具有對(duì)稱撲動(dòng)的四翼布局的相關(guān)系數(shù)比傳統(tǒng)布局的兩翼布局的要高，說明單剛體動(dòng)力學(xué)模型在四翼布局的適用范圍要廣，這是因?yàn)樗囊聿季值娘w行器在撲動(dòng)過程中會(huì)產(chǎn)生對(duì)稱的垂直于機(jī)體軸線方向的力，由于方向相反，會(huì)相互抵消，從而減少撲動(dòng)翼的慣性力對(duì)機(jī)體運(yùn)動(dòng)的影響。

圖7 不同布局的Z向氣動(dòng)力的相關(guān)系數(shù)隨著撲動(dòng)翼質(zhì)量占比的變化 [46]Fig.7 Coefficient of determination of the Z force with mass ratio for both configurations [46]

此外，隨著撲動(dòng)翼質(zhì)量占比的增加，兩種方法獲得的升力之間的相關(guān)系數(shù)逐漸下降，說明誤差越來越大。說明質(zhì)量占比越大時(shí)用單自由度的剛體運(yùn)動(dòng)模型得到的結(jié)果越不準(zhǔn)確，需要用較為復(fù)雜但更為精確的多體動(dòng)力學(xué)模型進(jìn)行建模分析。

3.3 撲翼飛行器動(dòng)力學(xué)特性分析

由(10)和式(11)可知，撲翼飛行器縱向的運(yùn)動(dòng)方程是一個(gè)非線性周期系統(tǒng)，目前有三種方法分析該系統(tǒng)的穩(wěn)定性。

第一種是直接平均法，常用于分析昆蟲或仿昆蟲撲翼飛行器的飛行動(dòng)力學(xué)特性[47-54]。

由于研究對(duì)象的撲動(dòng)頻率較高并且具有較小的撲動(dòng)翼質(zhì)量占比，在推導(dǎo)撲翼飛行器的動(dòng)力學(xué)方程時(shí)通常基于以下三個(gè)假設(shè)和方法[44]：1) 由于撲動(dòng)翼的撲動(dòng)周期遠(yuǎn)遠(yuǎn)小于機(jī)體運(yùn)動(dòng)的時(shí)間尺度，從而可以將非自治的周期系統(tǒng)在撲動(dòng)周期內(nèi)進(jìn)行平均，將系統(tǒng)方程變?yōu)樽灾蔚某Ｎ⒎址匠蹋?) 作用在機(jī)體上的撲動(dòng)翼慣性力引起的機(jī)體振動(dòng)量為小量，其高階量可以忽略；3) 撲動(dòng)翼產(chǎn)生慣性力關(guān)于撲動(dòng)周期的平均值近似為0，從而可以在方程中略去。

基于上述假設(shè)，縱向運(yùn)動(dòng)方程由非線性時(shí)間周期系統(tǒng)變?yōu)榱朔蔷€性的時(shí)不變系統(tǒng)，進(jìn)一步線化處理可以得到線性時(shí)不變系統(tǒng)，而后者的穩(wěn)定特性可以通過該線性系統(tǒng)矩陣的特征值獲得。

第二種分析方法是將非線性的周期系統(tǒng)變?yōu)榫€性周期系統(tǒng)，然后用Floquet方法去分析。Dietl和Garcia[55]第一次使用了該方法分析撲翼飛行器的穩(wěn)定性。他們的結(jié)果和Taylor等使用直接平均方法的結(jié)果類似，可以得到一個(gè)穩(wěn)定的振蕩特征模態(tài)、一個(gè)衰減特征模態(tài)和一個(gè)發(fā)散特征模態(tài)。

Su和Cesnik[56]同樣使用該方法分析了柔性撲翼飛行器的穩(wěn)定性。在撲翼飛行器的縱向方向，他們分別得到一個(gè)不穩(wěn)定的振蕩模態(tài)、一個(gè)發(fā)散模態(tài)以及一個(gè)衰減模態(tài)。他們還研究了系統(tǒng)的模態(tài)隨著撲動(dòng)翼的剛度和慣性力增加的變化趨勢(shì)，并發(fā)現(xiàn)撲動(dòng)翼的慣性力具有不穩(wěn)定效果，尤其在衰減模態(tài)中。

第三種是Taha等提出的高階平均方法[57]。相比于Floquet理論，平均理論最主要的優(yōu)勢(shì)是省去了繁瑣的配平過程。他們使用二階平均方法分析了撲翼飛行器的縱向穩(wěn)定性，發(fā)現(xiàn)對(duì)具有較低撲動(dòng)頻率系統(tǒng)，直接平均方法不能夠得到正確結(jié)論，當(dāng)撲動(dòng)頻率與機(jī)體的運(yùn)動(dòng)頻率之比大于100時(shí)才可以使用。

那么上述三種方法的適用范圍和對(duì)象究竟有什么區(qū)別呢？Wu和Sun[58]對(duì)此展開了研究，他們用直接平均方法和Floquet方法分別對(duì)兩種具有不同翅膀質(zhì)量占比的昆蟲(蚜蠅和鷹蛾)進(jìn)行了穩(wěn)定性分析。

由于蚜蠅具有相對(duì)高的撲動(dòng)頻率和小的翅膀質(zhì)量占比，因此在飛行過程中具有小的機(jī)體運(yùn)動(dòng)，鷹蛾則具有相對(duì)較低的撲動(dòng)頻率和較高的翅膀質(zhì)量占比，飛行過程中具有相對(duì)較大幅度的機(jī)體運(yùn)動(dòng)。最后他們采用N-S方程耦合多體動(dòng)力學(xué)方程進(jìn)行了驗(yàn)證。計(jì)算結(jié)果表明，對(duì)于蚜蠅這種研究對(duì)象該兩種方法得到的結(jié)果幾乎一致，而對(duì)于鷹蛾，直接平均方法則出現(xiàn)較大的誤差。說明對(duì)于撲動(dòng)頻率較低并且翅膀質(zhì)量占比較大的昆蟲，平均方法不能滿足要求，需要用考慮周期運(yùn)動(dòng)的Floquet方法。

以上方法針對(duì)的研究對(duì)象都為昆蟲或仿昆蟲型撲翼飛行器，由于仿鳥型撲翼飛行器的撲動(dòng)頻率更低并且翅膀質(zhì)量占比更大，因此在針對(duì)該類飛行器進(jìn)行飛行穩(wěn)定性分析時(shí)應(yīng)該使用高階平均方法或Floquet方法。

4 撲翼的柔性對(duì)飛行穩(wěn)定性的影響

有些研究在分析昆蟲或微型撲翼飛行器的飛行穩(wěn)定性時(shí), 為了分析方便，通常假設(shè)撲動(dòng)翼在氣動(dòng)力的作用下沒有變形，或者認(rèn)為變形不會(huì)對(duì)昆蟲或者撲翼飛行器的穩(wěn)定性產(chǎn)生重要影響[44-45,59-60]。但是最近的一些研究表明，撲動(dòng)翼的柔性不僅在氣動(dòng)分析中具有不可忽視的重要性[5,31]，而且柔性對(duì)撲翼飛行器的穩(wěn)定性也有重要影響。而對(duì)這一內(nèi)容的研究需要通過氣動(dòng)、結(jié)構(gòu)與飛行力學(xué)的耦合實(shí)現(xiàn)。

Richter和Patil[61]使用不同剛度的撲動(dòng)翼研究了撲動(dòng)翼柔性對(duì)撲翼飛行器穩(wěn)定性的影響。研究對(duì)象是Wood等設(shè)計(jì)的Robofly[2]。他們首先建立了撲翼飛行器氣動(dòng)結(jié)構(gòu)耦合飛行動(dòng)力學(xué)的模型。其中運(yùn)動(dòng)方程由拉格朗日方程推導(dǎo)，并使用低階的有限元模型來對(duì)撲動(dòng)翼的結(jié)構(gòu)進(jìn)行建模，使用準(zhǔn)定常氣動(dòng)理論計(jì)算氣動(dòng)力。流動(dòng)非定常效應(yīng)通過考慮簡(jiǎn)單的誘導(dǎo)入流理論。使用周期打靶法確定出三種不同柔性程度的撲動(dòng)翼和剛性翼在懸停狀態(tài)的周期解，然后使用Floquet理論對(duì)這四種情況進(jìn)行了穩(wěn)定性分析。

圖8為四種布局對(duì)應(yīng)系統(tǒng)的縱向特征值。圖中可以看出系統(tǒng)在縱向由四個(gè)非周期的兩個(gè)穩(wěn)定的和兩個(gè)不穩(wěn)定的指數(shù)型運(yùn)動(dòng)組成，沒有周期振蕩運(yùn)動(dòng)。其中特征模態(tài)1主要由俯仰運(yùn)動(dòng)和水平運(yùn)動(dòng)構(gòu)成。柔性布局對(duì)應(yīng)的半幅時(shí)間普遍要小于剛性翼的，且柔性程度越大，改變的程度越明顯。與特征模態(tài)1相反，四種布局的特征模態(tài)4是一個(gè)不穩(wěn)定的指數(shù)發(fā)散運(yùn)動(dòng)，而柔性布局比剛性翼布局更不穩(wěn)定。對(duì)于特征模態(tài)2和3，柔性的影響則不明顯。最后得出結(jié)論：柔性會(huì)放大穩(wěn)定和不穩(wěn)定特性，并且降低振蕩模態(tài)的阻尼作用。

圖8 四種結(jié)構(gòu)布局對(duì)應(yīng)系統(tǒng)的縱向特征值[61]Fig.8 Longitudinal eigenvalues of the flapping wing system with different wing structure[61]

Bluman、Sridhar和Kang[62]同樣研究了撲動(dòng)翼的柔性對(duì)懸停狀態(tài)下飛行器縱向穩(wěn)定性的影響。他們使用牛頓第二定律推導(dǎo)了運(yùn)動(dòng)方程。對(duì)剛性與柔性翼的氣動(dòng)計(jì)算，他們分別采用了準(zhǔn)定常方法以及二維N-S方程耦合歐拉-伯努利梁的程序?；谏鲜龇椒ㄋ麄?cè)u(píng)估了系統(tǒng)的縱向穩(wěn)定性。

為了更好地比較，剛性翼被設(shè)計(jì)為繞著前緣主動(dòng)旋轉(zhuǎn)，旋轉(zhuǎn)角度直接從柔性翼的結(jié)構(gòu)變形量得到。由于兩種布局的運(yùn)動(dòng)基本一致，剛性翼和柔性翼的運(yùn)動(dòng)差別主要來源于以下幾個(gè)方面：1)柔性翼的翼型彎度；2)流場(chǎng)不同的渦結(jié)構(gòu)；3)機(jī)體運(yùn)動(dòng)誤差的累積放大。

計(jì)算結(jié)果表明，具有剛性和柔性撲動(dòng)翼系統(tǒng)的縱向響應(yīng)有很大差別，如圖9所示。由圖可知，在機(jī)體x方向，柔性翼對(duì)應(yīng)的機(jī)體以相對(duì)較小的速度向著機(jī)尾方向運(yùn)動(dòng)，而剛性翼對(duì)應(yīng)的機(jī)體則以較大的速度向著機(jī)頭方向運(yùn)動(dòng)。在機(jī)體y方向上，運(yùn)動(dòng)方向雖然一致，但由于柔性翼比剛性翼產(chǎn)生更大的升力，所以柔性翼朝著y正方向的運(yùn)動(dòng)更為明顯。在俯仰方向，柔性翼的機(jī)體俯仰速度明顯小于剛性翼機(jī)體。

圖9 剛性翼與柔性翼的縱向穩(wěn)定性(u+、w+、q+分別代表無量綱的機(jī)體坐標(biāo)系下的x方向速度、y方向速度以及俯仰速率)[62]Fig.9 Longitudinal response of the FWMAV with rigid wing and flexible wing in a near-hover condition[62]

薛棟等[63]對(duì)撲動(dòng)翼柔性對(duì)撲翼飛行器的縱向穩(wěn)定性的影響也作了研究。通過對(duì)分別配置剛性撲動(dòng)翼和柔性撲動(dòng)翼的撲翼飛行器的縱向穩(wěn)定性進(jìn)行了比較，據(jù)此研究了撲動(dòng)翼的柔性對(duì)撲翼飛行器前飛過程中縱向穩(wěn)定性的影響。

具體方法是首先通過CFD對(duì)剛性撲動(dòng)翼的氣動(dòng)導(dǎo)數(shù)進(jìn)行計(jì)算，同時(shí)使用CFD/CSD耦合的方法對(duì)柔性撲動(dòng)翼的氣動(dòng)導(dǎo)數(shù)進(jìn)行計(jì)算，通過分析縱向運(yùn)動(dòng)小擾動(dòng)方程的特征根和特征向量，比較了兩種不同撲動(dòng)翼對(duì)應(yīng)系統(tǒng)的縱向穩(wěn)定性。結(jié)果為柔性撲動(dòng)翼對(duì)應(yīng)的系統(tǒng)為穩(wěn)定的，剛性翼對(duì)應(yīng)的系統(tǒng)則為不穩(wěn)定的，說明撲動(dòng)翼的柔性改變了系統(tǒng)的穩(wěn)定性。此外，研究結(jié)果還表明，撲翼飛行器的重心是影響縱向穩(wěn)定性的一個(gè)關(guān)鍵因素，通過設(shè)計(jì)重心的位置可以改變撲翼飛行器的縱向穩(wěn)定性。

以上三個(gè)研究分別通過不同方法對(duì)不同的研究對(duì)象進(jìn)行了撲動(dòng)翼柔性對(duì)飛行穩(wěn)定性的影響研究。研究結(jié)果表明，撲動(dòng)翼的柔性會(huì)影響系統(tǒng)的穩(wěn)定性，但具體是增穩(wěn)還是破壞穩(wěn)定，需要具體問題具體分析，不能一概而論。

5 總結(jié)與展望

雖然仿鳥型撲翼飛行器在飛行機(jī)動(dòng)性和飛行效率等方面具有巨大發(fā)展?jié)摿?，但目前仍面臨著很多尚未解決的難題。本文通過對(duì)仿鳥型撲翼飛行器的氣動(dòng)機(jī)理研究、氣動(dòng)結(jié)構(gòu)耦合研究、動(dòng)力學(xué)建模、飛行穩(wěn)定性分析方法以及撲動(dòng)翼柔性對(duì)飛行穩(wěn)定性的影響研究進(jìn)行了回顧和分析，建議今后針對(duì)仿鳥型撲翼飛行器的研究需關(guān)注以下幾個(gè)方面：

1) 目前缺乏對(duì)仿鳥型撲翼飛行器非定常氣動(dòng)機(jī)理的研究，尤其在雷諾數(shù)大于10000時(shí)前緣渦的穩(wěn)定條件和穩(wěn)定方法，這一方面值得研究人員繼續(xù)關(guān)注并作深入研究；

2) 由于仿鳥型撲翼飛行器的撲動(dòng)翼的氣動(dòng)和結(jié)構(gòu)之間存在強(qiáng)耦合，而非定常流場(chǎng)求解和非線性結(jié)構(gòu)動(dòng)力學(xué)計(jì)算問題都需要較大的計(jì)算代價(jià)，嚴(yán)重制約著撲動(dòng)翼氣動(dòng)性能評(píng)估和結(jié)構(gòu)設(shè)計(jì)速度。尺度律雖然是一種解決思路，但目前針對(duì)各向異性材料構(gòu)成的復(fù)雜撲動(dòng)翼結(jié)構(gòu)的尺度律分析理論仍為空白；

3) 仿鳥型撲翼飛行器的低撲動(dòng)頻率和高撲動(dòng)翼質(zhì)量占比，使得其柔性撲動(dòng)翼的撲動(dòng)運(yùn)動(dòng)與機(jī)體運(yùn)動(dòng)有較強(qiáng)的耦合，因此在分析氣動(dòng)特性和飛行穩(wěn)定性時(shí)需要將氣動(dòng)、結(jié)構(gòu)和飛行動(dòng)力學(xué)進(jìn)行耦合。由于每個(gè)問題的求解都需大量的時(shí)間成本和計(jì)算資源，三個(gè)問題的耦合計(jì)算會(huì)使計(jì)算量成倍增長(zhǎng)，如何利用較為精確的計(jì)算結(jié)果和實(shí)驗(yàn)結(jié)果構(gòu)造出滿足實(shí)用要求的高效和簡(jiǎn)化的多學(xué)科耦合模型是目前迫切需要解決的問題。

[1]Shyy W, Aono H, Kang C, et al. An introduction to flapping wing aerodynamics[M]. New York, NY: Cambridge University Press, 2013

[2]Wood R J. The first takeoff of a biologically inspired at-scale robotic insect[J]. IEEE Transactions on Robotics, 2013, 24(2): 341-347

[3]Send W, Fischer M, Jebens K, et al. Artificial hinged-wing bird with active torsion and partially linear kinematics[C]//ICAS2012. Brisbane, Australia, 2012: 1148-1157

[4]鄭祥明. 微型飛行器非線性飛行動(dòng)力學(xué)與智能控制研究[D]. 南京: 南京航空航天大學(xué), 2008

[5]Yang W Q, Song B F, Song W P, et al. The effects of span-wise and chord-wise flexibility on the aerodynamic performance of micro flapping-wing[J]. Chinese Science Bulletin, 2012, 57(22): 2887-2897

[6]Shyy W, Kang C K, Chirarattananon P, et al. Aerodynamics, sensing and control of insect-scale flapping-wing flight[J]. Proceedings Mathematical Physical & Engineering Sciences, 2016, 472(2186): 20150712

[7]Greenewalt C H. Dimensional relationships for flying animals[J]. Smithson Misc Collect, 1962, 144: 1-46

[8]Taylor G K, Nudds R L, Thomas A L R. Flying and swimming animals cruise at a Strouhal number tuned for high power efficiency[J]. Nature, 2003, 425: 707-711

[9]Lentink D, Dickinson M H. Rotational accelerations stabilize leading edge vortices on revolving fly wings[J]. J Exp Biol, 2009, 212: 2705-2719

[10]Knoller R. Die gesetze des luftwider-standes[J]. Flug-und Motortechnik (Wien), 1909, 3(21): 1-7

[11]Albert B. Ein beitrag zur erklaerung des segel-fluges[J]. Zeitschrift für Flugte-chnik und Motorluftschiffahrt, 1912, 3: 269-272

[12]Von Karman T, Burgers J M. Aerodynamic theory[M]. Berlin: Springer, 1934

[13]Bratt J B. Flow patterns in the wake of an oscillating airfoil[R]. Aeronautical Research Council R&M, 1953-2773

[14]Spedding G R. On the significance of unsteady effects in the aerodynamic performance of flying animals[J]. Contemp Math, 1993, 141: 401-419

[15]Shkarayev S, Silin D. Applications of actuator disk theory to membrane flapping wings[J]. AIAA J, 2010, 48(10): 2227-2234

[16]Peters D A, Karunamoorthy S, Cao W. Finite-state induced flow models, Part I: Two-dimensional thin airfoil[J]. Journal of Aircraft, 1995, 32: 313-322

[17]Liu H, Kawachi K. A numerical study of insect flight[J]. Journal of Computational Physics, 1998, 146(1): 124-156

[18]Shyy W, Lian Y, Tang J, et al. Aerodynamics of low reynolds number flyers[M]. Cambridge: Cambridge University Press, 2008

[19]Zeng R S, Ang H S. Aerodynamic computation of flapping wing simulating bird wings[J]. Journal of Nanjing University of Aeronautics & Astronautics, 2003, 35(1): 6-12. (in Chinese) 曾銳, 昂海松. 仿鳥復(fù)合振動(dòng)的撲翼氣動(dòng)分析[J]. 南京航空航天大學(xué)報(bào), 2003, 35(1): 6-12

[20]Gong K. Numerical simulation of finite-span flapping wing byeuler equations[D]. Xi’an: Northwestern Polytechnical University, 2003. (in Chinese)龔凱. 有限翼展撲動(dòng)翼的歐拉方程數(shù)值模擬[D]. 西安: 西北工業(yè)大學(xué), 2003

[21]Yang S L, Song W P, Song B F, et al. Achieving reliability and validity of predicted aerodynamic performance of flapping wings for micro air vehicle(MAV)[J]. Journal of Northwestern Polytechnical University, 2006, 24(6): 768-773. (in Chinese)) 楊淑利, 宋文萍, 宋筆鋒, 等. 微型撲翼飛行器機(jī)翼氣動(dòng)特性研究[J]. 西北工業(yè)大學(xué)學(xué)報(bào), 2006, 24(6): 768-773

[22]Xie H, Song W P, Song B F. Numerical solution of Navier-Stokes equations for flow over a flapping wing[J]. Journal of Northwestern Polytechnical University, 2008, 26(1): 104-109. (in Chinese)謝輝, 宋文萍, 宋筆鋒. 微型撲翼繞流的N-S方程數(shù)值模擬[J]. 西北工業(yè)大學(xué)學(xué)報(bào), 2008, 26(1): 104-109

[23]Yang W Q, Song B F, Song W P. Distance decreasing method for confirming corresponding cells of overset grids and its application[J]. Acta Aeronautica et Astronautica Sinica, 2009, 30(2): 205-212. (in Chinese)楊文青, 宋筆鋒, 宋文萍. 高效確定重疊網(wǎng)格對(duì)應(yīng)關(guān)系的距離減縮法及其應(yīng)用[J]. 航空學(xué)報(bào), 2009, 30(2): 205-212

[24]Ang H S, Zeng R, Duan W B, et al. Aerodynamic experimental investigation for mechanism of lift and thrust of flexible flapping-wing MAV[J]. Journal of Aerospace Power, 2007, 22(11): 1838-1845. (in Chinese) 昂海松, 曾銳, 段文博, 等. 柔性撲翼微型飛行器升力和推力機(jī)理的風(fēng)洞實(shí)驗(yàn)和飛行實(shí)驗(yàn)[J]. 航空動(dòng)力學(xué)報(bào), 2007, 22(11): 1838-1845

[25]Shao L M, Song B F, Xiong C, et al. Experimental investigation of flapping-wing MAV in wind tunnel[J]. Acta Aeronautica et Astronautica Sinica, 2007, 28(2): 275-280. (in Chinese)邵立民, 宋筆鋒, 熊超, 等. 微型撲翼飛行器風(fēng)洞試驗(yàn)初步研究[J]. 航空學(xué)報(bào), 2007, 28(2): 275-280

[26]Wang L. Design and optimization for the propulsion system of the flapping-wing micro air vehicle[D]. Xi’an: Northwestern Polytechnical University, 2008. (in Chinese)王利光. 微型撲翼飛行器動(dòng)力系統(tǒng)設(shè)計(jì)與優(yōu)化[D]. 西安: 西北工業(yè)大學(xué), 2008

[27]Wang L, Song B, Yang W, et al. Experimental characterization of a flexible membrane flapping-wing in forward flight[C]//ICAS2014. St. Petersburg. Russia, 2014

[28]Chin D D, Lentink D. Flapping wing aerodynamics: from insects to vertebrates[J]. Journal of Experimental Biology, 2016, 219(7): 920-932

[29]Muijres F T, Johansson L C, Hedenstrom A. Leading edge vortex in a slow-flying passerine[J]. Biol Lett, 2012, 8 (4): 554-557

[30]Hubel T Y, Riskin D K, Swartz S M, et al. Wake structure and wing kinematics: The flight of the lesser dog-faced fruit bat, cynopterus brachyotis[J]. J Exp Biol, 2010, 213: 34287-33440

[31]Kang C, Aono H, Cesnik C E S, et al. Effects of flexibility on the aerodynamic performance of flapping wings[J]. J Fluid Mech, 2011, 689: 32-74

[32]Shyy W, Aono H, Chimakurthi S K, et al. Recent progress in flapping wing aerodynamics and aeroelasticity[J]. Prog Aerosp Sci, 2010, 46: 284-327

[33]Heathcote S, Wang Z, Gursul I. Effect of spanwise flexibility on flapping wing propulsion[J]. Journal of Fluids and Structures, 2008, 24(2): 183-199

[34]Wu P, Ifju P, Stanford B, et al. A multidisciplinary experimental study of flapping wing aeroelasticity in thrust production[C]//AIAA Structures, Structural Dynamics, and Materials Conference, 2009: 4-7

[35]Kim D K, Lee J S, Lee J Y, et al. Anaeroelastic analysis of a flexible flapping wing using modified strip theory[C]//International Society for Optical Engineering(SPIE), 2008

[36]Unger R, Haupt M C, Horst P, et al. Structural design and aeroelastic analysis of an oscillating airfoil for flapping wing propulsion[R]. AIAA 2008-306

[37]Yang W, Song B, Wang L, et al. Dynamic fluid-structure coupling method of flexible flapping wing for MAV[J]. Journal of Aerospace Engineering, 2015, 28(6): 04015006

[38]Chen L L, Song B F, Song W P, et al. Research on aerodynamic-structural coupling of flexible flapping wing[J]. Acta Aeronautica et Astronautica Sinica, 2013, 34(12): 2668- 2681. (in Chinese ) 陳利麗, 宋筆鋒, 宋文萍, 等. 一種基于結(jié)構(gòu)動(dòng)力學(xué)的柔性撲翼氣動(dòng)結(jié)構(gòu)耦合方法研究[J]. 航空學(xué)報(bào), 2013, 34(12): 2668-2681

[39]Chen L L, Song B F, Song W P, et al. Numerical aerodynamic-structural coupling research for flexible flapping wing[J]. Acta Aerodynamica Sinica, 2015, 33(1): 125-133. (in Chinese ) 陳利麗, 宋筆鋒, 宋文萍, 等. 柔性撲翼氣動(dòng)結(jié)構(gòu)耦合特性數(shù)值研究[J]. 空氣動(dòng)力學(xué)學(xué)報(bào), 2015, 33(1): 125-133

[40]Kang C, Shyy W. Scaling law and enhancement of lift generation of an insect-size hovering flexible wing[J]. J R Soc Interface, 2013, 10: 20130361

[41]Combes S A, Daniel T L. Flexural stiffness in insect wings I: Scaling and the influence of wing venation[J]. J Exp Biol, 2003, 206: 2979-2987

[42]Combes S A, Daniel T L. Flexural stiffness in insect wings II: Spatial distribution and dynamic wing bending[J]. J Exp Biol, 2003, 206: 2989-2997

[43]Gebert G, Gallmeier P, Evers J. Equations of motion for flapping flight[C]//AIAA Atmospheric Flight Mechanics Conference. Monterey, CA, 2002. AIAA 2002-4872

[44]Sun M, Wang J, Xiong Y. Dynamic flight stability of hovering insects[J]. Acta Mech Sin, 2007, 23(3): 231-246

[45]Orlowski C T, Girard A R. Modeling and simulation of the nonlinear dynamics of flapping wing MAVs[J]. AIAA J, 2011, 49(5): 969-981

[46]Caetano J V, Weehuizen M B, Visser C C, et al. Rigid-body kinematics versus flapping kinematics of a flapping wing micro air vehicle[J]. Journal of Guidance, Control, and Dynamics, 2015, 38(12): 2257-2269

[47]Sun M, Xiong Y. Dynamic flight stability of a hovering bumblebee[J]. J Exp Biol, 2005, 208(3): 447-459

[48]Taylor G K, Bomphrey R J, Hoen J. Insect flight dynamics and control[C]//AIAA Aerospace Sciences Meeting and Exhibit. Reno, Nevada, 2006

[49]Taylor G K, Thomas A L R. Animal flight dynamics II: Longitudinal stability in flapping flight[J]. J Theor Biol, 2002, 214: 351-370

[50]Taylor G K, Thomas A L R. Dynamic flight stability in the desert locust[J]. J Theor Biol, 2003, 206: 2803-2829

[51]Taylor G K, Zbikowski R. Nonlinear time periodic models of the longitudinal flight dynamics of desert locusts[J]. J R Soc Interface, 2005, 2: 197-221

[52]Thomas A L R, Taylor G K. Animal flight dynamics I: Stability in gliding flight[J]. J Theor Biol, 2001, 212: 399-424

[53]Xiong Y, Sun M. Dynamic flight stability of a bumblebee in forward flight[J]. Acta Mechanica Sinica, 2008, 24(1): 25-36

[54]Sun M, Liu Y P, Wang J K. Dynamic flight stability of hovering insects: Theoretical analysis and numerical simulation[J]. Acta Aerodynamica Sinica, 2008, 26(z1): 6-13. (in Chinese ) 孫茂, 劉彥鵬, 王濟(jì)康. 昆蟲懸停飛行的動(dòng)穩(wěn)定性: 理論分析與數(shù)值模擬[J]. 空氣動(dòng)力學(xué)學(xué)報(bào), 2008, 26(z1): 6-13

[55]Dietl J M, Garcia E. Stability in ornitho-pter longitudinal flight dynamics[J]. J Guid Control Dyn, 2008, 31(4): 1157-1162

[56]Su W, Cesnik C E S. Flight dynamic stability of a flapping wing MAV in hover[C]//52nd AIAA/ASME/ASCE/AHS/ASC Struc-tures, Structural Dynamics, and Materials Conference. Denver, Colorado, 2011

[57]Taha H E, Tahmasian S, Woolsey C A, et al. The need for higher-order averaging in the stability analysis of hovering, flapping-wing flight[J]. Bioinspiration & Biomimetics, 2015, 10(1): 016002

[58]Wu J H, Sun M. Floquet stability analysis of the longitudinal dynamics of two hovering model insects[J]. Journal of the Royal Society Interface, 2012, 9(74): 2033-2046

[59]Faruque I, Sean H J. Dipteran insect flight dynamics. Part 1: Longitudinal motion about hover[J]. Journal of Theoretical Biology, 2010, 264(2): 538-552

[60]Wu J H, Zhang Y L, Sun M. Hovering of model insects: Simulation by coupling equations of motion with Navier-Stokes equations[J]. Journal of Experimental Biology, 2009, 212(Pt 20): 3313-3329

[61]Richter M, Patil M. Influence of Wing Flexibility on the Stability of Flapping Flight[C]//AIAA Atmospheric Flight Mechanics Conference, 2013: 2-5

[62]Bluman J E, Sridhar M K, Kang C K. The influence of wing flexibility on the stability of a biomimetic flapping wing Micro air vehicle in hover[C]//AIAA/ASCE/AHS/ASC Structures, Structural Dynamics, and Materials Conference. San Diego, California, USA, 2015

[63]Xue D, Song B, Song W, et al. Effect of wing flexibility on flight dynamics stability of flapping wing MAVs in forward flight[J]. International Journal of Micro Air Vehicles, 2016, 8(3): 170-180.