鄭 征，黎 輝，丁謙學(xué)，上官丹驊

(1.上海核工程研究設(shè)計院，上海 200233；2.北京應(yīng)用物理與計算數(shù)學(xué)研究所，北京 100190)

JMCT-S一次屏蔽計算源項生成功能開發(fā)

鄭 征1，黎 輝1，丁謙學(xué)1，上官丹驊2

(1.上海核工程研究設(shè)計院，上海 200233；2.北京應(yīng)用物理與計算數(shù)學(xué)研究所，北京 100190)

一次屏蔽計算對于評估工作人員的輻射劑量、確保反應(yīng)堆壓力容器(RPV)及堆內(nèi)構(gòu)件在整個反應(yīng)堆壽期內(nèi)的安全性以及防止混凝土屏蔽體及其外部部件和結(jié)構(gòu)被過度活化具有重要的意義。對于一次屏蔽源項的處理，JMCT-S程序自帶的源粒子抽樣功能無法完全滿足其計算需求。本文開發(fā)了JMCT-S程序的源項生成程序和源抽樣子程序，并在秦山一期和CAP1400一次屏蔽計算模型上進行了驗證和應(yīng)用。數(shù)值結(jié)果表明，推導(dǎo)的理論模型和開發(fā)的程序是正確的，從而為后續(xù)提高一次屏蔽設(shè)計精度提供了基礎(chǔ)。

JMCT-S；源項生成程序；一次屏蔽

一次屏蔽由圍繞堆芯的堆內(nèi)構(gòu)件(主要包括圍板、吊籃等)、水隙[包括圍板與吊籃、吊籃與壓力容器(RPV)之間的水隙等]、RPV和圍繞著RPV的混凝土結(jié)構(gòu)等組成。一次屏蔽具有降低堆芯和結(jié)構(gòu)材料產(chǎn)生的放射性對工作人員的輻射、確保RPV及堆內(nèi)構(gòu)件在整個反應(yīng)堆壽期內(nèi)安全性、確?；炷疗帘误w性能的穩(wěn)定性和完整性，以及防止混凝土一次屏蔽外的部件和結(jié)構(gòu)被過度活化等功能[1]。

國內(nèi)外通常采用蒙特卡洛方法和離散縱標(biāo)方法進行一次屏蔽設(shè)計。JMCT-S程序[2]是計算復(fù)雜三維幾何結(jié)構(gòu)中的粒子輸運的大型多功能蒙特卡洛程序。JMCT-S自帶的源描述功能可以設(shè)定多種類型的源粒子位置、能量和方向，常用的源包括點源、面源、體源、各向同性源、各向異性源、單能源和有能量分布的源等等。由于反應(yīng)堆堆芯幾何結(jié)構(gòu)復(fù)雜，以及堆芯源項在空間分布(PIN-BY-PIN分布)和能量分布(不同燃耗水平下由于重核素的增加導(dǎo)致不同組件裂變份額產(chǎn)生變化，從而使得不同組件具有不同的裂變譜)上的復(fù)雜性，當(dāng)前JMCT-S程序自帶通用源描述無法實現(xiàn)精確模擬，因此有必要開展JMCT-S的源項生成程序和源抽樣子程序的研究。

為了進一步完善JMCT-S程序的源項生成功能，提高JMCT-S程序的應(yīng)用范圍，本文研究了JMCT-S程序的源項理論模型，開發(fā)了相應(yīng)的源項生成程序和源抽樣子程序，實現(xiàn)了考慮裂變份額的三維全堆芯PIN-BY-PIN源分布，并且在秦山一期計算模型上進行了驗證，在CAP1400和CAP1700一次屏蔽設(shè)計中進行了應(yīng)用。

1 理論模型

對于堆芯源粒子抽樣，首先根據(jù)三維功率分布、裂變核素份額、裂變譜、每次裂變產(chǎn)生粒子數(shù)和能量等計算出源粒子的各個參數(shù)的累積分布函數(shù)(CDF)，然后基于CDF進行抽樣，確定源粒子的類型、權(quán)重、坐標(biāo)、方向、能量和所在柵元等參數(shù)。其中方向上服從裂變源各向同性分布，因此不需要計算CDF，采用各向同性抽樣。

1.1 累積分布函數(shù)計算方法

第m個燃料組件的能譜根據(jù)該組件的裂變核素份額、每次裂變產(chǎn)生粒子數(shù)和裂變譜等參數(shù)計算得到[3]：

(1)

式中：g,g′——能群編號；

m——組件編號；

n——核素編號；

p——粒子類型編號；

N——核素數(shù)目；

χ′(m,g,p)——第p種粒子、第g群、第m個組件粒子的產(chǎn)生概率；

χ(n,g,p)——第p種粒子、第g群、第n種核素粒子的產(chǎn)生概率。源粒子的能譜采用BUGLE[4]的能群結(jié)構(gòu)；

f(m,n)——第n種核素、第m個組件的粒子裂變份額；

ν(n,p)——第p種粒子、第n種核素每次裂變釋放粒子數(shù)。

根據(jù)各個燃料組件的能譜，可以計算出各個組件的能量CDF：

(2)

式中：E(m,g,p)——第p種粒子、第g群、第m個組件的CDF；

G(p)——第p種粒子的能群數(shù)目，p=1表示中子，G(1)=47，p=2表示光子，G(2)=20。

三維粒子源分布計算公式如下：

(3)

式中：S(q,p)——第p種粒子、第q個網(wǎng)格源強，粒子/(cm3·s)；

P(q)——第q個網(wǎng)格真實功率密度，MW/cm3，P(q)=P(i,j,k)，(i,j,k)表示網(wǎng)格編號；

C——單位轉(zhuǎn)換因子，6.24×1012MeV/(s·W)；

K(m)——第m個組件每次裂變釋放的能量，MeV。其中f(m,n)和K(m)根據(jù)組件燃耗和初始富集度通過六種核素的裂變份額以及每次裂變釋放的總能量插值得到。

根據(jù)公式(3)計算所得三維PIN-BY-PIN源分布，可以計算得到相應(yīng)的三維空間各個網(wǎng)格的CDF：

(4)

式中：q,q′,q″——網(wǎng)格編號；

Q——總網(wǎng)格數(shù)目；

D(q″,p)——第p種粒子、第q″個網(wǎng)格的CDF；

S(q,p)——第p種粒子、第q個網(wǎng)格的源強，粒子/(cm3·s)；

V(q)——第q個網(wǎng)格的體積，cm3。

堆芯裂變產(chǎn)生第p種粒子數(shù)目計算如下：

N′(p)=∑i,j,kS(q,p)·V(q)

(5)

式中：N′(p)——堆芯裂變產(chǎn)生第p種粒子的數(shù)目，粒子/s。

經(jīng)過5年的艱辛努力，西達里亞油田“外練筋骨皮、內(nèi)練一口氣”，各項指標(biāo)均出現(xiàn)逆轉(zhuǎn)，并實現(xiàn)質(zhì)的飛躍——原油日產(chǎn)量由2017年1月初的11噸升至目前的65噸，2018年折耗比2017年降低67%，2018年噸油操作成本比2015年降低63%，前三季度實現(xiàn)利潤915萬元。

根據(jù)堆芯裂變產(chǎn)生總的粒子數(shù)目，可以計算出粒子類型的CDF供JMCT-S抽樣：

(6)

式中：p,p′,p″——粒子類型編號；

P——粒子類型數(shù)目，P=1表示中子，P=2表示中子和光子；

F(p)——第p種粒子的CDF。

1.2 源粒子抽樣方法

考慮到變量之間的依賴關(guān)系，應(yīng)按照CDF維度從低到高的順序進行抽樣。首先確定粒子類型，再確定粒子空間位置，最后確定粒子能量。方向采用各向同性抽樣。

對于粒子類型，選取(0,1]區(qū)間內(nèi)的隨機數(shù)ξp，若F(p-1)<ξp≤F(p)，則表示源粒子為第p種粒子。

對于空間抽樣，選取(0,1]區(qū)間內(nèi)的隨機數(shù)ξq，若D(q-1,p)<ξq≤D(q,p)，則源粒子在第q個網(wǎng)格內(nèi)。源粒子的空間位置坐標(biāo)(x,y,z)計算公式如下：

(7)

式中：x,y,z——源粒子的空間位置坐標(biāo)，cm；

ξq1,ξq2,ξq3——(0,1]區(qū)間內(nèi)的隨機數(shù)。

對于能量抽樣，選取(0,1]區(qū)間內(nèi)的隨機數(shù)ξg，若E(g-1,q,p)<ξg≤E(g,q,p)，則源粒子在第g群內(nèi)。根據(jù)下式計算出源粒子的初始能量：

(8)

式中：E′——源粒子的初始能量，MeV；

ξg1——(0,1]區(qū)間內(nèi)的隨機數(shù)。

2 驗證及應(yīng)用

基于上文介紹的源項理論模型，本文編寫了JMCT-S源項生成程序及源抽樣子程序，并且在秦山一期計算模型進行了驗證，在CAP1400和CAP1700一次屏蔽設(shè)計中進行了應(yīng)用。

2.1 秦山一期輻照監(jiān)督管中子注量率計算

基于參考文獻[5]建立的秦山一期計算模型如圖1所示，本文采用新開發(fā)的源項程序進行了計算。

圖1 秦山一期計算模型示意圖Fig.1 Geometry model of Qinshan No.1 reactor

表1給出了在額定功率下第六根輻照監(jiān)督管中子注量率計算結(jié)果。從表1中可以看出，本文計算結(jié)果和測量值吻合良好，堆芯中平面結(jié)果相對誤差在5%以內(nèi)，上焊縫結(jié)果相對誤差在20%以內(nèi)。偏差的主要來源包括數(shù)據(jù)庫、計算模型和程序帶來的偏差等。

表1 秦山一期輻照監(jiān)督管中子注量率Table 1 Capsule neutron fluence rate of Qinshan No.1 reactor

2.2 CAP1400壓力容器快中子注量率計算

基于堆芯幾何結(jié)構(gòu)等參數(shù)，本文建立了JMCT-S的CAP1400計算模型如圖2所示。圖3給出了壓力容器內(nèi)表面、方位角方向從0°到45°、軸向上從堆芯活性區(qū)中平面到堆芯活性區(qū)頂部區(qū)域的快中子(>1.0MeV)注量率分布。

圖2 CAP1400計算模型示意圖Fig.2 Geometry model of CAP1400

圖3 CAP1400壓力容器內(nèi)表面快中子注量率分布Fig.3 Fast neutron fluence rate distribution on RPV inner surface of CAP1400

基于離散縱標(biāo)方法的2D/1D注量率綜合法在工程上被廣泛用于計算RPV中平面高度的快中子注量率[7]，因此選擇DORT[8]程序的計算結(jié)果作為參考。其中DORT計算采用BUGLE-96數(shù)據(jù)庫，中子劃分47群，求積組離散階數(shù)為S8，網(wǎng)格劃分數(shù)目分別為414×202(二維Rθ計算)，401×231(二維RZ計算)和401(一維R計算)，源項生成模型來自SORCERY程序[9]。

表2給出了壓力容器內(nèi)表面最大快中子注量率。從表中可以看出，JMCT-S和DORT計算結(jié)果吻合一致，相對誤差在2%以內(nèi)。

表2 CAP1400壓力容器內(nèi)表面最大快中子注量率Table 2 Maximum fast neutron fluence rate on RPV inner surface of CAP1400

2.3 CAP1700壓力容器快中子注量率計算

基于已有的CAP1700計算模型，本文采用新開發(fā)的源項程序進行了計算。表3給出了壓力容器內(nèi)表面最大快中子注量率(>1MeV)計算結(jié)果。從表中可以看出，當(dāng)采用相同的235U裂變譜時，本文計算結(jié)果和參考結(jié)果吻合良好，相對誤差在2%左右。當(dāng)采用混合裂變譜時，本文計算結(jié)果比采用235U裂變譜時增大了大約12%。其中TORT計算采用BUGLE-96數(shù)據(jù)庫，中子劃分47群，光子劃分20群，求積組離散階數(shù)為S8，網(wǎng)格劃分數(shù)目分別為199×90×99(三維1/8幾何RθZ計算)，源項生成模型來自SORCERY程序[9]。

表3 CAP1700壓力容器內(nèi)表面最大快中子注量率Table 3 Maximum fast neutron fluence rate on RPV inner surface of CAP1700

混合裂變譜包括235U、238U、239Pu、240Pu、241Pu和242Pu 6種核素。與235U每次裂變釋放的能量相比，238U、239Pu和241Pu每次釋放的能量更多；與235U每次裂變釋放的中子數(shù)目相比，238U、239Pu和240Pu和241Pu每次裂變釋放的中子數(shù)目更多；與235U裂變譜相比，238U裂變譜更軟，239Pu、240Pu和241Pu的裂變譜更硬[3]。這些原因?qū)е虏捎没旌狭炎冏V時計算結(jié)果較采用235U裂變譜時偏大。因此，本文的計算結(jié)果是合理的。

3 結(jié)論

本文推導(dǎo)了能量、空間和粒子類型等參數(shù)的累積分布函數(shù)計算公式，給出了源粒子能量、空間和粒子類型等參數(shù)的抽樣方法，編寫了相應(yīng)的程序模塊，并在秦山一期和CAP1400計算模型上進行了初步應(yīng)用驗證。數(shù)值結(jié)果表明，本文推導(dǎo)的理論模型和開發(fā)的程序是正確的，從而為后續(xù)提高一次屏蔽設(shè)計精度提供了基礎(chǔ)。

[1] 傅守信，胡建軍.秦山核電二期工程反應(yīng)堆主屏蔽設(shè)計[J]. 核動力工程，2003，24(S1)：64-67.

[2] JMCT-S 用戶使用手冊[R]. 北京應(yīng)用物理與計算數(shù)學(xué)研究所，2015.

[3] A. Vasiliev，H. Ferroukhi，M. A. Zimmermann，et al.，Development of a CASMO-4/SIMULATE-3/MCNPX Calculation Scheme for PWR Fast Neutron Fluence Analysis and Validation against RPV Scraping Test Data[J]. Annals of Nuclear Energy，2007，34：615-627.

[4] BUGLE-96，Coupled 47 Neutron，20 Gamma-Ray Group Cross Section Library Derived from ENDF/B-VI for LWR Shielding and Pressure Vessel Dosimetry Applications[R]. Oak Ridge National Laboratory，1999.

[5] 夏春梅，第六根輻照監(jiān)督管中子注量計算報告[R]. 上海核工程研究設(shè)計院，2014.

[6] 郭玥，秦山320MWe機組反應(yīng)堆壓力容器第六根輻照監(jiān)督管材料樣品中子能譜和快中子注量測定報告[R]. 中國原子能科學(xué)研究院反應(yīng)堆工程研究設(shè)計所，2014.

[7] J. A. Kulesza，Comparison of Three-Dimensional Flux Synthesis and Full Three-Dimensional Discrete Ordinates Methods for the Calculation of Reactor Cavity Bioshield Heat Generation Rates[J]. Radiation Transport and Protection，2011，175：228-237.

[8] W. A. Rhoades，R. L. Childs，The DORT Two-Dimensional Discrete Ordinates Transport Code[J]. Nuclear Science and Engineering，1988，99：88-89.

[9] SORCERY User Manual[R]. Westinghouse Electric Company，2006.

[10] W. A. Rhoades，D. B. Simpson，The TORT Three-Di-mensional Discrete Ordinates Neutron/Photon Transport Code[R]. Oak Ridge National Laboratory，1997.

DevelopmentonPrimaryShieldingCalculationSourceGenerationFunctionforJMCT-S

ZHENGZheng1，LIHui1，DINGQian-xue1，SHANGGUANDan-hua2

(1.Shanghai Nuclear Engineering Research and Design Institute，Shanghai，200233，China；2. Institute of Applied Physics and Computational Mathematics，Beijing，100190)

Primary shielding calculation plays an important role in evaluating the radiation dose of personnel，insuring the safety of reactor pressure vessel(RPV)and reactor internals during the whole lifetime，and preventing the concrete shielding as well as the outer components and structures from an excess of activation. The source particle sample function of the JMCT-S code cannot fully satisfy the calculation requirements for generating the source term of primary shielding calculation. In this paper，a source term generation code and a source sample code for JMCT-S are developed. Validation and application on primary shielding calculation model on Qinshan No.1 nuclear power plant and CAP1400 reactor are performed. Numerical results show that the theoretical model and the codes are both correct，and can be used to improve the accuracy of primary shielding calculation.

JMCT-S；source generation code；primary shielding

2017-03-20

鄭 征(1985—)，男，陜西人，工程師，博士，現(xiàn)從事核電廠屏蔽設(shè)計工作

TL328

A

0258-0918(2017)05-0805-05