李曉萍, 劉向農, 孫東方, 劉如佳, 董 婷

(合肥工業(yè)大學 汽車與交通工程學院,安徽 合肥 230009)

0 引 言

空氣源熱泵已經廣泛應用在生活中,但是當冬季蒸發(fā)器表面溫度低于0 ℃且低于空氣的露點溫度時,換熱器表面就會結霜。結霜會造成換熱器熱阻增加,從而導致?lián)Q熱器的性能下降。因此研究結霜機理對于換熱器的結構設計以及除霜具有一定的價值。

結霜過程不僅僅是一個相變過程,更是一個復雜的傳熱傳質過程。文獻[1]根據霜層生長初期與樹枝生長的相似性,基于有限擴散聚集(diffusion-limited aggregation,DLA)模型建立結霜分形模型,并利用該模型模擬空氣流速、冷表面接觸角對結霜的影響,但并未引入表面粗糙度對結霜的影響;文獻[2]基于DLA模型,在MATLAB平臺上建立霜層生長分形模型,模擬深冷表面霜層初期生長過程,但是并未討論霜層導熱系數(shù)隨時間的變化。隨著計算流體動力學(computational fluid dynamics,CFD)的發(fā)展,一些學者利用CFD來模擬霜層生長。文獻[3]利用Fluent建立霜層生長模型,基于經典成核理論將霜層生長分為冰晶成核和冰晶生長2個部分,建立水蒸氣向冰相的質量傳遞速率模型,這為模擬霜層的生長過程提供了另外一種方法;文獻[4]也利用Fluent模擬霜層的生長,從飽和水蒸氣濃度與其同溫度下的水蒸氣濃度差是霜層生長驅動力的角度出發(fā),建立了霜層生長速率的關系式,但是霜面水蒸氣采用飽和水蒸氣;文獻[5]提出霜面上水蒸氣應處于過飽和狀態(tài),通過與實驗數(shù)據的對比發(fā)現(xiàn),霜面水蒸氣過飽和模型優(yōu)于霜面水蒸氣飽和模型;文獻[6]指出,若霜面上濕空氣是飽和狀態(tài),則傳質速率會被高估,并提出了霜面水蒸氣過飽和度方程;文獻[7]指出,水蒸氣濃度差在結霜過程中占主導地位,利用CFD模擬結霜可以消除關于初始結霜厚度和密度的2個假設,同時還可以獲得霜導熱系數(shù)和濕空氣絕對濕度的分布。還有學者對結霜進行了實驗研究,提出了霜層生長速率的經驗關聯(lián)式。文獻[8]進行了霜層的實驗研究,對影響結霜的因素進行分析,提出了霜層生長的關聯(lián)式;文獻[9]針對水平冷板上霜層的生長過程,提出霜層密度的半經驗關聯(lián)式。但這些關聯(lián)式有一定的局限性。

雖然目前已經有大量關于結霜的研究,但是基于霜-氣界面水蒸氣處于過飽和狀態(tài)使用CFD模型對結霜形態(tài)的研究很少。CFD既可以用來獲得霜的形態(tài)、生長速率和霜物性參數(shù),也可以用于結構復雜的換熱器表面結霜的研究,以補充不同工況下結霜的實驗數(shù)據。

本文以霜-氣界面上過飽和水蒸氣與霜域內水蒸氣濃度差為結霜驅動勢建立相變率表達式,將相變率表達式導入到Fluent中,并耦合歐拉兩相流模型建立質量傳遞模型以研究水平冷板霜層生長與密實化過程。

1 物理模型

霜層生長的物理模型如圖1所示。數(shù)值模擬區(qū)域長75 mm,寬10 mm,濕空氣進口參數(shù)和低溫冷板溫度見表1所列。

在計算區(qū)域內,霜層的形態(tài)是時刻變化的,為了便于觀察不同時刻下霜層的形態(tài),引入霜-氣界面將計算域劃分為濕空氣域和霜域,采用冰晶體積分數(shù)作為閾值監(jiān)測霜-氣界面的更替變化。濕空氣域內冰晶體積分數(shù)為0;霜區(qū)域內冰晶體積分數(shù)小于臨界冰晶體積分數(shù),取10-6[3];霜氣界面上為臨界冰晶體積分數(shù)。

結霜是一個復雜的傳熱傳質過程,因此在計算過程中需要一些假設來簡化計算:將濕空氣視為不可壓縮的理想流體;鋁板溫度低,濕空氣直接凝固為冰;在結霜過程中忽略融化和升華。

表1 初始參數(shù)

2 霜層生長與致密化模型

在計算區(qū)域內,霜-氣界面的水蒸氣處于過飽和狀態(tài),這與霜域內的水蒸氣形成濃度差,結霜驅動勢形成,冷板上產生冰晶。隨著時間的增長,一部分濕空氣會增加霜層的厚度,另一部分濕空氣會增加霜的密實度。本文將重點討論霜厚度生長和密度生長規(guī)律。

本文采用Fluent中的歐拉兩相流模型。在模型中設置兩相物質,第1相為連續(xù)濕空氣(干空氣和水蒸氣),第2相為冰。將濕空氣-霜的相變率表達式以udf形式編程到Fluent中來實現(xiàn)兩相之間的質量傳遞。

2.1 質量傳遞模型

霜的形成需要經過冰晶成核、生長、塌陷回融等過程。冰晶始于冰晶核,隨著晶體生長條件的不同,晶體的形核可以分為均質成核或者異質成核。根據晶體成核理論,相變驅動勢主要與飽和水蒸氣壓力及其水蒸氣分壓力相關，冷凝驅動勢為同一溫度下水蒸氣濃度與飽和水蒸氣濃度之差。濕空氣冷凝為水的相變速率[10]如下：

(1)

其中:mwa為濕空氣與水之間的冷凝相變率,單位kg/(m2·s);ρa為濕空氣的密度,單位kg/m2;D為水蒸氣擴散率,單位m2/s;δ為單元網格中心至壁面的距離,單位m;Cphase為相變系數(shù)[11];ωsat為水蒸氣飽和質量分數(shù);ωwapor為水蒸氣質量分數(shù)。

霜氣界面上水蒸氣會越過飽和達到過飽和狀態(tài),霜面上過飽和水蒸氣與霜域內水蒸氣形成濃度差,霜域內不同位置處的濃度差導致了霜的不均勻分布。從結霜驅動勢ωsup-ωwapor角度出發(fā)修正冷凝相變速率,濕空氣中的水蒸氣以一定的擴散率進行質量傳遞,則單位時間內由濕空氣相向冰相的質量轉移速率為:

(2)

(3)

其中:mia為濕空氣與冰之間的相變率;Cphase為相變系數(shù),取106;ωsat為水蒸氣飽和質量分數(shù);ωsup為霜面上水蒸氣質量分數(shù);τ為修正系數(shù),取0.1。

過飽和度的定義如下：

(4)

文獻[6]通過對實驗數(shù)據分析,得到霜表面過飽和度的經驗式，即

(5)

本文使用的工況滿足該經驗式。在霜-氣界面上滿足S=Sfs,獲得霜面上的水蒸氣壓力，即

(6)

其中:Pv為水蒸氣分壓力;Pvs為飽和水蒸氣分壓力;Sfs為霜氣界面上濕空氣過飽和度;Pv,∞為來流水蒸氣分壓力;Pvs,∞為來流飽和水蒸氣分壓力;Pvs,fs為霜面飽和水蒸氣分壓力。

由水蒸氣壓力和理想氣體狀態(tài)方程確定水蒸氣的密度,進而確定水蒸氣的濕度,最后將濕度轉變?yōu)樗魵獾馁|量分數(shù)代入(2)式,便可計算濕空氣轉換為冰的相變率。

2.2 初始條件及邊界條件

開始時,計算域充滿了濕空氣,濕空氣的體積分數(shù)為1,冰相的體積分數(shù)為0,濕空氣密度設為不可壓縮理想氣體,定壓比熱容、熱導率、黏度均按照混合定律設置,水蒸氣的擴散率通過編程導入到Fluent中。其他物性參數(shù)均為軟件默認值。

濕空氣入口為速度邊界條件,出口為壓力邊界條件,冷板為恒溫邊界條件,上壁面為絕熱邊界條件。

3 結果與討論

3.1 模型驗證

霜層平均厚度和霜質量是衡量霜層生長的重要指標,將其模擬數(shù)據與實驗數(shù)據進行對比,完成本模型的驗證。模擬工況選取2種工況,模擬階段選取20～60 min,對比結果如圖2所示,所用的實驗數(shù)據均參考文獻[12]。

在模擬過程中,采用霜-氣界面的位置變化來追蹤霜層的增長,可以得到在不同時刻霜層的平均厚度,計算公式如下:

(7)

霜層質量的計算公式如下:

(8)

其中：j為濕空氣流動方向上的單元網格數(shù);N為濕空氣流動方向上全部網格數(shù);L為冷板長度,單位m;Δx為單元網格的長度,單位m;y(t)為霜面上網格縱坐標,單位m;k為霜層高度方向冷板至霜界面之間的網格數(shù);ρi為冰晶密度,單位kg/m2;A為每個單元網格的面積,單位m2;B為鋁板的寬度,單位m;S為霜區(qū)域的面積,單位m2。

實驗中提取t=60 min時霜的質量。在工況為ua=0.92 m/s、Ta=275 K、RH=85%、Tw=265 K條件下霜質量為0.302 g,模擬得到的霜質量為0.259 g;在工況為ua=0.92 m/s,Ta=275 K、RH=85%、Tw=260 K條件下霜質量為0.412 g,模擬得到的霜質量為0.404 g；2種工況下,誤差均在15%以內。由圖2可知，霜層平均厚度模擬與實驗結果的誤差均在30%以內，說明該模型可以有效預測結霜特征,可用于后續(xù)分析。

圖2 霜層平均厚度模擬結果和實驗結果對比

3.2 濕空氣溫度分布

霜的存在增加了濕空氣與冷板之間的熱阻,濕空氣溫度的分布制約著霜層形態(tài)。在t=2 h時不同截面上濕空氣溫度分布如圖3所示。霜-氣界面將計算域分為濕空氣域和霜域,傳熱方式不同導致2個區(qū)域內溫度變化速率的不同。

冷板前沿結霜驅動勢最大,故冷板前沿霜密集度較大,而霜的導熱系數(shù)是霜密度的函數(shù),因此冷板前沿溫度變化較大。

圖3 t=2 h時不同截面上濕空氣溫度分布

在冷板截面x=20 mm上不同時間的濕空氣溫度分布如圖4所示。從圖4可以看出,隨著時間增加,霜氣界面上濕空氣溫度升高,這會引起界面上過飽和水蒸氣質量分數(shù)降低,結霜驅動勢下降,最終導致霜層形態(tài)發(fā)生改變,霜層在后期增長緩慢。

抗菌藥說明書［適應癥］不符合抗菌藥說明書撰寫技術指導原則的主要表現(xiàn)有：適應癥沒有按照“本品適用于治療由對本品敏感的XXX、XXX和XXX菌引起的YYY病。”的規(guī)范描述；沒有遵循“如果獲得的證據僅僅支持用于較大人群的亞群（例如，疾病輕微的患者或特殊年齡組的患者）應予說明”的規(guī)定；沒有遵循“在某些情況下有理由限制適應癥，例如，建議藥品不作為某種感染的一線治療”應予描述的規(guī)定；遺漏使用限制的內容。

圖4 x=20 mm截面上濕空氣溫度分布

3.3 霜層形態(tài)

霜層形態(tài)反映霜的時間和空間特性。霜-氣界面上冰晶體積分數(shù)的變化率反映霜厚度增長速率;霜域內冰晶體積分數(shù)增長速度反映霜密實度增長速率。因此冰晶體積分數(shù)是表征霜層形態(tài)的重要指標。

冷板不同位置冰晶體積分數(shù)的變化如圖5所示，選取了距冷板0.1 mm,并且冷板截面在10、20、30、40、50、60、70 mm的7個位置。從圖5可以看出：冷板前端冰晶體積分數(shù)較大；在2 h之內,冰晶體積分數(shù)增長緩慢；在冷板截面10 mm位置處,冰晶體積分數(shù)在2～2.5 h內從0.086 32增長到0.150 40,增加了74.2%,此階段冰晶體積分數(shù)增長迅速。

圖5 冷板不同位置冰晶體積分數(shù)隨時間的變化

水蒸氣濃度差在結霜過程中占據主導地位。較高的溫差和較大的水蒸氣濃度差均有利于霜層生長。冷板前沿霜層較厚,霜層較密集,這一現(xiàn)象與結霜驅動勢有關,如圖5所示,冰晶體積分數(shù)在冷板截面20～30 mm之間下降較快,這是因為結霜驅動勢在此階段內減小速率最大。

2.5 h冷板不同位置結霜驅動勢的變化如圖6所示。從圖6可以看出,t=2.5 h時霜氣界面過飽和水蒸氣的質量分數(shù)以及在距冷板高度0.1 mm各截面的水蒸氣質量分數(shù)在截面10 mm處水蒸氣濃度差最大,霜層較厚,霜層較密；之后水蒸氣濃度差開始下降,且在截面20～25 mm之間濃度差下降15.6%,下降幅度最大,從而導致截面20～30 mm之間冰晶體積分數(shù)差異較大。

圖6 t=2.5 h冷板不同位置處結霜驅動勢的變化

3.4 霜物性參數(shù)

在不同的溫度和過飽和度區(qū)域內,霜的構造不同,霜的密度和導熱系數(shù)也不同。文獻[13]、文獻[14]結合實驗數(shù)據分別提出了在霜密度小于500 kg/m3霜層導熱系數(shù)的實驗關聯(lián)式。這種實驗關聯(lián)式更加貼近實際情況。本文采用文獻[13]關聯(lián)式計算霜層導熱系數(shù)。霜是由濕空氣和冰組成的多孔物質,采用加權平均方法計算可得:

ρfr=ρiαi+ρa(1-αi)

(9)

(10)

隨著霜面上水蒸氣過飽和度的下降,霜域內傳質速率大于霜面上傳質速率,水蒸氣主要擴散在霜層內部,孔隙率下降,霜的密度增長加快。

霜層密度和導熱系數(shù)隨時間的變化如圖7所示。在t=2 h以內密度增長緩慢,在t為2～3 h之間霜密度增長速度加快,從58.8 kg/m3增長到116.08 kg/m3,導熱系數(shù)增加了92.4%。這符合霜導熱系數(shù)與霜結構形態(tài)有關的理論,如圖5所示,霜密度分布與冰晶體積分數(shù)分布類似,因此冰晶體積分數(shù)是影響霜密度和導熱系數(shù)的重要因素。

圖7 不同時刻下霜密度和導熱系數(shù)

4 結 論

本文基于結霜驅動力建立了濕空氣轉變?yōu)樗南嘧兟时磉_式,數(shù)值模擬了水平冷板上霜生長與致密化過程；將霜層的平均厚度和霜層質量模擬數(shù)據與實驗數(shù)據進行對比,驗證該模型的合理性。用該模型分析了霜層形態(tài)和霜物性參數(shù)的變化,得到以下主要結論:

(1) 將霜層平均厚度、霜質量模擬數(shù)據與實驗數(shù)據對比,霜層平均厚度誤差在30%內,霜層質量誤差在15%內,驗證了本文模型的可靠性。

(2) 冷板不同截面上均存在霜-氣界面點,該界面將計算域分為濕空氣域和霜域,2個區(qū)域內傳熱方式不同。隨著時間的增加,霜-氣界面上濕空氣溫度升高,這會引起界面上過飽和水蒸氣質量分數(shù)降低,從而導致結霜驅動勢下降,霜層形態(tài)發(fā)生改變。

(3) 結霜初期,霜面上傳質速率處于主導地位,達到時間臨界點后,霜域內傳質速率大于霜面上傳質速率，這決定了結霜初期水蒸氣主要用來增加霜層厚度,后期水蒸氣主要用來增加霜的密度。霜密度與冰晶體積分數(shù)分布類似,因此冰晶體積分數(shù)是影響霜密度和導熱系數(shù)的重要因素。

本文主要工作是提出霜層生長與密實化數(shù)學模型,使用冷板結霜驗證模型的合理性,并利用該模型分析結霜機理,為換熱器結霜和除霜研究提供一定的理論基礎。