張 斌，李云召，吳宏春，王冬勇，劉 勇，于穎銳

(1.中國(guó)核動(dòng)力研究設(shè)計(jì)院 核反應(yīng)堆系統(tǒng)設(shè)計(jì)技術(shù)重點(diǎn)實(shí)驗(yàn)室，四川 成都 610213；2.西安交通大學(xué) 核科學(xué)與技術(shù)學(xué)院，陜西 西安 710049)

20世紀(jì)70年代以后，得益于粗網(wǎng)節(jié)塊法[1-2]的迅速發(fā)展，以組件均勻化理論[3-4]和粗網(wǎng)節(jié)塊方法為理論框架的兩步法計(jì)算方案逐漸成為壓水堆工程計(jì)算中普遍采用的燃料管理中子學(xué)數(shù)值計(jì)算方法。伴隨著粗網(wǎng)節(jié)塊法的大量工程應(yīng)用，適用于組件均勻化的均勻化理論得到了長(zhǎng)足的發(fā)展?；诰鶆蚧话阍韀5]，結(jié)合早期有限差分堆芯計(jì)算的傳統(tǒng)均勻化方法，Smith基于等效均勻化理論提出了能用于工程計(jì)算的廣義等效均勻化方法[6-7]，即在等效均勻化參數(shù)中引入不連續(xù)因子以達(dá)到各節(jié)塊的守恒條件，此方法很好地滿足了絕大部分商用壓水堆的工程應(yīng)用需求。

隨著科學(xué)研究的不斷深入以及計(jì)算條件和對(duì)核設(shè)計(jì)計(jì)算精度要求的不斷提高，傳統(tǒng)的兩步法計(jì)算方案面臨著越來越嚴(yán)峻的挑戰(zhàn)。為提高堆芯設(shè)計(jì)計(jì)算精度且滿足對(duì)新型反應(yīng)堆堆芯計(jì)算的精度要求，基于均勻化理論的改進(jìn)型兩步法計(jì)算方案全堆芯Pin-by-pin計(jì)算[8]成為下一代堆芯數(shù)值計(jì)算方法的研究熱點(diǎn)。全堆芯Pin-by-pin計(jì)算減少了堆芯計(jì)算過程中的近似與假設(shè)，能更加精細(xì)地考慮堆芯布置的非均勻性，直接求出單棒功率分布，便于堆芯燃料管理計(jì)算和相應(yīng)的堆芯安全分析。

有別于不連續(xù)因子在傳統(tǒng)組件均勻化計(jì)算中的廣泛應(yīng)用，壓水堆堆芯Pin-by-pin均勻化計(jì)算中除廣義等效均勻化方法外，超級(jí)均勻化(SPH)方法亦是其主流的均勻化技術(shù)之一。SPH方法由Kavenoky[9]于20世紀(jì)80年代提出。Herbet基于均勻化原理中特征值、各能群反應(yīng)率及中子泄漏率3個(gè)守恒關(guān)系不斷對(duì)SPH方法進(jìn)行完善，最終使其成為Pin-by-pin計(jì)算中的主流均勻化方法[10-11]。SPH方法的核心思想是放寬中子泄漏率守恒這一約束條件，通過SPH因子直接調(diào)整柵元均勻化截面使得均勻化前后各能群反應(yīng)率保持守恒。這種方法不需要均勻化少群常數(shù)中保存額外的均勻化常數(shù)(不連續(xù)因子)，在有很大內(nèi)存需求量的全堆芯Pin-by-pin計(jì)算中具有明顯優(yōu)勢(shì)。近年來，國(guó)內(nèi)外學(xué)者對(duì)SPH方法在程序?qū)崿F(xiàn)、能譜干涉修正、中子泄漏效應(yīng)處理等方面進(jìn)行了大量研究[12-19]，日本AEGIS/SCOPE2、德國(guó)DYN3D等程序已采用SPH方法作為其Pin-by-pin計(jì)算的均勻化方法。SPH方法不受具體堆芯計(jì)算模型限制，既可保證均勻化前后輸運(yùn)-擴(kuò)散守恒，也可保證輸運(yùn)-輸運(yùn)守恒，對(duì)于不同的均勻化過程體現(xiàn)在SPH因子修正少群常數(shù)時(shí)的不同，是壓水堆堆芯Pin-by-pin計(jì)算中首選的均勻化技術(shù)之一。

本文采用SPH方法計(jì)算燃料組件的Pin-by-pin等效均勻化參數(shù)，針對(duì)存在中子泄漏現(xiàn)象的反射層參數(shù)計(jì)算研究與空間泄漏相關(guān)的SPH方法。在保證反應(yīng)率守恒的基礎(chǔ)上，同時(shí)保證各柵元各能群的中子泄漏率守恒，解決存在中子泄漏時(shí)SPH因子迭代計(jì)算的不收斂問題，并分析壓水堆堆芯Pin-by-pin計(jì)算中應(yīng)用SPH因子的堆芯計(jì)算精度。

1 方法

1.1 SPH方法

均勻化少群截面的計(jì)算通常采用通量體積權(quán)重方法，如式(1)所示。

x=a,f,s,…

(1)

式中：g、h為能群編號(hào)；i為柵元編號(hào)；V為柵元體積；φ為反應(yīng)堆中子通量密度，cm-2·s-1；Σ為宏觀截面，cm-1；上標(biāo)hom表示均勻化后，het表示均勻化前。

式(1)假設(shè)柵格計(jì)算得到的非均勻中子通量密度等于均勻中子通量密度。此假設(shè)在只有均勻化少群截面與擴(kuò)散系數(shù)作為等效均勻化常數(shù)的情況下是不成立的。SPH方法通過調(diào)整均勻化少群截面使得式(1)成立。

(2)

在保證反應(yīng)率守恒的條件下，有：

(3)

整理可得：

(4)

聯(lián)立式(2)、(4)可得SPH因子的計(jì)算公式：

(5)

式(5)中非均勻柵元平均中子通量密度由柵格高階輸運(yùn)計(jì)算直接獲得；均勻柵元平均中子通量密度的求解須與下游堆芯計(jì)算中所使用的中子學(xué)求解方法保持一致。SPH因子通過均勻柵元平均中子通量密度計(jì)算得到，而獲得中子通量密度計(jì)算所需的均勻化少群常數(shù)與擴(kuò)散系數(shù)又需通過SPH因子的修正，因此SPH因子的計(jì)算是一反復(fù)迭代的計(jì)算過程，具體計(jì)算步驟如下。

1) 在單組件全反射邊界條件下進(jìn)行高階輸運(yùn)計(jì)算，產(chǎn)生非均勻柵元的平均中子通量密度、各柵元的少群均勻化截面和擴(kuò)散系數(shù)。

2) 初始化SPH因子。

3) 利用SPH因子修正步驟1中產(chǎn)生的少群均勻化截面與擴(kuò)散系數(shù)，并用下游堆芯計(jì)算中所使用的中子學(xué)求解方法對(duì)全反射邊界條件下柵元均勻化的組件問題進(jìn)行計(jì)算，產(chǎn)生均勻化柵元的平均中子通量密度。

4) 根據(jù)式(5)，利用步驟1、3中的中子通量密度計(jì)算結(jié)果更新SPH因子。在更新SPH因子前，需對(duì)步驟3中產(chǎn)生的平均中子通量密度進(jìn)行歸一化處理，使得整個(gè)組件每一能群下的平均中子通量密度與步驟1中的中子通量密度保持一致，歸一化公式如下：

(6)

5) 判斷SPH因子是否收斂，收斂條件如式(7)所示，若滿足上述收斂條件，則計(jì)算停止；否則，進(jìn)入下一次迭代計(jì)算，返回步驟3。

(7)

式中，n為迭代次數(shù)。

1.2 與空間泄漏相關(guān)的SPH方法

在柵格計(jì)算中，對(duì)單燃料組件采用全反射邊界條件進(jìn)行計(jì)算，整個(gè)問題在均勻化前后不存在中子泄漏，即利用SPH因子保證反應(yīng)率守恒的同時(shí)還保證了中子泄漏率的守恒，燃料組件特征值也因此而守恒。當(dāng)SPH方法應(yīng)用于真空邊界，即反射層參數(shù)計(jì)算時(shí)，SPH因子迭代計(jì)算的收斂性無法保證。這是因?yàn)樵谡婵者吔缦氯肷渲凶恿鳛?的邊界條件無法保證出射中子流密度的守恒，迭代過程中改變均勻化少群截面的同時(shí)，影響著整個(gè)問題的中子泄漏率及其特征值。因此，此類邊界條件存在時(shí)，SPH因子迭代計(jì)算的收斂性無法保證，且在收斂的情況下依舊無法保證整個(gè)問題中子泄漏率和特征值的守恒。

針對(duì)此類邊界下SPH因子的迭代計(jì)算，與空間泄漏相關(guān)的SPH方法能使迭代計(jì)算很好地收斂并保證各柵元中子泄漏率的守恒，其核心思想是以柵元級(jí)少群反照率邊界條件保證均勻化前后各柵元的中子泄漏率守恒。與空間泄漏相關(guān)的SPH方法將柵元級(jí)少群反照率的計(jì)算更新加入到SPH因子的迭代求解過程中，具體計(jì)算流程如下。

1) 真空邊界條件下柵格高階輸運(yùn)計(jì)算，產(chǎn)生非均勻柵元的平均中子通量密度、各柵元的少群均勻化截面和擴(kuò)散系數(shù)、組件在真空邊界上各柵元的入射和出射中子流密度和少群反照率。反照率計(jì)算公式如下：

(8)

式中：βi,g為第i個(gè)柵元第g群反照率；h為非均勻輸運(yùn)計(jì)算的細(xì)群編號(hào)；J為中子流密度；角標(biāo)in、out分別表示入射、出射。

2) 初始化SPH因子并以式(8)計(jì)算的柵元級(jí)少群反照率作為Pin-by-pin均勻化組件的首次邊界條件。

3) 利用SPH因子修正少群均勻化截面與擴(kuò)散系數(shù)，并用下游堆芯計(jì)算中所使用的中子學(xué)求解方法計(jì)算柵元級(jí)少群反照率邊界下的柵元均勻化組件問題，產(chǎn)生均勻化柵元的平均中子通量密度及邊界上各柵元的入射和出射中子流密度。

4) 按照式(6)歸一化處理平均中子通量密度，根據(jù)式(5)更新SPH因子，迭代過程中真空邊界上的柵元級(jí)少群反照率的更新如式(9)、(10)所示。

(9)

(10)

5) 判斷SPH因子是否收斂，收斂條件為式(7)。若滿足上述收斂條件，則計(jì)算停止；否則，進(jìn)入下一次迭代計(jì)算，返回步驟3。

與空間泄漏相關(guān)的SPH方法通過SPH因子保證各能群各柵元的反應(yīng)率守恒的同時(shí)，通過柵元級(jí)少群反照率邊界保證各群中子泄漏率的守恒，燃料組件特征值也因此而守恒，計(jì)算流程如圖1所示。

圖1 與空間泄漏相關(guān)的SPH方法計(jì)算流程圖Fig.1 Calculation flow of spatial neutron-leakage dependent SPH method

2 數(shù)值驗(yàn)證

為驗(yàn)證SPH方法在壓水堆堆芯Pin-by-pin計(jì)算中的計(jì)算精度，以KAIST基準(zhǔn)題2A問題下的二維全堆芯問題進(jìn)行驗(yàn)證與分析，堆芯布置和幾何結(jié)構(gòu)如圖2所示。

圖2 KAIST基準(zhǔn)題二維全堆芯布置和幾何結(jié)構(gòu)Fig.2 Layout and geometry of KAIST benchmark 2D whole core problem

采用堆芯Pin-by-pin計(jì)算和堆芯組件均勻化粗網(wǎng)計(jì)算兩種方案進(jìn)行堆芯計(jì)算并加以分析和比較。

1) 依據(jù)SPH方法及改進(jìn)型SPH方法的理論模型，針對(duì)擴(kuò)散方法(P1)和簡(jiǎn)化球諧函數(shù)方法(SP3)利用SPH因子分別產(chǎn)生燃料組件及反射層組件的兩群和七群的等效均勻化少群常數(shù)，并在指數(shù)函數(shù)展開節(jié)塊法中進(jìn)行全堆芯Pin-by-pin計(jì)算。

2) 根據(jù)傳統(tǒng)組件均勻化計(jì)算方法，產(chǎn)生兩群帶組件不連續(xù)因子(ADF)的組件等效均勻化常數(shù)并采用先進(jìn)變分節(jié)塊方法進(jìn)行堆芯計(jì)算。變分節(jié)塊方法在空間處理上直接針對(duì)中子通量密度分布進(jìn)行離散，因此不需要功率重構(gòu)過程而能直接得到中子通量密度分布，通過形狀因子的求解即可得到堆芯棒功率分布。

二維堆芯參考解由特征線一步法計(jì)算獲得，特征值參考解為0.979 48，棒功率分布如圖3所示。堆芯兩步法計(jì)算的特征值誤差和棒功率誤差列于表1。由表1可知：堆芯Pin-by-pin簡(jiǎn)化球諧函數(shù)的計(jì)算精度普遍高于組件均勻化的堆芯粗網(wǎng)節(jié)塊的計(jì)算精度；除因誤差抵消導(dǎo)致計(jì)算精度提高的兩群擴(kuò)散計(jì)算外，堆芯Pin-by-pin擴(kuò)散計(jì)算的計(jì)算精度均低于組件均勻化的堆芯粗網(wǎng)節(jié)塊計(jì)算。因此，簡(jiǎn)化球諧函數(shù)方法是堆芯Pin-by-pin計(jì)算的必要條件；隨著能群數(shù)的增加，堆芯Pin-by-pin簡(jiǎn)化球諧函數(shù)的計(jì)算精度顯著提高。七群簡(jiǎn)化球諧函數(shù)計(jì)算中特征值誤差小于50 pcm，棒功率相對(duì)誤差的最大值(MAX)低于4%，棒功率相對(duì)誤差的平均值(AVG)、均方根值(RMS)、平均相對(duì)值(MRE)都在1%以下；而在組件均勻化的粗網(wǎng)節(jié)塊計(jì)算中，特征值誤差為219 pcm，最大棒功率相對(duì)誤差接近10%。因此，七群是適用于堆芯Pin-by-pin計(jì)算精度的能群數(shù)目。

組件均勻化粗網(wǎng)節(jié)塊計(jì)算的棒功率誤差分布如圖4所示，堆芯Pin-by-pin擴(kuò)散計(jì)算和簡(jiǎn)化球諧函數(shù)計(jì)算的棒功率相對(duì)誤差分布分別如圖5、6所示。由圖5、6可發(fā)現(xiàn)，最大棒功率相對(duì)誤差出現(xiàn)在組件與組件交界面以及燃料柵元與反射層交界面附近。這是堆芯內(nèi)組件真實(shí)環(huán)境與單組件計(jì)算中全反射邊界條件之間的不同所致。

圖3 KAIST基準(zhǔn)題二維全堆芯棒功率分布Fig.3 Pin-power distribution of KAIST benchmark 2D whole core problem

表1 堆芯兩步法計(jì)算的特征值誤差和棒功率相對(duì)誤差Table 1 Two-step calculation results of kinf and pin-power relative error

圖4 組件均勻化粗網(wǎng)節(jié)塊計(jì)算的棒功率相對(duì)誤差分布Fig.4 Pin-power relative error distribution of assembly-homogenized nodal calculation

3 結(jié)論

本文在壓水堆堆芯Pin-by-pin計(jì)算中采用SPH方法作為均勻化技術(shù)，對(duì)燃料組件計(jì)算了SPH因子，產(chǎn)生了Pin-by-pin等效均勻化參數(shù)；針對(duì)存在中子泄漏現(xiàn)象的反射層組件，采用與空間泄漏相關(guān)的SPH方法產(chǎn)生了反射層組件的等效均勻化參數(shù)?；贙AIST基準(zhǔn)題，分析了壓水堆堆芯Pin-by-pin計(jì)算中應(yīng)用SPH因子的堆芯計(jì)算精度，與傳統(tǒng)組件均勻化計(jì)算方法相比，應(yīng)用SPH方法的壓水堆堆芯Pin-by-pin的計(jì)算精度更高?；谄呷旱暮?jiǎn)化球諧函數(shù)方法是高精度堆芯Pin-by-pin計(jì)算方案的實(shí)施策略，與KAIST基準(zhǔn)題一步法計(jì)算結(jié)果相比，其特征值誤差小于50 pcm，最大棒功率相對(duì)誤差小于4%，棒功率相對(duì)誤差的AVG、MRS、MRE都在1%以下。

圖5 堆芯Pin-by-pin擴(kuò)散計(jì)算的棒功率相對(duì)誤差分布Fig.5 Pin-power relative error distribution of Pin-by-pin diffusion calculation

圖6 堆芯Pin-by-pin SP3計(jì)算的棒功率相對(duì)誤差分布Fig.6 Pin-power relative error distribution of Pin-by-pin SP3 calculation