CO2在微通道內高速流動特性的研究

王亞賓,翁建華,崔曉鈺

(1.上海電力大學 能源與機械工程學院,上海 200090;2.上海理工大學 能源與動力工程學院,上海 200093)

1 引言

微通道的新型冷卻技術正在不斷發展。微系統技術通過對電子器件的微型化、集成化的方式實現了在有限空間內完成一系列復雜繁瑣的設備運行程序,如芯片級集成微系統等。從電子器件發展歷程來看,隨著電子器件的不斷微型化,散熱問題已成為電子器件研究制造環節所要解決的必要問題。設備尺寸減小以及功率的增加,當前一些用于電子設備的冷卻裝置已經不能滿足人們的需求,如何在系統整體功能不變的前提下,使得存儲空間更小、散熱更快是這個領域目前所要突破的關鍵[1]。

微型焦耳-湯姆孫效應制冷器因其具有結構簡單、重量輕、價格低、制冷速度快等優勢,被廣泛應用于電子器件的冷卻。考慮到焦湯效應制冷器的冷卻能力小,She Hailong等[2]研究了一種具有平行微通道的多層印刷電路板J-T,以氬氣為制冷氣體,入口壓力為8.02 MPa時的冷端溫度可達155.9 K,總制冷量3.66 W,忽略軸向熱傳導時冷端溫度可達154.7 K;S. P. Narayanan等[3]從制冷器設計的角度對實際微型制冷機熱交換器的控制方程進行數值求解,設計出一種回熱段采用曲折型設計的單層微通道制冷器來增強其換熱能力。

此外,很多學者對微通道流動特性展開研究。J. C. Harley等[4]使用硅晶片刻蝕的微通道對氮氣、氦氣和氬氣的流動進行了研究,所測得的摩擦因數與等溫假設下一階滑移邊界條件局部完全發展段的理論預測結果非常一致;M. J. Kohl等[5]將入口壓力分別為49.8、97.5和153 psig的空氣經過用3個硅芯片組成的矩形微通道,測得的壓力分布與數值預測吻合;C. Hong等[6]通過使用停滯壓力、局部壓力及質量流量獲得氮氣通過微通道時的局部馬赫數、溫度和摩擦因數,在范諾流假設的前提下局部范寧摩擦因子幾乎與布拉休斯的相關性重合。還有一些學者表明微尺度管段計算方法與宏觀管段計算不同。W. Peiyi等[7]使用氮氣、氬氣以及氦氣分別測量了硅和玻璃微通道內氣體流動時各管段截面的摩擦因數,結果表明,與宏觀管段相比微通道內的摩擦因數要高出10%～30%;R. Raju等[8]以N2和He為工作流體建立流體力學模型,研究了微通道內兩種不同情況的高速可壓縮氣體的流動傳熱特性,其結果與Monte Carlo模擬的結果存在一些明顯的差異;G. L. Morini等[9]研究相對粗糙度、可壓縮性和通道長徑比對轉捩雷諾數的影響時發現在沒計算失誤的情況下得到的臨界雷諾數與期望值不符。

上述學者研究的差異表明,對微通道內流體的流動特性研究是必要的。環境友好型工質CO2與制冷劑R134a和R1234yf相比,具有無毒、化學性質穩定、蒸發潛熱大、單位容積制冷量高等特點[10]。此外,CO2與微通道的結合提供了重量輕、壓實度高和傳熱系數高的CO2熱交換器[11]。以CO2為工質,建立矩形微通道流動模型,在絕熱流動與非絕熱流動的假設下對微通道的流動特性進行研究。

2 流動模型的確定

在微尺度流動中,當微通道的幾何尺寸遠大于流體分子的平均自由程時,常規尺度下的模型和方程仍然可以使用,但當流體流動的微通道的幾何尺寸很小時,不同的影響因素對其內部流體流動影響發生改變,從而導致在微通道中的流動狀況產生變化。在微尺度下的流動,還存在壁面的滑移效應、稀薄效應以及可壓縮性的影響。引入一個無量綱參數努森數(Kn)來劃分微尺度流動與常規尺度流動[12]。Kn數定義為分子平均自由程和流動特征尺寸之比,可由式1計算:

(1)

式中,λ是氣體分子平均自由程,單位為m;d是水力直徑,單位為m。

氣體分子平均自由程λ可由式2計算:

(2)

式中,μ是線性系數;R是為氣體常數,單位為J/(kg·K);T是流體局部溫度,單位為K;ρ是流體局部密度,單位為kg/m3。

根據Kn數的大小將氣體劃分為4個流動區域,對應的流動模型及流動尺度見表1。

表1 不同Kn數下流動尺度及流動模型

當Kn數小于10-3時,氣體分子之間的碰撞頻率遠高于氣體分子與壁面原子之間的碰撞頻率。經計算,研究的5種工況下,CO2氣體在微通道凹槽不同位置處的Kn數均遠小于10-3,因此流動區域可視為連續介質區,N-S方程和Fourier’s Law均適用,氣體與壁面交界處不存在速度滑移。

3 理論計算及數據分析

3.1 理論計算方法

對微通道進行分段計算,將長為L的微通道均勻分段,每一小段的微通道長度為l=L/n,n為段數。給定入口壓力和溫度時根據熱物性表可以查得入口截面處CO2的密度ρ0、比焓h0、氣體運動粘度v0、當地音速c0以及比熱比γ等。

CO2通過矩形微通道的速度可由式3計算得出:

(3)

式中,u是氣流速度,單位為m/s;s是質量流量,單位為g/s;ρ是密度,單位為kg/m3;A是截面面積,單位為m2。

馬赫數Ma由式4計算得出:

(4)

式中,γ是比熱比;R是氣體常數,單位為J/(kg·K);T是溫度,單位為K。

當Re<105時,沿程阻力系數f′可選用布拉休斯公式f′=0.316 4Re-0.25來計算,但是隨著Ma的增大,等溫流假設下的得到的摩擦因數f越來越偏離布拉休斯的表述,于是要引入對摩擦因數的修正,根據C. Hong等[13]的研究,修正后的摩擦因數fd可由式5計算得出:

(5)

假設第1段管段出口密度為ρ1,通過式6計算管段壓降:

(6)

第1段出口氣體壓力p1可由式7計算:

p1=p0-Δp1·10-6

(7)

式中,p0是管段進口壓力,單位為MPa;p1是管段出口壓力,單位為MPa。

第1段出口比焓可由式8計算:

(8)

式中,h0是入口CO2比焓,單位為kJ/kg;h1是出口CO2比焓,單位為kJ/kg;u0是入口速度,單位為m/s;u1是出口速度,單位為m/s。

圖1 計算流程框架

根據上述計算得出第1段出口參數,計算第2段、第3段直至第n段的出口參數,得到最終微通道出口端的溫度與壓力,進而分析各微管段流動的影響因素,迭代過程由MATLAB完成。

對于換熱量已知時的非絕熱壁流動,從能量平衡角度分析流動過程, 整個微通道視為一個整體取微管段進行分析。水平管道不考慮位能損失,即g(z2-z1)=0。根據式5和式8,絕熱壁q=0時,截面處的摩擦因數和比焓均是截面流速的函數。對于非絕熱壁面,外界對系統做功使得氣體焓值增加或降低,當槽道分段很多時,出口焓值就可近似由式9計算得出。假設凹槽受熱均勻,相對于絕熱壁流動,每個微管段焓值改變,其他流動參數隨之改變,通過迭代計算可得到換熱量已知時的非絕熱流動的出口參數。

(9)

式中,q是換熱量,單位為kJ/kg。

3.2 計算結果與分析

3.2.1 絕熱流動的計算結果與分析

以C. Hong等的研究為例驗證上述方法的準確性,其使用通過30%的KOH溶液刻蝕硅片與玻璃板結合的微通道,對氮氣在微通道的流動特性進行了研究,所采用的一種微通道長26.9 mm、寬1 020 μm、高112.7 μm,水力直徑為203 μm。取入口壓力Pin分別為0.19、0.56和1.01 MPa,入口處馬赫數Ma分別為0.27、0.32和0.34,入口溫度t分別為292.5、289和288 K。

將微通道劃分為2 690段,用該理論計算方法計算得到的Ma如圖2所示。將其通過實驗測得的數據與圖2進行對比,從圖2中可以看出,數據與實驗得到的數據十分吻合。

圖2 Ma沿流動方向的變化

調整微通道長度后對壓力和溫度進行計算,壓力與溫度的變化如圖3和圖4所示。在距入口端23 mm處,進口壓力為0.19 MPa的N2工況壓力降低至0.116 MPa,溫度從292.5 K降至285.66 K,降低了6.84 K;進口壓力為0.56 MPa的N2工況壓力降低至0.318 4 MPa,溫度從289 K降至277.41 K,降低了11.59 K;進口壓力為1.01 MPa的N2工況壓力降低至0.573 2 MPa,溫度從286 K降至272.84 K,降低了13.16 K,所得數據與試驗結果十分吻合。

圖3 壓力沿流動方向的變化

圖4 溫度沿流動方向的變化

對CO2在矩形微通道內的流動特性進行計算,槽道長L=40 mm,寬w=0.7 mm,高H=0.7 mm,水力直徑D=0.7 mm,入口壓力分別為0.3～0.7 MPa,質量流量s=0.5 g/s,溫度Tm=14 ℃,分析不同入口壓力對CO2流動參數的影響。

通過改變入口壓力,保持溫度及質量流量不變,5種狀態的馬赫數變化如圖5所示。可以看出,在進口低壓時相應的馬赫數反而越大,這是由于溫度不變,壓力降低,流體密度較低,從而質量流量恒定時流速較高,相應的入口處馬赫數越大。沿流動方向的摩擦因數越來越大,到達某一位置后流動達到臨界狀態發生壅塞。對于每個給定入口狀態,馬赫數都有確定的極限長度Lcr。馬赫數沿流動方向逐漸增加。槽道入口壓力和溫度恒定,質量流量相同,入口壓力越低,馬赫數在流動過程中增加越快。

圖5 五種狀態的Ma沿流動方向的變化

管段壓力沿流動方向的變化曲線如圖6所示。沿流動方向壓力不斷降低。入口壓力較低時的氣流速度較高,由式5可得,此時摩擦力相對于入口壓力較高的管段更大,使得在相同管段長度下的壓降更大,同時更容易速度達到臨界,管段發生壅塞。

圖6 壓力沿流動方向的變化

溫度沿流動方向的變化曲線如圖7所示。氣體溫度沿通道長度逐漸降低。在槽道未發生壅塞時,入口壓力為0.3和0.4 MPa的冷端溫度分別可達-9.91和-16.45 ℃,管段長度相同時,相對于另外3種工況制冷效果更為顯著。經計算,適當增加管段長度,使得入口壓力為0.5、0.6和0.7 MPa的流動接近臨界狀態,這時的冷端溫度分別可達-19.7、-22.96和-25.82 ℃,此時對應的極限長度分別為0.078、0.135和0.207 m。可見高壓工況的冷端溫度更低,但需要更長的管段去實現。低壓工況需要考慮管段壅塞,在管段長度一定時,低壓工況的溫降變化更為顯著。

圖7 溫度沿流動方向的變化

微通道高度H分別為0.7、0.35和0.175 mm時,質量流量調整為s=0.15 g/s,其他參數不變。此時不同凹槽高度時溫度的變化曲線如圖8所示。凹槽出口冷端溫度分別為13.85、13.39和-5.86 ℃,可以看到對于尺寸越小的微通道凹槽,溫降變化越大,相對于尺寸較大的凹槽,制冷效果更為顯著。

圖8 不同槽道高度時溫度的變化

3.2.2 非絕熱流動的計算結果與分析

非絕熱壁以CO2入口壓力0.5 MPa為例,槽道長L=40 mm,寬w=0.7 mm,高H=0.7 mm,水力直徑D=0.7 mm,質量流量s=0.5 g/s,溫度Tm=14 ℃。相對于絕熱壁工況,出口參數見表2。可以看出,外界對系統輸入熱量為正時,出口壓力降低,比焓增加,密度降低,流速增加,馬赫數增加;反之,輸入熱量為負時,出口壓力升高,比焓降低,密度增加,流速降低,馬赫數降低。

表2 不同壁面換熱功率時的出口流動參數

如果以CO2氣體進口壓力為0.3和0.4 MPa的工況進行非絕熱壁計算,槽道內流體達到臨界狀態時的位置有所變化。要使槽道受熱全部作用于有效管段,就需要提前設計槽道長度,使得氣體在不發生壅塞的前提下通過整個槽道。此時入口壓力為0.3和0.4 MPa時的非絕熱壁換熱極限管長見表3。由表3可知,加熱時極限長度減小,冷卻時極限長度增加,換熱條件下的極限管長長度隨換熱量而改變。

表3 不同壁面換熱功率時的極限管長

4 結語

本文針對不同入口壓力時的CO2氣體高速通過矩形微通道的流動特性進行研究,使用迭代計算得到冷端溫度,進而分析影響冷端溫度的因素。根據熱物性表,使用壓力與比焓擬合流動參數并用實驗數據驗證了理論計算方法的合理性。然后對給定質量流量的工況進行分析,結果表明,入口低壓時溫降變化更為顯著。微通道凹槽尺寸對冷端溫度的影響也非常重要。在換熱條件下的極限管段長度隨換熱量的不同而發生改變,在避免發生壅塞的前提下適當增加管段長度可以得到更低冷端溫度。所得結論總結如下:1)理論計算結果與試驗文獻結果十分吻合,驗證了理論計算方法的合理性;2)質量流量相同時,低壓工況溫降更大,但容易發生壅塞;3)槽道尺寸越小,溫降變化越為顯著;4)微通道內的高速流動,可以通過設計槽道尺寸、調節換熱量的方式得到理想的冷端溫度。