劉仕倡，王 侃，陳義學

(1.華北電力大學 核科學與工程學院，北京 102206；2.清華大學 工程物理系，北京 100084)

蒙特卡羅方法可根據原子核的反應截面來隨機模擬中子與不同核素相互作用時的反應概率。材料溫度將影響靶核熱運動，從而對反應截面產生影響，這就是所謂的多普勒效應。在核反應堆的高保真模擬中，燃料和慢化劑有非常精細的溫度變化和分布。因此，對于不同溫度下的截面必須采用特殊處理。

傳統方法是預先產生不同溫度點的截面用于插值。基于線性插值點數據卷積的Sigma1方法是預處理工具中最常用的方法，如NJOY[1]。然而，當溫度點數較多時，預生成的方法會帶來內存占用方面的問題。近年來，為減少內存的占用，在多物理耦合計算中提出了幾種在線溫度截面處理技術。MCNP6[2]提出了基于級數展開的溫度擬合方法。Becker等[3]提出了利用隨機抽樣算法計算有效多普勒展寬截面的方法。為提高隨機抽樣的效率，Serpent程序還提出了另一種拒絕抽樣算法(TMS方法)[4]。在OpenMC[5]中使用了另一種稱為多極點(multipole)表示的方法。Dean等[6]提出了應用Gauss-Hermite積分來計算多普勒展寬方程的方法。

然而，現有的方法均有其局限性。例如，MCNP6中的溫度擬合方法依賴于預先生成的系數。OpenMC中的多極點表示方法依賴于預先生成的多極點數據。拒絕抽樣方法不能使用徑跡長度估計器，這將限制該方法的應用。其他方法如Sigma1法和傳統的Gauss-Hermite法存在效率低的問題。因此，有必要提出一種新的在線多普勒展寬方法，該方法不需進行溫度擬合和多極點數據生成等額外的處理，從而提高計算精度和效率。

本文提出一種新的Gauss-Hermite求積方法，并在蒙特卡羅程序RMC[7]中實現，以提高傳統Gauss-Hermite求積方法的效率。

1 基于積分的多普勒展寬方法

σ(|v-v′|,0)P(v′,T)dv′

(1)

1.1 Sigma1方法

Sigma1方法是利用點截面的線性化特性進行卷積，將式(1)轉化為式(2)：

[e-α(v-V)2-e-α(v+V)2]dV

(2)

α=M/2kT

(3)

式中，v和V分別為中子速度和中子相對于靶核的速度。為簡化式(2)，定義了如下變量：

x2=αV2

(4)

y2=αv2

(5)

將式(2)轉化為：

(6)

再定義：

(7)

式(6)可拆分為兩部分：

式(6)最終可轉化為一系列Hn的疊加，該函數定義如下：

(8)

式中，n為次數。

在實際應用中，式(8)的求解無論采用誤差余函數還是Taylor展開均非常耗時。

1.2 Gauss-Hermite正交基函數法

Gauss-Hermite正交基函數定義為：

(9)

為使式(7)形式類似式(9)的形式，定義變量z=x-y，將式(7)寫為：

(z+y)2e-z2dz

(10)

式中：z為求積組的取值點；n為求積組取值點的數目，如圖1所示；y代表中子速度，x=z+y代表相對速度。對于式(10)，注意到指數項的存在，z的取值范圍為-4≤z≤4即可。需注意的是，式(10)的積分下界為負無窮。因此當y>4時，式(10)可用Gauss-Hermite多項式表示，而y<4時則不能。因此，在應用Gauss-Hermite方法時，必須考慮適用范圍為y>4。y<4時將采用高斯-勒讓德求積組或Sigma1算法。在Romano的實現[8]中也采用了相同的處理方法。在Romano的實現中，使用了16點正交求積組。

圖1 Gauss-Hermite多項式求積組Fig.1 Diagram of Gauss-Hermite quadrature

然而，退回到高斯-勒讓德求積組或Sigma1算法將降低Gauss-Hermite方法的效率。此外，Gauss-Hermite方法在能量網格中搜索相對能量x，以及插值相應的截面時效率較低。

2 改進Gauss-Hermite方法

針對傳統Gauss-Hermite方法存在的問題，提出了改進Gauss-Hermite方法。解決這些問題的一些關鍵技術如下。

2.1 低能段的修正

眾所周知，低能段截面是相對比較光滑的，如圖2所示。一般地，低能散射截面變化很小，接近1個常數；而低能吸收、裂變截面服從1/v率。因此，由式(1)可知，經過多普勒展寬后，吸收截面保持不變，而散射截面應乘以修正項，如式(11)所示。

圖2 900 K時238U總截面Fig.2 Total cross section of 238U at 900 K

(11)

erf是誤差函數，并且：

(12)

式中，n為核素。

通過上述處理，對于低能光滑段，利用式(11)修正散射截面時，不使用Gauss-Hermite方法，這可減少Gauss-Hermite方法的使用頻率，代之以計算成本較低的修正項。從而避免了Gauss-Hermite方法在低能段y<4不適用的問題。

低能區(qū)與共振能區(qū)的分界點很重要，必須采用合理的方法將低能區(qū)與共振能區(qū)劃分開。

2.2 共振段的開始和結束能量

如上所述，低能區(qū)與共振能區(qū)的分界點可稱為共振的起始能量或低能光滑段的結束點。低能光滑段的結束點可根據每個核素所在的最大溫度下多普勒展寬的作用范圍來確定，從而得到低能光滑段的結束點。對于某個中子能量點，根據中子與靶核之間的最大相對能量和最小相對能量，可找到上界和下界。然后計算出上、下界的差值，將低能光滑段結束點作為共振的起始能量，如圖3所示。

圖3 低能光滑段的結束點Fig.3 End point of low energy smooth segmen

利用上述算法，圖3中低能光滑段結束點為1.25 eV。而從圖2可看出在4.41 eV附近有1個很小的共振峰?？梢?，提出的上界-下界算法可有效地找出低能光滑段結束點。

共振段結束點的查找方法與NJOY的處理一致，即可分辨共振區(qū)的上界、閾能反應對應的最小能量以及1 MeV三者的最小能量。

2.3 多普勒展寬哈希表

另一個重要的改進是多普勒展寬哈希表的使用。在Gauss-Hermite方法中，尋找相對能量的能量網格序號非常耗時。RMC中采用了哈希表方法，來加速能量網格查找的速度[9]。例如，當前中子能量為E，哈希表先根據E找到對應的能量分箱i，然后再在能量分箱i中使用二分法查找能量E對應的網格序號。

當哈希表確定了能量分箱i后，可得到分箱i的能量上、下界Eupper和Elower，而根據Eupper和靶核速度可得到最大相對能量Emax，根據Elower和靶核速度可得到最小相對能量Emin。當中子能量在哈希表的分箱i內時，多普勒展寬后的相對能量將位于多普勒展寬哈希表的分箱i內，如圖4所示。

圖4 多普勒展寬哈希表Fig.4 Diagram of Doppler broadening Hash map

3 結果與分析

采用典型17×17壓水堆組件模型，如圖5所示，組件模型參數列于表1。對多普勒展寬的精度和效率進行了測試。在每個柵元中均有燃料、氣隙(16O)、輕水(1H、16O)。燃料溫度為900 K，其他材料溫度為600 K。蒙特卡羅計算中，每代使用100 000個中子，共1 000代，其中300個非活躍代。選取了采用基于ENDF/B-Ⅶ.0的精確溫度庫的RMC結果作為參考。

圖5 典型壓水堆組件模型Fig.5 Typical PWR assembly model

參數數值燃料密度,g/cm310.196 燃料富集度,%3燃料芯塊半徑,cm0.409 6燃料包殼外徑,cm0.475 0柵距,cm1.26冷卻劑密度,g/cm31.0

改進Gauss-Hermite方法采用300 K的ACE庫作為基礎數據庫。在輸運及燃耗計算中，均采用了哈希表加速，并對燃料棒中的中子能譜進行統計。

3.1 輸運計算

對于輸運計算，燃料中僅有235U、238U、16O。表2對kinf進行了比較?？煽闯鰝鹘yGauss-Hermite方法和改進Gauss-Hermite方法與精確庫的差別均在1倍標準偏差左右。傳統Gauss-Hermite方法的計算時間是參考值的9.3倍，而改進Gauss-Hermite方法時間僅增加了10%。

表2 壓水堆組件算例的kinf和計算時間Table 2 kinf and calculation time of PWR assembly case

圖6比較了3種方法的中子能譜，可看出傳統和改進Gauss-Hermite方法的能譜均與參考值很接近。大部分相對誤差在5%以內。較大的相對誤差僅出現在通量本身很小的能量點。

圖6 壓水堆組件算例的能譜相對誤差Fig.6 Relative error of energy sprectrum for PWR assembly case

3.2 輸運-燃耗耦合計算

第2個算例為燃耗后(58.065 MW·d/kg(HM))的PWR組件，燃料中共134種核素。此次計算還增加了燃料中核素的燃耗相關單群截面統計，因此相當于輸運-燃耗耦合計算進行到58.065 MW·d/kg(HM)的一步輸運計算。

從表3可看出，改進Gauss-Hermite方法的時間增加了29%，而kinf僅相差6.5 pcm。而傳統Gauss-Hermite方法的時間為精確庫的3.45倍。改進Gauss-Hermite方法在帶燃耗相關單群截面統計的輸運計算的效率較TMS方法(1.98倍)和OpenMC的Multipole方法(1.59倍)的高，因此更適合應用于輸運-燃耗耦合計算。

表3 帶單群截面統計的燃耗后PWR組件算例的kinf和計算時間Table 3 kinf and calculation time of burned PWR assembly case with one group cross section

4 結論

本文提出了一種改進Gauss-Hermite方法，并在蒙特卡羅程序RMC中進行了編程實現。壓水堆燃料組件算例的輸運和燃耗計算結果表明，改進Gauss-Hermite方法精度與傳統Gauss-Hermite方法相當，而且可大幅提高計算效率，特別是輸運-燃耗耦合計算的效率。該方法簡潔而直接，避免了額外的預處理過程，為在線多普勒展寬提供了一種新的有效方法。因此，該方法能適用于各種蒙特卡羅粒子輸運程序，具有良好的應用前景，為高保真多物理耦合計算奠定基礎。