謝 凱，李林森，沈 峰，費(fèi)立凱，王納秀

（1.國(guó)核（北京）科學(xué)技術(shù)研究院，北京 100029；2.中國(guó)科學(xué)院 上海應(yīng)用物理研究所，上海 201800）

釷基熔鹽堆（TMSR）核能系統(tǒng)是第4代核能系統(tǒng)，是一種使用氟化鹽作為冷卻劑，工作在常壓下的高溫反應(yīng)堆。相較于一般壓水反應(yīng)堆，釷基熔鹽堆具有固有安全性，并可使釷轉(zhuǎn)化為鈾實(shí)現(xiàn)增殖等優(yōu)點(diǎn)［1］。

TMSR-SF1是一種使用石墨包覆顆粒球形燃料的釷基熔鹽堆，其堆芯結(jié)構(gòu)復(fù)雜。傳統(tǒng)的熱工水力系統(tǒng)程序可模擬整個(gè)堆芯回路不同狀態(tài)下的工況，但無(wú)法對(duì)堆芯結(jié)構(gòu)進(jìn)行細(xì)致的模擬。ANSYS CFX 是一個(gè)采用有限容積法的大型商用CFD 軟件，有豐富的物理模型，可進(jìn)行流動(dòng)、耦合傳熱、多相流等問(wèn)題的模擬分析［2］。本文使用ANSYS DM、ANSYS Meshing建立幾何模型，在CFX-Pre中添加相關(guān)材料的物性，應(yīng)用CFX 的多孔介質(zhì)模型，初步模擬堆芯在穩(wěn)態(tài)運(yùn)行和喪失強(qiáng)迫循環(huán)后，堆芯的溫度分布情況，為堆芯的設(shè)計(jì)優(yōu)化提供一定的參考依據(jù)。

1 TMSR-SF1堆芯結(jié)構(gòu)

本文基于10 MW TMSR-SF1堆芯幾何結(jié)構(gòu)，堆芯整體結(jié)構(gòu)示意圖如圖1所示。

圖1 TMSR-SF1堆芯整體結(jié)構(gòu)Fig.1 TMSR-SF1core configuration

從石墨反射層向外側(cè)依次為石墨反射層的哈氏合金圍筒、下降環(huán)腔熔鹽、下降環(huán)腔外哈氏合金圍筒、氬氣層、不銹鋼、保溫層、保溫層外壁和大氣環(huán)境。

堆芯區(qū)域?yàn)橐粋€(gè)準(zhǔn)正八邊形結(jié)構(gòu)，內(nèi)部除了最上層球與最下層球?yàn)槭蛲?，均為燃料球，球與球?qū)訉佣逊e起來(lái)并接觸相切。上下為開(kāi)有熔鹽流道的石墨反射層，四周為開(kāi)有控制棒孔道和測(cè)量孔道的石墨反射層。

2 幾何建模和網(wǎng)格劃分

本文幾何建模使用ANSYS DM 軟件。為了簡(jiǎn)化建模，建立了1／4模型代替整個(gè)模型并設(shè)置旋轉(zhuǎn)對(duì)稱面進(jìn)行模擬。為了避免設(shè)置對(duì)稱后進(jìn)出口位置引起的問(wèn)題，本模型忽略了熔鹽上下腔室的區(qū)域，入口設(shè)為下反射層流道入口，出口為上反射層流道出口，建立了反射層高度的堆芯結(jié)構(gòu)模型，如圖2所示。

圖2 DM 建模俯視圖（a）和正視圖（b）Fig.2 Vertical view（a）and front view（b）of model

使用ANSYS Meshing軟件進(jìn)行幾何結(jié)構(gòu)的網(wǎng)格劃分。網(wǎng)格完成后，數(shù)量為198萬(wàn)，最低網(wǎng)格質(zhì)量為0.001 7，質(zhì)量小于0.1的網(wǎng)格均不位于交界面處，對(duì)模擬影響不大，符合建模模擬要求。

3 多孔介質(zhì)模型

ANSYS CFX 的多孔介質(zhì)模型中，認(rèn)為在每一小控制體內(nèi)，流體的體積［3］為：

式中：ε為多孔介質(zhì)孔隙率；V 為控制體的實(shí)際體積。

假定流體通過(guò)控制體的一個(gè)外邊界面An時(shí)，可流通的面積為：

式中：An＝nds，n 為面An的外側(cè)單位法向量，ds為面An的面積；K＝Kij為一個(gè)二階對(duì)稱張量，稱為面積孔隙率張量。在ANSYS 14.0版本中的ANSYS CFX，Kij只允許為各向同性，Kij＝εδij。

流體在ANSYS CFX 的多孔介質(zhì)模型中，遵守以下3個(gè)基本微分方程［3］：

質(zhì)量守恒方程：

動(dòng)量守恒方程：

能量守恒方程：

式中：μe 為有效黏度；SM為動(dòng)量源，在多孔介質(zhì)中可認(rèn)為是一個(gè)阻力梯度；H 為流體的焓；Γe為有效熱導(dǎo)率；Sh為體積內(nèi)熱源，本次模擬中將直接給定一個(gè)數(shù)值。

多孔介質(zhì)模型中，流體主要在流動(dòng)方向上產(chǎn)生動(dòng)量損失，式（4）中動(dòng)量源SM在一個(gè)方向上的表達(dá)式［3］為：

式中：k為滲透率；U 為流體在該方向上的速度；C 為阻力損失系數(shù)。滲透率k表示了多孔介質(zhì)允許流體穿過(guò)孔隙流動(dòng)性質(zhì)的量。本次模擬采用經(jīng)驗(yàn)關(guān)系式計(jì)算滲透率。對(duì)于由大小相同的球組成的多孔介質(zhì)，滲透率的計(jì)算公式［4］為：

式中：d 為水力直徑；KK為Kozeny常數(shù)，一般情況下取值為5。

對(duì)于由大小相同的球組成的多孔介質(zhì)，其水力直徑可采用如下公式［4］：

式中，dp為球的直徑。由式（7）、（8）可得滲透率k。

由式（6）、（7）、（9），得阻力損失系數(shù)C。

穩(wěn)態(tài)模擬時(shí)堆芯熔鹽區(qū)域的溫度分布可與系統(tǒng)程序RELAP5-3D 模擬結(jié)果進(jìn)行比較。使用RELAP5-3D 將堆芯區(qū)域劃分為多個(gè)水力部件并添加相應(yīng)熱構(gòu)件，僅建立堆芯球床區(qū)域的等效模型，模型示意圖如圖4a所示。在等效模型中添加堆芯區(qū)域相應(yīng)的進(jìn)出口條件，保持與CFX 相同的邊界條件，進(jìn)行穩(wěn)態(tài)模擬。模擬結(jié)果的不同高度溫度比較如圖4b所示。從圖中可看出，兩種建模的堆芯區(qū)域不同高度的平均溫度基本相同，證明了CFX 的模擬結(jié)果基本可靠。

4 邊界條件設(shè)置

本次建模是要模擬在堆芯區(qū)域的熔鹽喪失強(qiáng)迫循環(huán)后，堆芯能否依靠本身的散熱能力將衰變熱排出，也就是模擬事故后堆芯內(nèi)部的溫度分布情況。域的類型有多孔介質(zhì)域、流體域、固體域3 種。其中氬氣層是流體域，材料為氬氣，基本設(shè)置和流體域設(shè)置相同，并需考慮輻射換熱，增加了分離傳輸?shù)幕殷w輻射換熱模型［6］。

首先進(jìn)行的是熔鹽通過(guò)堆芯的穩(wěn)態(tài)模擬，堆芯處于熱平衡狀態(tài)，堆芯功率設(shè)為均勻的平均熱功率，并添加進(jìn)出口參數(shù)和其他一些參數(shù)，參數(shù)設(shè)置列于表1。

表1 邊界條件基本參數(shù)Table 1 Boundary condition parameters

5 初步計(jì)算結(jié)果

5.1 穩(wěn)態(tài)模擬初步結(jié)果

本次模擬在40步以內(nèi)收斂。CFX 穩(wěn)態(tài)模擬分析的溫度分布如圖3所示，各主要區(qū)域的溫度列于表2。

由圖3和表2可看出，整個(gè)堆芯區(qū)域溫度基本在873K 左右，在厚度為30cm 的保溫層有一個(gè)較大的溫度梯度，內(nèi)外溫差300K 左右。

圖3 穩(wěn)態(tài)模擬的堆芯縱截面溫度分布（a）和橫截面溫度分布（b）Fig.3 Temperature distributions of longitudinal section（a）and cross section（b）of core at steady state

表2 堆芯區(qū)域平均溫度Table 2 Average temperature in core

5.2 瞬態(tài)模擬初步結(jié)果

使用穩(wěn)態(tài)模擬的結(jié)果作為初始值進(jìn)行瞬態(tài)模擬，計(jì)算堆芯失去強(qiáng)迫循環(huán)后滿功率運(yùn)行10s的溫度情況，其溫度分布形式與穩(wěn)態(tài)結(jié)果無(wú)明顯變化，但熔鹽溫度有明顯上升。熔鹽的最高溫度達(dá)981K，平均溫度達(dá)935K，堆芯最高溫度上升了18K。堆芯熔鹽最高溫度按時(shí)間變化曲線如圖5所示，堆芯上部熔鹽流動(dòng)速度分布如圖6所示。

圖4 RELAP5-3D建模示意圖（a）和堆芯熔鹽不同高度平均溫度與RELAP5-3D模擬結(jié)果（b）Fig.4 RELAP5-3D modeling schematic diagram（a）and axial average temperature distribution comparison between CFX and RELAP5-3D（b）

圖5 瞬態(tài)前10s堆芯熔鹽最高溫度變化曲線Fig.5 Transient evolution of maximum temperature of molten salt at first 10s

從圖6可看出，由于密度差，堆芯熔鹽受熱上升，堆芯上部與上流道交界處熔鹽有相對(duì)較大的流速，在上流道內(nèi)產(chǎn)生了自然對(duì)流換熱［7］。

然后模擬堆芯功率按指數(shù)衰減、10～60s時(shí)的堆芯溫度分布情況。堆芯功率按式（10）衰減：

計(jì)算完成后得堆芯60s時(shí)的溫度分布，如圖7所示。

圖6 堆芯熔鹽出口速度分布云圖（a）和堆芯熔鹽出口速度分布矢量圖（b）Fig.6 Core velocity distribution contour（a）and velocity distribution vector（b）

圖7 瞬態(tài)模擬堆芯60s時(shí)的縱截面溫度分布（a）和橫截面溫度分布（b）Fig.7 Transient temperature distributions of longitudinal section（a）and cross section（b）of core at 60s

從圖7可看出，堆芯熔鹽的溫度分布較為均勻，圖7a中的最高溫度為957 K。整個(gè)堆芯熔鹽的最高溫度為968K，較10s時(shí)略有降低。10～60s中，堆芯最高溫度降低了13K，堆芯熔鹽最高溫度按時(shí)間變化曲線如圖8所示。

5.3 增加絕熱層壁面厚度的模擬結(jié)果

從以上結(jié)果可看出，堆芯溫度梯度大的地方發(fā)生在保溫層中。在以上的建模模擬中保溫層的厚度均為30cm，不一定是最佳的保溫層厚度。本次模擬將保溫層厚度由30m 增加到40cm，邊界條件設(shè)置保持不變，進(jìn)一步研究絕熱層壁面厚度對(duì)堆芯溫度分布的影響。絕熱層厚度為40cm 時(shí)，穩(wěn)態(tài)模擬中，各主要區(qū)域的溫度列于表3。

圖8 瞬態(tài)模擬10～60s堆芯熔鹽最高溫度變化曲線Fig.8 Transient evolution of maximum temperature of molten salt within 10-60s

對(duì)比表2和表3可看出，在堆芯穩(wěn)態(tài)運(yùn)行時(shí)，保溫層厚度由30cm 增加到40cm 后，堆芯溫度基本不變，保溫層的保溫效果并沒(méi)有明顯提升，外壁面溫度由508K 降到488K，降低了20K，對(duì)于降低由壁面產(chǎn)生的熱量損失的效果較小。

表3 保溫層厚度為40cm 時(shí)堆芯區(qū)域平均溫度Table 3 Average temperature in core for insulation layer thickness 40cm

進(jìn)行瞬態(tài)模擬中，邊界條件也保持不變，同樣模擬了在失去強(qiáng)迫循環(huán)后，堆芯功率前10s滿功率，然后按指數(shù)衰減的情況。堆芯熔鹽最高溫度按時(shí)間變化曲線如圖9所示。

圖9 保溫層厚度為40cm 時(shí)瞬態(tài)堆芯熔鹽最高溫度變化曲線Fig.9 Transient evolution of maximum temperature of molten salt with insulation layer thickness 40cm

從圖9可看出，堆芯最高溫度在10s內(nèi)從958K 增加到979K，增加了21K，然后在隨后的50s內(nèi)溫度下降到968K，較溫度峰值下降了11 K。本次模擬中峰值溫度較絕熱層為30cm時(shí)還低，是由于本次模擬中穩(wěn)態(tài)時(shí)堆芯最高溫度較低，其誤差可能來(lái)源于兩次網(wǎng)格的差異。

比較保溫層厚度由30cm 增加到40cm后，堆芯溫度變化情況可看出，雖然壁厚增加了1／3，但穩(wěn)態(tài)運(yùn)行時(shí)外壁面溫度并未下降很多。而在瞬態(tài)分析中堆芯主要靠壁面散熱，所以結(jié)果差異不大。

6 結(jié)論

本文基于10 MW TMSR-SF1 堆芯，建立了簡(jiǎn)化的堆芯多孔介質(zhì)模型，在主要考慮導(dǎo)熱換熱和浮力作用的情況下，對(duì)堆芯穩(wěn)態(tài)運(yùn)行時(shí)的溫度分布和發(fā)生事故后60s的瞬態(tài)溫度分布進(jìn)行了初步分析。

最終模擬了堆芯在喪失強(qiáng)迫循環(huán)，10s后停堆的60s內(nèi)，以外壁面為熱阱的散熱情況。在忽略了熔鹽輻射換熱影響的情況下，堆芯主要依靠導(dǎo)熱和自然對(duì)流進(jìn)行散熱，在保溫層厚度為30cm 和40cm 的情況下，堆芯除熔鹽區(qū)域外其他區(qū)域溫度沒(méi)有大的變化，堆芯熔鹽能達(dá)到的最高溫度為981K，低于熔鹽的沸點(diǎn)。

下一步研究將基于目前建立的多孔介質(zhì)模型，主要考慮對(duì)熔鹽加入輻射換熱模型，并進(jìn)一步改變保溫層厚度，對(duì)堆芯的溫度分布進(jìn)行更加準(zhǔn)確的模擬分析，并進(jìn)一步比較不同保溫層厚度對(duì)堆芯溫度的影響。

［1］ 江綿恒，徐洪杰，戴志敏.未來(lái)先進(jìn)核裂變能——TMSR 核能系統(tǒng)［J］.中國(guó)科學(xué)院院刊，2012（3）：366-372.JIANG Mianheng，XU Hongjie，DAI Zhimin.Advanced fission energy program：TMSR nuclear energy system［J］.Bulletin of Chinese Academy of Sciences，2012（3）：366-372（in Chinese）.

［2］ 謝龍漢，趙新宇，張烔明.ANSYS CFX 流體分析及仿真［M］.北京：電子工業(yè)出版社，2012.

［3］ ANSYS Inc.CFX-solver theory guide［M］.PA：ANSYS Inc.，2011.

［4］ LI Linsen，YUAN Haomin，WANG Kan.Coupling of RMC and CFX for analysis of pebble bed-advanced high temperature reactor core［J］.Nuclear Engineering and Design，2012，250：385-391.

［5］ KAVIANY M.Principles of heat transfer in porous media［M］.2nd ed.Berlin：Springer，1995.

［6］ ANSYS Inc.ANSYS CFX-pre user's guide［M］.PA：ANSYS Inc.，2011.

［7］ 楊世銘，陶文銓.傳熱學(xué)［M］.北京：高等教育出版社，2006.