結(jié)冰風(fēng)洞中液滴相變效應(yīng)數(shù)值模擬

郭向東，王梓旭，李明，劉蓓

1.中國(guó)空氣動(dòng)力研究與發(fā)展中心 空氣動(dòng)力學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，綿陽(yáng) 621000 2. 中國(guó)空氣動(dòng)力研究與發(fā)展中心 飛行器結(jié)冰與防除冰重點(diǎn)實(shí)驗(yàn)室，綿陽(yáng) 621000

飛機(jī)結(jié)冰廣泛存在于飛行實(shí)踐中，并嚴(yán)重威脅飛行安全[1-2]。研究飛機(jī)結(jié)冰的主要途徑包括數(shù)值模擬、結(jié)冰風(fēng)洞試驗(yàn)和飛行試驗(yàn)3種，其中結(jié)冰風(fēng)洞試驗(yàn)顧名思義是利用結(jié)冰風(fēng)洞在地面模擬飛機(jī)結(jié)冰過(guò)程，驗(yàn)證飛機(jī)防除冰系統(tǒng)，開(kāi)展飛機(jī)結(jié)冰問(wèn)題研究[3]。相對(duì)于數(shù)值模擬和飛行試驗(yàn)，結(jié)冰風(fēng)洞試驗(yàn)具有結(jié)冰條件易控、試驗(yàn)成本相對(duì)較低、結(jié)果可靠等優(yōu)點(diǎn)，是目前研究飛機(jī)結(jié)冰的主要手段。為模擬高空低溫低壓的云霧環(huán)境，結(jié)冰風(fēng)洞主要包括制冷系統(tǒng)、高度模擬系統(tǒng)和噴霧系統(tǒng)三大部分，其中噴霧系統(tǒng)利用噴嘴在風(fēng)洞穩(wěn)定段內(nèi)生成大量懸浮液滴，隨后液滴由低溫氣流攜帶進(jìn)入試驗(yàn)段，形成試驗(yàn)所需的過(guò)冷水滴云霧[4]。結(jié)冰試驗(yàn)中，氣流靜溫低于冰點(diǎn)，為防止噴霧系統(tǒng)中液態(tài)水結(jié)冰，須采用加熱方式提高液態(tài)水溫度，這樣便導(dǎo)致穩(wěn)定段內(nèi)生成的液滴溫度高于氣流靜溫，通常情況下，液滴會(huì)與氣流發(fā)生充分的熱交換，進(jìn)而在試驗(yàn)段內(nèi)達(dá)到過(guò)冷狀態(tài)(液滴溫度與試驗(yàn)段靜溫一致)，形成試驗(yàn)所需的過(guò)冷水滴[5]。但是，由于水滴與氣流間的熱交換過(guò)程十分復(fù)雜，涉及對(duì)流傳熱、相變傳質(zhì)和三維收縮效應(yīng)等多物理過(guò)程的耦合作用，試驗(yàn)參數(shù)的變化可能導(dǎo)致試驗(yàn)段內(nèi)液滴偏離過(guò)冷狀態(tài)，進(jìn)而改變積冰生長(zhǎng)過(guò)程和冰形特征，影響結(jié)冰風(fēng)洞試驗(yàn)結(jié)果，因此探索結(jié)冰風(fēng)洞中液滴傳熱過(guò)程，揭示試驗(yàn)參數(shù)的影響規(guī)律，評(píng)估液滴過(guò)冷狀態(tài)，對(duì)結(jié)冰風(fēng)洞試驗(yàn)準(zhǔn)確模擬結(jié)冰過(guò)程具有重要意義。

針對(duì)結(jié)冰風(fēng)洞液滴過(guò)冷問(wèn)題，Willbanks和Schulzt[6]率先發(fā)展了基于拉格朗日法的運(yùn)動(dòng)液滴傳質(zhì)傳熱計(jì)算方法，并對(duì)發(fā)動(dòng)機(jī)高空試驗(yàn)臺(tái)內(nèi)的液滴傳熱過(guò)程開(kāi)展了研究，隨后Miller等[7]采用這一方法，研究了NASA Glenn IRT結(jié)冰風(fēng)洞內(nèi)的液滴過(guò)冷問(wèn)題，發(fā)現(xiàn)初始液滴溫度對(duì)液滴過(guò)冷狀態(tài)無(wú)顯著影響，此外Bellucci等[8]采用試驗(yàn)方法對(duì)CIRA結(jié)冰風(fēng)洞高速試驗(yàn)段和低速試驗(yàn)段構(gòu)型的液滴過(guò)冷狀態(tài)進(jìn)行了評(píng)估。但是這些研究主要關(guān)注工程層面上液滴過(guò)冷狀態(tài)的評(píng)估，缺乏對(duì)結(jié)冰風(fēng)洞中液滴傳熱傳質(zhì)過(guò)程的規(guī)律性認(rèn)識(shí)。

近年來(lái)，隨著國(guó)內(nèi)大型結(jié)冰風(fēng)洞——3 m×2 m結(jié)冰風(fēng)洞的建成，液滴過(guò)冷問(wèn)題逐漸引起了一些國(guó)內(nèi)學(xué)者的關(guān)注[9-10]。其中郭向東等[10]發(fā)展了基于Euler法的氣液兩相耦合流動(dòng)計(jì)算方法，研究了3 m×2 m結(jié)冰風(fēng)洞主試驗(yàn)段構(gòu)型液滴傳熱過(guò)程，并評(píng)估了液滴過(guò)冷狀態(tài)，研究結(jié)果表明結(jié)冰風(fēng)洞中液滴傳熱過(guò)程可以分為準(zhǔn)一維傳熱和三維收縮傳熱兩個(gè)階段，其中三維收縮效應(yīng)對(duì)液滴過(guò)冷狀態(tài)起決定性作用，但是他們的研究中未考慮液滴相變傳質(zhì)對(duì)液滴傳熱過(guò)程的影響，導(dǎo)致對(duì)結(jié)冰風(fēng)洞液滴過(guò)冷狀態(tài)的評(píng)估結(jié)果較為保守。

因此，本文在郭向東等的計(jì)算方法基礎(chǔ)上，進(jìn)一步引入Hill矩方法[11]，發(fā)展了基于Euler法的氣液兩相傳質(zhì)傳熱耦合流動(dòng)計(jì)算方法，實(shí)現(xiàn)了液滴相變傳質(zhì)過(guò)程的模擬。進(jìn)而利用該方法，針對(duì)3 m×2 m結(jié)冰風(fēng)洞主試驗(yàn)段構(gòu)型，考察了液滴相變效應(yīng)對(duì)液滴傳熱過(guò)程的影響，揭示了相變效應(yīng)的影響規(guī)律，開(kāi)展了參數(shù)影響研究，分析了相對(duì)濕度、試驗(yàn)段氣流速度和液滴尺寸對(duì)液滴傳質(zhì)傳熱過(guò)程的影響，評(píng)估了液滴過(guò)冷狀態(tài)，為國(guó)內(nèi)飛機(jī)結(jié)冰地面模擬試驗(yàn)數(shù)據(jù)應(yīng)用評(píng)估提供支撐。

1 計(jì)算方法

1.1 控制方程和物理模型

本文借鑒了文獻(xiàn)[10]中基于Euler法的氣液兩相耦合流動(dòng)計(jì)算方法，進(jìn)一步耦合Hill矩函數(shù)方程組[12-14]，實(shí)現(xiàn)了相變傳質(zhì)過(guò)程的模擬。為獲得簡(jiǎn)化物理模型，對(duì)氣液兩相進(jìn)行以下假設(shè)[15]：①氣相為理想氣體，遵循理想氣體法則；②液滴為球形，不發(fā)生變形和破碎，液滴間不發(fā)生碰撞；③液滴內(nèi)溫度均勻分布；④忽略重力效應(yīng)；⑤ 黏性和熱傳導(dǎo)效應(yīng)僅發(fā)生在氣液兩相間；⑥不考慮液滴凍結(jié)過(guò)程以及水蒸氣凝結(jié)成核過(guò)程。

基于以上假設(shè)，簡(jiǎn)化后的控制方程組為

1) 氣相方程組

(1)

(2)

(3)

(4)

2) 液相方程組

(5)

(6)

(7)

3) 矩函數(shù)方程組

(8)

(9)

(10)

式中：下標(biāo)g和d分別代表氣相和液相;ρ、p、V和E分別為密度、壓力、速度和單位質(zhì)量總能量；氣相能量Eg和液相能量Ed的表達(dá)式分別為

(11)

式中：cv為氣相定容比熱；cl為液滴比熱；T為溫度。

σ和α分別為液滴有效密度和水蒸氣質(zhì)量分?jǐn)?shù)，其中σ對(duì)應(yīng)液態(tài)水含量(LWC)；n階矩函數(shù)定義為

(12)

其中：r為液滴半徑；f為粒子尺寸分布函數(shù)；Q0、Q1和Q2為零階、一階和二階矩函數(shù)，根據(jù)定義可見(jiàn)Q0表示單位體積內(nèi)液滴總個(gè)數(shù)。

S0為質(zhì)量源項(xiàng)，表示為

(13)

式中：ρl為液滴物理密度；dr/dt為液滴半徑變化率。

SD和SM分別為兩相動(dòng)量交換源項(xiàng)和兩相間傳質(zhì)引起的動(dòng)量傳遞源項(xiàng)，表示為

SD=Q0FD,SM=S0Vd

(14)

式中：FD為單個(gè)液滴受到的氣相作用力。

SEg和SEd分別為氣相和液相能量源項(xiàng)，表示為

(15)

式中：qheat與qrad分別為單個(gè)液滴與氣相間的對(duì)流傳熱率和輻射傳熱率；qgas,mass和qliquid,mass分別為單位質(zhì)量傳質(zhì)引起的氣相和液相能量變化率，表示為[16-17]

(16)

其中：cp為水蒸氣定壓比熱；L為汽化潛熱。

S1和S2分別為一階和二階矩函數(shù)方程源項(xiàng)，表示為

(17)

單個(gè)液滴受到的氣相作用力FD和單個(gè)液滴與氣相間的對(duì)流傳熱率qheat表示為[18]

(18)

qheat=Nukgπd(Tg-Td)

(19)

式中：d為液滴直徑；μg為氣相黏性系數(shù)；kg為氣相熱傳導(dǎo)系數(shù)；Red為液滴相對(duì)雷諾數(shù)，表示為

(20)

CD為液滴阻力系數(shù)，針對(duì)球型液滴，采用Schiller-Naumann圓球阻力模型，表示為

(21)

針對(duì)球型液滴，Nu采用Ranz-Marshall半經(jīng)驗(yàn)關(guān)系式，表示為

(22)

式中：Pr為Prandtl數(shù)。

單個(gè)液滴與氣流間的輻射傳熱率qrad表示為[19]

(23)

式中：ε為熱輻射發(fā)射率，根據(jù)文獻(xiàn)[19]，ε取0.95;σSB為Stefan-Boltzmann常數(shù)。

(24)

式中：DAB為二元質(zhì)量擴(kuò)散系數(shù)；Rv為水蒸氣氣體常數(shù)；pv,s為液滴表面飽和水氣壓；Hu為相對(duì)濕度，表示為

(25)

式中：pv,∞為氣相環(huán)境中水蒸氣分壓。

Sh為Sherwood數(shù)，針對(duì)球型液滴，Sh采用Ranz-Marshall半經(jīng)驗(yàn)關(guān)系式，表示為

(26)

式中：Sc為Schmidt數(shù)。根據(jù)文獻(xiàn)[16，21]，Ranz-Marshall半經(jīng)驗(yàn)關(guān)系式適用于Red<137、溫度低于0 ℃時(shí)球型液滴的傳質(zhì)傳熱過(guò)程，匹配本文的研究范圍(尤其對(duì)于高風(fēng)速工況下依然成立)，因此本文計(jì)算模型中Nu和Sh采用該關(guān)系式。

1.2 數(shù)值方法

本文采用與文獻(xiàn)[10]相同的數(shù)值方法，控制方程組采用有限體積方法離散，非定常項(xiàng)采用二階Euler格式離散，空間輸運(yùn)項(xiàng)采用二階迎風(fēng)型的Roe格式離散，源項(xiàng)采用隱式格式[10]。

氣液兩相邊界條件采用不同形式處理：氣相入口條件采用壓力入口，出口條件采用壓力出口，壁面采用無(wú)黏滑移壁面條件(1.1節(jié)假設(shè)⑤)；液相入口條件采用速度入口，出口采用出流條件，壁面處認(rèn)為液滴可以直接穿過(guò)，不考慮壁面對(duì)液相的作用。

1.3 驗(yàn)證算例

選取文獻(xiàn)[21]中的試驗(yàn)結(jié)果驗(yàn)證本文計(jì)算方法。文獻(xiàn)[21]采用懸掛液滴試驗(yàn)方法研究了液滴凍結(jié)過(guò)程中的液滴溫度變化特性，文中將液滴凍結(jié)過(guò)程分為過(guò)冷(Supercooling Stage)、復(fù)輝(Recalescence Stage)、凍結(jié)(Freezing Stage)和冷卻(Cooling Stage)4個(gè)階段，其中過(guò)冷階段對(duì)應(yīng)本文研究的內(nèi)容，因此選取文獻(xiàn)中過(guò)冷階段的試驗(yàn)結(jié)果驗(yàn)證本文計(jì)算方法。驗(yàn)證結(jié)果如圖1所示，圖中：Va和Ta分別表示氣流速率和氣流靜溫。從圖中可以看出：本文計(jì)算結(jié)果與試驗(yàn)結(jié)果吻合得較好，不同工況下液滴溫度誤差均小于1%，因此本文計(jì)算方法可以合理模擬液滴過(guò)冷階段的傳熱過(guò)程。針對(duì)以上的驗(yàn)證過(guò)程，應(yīng)該指出的是考慮到懸掛液滴法對(duì)液滴尺寸和試驗(yàn)風(fēng)速的限制，該方法僅獲得了大液滴(毫米量級(jí))低風(fēng)速(小于1 m/s)下的試驗(yàn)結(jié)果，而缺乏小液滴(微米量級(jí))高風(fēng)速的試驗(yàn)結(jié)果。進(jìn)一步，文獻(xiàn)[21]指出，盡管缺乏真實(shí)噴霧的試驗(yàn)數(shù)據(jù)，但是試驗(yàn)數(shù)據(jù)中液滴相對(duì)雷諾數(shù)范圍(36

2 結(jié)果與分析

本節(jié)針對(duì)3 m×2 m結(jié)冰風(fēng)洞主試驗(yàn)段收縮構(gòu)型，采用本文發(fā)展的計(jì)算方法，考察了基準(zhǔn)工況下液滴相變效應(yīng)對(duì)液滴傳熱過(guò)程的影響，揭示了相變效應(yīng)的影響規(guī)律，開(kāi)展了參數(shù)影響研究，分析了相對(duì)濕度、試驗(yàn)段氣流速度和液滴尺寸對(duì)液滴傳質(zhì)傳熱過(guò)程的影響，評(píng)估了試驗(yàn)段中心處液滴過(guò)冷狀態(tài)。

3 m×2 m結(jié)冰風(fēng)洞主試驗(yàn)段構(gòu)型如圖2所示，該構(gòu)型包括穩(wěn)定段、收縮段和試驗(yàn)段3部分，其中構(gòu)型入口截面為高8 m×寬11 m，位于穩(wěn)定段噴霧耙處，出口截面為高2 m×寬3 m，位于試驗(yàn)段中心處，構(gòu)型全長(zhǎng)16 m，收縮比為14.67。

根據(jù)文獻(xiàn)[10]的研究結(jié)果可知，試驗(yàn)段氣流速度和液滴尺寸是影響液滴過(guò)冷狀態(tài)的主要參數(shù)，因此本文選取相對(duì)濕度、試驗(yàn)段氣流速度和液滴尺寸3個(gè)參數(shù)開(kāi)展參數(shù)研究，計(jì)算工況矩陣如表1所示，表中LWCi、di、Tdi、T0i、Udi、UTS和Hui分別代表初始液態(tài)水含量、初始液滴直徑、初始液滴溫度、初始?xì)饬骺倻亍⒊跏家旱嗡俾?、試?yàn)段氣流速度、初始相對(duì)濕度。Case 1為基準(zhǔn)工況，考察典型工況下結(jié)冰風(fēng)洞中相變對(duì)液滴傳熱過(guò)程的影響；Case 2選取了70%、80%、90%和100%這4個(gè)典型的相對(duì)濕度狀態(tài)，考察相對(duì)濕度的影響；Case 3選取了55 、77、130、164 m/s這4個(gè)典型試驗(yàn)段風(fēng)速，考察氣流速度的影響，其中各風(fēng)速對(duì)應(yīng)的入口總壓分別為2 000、4 000、12 000、20 000 Pa，參考?jí)毫?05Pa；Case 4選取了40、60、80、100 μm這4個(gè)典型液滴直徑，考察液滴尺寸的影響，值得指出的是由于假設(shè)液滴為球型，不考慮液滴變形和破碎等動(dòng)力學(xué)特性，因此僅選取了小尺寸液滴開(kāi)展研究。

表1 計(jì)算工況矩陣Table 1 Calculation condition matrix

2.1 基準(zhǔn)工況

基準(zhǔn)工況選取FAR Part25 附錄C中典型結(jié)冰條件[22]，為考察相變對(duì)液滴傳熱過(guò)程的影響，下文通過(guò)對(duì)比構(gòu)型內(nèi)液滴和顆粒的傳熱過(guò)程，進(jìn)而揭示液滴相變效應(yīng)的影響規(guī)律，其中針對(duì)顆粒傳熱過(guò)程的模擬，本文采用文獻(xiàn)[10]中的計(jì)算方法來(lái)實(shí)現(xiàn)。

圖3給出了基準(zhǔn)工況下構(gòu)型中心線處氣流靜溫和液滴溫度變化曲線，從圖中可以看出，類(lèi)似于文獻(xiàn)[10]中顆粒傳熱過(guò)程，整個(gè)液滴傳熱過(guò)程近似以構(gòu)型6 m處截面為界分成準(zhǔn)一維傳熱和三維收縮傳熱兩個(gè)階段。在準(zhǔn)一維傳熱階段內(nèi)(入口至6 m處截面，包括穩(wěn)定段和收縮段入口區(qū)域)，風(fēng)洞型面變化較緩，氣流以準(zhǔn)一維方式低速流動(dòng)，氣流速度和靜溫變化較小，接近入口流速和總溫，液滴以準(zhǔn)一維穩(wěn)態(tài)方式傳熱；在三維收縮傳熱階段內(nèi)(6 m處截面至出口，包括收縮段后部區(qū)域和試驗(yàn)段)，風(fēng)洞型面變化較大，三維收縮效應(yīng)較強(qiáng)，氣流場(chǎng)在此階段內(nèi)發(fā)生了顯著變化，氣流速度加速增大，氣流靜溫則等熵下降，出口氣流速度越大，出口氣流靜溫越低，三維收縮效應(yīng)(氣流可壓縮性)顯著影響液滴傳熱過(guò)程，液滴以三維動(dòng)態(tài)方式傳熱。對(duì)比液滴與顆粒的溫度變化曲線可見(jiàn)：在準(zhǔn)一維傳熱階段，液滴溫度下降趨勢(shì)快于顆粒溫度，當(dāng)液滴溫度和顆粒溫度趨于穩(wěn)定時(shí)(約在0.6 m處)，液滴溫度(約-5.4 ℃)低于顆粒溫度(-5 ℃)，此時(shí)液滴溫度稱(chēng)為液滴濕球溫度[16]；在三維收縮傳熱階段，液滴溫度下降趨勢(shì)慢于顆粒，最終在構(gòu)型出口處(試驗(yàn)段中心)，液滴溫度(-7.41 ℃)高于顆粒溫度(-7.77 ℃)。

圖4給出了構(gòu)型中心線處液滴直徑和相對(duì)濕度變化曲線，從圖中可以看出：結(jié)冰風(fēng)洞中液滴經(jīng)歷了先蒸發(fā)后凝結(jié)兩個(gè)階段，在準(zhǔn)一維傳熱階段，氣流中水蒸氣未達(dá)到飽和狀態(tài)，液滴處于蒸發(fā)狀態(tài)，液滴直徑不斷減小，同時(shí)蒸發(fā)產(chǎn)生的水蒸氣則使氣流的相對(duì)濕度不斷增大；在三維收縮階段，氣流靜溫迅速下降，氣流飽和水汽壓降低，而水蒸氣分壓變化較小(適用于低速流動(dòng))，導(dǎo)致氣流的相對(duì)濕度迅速增大(對(duì)應(yīng)式(25))，氣流中水蒸氣從未飽和狀態(tài)轉(zhuǎn)換為過(guò)飽和狀態(tài)，液滴則從蒸發(fā)狀態(tài)轉(zhuǎn)換為凝結(jié)狀態(tài)。

進(jìn)一步理論分析相變效應(yīng)對(duì)液滴傳熱的影響。首先液相能量方程式(7)可以化簡(jiǎn)為

(27)

式中:由于液滴過(guò)冷階段輻射傳熱遠(yuǎn)小于對(duì)流傳熱和相變傳熱，因此忽略輻射傳熱項(xiàng)，同時(shí)聯(lián)立式(19)和式(24)，則式(27)可化簡(jiǎn)為

(28)

ξheat+ξmass

(29)

式中：ξheat為對(duì)流傳熱引起的液滴溫度空間變化率；ξmass為液滴相變引起的液滴溫度空間變化率。需要指出的是，式(29)是在一維定常流動(dòng)條件下得到的，但是對(duì)于三維收縮階段，構(gòu)型中心線處的氣流沿流向加速流動(dòng)，仍可以近似認(rèn)為氣流運(yùn)動(dòng)為準(zhǔn)一維，因此式(29)適用于構(gòu)型中心線處整個(gè)液滴傳熱過(guò)程的分析。

圖5給出了構(gòu)型中心線處ξmass和ξheat的變化曲線，從圖中可以看出，在準(zhǔn)一維傳熱階段，當(dāng)溫度較高的液滴進(jìn)入未飽和的低溫氣流(Hui<100%)中時(shí)，液滴與氣流間發(fā)生對(duì)流傳熱，液滴向氣流傳遞熱能(ξheat<0)，同時(shí)蒸發(fā)效應(yīng)從液滴吸收所需的汽化潛熱(ξmass<0)，而顆粒與氣流間僅發(fā)生對(duì)流傳熱過(guò)程，因此液滴溫度下降速度快于顆粒的速度(0～0.6 m)；當(dāng)液滴溫度下降至氣流靜溫時(shí)[23]，由于蒸發(fā)效應(yīng)仍持續(xù)從液滴吸收熱量，因此液滴溫度會(huì)持續(xù)下降，但對(duì)于顆粒而言，其溫度此時(shí)將會(huì)與氣流靜溫保持一致，顆粒與氣流間的對(duì)流傳熱過(guò)程停止，顆粒達(dá)到穩(wěn)定狀態(tài)；進(jìn)一步，當(dāng)液滴溫度下降至液滴濕球溫度時(shí)，蒸發(fā)效應(yīng)所需的汽化潛熱(ξmass<0)與氣流向液滴傳遞的熱能(ξheat>0)達(dá)到平衡(|ξmass|=|ξheat|)，液滴溫度趨于穩(wěn)定，此時(shí)ξmass與ξheat沿流向(x軸)對(duì)稱(chēng)分布(0.6～6 m)，液滴進(jìn)入穩(wěn)定蒸發(fā)階段。在三維收縮階段，三維收縮效應(yīng)導(dǎo)致氣流靜溫迅速下降，液滴穩(wěn)定蒸發(fā)過(guò)程被破壞，ξmass與ξheat不再對(duì)稱(chēng)分布，液滴溫度下降(ξmass+ξheat<0)；此階段內(nèi)對(duì)流傳熱的能量傳遞方向發(fā)生了轉(zhuǎn)換，即ξheat>0轉(zhuǎn)換為ξheat<0，促進(jìn)了液滴溫度下降，同時(shí)氣流中水蒸氣從未飽和狀態(tài)轉(zhuǎn)換為過(guò)飽和狀態(tài)，液滴則從蒸發(fā)狀態(tài)轉(zhuǎn)換為凝結(jié)狀態(tài)，即ξmass<0轉(zhuǎn)換為ξmass>0，則抑制了液滴溫度下降，而顆粒與氣流間僅存在對(duì)流換熱過(guò)程，缺少了凝結(jié)過(guò)程的抑制效應(yīng)，因此可見(jiàn)凝結(jié)效應(yīng)是三維收縮階段內(nèi)液滴溫度下降趨勢(shì)慢于顆粒的主要原因，最終導(dǎo)致構(gòu)型出口處(試驗(yàn)段中心)，液滴溫度高于顆粒溫度。

2.2 相對(duì)濕度影響

首先考察相對(duì)濕度的影響，在結(jié)冰風(fēng)洞試驗(yàn)中，為避免氣流濕度過(guò)低導(dǎo)致液滴大量蒸發(fā)，會(huì)采用加濕器提高風(fēng)洞回路內(nèi)氣流濕度，因此此處選取較高的相對(duì)濕度參數(shù)開(kāi)展研究。

圖6(a)和圖6(b)給出了不同初始相對(duì)濕度下構(gòu)型中心線處液滴直徑和相對(duì)濕度變化曲線，從圖中可以看出：在準(zhǔn)一維傳熱階段，隨著初始相對(duì)濕度的增加，液滴尺寸的減小趨勢(shì)減慢，表明液滴蒸發(fā)效應(yīng)減弱，尤其當(dāng)濕度100%時(shí)，液滴尺寸在溫度穩(wěn)定階段不再發(fā)生變化，蒸發(fā)效應(yīng)消失；與此相反，在三維收縮階段，隨著初始相對(duì)濕度的增加，水蒸氣過(guò)飽和程度增強(qiáng)，液滴尺寸的增大趨勢(shì)加快，表明凝結(jié)效應(yīng)增強(qiáng)。

圖6(c)給出了不同初始相對(duì)濕度下構(gòu)型中心線處液滴溫度變化曲線，從圖中可以看出：在準(zhǔn)一維傳熱階段，隨著初始相對(duì)濕度增加，液滴蒸發(fā)效應(yīng)減弱，液滴溫度下降趨勢(shì)減慢(如放大圖所示)，當(dāng)液滴趨于穩(wěn)定蒸發(fā)狀態(tài)時(shí)(約在0.6 m處)，初始相對(duì)濕度越大則液滴濕球溫度越高，尤其當(dāng)相對(duì)濕度為100%時(shí)，蒸發(fā)效應(yīng)消失，液滴溫度與顆粒溫度趨于一致；在三維收縮傳熱階段，隨著初始相對(duì)濕度增加，液滴凝結(jié)效應(yīng)增強(qiáng)，導(dǎo)致構(gòu)型出口處(試驗(yàn)段中心)液滴溫度不斷增大，尤其當(dāng)相對(duì)濕度100%時(shí)，液滴凝結(jié)效應(yīng)最強(qiáng)，出口液滴溫度最高，比顆粒溫度高近1 ℃。

圖7給出了不同初始相對(duì)濕度下構(gòu)型出口中心處液氣溫度差，從圖中可以看出：隨著初始相對(duì)濕度增加，液氣溫差不斷增大，其中初始相對(duì)濕度為70%時(shí)，液滴溫差與顆粒的一致，約為0.1 ℃，而當(dāng)初始相對(duì)濕度增加至100%時(shí)，液滴溫差增至約0.87 ℃，比顆粒溫度差高近0.8 ℃。由此可見(jiàn)，相對(duì)濕度會(huì)影響試驗(yàn)段內(nèi)液滴過(guò)冷狀態(tài)，但影響程度較弱，本文計(jì)算工況下最大液氣溫差(相對(duì)濕度為100%)小于1 ℃，液滴仍認(rèn)為處于過(guò)冷狀態(tài)。

2.3 試驗(yàn)段氣流速度影響

考察試驗(yàn)段氣流速度的影響，根據(jù)3 m×2 m結(jié)冰風(fēng)洞典型試驗(yàn)工況，選取55、77、130、164 m/s這4個(gè)試驗(yàn)段氣流速度開(kāi)展研究，其中各試驗(yàn)段氣流速度(構(gòu)型出口氣流速度)對(duì)應(yīng)的穩(wěn)定段氣流速度(構(gòu)型入口氣流速度)分別為4.11、5.70、9.15及11.00 m/s。

圖8給出了不同試驗(yàn)段氣流速度下構(gòu)型中心線處液滴直徑和相對(duì)濕度變化曲線，從圖中可以看出：在準(zhǔn)一維傳熱階段，隨著試驗(yàn)段氣流速度的增加，穩(wěn)定段氣流速度增大，液滴相變時(shí)間減少，液滴尺寸的減小趨勢(shì)減慢，表明蒸發(fā)效應(yīng)減弱；在三維收縮階段，隨著試驗(yàn)段氣流速度的增加，盡管液滴相變時(shí)間不斷減少，但是在三維收縮效應(yīng)的影響下，氣流靜溫顯著降低，水蒸氣過(guò)飽和度顯著增大，導(dǎo)致液滴尺寸的增大趨勢(shì)加快，進(jìn)而表明凝結(jié)效應(yīng)增強(qiáng)。

圖9給出了不同試驗(yàn)段氣流速度下構(gòu)型中心線處氣流靜溫和液滴溫度變化曲線，從圖中可以看出，在準(zhǔn)一維傳熱階段，隨著試驗(yàn)段氣流速度的增加，穩(wěn)定段氣流速度增大，但是由于穩(wěn)定段速度增大幅度較小(從約4 m/s增至11 m/s)，氣流可壓縮性較弱，因此各工況下的氣流靜溫和液滴溫度變化趨勢(shì)無(wú)顯著差異，氣流靜溫接近氣流總溫，并且液滴溫度穩(wěn)定后，液滴濕球溫度也基本一致。在三維收縮階段，隨著試驗(yàn)段氣流速度的增加，氣流可壓縮性顯著增強(qiáng)，試驗(yàn)段氣流靜溫顯著下降，具體的當(dāng)試驗(yàn)段氣流速度從55 m/s增至164 m/s時(shí)，試驗(yàn)段氣流靜溫約從-7 ℃降至-18 ℃，因此在結(jié)冰風(fēng)洞試驗(yàn)中，氣流可壓縮性導(dǎo)致的氣流靜溫變化是不可忽略的；同時(shí)，在這一階段隨著試驗(yàn)段氣流速度的增加，液滴溫度的下降趨勢(shì)加快，但是相對(duì)于顆粒，由于凝結(jié)效應(yīng)的不斷增強(qiáng)，抑制了液滴溫度的下降趨勢(shì)，導(dǎo)致出口處(試驗(yàn)段中心)液滴與顆粒溫差不斷增大。

圖10給出了不同試驗(yàn)段氣流速度下構(gòu)型出口中心處液氣溫度差，從圖中可以看出：隨著試驗(yàn)段氣流速度的增加，液氣溫差不斷增大，并且在相變效應(yīng)的影響下，其增大趨勢(shì)強(qiáng)于顆粒，其中試驗(yàn)段氣流速度為55 m/s時(shí)，液滴溫差約為0.3 ℃，而顆粒的溫差則趨近于0 ℃，當(dāng)試驗(yàn)段氣流速度增加至164 m/s時(shí)，液滴溫差增至約3.3 ℃，比顆粒的溫差高近2 ℃。由此可見(jiàn)，氣流速度會(huì)顯著影響試驗(yàn)段內(nèi)液滴過(guò)冷狀態(tài)，同時(shí)相變效應(yīng)則會(huì)顯著增強(qiáng)影響程度，其中小尺寸液滴(40 μm≤di≤100 μm)在高風(fēng)速時(shí)(UTS≥164 m/s)將偏離過(guò)冷狀態(tài)(液氣溫差超過(guò)3 ℃)。

2.4 液滴尺寸影響

最后，考察液滴尺寸的影響，由于采用球型液滴阻力模型，不考慮液滴變形和破碎效應(yīng)，因此在We數(shù)的約束下[24]，僅選擇尺寸較小的液滴，以較好地滿(mǎn)足球型液滴假設(shè)。

圖12給出了不同初始液滴尺寸下構(gòu)型中心線處液滴溫度變化曲線，從圖中可以看出：在準(zhǔn)一維傳熱階段，隨著初始液滴尺寸增加，液滴熱容量增大，液滴溫度下降趨勢(shì)減慢，液滴溫度穩(wěn)定后，液滴濕球溫度無(wú)顯著變化，液滴與顆粒溫差基本一致；在三維收縮階段，隨著初始液滴尺寸增加，液滴溫度下降趨勢(shì)減慢(與文獻(xiàn)[10]的結(jié)論一致)，但是相對(duì)于顆粒，凝結(jié)效應(yīng)不斷減弱，促進(jìn)了液滴溫度的下降趨勢(shì)，導(dǎo)致出口處(試驗(yàn)段中心)液滴與顆粒溫差不斷減小。

圖13給出了不同初始液滴尺寸下構(gòu)型出口中心處液氣溫度差，從圖中可以看出：隨著初始液滴尺寸增加，液氣溫差不斷增大，但是在相變效應(yīng)的影響下，其增大趨勢(shì)弱于顆粒的溫度差，其中液滴直徑為40 μm時(shí)，液滴溫度差約為0.6 ℃，而顆粒溫度差約為0.1 ℃，當(dāng)液滴直徑增加至100 μm時(shí)，液滴溫度差增至約0.85 ℃，比顆粒的溫差低近0.05 ℃。由此可見(jiàn)，液滴尺寸會(huì)影響試驗(yàn)段液滴過(guò)冷狀態(tài)，但相變效應(yīng)卻會(huì)減弱影響程差，其中在低風(fēng)速(UTS≤77 m/s)試驗(yàn)工況下，小尺寸液滴(di≤100 μm)仍可以認(rèn)為處于過(guò)冷狀態(tài)(液氣溫差小于1 ℃)，這與文獻(xiàn)[10]的評(píng)估結(jié)果一致。

3 結(jié) 論

本文發(fā)展了基于Euler法的氣液兩相傳質(zhì)傳熱耦合流動(dòng)計(jì)算方法，模擬了結(jié)冰風(fēng)洞中氣液兩相傳質(zhì)傳熱耦合流動(dòng)過(guò)程，研究了結(jié)冰風(fēng)洞中液滴相變效應(yīng)。

1) 結(jié)冰風(fēng)洞中液滴經(jīng)歷了先蒸發(fā)后凝結(jié)兩個(gè)階段。在準(zhǔn)一維傳熱階段，液滴為蒸發(fā)狀態(tài)，蒸發(fā)效應(yīng)加快了液滴溫度下降趨勢(shì)，并且促使液滴溫度趨于濕球溫度；在三維收縮階段，液滴則從蒸發(fā)狀態(tài)轉(zhuǎn)換為凝結(jié)狀態(tài)，凝結(jié)效應(yīng)減慢了液滴溫度下降趨勢(shì)，導(dǎo)致構(gòu)型出口處(試驗(yàn)段中心)液氣溫差增大，影響液滴過(guò)冷狀態(tài)。

2) 初始相對(duì)濕度增加，會(huì)減弱蒸發(fā)效應(yīng)而增強(qiáng)凝結(jié)效應(yīng)，進(jìn)而影響試驗(yàn)段內(nèi)液滴過(guò)冷狀態(tài)。但其影響程度較弱，本文計(jì)算工況下最大液氣溫差(相對(duì)濕度為100%)小于1 ℃，液滴仍認(rèn)為處于過(guò)冷狀態(tài)。

3) 氣流速度會(huì)顯著影響試驗(yàn)段內(nèi)液滴過(guò)冷狀態(tài)，隨著氣流速度的增加，蒸發(fā)效應(yīng)減弱而凝結(jié)效應(yīng)增強(qiáng)，進(jìn)而顯著增強(qiáng)了影響程度，其中小尺寸液滴(40 μm≤di≤100 μm)在高風(fēng)速時(shí)(UTS≥164 m/s)將偏離過(guò)冷狀態(tài)(液氣溫差超過(guò)3 ℃)。

4) 初始液滴尺寸增加，會(huì)導(dǎo)致蒸發(fā)效應(yīng)和凝結(jié)效應(yīng)均減弱，進(jìn)而減弱了液滴尺寸對(duì)試驗(yàn)段內(nèi)液滴過(guò)冷狀態(tài)的影響程度，對(duì)于小尺寸液滴(di≤100 μm)在低風(fēng)速(UTS≤77 m/s)試驗(yàn)工況下，仍可以認(rèn)為處于過(guò)冷狀態(tài)(液氣溫差小于1 ℃)。

需要指出的是，實(shí)際結(jié)冰風(fēng)洞中，氣液兩相流動(dòng)過(guò)程十分復(fù)雜，涉及動(dòng)量傳遞、質(zhì)量傳遞和能量傳遞等多物理過(guò)程的強(qiáng)耦合作用，而本文發(fā)展的基于Euler法的氣液兩相傳質(zhì)傳熱耦合流動(dòng)計(jì)算方法是基于一定假設(shè)下的簡(jiǎn)化分析方法，因此開(kāi)展計(jì)算方法驗(yàn)證試驗(yàn)則是下一步急需開(kāi)展的工作。

此外，針對(duì)液滴直徑超過(guò)100 μm的大尺寸液滴，由于其存在變形和破碎等動(dòng)力學(xué)效應(yīng)，超出了本文采用的球形液滴數(shù)學(xué)模型適用范圍，因此建立大尺寸液滴傳質(zhì)傳熱模型，則是今后需要開(kāi)展的研究。

參 考 文 獻(xiàn)

[1] 杜雁霞， 李明， 桂業(yè)偉， 等. 飛機(jī)結(jié)冰熱力學(xué)行為研究綜述[J]. 航空學(xué)報(bào)， 2017， 38(2)： 520706.

DU Y X，LI M, GUI Y W, et al. Review of thermodynamic behaviors in aircraft icing process[J]. Acta Aeronautica et Astronautica Sinica, 2017, 38(2): 520706 (in Chinese).

[2] 桂業(yè)偉， 周志宏， 李穎暉， 等. 關(guān)于飛機(jī)結(jié)冰的多重安全邊界問(wèn)題[J]. 航空學(xué)報(bào)， 2017， 38(2): 520723.

GUI Y W, ZHOU Z H, LI Y H, et al. Multiple safety boundaries protection on aircraft icing[J]. Acta Aeronautica et Astronautica Sinica, 2017, 38(2): 520723 (in Chinese).

[3] 易賢， 王斌， 李偉斌， 等. 飛機(jī)結(jié)冰冰形測(cè)量方法研究進(jìn)展[J]. 航空學(xué)報(bào)， 2017, 38(2): 520700.

YI X, WANG B, LI W B, et al. Research progress on ice shape measurement approaches for aircraft icing[J]. Acta Aeronautica et Astronautica Sinica, 2017, 38(2): 520700 (in Chinese).

[4] RAGNI A, ESPOSITO B, MARRAZZO M, et al. Calibration of the CIRA IWT in the high speed configuration[C]//43th AIAA Aerospace Science Meeting and Exhibit. Reston, VA: AIAA, 2005.

[5] STEEN L E, IDE R F, van ZANTE J F, et al. NASA Glenn icing research tunnel: 2014 and 2015 cloud calibration procedures and results: NASA/TM-2015-218758 [R]. Washington, D. C. : NASA, 2015.

[6] WILLBANKS C E, SCHULZT R J. Analytical study of icing simulation for turbine engines in altitude test cells[J]. Journal of Aircraft, 1975, 12(12): 960-967.

[7] MILLER D R, ADDY H E, IDE R F. A study of large droplet ice accretions in the nasa glenn irt at near-freezing conditions: NASA TM-1996-107142-REV1 [R]. Washington, D. C. : NASA, 1996.

[8] BELLUCCI M, ESPOSITO B M, MARRAZZO M, et al. Calibration of the CIRA IWT in the low speed configuration[C]∥45th AIAA Aerospace Science Meeting and Exhibit. Reston, VA： AIAA，2007.

[9] 易賢, 馬洪林, 王開(kāi)春, 等. 結(jié)冰風(fēng)洞液滴運(yùn)動(dòng)及傳質(zhì)傳熱特性分析[J]. 四川大學(xué)學(xué)報(bào), 2012, 44(2): 132-135.

YI X, MA H L, WANG K C, et al. Analysis of water droplet movement and heat / mass transfer in an icing wind tunnel[J]. Journal of Sichuan University, 2012, 44(2): 132-135 (in Chinese).

[10] 郭向東，王梓旭，李明，等．結(jié)冰風(fēng)洞中液滴過(guò)冷特性數(shù)值研究[J]．航空學(xué)報(bào), 2018, 39(1): 121254.

GUO X D, WANG Z X, LI M, et al. Numerical study of droplet supercooling in an icing wind tunnel [J]. Acta Aeronautica et Astronautica Sinica, 2018, 39(1): 121254 (in Chinese).

[11] HILL P G. Condensation of water vapor during supersonic expansion in nozzles[J]. Journal of Fluid Mechanics, 1966, 25(3): 593-620.

[12] CHENG W, LUO X S, van DONGEN M E H. On condensation-induced waves[J]. Journal of Fluid Mechanics, 2010, 651(1): 145-164.

[13] 郭向東, 黃生洪, 吳穎川, 等. 氫燃料燃燒加熱風(fēng)洞中水蒸氣相變效應(yīng)的數(shù)值研究[J]. 推進(jìn)技術(shù), 2017, 38(4):932-941.

GUO X D, HUANG S H, WU Y C, et al. Numerical evaluation of water vapor phase transition effects in a hydrogen-fueled combustion-heated wind tunnel[J]. Journal of Propulsion Technology, 2017, 38(4): 932-941 (in Chinese).

[14] CAO Y, LUO X S. Numerical study on Homogeneous condensation in a vortex[J]. Procedia Engineering, 2015, 126: 607-611.

[15] LUO X S, WANG G, OLIVIER H. Parametric investigation of particle acceleration in high enthalpy conical nozzle flows for coating applications[J]. Shock Waves, 2008, 17(5): 351-362.

[16] KERSEY J, LOTH E, LANKFORD D. Effects of Evaporating Droplets on Shock Waves[J]. AIAA Journal, 2010, 48(9): 1975-1986.

[17] RUSSO E, KUERTEN J G M, van der Geld C W M, et al. Water droplet condensation and evaporation in turbulent channel flow[J]. Journal of Fluid Mechanics, 2014, 749: 666-700.

[18] SAITO T. Numerical analysis of dusty-gas flows[J]. Journal of Computational Physics, 2002, 176(1): 129-144.

[19] FUJITA A, KUROSE R, KOMORI S. Experimental study on effect of relative humidity on heat transfer of an evaporating water droplet in air flow[J]. International Journal of Heat and Mass Transfer, 2010, 36(3): 244-247.

[20] MILLER R S, HARSTAD K, BELLAN J. Evaluation of equilibrium and non-equilibrium evaporation models for many-droplet gas-liquid flow simulations[J]. International Journal of Multiphase Flow, 1998, 24(6): 1025-1055.

[21] HINDMARSH J P, RUSSEL A B, CHEN X D. Experimental and numerical analysis of the temperature transition of a suspended freezing water droplet[J]. International Journal of Heat and Mass Transfer, 2003, 46(7): 1199-1213.

[22] COBER S, BERNSTEIN B, JECK R, et al. Data and analysis for the development of an engineering standard for supercooled large drop conditions: DOT/FAA/AR-09/10[R]. Washington, D. C.: FAA, 2009.

[23] KUROSE R, FUJITA A, KOMORI S. Effect of relative humidity on heat transfer across the surface of an evaporating water droplet in air flow[J]. Journal of Fluid Mechanics, 2009, 624: 57-67.

[24] LUXFORD G. Experimental and modelling investigation of the deformation, drag and break-up of drizzle droplets subjected to strong aerodynamic forces in relation to SLD aircraft icing[D]. Bedfordshire: Cranfield University, 2005: 5-7 .