楊 璞,嚴 睿,鄒 楊,于世和,朱貴鳳,周 波, 馬玉雯,*
(1.中國科學院 上海應用物理研究所,上海 201800;2.中國科學院大學,北京 100049)
TMSR-SF1是熱功率為10 MW的熔鹽冷卻球床實驗堆,燃料元件采用球形燃料,冷卻劑采用2LiF-BeF2(7Li豐度>99.99%),相比于傳統水堆,TMSR-SF1在堆芯結構、冷卻劑、燃料組件等方面有很大的不同,而準確計算堆內釋熱率分布對于這種新型反應堆的熱工水力設計、瞬態分析、結構力學設計等都有重要意義。本文使用蒙特卡羅計算程序MCNP(版本5-1.51)對TMSR-SF1壽期初(BOL)及壽期末(EOL)堆內部件(包括燃料球、冷卻劑、控制棒、石墨反射層及合金結構材料等)的釋熱率進行詳細計算。MCNP程序可通過F6及F7能量計數卡來統計裂變碎片、中子、瞬發γ和俘獲γ的能量沉積,但該程序無法直接計算緩發β及緩發γ的能量沉積,因此一般通過統一乘1個系數來計入這兩者對堆內部件釋熱率的貢獻,但這種方法忽略了緩發β與緩發γ在堆內能量沉積分布的不同。本文通過使用光子產生偏倚卡(pikmt),經過3次MCNP輸運計算,分別統計裂變碎片、中子、瞬發γ、俘獲γ、緩發β及緩發γ在TMSR-SF1堆內不同部件內的能量沉積,并最終計算出部件內的釋熱率。
TMSR-SF1堆內釋熱率分布計算需要分別進行3次MCNP輸運計算,即:1) 柵元體積計算;2) 使用F6:n、F6:p和F7:n統計卡的常規MCNP計算;3) 使用光子產生偏倚卡及F6:p統計卡的MCNP計算。計算結果處理需要兩步,首先需通過常規MCNP計算F6:n、F6:p和F7:n卡給出的統計結果來求得平均單次裂變放出的可利用能量(Q),然后根據Q值來計算堆內不同部件的總釋熱率及體積釋熱率。
根據反應堆總熱功率和平均單次裂變放出的可利用能量(即Q)可計算反應堆單位時間的裂變次數,因此對于反應堆釋熱率計算,Q值的計算是非常重要的一步。反應堆Q值由以下幾部分組成[9]:裂變碎片的能量貢獻(Qfp)、中子的能量貢獻(Qn)、瞬發γ的能量貢獻(Qγp)、俘獲γ的能量貢獻(Qγc)、裂變產物衰變放出的緩發β粒子能量貢獻(Qβ)和緩發γ的能量貢獻(Qγd)。MCNP程序中的F6和F7卡可用來統計裂變碎片、中子、緩發和俘獲γ的能量。其中F6:n用于統計裂變碎片和中子的能量,F6:p用于統計瞬發和俘獲γ的能量,F7用于統計裂變碎片、中子、瞬發γ的能量且與F6卡統計原理不同,F7卡假設所有能量就地沉積。因此,可用下式計算Q值[10]:
(1)
(2)
(3)

由于MCNP本身無法直接計算緩發γ的能量貢獻Qγd和緩發β的能量貢獻Qβ,因此這兩部分直接使用ENDF/B-Ⅶ.0數據庫中235U的緩發γ和緩發β的能量數據,分別取6.33 MeV和6.50 MeV[11]。一般反應堆內的緩發γ主要來源于235U裂變產物的衰變,因此使用235U的緩發γ的能量作為整個反應堆的Qγd值是一個合理的假設[10,12-13]。β粒子穿透力很低,可認為β粒子的全部能量就地沉積,因此堆內緩發β的能量沉積近似和反應堆的裂變率呈正比,對于使用U燃料的反應堆,這部分能量可直接使用235U的緩發β的能量數據[10]。
反應堆堆內部件的釋熱率主要來源于裂變碎片、中子、緩發β、瞬發γ、俘獲γ和緩發γ的能量沉積。MCNP通過柵元來對反應堆內的實際部件進行建模,某個柵元i內裂變碎片、中子、瞬發和俘獲γ的能量沉積可使用下式計算:
(4)

(5)

采用SPSS15.0統計學軟件,計量資料采用()的形式表示,用t檢驗;計數資料采用x2檢驗,P<0.05為差異具有統計學意義。
(6)

柵元i的釋熱率為上述所有粒子在柵元i內能量沉積的總和(式(7))。柵元i的釋熱率除以其柵元體積便得到柵元i的體積釋熱率(式(8))。反應堆的總釋熱率(即總功率)為所有柵元釋熱率的總和(式(9))。
(7)
(8)
(9)

圖1 釋熱率計算流程Fig.1 Computational step for heat rate
為簡化上述計算過程中數據傳輸和消除容易出錯的手動數據傳輸過程,本文使用Python語言編寫用于輸入輸出處理的程序。圖1為整個計算流程。計算過程包括3次MCNP運行:第1次運行采用固定源模式并配合材料無效卡(void),用于計算柵元的體積;第2次運行則采用臨界模式(使用kcode卡)配合F6:n、F6:p和F7:n計數卡,用來計算裂變碎片、中子、緩發β、瞬發γ和俘獲γ的能量沉積;第3次運行在臨界模式使用pikmt卡且配合F6:p計數卡,用以計算緩發γ的能量沉積。最后,采用自編程序自動讀取輸出文件中的統計結果,通過數據處理得到Q值及各柵元的釋熱率及平均體積釋熱率,并將這些計算結果整理和輸出到excel文件中。
TMSR-SF1堆芯主要由堆芯活性區、石墨反射層和哈氏合金圍筒組成。堆芯活性區為燃料球和石墨球的隨機堆積區,燃料球及石墨球間空隙形成熔鹽隨機不規則流道,供其自下往上流動帶走熱量,活性區由中心圓柱體和上下圓臺組成,其中圓柱體直徑為135 cm、高度為180 cm,上下圓臺直徑為30 cm、高度為30 cm。活性區外為石墨反射層,其內布置16個控制棒通道、2個燃料球注入通道、1個熔鹽裝卸通道、1個中子源通道、2個物理啟動用測量通道、7個中子通量與能譜測量通道、6個溫度測量通道、3個隨堆監督通道及2個備用通道。每個通道內還包含套管,其中溫度測量通道、中子通量與能譜測量通道、燃料球注入通道及中子源通道內為哈氏合金套管,其他通道內則為碳-碳復合材料(C-C)套管。反射層外為堆芯圍筒,其厚度為3 cm。整個堆芯直徑為292.6 cm、高度為308.4 cm。堆芯具體結構如圖2所示。

圖2 TMSR-SF1堆芯結構Fig.2 Schematic of TMSR-SF1 core
TMSR-SF1活性區裝載的燃料球與石墨球直徑為6 cm,其中燃料球具體結構示于圖3,主要參數列于表1。

圖3 燃料球及TRISO包覆顆粒結構Fig.3 Fuel pebble and TRISO coated particle structure

參數數值直徑,cm6燃料區域直徑,cm5石墨基體和石墨球殼密度,g/cm31.73235U富集度,%17.0單個燃料球U裝量(純鈾),g7.00燃料顆粒直徑,mm0.5UO2核芯密度,g/cm310.4包覆層材料(從內到外)PrC/PyC/SiC/PyC包覆層厚度,μm95/40/35/40包覆層密度,g/cm31.1/1.90/3.18/1.90
考慮到TRISO包覆顆粒在燃料球內的填充比一定的情況下,其隨機分布與規則分布計算的燃料球無限介質增殖因數并無太大差別[15],因此為簡化建模過程,本文采用簡單立方建模方式來模擬TRISO顆粒在燃料球內的排布,為保證單顆燃料球鈾裝量為7.00 g,簡單立方的邊長被設為0.178 cm。TRISO及燃料球模型結構如圖4所示。
TMSR-SF1堆內燃料球采用隨機堆積方式,燃料球在堆芯內的占空比為60%,但在MCNP建模過程中由于隨機堆積模型存在較大的困難和不確定性[7],因此,本文采用圖5所示的體心立方排布方式(占空比60%)對堆芯內的燃料球進行等效建模,體心立方邊長為7.224 cm。

圖4 TRISO顆粒及燃料球模型Fig.4 Geometry model of TRISO and fuel pebble

圖5 隨機球床堆的燃料球堆積等效Fig.5 Equivalent of random packing of fuel pebble
圖6為最終的MCNP模型。模型主要包括堆芯活性區、石墨反射層、堆內通道、通道套管、控制棒及堆芯哈氏合金圍筒等結構。

圖6 TMSR-SF1的MCNP模型Fig.6 MCNP model of TMSR-SF1
表2為TMSR-SF1在壽期初和壽期末時的Q值計算結果。從表2可見,壽期末的Q值相較于壽期初略有增加,這主要是由于壽期末有一定量的239Pu產生,而239Pu的單次裂變能大于235U的。

表2 TMSR-SF1 Q值計算結果Table 2 Fission Q value of TMSR-SF1 model
表3為TMSR-SF1在壽期初和壽期末時堆內不同部件的釋熱率計算結果,其中合金構件包括堆芯圍桶及通道合金套管。由表3可見,TMSR-SF1的堆內釋熱主要集中在燃料球內,其釋熱率占堆內總釋熱率的94%以上,熔鹽和反射層由于在堆內體積較大,因此釋熱率也較高,分別占到總釋熱率的3%和1%以上,其他堆內部件占總釋熱率的比例都小于1%。與壽期初相比,壽期末燃料球釋熱率略有減小,這主要歸因于燃料球的釋熱率主要來源于裂變碎片和裂變中子的能量沉積,且壽期末單位時間的裂變次數小于壽期初的。此外,控制棒與石墨球的釋熱率也有所減少,前者主要是壽期末控制棒移動到堆頂,導致能入射到控制棒的中子及γ減少造成的。而石墨球的釋熱率減少,主要原因是堆芯活性區中子和γ軸向通量密度分布的峰值上移,而石墨球位于堆芯活性底部,所以入射到石墨球的中子及γ也有所減少。而由于壽期末控制棒上移,堆內產生的裂變中子被控制棒吸收的概率減小,中子及中子被吸收后產生的俘獲γ在堆芯活性區外的部件內沉積的能量增加,導致反射層、冷卻劑熔鹽、孔道C-C套管及堆內合金構件的釋熱率增加。

表3 TMSR-SF1堆內釋熱率分布計算結果Table 3 Heat rate distribution results of TMSR-SF1
圖7為壽期初和壽期末堆內主要結構件(非燃料球)的釋熱來源。由圖7可見,結構件的釋熱率主要來源于γ的能量沉積,其中包含合金材料的控制棒及合金構件的釋熱主要是來源于俘獲γ的能量沉積,這是因為合金材料會俘獲中子發生(n,γ)反應并放出大量的俘獲γ,而俘獲γ的大部分能量會就地在合金材料內沉積。
表4為不同部件內的平均和最大體積釋熱率。圖8為壽期初和壽期末全堆體積釋熱率的精細分布。可見,燃料球的體積釋熱率遠大于其他堆內部件,壽期初的最大體積釋熱率為22.14 W/cm3,壽期末的最大體積釋熱率為20.26 W/cm3。反射層靠近活性區的部分體積釋熱率較大,最大值達到約0.7 W/cm3,而遠離活性區的部分體積釋熱率僅約10-5W/cm3量級。此外,相比于壽期初,壽期末反射層周向的體積釋熱率分布更加均勻。控制棒的最大體積釋熱率位于控制棒底端的配重合金棒頭,壽期初和壽期末的最大體積釋熱率分別可達1.408 W/cm3和0.9502 W/cm3。堆內合金構件的最大體積釋熱率位于中子源通道的合金套管,壽期初為0.792 1 W/cm3,壽期末增加到0.822 4 W/cm3。

圖7 TMSR-SF1堆內結構件釋熱來源Fig.7 Heat source of component in TMSR-SF1

堆內構件平均體積釋熱率/(W·cm-3)最大體積釋熱率/(W·cm-3)BOLEOLBOLEOL燃料球6.3636.34922.14020.260石墨球0.0470.0380.2230.196反射層0.0090.0110.7130.751冷卻劑熔鹽0.1290.1281.5341.464孔道C-C套管0.0260.0280.1550.125控制棒0.0740.0541.4080.950合金構件0.0300.0400.7920.822
根據經驗,相對誤差小于5%的結果通常認為是可靠的[16]。表5為上述計算結果的相對誤差。體積計算投入粒子數為5×1010,體積計算結果的相對誤差不超過0.02%。能量沉積計算中單代投入粒子數為100 000,有效迭代次數為4 550,所有F6、F7計數卡結果的相對誤差都小于1%。

圖8 TMSR-SF1體積釋熱率分布Fig.8 Volumetric heat rate distribution of TMSR-SF1

堆內構件相對誤差/%體積計算BOLF6:n、F6:p、F7:n卡統計使用pikmt卡的F6:p統計體積計算EOLF6:n、F6:p、F7:n卡統計使用pikmt卡的F6:p統計燃料球0~0.010.01~0.020.01~0.020~0.010.01~0.020.01~0.02石墨球0~0.020.01~0.020.02~0.270~0.020.01~0.020.01~0.14反射層0~0.010.02~0.820.04~0.890~0.010.02~0.780.04~0.81冷卻劑熔鹽0~0.010.03~0.540.01~0.650~0.010.03~0.480.01~0.62孔道C-C套管0~0.020.02~0.460.02~0.570~0.020.02~0.450.02~0.55控制棒0~0.020.04~0.800.05~0.910~0.020.04~0.850.05~0.99合金構件0~0.020.05~0.760.04~0.81 0~0.020.05~0.720.03~0.78
通過使用光子產生偏倚卡,經過3次MCNP輸運計算,對TMSR-SF1壽期初和壽期末堆內的釋熱率進行了詳細計算研究,得到了堆內不同部件的總釋熱率、體積釋熱率分布和最大體積釋熱率,并得到以下結論。
1) 燃料球釋熱率占堆內總釋熱率的94%以上,熔鹽和反射層由于在堆內體積較大,因此釋熱率也較高,分別占總釋熱率的3%和1%以上,其他堆內部件釋熱率占總釋熱率的比例都小于1%。
2) 與壽期初相比,壽期末燃料球、控制棒與石墨球的釋熱率有所減少,而其他構件包括反射層、冷卻劑熔鹽、孔道C-C套管及合金構件的釋熱率都有所增加。
3) 包含合金材料的控制棒及堆內合金構件的釋熱率主要來源是俘獲γ的能量沉積。控制棒的最大體積釋熱率位于控制棒底端的配重合金棒頭,合金構件的最大體積釋熱率位于中子源通道的合金套管。