吳卓恒，余 莉,2，張思宇，賈 賀

(1.飛行器環(huán)境控制與生命保障工業(yè)和信息化部實(shí)驗(yàn)室(南京航空航天大學(xué))，南京 210016；2.南京航空航天大學(xué) 航空學(xué)院，南京 210016；3.北京空間機(jī)電研究所，北京 100094)

自旋翼無(wú)人機(jī)因其飛行性能好、不受地形限制等優(yōu)點(diǎn)，被廣泛用于地質(zhì)探測(cè)、搜救偵查以及科學(xué)研究等領(lǐng)域。但是，當(dāng)在飛行過(guò)程中出現(xiàn)故障時(shí)，可能會(huì)發(fā)生嚴(yán)重的安全事故甚至墜毀。目前，已經(jīng)廣泛采用降落傘對(duì)固定翼無(wú)人機(jī)進(jìn)行減速回收[1]，而對(duì)旋翼機(jī)回收的研究工作較少。當(dāng)失控狀態(tài)的旋翼機(jī)處于垂直下降時(shí)，旋翼產(chǎn)生的旋轉(zhuǎn)流場(chǎng)[2]加劇了前體尾流的非定常變化，導(dǎo)致流向降落傘的流場(chǎng)變化紊亂.為探究旋翼擾動(dòng)下物傘系統(tǒng)的流場(chǎng)特性以及傘衣氣動(dòng)性能變化，需要對(duì)旋翼機(jī)/降落傘的非定常復(fù)合流場(chǎng)展開(kāi)深入的研究。

目前，旋翼流場(chǎng)的研究主要通過(guò)風(fēng)洞試驗(yàn)和數(shù)值模擬，風(fēng)洞試驗(yàn)?zāi)軌蛘鎸?shí)反映旋翼的工作狀態(tài)[3-4]，通過(guò)控制旋翼的狀態(tài)參數(shù)獲得準(zhǔn)確的試驗(yàn)數(shù)據(jù)，但存在洞壁效應(yīng)[5]、實(shí)驗(yàn)數(shù)據(jù)可重復(fù)性差[6]、成本代價(jià)高等問(wèn)題。隨著計(jì)算流體力學(xué)(computational fluid dynamics，CFD)的發(fā)展，數(shù)值模擬逐漸成為研究旋翼流場(chǎng)的重要手段，該方法通過(guò)直接求解流場(chǎng)控制方程(Euler方程或Navier-Stokes方程)[7-8]，能真實(shí)還原流場(chǎng)的細(xì)節(jié)變化，具備較高的求解精度。如文獻(xiàn)[9-11]采用CFD方法對(duì)懸停狀態(tài)下的旋翼流場(chǎng)進(jìn)行了數(shù)值模擬，得到了槳葉附近真實(shí)的流場(chǎng)細(xì)節(jié)變化；文獻(xiàn)[12-13]觀察到前飛狀態(tài)下的旋翼槳尖渦，分析了旋翼尾跡渦和下洗流的非定常流場(chǎng)特征變化；文獻(xiàn)[14]發(fā)現(xiàn)旋翼垂直下降流場(chǎng)會(huì)出現(xiàn)渦環(huán)狀態(tài)，隨下降速度增加旋翼下洗流逐漸上移直至被來(lái)流抵消。

眾多學(xué)者也廣泛開(kāi)展了旋翼復(fù)合流場(chǎng)的數(shù)值研究，如文獻(xiàn)[15-17]對(duì)傾轉(zhuǎn)旋翼機(jī)/機(jī)翼的復(fù)合流場(chǎng)進(jìn)行了CFD數(shù)值模擬，精確捕捉到了流場(chǎng)的旋渦結(jié)構(gòu)，觀察到旋翼/機(jī)身產(chǎn)生的“噴泉效應(yīng)”[15-16]現(xiàn)象；文獻(xiàn)[18-20]對(duì)旋翼與艦船的復(fù)合流場(chǎng)進(jìn)行數(shù)值模擬，發(fā)現(xiàn)旋翼繞流與艦船尾流發(fā)生疊加[18]，會(huì)產(chǎn)生旋渦、回流、分離[19]等多種復(fù)雜的流動(dòng)形式，在側(cè)風(fēng)的摻混下進(jìn)一步形成渦流紊亂區(qū)，渦流區(qū)范圍和氣流下洗明顯增加，甲板流場(chǎng)變化極不穩(wěn)定[20]。目前，對(duì)旋翼復(fù)合流場(chǎng)的數(shù)值研究還有很多[21-22]，但旋翼機(jī)/降落傘非定常復(fù)合流場(chǎng)的研究工作卻很少。

本文基于PISO(pressure implicit split operator)算法對(duì)旋翼機(jī)/降落傘復(fù)合流場(chǎng)展開(kāi)非定常數(shù)值計(jì)算，為提高網(wǎng)格更新效率，采用了不同變形尺度網(wǎng)格分類處理的思路。隨后，研究了旋翼轉(zhuǎn)動(dòng)對(duì)物傘系統(tǒng)流場(chǎng)分布和旋渦演變的非定常影響，分析了不同旋翼轉(zhuǎn)速下流場(chǎng)特性參數(shù)和降落傘氣動(dòng)特性的變化情況，為進(jìn)一步探究旋翼機(jī)與降落傘的氣動(dòng)干擾提供一定參考。

1 計(jì)算方法

1.1 數(shù)值方法及湍流模型

本文基于有限體積法采用MUSCL(monotone Upstream-centered schemes for conservation laws)三階格式進(jìn)行流場(chǎng)方程離散，為提高瞬態(tài)計(jì)算精度并加速收斂采用了PISO算法，壓力插值選擇Standard格式，為有效處理瞬態(tài)計(jì)算偽擴(kuò)散問(wèn)題選擇Green-gauss Node-based格式梯度插值。邊界條件采用速度進(jìn)口和壓力出口，遠(yuǎn)場(chǎng)為對(duì)稱邊界條件，傘衣面以及旋翼機(jī)(包括前體和旋翼)為無(wú)滑移邊界條件，旋轉(zhuǎn)區(qū)域Ωr和自由流區(qū)域Ωf的交界面Γinterface為Ineterface邊界條件。

為獲得精細(xì)的流場(chǎng)結(jié)構(gòu)、有效求解固壁繞流問(wèn)題，本文計(jì)算采用了標(biāo)準(zhǔn)壁面函數(shù)法(standard wall functions)，湍流模型選擇Realizablek-ε兩方程模型，其控制方程為

(1)

1.2 網(wǎng)格更新模型

迭代格式的時(shí)間推進(jìn)方法在每個(gè)時(shí)間步內(nèi)都要對(duì)流場(chǎng)網(wǎng)格進(jìn)行更新，以確保網(wǎng)格變形后能夠繼續(xù)滿足計(jì)算要求。為了提高網(wǎng)格更新效率，本文將Diffusion Smoothing和Remeshing兩種網(wǎng)格更新方法相結(jié)合，對(duì)不同變形尺度的網(wǎng)格進(jìn)行分類處理：對(duì)小幅運(yùn)動(dòng)的網(wǎng)格采用Diffusion Smoothing[23]方法通過(guò)求解擴(kuò)散方程得到位移后網(wǎng)格節(jié)點(diǎn)的位置，但該方法未改變網(wǎng)格拓?fù)浣Y(jié)構(gòu)，在大變形區(qū)域容易產(chǎn)生負(fù)體積，此時(shí)引入Remeshing[23-24]更新方法進(jìn)行彌補(bǔ)，對(duì)超過(guò)設(shè)定的歪斜率或尺寸閾值的網(wǎng)格進(jìn)行標(biāo)記，并局部重新劃分這些網(wǎng)格單元，技術(shù)途徑如圖1所示。

Diffusion Smoothing方法通過(guò)求解以下控制方程來(lái)獲得網(wǎng)格節(jié)點(diǎn)的位置：

(2)

式中：u為網(wǎng)格運(yùn)動(dòng)速度，γ=(1/dα)為擴(kuò)散系數(shù)(α取1.5)，d為正則化后的網(wǎng)格節(jié)點(diǎn)與邊界的距離，x2為下一個(gè)時(shí)間步網(wǎng)格節(jié)點(diǎn)的位置。

圖1 本文網(wǎng)格更新技術(shù)途徑

2 研究對(duì)象及數(shù)值模型

2.1 研究對(duì)象及計(jì)算工況

本文的研究對(duì)象由三槳葉圓錐形旋翼機(jī)(槳葉角為0°)和半球形降落傘組成，其外形及幾何參數(shù)關(guān)系如圖2所示。不考慮傘衣的織物透氣性和傘繩對(duì)流場(chǎng)的影響，物傘系統(tǒng)結(jié)構(gòu)參數(shù)、來(lái)流及計(jì)算工況見(jiàn)表1、2。

圖2 旋翼機(jī)-降落傘系統(tǒng)三維模型

表1 物傘系統(tǒng)結(jié)構(gòu)參數(shù)

表2 來(lái)流條件及計(jì)算工況

2.2 數(shù)值模型

本文建立直徑6Dt、高9Dt(Dt為傘衣投影直徑)的圓柱形流場(chǎng)模型，如圖3所示。流場(chǎng)計(jì)算域分為兩個(gè)部分：旋轉(zhuǎn)流域Ωr和自由流域Ωf。旋翼位于旋轉(zhuǎn)流域，其尺寸為直徑2.2l、高0.15l的圓柱形流場(chǎng)，其他部分均處于自由流域，兩個(gè)流域交界面為Γinterface。

旋轉(zhuǎn)流域在每個(gè)時(shí)間步內(nèi)均涉及網(wǎng)格的移動(dòng)和更新，分別對(duì)旋轉(zhuǎn)流域和自由流域劃分網(wǎng)格，并將兩流域交界面Γinterface上的網(wǎng)格節(jié)點(diǎn)(如Q1、Q2和Q3)進(jìn)行連接，保證了兩個(gè)流域邊界信息傳遞的連續(xù)性，如圖4所示。計(jì)算采用了適應(yīng)性良好的三角形網(wǎng)格，為了增加計(jì)算精度和捕捉流場(chǎng)細(xì)節(jié)，對(duì)旋翼機(jī)和降落傘大梯度復(fù)雜流域的網(wǎng)格進(jìn)行加密，網(wǎng)格總數(shù)為118萬(wàn)。

圖3 計(jì)算域和物傘系統(tǒng)

圖4 物傘系統(tǒng)數(shù)值模型

2.3 網(wǎng)格適應(yīng)性驗(yàn)證

本文首先對(duì)旋翼轉(zhuǎn)速ω=0 rad/s時(shí)物傘系統(tǒng)的繞流流場(chǎng)進(jìn)行了數(shù)值計(jì)算，穩(wěn)定時(shí)傘衣的平均阻力系數(shù)為0.71，與文獻(xiàn)[25]中半球形傘阻力系數(shù)0.75的相對(duì)誤差僅為5.33%，說(shuō)明本文計(jì)算方法具有較好的準(zhǔn)確性。

表3 數(shù)值結(jié)果對(duì)比

3 物傘系統(tǒng)非定常尾流特征

3.1 流場(chǎng)分布規(guī)律

本文選取旋翼轉(zhuǎn)速ω=0 rad/s和50 rad/s兩種工況下的流場(chǎng)分布和旋渦演變進(jìn)行分析。圖5為計(jì)算穩(wěn)定(t=15 s)旋翼機(jī)-降落傘系統(tǒng)的流場(chǎng)分布，可以看出：高速轉(zhuǎn)動(dòng)的槳葉直接改變了旋翼機(jī)的尾流結(jié)構(gòu)，將其分割成兩個(gè)擾動(dòng)區(qū)域，旋渦中心由Q分裂為Q1和Q2，尾流低壓區(qū)P1也以槳葉平面P為分界面形成兩個(gè)低壓區(qū)；槳葉運(yùn)動(dòng)迫使旋翼機(jī)尾流向傘衣入口延伸，旋渦中心Q1較Q向上移動(dòng)；旋翼機(jī)尾流低壓區(qū)P1長(zhǎng)度增加，傘衣入口正壓區(qū)P2長(zhǎng)度減小且呈不對(duì)稱分布；旋翼機(jī)尾流的低動(dòng)壓區(qū)V完全移至旋翼上方，縮短了到傘衣入口的距離。

圖5 旋翼機(jī)-降落傘系統(tǒng)流場(chǎng)分布(t=15 s)

3.2 尾渦演變分析

圖6為計(jì)算穩(wěn)定階段旋翼機(jī)-降落傘系統(tǒng)的尾渦變化情況，可以看出：旋翼機(jī)和傘衣表面時(shí)刻出現(xiàn)旋渦的生成與脫離，旋翼轉(zhuǎn)動(dòng)導(dǎo)致前體尾流區(qū)的渦量呈不對(duì)稱分布，干擾了前體尾渦的形成，前體尾流中的旋渦數(shù)量明顯減少。同時(shí)，前體尾流的負(fù)渦量區(qū)在旋翼擾動(dòng)下逐漸后移，與傘衣入口的負(fù)渦量區(qū)相連通，促進(jìn)了傘衣尾渦的脫離，傘衣尾流中的旋渦數(shù)量明顯增多。

圖6 物傘系統(tǒng)旋渦變化(上：ω=0 rad/s；下：ω=50 rad/s)

圖7為物傘系統(tǒng)渦量Q的等值面分布，從圖7(a)可以看出，由于旋翼機(jī)前體和傘衣的鈍體效應(yīng)[20]，來(lái)流繞流物傘表面后形成表面渦流區(qū)VS-1和VS-2，旋翼轉(zhuǎn)動(dòng)擾亂了旋翼機(jī)表面渦流區(qū)VS-1的渦量分布，形成旋轉(zhuǎn)渦流區(qū)VS-1R和表面渦流區(qū)VS-1S。旋翼機(jī)表面渦流區(qū)VS-1繼續(xù)向后發(fā)展，生成脫落渦流區(qū)VB進(jìn)入流場(chǎng)，如圖7(b)所示，而旋翼轉(zhuǎn)動(dòng)時(shí)的脫落渦流區(qū)VB明顯小于無(wú)旋翼轉(zhuǎn)動(dòng)，旋翼擾亂了前體表面渦流分布后形成旋轉(zhuǎn)渦流區(qū)，將前體尾渦“打碎”，導(dǎo)致前體尾流中的脫落渦流區(qū)范圍變小，進(jìn)入傘衣的尾渦渦流隨之減少，傘衣表面渦流區(qū)VS-2范圍變小。

圖7 物傘系統(tǒng)渦量等值面分布(t=15 s)

4 旋翼機(jī)-降落傘氣動(dòng)干擾分析

4.1 流場(chǎng)特性參數(shù)分析

圖8為不同旋翼轉(zhuǎn)速下傘衣入口處旋渦強(qiáng)度的徑向分布，可以看出：由于分離流動(dòng)，旋渦強(qiáng)度在傘衣裙邊(R=±2.5 m)處的變化最為劇烈，此時(shí)曲線峰值較高且變化梯度較大，在傘衣內(nèi)部變化則較為平緩，隨旋翼轉(zhuǎn)速增加，旋渦強(qiáng)度曲線的峰值減小，徑向方向的旋渦強(qiáng)度逐漸減弱。

圖8 傘衣入口處旋渦強(qiáng)度的徑向分布(t=15 s)

圖9為傘衣入口處旋渦強(qiáng)度隨時(shí)間的變化，可以看出：由于旋翼機(jī)尾渦的周期性影響，傘衣入口處的旋渦強(qiáng)度隨時(shí)間發(fā)生著波動(dòng)，隨旋翼轉(zhuǎn)速增加，旋翼對(duì)前體尾渦的擾亂作用加強(qiáng)，加速了渦流區(qū)的黏性耗散，進(jìn)入流場(chǎng)的尾渦渦量減少，傘衣入口處的旋渦強(qiáng)度逐漸減弱。

圖9 傘衣入口旋渦強(qiáng)度隨時(shí)間的變化

為探究旋翼轉(zhuǎn)速對(duì)旋翼機(jī)和傘衣周圍流場(chǎng)壓力的影響情況，在對(duì)稱軸旋翼機(jī)后方0.5 m處(A點(diǎn))和傘衣前方0.5 m處(B點(diǎn))進(jìn)行流場(chǎng)壓力監(jiān)測(cè)，得到兩特征位置的壓力變化情況(如圖10所示)，從圖10中可以看出：旋翼運(yùn)動(dòng)擾亂了旋翼機(jī)的尾流負(fù)壓區(qū)，隨旋翼轉(zhuǎn)速增加，旋翼后方的流場(chǎng)壓力逐漸增加，而傘衣入口的流場(chǎng)壓力逐漸減小，當(dāng)ω≥100 rad/s時(shí)，旋翼機(jī)后方和傘衣入口前方的流場(chǎng)壓力變化逐漸延緩。

圖10 A、B點(diǎn)流場(chǎng)壓力變化

4.2 降落傘氣動(dòng)特性分析

圖11為傘衣沿子午線上的內(nèi)外壓力及壓強(qiáng)系數(shù)分布，可以看出：旋翼轉(zhuǎn)動(dòng)對(duì)傘衣面壓強(qiáng)的影響主要在傘衣內(nèi)側(cè)，隨著旋翼轉(zhuǎn)速增加，傘衣面外側(cè)壓力幾乎不變，而內(nèi)側(cè)壓力和壓強(qiáng)系數(shù)均逐漸減小，內(nèi)外壓差減小，傘衣氣動(dòng)阻力減小。

圖11 傘衣沿子午線上內(nèi)、外壓力及壓強(qiáng)系數(shù)分布

圖12 阻力系數(shù)CD隨時(shí)間的變化

圖13 平均阻力系數(shù)隨旋翼轉(zhuǎn)速的變化

5 結(jié) 論

1)本文建立的數(shù)值方法對(duì)旋翼機(jī)/降落傘非定常復(fù)合流場(chǎng)的計(jì)算具有較好的準(zhǔn)確性，同時(shí)能精確捕捉流場(chǎng)的尾渦細(xì)節(jié)變化。

2)旋翼轉(zhuǎn)動(dòng)將前體尾流分割成兩個(gè)擾動(dòng)區(qū)域，前體尾流區(qū)長(zhǎng)度增加，對(duì)降落傘的影響逐漸增強(qiáng)，導(dǎo)致傘衣入口的流場(chǎng)結(jié)構(gòu)呈不對(duì)稱分布。

3)旋翼機(jī)尾部負(fù)渦量區(qū)受旋翼影響逐漸后移，與傘衣入口處的負(fù)渦量區(qū)相連通，促進(jìn)了降落傘尾渦的脫離，傘衣尾流中的旋渦數(shù)量明顯增加。

4)旋翼轉(zhuǎn)動(dòng)擾亂了前體表面的渦流分布，在旋翼周圍形成旋轉(zhuǎn)渦流區(qū)，前體尾流脫落渦流區(qū)范圍變小，進(jìn)入傘衣的旋渦強(qiáng)度逐漸減弱.旋渦強(qiáng)度在傘衣裙邊處變化最為劇烈，在傘衣入口內(nèi)側(cè)相對(duì)平緩。

5)隨旋翼轉(zhuǎn)速增加，傘衣入口正壓區(qū)長(zhǎng)度減小，傘衣外側(cè)壓力幾乎不變，內(nèi)側(cè)壓力逐漸減小，降落傘的阻力系數(shù)隨之減小，當(dāng)ω≥100 rad/s時(shí)，平均阻力系數(shù)和壓強(qiáng)系數(shù)減小趨勢(shì)逐漸延緩。