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

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

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

1 引言

水滴撞擊特性研究是飛機(jī)及發(fā)動機(jī)結(jié)冰、防冰數(shù)值模擬的前提和基礎(chǔ)。對水滴撞擊到結(jié)冰表面的質(zhì)量和分布的準(zhǔn)確模擬在很大程度上決定了結(jié)冰、防冰計算的準(zhǔn)確性，為此準(zhǔn)確的水滴撞擊特性計算一直是結(jié)冰研究中的一個熱點(diǎn)。

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

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

2 控制方程

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

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

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

式中：dp為水滴直徑；ρp為水滴密度；CD為阻力系數(shù)[18]，因水滴尺寸較小，可假設(shè)為剛性球體，采用球形顆粒的阻力系數(shù)；ua為當(dāng)?shù)乜諝馑矔r速度；Re為相對雷諾數(shù)，定義為為空氣密度。

3 隨機(jī)軌道模型

采用隨機(jī)軌道模型[19]模擬空氣湍流脈動。方程(2)中當(dāng)?shù)乜諝馑矔r速度可寫為時均速度與脈動速度之和，即ua=+u′?；贐ussinise 假設(shè)的各向同性湍流僅可以得到流場的時均速度，脈動速度需要由隨機(jī)軌道模型產(chǎn)生。隨機(jī)軌道模型假設(shè)氣體脈動速度服從高斯分布，由如下公式計算：

式中：ζ為均值為0、方差為1的正態(tài)分布隨機(jī)數(shù)；k為計算的湍動能。

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

式中：CL為常數(shù)，對k-ε湍流模型以及改進(jìn)形式的湍流模型，CL約為0.15。

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

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

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

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

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

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

4 水滴收集系數(shù)的計算

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

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

5 計算分析

5.1 計算模型

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

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

5.2 計算工況

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

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

5.3 計算結(jié)果分析

5.3.1 空氣場計算結(jié)果分析

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

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

翼型表面壓力系數(shù)計算結(jié)果與實(shí)驗(yàn)結(jié)果的對比如圖2 所示?？煽闯?，壓力面上計算結(jié)果與實(shí)驗(yàn)結(jié)果吻合非常好；吸力面上氣流先加速后減速，使得過駐點(diǎn)后的吸力面上有個低壓區(qū)，此處壓力系數(shù)與實(shí)驗(yàn)值稍有差別?？傮w看，計算結(jié)果與實(shí)驗(yàn)結(jié)果吻合良好，空氣場計算結(jié)果準(zhǔn)確可靠。

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

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

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

圖4給出了空氣為時均速度和瞬時速度時水滴的運(yùn)動軌跡?？梢钥闯觯諝鉃闀r均速度時，水滴軌跡非常整齊且擾流過翼型后有明顯的遮蔽區(qū)。在駐點(diǎn)附近，由于y方向速度小，水滴主要受x方向速度影響，水滴撞擊在翼型表面。離開駐點(diǎn)向后，由于y方向速度增加到與x方向速度相同量級，水滴軌跡發(fā)生偏轉(zhuǎn)，繞過翼型前緣；之后y方向速度又迅速降低，x方向速度增加，水滴又主要受x方向速度影響，形成遮蔽區(qū)。由于駐點(diǎn)后吸力面y方向速度比壓力面的高，因此吸力面形成的遮蔽區(qū)更大。

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

空氣為瞬時速度時，水滴軌跡無序且遮蔽區(qū)不明顯。這是因?yàn)樗芜\(yùn)動受到了湍流脈動的影響，而脈動速度是隨機(jī)的，使得同一截面上水滴的運(yùn)動速度大小和方向不同所致。

5.3.3 水滴直徑對水滴收集系數(shù)的影響分析

圖5 給出了水滴等效平均直徑為20，30，40 μm時湍流脈動對水滴收集系數(shù)的影響?？煽闯?，當(dāng)MVD為20 μm時，駐點(diǎn)附近湍流脈動的影響較明顯；MVD為40 μm 時，駐點(diǎn)附近湍流脈動的影響已變得非常小。說明隨著水滴粒徑的增大，湍流脈動的影響逐漸減小。這是因?yàn)榱皆龃?，水滴保持自身慣性的能力增強(qiáng)，受空氣速度的影響減小，即使考慮了空氣脈動速度，結(jié)果變化也非常小。另外，撞擊極限附近的差異雖然也隨著粒徑的增大而逐漸減小，但是依然較明顯。其原因是在撞擊極限附近，受空氣脈動速度的影響，水滴逐漸撞擊在翼型表面，局部收集系數(shù)曲線緩慢變?yōu)榱?；而不考慮脈動速度時，水滴達(dá)到最遠(yuǎn)撞擊點(diǎn)之后，便不再有水滴撞擊到翼型表面，局部收集系數(shù)曲線在撞擊極限位置會突然變?yōu)榱?。?jù)此，雖然水滴直徑增大，湍流脈動的影響減小，但撞擊極限附近的影響依然較明顯。

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

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

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

圖6 水滴直徑為Langmiur D分布時水收集系數(shù)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 結(jié)論

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

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

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

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

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