翼型表面湍流脈動對水滴撞擊特性的影響

周 琳，曹逸韜，，郭垠昊，朱鵬飛

(1.中國航發湖南動力機械研究所，湖南 株洲 412002；2.西北工業大學動力與能源學院，西安 710072；3.陸裝駐洛陽地區航空軍事代表室，河南 新鄉 453000)

1 引言

水滴撞擊特性研究是飛機及發動機結冰、防冰數值模擬的前提和基礎。對水滴撞擊到結冰表面的質量和分布的準確模擬在很大程度上決定了結冰、防冰計算的準確性，為此準確的水滴撞擊特性計算一直是結冰研究中的一個熱點。

水滴撞擊特性計算一般有歐拉-拉格朗日法和歐拉-歐拉法兩種方法。近年來，采用歐拉-拉格朗日法的研究主要針對計算效率[1]、在非結構網格中的應用[2-3]和在三維數值計算中的應用[4-5]開展，而采用歐拉-歐拉法的研究則是針對邊界條件設置[6-7]、歐拉-歐拉法在三維數值計算中的應用[8-12]以及兩相流計算中的奇異性問題[13-14]展開。目前，兩種方法都可以用于迎風部件(如后掠機翼、尾翼、進氣唇口等)的水滴撞擊特性計算。由于已有研究對象大多為飛機的迎風面，氣流流過這些表面時湍流脈動不是很強烈，故研究中均未考慮湍流脈動對水滴運動軌跡以及水滴撞擊的影響，進而造成駐點及撞擊極限附近的數值結果與實驗結果存在一定差異。

本文以二維NACA23012翼型為研究對象，開展了湍流脈動對水滴撞擊特性影響的數值研究。采用歐拉-拉格朗日法求解空氣-水滴兩相流動，相間作用采用單向耦合；運用隨機軌道模型計算湍流脈動及空氣的瞬時速度。分析湍流脈動對水滴軌跡、撞擊特性以及水滴局部收集系數的影響，研究水滴直徑變化時湍流脈動對水滴撞擊特性的影響規律。

2 控制方程

在水滴撞擊特性數值模擬中，流場是空氣-水滴共存的兩相流場。本文研究的水滴直徑不大于40 μm 且液態水含量最大不超過3 g/m3。在此條件下，液態水的體積分數在10-6量級，遠小于兩相流中稠密與稀疏流動的劃分界限0.004[15]。因此本文研究的問題屬于稀疏兩相流，在求解中只考慮空氣相對水滴的影響，忽略水滴對空氣的影響。采用歐拉-拉格朗日法進行兩相流計算，相間采用單向耦合。對空氣場采用歐拉法描述，控制方程為雷諾時均Navier-Stokes方程，湍流模型采用標準k-ε模型；對水滴相采用拉格朗日法描述，通過分析作用于水滴上的各種作用力，并由牛頓第二定律給出：

式中：F為作用于水滴上的各種力之和，包括阻力、重力、Basset 力、壓強梯度力、Magnus 力、Saffman 力等；m為水滴質量；up為水滴速度。

因研究的水滴直徑較小，根據相關文獻[16-18]，Basset力、壓強梯度力、重力、Magnus力和Saffman力均可以忽略，故只考慮作用于水滴上的阻力。據此，方程(1)可改寫為：

式中：dp為水滴直徑；ρp為水滴密度；CD為阻力系數[18]，因水滴尺寸較小，可假設為剛性球體，采用球形顆粒的阻力系數；ua為當地空氣瞬時速度；Re為相對雷諾數，定義為為空氣密度。

3 隨機軌道模型

采用隨機軌道模型[19]模擬空氣湍流脈動。方程(2)中當地空氣瞬時速度可寫為時均速度與脈動速度之和，即ua=+u′。基于Bussinise 假設的各向同性湍流僅可以得到流場的時均速度，脈動速度需要由隨機軌道模型產生。隨機軌道模型假設氣體脈動速度服從高斯分布，由如下公式計算：

式中：ζ為均值為0、方差為1的正態分布隨機數；k為計算的湍動能。

水滴的拉格朗日積分時間TL表示為：

式中：CL為常數，對k-ε湍流模型以及改進形式的湍流模型，CL約為0.15。

渦的特征時間τe表示為：

式中：r為0到1均勻分布隨機數，也即渦的特征時間是拉格朗日積分時間的一個隨機變量。

水滴穿過渦的時間τcross為：

式中：τ為水滴的弛豫時間；Le為渦的特征尺度；|u-up|為空氣與水滴的相對速度。

水滴與渦的作用時間τint取為渦的特征時間與水滴穿過渦的時間中的最小值，即：

隨機軌道模型假設水滴連續地與一系列渦作用，每個渦具有正態分布的脈動速度及特征時間。當渦的特征時間結束或水滴穿過渦時，水滴將會與下一個渦作用。即當達到τint這一時間時，由式(3)獲得新的瞬時速度。空氣采用瞬時速度后，由于脈動速度是隨機量，水滴的軌跡勢必會表現出無序性，這給水滴收集系數的計算帶來了困難。下文從水滴收集系數定義出發，推導空氣采用瞬時速度時水滴收集系數的計算公式。

4 水滴收集系數的計算

水滴收集系數是研究水滴撞擊特性的一個重要參數，有局部收集系數與總收集系數兩種。為反映部件表面撞擊水量的分布，一般采用局部收集系數，其定義為物面某局部區域實際所收集的水量與該物面可能收集的最大水量之比[20-21]：

式中：Wβ為局部的實際水收集量，表示撞擊到結構表面微元面積Am內水滴所攜帶的質量流量m˙i之和；Wβmax為局部可能收集的最大水量，可表示為Wβmax=LWC·U∞·Am(LWC為大氣中的液態水含量，U∞為來流速度)。因此，局部收集系數計算公式可寫為：

5 計算分析

5.1 計算模型

計算模型如圖1 所示。NACA23012 翼型弦長0.914 4 m，遠場邊界設置在約3倍弦長處，翼型表面為壁面邊界。

圖1 NACA23012計算模型Fig.1 Calculation model of NACA23012

5.2 計算工況

NACA23012 翼型的計算工況如表1 所示。表中，MVD為水滴等效平均直徑。

表1 NACA23012翼型計算工況Table 1 NACA23012 aerofoil computational condition

5.3 計算結果分析

5.3.1 空氣場計算結果分析

空氣-水滴兩相流數值計算采用單向耦合，首先計算空氣場，然后基于空氣場計算結果進行兩相流場內水滴軌跡和撞擊特性的計算。為考察空氣場計算結果的準確性，選取翼型表面壓力系數的計算結果與實驗結果[22]進行對比。壓力系數定義為：

式中：p為翼型表面壓力；p∞為來流壓力；ρ∞為來流密度；V∞為來流速度。

翼型表面壓力系數計算結果與實驗結果的對比如圖2 所示。可看出，壓力面上計算結果與實驗結果吻合非常好；吸力面上氣流先加速后減速，使得過駐點后的吸力面上有個低壓區，此處壓力系數與實驗值稍有差別。總體看，計算結果與實驗結果吻合良好，空氣場計算結果準確可靠。

5.3.2 湍流脈動對水滴軌跡的影響分析

圖3給出了空氣場x方向和y方向的時均速度分布計算結果。可以看出，x方向速度相對較高，除駐點外都在10 m/s量級；駐點后吸力面的速度增加，最高達到了100 m/s。而y方向速度總體相對較低，大部分區域在1 m/s量級；駐點后y方向速度迅速增大，在貼近駐點的吸力面和壓力面的一小部分區域，速度也達到了10 m/s量級。據此，空氣為時均速度時，大部分區域的水滴速度主要受空氣場x方向速度影響，駐點后一小部分區域的水滴速度將受到空氣場x方向和y方向速度的共同作用。

圖3 空氣場x方向和y方向的時均速度分布Fig.3 Time-averaged velocity distribution along x and y direction

圖4給出了空氣為時均速度和瞬時速度時水滴的運動軌跡。可以看出，空氣為時均速度時，水滴軌跡非常整齊且擾流過翼型后有明顯的遮蔽區。在駐點附近，由于y方向速度小，水滴主要受x方向速度影響，水滴撞擊在翼型表面。離開駐點向后，由于y方向速度增加到與x方向速度相同量級，水滴軌跡發生偏轉，繞過翼型前緣；之后y方向速度又迅速降低，x方向速度增加，水滴又主要受x方向速度影響，形成遮蔽區。由于駐點后吸力面y方向速度比壓力面的高，因此吸力面形成的遮蔽區更大。

圖4 空氣為時均速度和瞬時速度時水滴的運動軌跡Fig.4 Trajectory of the droplets at time-averaged and transient speed

空氣為瞬時速度時，水滴軌跡無序且遮蔽區不明顯。這是因為水滴運動受到了湍流脈動的影響，而脈動速度是隨機的，使得同一截面上水滴的運動速度大小和方向不同所致。

5.3.3 水滴直徑對水滴收集系數的影響分析

圖5 給出了水滴等效平均直徑為20，30，40 μm時湍流脈動對水滴收集系數的影響。可看出，當MVD為20 μm時，駐點附近湍流脈動的影響較明顯；MVD為40 μm 時，駐點附近湍流脈動的影響已變得非常小。說明隨著水滴粒徑的增大，湍流脈動的影響逐漸減小。這是因為粒徑增大，水滴保持自身慣性的能力增強，受空氣速度的影響減小，即使考慮了空氣脈動速度，結果變化也非常小。另外，撞擊極限附近的差異雖然也隨著粒徑的增大而逐漸減小，但是依然較明顯。其原因是在撞擊極限附近，受空氣脈動速度的影響，水滴逐漸撞擊在翼型表面，局部收集系數曲線緩慢變為零；而不考慮脈動速度時，水滴達到最遠撞擊點之后，便不再有水滴撞擊到翼型表面，局部收集系數曲線在撞擊極限位置會突然變為零。據此，雖然水滴直徑增大，湍流脈動的影響減小，但撞擊極限附近的影響依然較明顯。

圖5 湍流脈動對水滴局部收集系數的影響Fig.5 Influence of fluctuation for local collected coefficient

5.3.4 水滴直徑Langmiur D 分布時脈動速度的影響分析

實際飛行時，大氣中的水滴直徑分布是不均勻的，一般認為符合Langmiur D 分布。選取最常用的7 種直徑的Langmiur D 分布，如表2 所示，計算分析了湍流脈動對水滴局部收集系數的影響，結果見圖6。可以看出，MVD=20 μm 時，在駐點附近，受湍流脈動影響，水滴收集系數略有降低，而在撞擊極限附近，脈動速度影響不明顯。其原因是不同粒徑水滴的撞擊極限不同，因此在撞擊極限附近水滴收集系數逐漸降低，而水滴受脈動速度影響時，水滴仍然會依次沉降在翼型表面，但受撞擊極限附近y方向脈動速度的影響，水滴局部收集系數略偏高。當MVD由30 μm 增加至40 μm 時，湍流脈動影響逐漸減小。因此，采用7種粒徑的Langmiur D分布時，在一定程度上能夠降低忽略湍流脈動帶來的影響。

圖6 水滴直徑為Langmiur D分布時水收集系數Fig.6 Local collection efficiency comparison(Langmiur D distribution)

表2 7種水滴直徑的Langmiur D分布Table 2 Langmiur D distribution of seven kinds of droplet diameters

6 結論

研究了翼型表面湍流脈動對水滴撞擊特性以及水滴收集系數的影響，得到如下結論：

(1) 湍流脈動使水滴的運動軌跡變得無序且遮蔽區不明顯。

(2) 湍流脈動對水滴撞擊特性的影響主要體現在駐點附近和撞擊極限附近，其他區域脈動影響不明顯。

(3) 隨著水滴直徑的增加，湍流脈動的影響逐漸減小，但撞擊極限附近的影響依然較明顯。

(4) 采用7 種粒徑的Langmiur D 分布時，湍流脈動的影響非常小。