壓水堆堆芯Pin-by-pin計算廣義等效均勻化方法研究①

張斌 李云召 吳宏春 劉勇 王冬勇 王星博 鐘旻霄

(1.中國核動力研究設計院核反應堆系統設計技術重點實驗室 四川成都 610213；2.西安交通大學核科學與技術學院 陜西西安 710049)

20世紀70年代以后，得益于粗網節塊法[1,2]的迅速發展，以組件均勻化理論[3,4]和粗網節塊方法為理論框架的兩步法計算方案逐漸成為了壓水堆工程計算中普遍采用的燃料管理中子學數值計算方法。伴隨著粗網節塊法的大量工程應用，適用于組件均勻化的均勻化理論得到了長足的發展。基于均勻化一般原理，結合早期有限差分堆芯計算的傳統均勻化方法，Kord Smith基于等效均勻化理論提出了能夠用于工程計算的廣義等效均勻化方法[5]，即在等效均勻化參數中引入不連續因子以達到各節塊的守恒條件，此化方法很好地滿足了絕大部分商用壓水堆的工程應用需求。

隨著科學研究的不斷深入、計算條件和對核設計計算精度要求的不斷提高，傳統的兩步法計算方案面臨著越來越嚴峻的挑戰。為提高堆芯設計計算精度且滿足對新型反應堆堆芯計算的精度要求，基于均勻化理論的改進型兩步法計算方案全堆芯Pin-by-pin計算成為了下一代堆芯數值計算方法的研究熱點。全堆芯Pin-by-pin計算減少了堆芯計算過程中的近似與假設，能夠更加精細地考慮堆芯布置的非均勻性，直接求出單棒功率分布，便于堆芯燃料管理計算和相應的堆芯安全分析。

區別于組件均勻化的廣義等效均勻化方法，Pin-bypin均勻化計算中，不連續因子的計算[6]不僅與非均勻介質本身及其所處環境有關，也與堆芯計算方法相關。本文研究了廣義等效均勻化方法在壓水堆堆芯Pin-by-pin計算中的理論模型，針對矩形指數函數展開擴散方法和SP3方法推導了零階不連續因子和二階不連續因子的計算公式，通過SP3方程形式構造偽固定源問題解決了高階輸運方程無法提供SP3方程中二階中子通量密度的問題，并分析了壓水堆堆芯Pin-by-pin計算中應用Pin-by-pin不連續因子的堆芯計算精度。

1 不連續因子的計算方法研究

廣義等效均勻化方法通過放寬節塊表面中子通量密度這一約束條件，來實現反應率及中子泄漏率的守恒，并提出了新的節塊與節塊之間的邊界條件，即均勻化后堆芯內界面上中子通量密度乘上不連續因子后保持連續，第i節塊第g能群的不連續因子定義式如下：

式中：f——不連續因子；s——節塊表面；——節塊非均勻中子面通量，單位（cm2·s）-1；——節塊均勻中子面通量，單位（cm2·s）-1。

在柵格計算中，由高階輸運計算可以直接求得各個柵元各表面的非均勻中子面通量：

而柵元均勻中子面通量的計算方法則必須要與下游堆芯計算中所使用的中子學求解方法保持一致。

1.1 矩形指數函數展開擴散方法

本文堆芯低階輸運計算采用矩形指數函數展開節塊方法，在矩形節塊中，中子通量密度解析解取六階近似，展開方向為各坐標軸的正負方向，則：

式中：Ci——展開系數；S——中子源項；D——擴散系數；Σr——移出截面；在保證右手系的前提下，。

公式(3)中的六個積分常數可由包含節塊平均中子通量密度和節塊五個邊界面上凈中子流密度（以出射流為正）在內的六個定解條件確定。

當堆芯中子學計算采用指數函數展開擴散方法求解時，根據公式(3)和表面出射中子流密度的定義式，可以得到指數函數展開節塊方法內第i節塊第g能群入射中子流密度與出射中子流密度之間的響應關系，以u+方向為例：

式中：h——柵元長度，單位cm；out——出射中子流；in——入射中子流；——系數；S——中子源。

圖1 KAIST基準題二維全堆芯布置和幾何結構

圖2 KAIST基準題二維全堆芯棒功率分布

根據菲克定律中柵元中子面通量、柵元凈中子流密度與出射入射中子流密度的關系可以得到如下關系式：

根據均勻化前后特征值、各能群反應率及中子泄漏率守恒，則在單組件柵格計算中有如下等式：

聯立公式(5)和公式(6)：

聯立公式(1)、公式(2)和公式(7)即可求出指數函數展開擴散方法下柵元的不連續因子。

1.2 矩形指數函數展開SP3方法

當堆芯中子學計算采用矩形指數函數展開SP3方法求解時，通過與擴散方法中相似的推導，根據均勻化前后特征值、各能群反應率及中子泄漏率守恒，零階均勻化面通量與二階均勻化面通量表示如下：

圖3 組件均勻化粗網節塊計算的棒功率誤差分布

柵格非均勻一步法計算通常采用高階輸運方法進行計算，公式(8)和公式（9）中非均勻零階中子通量密度和零階中子流密度即為柵格輸運計算直接獲得，但非均勻二階中子通量密度和二階中子流密度無法無法由柵格計算直接提供。根據SP3方程的形式，非均勻二階中子通量密度矩可由公式(10)計算得到。通過高階輸運的柵格計算得到非均勻零階中子通量密度和零階中子流密度，根據公式(10)右端項構造零階中子通量密度的泄漏項作為該方程的固定源項，通過求解此偽固定源問題即可近似求得非均勻二階中子通量密度和二階中子流密度。

通過偽固定源方法近似求解出柵格計算中的非均勻二階中子通量密度和二階中子流密度后，根據公式(8)和公式(9)可以求出各階柵元均勻中子面通量。需要注意的是二階面通量的數值很小，且數值有正有負，而偽固定源方法是種近似方法，若按照不連續因子的定義式(1)進行SP3方程中二階通量的不連續因子的求解會出現很大的數值甚至是負值。不合理的不連續因子會導致全堆芯SP3方程迭代求解時的不收斂。為避免出現這種情況，對二階中子通量密度的不連續因子的求解關系式重新定義如下：

式中：s——節塊表面。

根據公式(1)和公式(11)可以得到SP3方程中零階和二階不連續因子，相應地，在SP3方法進行全堆芯Pin-by-pin計算求解時，節塊內部邊界的耦合條件應為（以u方向為例）：

2 計算結果分析

為驗證超級均勻化方法在壓水堆堆芯Pin-by-pin計算中的計算精度，本節以KAIST基準題2A問題[7]下的二維全堆芯問題進行驗證與分析，堆芯布置和幾何結構如圖1所示。

圖4 堆芯Pin-by-pin擴散計算的棒功率相對誤差分布

圖5 堆芯Pin-by-pin SP3計算的棒功率相對誤差分布

表1 堆芯兩步法計算的特征值誤差和棒功率誤差

采用堆芯Pin-by-pin計算和堆芯組件均勻化粗網計算兩種方案進行堆芯計算并加以分析和比較：

（1）依據廣義等效均勻化理論模型，針對擴散方法（P1）和簡化球諧函數方法（SP3）利用柵元不連續因子（PDF）分別產生兩群和七群的等效均勻化少群常數，并在指數函數展開節塊法中進行全堆芯Pin-by-pin計算；

（2）根據傳統組件均勻化計算方法，產生兩群帶組件不連續因子（ADF）的組件等效均勻化常數并采用先進變分節塊方法進行堆芯計算。變分節塊方法在空間處理上直接針對中子通量密度分布進行離散，因此不需要功率重構過程而能直接得到中子通量密度分布，通過形狀因子的求解即可得到堆芯棒功率分布。

二維堆芯參考解由Bamboo-Lattice2.0一步法計算獲得，特征值參考解為0.97948，棒功率分布如圖2所示。堆芯兩步法計算的特征值誤差和棒功率誤差如表1所示。由表1的結果可知：堆芯Pin-by-pin SP3計算的計算精度普遍高于組件均勻化的堆芯粗網節塊計算；除由于誤差抵消導致計算精度提高的兩群擴散計算外，堆芯Pin-by-pin擴散計算的計算精度都低于組件均勻化的堆芯粗網節塊計算。因此，SP3方法是堆芯Pin-by-pin計算的必要條件；隨著能群數的增加，堆芯Pin-by-pin簡化球諧函數計算精度顯著提高。七群SP3計算中特征值誤差小于50pcm，最大棒功率相對誤差低于2%，相對棒功率誤差的AVG，MRS，MRE在1%以下；而在組件均勻化的粗網節塊計算中，特征值誤差為219pcm，最大棒功率相對誤差接近10%。因此，七群是適用于堆芯Pin-by-pin計算精度的能群數目。

組件均勻化粗網節塊計算的棒功率誤差分布如圖3所示，堆芯Pin-by-pin擴散計算和SP3計算的棒功率相對誤差分布分別如圖4和圖5所示。可以發現，最大棒功率相對誤差出現在組件與組件交界面以及燃料柵元與反射層交界面附近。這是堆芯內組件真實環境與單組件計算中全反射邊界條件之間的不同所導致的。

3 結語

本文在壓水堆堆芯Pin-by-pin計算中采用廣義等效均勻化方法作為均勻化技術，對矩形指數函數展開擴散方法和SP3方法推導了零階不連續因子和二階不連續因子的計算公式，產生了適用于擴散計算和SP3計算的Pin-bypin等效均勻化參數。基于KAIST基準題，分析了壓水堆堆芯Pin-by-pin計算中應用Pin-by-pin不連續因子的堆芯計算精度，與傳統組件均勻化計算方法相比，應用了廣義等效均勻化方法的壓水堆堆芯Pin-by-pin計算的計算精度更高。7群簡化球諧函數方法是高精度堆芯Pin-by-pin計算方案的實施策略，與KAIST基準題一步法計算結果相比，其特征值誤差小于50pcm，最大棒功率相對誤差小于2%，相對棒功率誤差的AVG，MRS，MRE都在1%以下。