液態金屬冷卻快堆堆芯物理分析軟件LoongSARAX的驗證與確認

路 瑤,杜夏楠,李愛鑫,高杰豪,陳文杰,鄭友琦,王永平,吳宏春

(西安交通大學 能源與動力工程學院,陜西 西安 710049)

在第四代反應堆中,快堆憑借高的平均中子通量具有增殖裂變核燃料和嬗變長衰變期錒系元素兩大特點和優勢[1-2],是我國核能三步走戰略的重要環節。與傳統的核電廠熱堆相比,快堆通常設計得較為緊湊,以提高其核燃料增殖性能,因堆芯內功率密度高,需要采用換熱能力足夠強的冷卻劑帶走堆芯的熱量。目前在建或在運行的快堆均采用液態金屬鈉或鉛作為冷卻劑,這得益于液態金屬強大的載熱能力。同時采用液態金屬冷卻的快堆設計容易實現小型化,在海、陸、空、天以及特種同位素生產等領域均有迫切需求。

為提高堆芯性能,液態金屬冷卻快堆設計與傳統壓水堆相比,在燃料形式、能譜復雜程度以及堆芯尺寸規模等方面都具有非常顯著的差別,其表現出的特征使得適用于傳統壓水堆設計開發的堆芯設計軟件和方法不再適用。因此,快堆軟件的自主化研發對于我國實現快堆穩步發展、技術彎道超車具有重要意義。

1 核電軟件驗證與確認方法

世界上很多國家有不同的核電軟件,相對于國內,國外的核電起步較早,不同的大型核電技術公司有針對其核電軟件獨特且完善的驗證確認方法[3],如美國西屋公司將核電軟件的驗證與確認稱為合格性測試,法國阿海琺公司的核電軟件——堆芯燃料管理軟件分為兩部分,分別是組件程序APOLLO2和堆芯程序ARTEMIS,分別用數值驗證和基準題的方式進行驗證。

在核工業領域,有6種常用的軟件確認方法[4],分別是解析方法(BASICS)、程序對標(CODE)、電廠測量數據比對(PLANT)、實驗測量(TEST)、審查測試(LICENSING)和回歸測試(REGRESSION), 針對以上6種方法的使用頻率做統計,可知審查測試、電廠測量數據比對與實驗測量這3種確認方式應用相較廣泛[5]。因此,本文使用實驗測量的方式對LoongSARAX進行確認。

2 LoongSARAX的驗證與確認策略

2.1 LoongSARAX程序簡介

NECP-SARAX是西安交通大學核工程計算物理實驗室自主開發的先進反應堆中子學分析計算系統[6-8],LoongSARAX程序則是在NECP-SARAX程序的基礎上,針對液態金屬冷卻快堆的設計計算與安全審評需求,定制化開發的版本。該版本適用于裝載二氧化鈾、MOX以及金屬燃料,采用液態金屬鈉、鉛或鉛鉍作為冷卻劑的快中子反應堆堆芯物理分析計算,主要包括截面產生模塊TULIP[9]和堆芯穩態分析模塊LAVENDER[10]。

2.2 LoongSARAX程序的驗證與確認方法

LoongSARAX軟件包括組件程序和堆芯計算程序,若將軟件稱為系統,則組件程序和堆芯計算程序稱為子系統,組件程序包括組件共振計算模塊、組件輸運計算模塊等,堆芯計算程序包括中子輸運計算模塊和反應性系數計算模塊。

LoongSARAX的驗證流程如下:將整個程序分成不同的模塊,確定每個模塊驗證所需要的基準題并生成驗證矩陣,在此基礎上,再分為3個階段分別進行模塊驗證、子系統驗證和系統確認,按照驗證計劃完成驗證。LoongSARAX子系統各模塊及其驗證方法列于表1。

表1 LoongSARAX子系統驗證Table 1 LoongSARAX subsystem verification

本文對子系統中每個模塊進行單獨的驗證后,將所有模塊集成一個完整的計算程序,并對其進行軟件的系統確認。本文采用公開文獻以及評估反應堆物理基準實驗的國際手冊中的相關零功率實驗對LoongSARAX程序進行系統確認,其中臨界反應性、能譜特征、反應性效應、反應性系數、動力學參數以及反應率分布所采用的確認方法均為實驗測量。

本文在確認子系統驗證和系統確認的具體內容后,搜集每個模塊所需基準題,并建立驗證矩陣和確認矩陣,如表2、3所列。

表2 LoongSARAX驗證矩陣Table 2 LoongSARAX verification matrix

表3 LoongSARAX 確認矩陣Table 3 LoongSARAX validation matrix

3 LoongSARAX的驗證、確認及不確定度量化

3.1 組件程序驗證

TULIP組件程序主要驗證3個模塊,分別為共振計算模塊、輸運計算模塊、少群截面均勻計算模塊[11]。3個模塊的驗證均是通過對例題進行校算來實現的,例題的參考解kinf、能群截面及通量由MCNP程序計算產生。為驗證共振計算模塊和少群截面均勻計算模塊,利用MCNP程序統計1 968群以及均勻化33群的核素截面,并展開對比分析。

本文在驗證時選取了不同液態金屬冷卻劑(鈉、鉛)、不同富集度情況下均勻模型組件、非均勻模型的一維平板組件及六邊形組件。分別用TULIP程序和MCNP程序對均勻模型下各燃料組件進行計算并進行對比,包括kinf、總截面及裂變截面(235U、238U、239Pu、240Pu)、冷卻劑Na或Pb的總截面、56Fe的總截面以及歸一化能譜,截面及能譜都包括1 968群和33群。以ZPPR17A燃料組件均勻模型組件為例,部分截面計算結果如圖1所示,能譜計算結果如圖2所示。均勻模型所有組件的kinf計算結果如表4所列。C/E為計算值與實驗值的比值。

圖1 均勻燃料組件截面計算結果Fig.1 Calculation result of cross section for homogeneous fuel assembly

圖2 均勻燃料組件能譜計算結果Fig.2 Calculation result of energy spectrum for homogeneous fuel assembly

表4 均勻燃料組件kinf計算結果Table 4 kinf calculation result of homogeneous fuel assembly

從圖1可看出,中高能量段時兩個程序的計算結果吻合較好、偏差較小,低能量段時偏差較大。從圖2可發現,低能量段時通量接近于0,這是導致粒子輸運程序統計偏差較大的原因,因此低能量段的結果不作為參考。由于MCNP和LoongSARAX在處理高能區閾能反應時的方式不同,因此導致圖2在高能群位置相對偏差開始下降。其余核素的1 968群截面計算結果偏差的均方根在[0,2.1%]范圍內,33群截面計算結果偏差的均方根在[0,7%]范圍內。

在計算一維非均勻組件時,若是六邊形組件,則采用等效一維圓環模型,即根據體積內核子密度守恒,將不同材料等效成半徑不同的一維圓環在徑向上進行堆疊,如圖3所示,以從內到外第3圈為例,在面積不變的情況下,根據圓環面積的計算公式得出黃色區域的半徑。以鈉冷快堆MET1000一維模型ZPPR17A燃料組件為例,僅列出240Pu、23Na總截面計算結果,如圖4所示。能譜計算結果如圖5所示,一維模型所有組件的kinf列于表5。

圖3 等效一維圓環模型Fig.3 Equivalent one-dimensional ring model

圖4 非均勻組件截面計算結果Fig.4 Calculation results of cross section for heterogeneous fuel assembly

圖5 非均勻組件能譜計算結果Fig.5 Calculation results of energy spectrum for heterogeneous fuel assembly

表5 非均勻燃料組件kinf計算結果Table 5 kinf calculation result of heterogeneous fuel assembly

由圖4可看出,中高能量段TULIP程序與MCNP程序計算截面吻合程度較高、偏差較小,低能量段偏差過大。從圖5可看出,由于典型快譜的低能量段的中子通量非常低,導致低能量段的統計偏差很大,因此,低能量段的結果不作為參考解。

3.2 堆芯程序驗證

針對堆芯程序的驗證,重點關注其中子輸運模塊的計算精度,為此,本文將基于TAKEDA基準題對堆芯程序中的中子輸運模塊展開驗證。

以TAKEDA基準題2為例,采用程序對標的方法驗證輸運計算模塊。TAKEDA基準題2模擬的是一小型快中子增殖反應堆(FBR)[12]。小型快中子增殖反應堆含有4個材料區域,分別為燃料區、軸向增殖區、徑向增殖區和控制棒區。當控制棒提出時,采用鈉來填充空出的區域。

計算時考慮兩種情況:控制棒全提(情況1)和控制棒插入一半(情況2)。有效增殖因數及控制棒價值計算結果列于表6,平均中子通量密度偏差計算結果示于圖6,其中G表示能群。可看出,LAVENDER與MCNP程序所提供的參考值吻合很好,keff偏差小于40 pcm,控制棒價值偏差小于2%,各區域、能群的中子通量密度偏差均小于5%。

圖6 平均中子通量密度偏差計算結果Fig.6 Calculation bias results of average neutron flux density

表6 有效增殖因數及控制棒價值計算結果Table 6 Calculation results of keff and control rod value

3.3 LoongSARAX程序驗證

以OECD循環快堆和超鳳凰堆基準題為主要計算對象,在輸運、功率、燃耗、反應性效應等計算功能方面進行LoongSARAX程序驗證。

3.3.1OECD循環快堆 钚循環物理工作組(WPPR)指定了兩種快速燃燒器基準設計(1種氧化物和1種金屬)。本文通過計算金屬燃料反應堆進行驗證。金屬燃料反應堆堆芯的功率是1 575 MWth,周期長度為365 d,容量系數為85%。反應堆包含420個燃料組件和30個控制組件,并排列成一堆芯為1/6對稱的配置,燃料組件由鈾、钚組成[13]。

本文分為5個燃耗步且每個燃耗步時長為73 d,進行1個周期的循環。計算結果列于表7,其中BOL表示壽期初燃耗重核質量,EOL表示壽期末燃耗重核質量。從表7可知,LoongSARAX的計算結果和其他單位的計算結果相比,偏差較小。

表7 燃耗重核質量及其變化計算結果Table 7 Calculation results of burnup heavy core mass and its variation

3.3.2超鳳凰堆 超鳳凰堆是世界上運行的第一座商業規模的快中子增殖堆,也是至今已建成的最大的鈉冷卻的快中子反應堆,該反應堆裝載358個燃料組件,反應堆額定熱功率為3 000 MW,堆芯燃料區分區裝載有兩種不同的MOX燃料,增殖區裝載低富集度的二氧化鈾增殖燃料[14]。LoongSARAX通過超鳳凰堆計算驗證了輸運、功率及反應性效應。其中,參考解來自于基準題報告中Serpent的計算結果。

本文完成了13種堆芯狀態下的反應堆堆芯臨界計算。這13種堆芯狀態分別在幾何溫度、截面溫度、控制棒插入深度上存在差別。13種堆芯狀態的穩態計算結果與參考解之間的偏差列于表8。

表8 堆芯穩態計算結果Table 8 Core steady state calculation results

從表8可看出,LoongSARAX的計算結果和參考解相比不存在較大的偏差。除堆芯狀態較為復雜的算例9外,其他所有算例和參考解的偏差都在400 pcm以內。表明LoongSARAX程序在堆芯臨界計算上存在較高的計算精度。

堆芯功率分布計算基于的堆芯狀態為幾何膨脹溫度673 K,調節棒吸收體插入活性區高度的40%。計算結果詳見文獻[15],計算結果和參考解的相對偏差如圖7所示。

圖7 功率分布計算偏差Fig.7 Calculation bias of power distribution

如圖7所示,功率分布計算的計算結果和參考解相比在活性區不存在明顯的偏差。增殖區相對偏差較大,絕對偏差不大,可認為是由于增殖區功率較低導致的。因此,可證明LoongSARAX程序對功率分布的計算結果具有較高的計算精度。

本文將基于與堆芯臨界狀態相同的堆芯狀態計算堆芯的反應性常數和反應性系數,其中包括燃料多普勒常數、鈉密度系數、鈉空泡系數、燃料軸向膨脹系數、包殼膨脹系數、組件壁膨脹系數、柵格膨脹系數。具體計算假設和計算方法參考基準題報告[14],反應性系數和反應性常數的計算結果和計算偏差如表9所列。

表9 反應性系數和反應性常數計算結果Table 9 Calculation results of reactivity coefficient and reactivity constant

整體來看,LoongSARAX的計算結果和參考結果不存在較大的偏差。因此,可證明LoongSARAX程序在計算反應性系數和反應性常數時存在較高的精度。

3.4 LoongSARAX程序確認

本文以JOYO以及ZPPR17A兩個零功率實驗裝置為主要計算對象,對程序在確認過程中的對比參數進行介紹。

3.4.1JOYO反應堆 JOYO是日本的第一個實驗反應快堆,并在1977年達到了初始臨界[16]。JOYO的主要作用是改進快堆技術、進行輻照燃料,并創新將快堆投入實際使用的技術。JOYO是鈉冷快堆,采用六邊形組件,以鈾钚混合氧化物為燃料。在穩定的溫度和恒定鈉流量條件下,該反應堆主要通過調整控制棒的位置來實現臨界。當堆芯裝載70個燃料組件時,簡稱為JOYO-70,其中有兩根調節棒組件半插入堆芯。本文主要進行了臨界反應性、控制棒價值、鈉空泡反應性以及燃料替換反應性計算,結果列于表10～13。

表10 JOYO-70堆芯臨界計算結果Table 10 JOYO-70 core critical calculation results

從表10可看出,SARAX對臨界測量實驗的建模計算與實驗結果吻合較好,有效增殖因數偏差為365 pcm。

在本次計算中,采用反應性差值法進行計算控制棒價值。選取JOYO-70所有控制棒均為半插入狀態時的βeff作為本次計算采用的βeff,所有控制棒都撤出(ARO)時的有效增殖因數的測量結果為1.024 95,再分別計算每根控制棒插入時的keff,與ARO狀態的有效增殖因數作差,即可得控制棒價值,最終測量結果如表11所列。可看出,LoongSARAX計算控制棒價值的結果較好,偏差較小。

表11 JOYO-70控制棒價值計算結果Table 11 JOYO-70 control rod worth calculation results

通過將燃料組件替換為鈉空泡實驗組件的方式得到鈉空泡反應性,如表12所列。可看出,鈉空泡反應性數值的絕對值減小。由于[000]處于堆芯,燃料數量多,中子通量最多,由堆芯徑向往外中子通量減少,引入負反應性減少,[4f1]、[6F1]均在燃料區外,處于增殖區,移除鈉引起正反應性。實驗值與測量值偏差較小,所計算偏差均小于實驗不確定度。

表12 JOYO-70鈉空泡反應性計算結果Table 12 JOYO-70 calculation results of sodium void reactivity

燃料替換反應性的所有測量都是在低功率運行(溫度分布均勻)下進行的,平均冷卻劑溫度范圍為240～250 ℃。計算結果列于表13,可看出,兩者結果較為吻合。

表13 JOYO-70燃料替換反應性計算結果Table 13 JOYO-70 fuel replacement reactivity calculation results

3.4.2ZPPR17A反應堆 ZPPR17A是零功率物理反應堆[17],其特征是軸向不均勻,臨界堆芯狀態沒有控制棒位置。它被認為是JUPITER-Ⅲ系列實驗的參考堆型之一。ZPPR-17A的燃料是Pu-U-Mo板,鈉為冷卻劑,反射層材料為不銹鋼。本文采用LoongSARAX計算的臨界反應性為1.003 43,與實驗參考值的偏差為333 pcm,本文還計算了控制棒價值以及反應率分布,結果如表14及圖8所示。可看出,計算結果與參考值較為吻合。

a——239Pu在z=5.16 cm處沿x軸總反應率分布;b——238U在z=28.02 cm處沿x軸裂變反應率分布;c——235U在z=28 cm處沿y軸裂變反應率分布;d——238U在148-70位置處沿z軸俘獲反應率分布圖8 部分反應率分布計算結果Fig.8 Calculation results of partial reaction rate distribution

表14 堆芯控制棒價值計算結果Table 14 Core control rod worth calculation results

基準題報告中,在燃料組件或增殖組件中裝入金屬探測箔片,通過測量箔片的感生放射性活度,并結合已知的探測箔片位置和材料特性,可得到相應區域的反應率,包括軸向高度z=5.16 cm、z=28.02 cm處沿x軸反應率分布、計算組件位置位于148-50(堆芯中心位置)和148-70(外燃料區內,具體位置見文獻[17])的z軸反應率分布以及軸向高度z=28 cm處沿y軸的反應率分布(包括總反應率、裂變反應率、俘獲反應率)。利用LoongSARAX計算反應率,所得結果已歸一化。圖8中僅展示其部分反應率分布。可看出,反應率分布的計算值與實驗測量值的偏差較小。

3.5 程序確認結果分析

在3.3節的基礎上,LoongSARAX針對確認矩陣中的算例及其他鈉冷快堆如中國實驗快堆(CEFR)均開展了計算,結果表明,與實驗測量值相比,堆芯有效增殖因數計算偏差均小于400 pcm,其他反應性與參考解的偏差均小于15%。受限于篇幅,其他計算結果可參考文獻[18-19]。

特別地,針對CEFR的啟動物理試驗,在臨界反應性、控制棒價值、鈉空泡系數、等溫溫度系數以及燃料組件替代反應性測量過程中共產生了60余組臨界狀態。具體結果如圖9所示。

圖9 中國實驗快堆各臨界狀態計算結果偏差Fig.9 Calculation bias results of CEFR core critical reactivity

基于以上結果,在不依賴于總體分布類型的情況下,本文應用非參數統計法對CEFR的臨界反應性進行了不確定度量化[20]。

為了提高結果的準確性和可靠性,本文對偏差數據進行篩選,剔除可能由于失誤或不可控因素導致偏差數據存在異常值或數據誤差較大的情況。

如式(1)[21]所示,假設從具有連續累積分布函數F(x)的總體中抽取大小為n的樣本,將樣本值按增序排列為x1、x2、…、xn,對于包含在樣本中從第r個最小值xr到第s個最大值xn-s+1之間的比例可表示為F(xn-s+1)-F(xr)。這個比例稱為區間(xr,xn-s+1)的覆蓋率。如式(2)[21]所示,計算這個覆蓋率至少達到某個給定數值α的概率,α表示100%的總體被包含在區間xr和xn-s+1之間的概率,這個概率僅依賴于n和m。因此,在樣本量n一定的情況下,篩選的具體個數m根據置信度、概率查表[21]獲得。

un-m(1-u)m-1du

(1)

(2)

在99%的置信水平下,通過查表[21]可知,經剔除篩選偏差數據后,LoongSARAX對于鈉冷快堆CEFR的計算偏差有90%的概率位于[-389 pcm,300 pcm]以內。

4 總結

針對液態金屬冷卻快堆的研發需求,本文基于LoongSARAX軟件提出了堆芯核設計程序的驗證與確認策略,通過搜集國際上關于液態金屬冷卻快堆物理計算基準題建立了相應的驗證與確認矩陣,在針對LoongSARAX制定驗證確認流程后,開展了組件程序、堆芯程序各個模塊的驗證以及LoongSARAX程序確認的計算。

本文主要獲得了各堆芯及其在實驗過程中不同狀態點下的有效增殖因數,以及其他反應性的結果與偏差,并針對有效增殖因數開展了計算結果的不確定度量化。由計算分析可得,與實驗測量值相比,LoongSARAX對于臨界反應性、控制棒價值、鈉空泡價值和反應率分布等結果的計算精度較好,可說明LoongSARAX在液態金屬冷卻快堆的計算分析中能夠滿足多種中子學參數的分析需求且具備較高的計算精度。