VVER堆芯捕集器換熱特性數(shù)值模擬研究

朱光昱 郭 超 劉巧鳳 李 春 依 巖

1（生態(tài)環(huán)境部核與輻射安全中心 北京 100082）

2（中國核電工程有限公司 核電安全研究中心 北京 100840）

堆芯捕集器是核電廠嚴(yán)重事故緩解手段之一，該技術(shù)通過在反應(yīng)堆壓力容器（Reactor Pressure Vessel，RPV）外部設(shè)置的收集裝置對熔融物進(jìn)行冷卻并最終實現(xiàn)滯留。目前，國內(nèi)外采用堆芯捕集器設(shè)計包括：俄羅斯OKB Gidropress 設(shè)計開發(fā)的VVER［1-2］（Vodo-Vodyanoi Energetichesky Reactor）；由法國和德國聯(lián)合開發(fā)的改進(jìn)型壓水堆EPR［3-4］（Europe Pressurized Water Reactor）；由通用電氣公司和東芝聯(lián)合設(shè)計的新型沸水堆ESBWR（The Economic Simplified Boiling Water Reactor）［5］；由韓國電力公司設(shè)計的APR1400 及其歐洲版本EUAPR1400（Europe-Advanced Power Reactor）［6］。

上述設(shè)計中，VVER和ESBWR采用了坩堝式設(shè)計而EPR和APR1400采用了熔融物擴(kuò)展冷卻方式。其中，VVER 的堆芯捕集器［1-2］主體布置在壓力容器下方的堆坑中，坩堝外壁面設(shè)置有冷卻水通道。在事故發(fā)生后，首先通過犧牲材料與堆芯熔融物發(fā)生反應(yīng)，降低熔融物整體的體積釋熱率并改善熔融池流動性。我國學(xué)者結(jié)合上述設(shè)計提出了一種雙層坩堝形式的堆芯捕集器，通過增設(shè)冷卻管提升熔融物冷卻速率［7］。

EPR 堆芯捕集器［3-4］設(shè)計采用熔融物擴(kuò)展的方式增加換熱面積，從而使得堆芯熔融物快速冷卻。事故發(fā)生后首先通過堆腔內(nèi)的熔融物臨時滯留裝置，將分批轉(zhuǎn)移的熔融物進(jìn)行收集，再統(tǒng)一將熔融物鋪展在面積約為170m2的捕集器擴(kuò)展室中，隨后通過捕集器底面、側(cè)壁面，以及熔融物頂部淹沒的方式對熔融物進(jìn)行冷卻。EPR堆芯捕集器預(yù)計可將熔融物的時間縮短至幾天之內(nèi)，但系統(tǒng)運(yùn)行流程較多，熔融物擴(kuò)展流動的不確定性對系統(tǒng)設(shè)計功能的可靠性影響較大，而且堆芯捕集器擴(kuò)展室占據(jù)了大量的空間不利于安全殼內(nèi)的設(shè)備布置。

為了研究坩堝式堆芯捕集器中熔融物冷卻過程，本文采用多物理場耦合計算軟件COMSOL，結(jié)合VVER 堆芯捕集器設(shè)計資料搭建計算模型，對熔融物冷卻過程中的熔融池的溫度、結(jié)殼分布等參數(shù)進(jìn)行了數(shù)值模擬研究，探討了事故緩解進(jìn)程，為該系統(tǒng)的安全性評價提供了一定的依據(jù)。

2 計算模型設(shè)置

2.1 幾何與邊界條件

VVER 堆芯捕集器整體結(jié)構(gòu)如圖1所示，主要由位于RPV 下部的托盤和位于堆腔內(nèi)的換熱器組成。托盤內(nèi)設(shè)置有肋片以協(xié)助收集RPV 破口流出的熔融物，并引導(dǎo)熔融物向下部的換熱器遷移。換熱器的內(nèi)部設(shè)置有填料籃來裝載犧牲材料，其外部設(shè)置有冷卻水通道。

圖1 VVER堆芯捕集器結(jié)構(gòu)示意圖Fig.1 Structure diagram of VVER core catcher

圖2 為參考VVER堆芯捕集器的換熱器部分建立的軸對稱計算域模型。坩堝式換熱器主體為圓柱體，半徑約2.355m，底部為圓錐體，錐段與水平方向的傾角為12°。根據(jù)鄭濱［8］對VVER 堆芯捕集器內(nèi)熔融池的計算結(jié)果，在分層階段氧化物層和金屬層體積分別為26.1m3和4.6m3，由此得到模型中氧化物層高度為1.61m，金屬層在中軸處的高度為0.5m。熔融物的衰變熱全部集中在氧化物層中，事故發(fā)生后熔融池內(nèi)的衰變熱變化在文獻(xiàn)［8］中有詳細(xì)介紹。

圖2 VVER堆芯捕集器的換熱器軸對稱計算域模型示意圖Fig.2 Schematic of axisymmetric computational domain model for heat exchanger of VVER core catcher

計算過程中，考慮安全殼內(nèi)壓力維持在0.12MPa，在系統(tǒng)長時間運(yùn)行后認(rèn)為冷卻水過冷度很低，因此將坩堝外壁面溫度設(shè)置為接近飽和溫度的378K。通過附加薄層的方式考慮不銹鋼坩堝對導(dǎo)熱的影響，根據(jù)設(shè)計參數(shù)薄層厚度可設(shè)置為60mm。氧化物層頂部與冷卻水直接接觸可認(rèn)為長期處于沸騰換熱狀態(tài)，因此頂部設(shè)置為對流換熱條件。由于缺乏實際參數(shù)，換熱系數(shù)根據(jù)大容積沸騰米海耶夫公式估算［9］：

式中：h為沸騰換熱系數(shù)，W·(m2·K)-1；Tup和Tsat分別表示頂部邊界網(wǎng)格溫度和飽和溫度，K；P為絕對壓力，Pa。

2.2 模型設(shè)置

通過COMSOL 軟件相變模型計算熔融物的熔化凝固過程，該模型考慮物質(zhì)相變過程中固液相線溫度之間為糊狀區(qū)，通過插值方式計算糊狀區(qū)的密度、比熱容、導(dǎo)熱系數(shù)。在VVER堆芯捕集器熔融池中，金屬層中98%的質(zhì)量為Fe，其他為少量的金屬Al、U 和Zr［8］，因此，層內(nèi)物質(zhì)的相變潛熱和膨脹系數(shù)均設(shè)置為Fe 的物性。金屬層相變溫度區(qū)間也參照Fe的熔點設(shè)置，同時考慮了混合物熔點會適當(dāng)降低。氧化物層中熔融物的固液相線間溫差較大，在凝固過程中黏度會逐漸變化，因此對流動過程影響較大。EPR堆芯捕集器設(shè)計過程中對熔融物擴(kuò)展流動特性進(jìn)行了大量實驗研究［10-13］，實驗中各元素質(zhì)量分?jǐn)?shù)和氧化物相的固液相線溫度如表1所示，表中同時給出了VVER氧化物層各元素質(zhì)量分?jǐn)?shù)。

表1 熔融物成分和固液相線溫度Table 1Melt composition and solidus-liquidus temperature

通過對比，本文選擇與VVER 氧化物層主要成分接近的VULCANO 實驗參數(shù)設(shè)置氧化物層的相變潛熱和固液相線溫度。湍流場計算采用低雷諾數(shù)k-ε模型，固液相線間氧化物黏度μ的變化可根據(jù)懸浮液黏度公式（2）計算：

式中：μl為熔融物液相黏度，Pa·s；θs為固相體積分?jǐn)?shù)，可直接在計算過程中調(diào)用；系數(shù)C取值范圍為3.5～8。

計算域內(nèi)固相體積分?jǐn)?shù)可以直接從相變模型中實時調(diào)取。Journeau等［13］為了簡化計算提出了將黏度與溫度擬合的估算方式，同時推薦系數(shù)C取6，據(jù)此可以將式（2）推導(dǎo)為式（3）：

式中：T、Tliquidus和Tsolidus分別為熔融物溫度、液相線溫度和固相線溫度，K。

文獻(xiàn)［8］中給出了其他計算所需的參數(shù)，其中采用氧化物層膨脹系數(shù)為0.000 1K-1。在本文中，通過Boussinesq 假設(shè)計算熔融池中的自然流動過程，氧化層膨脹系數(shù)參考一般陶瓷類材料的膨脹系數(shù)設(shè)置，由此計算得到的氧化物層自然循環(huán)速度相對文獻(xiàn)［8］更加保守。綜合以上內(nèi)容，模擬過程中使用的物性參數(shù)具體如表2所示。

表2 物性參數(shù)Table 2 Physical properties

2.3 模型參數(shù)敏感性分析

COMSOL Multiphysics 默認(rèn)采用自適應(yīng)時間步長來求解瞬態(tài)問題，可以通過調(diào)整容差因子（Tolerance Factor）可以改變軟件對于步長的限制，從而提升計算效率。模型敏感性分析過程中，采用事故后2～5 h的衰變熱作為氧化物層加熱功率。

圖3所示為不同網(wǎng)格數(shù)和容差因子下的氧化物層中軸溫度分布，網(wǎng)格數(shù)18 276 和14 472的計算溫度偏差已經(jīng)縮小至1%以內(nèi)，因此，后續(xù)長期冷卻模擬采用的網(wǎng)格數(shù)為14 472。通過在同一臺計算機(jī)上對比研究發(fā)現(xiàn)，網(wǎng)格數(shù)18 276、14 472和11 790的模型中，將容差因子由0.5提高至2可分別將數(shù)值模擬時間縮短10.86%、10.05%和25.15%，可見，在計算模型達(dá)到網(wǎng)格無關(guān)性的程度后，容差因子對數(shù)值模擬時間的影響接近一致。然而研究過程中也發(fā)現(xiàn)過高的容差因子不會縮短數(shù)值模擬時間，對于網(wǎng)格數(shù)為14 472的模型，容差因子增加至4 反而增加了數(shù)值模擬時間，對于網(wǎng)格數(shù)為18 276的模型，容差因子設(shè)置為4 會導(dǎo)致計算過程不收斂。綜上所述，長期冷卻計算過程中采用的容差因子為2。

圖3 不同網(wǎng)格數(shù)和容差因子數(shù)下的模型中軸溫度分布Fig.3 Temperature distribution on the axis of calculation domain under different mesh number and tolerance factor

3 計算結(jié)果分析

3.1 速度場和溫度場

圖4和圖5分別為計算時長24 h 和96 h 時堆芯捕集器內(nèi)的速度場和溫度場，計算時氧化物層和金屬層的初始溫度分別設(shè)置為2 400K和2 200K。由于金屬層內(nèi)沒有衰變熱，其溫度迅速降低900K 以下，并在3 h左右就已經(jīng)完全凝固。氧化物層中上部在衰變熱和冷卻邊界作用下形成了不斷變化的自然對流，雖然自然對流的強(qiáng)度不高，但在其攪拌作用下該區(qū)域未發(fā)生明顯的熱分層。由于熔融物黏度隨著溫度降低以指數(shù)形式增加，在氧化物層下部區(qū)域的熔融物溫度較低，因此黏度很大流動緩慢，導(dǎo)致熱分層現(xiàn)象十分明顯。

圖4 熔融池內(nèi)速度分布 (a) 24 h，(b) 96 hFig.4 Distribution of velocity field in corium pool (a) 24 h, (b) 96 h

圖5 熔融池內(nèi)溫度分布 (a) 24 h，(b) 96 hFig.5 Distribution of temperature field in corium pool (a) 24 h, (b) 96 h

圖6 為不同計算時長下與氧化物層接觸的堆芯捕集器側(cè)壁面導(dǎo)出的熱流密度。側(cè)壁面受熔融池沖刷的位置熱流密度普遍在0.56MW·m-2以下，隨著側(cè)壁面處熔融物凝固結(jié)殼熱流密度降低至0.2MW·m-2以下，在計算時長2 h 時，局部最大熱流密度也僅達(dá)到約0.67MW·m-2。坩堝式堆芯捕集器的側(cè)壁面的導(dǎo)出熱流遠(yuǎn)低于大容積沸騰臨界熱流密度，可以認(rèn)為側(cè)壁面不存在燒毀風(fēng)險，堆芯捕集器的底部壁面與氧化物層之間由金屬層隔離，因此，導(dǎo)熱量很小也不存在燒毀風(fēng)險。

圖6 不同計算時長下的熔融池平均速度和平均溫度Fig.6Heat flux density of corium pool under different calculation durations

圖7為在計算時長12～96 h的范圍內(nèi)，堆芯捕集器內(nèi)熔融池的平均速度和平均溫度的變化。隨著計算時間增加，熔融池平均速度呈現(xiàn)波動降低的趨勢，其中的波動過程可能是氧化物層內(nèi)的自然對流不斷變化導(dǎo)致的。在當(dāng)前計算時間范圍內(nèi)，衰變熱由約13.3MW下降至約7.8MW，而捕集器內(nèi)平均溫度下降速率由12 h 時約5K·h-1迅速降低至24 h 時的約2K·h-1，到96 h 時已經(jīng)降低至1.1K·h-1左右，可見當(dāng)前事故緩解策略的降溫過程十分緩慢。在事故末期氧化物層完全凝固后，系統(tǒng)的熱量導(dǎo)出能力將更低，因此采用坩堝式堆芯捕集器作為嚴(yán)重事故緩解策略時，需要相關(guān)設(shè)備和支持輔助系統(tǒng)在很長時間內(nèi)保持可運(yùn)行性。

圖7 熔融池平均速度和平均溫度Fig.7 Average velocity and average temperature of corium pool

3.2 熔融池結(jié)殼形式

圖8 為計算時長4 h、12 h、24 h和96 h時熔融池內(nèi)的固相分布。在氧化物層冷卻過程中，上部、側(cè)壁和下部的結(jié)殼厚度均逐漸增加，同時熔融池內(nèi)的最小固相分?jǐn)?shù)也會逐漸增加。

圖8 熔融池內(nèi)固相體積分?jǐn)?shù) (a) 4 h，(b) 12 h，(c) 24 h，(d) 96 hFig.8 Soild volume fraction of corium in corium pool (a) 4 h, (b) 12 h, (c) 24 h, (d) 96 h

值得注意的是，在氧化物層上部被冷卻的大部分熔融物在完全凝固之前會受重力影響并隨著自然對流向底部轉(zhuǎn)移，這導(dǎo)致上部的結(jié)殼厚度的增加不大，而高固相份額的熔融物聚集在氧化物層下部導(dǎo)致該區(qū)域的結(jié)殼厚度顯著增加。側(cè)壁上凝固的熔融物也體現(xiàn)出向下遷移趨勢，導(dǎo)致側(cè)壁上形成了上薄下厚的結(jié)殼形式。氧化物層內(nèi)硬殼呈現(xiàn)上述形態(tài)是因為當(dāng)前熔融物固液相線間的相變溫差很大，遇冷后不會迅速完全凝固，因此存在受重力和流動影響而遷移的現(xiàn)象。在當(dāng)前結(jié)殼方式影響下，氧化物層向側(cè)壁的導(dǎo)出熱量占與總釋熱率的比值，由12 h 時的約60%逐漸降低至96 h時的約41%。考慮到頂部與水直接接觸部分的換熱效果更好，這種結(jié)殼方式對提升熔融池的冷卻速率是有利的。

4 結(jié)語

采用多物理場耦合計算軟件COMSOL 搭建了VVER 堆芯捕集器計算模型，對捕集器中熔融池的冷卻過程進(jìn)行了數(shù)值模擬研究，主要結(jié)論如下：

1）在熔融池分層結(jié)構(gòu)下，金屬層會迅速凝固，氧化物層中上部會形成不斷變化的自然對流流動且溫度分布較為均勻，氧化物層下部流動緩慢并出現(xiàn)了明顯的熱分層。

2）在氧化物層上部被冷卻的大部分熔融物在完全凝固之前會受重力和自然對流影響向下部轉(zhuǎn)移，導(dǎo)致上部結(jié)殼厚度增加緩慢，而底部結(jié)殼厚度增加迅速。

3）坩堝式堆芯捕集器的安全裕量較高，冷卻壁面不存在燒毀風(fēng)險，但其中的熔融物冷卻過程十分緩慢，需要相關(guān)設(shè)備和支持輔助系統(tǒng)在很長時間內(nèi)保持可運(yùn)行性。

作者貢獻(xiàn)聲明朱光昱負(fù)責(zé)數(shù)值模擬和文章撰寫；郭超、劉巧鳳負(fù)責(zé)文章撰寫；李春、依巖負(fù)責(zé)工作支持。