核反應堆物理確定論計算程序系統NECP軟件包的研發與應用

吳宏春，曹良志，鄭友琦，李云召，祖鐵軍，劉宙宇

(西安交通大學 核科學與技術學院，陜西 西安 710049)

近年來，隨著核反應堆用途的不斷拓展，人們對核能安全性和經濟性的要求不斷提高，核反應堆的設計越來越趨于復雜化和精細化，對傳統的核反應堆設計技術提出了諸多新的挑戰，特別是對核反應堆物理計算方法和數值工具提出了更高的要求。在上述需求的驅動下，國際核強國在過去20年的時間中對其物理計算方法和工具進行了深入、廣泛的改進，甚至重新研發了新一代的反應堆物理計算程序系統以應對新堆型研發和設計要求提高對物理計算所提出的挑戰，具體包括以下兩方面。

1) 物理計算方法和程序的通用化。目前反應堆物理程序開發的進展表明，物理計算程序具有更寬的能譜適用范圍，滿足不同反應堆的物理計算需求已是大勢所趨。國際知名的核設計程序如美國西屋公司開發的APA、法國CEA開發的APOLLO-3等均已發展為能同時用于熱堆/快堆的物理計算。

2) 物理計算建模的精細化、精確化。隨著數值反應堆概念的提出并被廣泛接受，基于高性能計算條件的“一步法”反應堆計算方法成為另一重要的研究方向，為實現對核反應堆堆芯的高分辨率建模和高精度計算，產生了一批國際知名的“高保真”物理計算程序，如美國阿貢國家實驗室的PROTEUS及密歇根大學的MPACT等。

國內在核反應堆設計方法及相應的計算程序開發上一直受到國外的技術封鎖，為打破封鎖、實現核電軟件的完全自主化，我國在最近10年間開展了大量針對國內在建、在研壓水堆的自主化物理程序研發，推出了TORCH、PCM、cosNU等自主化程序[1-2]。然而，上述程序在通用化、精細化或精確化等方面仍有待進一步加強。

西安交通大學反應堆物理團隊在核數據及共振理論、中子輸運理論及數值方法等方面開展了大量的基礎研究，發展了子群方法[3]、超細群方法[4]，提出了小波展開方法[5]、小波-子群耦合方法[6]、全局-局部耦合方法[7]等共振計算方法；發展了模塊化特征線方法[8]、2D/1D特征線方法[9]，提出了非結構幾何SN節塊方法[10]等中子輸運數值方法，并最終形成了一整套具有強通用性和高精度的物理計算程序系統，稱之為NECP軟件包。NECP軟件包主體包含NECP-Atlas、NECP-Bamboo、NECP-X和NECP-SARAX(簡稱為Atlas、Bamboo、X和SARAX)等自主化程序。Atlas程序采用改進的截面處理方法制作高精度的堆用多群核數據；Bamboo程序基于組件均勻化或柵元均勻化，采用改進的子群方法、特征線方法和變分節塊方法進行壓水堆燃料管理計算；X程序基于全堆一步計算，采用全局-局部耦合共振計算方法、2D/1D特征線方法進行輕水堆高保真計算；SARAX程序基于快譜均勻化思想，采用超細群方法、輸運節塊方法進行快譜反應堆計算。NECP軟件包可實現對各類壓水堆和快堆的高精度核設計和燃料管理計算，并具備了進行高保真物理計算的能力。

本文從NECP軟件包的核心算法模型出發，對該軟件包內主體程序的特點和驗證與確認情況進行介紹，并對軟件包在實際工程計算中的應用情況進行代表性的闡述。

1 核數據處理程序Atlas

1.1 程序簡介

Atlas程序[11]是為滿足物理計算程序的堆用核數據需求而開發的核數據處理程序，可將評價核數據庫中的數據加工為Bamboo、X、SARAX所需的核截面。該程序可處理ENDF/B-Ⅷ.0、ENDF/B-Ⅶ.1、ENDF/B-Ⅶ.0、ENDF/B-Ⅵ.8、CENDL-3.1、JEFF-3.2及JENDL-4.0等評價核數據，輸出WIMS-D、MATXS和ACE格式數據庫。

1.2 理論介紹

Atlas程序可將評價核數據庫中0 K下的數據轉化到所需溫度下的連續截面或多群截面，通過對傳統核數據處理程序中使用的近似模型進行改進優化，提高了堆用核截面數據的精度。

對于熱譜反應堆，由于重核素的熱運動及其彈性散射截面的共振，使得中子散射后有很大的概率獲得能量發生上散射。傳統核數據處理程序中并未對共振彈性散射進行專門的處理，而是假設超熱能區靶核處于靜止狀態并且其彈性散射截面為常數，這樣處理后中子發生散射后總是損失能量。針對該問題，Atlas程序采用不同階數的散射矩進行多普勒展寬，精確考慮靶核熱運動及散射截面的共振對散射中子能量分布的影響。在超熱能區可認為靶核不受化學鍵束縛，將靶核的運動按自由氣體模型處理，可得到不同溫度T下各階散射矩的計算式：

(1)

對于快譜或中間能譜反應堆，不可分辨共振區核數據的處理對后續計算精度的影響顯著。Atlas程序使用Multi-Level Breit-Wigner公式。該公式可避免采用Single-Level Breit-Wigner公式時出現的負彈性散射截面現象，提高了結果的魯棒性，且可實現對競爭反應概率表的計算，方便堆物理計算程序實現閾能反應的自屏計算。此外，對于積分統計方法，Atlas程序對波動積分中的孤立項、序列內干涉項和序列間干涉項進行了精確的計算，提高了有效自屏截面的精度。

1.3 程序驗證

對Atlas程序進行了廣泛的基準題驗證，包括國際原子能機構的WLUP基準題、OECD的ICSBEP基準題、IRPhEP基準題及VERA基準題等。反應堆的能譜涉及快譜、中間能譜和熱譜。圖1示出采用Atlas程序制作的數據庫和NJOY2016制作的數據庫對ICSBEP基準題的計算結果，計算程序為MCNP。由圖1可見，MCNP使用兩個數據庫的計算結果偏差小于25 pcm。采用BEAVRS基準題驗證了Atlas程序中共振彈性散射效應的處理精度，表1列出3種計算結果的對比，包括采用Atlas程序制作的未考慮共振彈性散射的數據庫(方法1)、考慮共振彈性散射的數據庫(方法2)和直接采用多普勒舍棄修正方法制作的數據庫(方法3)，計算程序均采用OpenMC程序。表1中，偏差1表示方法1、2間的偏差(方法2的值-方法1的值)，偏差2表示方法2、3間的偏差(方法3的值-方法2的值)。結果表明，對共振彈性散射的準確處理將使壓水堆keff計算誤差降低100～200 pcm。對不可分辨共振區截面處理的改進相比于傳統方法可得到更準確的多普勒常數，計算結果列于表2。

圖1 ICSBEP基準題計算結果Fig.1 Calculation result of ICSBEP benchmark

表1 BEAVRS基準題 keff計算結果Table 1 keff calculation result of BEAVRS benchmark

注：括號內為蒙特卡羅程序OpenMC的統計偏差

表2 鈉冷快堆燃料組件無限增殖因子kinf和多普勒常數KD的結果Table 2 kinf and KD results of subassembly in sodium-cooled fast reactor

2 壓水堆堆芯燃料管理兩步法計算程序Bamboo

2.1 程序簡介

Bamboo程序[12-13]由二維柵格計算程序Bamboo-Lattice、三維穩態分析程序Bamboo-Core和三維瞬態分析程序Bamboo-Transient組成。程序可用于壓水堆核設計和燃料管理計算，具有臨界和目標軸向功率偏移搜索、啟動物理試驗模擬、基負荷運行模擬、功率歷史跟蹤模擬及標簽換料計算等主要功能。目前，Bamboo程序已開發了兩個版本，其中，Bamboo1.0是基于組件均勻化的計算程序，Bamboo2.0是基于更加先進的柵元均勻化的pin-by-pin計算程序。

2.2 理論模型

Bamboo程序基于由ENDF/B-Ⅶ評價核數據庫經由Atlas程序加工制作的69群多群截面數據庫、核素燃耗數據庫和動力學參數庫，包含了258種核素，可直接分析硼、銀-銦-鎘、鎢、鉿、釓等多種中子吸收體材料。采用全組件一步子群共振計算方法，通過式(2)給出有效自屏截面：

(2)

組件輸運計算采用組件模塊化特征線方法，顯式考慮涂硼可燃毒物(IFBA)和水隙等重要精細結構的影響。采用顯式中子俘獲釋熱方法，如式(3)所示，以提高對可燃毒物棒功率計算的精度。同時采用高階燃耗方法，實現對各種可燃毒物燃耗效應的精細計算。

(3)

其中：P′m，g為第g群中子在第m個區域內核反應發出的功率；m為空間區域索引號；κ′i為第i個核素單次裂變釋放的能量；Ni為第i個核素的原子核密度；σf，i，g為第g群第i個核素的微觀裂變截面；qi為第i個核素的中子俘獲釋能；σc，i，g為第g群第i個核素的微觀俘獲截面。

堆芯擴散計算采用變分節塊方法，直接在三維空間內將中子通量密度進行正交基函數展開，如式(4)所示。利用變分原理獲得相應展開系數矩的離散代數方程組，避免了橫向積分，可完全消除控制棒尖齒效應，能直接給出節塊內連續的中子通量密度分布。

(4)

其中：φg(r)為在空間位置r處第g群的中子通量密度；φi，g為第g群中子通量密度的第i階展開矩；fi(r)為正交基函數系的第i階基函數；f為fi(r)構成的向量函數；φg為φi，g構成的向量。

2.3 驗證與應用

Bamboo程序經過了大量模塊驗證、子程序驗證和系統驗證，已成功應用于BEAVRS、Wats-Bar、M310、CNP300、華龍一號、AP1000和EPR等多種型號壓水堆的堆芯物理計算。以AP1000堆芯為例，表3列出首次臨界試驗參數的計算結果，表4列出控制棒組積分價值的計算結果。與實際測量值相比，計算結果均滿足核安全審評的要求，表明Bamboo達到了與AP1000核設計程序APA系統相當的計算精度。

表3 啟動物理試驗參數的計算結果Table 3 Calculation result of startup physical test parameter

表4 控制棒組積分價值的計算結果Table 4 Calculation result of control rod integral worth

3 輕水堆全堆一步法高保真物理程序X

3.1 程序簡介

X程序[15]是面向輕水堆開發的基于一步法計算的反應堆物理計算程序。程序基于高分辨率精細化建模、全堆共振計算、全堆大規模并行輸運計算及pin-by-pin物理熱工耦合實現了全堆芯高保真物理計算，給出了反應堆全三維精細的功率分布。

3.2 理論方法

X程序包括基于ENDF/B-Ⅵ.8和ENDF/B-Ⅶ.0評價數據庫經由Atlas程序加工制作的多群數據庫和連續能量數據庫。為能精細處理復雜的反應堆堆芯結構，采用構造實體幾何(CSG)建模方法，可直接描述包括格架、堆芯圍板、吊籃等堆內構件。

為實現準確高效的全堆芯共振計算，提出了全局-局部耦合共振計算方法，利用中子流方法獲得每個燃料棒的丹可夫因子Cb：

(5)

其中：φf，1為真實幾何通過特征線方法(MOC)計算的中子通量密度；φf，0為利用碰撞概率方法求解的孤立棒問題中子通量密度。利用丹可夫因子守恒，將全堆芯共振問題等效為一系列一維圓柱問題，通過共振偽核素子群方法求解等效后的一維圓柱問題，獲得每個燃料棒的共振自屏截面。

考慮到堆芯的非均勻性主要集中在徑向，輸運計算采取了基于泄漏項分割技術的2D/1D耦合方法，其中徑向和軸向方程如式(6)、(7)所示。利用MOC求解二維方程、SN方法求解一維方程。

Σt,g,k(x,y)ψg,m,k(x,y)=

(6)

(7)

為提高計算效率，采用空間并行、角度并行和特征線并行三重并行技術，實現大規模并行計算。采用無矩陣粗網有限差分(CMFD)加速方法，通過掃描方式計算矩陣向量積，實現CMFD方程的高效求解和有效加速。

X程序擁有精細燃耗數據庫(含有1 547個核素)和經過優化壓縮的燃耗數據庫(包含233個核素)以滿足不同應用需求。利用切比雪夫有理近似方法(CRAM)進行點燃耗方程的求解。在此基礎上通過預估校正、子步方法實現輸運計算和燃耗計算之間的耦合?；诎魞茸訁^尺度實現了精細的物理計算和pin-by-pin子通道熱工計算的耦合，結合內部耦合和外部耦合的混合耦合方式實現了耦合數據的高效傳遞。

3.3 驗證與應用

利用CASL計劃的VERA系列基準題對X程序進行了詳細的對比驗證。計算基于全堆精細幾何建模，如圖2、3所示，直接給出了每根燃料棒的功率。結果表明：棒功率最大相對偏差和均方根偏差為：二維組件分別在1.5%和1.0%以內，二維堆芯分別在2.5%和1.0%以內。表5列出VERA系列基準題中代表性的VERA5三維堆芯算例的詳細結果。由表5可見，三維堆芯特征值偏差在100 pcm以內。組件徑向歸一化功率最大相對偏差為1.38%，均方根偏差為0.69%。

圖2 堆芯布置Fig.2 Reactor core layout

圖3 精細幾何建模Fig.3 High-fidelity geometry modeling

表5 VERA5三維堆芯計算結果比較Table 5 Calculation result comparison of 3D core of VERA5

將X程序應用到大型商用壓水堆Watts Bar熱態滿功率的計算分析中，獲得了pin-by-pin分辨率的燃料棒溫度分布、冷卻劑溫度分布，如圖4所示。與基于蒙特卡羅的耦合程序系統MC21/CTF相比，組件徑向歸一化功率的最大相對偏差和均方根偏差分別為1.7%和0.8%。

圖4 燃料棒和冷卻劑溫度的三維分布Fig.4 3D temperature distribution of fuel rod and coolant

4 快譜反應堆物理計算程序SARAX

4.1 程序簡介

SARAX程序[16-17]是針對快譜反應堆和局部裝載固體慢化劑的中間能譜反應堆研發的物理計算程序，由快譜/中間能譜專用核數據庫、截面產生程序TULIP、三維穩態堆芯程序LAVENDER組成。

4.2 理論介紹

SARAX程序基于ENDF/B-Ⅶ.0評價數據庫經由Atlas程序加工制作的點截面數據庫和超細群數據庫。為處理快堆復雜的共振效應，采用基于點截面的方法計算有效自屏截面。超細群采用1 968群的能群結構，在每個能群中采用插值的方式產生200個能量點的截面數據。鑒于在中高能量區間，中子能譜很好地滿足窄共振近似，因此超細群的有效自屏截面可由式(8)獲得：

(8)

為處理組件的非均勻效應，基于Wigner單項有理近似，假設在1 968群各能群內的逃脫截面為常數，則非均勻系統下的有效自屏截面通過式(9)計算：

(9)

為準確計算快堆堆芯強烈的輸運效應，本文采用非結構幾何的SN節塊方法求解中子輸運方程以獲得堆芯中子通量密度分布：

Σtψm(x,y,z)=Qm(x,y,z)

(10)

4.3 驗證與應用

為驗證SARAX程序的正確性及對于鈉冷快堆計算的適用性，分別采用數值基準題和反應堆物理實驗基準算例的測量數據對程序進行了驗證與確認。在驗證方面，主要基于IAEA及OECD發布的各類快堆基準算例。數值結果表明，對于各基準題算例，與蒙特卡羅程序計算的參考解相比，keff計算偏差在100～300 pcm范圍內，功率相對偏差在3%范圍內，反應性系數相對偏差不超過7%。在確認方面，針對日本長陽實驗堆(JOYO)、美國零功率物理試驗堆(ZPPR)及俄羅斯大科學裝置(BFS)等物理試驗進行了建模計算。計算結果表明，各臨界試驗的計算結果標準差為170 pcm，各類反應性計算結果計算的標準差小于5%，燃耗效應計算相對誤差小于10%。

在此基礎上，將SARAX程序應用于法國鳳凰快堆壽期末物理試驗與中國實驗快堆(CEFR)啟動物理試驗的模擬。表6、7分別列出鳳凰快堆提棒臨界試驗和CEFR首次臨界試驗的模擬結果。數值結果表明，SARAX程序對鳳凰快堆的堆芯臨界計算誤差小于200 pcm，達到甚至優于法、俄專用核設計程序的精度。對CEFR的臨界計算誤差小于100 pcm。

表6 鳳凰快堆控制棒提棒臨界的模擬結果 Table 6 Simulation result of critical state at control rod withdrawal of Phoenix reactor

表7 CEFR凈堆臨界的模擬結果Table 7 Simulation result of critical state of CEFR

5 結論

本文介紹了西安交通大學反應堆物理團隊開發的、具有完全自主知識產權的核反應堆物理確定論計算程序系統NECP軟件包的核心理論模型和主要驗證及應用情況。該軟件包是在西安交通大學20多年來堆用核截面處理、共振理論、中子輸運計算方法等方面研究的基礎上，面向現有壓水堆技術改造和新型反應堆研發需求開發的通用化確定論程序系統。軟件包采用先進的理論模型和現代軟件工程開發架構，形成了標準化、模塊化的計算體系。所開展的大量驗證和確認結果表明，該軟件包功能齊全、計算精度高，可滿足以壓水堆為代表的熱譜堆和以鈉冷快堆為代表的快譜堆的計算要求，同時對于壓水堆具備了高保真物理計算的能力。

NECP軟件包已應用于國內大型壓水堆項目、中國示范快堆項目等，國內用戶基本包括了各大核電公司和研究院所，國外用戶已達9個。已成功應用于AP1000啟動物理試驗的計算，預測的臨界硼濃度與測量值相比誤差僅為30 ppm；應用于CEFR啟動物理試驗的計算，預測的臨界堆芯裝載與實際測量值相比小于30 pcm。應用情況表明，NECP軟件包的計算精度達到甚至優于國外先進的核設計專用程序。