李精精,黃政,陳巧艷
中國核電工程有限公司,北京 100840
“ 華龍一號” ( Hua-long pressurized reactor 1000,HPR1000)非能動安全殼熱量導出系統(passive containment heat removal system,PCS)用于在設計擴展工況下安全殼的長期排熱,包括與全廠斷電、噴淋系統故障相關的事故。在電站發生設計擴展工況時,將安全殼壓力和溫度降低至可以接受的水平,并保持安全殼的完整性。在設計擴展工況,特別是嚴重事故情況下,PCS 作用可能會造成安全殼內呈現局部熱工參數分布不均勻的情形,此時傳統的集總參數程序無法進行準確模擬。
近年來隨著兩相計算流體力學(computational fluid dynamics,CFD)技術的發展,采用CFD 技術開展含不可凝氣體的蒸汽冷凝方法的研究日益增多[1?8]。然而,對于較為復雜的核動力系統或設備,采用計算流體力學程序進行計算時會消耗大量的計算資源。因此,許多學者開展了1D 系統程序與3D 流體力學程序的耦合研究[9?14]。Bezlepkin[15]針對VVER-1 200 電站的PCS 開展內部流場計算,其研究中將PCS 等效為定溫邊界,研究發現將傳熱管等效為等面積的平板不會導致換熱表面明顯的流場變化。邊浩志[16]建立了擴散邊界層冷凝模型,并通過修正Suction 效應,完成了豎直管束周圍冷凝傳熱分析。目前,關于含不凝性氣體蒸汽冷凝研究,多以定壁溫作為邊界條件,這與PCS 運行的實際工程情況不符,且目前關于PCS 程序與CFD 程序耦合計算研究的內容還較少,而PCS 運行狀態與殼內熱工水力狀態又是相互影響的耦合效應。因此,本文基于兩相流漂移模型開發HPR1 000 PCS 一維自然循環瞬態程序,開展PCS 程序與STAR-CCM+程序的耦合計算方法研究,分析PCS 程序作用下的局部熱工水力行為。
PCS 內流體熱工水力參數的變化遵循質量、動量和能量守恒的基本規律,再輔以流體狀態方程及相應的輔助方程可構成封閉方程。在漂移流模型中動量方程表示的壓降除摩擦阻力壓降、重力壓降和加速壓降外,還要考慮漂移流壓降梯度。在PCS 程序開發過程中所采用的數學物理模型參見文獻[17]。
在PCS 程序中,使用的PCS 換熱管管外換熱系數公式為
進而可以計算冷凝率為
式中:a、b、c、d為實驗獲得的系數;Ns為水蒸汽的體積分數;Ci為冷凝量修正系數;Tb為流體溫度,℃;Tw為壁面溫度, ℃;Pt為總壓,Pa; ?tw為壁面過冷度 ,℃;A為PCS 換熱管面積,m2;Hlg為汽化潛熱。
雙組份混合物的擴散系數[18]計算為
式中:T0、P0分別為標準狀態下的溫度和壓力,分別為273.15 K 和101 325 Pa;a的取值范圍為1.5~2.0,該處取1.81;D0為標準狀態下的擴散系數,取值為2.56×10?5m2/s。
在1∶1 比例的PCS 功能驗證實驗中,PCS 換熱器采用雙排豎管布置方式,C 形傳熱管由上、下聯箱連接,如圖1 所示。

圖1 實驗裝置示意
圖1中,冷凝罐組件直徑4.5 m、直筒段高度6 m、自由容積約127 m3。實驗開展了不同壓力和氣體(蒸汽、空氣和氦氣)配比組合條件下的準穩態實驗工況,根據不同工況實驗情況,對PCS 程序開展驗證工作,計算值與實驗值之間的誤差如圖2 所示。

圖2 PCS 程序準穩態工況驗證
由圖2 可知,計算值與實驗值之間的絕大部分誤差落在±20%范圍內。其中,誤差落在10%以內的有81 個,占總數據的62%;誤差落在10%~20%的有33 個,占總數據的25%;誤差超過20%的有17 個,占總數據的13%。對實驗數據進行分析發現,誤差超過20%的工況均為低壓工況。當冷凝罐內壓力較高時,由于水蒸汽的份額和混合氣體的溫度較高,換熱器管外具有較高的冷凝換熱系數,使系統有能力在回路出口附近維持穩定的閃蒸進程,因而系統的流量比較大,且流動穩定;而隨著冷凝罐壓力的降低,不可凝性氣體的相對份額升高,換熱器的冷凝換熱系數顯著下降,導致自然循環流動逐步發生周期性波動,且波動周期越來越長,回路中呈現兩相流–單相流交替流動狀態,從而增加了實驗測量誤差和計算誤差。事故狀態下安全殼內壓力、溫度迅速上升,并維持在較高的溫度、壓力范圍內。HPR1000 PCS 啟動壓力為0.24 MPa,所開發的PCS 程序能夠較為準確地模擬事故狀態下PCS 的運行。
引導基金是由政府設立并按照市場化方式運作的政策性基金,不以營利為目的,通過財政性資金投入,引導社會資本支持科技型企業發展,促進科技成果轉化和產業化,全面提升科技型中小企業的創新能力。引導基金的資金來源包括市級財政專項資金,引導基金資金存放銀行或者購買國債所得收益,引導基金投資退出返回的本金及收益,個人、企業或者社會機構無償捐贈的資金等。
在CFD 計算中,含不凝性氣體蒸汽冷凝模型的處理是極其重要的。一般認為,冷凝發生在近壁面第一層網格內,因此,可以通過在壁面第一層網格內或者壁面邊界設置冷凝質量通量的方式完成冷凝的計算。需要注意的是,在計算求解時,需判斷壁面溫度對應的飽和蒸汽壓力是否小于近壁面水蒸汽分壓,同時需判斷壁面溫度是否小于近壁面水蒸汽溫度,如滿足以上條件則判定會發生冷凝。King[19]對平板冷凝模型開展了計算分析,并給出了理論值,分析模型如圖3 所示。本文即采用在壁面第一層網格設置冷凝通量的方法對Sparrow 等的計算開展CFD 計算分析。

圖3 平板冷凝模型
對圖3 所示的模型開展二維計算,取平板溫度為定壁溫330 K,入口流速為0.5 m/s,壓力為101 325 Pa,入口水蒸汽體積分數分別為0.5 和0.75、水蒸汽溫度分別為360 和369 K 這2 種工況。本文通過CFD 計算得到的值與Sparrow 等給出的參考值對比如圖4 所示。

圖4 CFD 計算值與參考值對比
由圖4 可知,CFD 計算的壁面冷凝通量與參考值符合較好,只在入口1 mm 位置附近出現了計算值略低于參考值的情況。因此,在壁面第一層網格采用冷凝通量的方法可以很好地模擬壁面附近水蒸汽冷凝過程。
PCS 程序與STAR-CCM+程序之間需要完成數據的傳遞。首先,需要將STAR-CCM+計算得到的域內的溫度、壓力、氣體組分等信息傳遞給PCS 程序;PCS 程序計算完成后再將結果傳遞給STAR-CCM+程序。在STAR-CCM+程序PCS 壁面第一層網格或PCS 壁面,將PCS 程序計算得到的能量源項和質量源項以熱阱的形式減掉,具體如圖5 所示。也可只將PCS 程序計算得到的壁面溫度賦值給STAR-CCM+程序,PCS 的冷凝換量通過式(1)獲得。

圖5 PCS 系統程序與STAR-CCM+程序數據交換
為了完成以上過程的數據傳遞,需要使用C++程序編寫用戶自定義程序,該程序會從STAR-CCM+讀取所需的環境參數信息,并傳遞給PCS 程序,待PCS 程序計算完成后將計算結果傳遞給STAR-CCM+開始下一步迭代,從而完成數據交換。2 個程序的耦合計算發生在STARCCM+程序的每個時間步長內。
在耦合計算研究中,以換熱面積相當的原則,將PCS 換熱器等效為沿圓周分布的一系列管束。為了節省計算資源,取罐體的1/4 體積作為計算域,具體如圖6 所示。

圖6 耦合計算模型
假設罐體內初始全部為空氣,溫度為35 ℃。155 ℃的水蒸汽以0.25 g/s 的質量流量噴入罐體內,為了兼顧計算速度,STAR-CCM+程序每迭代10 步調用1 次PCS 程序。
PCS 冷凝速率和罐體內氣體體積分數如圖7所示。

圖7 PCS 冷凝速率及罐體內氣體組分
由圖7 可知,在初始時刻罐體內只有空氣,隨著水蒸汽的不斷注入,罐體內水蒸汽體積分數逐漸升高,最終罐體內空氣和水蒸汽體積分數達到平衡狀態。與此同時,PCS 的冷凝速率在開始時刻為0,約20 s 后PCS 開始出現冷凝速率并逐漸增大,最終達到比蒸汽注入速率略低的水平并維持在動態平衡狀態。
分別對有PCS 作用和無PCS 作用下罐體內的溫度、壓力開展計算如圖8 所示。

圖8 殼內平均壓力和溫度變化曲線
由圖8 可知,隨著水蒸汽的噴入,罐體內壓力、溫度迅速上升。考慮PCS 的工況下,罐體內的溫度和壓力可以得到有效控制,而對比工況中,PCS 不動作時罐體內的壓力、溫度維持持續上升趨勢。考慮PCS 作用情況下,罐體內的溫度最終維持在約142 ℃,罐體內的壓力維持在約0.33 Pa。
輸出PCS 周圍的溫度場和水蒸氣體積分數場如圖9 所示。

圖9 PCS 周圍水蒸汽濃度場和溫度場
由圖9 可知,在水蒸汽冷凝和壁面對流換熱雙重作用下,近壁面網格內溫度最低,向外擴展溫度慢慢升高。PCS 近壁面的水蒸汽濃度較低,且靠近噴放口的一側PCS 近壁面水蒸汽體積分數低于靠近容器壁一側。這主要是由于自罐體中心線位置噴入的蒸汽在上升過程中不斷卷吸周圍的氣體,待水蒸汽達到頂部并與頂部接觸后折返到達PCS 壁面周圍,再向下運動,促使罐體內的空氣在罐體下部沉積。
本研究針對HPR1 000 PCS 程序開展與STARCCM+程序的耦合計算方法研究,得到如下結論:
1)所開發的PCS 程序能夠較為準確地模擬事故狀態下PCS 的運行,計算值與實驗值的絕大部分誤差在±20%范圍內。
2)通過在傳熱管壁面附近第一層網格內設置冷凝通量的方案是可行的,且CFD 計算值與分析值符合較好。
3)本研究所建立的PCS 程序與STAR-CCM+程序耦合計算方法能夠用于后續HPR1000 及華龍后續機型的PCS 研究。