石墨慢化通道式熔鹽堆有效緩發中子份額計算方法研究

崔 勇，陳金根,*，何 帆，李曉曉，蔡翔舟,*

(1.中國科學院 上海應用物理研究所，上海 201800； 2.中國科學院 先進核能創新研究院，上海 201800；3.中國科學院大學，北京 100049)

熔鹽堆是6種第4代核能系統候選堆型中唯一使用液態燃料的反應堆，與固態燃料堆不同，熔融的氟化鹽混合物(LiF，BeF2，ZrF4，UF4或ThF4)兼作冷卻劑與燃料，充滿整個一回路(含堆芯)。燃料流動使得緩發中子先驅核(DNP)在堆內的分布受到流場的影響。相應地，緩發中子的價值也會發生改變。部分具有較長半衰期的DNP會隨燃料流出堆外并衰變[1]，從而導致有效緩發中子份額(βeff)降低。其次，裂變能直接釋放于冷卻劑中，因此傳統的流動-傳熱計算方法不再適用。最后，熔鹽堆適合采用釷鈾燃料循環[2]，因而易裂變材料中可能包含233U和錒系核素(MA)[3]，這些重核的緩發中子數據與采用鈾钚燃料循環的固態堆不同。這些特殊性使得熔鹽堆βeff的計算較為復雜。

βeff的精確計算對研究熔鹽堆的瞬態安全特性具有重要意義。如在失流事故中，緩發中子源項的再分布會向堆內引入正反應性，該反應性的量化需對不同流動工況下的βeff進行評估[3]；在反應性引入事故中，需確定所引入反應性與βeff的關系以得到反應堆的響應特性[4-5]。近年來，國際國內針對DNP的對流[6]、擴散[7]和湍流輸運效應[8]對熔鹽堆穩態與瞬態特性的影響開展了廣泛研究。對于罐式熔鹽快堆，由于其堆內無慢化材料，可采用解析方法[3]、數值方法[3,8-9]和蒙特卡羅方法[3]計算其βeff與DNP分布。對于通道式熔鹽堆，由于存在石墨、氫化鋯等慢化材料，其流動-傳熱方式具有特殊性[5,10]，對快堆提出的βeff計算方法不再適用。

本文利用解析方法和基于四階多項式節塊展開法與并聯多通道模型的三維核熱耦合程序TMSR3D[11]，對熔鹽實驗堆(MSRE)進行穩態計算，得到堆內DNP分布情況，研究入口流量、燃料的堆外流動時間等參數對βeff的影響。

1 計算模型與方法

1.1 計算模型

MSRE堆芯與組件幾何結構[12]如圖1所示。堆芯由正方形石墨組件構成，組件四周開孔作為燃料鹽流道，所有孔道與上下腔室連通。堆芯四周為石墨反射層。表1、2分別列出MSRE的主要設計參數[4,12]和裝載235U燃料時的緩發中子數據[4]。詳盡的堆芯設計與熱工參數參見文獻[4]。

表1 MSRE設計參數Table 1 Design parameter of MSRE

a——堆芯縱截面；b——堆芯橫截面；c——組件幾何圖1 MSRE堆芯與組件幾何結構Fig.1 Core and assembly configuration of MSRE

1.2 計算方法

熔鹽堆穩態下中子學計算模型為：

(1)

(2)

1.2.1解析方法 利用圓柱裸堆幾何對MSRE的臨界通量進行近似，記MSRE的堆芯半徑和高度分別為R和H，可得到圓柱坐標下的通量[13]為：

(3)

其中：J0為零階貝塞爾函數；φ0為待定系數，由反應堆功率決定。

僅考慮軸向流動，記流速為u，將式(3)代入(2)中得到：

(4)

對上式在(0，z)區間積分得到：

c0exp(-λiz/u)

(5)

其中，c0為堆芯入口z=0處的DNP濃度。

設堆芯中熔鹽占堆芯總體積的比值恒為f，且堆芯出口處先驅核濃度均勻混合，其出口濃度ci，out可表示為：

(6)

設熔鹽在堆芯外回路中時間為τ，先驅核在堆外回路僅發生衰變，因此入口處先驅核濃度ci，in可表示為：

ci,in=ci,outexp(-λiτ)

(7)

由式(3)可得堆芯內中子產生率N為：

sin(πz/H)2πrfdrdz

(8)

從而得到歸一化緩發中子源項Si，d為：

Si,d(r,z)=λici(r,z)/N

(9)

由于熔鹽流速u不隨半徑變化，則得到熔鹽在堆芯中的流動時間tc為：

tc=H/u

(10)

將式(5)、(7)和(8)代入式(9)中得到緩發中子源項分布：

Si,d(r,z)=α0,iexp(-λitcz/H)+

(α1,isin(πz/H)-α2,icos(πz/H)+

α2,iexp(-λitcz/H))J0(2.405r/R)

(11)

α0,i=βi·

(12)

(13)

(14)

其中，J1為一階貝塞爾函數。第i組緩發中子有效份額βi，eff為：

(15)

1.2.2數值方法 熔鹽堆穩態G群中子擴散方程和I組DNP守恒方程[11]如下：

g=1,2,…,G;i=1,2,…,I

(16)

(17)

其中，Dg、Σt，g和Σg′g分別為g群擴散系數、總截面和散射矩陣?；诠潐K展開法和有限差分方法[11]可求解得到穩態下的中子通量分布和DNP分布。

熱工水力學計算采用并聯多通道模型[5,11,14]，包括燃料鹽的質量、動量和能量守恒方程，即：

(18)

(19)

(20)

其中：i和j分別為通道編號和軸向節塊編號；M和N分別為軸向總節塊數和堆芯總通道數；G、A和W分別為質量流速、等效通道流通面積和質量流量；ρ、f和h為熔鹽的密度、摩阻系數和焓；g為重力加速度；D為通道水力直徑；Qf和Qg分別為熔鹽和石墨體積熱源；Δz和Δp分別為節塊高度和壓降。

對式(18)和(19)耦合求解[14]可得到各通道內的流場和壓力場分布。對于單個燃料通道，石墨慢化劑可簡化為空心圓柱，只考慮徑向導熱，其能量守恒方程如下：

(21)

其中，Tg和λg分別為石墨溫度和其導熱系數。石墨內邊界采用熔鹽-石墨對流換熱條件，外邊界為絕熱邊界條件[5]。根據式(20)得到熔鹽溫度分布后，可解析求解式(21)得到石墨溫度分布[5]。

(22)

其中，V為積分域。

基于TMSR3D的有效緩發中子份額計算流程如圖2所示。具體為：1) 采用DRAGON5[15]計算各種組件的少群參數，并進行擬合；2) 物理熱工耦合計算；3) 共軛通量計算；4) 根據式(22)計算βeff。

圖2 基于數值方法的有效緩發中子份額計算流程Fig.2 Calculation scheme of effective delayed neutron fraction based on numerical method

2 計算結果與分析

2.1 緩發中子先驅核損失份額及其分布

額定工況下，MSRE裝載235U時DNP損失份額列于表3?？梢钥吹剑疚幕诮馕龇椒ê蛿抵捣椒ㄓ嬎闼媒Y果與其他程序[16]均符合較好。對于第2組DNP，由于其具有較大的裂變份額和較小的衰變常量，故損失份額最大。

表3 熔鹽實驗堆額定工況下緩發中子先驅核損失份額Table 3 Loss of delayed neutron precursor fraction of molten salt reactor experiment under rated condition

基于數值方法計算得到額定工況下DNP濃度的分布，如圖3所示，由于熔鹽的流動效應，每組DNP濃度的峰值均沿軸向有一定偏移，且流動對半衰期較長的DNP(如第2組)影響較大。這是因為半衰期越長，DNP在衰變前隨熔鹽行進的距離越長，甚至流出堆芯，造成反應性的損失。

圖3 緩發中子先驅核濃度分布Fig.3 Concentration distribution of delayed neutron precursor

圖4示出額定工況下最熱通道內DNP與緩發中子軸向的分布。其中，DNP濃度由數值程序TMSR3D計算得到，緩發中子分布由解析方法得到?？梢钥吹剑?組DNP歸一化濃度高于其他各組DNP。這是由于DNP的穩定分布不僅與裂變份額(β)有關，也與半衰期(1/λ)有關。裂變份額越大，半衰期越長，其穩態濃度越高，即DNP濃度正比于β/λ。由于第2組DNP具有最大的β/λ值，因此其濃度高于其他各組。緩發中子由DNP經β衰變得到，如無燃料流動效應，穩態下各組緩發中子的分布次序應與DNP份額的大小順序一致。對于液態燃料熔鹽堆，還需考慮流動對緩發中子濃度的影響，如對于第5組，盡管其DNP濃度較小(圖4a)，然而其具有較短的半衰期，DNP在堆內產生后短期內即衰變為緩發中子，因流動導致流出堆芯損失的緩發中子較少，因而第5組緩發中子具有較大的濃度(圖4b)。

圖4 額定工況下最熱通道內DNP與緩發中子軸向分布Fig.4 Axial distribution of DNP and delayed neutron in the hottest channel at rated condition

2.2 入口流量對βeff的影響

液態燃料熔鹽堆中燃料的流動導致引起DNP的空間再分布，部分DNP流出堆芯并衰變，從而對βeff產生影響。圖5示出采用解析方法和數值方法計算得到的不同額定流量百分比下堆芯βeff的變化?？梢钥吹?，兩種方法計算得到的βeff變化趨勢基本一致，即入口流量增加，βeff減小。這是由于流量增加，更多的DNP流出堆芯且在堆外衰變，導致了反應性損失，從而引起了βeff的減小。除0%額定流量外，解析方法得到的計算結果始終高于數值程序計算結果，二者最大相對誤差σmax為6.26%，額定流量下相對誤差σnormal為4.18%。

圖5 不同流量下的有效緩發中子份額Fig.5 Effective delayed neutron fraction at different percentages of flow rate

2.3 堆外流動時間對βeff的影響

對于半衰期較長的DNP，其隨燃料流出堆芯后可能未發生衰變重新進入堆芯，從而使堆芯反應性及βeff產生波動，這一波動主要與燃料在堆外回路的流動時間有關[17]。圖6示出采用解析方法和數值方法計算得到的堆芯βeff隨燃料在堆外流動時間的變化，此過程維持入口流量為額定值?？梢钥吹?，βeff隨燃料在堆外的流動時間的增加而減小，80 s后趨于穩定。這是因為隨著燃料在堆外流動時間的增加，更多DNP將在堆外衰變，返回堆芯的DNP減少，從而使βeff減小。由表2可知，第1組DNP的平均壽命最大，約為80.65 s，因此當堆外流動時間超過80.65 s時，所有DNP均將在堆外衰變，因而無任何DNP返回堆芯，βeff降低至最小值，趨于穩定。可以看到，60 s前，兩種方法計算結果存在一定偏差，解析方法仍高于數值方法計算結果，二者最大相對誤差為10.61%，額定工況下相對誤差為4.33%。結合2.2節討論可發現，解析方法計算得到的βeff相對于數值方法偏大。這是由于，解析方法忽略了堆芯的非均勻性和反射層的作用，采用均勻的流速場，認為堆內所有緩發中子具有相同的價值，即為1，而堆外所有緩發中子價值為0?？偟膩碚f，由于燃料的流動性，相對于先驅核衰變產生緩發中子的位置，先驅核產生的位置中子價值是比較高的[3]，這也就是說，解析方法總體上高估了堆內緩發中子的價值，導致βeff計算結果高于數值方法的計算結果。

圖6 有效緩發中子份額隨燃料在堆外回路流動時間的變化Fig.6 Effective delayed neutron fraction change with fuel residence time out of core

表2 緩發中子先驅核份額及其衰變常量Table 2 Delayed neutron precursor fraction and precursor decay constant

3 結論

有效緩發中子份額βeff的精確計算可為反應堆設計和安全分析提供合理準則和依據。燃料的流動性、特殊的傳熱方式和釷鈾燃料循環的使用導致液態燃料熔鹽堆βeff的計算方法與固態燃料反應堆不同。本文針對石墨慢化通道式熔鹽堆，分別基于解析方法和數值方法提出了計算βeff的數學模型，計算了MSRE在額定工況下的DNP損失份額和堆內DNP濃度分布，并分析了燃料在堆外流動時間和入口流量對βeff的影響，可得到以下結論。

1) 兩種方法均可對DNP行為提供合理描述。解析方法忽略了堆芯的非均勻性和反射層的作用，總體上高估了堆內緩發中子的價值，導致βeff計算結果相對于數值方法偏大。

2) 固定燃料在堆外流動時間，βeff隨入口流量的增加而減小。

3) 固定入口流量，βeff隨燃料在堆外流動時間的增加而減小，80 s后趨于穩定。