四旋翼無人機旋翼對機身非定常氣動干擾特性

朱毅飛 ,林德福 ,莫靂 ,葉建川

(1.北京理工大學 宇航學院,北京 100081;2.北京理工大學 無人機自主控制技術北京市重點實驗室,北京 100081)

0 引言

目前,隨著科技的發展和需求的不斷增長,無人機產業迅速發展,在軍用和民用方面,無人機的地位日益提高。多旋翼無人機作為無人機中最為廣泛的一個種類,被廣泛應用于軍事偵查、軍事打擊[1-3]、災害救援[4-6]、植物保護、物流運輸、勘探等領域。由于這類飛行器旋轉升力面多,耦合干擾強烈,非定常流動明顯,亟需深入分析其氣動特性的機理。然而當前多旋翼無人機雖然工程化樣機及產品眾多,但是機理研究尚有待進一步深化。

Hwang 等[7]采用數值計算方法,對不同構型無人機的氣動特性和旋翼之間的氣動干擾作了研究,采用滑移網格的方法來建立旋轉的計算域,所得結果表明,氣動干擾存在于懸停工況和前飛工況下,對于“+”構型和“×”構型,由于旋翼尾流誘導出的上洗氣流和下洗氣流特性不同而導致氣動干擾的程度不同。進而對兩種不同構型無人機旋翼之間的氣動干擾進行了分析與對比,卻未對旋翼與機身的干擾進行探究。

Barcelos 等[8]用勢流法通過對比單獨的旋翼,發現有關“+”構型和“×”構型之間力和力矩存在差異。結果表明,“+”構型無人機在順時針旋轉的前旋翼下有向右滾轉的趨勢,在逆時針旋轉的前旋翼下有向左滾轉的趨勢。對于高前進比的“×”構型無人機,Breaststroke 布局(無人機對稱面附近的前槳槳尖為后行)的旋翼效率略高于Bear-hug 布局(無人機對稱面附近的前槳槳尖為前行),俯仰力矩配平需求較小。這種現象是由前行槳葉和后行槳葉誘導的上洗氣流和下洗氣流不對稱導致的?！?”構型無人機比“×”構型無人機有更高的旋翼效率;“+”構型無人機比“×”構型無人機需要更多的滾轉配平力矩,而“×”構型無人機則需要更多的俯仰配平力矩。

Luo 等[9]通過理論推導和計算流體動力學(CFD)仿真,得出了一種考慮尾跡干擾的四旋翼無人機前飛的數學模型,可用于飛行控制。研究表明,高速前飛模型可以通過Glauert 高速近似法[10]模擬為圓形固定翼,簡化了旋翼在前飛過程中誘導處的復雜渦系,推導出了前飛狀態下推力系數的計算公式,建立了考慮相互干擾的四旋翼前飛模型。但該研究未考慮非定常的氣動干擾影響。

Fu 等[11]對旋翼進行補充實驗,探究了滑流的氣動特性,同時將拼接結構網格、滑移網格和雷諾平均納維-斯托克斯(RANS)方程的計算研究作為依據來證明結論。結果表明,旋翼氣動計算結果與實驗結果吻合較好,能正確預測旋翼的滑流特性,旋翼之間的干擾對飛行器的性能有重要影響。

Misiorowaski 等[12]采用包含分離渦仿真模型和滑移網格的求解器AcuSolve 來求解Navier-Stokes 方程,以驗證四旋翼無人機在側向飛行條件下的性能。在10 m/s 前飛和低頭5°俯仰角的工況下,開展2 個不同構型無人機——“×”構型無人機和“+”構型無人機的實驗。在“×”構型無人機中,由于前槳的下洗作用,后槳存在拉力和扭矩損失,同時俯仰力矩也有所減小。在“+”構型無人機中,相比于單旋翼而言,東側和西側的旋翼由于北側旋翼的上洗作用產生額外的拉力。由此可見,由于“×”構型無人機中南部旋翼的前部拉力損失較大,后槳需要更高的轉速。

綜合上述研究結果,氣動干擾對無人機的氣動特性將產生重要影響,其中作用在機身上的力和力矩是影響無人機飛行性能優劣的重要因素[13-16],然而目前,對無人機旋翼與機身瞬態氣動干擾的研究仍顯不足,無法獲得無人機旋翼與機身氣動干擾在時間上的分布規律及干擾的產生機理。

本文通過CFD 計算方法[15]對無人機旋翼和機身的氣動干擾的瞬態變化進行分析,總結無人機旋翼和機身氣動干擾的基本規律,并獲得適用于不同需求下不同精度的擬合方程。

在計算方法上,一是采用四面體非結構網格劃分計算域;二是使用滑移網格方法建立旋轉域[9,17-19,23],模擬旋翼的旋轉;三是在壁面附近采用棱柱網格,將第1 層網格節點布置在黏性底層之內,令y+＜5(y+為無量綱壁面距離,選取合適的y+值可以提高邊界層的求解精度),來提高湍流的計算精度[21],捕捉流動細節。其次假設無人機機身及旋翼為剛體,忽略槳葉的形變及揮舞[12,20]。另外采用非定常計算方法[9,19],計算機身在旋翼旋轉的每個周期內所受的力和力矩。最后,通過分析CFD 方法計算所得到的數據,獲得了機身所受力和力矩的變化規律,并給出機身升力系數、阻力系數和俯仰力矩系數的擬合方法。本文結論是基于某一典型前飛工況得出的,但在不同飛行工況下無人機機身氣動干擾特性的分析方法是一致的。本文可以為相似構型和工況的無人機機身氣動干擾分析提供基本方法和參考。

1 模型及坐標系定義

無人機模型使用三維計算機輔助設計軟件Solidworks 生成,其三視圖以及軸測圖如圖1 所示。圖1(c)中,R為旋翼半徑,R=259.08 mm,d為旋翼中心到無人機縱向對稱面的距離,d=313.85 mm,yh為機身寬度的一半,yh=45.00 mm,Ωf為前旋翼轉速,Ωr為后旋翼轉速。

圖1 給出了無人機模型的部分尺寸和旋翼編號及旋轉方向。在該模型中,螺旋槳選用中國TMOTOR 公司生產的MF2009 型多軸旋翼碳纖維塑料聚合物正反折疊槳,直徑為20.4 in(518.16 mm)。螺旋槳旋轉域直徑為550 mm,厚度為30 mm。在圖1(c)中,旋翼1、旋翼3 為逆時針旋轉,旋翼2、旋翼4 為順時針旋轉[8]。無人機前飛時,為了保持其前飛姿態,需要額外的低頭力矩來抵消恢復力矩,因此前飛時的無人機后旋翼轉速一般高于前旋翼,而由于無人機的對稱性,旋翼2 和旋翼3的轉速相同,旋翼1 和旋翼4 的轉速相同。需要特別說明的是,本文的研究是基于圖1(c)中所標注的旋翼轉向進行的。在該轉向下,前旋翼靠近縱向對稱面的一側為前行槳,槳尖的局部速度大;后旋翼靠近縱向對稱面的一側為后行槳,槳尖的局部速度小。然而,若無人機的旋翼轉向與本文所規定的相反,則前旋翼靠近縱向對稱面的一側為后行槳,后旋翼靠近縱向對稱面的一側為前行槳。在這種情況下分析出的最終結果有可能與本文所述結論有所不同,但計算和分析方法通用。因此,對于旋翼反轉的情況,本文的內容仍然具有參考意義。

圖1 無人機模型Fig.1 Quadcopter model

無人機機體坐標系Obxbybzb的定義如圖2 所示:坐標系原點Ob位于無人機的質心處,xb軸指向無人機尾部,zb軸垂直于槳盤平面指向無人機上方,yb軸方向遵循右手定則。

圖2 無人機機體坐標系Fig.2 Body fixed coordinate

假設無人機無側滑和滾轉,速度坐標系Ovxvyvzv可以由機體坐標系Obxbybzb通過繞機體坐標系yb軸旋轉α得到,α為攻角。

無人機在機體坐標系Obxbybzb和速度坐標系Ovxvyvzv下的力、力矩、速度、角度定義如圖3 所示。圖3 中:分別為無人機機身受力在機體坐標系xb軸、yb軸、zb軸上的分量;分別為無人機機身所受力矩在機體坐標系xb、yb、zb軸上的分量;Lf為無人機機身所受升力;Df為無人機機身所受阻力;分別為旋翼受力在機體坐標系xb軸、yb軸、zb軸上的分量;分別為旋翼所受力矩在機體坐標系xb軸、yb軸、zb軸上的分量;Lb為無人機旋翼所受升力;Db為無人機旋翼所受阻力;T為旋翼拉力;H為旋翼Hub 力;Q為旋翼扭矩;v∞為來流速度。

圖3 無人機力、力矩、角度定義Fig.3 Definition of force,moment and angle

圖3(a)、圖3(c)為無人機機身的力、力矩、角度定義,圖3(b)、圖3(d)為無人機某一槳盤的力、力矩、角度定義。對于機身,機身所受力和力矩的方向與機體系正方向一致,阻力方向與來流方向一致,升力與阻力垂直,并且指向無人機上方。攻角為來流方向與機體坐標系xb軸之間的夾角,無人機低頭時為正。對于槳盤,由于旋翼與機身之間沒有安裝角,機身攻角也是槳盤攻角,槳盤的力的方向定義與機身相同。

因此,可得機體坐標系與速度坐標系下相關力和力矩與攻角的轉換關系:

正升力提供的是使無人機上升的力,正俯仰力矩的作用效果使無人機抬頭。由于無人機前飛時始終處于低頭狀態,單獨機身的升力一般為負值。

由于本文是對無人機的非定常氣動干擾進行探究,需要對于初始狀態,即t=0 s 時刻下的狀態進行特別說明。

圖4 為無人機的初始狀態,其中φf為前旋翼的方位角,φr為后旋翼的方位角,前、后旋翼的初始角度分別為φf0=55°、φr0=55°。

圖4 無人機初始狀態Fig.4 Initial state of UAV

2 CFD 計算方法

下面首先介紹本文采用的CFD 計算工況、求解模型以及網格生成策略,并進行網格敏感性分析,然后驗證計算域劃分的合理性,最后進行風洞實驗驗證以驗證CFD 計算的準確性。

2.1 CFD 計算工況

無人機所處的飛行工況如表1 所示。表1 中,ff為無人機前旋翼的旋轉頻率,fr為無人機后旋翼的旋轉頻率。

表1 無人機飛行工況Tab.1 Flight conditions

定義f為無人機旋翼的旋轉頻率:

由于每片旋翼由2 片槳葉組成,旋翼每旋轉半周即可被視作一個周期,旋轉一周則對應2 個周期。

CFD 計算方法采用非定常流動模型,流動介質為不可壓縮假設下的空氣,氣體密度為1.225 kg/m3,采用Pressure-Based Solver 求解器進行求解,湍流模型為可實現k-ε二方程湍流模型[9],該模型在標準k-ε模型的基礎上提供了旋流修正,對旋轉流動、流動分離計算有更好的性能,并采用半隱式壓力速度耦合(SIMPLE)算法和2 階迎風格式進行計算[23]。采用四面體的非結構網格來劃分計算域,在靠近壁面處以及流動較為復雜的區域進行了加密處理,邊界層通過棱柱網格來進行較為精確的模擬[21],保證近壁面處的y+＜5。采用滑移網格[9],用真實的運動對旋翼旋轉進行相對精確的模擬。圖5 通過對計算域中某一切面的展示說明了網格劃分方法,展示了網格劃分的效果。

圖5 計算域網格劃分Fig.5 Mesh generation of computational domain

由于需要對無人機旋翼對機身所產生的影響進行分析,基于圖6(a)有旋翼模型建立了如圖6(b)所示的純機身模型,并在相同的計算條件下進行計算。

圖6 有旋翼、無旋翼無人機模型Fig.6 UAV models with/without rotors

2.2 網格敏感性分析

為了提高計算效率,要求在保證精度的情況下網格數量盡可能少。因此,將2 種不同尺寸的網格在相同計算條件下進行結果對比,計算條件采用表1 中所示的速度和攻角。為了縮短計算時間,采用無旋翼的純機身模型進行對比,計算結果如表2 所示。其中,CL,f、CD,f、CM,f分別為無旋翼無人機的機身升力系數、阻力系數、俯仰力矩系數。

表2 網格敏感性分析結果對比Tab.2 Comparison of mesh sensitivity analysis results

由表2 可知,2 組網格數量分別為310 萬和510 萬,計算所得升力系數誤差的絕對值為3.62%,阻力系數誤差的絕對值為10.38%,俯仰力矩系數誤差的絕對值為4.11%,因此可認為310 萬網格計算的結果精度滿足要求。

2.3 風洞實驗驗證

風洞實驗在天津大學航空航天研究院進行。該風洞的截面尺寸為1.2 m×1.2 m,能提供的最大風速為30 m/s。由于風洞尺寸的限制,全尺寸無人機模型無法進入風洞進行實驗,因此采用無人機所配螺旋槳進行單個旋翼的風洞測試,通過記錄并對比單個旋翼的受力情況來驗證CFD 計算方法的準確性[22]。

圖7 為風洞和實驗用旋翼。單個旋翼的測試條件如表3 所示,風洞實驗時旋翼攻角與圖3(d)中槳盤攻角的定義相同。

圖7 風洞和實驗用旋翼Fig.7 Wind tunnel and test rotor

表3 風洞實驗測試條件Tab.3 Wind tunnel test conditions

采用德國ME-Meβsysteme 公司生產的MEK6D40 傳感器來測量單個旋翼所受的拉力T,采樣率為1 000 Hz,并對風洞實驗數據和CFD 計算方法得出的數據進行對比,結果如表4 和圖8 所示。

表4 風洞實驗拉力對比Tab.4 Comparison of tensile forces T

圖8 不同轉速下風洞實驗與CFD 計算結果對比Fig.8 Comparison between wind tunnel test result and CFD calculated result at different rotation speeds

由表4 和圖8 可知:CFD 計算結果與風洞實驗數據誤差在15%～20%之間,并且隨著轉速升高誤差逐漸增大。本文中無人機旋翼的轉速在3 200 r/min 左右,在該轉速下誤差大約控制在15%以內,滿足精度要求,因此可以認為CFD 計算方法結果可靠。

3 氣動特性分析

對無人機旋翼-機身氣動干擾特性進行分析,包括機身的升力特性、阻力特性和俯仰力矩特性。對機身而言,若忽略旋翼所產生的誘導速度,則升力系數、阻力系數、俯仰力矩系數可以由(4)式～(6)式求得:

式中:ρ為空氣密度;S為參考面積,取無人機機身前視圖的投影面積;l為參考長度,取無人機機身的長度;Mf為俯仰力矩。對于CFD 計算結果,氣動干擾產生的力和力矩系數可以通過將兩個無人機模型各自的CFD 計算數據相減獲得:

式中:CL,u、CD,u和CM,u分別為有旋翼無人機的升力系數、阻力系數、俯仰力矩系數。

3.1 升力特性分析

圖9 為ΔCL快速傅里葉變換(FFT)后的結果。由圖9 可知,頻率最低的峰值出現在103.8 Hz 和111.1 Hz 處,且頻率峰值都是成對出現,高頻分量與最低頻率呈倍數關系,且高于2 階的模態其對ΔCL周期變化的貢獻以近似線性的規律逐漸減小。對比表1 的結果,可以發現2ff=104.30 Hz≈103.8 Hz,2fr=110.77 Hz≈111.1 Hz。

圖9 ΔCL FFT 結果Fig.9 FFT results of ΔCL

FFT 結果的最低頻分量與無人機前后旋翼的轉速頻率的2 倍近似相等,即ΔCL變化的1 階模態與無人機旋翼的轉速有直接關系。從結果的幅值來看,前旋翼對于機身升力的影響顯著大于后旋翼;頻率越高,對應分量的影響越小。因此,旋翼對升力的影響可以用傅里葉級數來表示如下:

表5 和表6 對FFT 所得到的幅值結果與擬合方程中的系數做了對比,為了使結果更直觀,取擬合方程系數的絕對值進行比較,計算出的誤差為兩項變化的絕對值。結果表明,擬合系數與FFT 獲得的幅值誤差普遍在7%以下,最大不超過15%,吻合情況良好。故推測,前5 階模態的信息能夠反映ΔCL變化絕大部分信息,忽略高于5 階模態不會對擬合精度造成大的影響,前5 階模態能夠精確地體現出旋翼對機身的氣動干擾特點。

表5 前旋翼擬合系數對比Tab.5 Comparison of fitting coefficients of front rotor

表6 后旋翼擬合系數對比Tab.6 Comparison of fitting coefficients of rear rotor

根據時間t和前后旋翼轉速的關系

可得擬合后的結果,如圖10 所示。

圖10(a)為前旋翼方位角在2 000°～8 000°之間ΔCL的傅里葉級數擬合結果以及原始數據。由圖10 可知,1～5 階模態傅里葉級數能夠很好地描述無人機機身ΔCL的影響,擬合曲線能反映原始數據的大部分信息,并且描述原始數據的變化細節。

為了進一步簡化方程的形式,在只保留單獨1 階模態信息的情況下,ΔCL可以寫為

因此,擬合方程可寫為

在忽略2 階及更高階模態以后,ΔCL的擬合結果如圖11 所示。與圖10 相比,圖11 中擬合曲線的變化幅值有所減小,這是因為高階模態信息被剔除后,擬合曲線失去了這一部分對幅值的貢獻,從而導致了變化幅值的減小。同時可以發現,擬合曲線中的諸多細節被忽略。例如圖11(b)中紅圈所示的部分,其變化趨勢在圖10(b)中能夠通過高階模態項在一定程度上被表達出來,但單1 階模態無法描述該變化情況。然而,雖然單1 階模態相比1～5 階模態的擬合結果而言無法描述出原始數據的一些細節變化情況,但仍能反映出ΔCL的總體變化趨勢,即ΔCL總體上由2 個周期性變化構成,其中,短周期描述了ΔCL隨旋翼旋轉的變化規律,長周期描述了前、后旋翼轉速差的影響,ΔCL的長、短周期均以正弦規律隨時間變化。

圖10 ΔCL的1～5 階傅里葉級數擬合結果Fig.10 Fitting result of the first-fifth-order Fourier series of ΔCL

圖11 ΔCL的單1 階傅里葉級數擬合結果Fig.11 Fitting result of single first-order Fourier series of ΔCL

3.2 阻力特性分析

ΔCD的擬合結果如圖12 所示。由圖12 可以看出,ΔCD與ΔCL有相似的頻率分布,1 階模態和2 階模態對ΔCD的貢獻相似,高于2 階的模態其對ΔCD周期變化的貢獻以近似線性的規律逐漸減小。ΔCD的傅里葉級數可以寫為

圖12 ΔCD的FFT 結果Fig.12 FFT results of ΔCD

同樣,保留前5 階模態的信息,(16)式可以寫為

經擬合,可得(17)式的具體表達式為

圖13 為前旋翼方位角在2 000°～8 000°之間ΔCD的傅里葉級數擬合結果以及原始數據。由圖13中可知,包含1～5 階模態的擬合曲線能夠很好地與CFD 計算方法結果吻合,能反映出ΔCD變化的細節部分。從原始數據的走勢分析,若采用單1 階模態的方程擬合,則擬合曲線與原始數據之間會存在較大的差別。

圖13 ΔCD的1～5 階傅里葉級數擬合結果Fig.13 Fitting result of the first-fifth-order Fourier series of ΔCD

表7 計算并歸納了1～5 階模態傅里葉級數擬合結果與原始數據的R2值。從表7 中數據來看,1～3 階模態下的擬合結果所對應的R2值為0.816 0,大于0.8。因此,為簡化擬合方程形式同時保證擬合精度,保留1～3 階模態信息,經擬合,包含1～3 階模態信息的ΔCD傅里葉級數擬合結果可以寫為

表7 ΔCD各階模態傅里葉級數擬合結果R2值Tab.7 R2 of fitting results of ΔCD with different Fourier series

圖14 為包含1～3 階模態信息的ΔCD傅里葉級數擬合結果,與ΔCL相比,結果更為復雜。圖14(b)與圖13(b)的周期相同,圖14(c)為另一周期內的結果。由圖14 可知,相比于1～5 階模態的傅里葉級數的擬合結果,1～3 階模態能夠在一定程度上描述ΔCD的變化規律,但在幅值上有所損失,并且1～3 階模態擬合曲線的峰值所對應的相位與原始數據之間存在少量偏差。通過對比圖14(b)和圖14(c)不難發現,隨著時間的推移,ΔCD在不同周期內的變化規律也有所不同。圖14(b)中,ΔCD的第1 個突變出現在φf=105.1°,第1 個峰值出現在φf=138.9°;在圖14(c)中,ΔCD的第1 個突變出現在φf=65.71°,第1 個峰值出現在φf=140.8°。圖14(b)、圖14(c)所對應的周期中,突變出現的位置存在大約40°的相位差,但峰值與谷值的相位差距很小。表明ΔCD的變化中,波峰和波谷的變化是由于前旋翼的旋轉所造成的,而突變的產生與其他高階模態的信息有關。因此,1～3 階模態能夠反映出ΔCD的總體變化,但無法準確地描述突變的產生。根據圖13 中的結果,對于產生突變的高階模態,至少需要1～5 階傅里葉級數才能進行較為準確的擬合。

圖14 ΔCD的1～3 階傅里葉級數擬合結果Fig.14 Fitting result of the firsst-third-order Fourier series of ΔCD

3.3 俯仰力矩特性分析

圖15 為ΔCM在FFT 后得到的頻率分布。由圖15 可知,ΔCM的頻率分布與ΔCL和ΔCM的規律一致,最低頻分量都為103.8 Hz 和111.1 Hz,且高頻分量為最低頻分量的倍數,高于2 階的模態其對ΔCM周期變化的貢獻以近似線性的規律逐漸減小。但不同的是,1 階模態對ΔCM產生的影響小于2～5 階,且對于1 階模態而言,后旋翼的貢獻略大于前旋翼。同樣,ΔCM可以用形如(20)式的傅里葉級數描述出來:

圖15 ΔCM的FFT 結果Fig.15 FFT result of ΔCM

圖16(a)為前旋翼方位角在2 000°～8 000°之間ΔCM的傅里葉級數擬合結果以及原始數據,圖16(b)為圖16(a)中某一個周期結果的局部圖,該周期與圖10(b)中的周期相同。由圖16 可知,除最小值附近的某些數據點外,1～5 階傅里葉級數能夠較好地反映出ΔCM的變化規律,包括一些變化細節。與圖8 中的ΔCL對比,發現圖16(b)中ΔCM第1 個波峰和第1 個波谷出現的相位與ΔCL相同。由此可以推斷,ΔCM在φf=97.62°和φf=146.4°出現的突變是由于機身升力的突變所引起的,在φf=97.62°時,造成機身升力減小的因素出現在無人機質心后部;在φf=146.4°時,造成機身升力減小的因素出現在無人機質心前部,因此造成了ΔCM先突然增大、后突然減小的現象。

圖16 ΔCM的1～5 階傅里葉級數擬合結果Fig.16 Fitting result of the first-fifth order Fourier series of ΔCM

由圖15 所得結論,ΔCM受2 階模態的影響顯著大于1 階模態,因此為簡化擬合方程,在只保留單2 階模態信息的情況下,ΔCM的傅里葉級數擬合結果為

圖17 所示的擬合結果表明,單2 階模態傅里葉級數的擬合結果能夠反映ΔCM的變化趨勢,但無法準確描述出波峰和波谷的幅值,并且在相位上也有誤差?？傮w上來說,ΔCM與正弦函數的變化趨勢相似,并且變化頻率與基頻的2 倍相關性較高。

圖17 ΔCM的2 階傅里葉級數擬合結果Fig.17 Fitting result of the second-order Fourier series of ΔCM

為了更精確地描述出ΔCM的變化,總結了包含各個模態方程擬合結果的R2值,如表8 所示。

表8 ΔCM各階模態傅里葉級數擬合結果R2值Tab.8 R2 of fitting results of ΔCM with different orders of Fourier series

從表8 中數據來看,1～4 階模態的擬合結果所對應的R2值0.818 7 大于0.8。這說明方程至少包含1～4 階的模態信息才能夠以較高的精度反映原始數據的變化情況。保留1～4 階模態信息,經過擬合,ΔCM的傅里葉級數擬合結果為

圖18 為1～4 階模態傅里葉級數的擬合結果。從圖18 中可以看到,擬合結果能夠較為精確地描述ΔCM的變化趨勢。與圖17 中的擬合曲線相比,1～4 階模態包含的信息能夠更好地描述波峰和波谷的幅值,反映出了無人機機身升力突變對俯仰力矩造成的突變影響,并且擬合曲線在相位上的誤差也有所減小。因此,對于ΔCM而言,最少需要1～4 階模態的信息才能反映變化的細節,而單2 階模態只能夠描述出ΔCM在周期變化上的規律。若只對ΔCM的周期變化感興趣,則2 階模態可以滿足要求;若需要獲得更多關于ΔCM的變化信息,則最少需要1～4 階模態。

圖18 ΔCM的1～4 階傅里葉級數擬合結果Fig.18 Fitting result of the first-fourth-order Fourier series of ΔCM

4 結論

本文采用數值計算方法,對典型構型的四旋翼無人機旋翼與機身之間的非定常氣動干擾進行了分析。通過對CFD 計算結果的擬合與分析,得出如下主要結論:

1)無人機機身升力系數的變化ΔCL與無人機旋翼的旋轉周期密切相關,并且前旋翼的影響程度顯著大于后旋翼。ΔCL變化的1 階模態頻率為無人機前、后旋翼的旋轉頻率的2 倍,而更高階模態的頻率則是1 階模態頻率的整數倍頻。

2)無人機機身升力系數變化ΔCL可以用1～5 階傅里葉級數較為精確地描述。單1 階模態的傅里葉級數能夠描述出ΔCL的總體變化趨勢?？傮w上說,ΔCL由2 個周期變化構成,其中,短周期描述了ΔCL隨旋翼旋轉方位角的變化規律,長周期描述了前后旋翼轉速差的影響。ΔCL的長短周期均以正弦規律隨時間變化。

3)無人機機身阻力系數的變化ΔCD可以用1～5 階傅里葉級數較為精確地描述出來。其中1～3 階傅里葉級數能夠描述出ΔCD的總體變化情況。前旋翼對于ΔCD的影響較大,影響ΔCD曲線的波峰和波谷,但其中的突變需要更高階模態的信息來體現。

4)無人機機身俯仰力矩系數的變化ΔCM可以用1～5 階傅里葉級數較為精確地描述出來。ΔCM的突變是由于機身升力的突變所引起的。ΔCM受2～5 階模態的影響較大。其中2 階模態能夠描述出ΔCM在周期變化上的規律,而1～4 階模態包含的信息能夠更好地描述波峰和波谷的幅值,反映出了無人機機身升力突變對俯仰力矩造成的突變影響?？傮w而言,ΔCM與正弦函數的變化趨勢相似,并且變化頻率與基頻的2 倍相關性較高。

5)對于無人機機身升力、阻力、俯仰力矩受旋翼氣動干擾影響的規律,可以根據實際應用中所關注的重點來靈活選取擬合模型。若重點關注其細節變化,則需要采用包含更多模態信息的高階傅里葉級數擬合結果;若只關心其周期變化的規律,不著重考察變化中的細節,則可以采用相對簡化的模型來代替。