劉 聰，諶哲宇，儲(chǔ)天羽，高晟智，董宇煒

(中國民航大學(xué) 空中交通管理學(xué)院，天津 300300)

0 引言

小型旋翼無人機(jī)具有低成本性、高適應(yīng)性、可懸停及垂直起降的高機(jī)動(dòng)性等優(yōu)點(diǎn)，在物流、植保、安防等行業(yè)得到大量應(yīng)用[1-2]，其運(yùn)行安全性也逐漸引起重視[3-5]。但單機(jī)執(zhí)行任務(wù)易受限，隨著無人飛行技術(shù)的發(fā)展和任務(wù)的復(fù)雜化，多無人機(jī)編隊(duì)飛行具有廣闊發(fā)展前景[6-7]。與有人機(jī)類似，編隊(duì)飛行提高任務(wù)執(zhí)行效率的同時(shí)，多機(jī)之間的氣動(dòng)干擾效應(yīng)不容忽視。當(dāng)后機(jī)處于前機(jī)的渦流場(chǎng)中，后機(jī)氣動(dòng)力變化，嚴(yán)重影響后機(jī)的飛行穩(wěn)定性，若無法操控則有墜機(jī)的風(fēng)險(xiǎn)[8-9]。旋翼為無人機(jī)主要?jiǎng)恿Σ考?，是產(chǎn)生氣動(dòng)力和渦流場(chǎng)的主要來源。因此，有必要針對(duì)不同間隔的布局形式，分析無人機(jī)旋翼間的流場(chǎng)干擾特征，并定量給出間隔與氣動(dòng)力參數(shù)的對(duì)應(yīng)關(guān)系，這有助于近距密集編隊(duì)飛行時(shí)旋翼無人機(jī)間的安全間隔設(shè)置并制定安全的飛行控制策略。

在無人機(jī)旋翼間的氣動(dòng)干擾研究方面，郭瑞雪等[10]總結(jié)了國內(nèi)外研究現(xiàn)狀，多旋翼飛行器氣動(dòng)干擾特性受到自身翼間干擾和非定常來流影響，并且指出使用計(jì)算流體力學(xué)CFD方法模擬分析相鄰旋翼流場(chǎng)的疊加和干擾的細(xì)節(jié)特性是可行且必要的。沈奧等[11]使用CFD方法模擬了旋翼無人機(jī)單機(jī)懸停狀態(tài)和飛行狀態(tài)的流場(chǎng)，根據(jù)單機(jī)流場(chǎng)定性劃分了不同區(qū)域，指出了區(qū)域不同會(huì)對(duì)密集編隊(duì)時(shí)的氣動(dòng)產(chǎn)生影響，但未針對(duì)雙機(jī)或雙旋翼進(jìn)行模擬并給出定量的結(jié)果。縱向間隔的研究以共軸雙旋翼為代表，因其優(yōu)良的升力特性被直升機(jī)和無人機(jī)等廣泛采用，共軸雙旋翼的氣動(dòng)干擾機(jī)理[12-14]進(jìn)行較深入研究。但共軸雙旋翼的共軸結(jié)構(gòu)決定了兩旋翼僅有縱向間隔。另外，雷瑤等[15]和李沛等[16]等從無人機(jī)的設(shè)計(jì)布局角度研究了旋翼間距對(duì)旋翼飛行器氣動(dòng)性能的影響，試圖找到最優(yōu)的旋翼間距。此間距為橫向間隔，即旋翼的旋轉(zhuǎn)面在同一平面內(nèi)。以上關(guān)于雙旋翼間氣動(dòng)的干擾研究，增強(qiáng)了對(duì)其流場(chǎng)作用機(jī)理的了解和認(rèn)知；但是，雙旋翼的位置關(guān)系大多局限在僅有縱向或橫向的單一間隔，同時(shí)存在橫向和縱向間隔和氣動(dòng)干擾研究則很少涉及，而這在無人機(jī)編隊(duì)飛行時(shí)為常見形態(tài)。

因此本文擬以雙旋翼單元為研究對(duì)象，主要采用計(jì)算流體力學(xué)CFD仿真的方法，分析懸停狀態(tài)時(shí)雙旋翼處在不同橫向和縱向間隔位置時(shí)的流場(chǎng)特征和氣動(dòng)參數(shù)變化情況，為無人機(jī)編隊(duì)時(shí)的安全間隔布局提供理論參考。

1 計(jì)算模型和方法

1.1 幾何模型

無人機(jī)旋翼為其主要?jiǎng)恿Σ考?，旋翼之間的氣動(dòng)干擾直接影響整機(jī)的穩(wěn)定性，為簡化研究，暫不研究機(jī)身與旋翼之間的干擾。本文研究采用的無人機(jī)旋翼模型如圖1所示，槳半徑11.94 cm，槳距12.70 cm，此類槳廣泛應(yīng)用于大疆精靈系列無人機(jī)。沿半徑方向其弦長和扭轉(zhuǎn)角度變化較大，分布如圖2所示。在29%相對(duì)半徑位置附近扭轉(zhuǎn)角和弦長取得最大值，最大扭轉(zhuǎn)角和平均扭轉(zhuǎn)角分別為19.04°和12.75°，最大弦長和平均弦長分別為6.06 cm和2.00 cm。

圖1 旋翼模型Fig.1 Rotor model

圖2 沿半徑方向弦長和扭轉(zhuǎn)角分布Fig.2 Distribution of chord length and twist angle along radial direction

1.2 計(jì)算域和網(wǎng)格

對(duì)旋翼流場(chǎng)的仿真采用多重參考系方法(Multiple reference frame，MRF)，經(jīng)與實(shí)驗(yàn)對(duì)比該方法使用于旋翼的流場(chǎng)仿真[17]。此方法需要將計(jì)算域劃分為包含旋翼在內(nèi)的旋轉(zhuǎn)域和外部靜止域，如圖3所示：靜止域尺寸為5D×2D×10D，旋轉(zhuǎn)域直徑1.1D，厚度約0.13D。旋轉(zhuǎn)域網(wǎng)格保持靜止，設(shè)置其旋轉(zhuǎn)坐標(biāo)系并以作用的科氏力和向心力來定常計(jì)算；靜止域按慣性系坐標(biāo)進(jìn)行計(jì)算；旋轉(zhuǎn)域與靜止域間通過交界面?zhèn)鬟f流動(dòng)參數(shù)。當(dāng)分析2個(gè)旋翼間隔影響時(shí)，固定1#旋翼，移動(dòng)2#旋翼相對(duì)水平位置X/D和相對(duì)豎直位置Z/D。

圖3 計(jì)算域Fig.3 Computational domain

計(jì)算域流場(chǎng)采用結(jié)構(gòu)-非結(jié)構(gòu)混合網(wǎng)格進(jìn)行劃分，網(wǎng)格單元總共約980萬。其中靜止域采用結(jié)構(gòu)網(wǎng)格，約400萬網(wǎng)格單元；2個(gè)旋轉(zhuǎn)域采用非結(jié)構(gòu)網(wǎng)格，各290萬。過旋翼豎直截面的網(wǎng)格如圖4所示。為更好捕捉旋翼流場(chǎng)特征，對(duì)旋翼壁面及其前后緣位置進(jìn)行網(wǎng)格加密，如圖5。另外，旋翼表面采用棱柱網(wǎng)格劃分10層邊界層網(wǎng)格，并保證第1層網(wǎng)格對(duì)應(yīng)y+在1附近。

圖4 中截面網(wǎng)格Fig.4 Mesh of middle section

圖5 旋翼壁面網(wǎng)格Fig.5 Mesh of rotor wall

靜止域外流場(chǎng)6個(gè)邊界設(shè)置為壓力出口邊界條件，模擬靜止無窮遠(yuǎn)的環(huán)境；靜動(dòng)域的交界面設(shè)置為interface邊界；旋翼表面設(shè)為無滑移壁面邊界，并指定旋轉(zhuǎn)域的旋轉(zhuǎn)中心和轉(zhuǎn)速。

1.3 控制方程和計(jì)算方法

小尺寸旋翼根據(jù)其直徑D(0.239 m)和設(shè)計(jì)轉(zhuǎn)速(5 000 r/min)，槳尖線速度約60 m/s，馬赫數(shù)Ma<0.3，雷諾數(shù)量級(jí)105，可認(rèn)為不可壓縮流。故本文采用不可壓且有黏的三維N-S方程模擬流場(chǎng)，采用有限體積法對(duì)方程進(jìn)行空間離散，控制方程如式(1)所示：

(1)

矩陣向量U、F、G、H和J如式(2)～(6)所示：

(2)

(3)

(4)

(5)

J=[0ρfxρfyρfzρw(ufx+vfy+wfz)]Τ

(6)

式中：x,y,z分別為三維坐標(biāo)方向；t為時(shí)間;ρ為密度，kg/m3;V為合速度;u,v,w分別為x,y,z方向的分速度，m/s;p為靜壓，Pa;e為比內(nèi)能，J/kg;fx，fy，fz分別為x，y，z方向的粘性切應(yīng)力，N/m2。

采取SIMPLE算法進(jìn)行速度場(chǎng)和壓力場(chǎng)的耦合。湍流模型選用SSTk-ω模型，該模型結(jié)合了自由流中的k-ε和近壁的k-ω模型優(yōu)點(diǎn)，適用于旋翼氣動(dòng)的仿真模擬，并采用二階格式對(duì)壓力、動(dòng)量、能量和湍流項(xiàng)進(jìn)行差分求解。以計(jì)算殘差小于1×10-4作為收斂判斷標(biāo)準(zhǔn)。在2.1節(jié)中將對(duì)此計(jì)算方法進(jìn)行實(shí)驗(yàn)驗(yàn)證。

1.4 計(jì)算工況和參數(shù)定義

本文研究主要分析旋翼不同位置間隔下的干擾情況，計(jì)算工況主要改變旋翼間的橫向間隔X/D和縱向間隔Z/D，具體間隔設(shè)置如表1所示。

表1 計(jì)算工況Table 1 Calculation conditions

結(jié)果分析時(shí)使用無量綱系數(shù)評(píng)估流場(chǎng)干涉下旋翼的氣動(dòng)性能，即拉力系數(shù)CT、扭矩系數(shù)CQ、俯仰力矩系數(shù)CM。其定義式如式(7)～(9)所示：

(7)

(8)

(9)

式中：T為旋翼拉力，N；Q為旋翼扭矩，N·M；M為水平內(nèi)的合力矩，N·M；R為旋翼半徑，m；Ω為旋翼角速度，rad/s；ρ為空氣密度，取1.225 kg/m3。

2 結(jié)果與分析

2.1 單旋翼流場(chǎng)仿真結(jié)果與實(shí)驗(yàn)驗(yàn)證

下文展示單個(gè)旋翼在不同轉(zhuǎn)速下的氣動(dòng)特征，作為雙旋翼干擾結(jié)果的對(duì)比基準(zhǔn)；并與實(shí)驗(yàn)結(jié)果對(duì)比，驗(yàn)證數(shù)值模擬方法的準(zhǔn)確性。圖6～7為2 000,5 000,8 000 RPM的中截面速度和渦量仿真結(jié)果。從速度云圖上看，每側(cè)槳葉下方形成1股高速區(qū)，且向外半徑方向向外擴(kuò)散；同時(shí)在中心區(qū)域，由于高速氣流誘導(dǎo)，形成1股較低速度區(qū)。隨著轉(zhuǎn)速升高，旋翼下方氣流速度增大同時(shí)，擴(kuò)散范圍逐漸增大。從圖7渦量云圖上看兩側(cè)槳尖形成的槳尖渦最為明顯，內(nèi)部則誘導(dǎo)出較為復(fù)雜的二次渦系。結(jié)合速度和渦量的仿真結(jié)果，當(dāng)轉(zhuǎn)速為5 000 RPM時(shí)旋翼下洗流場(chǎng)主要影響區(qū)域在X=1.5D、Z=2.0D范圍內(nèi)；當(dāng)轉(zhuǎn)速為8 000 RPM時(shí)，此區(qū)域擴(kuò)展為約X=2.0D、Z=3.0D。因此進(jìn)行2個(gè)旋翼干擾分析時(shí)，設(shè)置雙旋翼的相對(duì)位置間隔覆蓋單個(gè)旋翼流場(chǎng)影響范圍。

圖6 旋翼截面速度分布Fig.6 Velocity distribution of rotor section

圖7 旋翼截面渦量分布Fig.7 Vorticity distribution of rotor section

圖8為單旋翼氣動(dòng)性能測(cè)試臺(tái)。該實(shí)驗(yàn)臺(tái)用于測(cè)試定距螺旋槳或旋翼的拉扭曲線和力效實(shí)驗(yàn)等，也可置于風(fēng)洞中測(cè)試飛行狀態(tài)的性能數(shù)據(jù)。實(shí)驗(yàn)臺(tái)主要由動(dòng)力系統(tǒng)、傳感器系統(tǒng)和數(shù)據(jù)采集系統(tǒng)組成。動(dòng)力系統(tǒng)由電池、電調(diào)和電機(jī)組成，用于驅(qū)動(dòng)旋翼旋轉(zhuǎn)，控制信號(hào)由上位機(jī)發(fā)出，通過PWM信號(hào)控制電機(jī)轉(zhuǎn)速，實(shí)驗(yàn)轉(zhuǎn)速范圍為2 000～8 000 RPM；傳感器系統(tǒng)主要測(cè)量參數(shù)為輸入電壓、電流、轉(zhuǎn)速、電機(jī)溫度和旋翼拉力、扭矩。傳感器測(cè)量信號(hào)通過16位高速數(shù)據(jù)采集系統(tǒng)采集和處理數(shù)據(jù)。實(shí)驗(yàn)臺(tái)的測(cè)量誤差主要取決于各傳感器的測(cè)量精度，主要參數(shù)相對(duì)誤差如表2所示。

圖8 單旋翼氣動(dòng)性能測(cè)試臺(tái)Fig.8 Test system for aerodynamic performance of single rotor

表2 實(shí)驗(yàn)參數(shù)相對(duì)誤差Table 2 Relative errors of experimental parameters

數(shù)值仿真和實(shí)驗(yàn)結(jié)果對(duì)比如表3所示?？梢钥闯觯?dāng)轉(zhuǎn)速較低時(shí)，拉力系數(shù)和扭矩系數(shù)的計(jì)算誤差較大，這主要是由于低轉(zhuǎn)速下拉力和扭矩的絕對(duì)值均較小，放大了實(shí)驗(yàn)和模擬的相對(duì)誤差；在中高轉(zhuǎn)速時(shí)，拉力系數(shù)和扭矩系數(shù)的相對(duì)誤差較小，且比較穩(wěn)定，分別約5%和9%。上文的研究中旋翼轉(zhuǎn)速設(shè)定在懸停設(shè)計(jì)轉(zhuǎn)速5 000 RPM。

表3 單旋翼仿真與實(shí)驗(yàn)結(jié)果對(duì)比Table 3 Comparison on simulation and experimental results of single rotor

2.2 不同間隔下雙旋翼的流場(chǎng)分析

為分析在不同間隔距離下的雙旋翼流場(chǎng)干擾情況，進(jìn)行速度場(chǎng)和渦量場(chǎng)的對(duì)比，分別如圖9和圖10所示。從圖9～10中可以看出，無論是速度場(chǎng)還是渦量場(chǎng)，1#旋翼基本不受影響，而2#旋翼所受影響非常明顯。因此，著重分析2#旋翼的流場(chǎng)變化。

圖9 不同間隔雙旋翼截面速度分布Fig.9 Velocity distribution of twin rotors sections with different separations

圖10 不同間隔雙旋翼截面渦量分布Fig.10 Vorticity distribution of twin rotors sections with different separations

當(dāng)兩旋翼之間無橫向間隔，即2#旋翼位于1#旋翼的正下方時(shí)，2#旋翼的下洗氣流受1#旋翼氣流限制，未呈現(xiàn)出擴(kuò)散趨勢(shì)，速度和渦量場(chǎng)反應(yīng)出其下洗氣流局限在旋翼下方1.0D范圍內(nèi)，且高速區(qū)和高渦量區(qū)縱向范圍縮小。此時(shí)，流場(chǎng)分布仍保持對(duì)稱形態(tài)。

當(dāng)2#旋翼在橫向具有0.5D和1.0D間隔時(shí)，2#旋翼下洗流場(chǎng)出現(xiàn)非對(duì)稱特征，這主要由2#旋翼兩側(cè)槳葉上游來流的非均勻性決定。如圖9(b)～(c)和圖10(b)～(c)所示，2#旋翼主要處于1#旋翼右半側(cè)槳的下洗流中，此部分氣流除具有垂直速度外，還具有與右半槳旋轉(zhuǎn)線速度同向的水平分速度；2#旋翼左半側(cè)槳旋轉(zhuǎn)速度方向恰與1#旋翼右半側(cè)槳相反，即對(duì)于2#旋翼左半側(cè)槳來說，相當(dāng)于增大了其來流的相對(duì)速度。同理，2#旋翼右半側(cè)槳的相對(duì)來流速度相應(yīng)減小。另外豎直間隔不同，兩側(cè)來流的豎直速度分量也不同，同樣也會(huì)加劇此非對(duì)稱性?？傮w來說，2#旋翼相對(duì)1#旋翼向右側(cè)偏置，增強(qiáng)了2#旋翼下洗左半側(cè)氣流的速度和渦量，而右半側(cè)氣流的速度和渦量則相對(duì)減弱；右半側(cè)的下洗流在橫向的擴(kuò)散受到限制，這與無橫向間隔時(shí)相同。此非對(duì)稱性的流場(chǎng)分布，還會(huì)體現(xiàn)在對(duì)無量綱氣動(dòng)參數(shù)，尤其是對(duì)俯仰力矩系數(shù)的影響上，將在下文中討論。

當(dāng)雙旋翼的橫向間隔在1.5D及以上時(shí)，旋翼流場(chǎng)干涉現(xiàn)象基本消失?？v向間隔1.0D時(shí)，1#旋翼下洗流不經(jīng)過2#旋翼旋轉(zhuǎn)盤面，即使兩者的下洗流場(chǎng)有所干涉，但影響不大；當(dāng)縱向間隔在2.0D及以上時(shí)，2#旋翼來流即使受到1#旋翼下洗流影響，但由于距離較遠(yuǎn)，來流速度較低，渦強(qiáng)度也較弱，未對(duì)2#旋翼產(chǎn)生明顯影響。

2.3 不同間隔下雙旋翼的氣動(dòng)參數(shù)結(jié)果

當(dāng)旋翼間存在氣動(dòng)干擾時(shí)，改變流場(chǎng)分布同時(shí)會(huì)引起旋翼氣動(dòng)力和氣動(dòng)力矩的顯著變化。圖11～13分別給出了雙旋翼在不同間隔下的拉力系數(shù)、扭矩系數(shù)、俯仰力矩系數(shù)變化曲線。1#旋翼的氣動(dòng)力和力矩系數(shù)并不受間隔的影響，2#旋翼的氣動(dòng)參數(shù)則變化劇烈。

圖11 不同間隔的拉力系數(shù)曲線Fig.11 Tension coefficient curves under different separations

圖12 不同間隔的扭矩系數(shù)曲線Fig.12 Torque coefficient curves under different separations

圖13 不同間隔的俯仰力矩系數(shù)曲線Fig.13 Pitching moment coefficient curves under different separations

從拉力系數(shù)變化曲線可看出，2#旋翼處于1#旋翼的下洗流中時(shí)，造成其旋翼拉力下降；當(dāng)2#旋翼位于1#旋翼的正下方時(shí)，其拉力損失最大，且縱向間隔越小，拉力減小越明顯。間隔X=0.0D，Z=(0.5-1.5)D范圍內(nèi)，拉力系數(shù)減小約50%；當(dāng)橫向間隔X=0.5D和1.0D時(shí)分別減小約20%和10%；此種現(xiàn)象主要原因是上游旋翼下洗氣流產(chǎn)生向下的誘導(dǎo)速度，減小了旋翼剖面迎角，進(jìn)而減小了拉力。當(dāng)橫向間隔X≥1.5D以上時(shí)，2#旋翼受來流的向下誘導(dǎo)速度影響較小，其拉力則基本不受影響。

扭矩系數(shù)反映了旋翼氣動(dòng)阻力的大小。從圖12可以看出，2#旋翼扭矩系數(shù)受影響最嚴(yán)重區(qū)域?yàn)闄M向間隔X=0.0D，即正下方區(qū)域，其扭矩減??；且隨著縱向間隔越小，扭矩減小越明顯，當(dāng)Z=0.5D時(shí)減小約20%。而當(dāng)存在橫向間隔時(shí)，除了X=0.5D,Z≤1.0D區(qū)域，2#旋翼扭矩則相對(duì)增大，不過增大幅度在10%以內(nèi)。此主要原因如上文所述，當(dāng)2個(gè)旋翼橫向錯(cuò)置時(shí)，2#旋翼的半側(cè)旋翼主要位于1#旋翼下洗流中，其旋轉(zhuǎn)線速度方向與下洗流的水平分速度方向相反，因此增加了旋轉(zhuǎn)的氣動(dòng)阻力；當(dāng)無橫向間隔時(shí)，2#翼的旋轉(zhuǎn)方向與來流水平分速度方向相同，則減小了氣動(dòng)阻力。

俯仰力矩系數(shù)反映了旋翼當(dāng)前的水平穩(wěn)定狀態(tài)，當(dāng)俯仰力矩系數(shù)為零時(shí)，旋翼可保持水平狀態(tài)；當(dāng)此系數(shù)大于0時(shí)，則說明此旋翼在俯仰力矩作用下，有傾覆的危險(xiǎn)。從圖13中可以看出，當(dāng)無橫向間隔X/D=0時(shí)，雖然2#旋翼完全處于下洗流區(qū)，但由于流場(chǎng)仍可保持對(duì)稱性，兩側(cè)的半旋翼氣動(dòng)力的變化相差不大，其俯仰力矩系數(shù)很小。X=0.5D且1.0D≤Z≤3.0D區(qū)間為最危險(xiǎn)區(qū)域。從上文的流場(chǎng)分析可知，此時(shí)2#旋翼的流場(chǎng)出現(xiàn)顯著的非對(duì)稱特征，兩側(cè)旋翼的氣動(dòng)力不平衡。當(dāng)橫向間隔X/D繼續(xù)增大，2#旋翼逐漸離開下洗流區(qū)，非對(duì)稱流場(chǎng)現(xiàn)象減弱，俯仰力矩系數(shù)也逐漸減小。

3 結(jié)論

1)雙旋翼氣動(dòng)干擾，1#旋翼流場(chǎng)和氣動(dòng)參數(shù)基本不受影響，對(duì)2#旋翼的氣動(dòng)特性影響顯著。當(dāng)橫向間隔X≤1.0D時(shí)，2#旋翼速度場(chǎng)和渦量場(chǎng)擴(kuò)散范圍受限；X≥1.5D時(shí)，2#旋翼離開干擾區(qū)，流場(chǎng)和氣動(dòng)參數(shù)影響較小。

2)2#旋翼位于1#旋翼正下方即無橫向間隔時(shí)，其流場(chǎng)仍可保持對(duì)稱性；此范圍2#旋翼所受氣動(dòng)力變化主要為拉力和扭矩的降低，且縱向間隔越小表現(xiàn)越明顯，當(dāng)Z≥2.5D時(shí)此影響基本消失。

3)2#旋翼相對(duì)1#旋翼橫向間隔0.5D≤X≤1.0D時(shí)，其流場(chǎng)出現(xiàn)明顯非對(duì)稱特征；X=0.5D且1.0D≤Z≤3.0D區(qū)域，2#旋翼受較大俯仰力矩作用，為危險(xiǎn)區(qū)域，懸停時(shí)應(yīng)避免進(jìn)入。