葉茂菁, 吳天宇, 劉鎧瑞, 王利民, 車得福
(西安交通大學能源與動力工程學院, 陜西 西安 710100)
作為一種新型的動力循環技術, 超臨界二氧化碳(S-CO2) 布雷頓循環系統具有設備少、 結構緊湊、 占地面積小、 熱慣性小、 效率高和運行靈活等特點[1-3]。 因此, S-CO2布雷頓循環在很多領域具有廣泛的應用前景, 如高溫太陽能發電站和增強型地熱系統, 以及環境工程[4-8]。 作為S-CO2布雷頓循環重點關注的問題, S-CO2在管內的傳熱特性影響著系統整體效率。 因此, 充分了解并掌握S-CO2傳熱機理對于換熱裝置設計和安全穩定運行至關重要。
Hall 和Adebiyi[9]對水平圓管內S-CO2流動傳熱進行了實驗研究, 經觀察圓管底部的傳熱能力較強, 而頂部的傳熱能力較差。 Yan 等[10]研究了SCO2在水平圓管頂母線和底母線區域不同的流動傳熱行為, 解釋了不同熱流密度和質量流速下頂母線內壁溫度分布產生差異的傳熱機理。 Zhang 等[11]研究了S-CO2在水平半圓形通道的耦合傳熱特性,并詳細分析了浮力引起的二次流對局部熱性能的影響, 發現二次流會導致熱邊界層變薄, 從而增強壁面局部換熱。 Li 等[12]研究了超臨界水在側壁加熱腔內的自然對流和溫度分布, 發現與理想氣體的單渦流動不同, 超臨界水在腔內表現出更薄的邊界層和更高的速度; 強流動邊界層會進一步導致雙渦流型和過熱現象, 也會影響自然對流和溫度分布。
綜上可以看出, 關于浮升力方向與流動方向垂直的流動傳熱研究主要集中于S-CO2不同形狀截面流動分層[13-14]、 局部強化傳熱機理[15-17]和其他超臨界流體在腔體中的自然對流, 很少涉及S-CO2僅在浮升力驅動下的腔內自然對流流動發展過程。
為研究水平圓管垂直軸向方向上浮升力驅動下的流動傳熱行為, 本文基于二維圓腔的物理模型,在恒定壁溫條件下研究超臨界二氧化碳自然對流的發展過程, 以及壓力、 管徑、 溫度區間對發展過程的影響, 為超臨界二氧化碳水平圓管非穩態流動傳熱提供參考。
為研究垂直于軸向方向浮升力驅動下S-CO2流動傳熱行為, 以及壓力、 管徑和溫度區間對流動傳熱行為的影響, 建立二維圓腔物理模型, 如圖1所示。

圖1 物理模型
圓腔直徑分別為2 mm、 4 mm、 8 mm, 重力方向豎直向下, 管壁為恒定溫度且與流體之間無滑移, 初始時刻流體與壁溫存在一定溫差。
具體工況邊界條件如下: 流體涉及擬臨界和遠離擬臨界溫度區間分別為300~360 K 和360~600 K。對于擬臨界溫度區間的工況, 壁溫為360 K, 初始時刻流體溫度為300 K; 而遠離擬臨界溫度區間的工況, 壁溫為600 K, 初始時刻流體溫度為360 K, 流體壓力分別為8 MPa、 9 MPa、 10 MPa。
1.2.1 控制方程
非穩態對流換熱的控制方程組[12]在直角坐標系中表示如下。
1) 連續性方程:
2) 動量方程:
3) 能量方程:
受熱壁面設置為恒定壁溫并設置為無滑移, 以流體溫度作為初始化溫度場的溫度。 采用有限體積法離散各項控制方程, 選擇層流模型, 壓力-速度耦合方程采用SIMPLEC 算法進行求解, 動量、 能量方程均采用二階迎風格式。 由于二維圓腔中壓降遠小于操作壓力, 故忽略壓力對物性的影響, 認為S-CO2的物性僅隨溫度發生變化。
為保證計算結果的準確性, 通過選取NIST 數據庫[18]中的CO2物性數據在Fluent 中進行插值分段線性擬合, 因為Fluent 中僅支持50 個點的物性分段線性擬合, 故在300 ~360 K、 360 ~600 K 兩個溫度區間內選取50 個點進行曲線擬合。 計算中對殘差進行監視以判斷計算是否達到收斂。
1.2.2 網格劃分及無關性驗證
使用ICEM 軟件建立二維圓腔的模型, 建模的過程對邊界層進行局部加密, 并使用O 型切分以改善網格正交性, 如圖2 所示。

圖2 網格劃分
在保證節點相對比例不變的情況下, 做出6 套網格進行網格無關性驗證, 以確保計算結果與網格數量無關。 驗證工況如下: 圓腔直徑4 mm, 壓力8 MPa, 壁面溫度330 K, 初始時刻流體溫度300 K。如圖3 所示, 隨著網格數增加, 測點速度曲線逐漸趨于重合。 當網格數目為160 564 時, 與網格數目為198 404 的工況相比, 測點速度的變化幅度在0.5%以內, 此時已滿足網格無關性條件。 因此,選用網格數為160 564 的網格進行模擬。

圖3 網格無關性驗證
1.2.3 時間步長無關性驗證
選擇網格數目為160 564 的網格進行研究, 設計時間步長為0.000 5 s、 0.001 s、 0.002 s、 0.004 s、0.005 s、 0.01 s 的6 組算例進行時間步長無關性驗證, 以確保計算結果與時間步長無關。 對比測點速度曲線隨時間步長的變化, 如圖4 所示。 隨著時間步長逐漸減小, 曲線逐漸趨于重合。 當時間步長為0.001 s 時, 與時間步長為0.005 s 的工況相比,測點速度的變化幅度在0.5%以內, 此時已滿足時間步長無關性條件。 因此, 選取0.001 s 作為時間步長進行仿真計算。

圖4 時間步長無關性驗證
為保證數值模擬方法的可靠性, 需進行模型驗證。 目前關于超臨界二氧化碳自然對流的研究報道較少, 因此采用普通流體空氣進行驗證[12]。 參考Barakos 和Mitsoulis[19]關于方腔自然對流的案例,根據文獻建立由80×80 的網格組成的方腔, 如圖5 所示。 在文獻[19] 中, 模型的大小即方腔的邊長是由Ra數計算得出的。 對于文獻中給出Ra=1.89×105, 對應方腔邊長為363.875 mm。 邊界條件如下: 上下邊界絕熱, 左邊界溫度為303 K, 右邊界溫度為283 K。 其中, 工質為空氣, 初始溫度為293 K。

圖5 方腔模型
不同作者G.Barakos、 Krane 和Jessee 與本研究的驗證結果對比如圖6 所示, 圖中W為寬度,H為高度。 通過對比三個無量綱參數在數值模擬與實驗條件的差異, 可以得出本研究模型與Barakos[19]模型的偏差較小, 偏差最大值在5%以內, 所以該模型符合層流自然對流[20]。 因此, 本文采用的模型可以用于S-CO2自然對流傳熱流動特性研究。

圖6 模型驗證
以壓力8 MPa、 管徑4 mm 的工況為例, 關注各個時間下, S-CO2流體在二維圓腔內擬臨界點附近流動時, 速度場和溫度場隨傳熱過程的變化情況, 并對擬臨界點工況與遠離擬臨界點工況的速度場和溫度場隨時間變化進行對比分析。
為了清晰地描述S-CO2流體在受熱發展中的變化過程, 按時間尺度將整個對流受熱過程分為3 個階段, 分別為初始階段、 發展階段和穩定階段。
2.1.1 速度場對比
從圖7 (a) 中可以看出, 對于溫度區間處于擬臨界點的工況, 在流動傳熱的初始階段, 流場中大部分區域的速度為0 m/s, 高溫流體從圓腔底部沿著管壁向上流動, 產生薄薄的速度邊界層。 到達圓腔頂部附近時, 流體脫離壁面向下流動, 且向壁面方向彎曲。 此時在圓腔頂部會產生很多細小的渦流, 圓腔底部有不規則的羽流出現。

圖7 初始階段二維圓腔自然對流速度場
而對于溫度區間遠離擬臨界點的工況, 如圖7(b) 所示, 高溫流體同樣從圓腔底部沿壁面向上流, 而后向下流動并向壁面彎曲, 形成環狀流, 在圓腔底部并沒有羽流的產生。
對比圖7 (a)、 (b) 可以發現, 溫度區間處于擬臨界點的工況, 在流動傳熱初級階段產生的速度邊界層要遠薄于溫度區間遠離擬臨界點工況產生的速度邊界層。 此邊界層不同于圓管內流動所產生的邊界層, 本文對應工況產生的邊界層具有如下特點: 速度從壁面開始沿徑向增長到最大值, 而后速度開始減小直至速度為0 m/s。 這是因為貼壁處速度為0 m/s, 且在靠近主流區域溫度不均勻作用消失, 速度為0 m/s, 所以在接近壁面的中間區域產生峰值。
當流動傳熱處于發展階段時, 如圖8、 圖9 所示, 對于溫度區間處于擬臨界點的工況, 彎曲向下流動的流體在圓腔內部脫離, 形成一個個不規則的渦流, 布滿圓腔內部, 并且隨著時間的增長, 流體在沿壁面流動時形成斷流, 內部渦流減少, 流場內流動混亂。 對于溫度區間遠離擬臨界點的工況, 此時流動邊界層內流動特征保持不變, 向下流動的流體開始脫落, 逐漸移動到圓腔底部。 并且隨著時間的增長, 圓腔內整體速度減小, 大部分區域速度為0 m/s, 但仍保持邊界層的特征。

圖8 發展階段二維圓腔自然對流速度場

圖9 發展階段二維圓腔自然對流速度場
從圖10 中可以看出, 當流動傳熱處于穩定階段時, 對于溫度區間處于擬臨界點的工況, 此時壁面無流體脫落, 且邊界層內速度整體減小, 僅在圓腔底部存在少量混亂的渦流。 而對于溫度區間遠離擬臨界點的工況, 此時邊界層內速度再次降低, 整個圓腔內大部分區域速度為0 m/s, 邊界層特征逐漸模糊。

圖10 穩定階段二維圓腔自然對流速度場
從二維圓腔自然對流速度場的整個發展過程,即初始階段、 發展階段及穩定階段可以看出, 對于溫度區間處于擬臨界點的工況, 流場僅在初始階段呈現對稱的趨勢, 且隨著時間的增長, 流場逐漸混亂。 這是因為在擬臨界區域物性隨溫度劇烈變化,導致在流動發展過程中流場相對混亂。 而對于溫度區間遠離擬臨界點的工況, 流場在整個發展階段都呈現對稱的趨勢。 因為溫度區間遠離擬臨界點的工況物性隨溫度變化較為平緩, 所以流場在整個發展過程中相對穩定。
2.1.2 溫度場對比
如圖11 所示, 當溫度區間處在擬臨界點附近時, 在流動傳熱的初始階段, 大部分流體并未發生明顯溫升, 高溫流體主要集中于圓腔頂端和近壁面區域; 同時圓腔頂端有部分高溫流體脫落, 進入主流區, 圓腔底部也有部分羽流擺脫壁面的黏性作用進入主流區, 此時溫度場呈對稱分布。 而對于溫度區間遠離擬臨界點的工況, 在流動傳熱的初始階段, 已有大部分流體發生明顯溫升; 高溫流體同樣主要集中于圓腔頂端和近壁面區域, 此時高溫流體沿著壁面由圓腔頂部向下彎曲流動, 當流體再次接近壁面時, 流體旋向改變流向低溫主流區。

圖11 初始階段二維圓腔自然對流溫度場
對比圖11 (a)、 (b) 可以發現, 溫度區間處于擬臨界點的工況, 在流動傳熱初級階段產生的溫度邊界層遠薄于溫度區間遠離擬臨界點工況產生的溫度邊界層。 因為導熱系數在360 K 附近較低, 而在遠離擬臨界區600 K 附近較高, 所以高溫流體向內傳遞熱量的速度更快, 導致溫度區間遠離擬臨界點工況產生的溫度邊界層更厚。
如圖12、 圖13 所示, 流動傳熱處于發展階段時, 對于溫度區間處于擬臨界點的工況, 溫度場上下層分明, 熱量逐漸從上部高溫流體傳向下部低溫流體; 高溫流體所占的面積隨發展過程擴大, 主流區流體溫度緩慢上升, 在圓腔底端仍有部分羽流擺脫壁面的黏性作用進入主流區。 而對于溫度區間遠離擬臨界點的工況, 在0.5 s 時圓腔內大部分流體已出現明顯溫升, 而到1 s 時溫度場已基本穩定,僅有少部分低溫流體存在于圓腔底端。

圖12 發展階段二維圓腔自然對流溫度場

圖13 發展階段二維圓腔自然對流溫度場
從圖14 中可以看出, 當流動傳熱處于穩定階段時, 對于溫度區間處于擬臨界點的工況, 圓腔內所有流體均出現明顯溫升, 但是圓腔下半部分流體溫升較慢, 圓腔上、 下半部分平均溫差較大。 而對于溫度區間遠離擬臨界點的工況, 在2 s 時溫度場完全穩定, 整個圓腔內溫度分布均勻, 保持在600 K。

圖14 穩定階段二維圓腔自然對流溫度場
從二維圓腔自然對流溫度場的整個發展過程,即初始階段、 發展階段及穩定階段可以看出, 對于溫度區間處于擬臨界點的工況, 溫度場僅在初始階段呈現對稱的趨勢, 且隨著時間的增長, 溫度場逐漸分層; 相較于速度場, 溫度場在整個發展過程中更加穩定, 熱量自上而下逐層傳遞。 而對于溫度區間遠離擬臨界點的工況, 溫度場在整個發展階段都呈現對稱的趨勢。 因為遠離擬臨界點的工況物性隨溫度變化較為平緩, 所以溫度場在整個發展過程中相對穩定, 并且熱量傳遞速度高于溫度區間處于擬臨界點的工況。
此外, 溫度區間處于擬臨界點的工況溫升為60 K, 遠小于遠離擬臨界區域工況的溫升240 K;而溫度區間遠離擬臨界點的工況下達到溫度場穩定所需要的時間快于溫度區間處于擬臨界區域工況下的所需時間, 說明CO2在擬臨界區域內物性的劇烈變化不利于流動傳熱。
不同壓力下S-CO2的物性存在很大差異, 所以壓力對圓腔內流動傳熱的影響不容忽視。 因此, 為研究壓力對圓腔內S-CO2流動傳熱的影響, 分別設置8 MPa、 9 MPa、 10 MPa 三種壓力, 并選擇溫度區間處于擬臨界點的工況來研究壓力對圓腔S-CO2流動傳熱的影響。 為便于比較分析, 選取過圓腔圓心水平方向(Y=0) 軸線上速度和溫度的分布來分析壓力對圓腔內流動傳熱的影響。
2.2.1 壓力對圓腔流動傳熱速度分布的影響
在流動傳熱發展的初始階段, 如圖15 (a) 所示, 三種壓力邊界下的工況速度邊界層的厚度相差不大, 壓力為8 MPa 的工況相較于壓力為9 MPa、10 MPa 的工況對應速度變化最為劇烈, 且整體速度略高于9 MPa 和10 MPa 的工況。 這是因為SCO2在8 MPa 下的擬臨界溫度比在9 MPa、 10 MPa的低, 相同壁面溫度下, 8 MPa 對應的工況更容易到達擬臨界點, 故壓力為8 MPa 的工況在初始階段速度變化最為劇烈且整體溫度相對較高。 在流動傳熱的發展階段, 如圖15 (b)、 (c) 所示, 近壁面處的速度邊界層依然顯著, 并且浮升力已作用到主流區; 主流區域速度分布復雜并出現多個峰值, 是此階段圓腔底部上升的羽流及圓腔頂部脫落的渦流作用于主流區所致。 如圖15 (d) 所示, 在流動傳熱進入穩定階段時, 邊界層、 主流區速度降低, 整體主流區速度場趨于平穩。

圖15 不同時間下二維圓腔Y=0 處速度分布
從圖15 中還可以觀察到: 0.2 s 時, 邊界層內速度快速上升, 主流區速度較低; 0.5 s 時, 邊界層內速度保持不變, 主流區速度上升且波動劇烈;1 s 時, 邊界層內速度保持不變, 主流區速度降低且開始保持穩定; 2 s 時, 邊界層內速度開始大幅降低, 主流區速度不變, 保持較低水平。 因為邊界層內流速與溫差、 流體的黏性力相關, 所以從初始階段至穩定階段邊界層內流速始終保持較高水平,直到2 s 流動傳熱充分進行時, 邊界層內的流速才開始降低。
2.2.2 壓力對圓腔流動傳熱溫度分布的影響
如圖16 (a) 所示, 在流動傳熱發展初始階段, 不同壓力條件下的工況均產生溫度邊界層, 由于壁面溫度高于管內流體的溫度, 近壁面處流體先受到壁面傳熱影響開始升溫, 產生溫度邊界層。 對于Y=0的溫度場, 還未影響到主流區, 故主流區域溫度保持在300 K。

圖16 不同時間下二維圓腔Y=0 處溫度分布
進入發展階段, 如圖16 (b)、 (c) 所示,0.5 s時, 邊界層溫度變化的影響傳遞到了主流區,同時受到圓腔底部向上發展的不規則羽流影響, 主流區產生不規則波動的溫度場, 溫度場在此時開始整體上升且具有多個峰值。 溫度場邊界層仍然保持, 且在溫度場邊界層與主流區的交界區域開始出現溫度極小值。 這是因為主流區熱量自上而下傳遞, 在Y=0 軸線上, 熱量傳遞的過程是高溫主流區向邊界層附近傳熱, 同時壁面向邊界層傳熱, 所以出現溫度極小值。
隨后在1 ~2 s 的時間內, 溫度場逐漸趨于穩定, 呈現出溫度自壁面開始降低, 達到極小值后開始升高, 而后保持穩定。 此時可以看出10 MPa 工況的傳熱發展過程明顯快于8 MPa、 9 MPa 工況,且壓力越高傳熱發展過程越快。 這是因為對于SCO2, 當溫度大于擬臨界點時, 壓力越大, 導熱系數越高, 所以傳熱發展過程所需要的時間越短。
管徑大小也是影響圓腔內流體傳熱的重要因素, 對三種不同管徑下的圓腔對流傳熱情況進行分析, 重點關注工況達到穩定發展狀態的時間, 即主流區溫度不再波動, 開始穩定上升。
如圖17 所示, 2 mm 管徑的傳熱過程發展最快, 在0.3 s 時已達到穩定傳熱狀態。 4 mm 管徑的工況達到穩定發展狀態需要0.8 s, 而8 mm 管徑的工況達到穩定發展狀態則需要2 s。

圖17 不同管徑下不同時間的二維圓腔Y=0 處
此外, 還可以觀察到這樣的現象: 管徑越大,主流區在初始發展階段波動越劇烈。 這表明管徑越大, S-CO2圓腔自然對流越劇烈。
對于S-CO2圓腔內自然對流傳熱, 管徑越大,傳熱達到穩定的時間越長, 主流區在初始發展階段波動越劇烈。
本文對S-CO2在二維圓腔中自然對流流動傳熱過程, 以及壓力、 管徑、 溫度區間對流動傳熱過程的影響進行了研究, 得到以下結論:
1) 對于溫度區間處于擬臨界點的工況, 流場在初始階段呈現對稱的趨勢, 但隨著時間的增長,流場逐漸混亂; 對于溫度區間遠離擬臨界點的工況, 流場在整個發展階段都呈現對稱的趨勢。 相較于速度場, 溫度場在整個發展過程中更加穩定, 熱量自上而下進行傳遞。
2) 溫度區間處于擬臨界點的工況, 在流動傳熱初級階段產生的溫度邊界層遠薄于溫度區間遠離擬臨界點工況產生的溫度邊界層。 這是因為在360 K附近導熱系數較低, 傳熱熱阻較高; 而在遠離臨界區600 K 附近導熱系數較高, 傳熱熱阻較低, 導致溫度區間遠離擬臨界點工況產生的溫度邊界層更厚。
3) 溫度區間處于擬臨界點的工況, 溫升為60 K, 遠小于遠離臨界區域工況的溫升240 K, 而溫度區間遠離擬臨界點的工況下達到溫度場穩定所需要的時間小于溫度區間處于擬臨界區域的工況下所需的時間, 說明CO2在擬臨界區域內物性的劇烈變化不利于流動傳熱。
4) 對于S-CO2圓腔內自然對流傳熱, 壓力越大, 傳熱發展過程越快。 這是因為對于S-CO2, 當溫度大于擬臨界點時, 壓力越大, 導熱系數越高,故傳熱發展過程所需要的時間越短; 對于S-CO2圓腔內自然對流傳熱, 管徑越大, 傳熱達到穩定的時間越長, 主流區在初始發展階段波動越劇烈。