液態水含量對防冰表面水膜流動換熱的影響

鄭　梅，朱劍鋆，董　威

（上海交通大學機械與動力工程學院，上海200240）

液態水含量對防冰表面水膜流動換熱的影響

鄭梅，朱劍鋆，董威

（上海交通大學機械與動力工程學院，上海200240）

為了研究液態水含量對防冰表面水膜流動換熱的影響，基于機翼防冰表面水膜及空氣相互作用機理，并考慮水膜表面傳熱傳質過程，建立了水膜與空氣的流動換熱模型，得到溢流水膜及空氣邊界層流動換熱的積分控制方程，通過對比文獻試驗結果驗證了模型的準確性。在此基礎上，比較了不同液態水含量條件下防冰表面水膜厚度及主要熱流量的分布情況。結果表明：液態水含量對水膜沿表面厚度分布有明顯影響，而對換熱過程中各項熱流的影響主要集中在水滴撞擊區域，加熱熱流與散熱熱流隨液態水含量的增加呈現相反的變化趨勢。

水膜；防冰；流動換熱；液態水含量；飛機

0　引言

結冰是影響飛行安全的一大隱患。當飛機穿越由亞穩態過冷水滴組成的云團時，部件迎風表面就容易結冰［1］。氣動面結冰會破壞飛機的氣動外形，影響其穩定性和操縱性；發動機進口部件結冰會改變發動機進氣條件，造成進氣畸變和推力減小；近些年來，甚至還發生過由于結冰而直接引發的墜機事故［2-3］。

熱氣防冰系統是目前用于防止飛機結冰而使用最為廣泛且技術最為成熟的裝置。該系統從發動機壓氣機中引出帶有一定壓力的熱空氣，對防冰部件表面進行加熱，使得撞擊表面的過冷水滴受熱蒸發，從而避免結冰［4］。熱氣防冰系統能夠處于2種工作狀態，完全蒸發系統能夠完全蒸發部件表面的撞擊水，而不完全蒸發系統則允許部件表面溢流水存在［5］。對于不完全蒸發熱氣防冰系統，冷水滴撞擊防冰表面后，未完全蒸發的液態水在氣動力作用下沿表面向后流動，形成溢流水。溢流水膜的流動形態將直接決定可能的結冰區域和所需的防冰區域，同時對防冰部件所需的防冰熱流分布造成直接影響。因此，國內外對此開展了許多研究［6-11］。Messinger模型［12］是結冰與防冰表面流動換熱分析中最為經典的熱平衡模型。結合質量守恒和能量守恒定律，模型可以求解穩態情況下表面的溫度分布。此后，形成功能相對完善的結冰計算軟件，如：LEWICE［13］、ANTICE［14］、ONERA［15］、CANICE［16］和FENSAP-ICE［17］等。很明顯，表面水膜的流動形態對防冰表面的流動換熱過程影響很大，但現有的研究工作大多圍繞結冰表面展開，而對防冰表面溢流水膜流動形態及厚度分布的預測較少。

本文主要關注水膜在防冰表面的流動問題。建立了防冰表面水膜流動換熱的數學模型，并獲得了水膜沿防冰表面流向的厚度以及換熱過程中各項熱流的分布情況，以及不同液態水含量對防冰表面流動換熱的影響。

1　水膜流動的物理描述

防冰部件表面的水膜流動現象是伴隨水滴撞擊產生的。過冷水滴到達防冰部件表面，會在部件表面鋪展形成1層水膜。由于受到防冰系統的加熱，一部分水會在流動過程中蒸發，而剩余部分會在氣動力的作用下沿部件表面向后流動形成溢流水。相關研究表明，在實際飛行過程中部件表面形成的水膜是非常薄的，且在流動過程中其厚度隨著蒸發等作用不斷變化。

2維翼型表面水膜流動如圖1所示。在水膜流動過程中存在質量、動量及能量的交換。從圖中可見，其前方流入的水膜質量與撞擊質量都會導致控制體積內水膜質量的增大，而水膜的蒸發又導致控制體積質量的減小。水膜沿表面的流動受空氣剪切應力和壓力差驅動，并和壁面黏滯力共同對水膜沿壁面法向速度分布產生影響。此外，水膜還受到部件表面加熱和空氣對流換熱的作用，使其溫度分布改變。為了描述這些物理過程，將控制體積分為水膜層和空氣層，分別建立相應的控制方程，并給出對應的邊界條件。

圖1　防冰機翼表面水膜流動

2　數學模型

在翼型表面建立曲線坐標系（圖1），沿表面流向方向為s軸，法向方向為y軸。分別取翼型防冰表面上方水膜和空氣邊界層的控制體積為研究對象。δ1為水膜厚度，δ2為空氣邊界層厚度。控制體積沿s方向長度取為ds，展向方向取單位寬度。為了便于數學描述，采用如下假設：

（1）空氣邊界層和水膜流動均為2維不可壓縮層流且為定常層流流動。

（2）由于水膜非常薄，同時空氣中液態水含量很小，因此僅考慮空氣對水膜的單向耦合作用，同時由于水膜重力相對于空氣剪切力十分小，因此忽略重力影響，僅考慮空氣剪切應力與壓力差對水膜流動的驅動效應。

（3）僅考慮防冰系統運行足夠長時間后表面水膜的流動情況，即穩態情況下的水膜流動問題，且不考慮防冰失敗形成積冰的情況。

（4）控制方程中的所有物性參數均為常數。

2.1水膜流動控制方程

2.1.1連續性方程

取控制體積內的水膜為研究對象，進入控制體的質量包括前一控制體溢流到當前控制體的質量m.in和由于水滴碰撞進入控制體積的質量m.imp。流出控制體積的質量包括當前控制體積溢流到下一控制體積的質量m.out和控制體積表面由于蒸發損失的質量m.evap。根據質量守恒原理，水膜流動的連續性方程應滿足

控制體積內撞擊水量m.imp為

式中：LWC為液態水含量；β為局部水收集系數；u∞為無限遠處的來流速度。

控制體積內的蒸發水量m.evap［2］為

式中：hc為對流換熱系數；Cp，air為空氣定壓比熱容；Pv，w為溢流水表面飽和蒸汽壓；Pv，e為自由來流邊界層飽和蒸汽壓。

在一定條件下，飽和蒸汽壓可以表示為關于溫度的函數［8］

考慮水膜沿y方向的速度分布，控制體積進、出口的質量流量為

因此，水膜的連續性方程可以簡化為

2.1.2動量方程

防冰部件表面的水膜在流動方向上主要受到空氣剪切應力、沿流動方向的壓力差以及壁面黏滯力的作用。水膜沿表面流動方向受到的合力Fs-water［18］為

式中：p為機翼表面靜壓力；ε為壓力梯度增量系數，為極小值；τ0為壁面黏滯應力；τair為氣-液界面空氣對水膜的剪切應力，由τ=μ（du/dy）確定。將上式展開整理后并忽略2階小量可得

考慮控制體積各表面的動量通量。由此，總的動量通量為

式中：u∞，s為水滴撞擊表面時沿s方向的速度。

由動量定理可知，對于控制體積，動量變化率和儲存動量增量之和等于作用于控制體積上的合力［19］。運用定常流動的假設，水膜控制體積的動量方程為

式中：uδ1為氣-液界面水膜速度。

對于不可壓縮層流，速度分布可以用1個多項式表示，由邊界條件可以確定其各項系數。對于水膜流動，由于其厚度較薄，可以假設邊界層速度為2次方多項式分布。在沿流向方向的某一給定位置，水膜邊界層速度分布是關于y的拋物線函數。運用壁面無滑移假設，以及水膜在氣-液界面處的邊界條件，水膜內速度分布為

其中：uδ1、δ1、τair均為關于s的函數。

2.1.3能量方程

從圖1中可見，進入控制體積的能量包括水滴撞擊防冰表面帶入的動能Qimp、水膜由前一控制體積流入當前控制體積帶入的能量Qin和由防冰表面導入的能量Qw。離開控制體積的能量包括流出當前控制體積帶走的能量Qout、由于蒸發損失的能量Qevap、控制體積表面與空氣對流換熱損失的能量Qc和加熱碰撞水滴損失的能量Qd。根據熱力學第1定律，水膜控制體積的能量方程為

考慮換熱面積1·ds，單位面積控制體積進、出口熱流密度之差為

式中：Cp，water為液態水的定壓比熱容；Twater為水膜內的溫度分布；T∞為自由來流溫度。

單位面積各項熱流密度為

由于水膜厚度較薄，其邊界層內的溫度分布可以假設為線性分布，代入邊界條件求解后可以表示為

由此，水膜的能量積分方程可以簡化為

式中：λwater為液態水的導熱系數；Le為液態水的蒸發潛熱；Tw-s為氣-液界面水膜表面溫度；Td為碰撞水滴的溫度，可以視為與來流溫度相等。

2.2空氣邊界層流動控制方程

取空氣邊界層的控制體積為研究對象，與水膜控制體積的分析方法相類似，在s方向上，空氣控制體積的總動量通量為

其中：ue為空氣邊界層靠近主流位置的速度；ρair為自由來流空氣密度；uair為空氣邊界層的速度分布。

沿s方向，空氣控制體積所受合力為

其中：τ'air為氣-液界面水膜對空氣的黏滯力，與τair互為作用力和反作用力。簡化可得

根據動量定理可以寫出完整的動量積分方程，并簡化可得

同樣假設空氣邊界層速度為1個2次多項式分布，代入邊界條件后可求得

氣-液界面黏滯力由邊界層速度梯度決定

3　求解過程

表面水膜的流動需要通過求解水滴撞擊特性以獲得水滴撞擊量，而空氣流場求解是獲得水滴撞擊特性的前提。2維翼型的外流場數據通過商業軟件ANSYS Fluent獲得，而局部水收集系數采用歐拉法，基于Fluent用戶自定義標量方程求解。

圖2　計算流程

水膜流動控制方程采用Matlab編程求解。在每1個迭代步中，依次求解空氣邊界層動量方程、水膜連續性方程、動量方程和能量方程。在每1迭代步結束時將求解得到的水膜厚度、空氣邊界層厚度、水膜上表面溫度和厚度與迭代開始前的值進行對比，以判定迭代是否收斂。若迭代收斂則輸出計算結果，反之則將計算結果作為迭代初始值進行新1輪迭代直至收斂。整個計算過程如圖2所示。

4　算例驗證

采用文獻［20］中的試驗作為算例，對建立的2維水膜流動換熱模型進行驗證。文獻［20］中給出了利用DIP技術獲得的NACA 0012翼型表面液態水的流動形態和水膜厚度。因此，2維模型的有效性主要通過對比文獻中穩態情況下的連續水膜厚度測量結果來驗證。

對比試驗中的模型為截面NACA 0012翼型，其弦長為0.101 m，攻角為0°，試驗均在室溫條件下進行。試驗過程中液態水含量（LWC）以及平均水滴直徑（DMV）等參數見表1。

表1　模型計算狀態參數

在進行水膜流動計算前，先采用ANSYS Flu-ent 14.5軟件計算翼型表面空氣流場并迭代至收斂；在此基礎上，再采用歐拉法［21］計算水滴撞擊特性并獲得翼型表面的局部水收集系數。翼型表面的局部水收集系數分布如圖3所示。

在獲得空氣流場及水滴撞擊特性計算結果的基礎上，求解2維水膜流動換熱模型，可得到水膜在翼型表面各位置處的厚度分布。

圖3　局部水收集系數分布

翼型表面水膜厚度的計算結果與文獻中試驗測量結果的對比如圖4所示。從圖中可見，當表面液態水以穩態連續水膜形式存在時，即水膜未破裂之前，本文模型的計算結果與試驗結果在總體趨勢上吻合良好，水膜厚度均沿表面流向不斷增加，僅在駐點區域附近計算結果與試驗結果有所偏差。但從整體趨勢和水膜厚度的數量級的吻合度方面而言，提出的2維水膜流動換熱模型在一定程度上能較為準確地反映翼型表面的水膜流動情況。

圖4　水膜厚度分布

5　液態水含量對水膜流動換熱的影響

在不同氣象參數條件下，防冰部件表面的水膜流動會呈現不同的特性，從而影響表面換熱過程。液態水含量作為主要的氣象參數之一，對防冰表面水膜流動換熱過程會產生影響。

不同液態水含量下算例的計算狀態參數見表2。

表2　不同液態水含量下的計算狀態參數

圖5　不同液態水含量下水膜厚度的分布

5.1水膜厚度分布

不同液態水含量條件下水膜厚度沿翼型表面的分布如圖5所示。從圖中可見，機翼表面水膜厚度的數量級均在微米級，沿翼型表面的分布呈現先增后減的趨勢。在水滴撞擊區域，部件表面水膜質量的改變主要受到撞擊與蒸發的綜合作用。由于在駐點附近撞擊質量相比蒸發質量要大得多，因此水膜沿表面向后流動過程中先不斷變厚；當水膜經過一定流動后，其蒸發質量會超過撞擊量，此時表面水膜達到最大厚度，此后水膜的厚度會在流動過程中不斷減小，直至被完全蒸發。隨著液態水含量的增加，水膜沿翼型表面的鋪展區域及最大厚度均有所增加。由于液態水含量的增加會使翼型表面水滴撞擊質量增加，不同液態水含量條件下機翼表面水滴撞擊質量的分布如圖6所示；而在相同防冰狀態下由于水膜厚度很薄，水膜表面與氣流溫差幾乎保持不變，使表面液態水蒸發量近似相等。因此，根據水膜質量守恒方程可知，表面溢流水增加，即表現為水膜厚度增加。

圖6　不同液態水含量下機翼表面水滴撞擊質量的分布

圖7　不同液態水含量下機翼表面導熱熱流的分布

圖8　不同液態水含量下水膜表面蒸發散熱熱流的分布

圖9　不同液態水含量下水膜表面對流換熱熱流的分布

5.2主要熱流項分布

不同液態水含量下水膜換熱過程中主要熱流項的分布情況如圖7～9所示。從各圖的局部放大圖中均可見，液態水含量對于水膜流動過程中換熱熱流的影響主要集中在翼型表面的水滴撞擊區域。從圖7中可見，機翼壁面導熱熱流隨液態水含量的增加而增加。而從圖8、9中可見，隨著液態水含量的增加，在水滴撞擊區域內，水膜表面由于蒸發和對流換熱作用所損失的熱量反而略有減小。因為在水滴撞擊區域內，水滴撞擊質量隨液態水含量的增加而增加，使得溢流水量增加，導致水膜控制體積進、出口焓值之差增大，即qout-qin增加。因此，由水膜能量平衡方程可知，當液態水含量增加時，即使熱損失有所減少，但壁面導熱熱流仍有所增加。

6　結論

以防冰表面水膜的流動換熱作為研究對象，基于水膜及空氣的相互作用機理，建立了水膜流動的連續性方程、動量方程和能量方程以及空氣邊界層動量方程。通過對文獻算例的求解，本文模型的有效性得到驗證。通過改變算例的液態水含量，研究了液態水含量對表面水膜厚度及各項熱流分布的影響。主要結論如下：

（1）在駐點附近的水滴撞擊區域，水滴撞擊質量對水膜厚度影響較大；到達撞擊極限后，水膜表面蒸發作用開始占主導地位，由此，水膜厚度沿機翼表面的分布呈現先增后減的趨勢；

（2）液態水含量的增加使得表面撞擊區域內的收集水量增加，表現為水膜在表面的鋪展面積及最大厚度均隨之增加；

（3）液態水含量對壁面導熱熱流、蒸發及對流換熱作用產生的熱損失的影響主要集中在水滴撞擊區域，且加熱熱流和散熱熱流隨液態水含量的增加呈現相反的變化趨勢。

［1］Melody J W.In-flight characterization of aircraft icing［D］.Illinois：University of Illinois at Urbana-Champaign，2004.

［2］裘燮綱，韓鳳華.飛機防冰系統［M］.北京：航空工業出版社，1985：44-58. QIU Xiegang，HAN Fenghua.Aircraft anti-icing system［M］.Beijing：Aviation Industry Press，1985：44-58.（in Chinese）

［3］Addy H E.Ice accretions and icing effects for modern airfoils［R］. NASA-TP-2000-210031.

［4］常士楠.大型飛機的防/除冰問題［C］//大型飛機飛機關鍵技術高層論壇暨中國航空學會學術年會.北京：中國航空學會，2007：1-7. CHANG Shinan.Large aircraft anti-icing and deicing［C］//China Aviation Association Annual Meeting.Beijing：Chinese Society of Aeronautics and Astronautics，2007：1-7.（in Chinese）

［5］周莉，徐浩軍，龔勝科，等.飛機結冰特性及防除冰技術研究［J］.中國安全科學學報，2010，20（6）：105-110. ZHOU Li，XU Haojun，GONG Shengke，et al.Research of aircraft icing characteristics and antiicing and deicing technology［J］.China Safety Science Journal，2010，20（6）：105-110.（in Chinese）

［6］Silva G A L，Silvares O M，Zerbini G J.Numerical simulation of airfoil thermal anti-ice operation，part 1：mathematical modelling［J］.Journal of Aircraft，2007，44（2）：627-633.

［7］Silva G A L，Silvares O M，Zerbini G J.Numerical simulation of airfoil thermal anti-ice operation，part 2：implementation and results［J］.Journal of Aircraft，2007，44（2）：634-641.

［8］Al-Khalil K M.Numerical simulation of an aircraft anti-icing system incorporating a rivulet model for the runback water［D］.Toledo：The University of Toledo，1991.

［9］朱劍鋆，董葳.防冰部件表面流動換熱與溫度計算分析［J］.燃氣渦輪試驗與研究，2011，24（1）：15-18. ZHU Jianjun，DONG Wei.Heat transfer analyis and temperature prediction of the anti-icing airfoil surface［J］.Gas Turbine Experiment and Research，2011，24（1）：15-18.（in Chinese）

［10］易賢，朱國林.考慮傳質傳熱效應的翼型積冰計算［J］.空氣動力學學報，2004，22（4）：490-493. YI Xian，ZHU Guolin.Computation of glaze ice accretion on airfoil［J］. ACTA Aerodynamica Sinica，2004，22（4）：490-493.（in Chinese）

［11］Zheng M，Dong W，Lei G L，et al.Study of the flow and heat transfer of water film on hot air anti-icing airfoil surface［R］.IHTC15-8671.

［12］Messinger B L.Equilibrium temperature of an unheated icing surface as a function of air-speed［J］.Journal of the Aeronautical Sciences，1953，20（1）：29-42.

［13］Wright W B，Al-Khalil K M，Miller D R.Validation of NASA thermal ice protection computer codes，part 2-the validation of LEWICE/thermal［R］.AIAA-97-0050.

［14］Al-Khalil K M，Wright W B，Miller D R.Validation of NASA thermal ice protection computer codes，part 3-the validation of ANTICE［R］. AIAA-97-0051.

［15］Hedde T，Guffond D.ONERA three-dimensional icing model［J］.AIAA Journal，1995，33（6）：1038-1045.

［16］Morency F，Brahimi M，Tezok F，et al.Hotairanti-icing system modelization in the ice prediction code CANICE［R］.AIAA-98-0192.

［17］Honsek R，Habashi W G.FENSAP-ICE：eulerian modeling of droplet impingement in the SLD regime of aircraft icing［R］.AIAA-2006-465.

［18］Hughes W F，Brighton J A.Schanm’s outline of theory and problems of fluid dynamics［M］.New York：McGraw-Hill，1967：80-102.

［19］Kays W M，Crawford M E，Weigand B.Convective heat and mass transfer［M］.New York：McGraw-Hill Series in Mechanical Engineering，2005：427-452.

［20］Zhang K，Johnson B，Rothmayer A P，et al.An experimental investigation on wind-driven rivulet/film flows over a NACA0012 airfoil by using digital image projection technique［R］.AIAA-2014-0741.

［21］Dong W，Zhu J J，Chen Y，et al.Numerical analysis of oil anti-icing system performance of turboshaft engine inlet strut［R］.ASME 2014-GT-26564.

（編輯：張寶玲）

Influence of Liquid Water Content on Flow and Heat Transfer of Water Film on Anti-Icing Surfaces

ZHENG Mei，ZHU Jian-jun，DONG Wei
（School of Mechanical Engineering，Shanghai Jiao Tong University，Shanghai 200240，China）

In order to study the effect of liquid water content on flow and heat transfer of water film on anti-icing surface，a mathematical model was developed based on the interaction mechanism of the water film and air boundary layer and considering the mass and heat transfer of the water film.Governing equations of water film and air flow in the model were both in integral form.The validation of the model was carried out by the comparison between the computation of the model and the measurement data from the reference.On this basis，the film thickness and main heat flux distribution on the anti-icing surface under different liquid water contents were compared. Results show that the thickness distribution of the water film is affected obviously by liquid water content while the influence for the heat transfer mainly concentrates on the droplets impingement region.In addition，the opposite trend can be observed in the distributions of the heating heat fluxes and heat losses with liquid water content increasing.

water film；anti-icing；flow and heat transfer；liquid water content；aircraft

V 211.3

A

10.13477/j.cnki.aeroengine.2016.01.012

2014-11-19基金項目：國家自然科學基金（51076103、11272212）、國家重點基礎研究發展計劃（2015CB755800）資助

鄭梅（1991），女，在讀博士研究生，研究方向為飛機防冰；E-mail：may-zheng@sjtu.edu.cn。

引用格式：鄭梅，朱劍鋆，董威.液態水含量對防冰表面水膜流動換熱的影響［J］.航空發動機，2016，42（1）：59-64.ZHENG Mei ZHU Jianjun DONG Wei. Influence of liquid water content on flowand heattransfer of water film on anti-icing surfaces［J］.Aeroengine，2016，42（1）：59-64.