劉文君 蔡 毓

（廣西大學(xué)計(jì)算機(jī)與電子信息學(xué)院，廣西 南寧 530004）

自然界中的鳥類能夠在復(fù)雜和高度動(dòng)態(tài)的環(huán)境中飛行，具有穩(wěn)定性、機(jī)動(dòng)性和高效性，對(duì)環(huán)境干擾有很強(qiáng)的抵抗力，從而可以進(jìn)行長距離遷徙。然而，由生物飛行靈感啟發(fā)研制的各種規(guī)模結(jié)構(gòu)、傳感、穩(wěn)定性控制的仿生撲翼飛行器還沒有自主地表現(xiàn)出這些特性[1]。仿鳥型撲翼飛行器和仿昆蟲型撲翼飛行器的飛行機(jī)理有較大差別。仿昆蟲型撲翼飛行器主要靠撲翼撲動(dòng)產(chǎn)生升力和推力，由于其撲動(dòng)頻率高，且機(jī)翼質(zhì)量占比小，所以氣動(dòng)、結(jié)構(gòu)和飛行力學(xué)之間的耦合問題影響較小。對(duì)于仿鳥型撲翼飛行器而言，其具有低撲動(dòng)頻率、高機(jī)翼質(zhì)量占比和耦合關(guān)系性強(qiáng)等特點(diǎn)。

氣動(dòng)效率是仿生撲翼飛行器研究中備受關(guān)注的領(lǐng)域，是該飛行器在飛行狀態(tài)下受到的升力、阻力、推力、力的方向、力的大小以及穩(wěn)定性等影響的客觀反應(yīng)，也是檢驗(yàn)飛行器外形布局及氣動(dòng)關(guān)系設(shè)計(jì)是否合理的重要依據(jù)[2]。目前，關(guān)于撲翼的氣動(dòng)效率，在生物觀測(cè)、理論分析、數(shù)值模擬、飛行試驗(yàn)4 個(gè)層面都在逐步開展研究。在數(shù)值模擬方面，Young J[3]采用可壓縮二維Navier-Stokes 求解器，對(duì)斜向振蕩的NACA0012 翼型的流場(chǎng) 進(jìn)行了數(shù)值模擬，用數(shù)值粒子追蹤法對(duì)翼型的尾跡進(jìn)行了可視化。Lai JCS[4]研究了Strouhal 數(shù)推進(jìn)效率的峰值，得出僅用Strouhal 數(shù)不足以表征撲翼推進(jìn)效率的結(jié)論。Trizila P[5]使用代理建模技術(shù)對(duì)前緣渦、尾流及射流進(jìn)行了數(shù)值研究。肖天航等人[6]提出一種Delaunary 圖映射網(wǎng)格變形技術(shù)和非結(jié)構(gòu)嵌套網(wǎng)格的數(shù)值計(jì)算方法。

雖然針對(duì)撲翼運(yùn)動(dòng)的研究已經(jīng)初步揭示了產(chǎn)生升推力的飛行機(jī)制，然而定性的認(rèn)識(shí)還不能完全滿足撲翼飛行器設(shè)計(jì)層面的需求。一個(gè)重要的俯仰運(yùn)動(dòng)學(xué)案例[7]在中沖程、中下沖程不同的迎角下，分析二維、三維CFD 模型流場(chǎng)非定常瞬態(tài)狀態(tài)，二者經(jīng)驗(yàn)系數(shù)的變化曲線基本一致。該案例表明二維模型為微型飛行器設(shè)計(jì)中的實(shí)際應(yīng)用提供了重要的仿真依據(jù)。該文在進(jìn)行二維數(shù)值模擬的同時(shí)，對(duì)結(jié)果進(jìn)行了積分計(jì)算，得到整個(gè)機(jī)翼平面的氣動(dòng)力和氣動(dòng)效率。

1 實(shí)驗(yàn)?zāi)Ｐ?/h2>

1.1 撲動(dòng)機(jī)構(gòu)模型

在研制仿生撲翼飛行器時(shí)，可借助仿生學(xué)原理，參考鳥類的翼型構(gòu)造、形狀[8]，來研制機(jī)翼結(jié)構(gòu)，如圖1 所示。飛行時(shí)產(chǎn)生升力克服其自身重量，產(chǎn)生推力克服飛行時(shí)的阻力。

該文在實(shí)驗(yàn)中的飛行器機(jī)構(gòu)簡(jiǎn)圖，如圖2（a）所示，并將其機(jī)翼撲動(dòng)基于葉素理論將機(jī)翼平面按展向分成若干翼帶，進(jìn)行三維數(shù)學(xué)建模，如圖2（b）所示。

圖1 鳥翼的構(gòu)造和不同階段翅膀的形狀變化

1.2 機(jī)翼運(yùn)動(dòng)分析

機(jī)構(gòu)簡(jiǎn)圖的運(yùn)動(dòng)角度位移是分析機(jī)構(gòu)運(yùn)動(dòng)性能、優(yōu)化機(jī)構(gòu)的依據(jù)。機(jī)翼的運(yùn)動(dòng)分析，是根據(jù)原動(dòng)件的運(yùn)動(dòng)規(guī)律，求解從動(dòng)件的運(yùn)動(dòng)規(guī)律。機(jī)構(gòu)運(yùn)動(dòng)如圖3 所示。

圖2 撲翼飛行器模型機(jī)構(gòu)簡(jiǎn)圖及機(jī)翼建模

在機(jī)構(gòu)簡(jiǎn)圖中的支鏈上引入閉環(huán)結(jié)構(gòu)，利用閉環(huán)矢量方程建立連桿機(jī)構(gòu)的約束方程。閉環(huán)方程為

式中：R1是曲柄，R2、R3、R4為機(jī)械連桿，R5為機(jī)構(gòu)中的輔助桿。R1繞固定鉸鏈點(diǎn)O 轉(zhuǎn)動(dòng)，即A 點(diǎn)是以O(shè) 為原點(diǎn)、R1為半徑的圓上運(yùn)動(dòng)的動(dòng)點(diǎn)。α1是連桿機(jī)構(gòu)的設(shè)定角度，α2是R1繞定點(diǎn)O 轉(zhuǎn)過的瞬時(shí)角，α3是翼?xiàng)UAB 繞動(dòng)的瞬時(shí)撲動(dòng)角，α4是平板翼型中的機(jī)翼機(jī)構(gòu)的位置角。沉浮運(yùn)動(dòng)參數(shù)示意圖如圖4 所示。沖程平面角為0。全局坐標(biāo)系（X，Y，Z）固定在沖程平面的中心，Z 方向垂直于沖程平面，Y 方向垂直于機(jī)翼軸，X 方向與沖程平面平行。撲翼運(yùn)動(dòng)規(guī)律接近于簡(jiǎn)諧運(yùn)動(dòng)，撲動(dòng)飛行可以同時(shí)產(chǎn)生升力和推力。平板翼型沉浮運(yùn)動(dòng)位移被定義為

式中：hm是拍打振幅，f 是拍打頻率，t 是拍打時(shí)間，h（t）是t 時(shí)間內(nèi)的位移。

2 數(shù)值計(jì)算方法

該文的數(shù)值模擬基于ANSYS 軟件平臺(tái)進(jìn)行。該平臺(tái)集成了N-S 方程流場(chǎng)解算器、動(dòng)態(tài)網(wǎng)格生成技術(shù)、空氣動(dòng)力學(xué)算法等多方面模塊，可用于理想和多種設(shè)置氣體/流體介質(zhì)的定常/非定常數(shù)值模擬仿真計(jì)算。

圖3 撲翼機(jī)構(gòu)運(yùn)動(dòng)簡(jiǎn)圖

流場(chǎng)計(jì)算的基本過程是在空間上用有限體積法，將計(jì)算域離散成許多小的體積單元，在每個(gè)體積單元上對(duì)離散后的控制方程組進(jìn)行求解。這也是流場(chǎng)計(jì)算方法的本質(zhì)。采用分離式求解器（Segregated solver），逐一、順序地求解各變量。

假定初始?jí)毫退俣龋_定離散方程的系數(shù)及常數(shù)項(xiàng)；求解連續(xù)方程、動(dòng)量方程、能量方程；求解湍流方程及其標(biāo)量方程；判斷當(dāng)前時(shí)間步上的計(jì)算是否收斂。

圖4 沉浮運(yùn)動(dòng)懸停參數(shù)示意圖

2.1 動(dòng)態(tài)網(wǎng)格技術(shù)

在非定常流場(chǎng)數(shù)值模擬的動(dòng)態(tài)邊界運(yùn)動(dòng)中，動(dòng)態(tài)網(wǎng)格生成是非常重要的組成部分。動(dòng)網(wǎng)格計(jì)算中，網(wǎng)格的動(dòng)態(tài)變化過程可以用3 種模型進(jìn)行計(jì)算，即彈簧近似光滑模型（spring-based smoothing）、動(dòng)態(tài)分層模型（dynamic layering）和局部重劃模型（local remeshing）。該文采用彈簧近似光滑法和局部重劃法。

彈簧近似法適用于該研究中的非結(jié)構(gòu)化網(wǎng)格，網(wǎng)格拓?fù)洳蛔?，通過近似彈簧的壓縮或拉伸，實(shí)現(xiàn)網(wǎng)格和模擬流場(chǎng)計(jì)算區(qū)域的變化。網(wǎng)格邊界處位移變化后，彈簧系統(tǒng)會(huì)經(jīng)過調(diào)整，使已破壞的原有平衡達(dá)到新的平衡，生成力的大小由胡克定律計(jì)算得到。在控制臺(tái)窗口掌握區(qū)域范圍、體積統(tǒng)計(jì)以及連通性信息。網(wǎng)格檢查沒有出現(xiàn)網(wǎng)格體積為負(fù)數(shù)的狀態(tài)。該次數(shù)值模擬實(shí)驗(yàn)需對(duì)模型進(jìn)行不同方向的簡(jiǎn)諧運(yùn)動(dòng)，多次模擬之后得到相應(yīng)良好的設(shè)定尺寸，根據(jù)預(yù)先設(shè)置的最大、最小網(wǎng)格規(guī)格，同時(shí)采用網(wǎng)格重劃法在網(wǎng)格更新計(jì)算區(qū)域進(jìn)行網(wǎng)格重劃。在數(shù)值模擬中，還要注意更新區(qū)域的完整性。

2.2 非定常數(shù)值模擬方法

數(shù)值模擬中的控制方程是具有恒定密度和黏度的非定常Navier-Stokes 方程[9]。基于動(dòng)態(tài)網(wǎng)格，采用二階迎風(fēng)格式（second order upwind）有限體積法，考慮量值在物理上的值點(diǎn)分布，即曲線的曲率影響，具有二階精度截差，同時(shí)也會(huì)考慮氣流的流動(dòng)方向，所以數(shù)值解在物理上是合理的。利用BLU-SGS 隱式算法進(jìn)行時(shí)間推進(jìn)，求解低速非定常雷諾平均N-S 方程。同一算例應(yīng)用不同的湍流模型，得到的數(shù)值結(jié)果會(huì)存在很大不同，與此同時(shí)，采用湍流模型中試驗(yàn)最成功的一種模型，Spalart-Allmaras 湍流模型[10]，可以合理地處理邊界層的黏性影響區(qū)域。

文中分別進(jìn)行了穩(wěn)態(tài)和瞬時(shí)模擬計(jì)算，對(duì)于穩(wěn)態(tài)問題，壓力修正法采用SIMPLEC 算法，對(duì)于瞬態(tài)問題，同時(shí)為了更好地滿足動(dòng)量方程、連續(xù)性方程，采用PISO 算法。

2.3 算例準(zhǔn)確性驗(yàn)證

該算例重點(diǎn)討論外部空氣動(dòng)力學(xué)應(yīng)用的無量綱系數(shù)，對(duì)氣動(dòng)進(jìn)行配置，計(jì)算典型的阻力系數(shù)?；谠撐难芯康牡退賳栴}，周圍流場(chǎng)中的空氣密度保持不變，為不可壓縮流，并且不會(huì)出現(xiàn)馬赫數(shù)。氣動(dòng)系數(shù)是關(guān)于雷諾數(shù)的函數(shù)CD=f（M，Re）。如圖5 所示，平板模型垂直于流動(dòng)氣流，這種機(jī)構(gòu)在任何常規(guī)配置中都會(huì)產(chǎn)生最大的阻力系數(shù)。

雷諾數(shù)計(jì)算公式為：

阻力系數(shù)計(jì)算公式為：

其中，S 是平板單位跨度的正面面積。文獻(xiàn)[7]中的阻力系數(shù)為2.0，該文數(shù)值模擬方法下，迭代收斂后的阻力系數(shù)為1.9986。數(shù)值模擬和理論方法對(duì)該算例進(jìn)行比較，數(shù)據(jù)吻合，所以將利用該數(shù)值模擬方法進(jìn)行進(jìn)一步研究。

2.4 數(shù)據(jù)處理方法

對(duì)不同參數(shù)的二維模型進(jìn)行非定常流場(chǎng)數(shù)值模擬，得到已設(shè)置的氣流條件下的升力系數(shù)和阻力系數(shù)，與此同時(shí)，將每個(gè)翼帶的二維氣動(dòng)力，通過積分計(jì)算形式，得到整個(gè)機(jī)翼平面氣動(dòng)力。后續(xù)通過如下數(shù)據(jù)處理方法做進(jìn)一步研究。

減縮頻率定義為

雷諾數(shù)表示遷移慣性力與黏性力的比值，定義為

經(jīng)過無量綱化，拍打運(yùn)動(dòng)的升力系數(shù)和推力系數(shù)定義為

式中：CL、CT分別是升力系數(shù)和推力系數(shù)定義，L、T 分別是升力和推力，是來流密度，S 是平板跨度的正面面積。

依據(jù)文獻(xiàn)[11]，平均推力系數(shù)和平均氣動(dòng)力輸入系數(shù)分別定義為

同時(shí)平均氣動(dòng)力系數(shù)定義為

依據(jù)文獻(xiàn)[12]，升力效率和推進(jìn)效率分別定義為

3 數(shù)值結(jié)果與分析

對(duì)于撲翼飛行，機(jī)翼產(chǎn)生顯著升力來抵消自身重力。氣動(dòng)效率是與給定功率的升力和推力有關(guān)的運(yùn)動(dòng)效率，能保持在空中的能源供應(yīng)。該文從上述2 個(gè)方面對(duì)撲翼氣動(dòng)效率展開研究。參數(shù)空間包括撲打頻率f，風(fēng)速，減縮頻率k，雷諾數(shù)Re。在二維撲翼無粘平板經(jīng)典理論中，對(duì)于低振幅、低減縮頻率的情況，撲翼在周期內(nèi)的平均推力不依賴初始攻角[13]，并且對(duì)于三維流動(dòng)的撲翼，只要運(yùn)動(dòng)不發(fā)生明顯的流動(dòng)分離現(xiàn)象，同樣適用該理論[13]。該文數(shù)值研究中平均攻角α 設(shè)置為零。

3.1 撲動(dòng)頻率f 和風(fēng)速對(duì)氣動(dòng)效率的影響

低頻振蕩是導(dǎo)致流動(dòng)失穩(wěn)的一個(gè)主要因素[13]，由飛行速度和撲打頻率共同作用。對(duì)于撲翼飛行，機(jī)翼產(chǎn)生顯著升力來抵消自身重力。在不同撲打頻率f 和風(fēng)速下，總氣動(dòng)力會(huì)隨撲打頻率的提高而增大。在飛行器起飛階段，提高撲打頻率，增加推力，升力會(huì)在達(dá)到一定撲打頻率后有增加的趨勢(shì)。結(jié)合如圖5 所示的氣動(dòng)效率實(shí)驗(yàn)結(jié)果，在撲打頻率為7 Hz、8 Hz 時(shí)，具有較高速度，達(dá)到最大升力效率。進(jìn)入巡航階段，降低撲打頻率，達(dá)到最大推進(jìn)效率。

3.2 減縮頻率k 和雷諾數(shù)Re對(duì)氣動(dòng)效率的影響

該節(jié)實(shí)驗(yàn)計(jì)算考慮有限翼展的拍打運(yùn)動(dòng)，針對(duì)k=0.2～1.5、Re=1.7×104～5.2×104的減縮頻率、低雷諾參數(shù)進(jìn)行大量數(shù)值計(jì)算，得到定性以及定量的結(jié)果。

隨著減縮頻率k 的增加，實(shí)驗(yàn)中雷諾數(shù)下的氣動(dòng)力均增加。升力效率和推進(jìn)效率分別如圖6 所示。計(jì)算出的氣動(dòng)效率表明，隨著雷諾數(shù)Re 的增加，實(shí)驗(yàn)中減縮頻率k 下的推進(jìn)效率逐步下降至接近平緩，而升力效率的結(jié)果變化更陡峭，出現(xiàn)較高的峰值。

從續(xù)航時(shí)間考慮，應(yīng)以較低的撲動(dòng)頻率飛行，較高的撲動(dòng)頻率對(duì)撲翼飛行沒有本質(zhì)優(yōu)化，卻會(huì)造成機(jī)翼架構(gòu)的損壞，正負(fù)升力抵消。

圖5 氣動(dòng)效率隨撲打頻率和風(fēng)速變化的比較

圖6 氣動(dòng)效率隨減縮頻率和雷諾數(shù)變化的比較

4 結(jié)論

為了探討仿生撲翼飛行器產(chǎn)生升力和推力的飛行機(jī)理，該文在研究考慮撲翼氣動(dòng)力和氣動(dòng)效率數(shù)值計(jì)算方法的基礎(chǔ)上，利用數(shù)值模擬仿真，針對(duì)不同撲動(dòng)頻率、風(fēng)速、減縮頻率、雷諾數(shù)、撲動(dòng)幅度等參數(shù)變量，研究揭示了該飛行器撲翼機(jī)構(gòu)數(shù)值模擬計(jì)算的總氣動(dòng)力和氣動(dòng)效率的變化，從而研究分析主要飛行規(guī)律。為仿生撲翼飛行器的進(jìn)一步研究提供了參考依據(jù)。