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

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

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

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

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

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

1 計算方法

1.1 數值方法及湍流模型

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

為獲得精細的流場結構、有效求解固壁繞流問題，本文計算采用了標準壁面函數法(standard wall functions)，湍流模型選擇Realizablek-ε兩方程模型，其控制方程為

(1)

1.2 網格更新模型

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

Diffusion Smoothing方法通過求解以下控制方程來獲得網格節(jié)點的位置：

(2)

式中：u為網格運動速度，γ=(1/dα)為擴散系數(α取1.5)，d為正則化后的網格節(jié)點與邊界的距離，x2為下一個時間步網格節(jié)點的位置。

圖1 本文網格更新技術途徑

2 研究對象及數值模型

2.1 研究對象及計算工況

本文的研究對象由三槳葉圓錐形旋翼機(槳葉角為0°)和半球形降落傘組成，其外形及幾何參數關系如圖2所示。不考慮傘衣的織物透氣性和傘繩對流場的影響，物傘系統(tǒng)結構參數、來流及計算工況見表1、2。

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

表1 物傘系統(tǒng)結構參數

表2 來流條件及計算工況

2.2 數值模型

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

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

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

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

2.3 網格適應性驗證

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

表3 數值結果對比

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

3.1 流場分布規(guī)律

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

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

3.2 尾渦演變分析

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

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

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

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

4 旋翼機-降落傘氣動干擾分析

4.1 流場特性參數分析

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

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

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

圖9 傘衣入口旋渦強度隨時間的變化

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

圖10 A、B點流場壓力變化

4.2 降落傘氣動特性分析

圖11為傘衣沿子午線上的內外壓力及壓強系數分布，可以看出：旋翼轉動對傘衣面壓強的影響主要在傘衣內側，隨著旋翼轉速增加，傘衣面外側壓力幾乎不變，而內側壓力和壓強系數均逐漸減小，內外壓差減小，傘衣氣動阻力減小。

圖11 傘衣沿子午線上內、外壓力及壓強系數分布

圖12 阻力系數CD隨時間的變化

圖13 平均阻力系數隨旋翼轉速的變化

5 結 論

1)本文建立的數值方法對旋翼機/降落傘非定常復合流場的計算具有較好的準確性，同時能精確捕捉流場的尾渦細節(jié)變化。

2)旋翼轉動將前體尾流分割成兩個擾動區(qū)域，前體尾流區(qū)長度增加，對降落傘的影響逐漸增強，導致傘衣入口的流場結構呈不對稱分布。

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

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

5)隨旋翼轉速增加，傘衣入口正壓區(qū)長度減小，傘衣外側壓力幾乎不變，內側壓力逐漸減小，降落傘的阻力系數隨之減小，當ω≥100 rad/s時，平均阻力系數和壓強系數減小趨勢逐漸延緩。