呂煒楓 陳明亮 黃倩倩 宋 譽(yù) 洪 浩

（中廣核工程有限公司核電安全監(jiān)控技術(shù)與裝備國家重點(diǎn)實(shí)驗(yàn)室 深圳 518172）

在壓水堆核電站功率運(yùn)行期間，來自堆芯的中子將會對堆外設(shè)備和墻體材料進(jìn)行活化，從而產(chǎn)生活化產(chǎn)物源項(xiàng)，進(jìn)而對停堆期間和退役期間核電站內(nèi)工作人員產(chǎn)生輻射影響，而且隨著核電站運(yùn)行年限的增長，中子活化產(chǎn)生的長壽命核素（如60Co）的活度將逐漸增加，因此有必要建立堆外設(shè)備和墻體材料內(nèi)活化產(chǎn)物源項(xiàng)的精細(xì)評估手段，從而為后端的屏蔽設(shè)計(jì)和廢物處理提供輸入。

美國橡樹嶺實(shí)驗(yàn)室基于矩陣指數(shù)法開發(fā)了求解燃耗方程的燃料組件中裂變產(chǎn)物核素和活化產(chǎn)物核素源項(xiàng)的計(jì)算程序 ORIGEN-S［1］。ORIGEN-S 程序基于點(diǎn)燃耗進(jìn)行計(jì)算，其已廣泛應(yīng)用于國內(nèi)外核電站設(shè)計(jì)。

歐洲原子能共同體和英國原子能管理局共同開發(fā)了歐洲壓水堆活化產(chǎn)物源項(xiàng)計(jì)算程序FISPACT-2007［2］。FISPACT-2007基于歐拉指數(shù)法開發(fā)，可計(jì)算高能中子裝置的材料活化以及壓水堆和聚變堆中材料的活化情況。其計(jì)算配套的核數(shù)據(jù)庫為EAF-2007。

華北電力大學(xué)韓文靜等［3-5］分別基于歐拉指數(shù)法和切比雪夫有理近似法（Chebyshev Rational Approximation Method，CRAM）對燃耗方程進(jìn)行求解，采用EAF-2007數(shù)據(jù)庫開發(fā)了活化產(chǎn)物源項(xiàng)計(jì)算程序EuACT和ABURN。其中，EuACT程序?yàn)榛跉W拉指數(shù)法開發(fā)，ABURN程序?yàn)榛贑RAM方法開發(fā)。其基于ABURN程序計(jì)算了聚變堆第一壁活化例題和UO2燃料燃耗例題，并將ABURN程序的計(jì)算結(jié)果與歐洲活化產(chǎn)物源項(xiàng)計(jì)算程序FISPACT-2007進(jìn)行對比。結(jié)果表明：ABURN程序可達(dá)到FISPACT-2007程序同等精度，并且由于采用了CRAM方法，程序在燃耗步設(shè)置方面具有高度的靈活性，初步驗(yàn)證了ABURN程序的可用性與準(zhǔn)確性。

中國核動力研究設(shè)計(jì)院胡建軍等［6］采用二維離散縱標(biāo)輸運(yùn)程序DORT及ENDF/B-VI庫計(jì)算壓水堆壓力容器內(nèi)中子注量率分布，用自行研制的活化產(chǎn)物源項(xiàng)計(jì)算程序計(jì)算冷卻劑16N和17N源項(xiàng)，并驗(yàn)證堆芯徑向與軸向功率分布對計(jì)算結(jié)果的影響。其自行研制的活化源項(xiàng)計(jì)算程序名稱為SCACT。

中國科學(xué)院核能安全技術(shù)研究所張彬航等［7］基于超級蒙特卡羅核模擬軟件系統(tǒng)SuperMC，采用CRAM方法發(fā)展了基于深度優(yōu)先搜索的活化鏈動態(tài)構(gòu)建方法和大規(guī)模矩陣自適應(yīng)降階方法，進(jìn)行了活化產(chǎn)物源項(xiàng)計(jì)算的初步研究。通過選用IAEA-ACB（InternationalAtomic Energy Agency-Activation Calculation Benchmark）國際基準(zhǔn)例題及壓水堆燃料包殼材料例題，初步驗(yàn)證了該活化產(chǎn)物源項(xiàng)計(jì)算方法的正確性，且測試結(jié)果表明發(fā)展的大規(guī)模矩陣自適應(yīng)降階方法能夠有效地提高活化系數(shù)矩陣的求解效率。

綜上所述，國內(nèi)外各研究機(jī)構(gòu)針對壓水堆核電站活化產(chǎn)物源項(xiàng)計(jì)算程序進(jìn)行了大量的研究工作，但因程序與工程設(shè)計(jì)脫節(jié)以及軟件認(rèn)證的問題，除了國外開發(fā)的相關(guān)程序，尚未有自主研發(fā)的活化產(chǎn)物源項(xiàng)計(jì)算程序在工程設(shè)計(jì)中大規(guī)模應(yīng)用。在華龍一號“走出去”的大背景下，有必要研發(fā)可應(yīng)用于工程實(shí)際的具有自主知識產(chǎn)權(quán)的壓水堆活化產(chǎn)物源項(xiàng)計(jì)算程序。

1 核電站內(nèi)中子活化產(chǎn)物源項(xiàng)計(jì)算模型

1.1 計(jì)算方程組的建立

核電站內(nèi)受中子輻照的非裂變材料內(nèi)放射性核素產(chǎn)生途徑有以下兩種：先驅(qū)核素衰變產(chǎn)生和先驅(qū)核素中子俘獲產(chǎn)生。核電站內(nèi)受中子輻照的非裂變材料內(nèi)放射性核素消失途徑有以下兩種：放射性衰變消失以及與中子反應(yīng)消失。可分別建立放射性核素i的產(chǎn)生率和消失率計(jì)算公式如下：

式中：Si(t)為t時(shí)刻放射性核素i的產(chǎn)生率，s-1；Pi(t)為t時(shí)刻放射性核素i的消失率，s-1；Φk為k群中子注量率，m-2·s-1；Ni(t)為t時(shí)刻中子輻照材料內(nèi)放射性核素i的積存量；V為中子輻照材料體積，m3；nj為先驅(qū)核素 j的核子密度，m-3；σɑ，i，k為 k群中子對先驅(qū)核素i的微觀吸收截面，m2；σx，j，k為k群中子對先驅(qū)核素j的微觀x反應(yīng)截面，x反應(yīng)指（n，r）、（n，2n）、（n，p）、（n，ɑ）、（n，t）等中子反應(yīng)，m2；λi為 i核素的衰變常數(shù)，s-1；Ki，j為衰變分支比，無量綱；Γx，i，j為先驅(qū)核素 j通過與中子的x反應(yīng)產(chǎn)生核素i的產(chǎn)額，無量綱。

根據(jù)以上放射性核素產(chǎn)生率和消失率的計(jì)算公式，可建立放射性核素i積存量的微分方程如下：

求解以上微分方程可得中子輻照材料內(nèi)放射性核素i的積存量。如需計(jì)算n種核素，則以上方程可寫成矩陣的形式：

其中：

在以上方程中，非放射性核素也在n種核素之內(nèi)，由此先驅(qū)核素中子俘獲產(chǎn)生項(xiàng)可包含放射性核素和非放射性核素中子俘獲活化的貢獻(xiàn)。

對于式（3）中方程，其解為：

式中：X(t)是核素積存量矩陣；A是包含衰變和中子反應(yīng)的轉(zhuǎn)換矩陣（活化矩陣）。X(0)是初始核素積存量矩陣。

1.2 矩陣指數(shù)法求解

對式（8）中指數(shù)函數(shù)部分進(jìn)行Taylor級數(shù)展開，有：

式中：I表示單位矩陣。矩陣指數(shù)法在理論上可以精確求解燃耗問題，但是在具體的例子中，求矩陣冪的精確值時(shí)會產(chǎn)生各種困難。在這種情形下，就需要借助Bateman方程或者使用Gauss-Seidel迭代法來解決。

燃耗矩陣A的規(guī)模通常在2000×2000量級，但是在實(shí)際的活化源項(xiàng)計(jì)算中，往往僅關(guān)注于幾個(gè)主要的核素，極少出現(xiàn)一次性計(jì)算2000×2000的情況，考慮到核電廠實(shí)際情況，這個(gè)規(guī)模將大大縮小。

1.3 CRAM方法求解

CRAM方法的核心是在[0，∞)上對e-x進(jìn)行有理逼近。定義πk，l表示所有如下形式實(shí)有理函數(shù)rk，l(x)的集合：

式中：pj是不高于j階的多項(xiàng)式。根據(jù)逼近理論可知存在唯一的 r^k，l∈ πk，l滿足：

對于給定的k和l，可通過Remez算法或Carathéodory-Fejér方法得到逼近表達(dá)式。由式（10）可知，有理函數(shù)r^k，k的極點(diǎn){θ1，θ2，...，θk}是分立極點(diǎn)，因此可以通過部分分式展開得到：

式中：α0是函數(shù)在無窮遠(yuǎn)的極限；αi表示極點(diǎn)θi處的留數(shù)。式（12）也可以通過留數(shù)定理得到。另外，注意到αi和θi共軛出現(xiàn)，因此，式（12）可以簡化成：

對于CRAM方法計(jì)算矩陣指數(shù)，主要的難點(diǎn)在于對給定的k值求出相應(yīng)的有理函數(shù)的系數(shù)。原則上來說，p^k和q^k的系數(shù)都可以通過Remez型方法得到，但是計(jì)算過程需要相當(dāng)有技巧的算法和非常高精度的計(jì)算過程。文獻(xiàn)調(diào)研顯示14階系數(shù)的計(jì)算結(jié)果已經(jīng)具有較高的精度［5，8］。基于此，CPAP程序采用16階系數(shù)以進(jìn)一步提高精度，采用的CRAM系數(shù)見表1。

2 活化產(chǎn)物源項(xiàng)計(jì)算程序CPAP的開發(fā)

基于以上中子活化產(chǎn)物源項(xiàng)計(jì)算模型，采用面向?qū)ο蟮木幊趟季S以及當(dāng)前主流的圖形程序界面，開發(fā)了具有良好人機(jī)界面的各堆型壓水堆核電站通用的中子活化源項(xiàng)計(jì)算程序CPAP。其中，程序圖形化界面和計(jì)算模塊均采用C++語言開發(fā)。CPAP程序的計(jì)算流程圖見圖1和圖2。

CPAP程序有如下特點(diǎn)：

1）各堆型壓水堆核電站在計(jì)算中子活化產(chǎn)物源項(xiàng)時(shí)，一般采用矩陣指數(shù)法或者CRAM方法，由此CPAP程序中提供了矩陣指數(shù)法和CRAM方法兩種選項(xiàng)；

2）CPAP程序中采用的CRAM方法基于16階系數(shù)，計(jì)算精度進(jìn)一步提高；

3）CPAP程序采用的中子反應(yīng)多群截面數(shù)據(jù)庫基于ENDF/B-VII.0中的中子反應(yīng)子庫（NSUB=10）及EAF-2010中子活化反應(yīng)數(shù)據(jù)庫加工得出，涉及860種核素；

4）CPAP程序的截面數(shù)據(jù)庫采用47群結(jié)構(gòu)，能群劃分參照BUGLE-96多群數(shù)據(jù)庫的能群邊界；加工多群截面數(shù)據(jù)庫時(shí)，采用傳統(tǒng)核截面加工處理軟件NJOY中的相關(guān)運(yùn)行模塊，最終制成多群數(shù)據(jù)文件時(shí)采用壓水堆堆芯外的中子能譜作為權(quán)重譜；截面數(shù)據(jù)庫涉及的反應(yīng)類型共11種，包括總反應(yīng)、彈性散射反應(yīng)、非彈性散射反應(yīng)、輻射俘獲反應(yīng)、（n，2n）、（n，3n）、（n，α）、（n，p）、（n，t）、（n，2p）和裂變反應(yīng)；

表116 階CRAM系數(shù)Table 1 16th order coefficient of CRAM method

圖1 CPAP程序計(jì)算流程圖1Fig.1 Calculation flowchart 1 of CPAP code

5）具有良好的人機(jī)界面，根據(jù)用戶需求預(yù)先嵌入大量常見材料數(shù)據(jù)庫。

CPAP程序界面見圖3。

3 活化產(chǎn)物源項(xiàng)計(jì)算程序CPAP的測試驗(yàn)證

對于活化產(chǎn)物源項(xiàng)計(jì)算程序CPAP的測試包括算法的測試以及數(shù)據(jù)庫的測試兩部分。以下分別說明。

3.1 算法的測試

選取了部分國際上通用的基準(zhǔn)題以及自主建立的典型例題對CPAP程序的兩種算法進(jìn)行測試。算法測試時(shí)CPAP程序所使用的核素豐度、反應(yīng)類型、截面數(shù)據(jù)和衰變常數(shù)均取自EAF-2007數(shù)據(jù)庫，與FISPACT-2007程序一致。本文以兩個(gè)典型的活化算例說明程序的算法測試結(jié)果。

3.1.1 AISI304活化計(jì)算測試

在本測試算例中，AISI304的初始質(zhì)量為1 kg，其密度為8.04 g·cm-3。元素質(zhì)量成分見表2。中子通量密度為 3.641×1015cm-2·s-1，中子輻照時(shí)間為1 a。中子能譜為100群的GAM-II格式能譜［9］。

AISI304活化例題的計(jì)算結(jié)果以及其與FISPACT-2007程序計(jì)算結(jié)果之間的對比如表3所示。由于計(jì)算所得核素較多，僅選取典型核素的結(jié)果展示。

圖2 CPAP程序計(jì)算流程圖2Fig.2 Calculation flowchart 2 of CPAP code

圖3 CPAP程序界面展示Fig.3 Interface display of CPAP code

表2 AISI304的元素成分Table 2 Elemental content ofAISI304

表3中，Nf表示使用FISPACT-2007程序計(jì)算出的核素的原子數(shù)量；Nc表示使用CPAP程序中的CRAM求解器計(jì)算出的核素的原子數(shù)量；Nm表示使用CPAP程序中的矩陣指數(shù)求解器計(jì)算出的核素的 原 子 數(shù) 量 ；delta_c=（Nc-Nf）/Nf，delta_m=（Nm-Nf）/Nf。

3.1.2 混凝土活化計(jì)算測試

在本測試算例中，混凝土的初始質(zhì)量為1 kg，其密度為2.35 g·cm-3。元素質(zhì)量成分見表4。中子通量密度為3.641×1015cm-2·s-1，中子輻照時(shí)間為1 a。中子能譜為100群的GAM-II格式能譜［9］。

混凝土活化例題的計(jì)算結(jié)果以及其與FISPACT-2007程序計(jì)算結(jié)果之間的對比如表5所示。由于計(jì)算所得核素較多，僅選取典型核素的結(jié)果進(jìn)行展示。

3.1.3 測試結(jié)論

經(jīng)測試，對于測試算例，CPAP程序的兩種算法在活化計(jì)算時(shí)大部分核素的計(jì)算結(jié)果與FISPACT-2007程序的計(jì)算結(jié)果之間的誤差范圍均在1%以下。但在H和He的計(jì)算時(shí)CPAP程序的計(jì)算結(jié)果與FISPACT-2007程序有較大誤差，這主要是因?yàn)镕ISPACT-2007程序在計(jì)算H和He時(shí)實(shí)際上是合并計(jì)算了（n，p）和（n，α）反應(yīng)產(chǎn)生的質(zhì)子和α粒子，而CPAP程序不把（n，p）和（n，α）反應(yīng)產(chǎn)生的質(zhì)子和α粒子計(jì)算為H和He元素，所以FISPACT-2007程序計(jì)算結(jié)果高于CPAP程序的計(jì)算結(jié)果。

表3 算法測試時(shí)AISI304活化例題的計(jì)算結(jié)果Table 3 Calculation result ofAISI304 activation during algorithm test

表4 混凝土的元素成分Table 4 Elemental content of concrete

3.2 截面庫和衰變數(shù)據(jù)庫的測試

3.2.1 測試方法

CPAP程序的截面庫和衰變數(shù)據(jù)庫的測試采用與ORIGEN-S程序?qū)Ρ确治龅姆绞健＞唧w方法為首先使用CPAP程序的截面庫和衰變數(shù)據(jù)庫替換ORIGEN-S程序自帶的截面庫和衰變數(shù)據(jù)庫，然后基于替換截面庫和衰變數(shù)據(jù)庫后的ORIGEN-S程序構(gòu)建典型算例，并與原版ORIGEN-S程序的計(jì)算結(jié)果進(jìn)行對比。

由于CPAP程序的截面數(shù)據(jù)庫為47群的多群截面庫，而ORIGEN-S程序最終讀取和使用的中子反應(yīng)截面為單群截面，由此需要對CPAP程序的截面數(shù)據(jù)庫進(jìn)行并群處理。并群所基于的中子能譜為使用蒙特卡羅輸運(yùn)計(jì)算程序?qū)Φ湫蛪核训亩研具M(jìn)行輸運(yùn)計(jì)算后按照BUGLE-96庫中能群邊界統(tǒng)計(jì)得出的分為47群的中子能譜。

表5 算法測試時(shí)混凝土活化例題的計(jì)算結(jié)果Table 5 Calculation result of concrete activation during algorithm test

3.2.2 測試算例

對于CPAP程序的截面庫和衰變數(shù)據(jù)庫的測試算例見表6。其中混凝土和AISI304的元素成分分別見表2和表4。

選取各算例中典型核素的計(jì)算結(jié)果與ORIGEN-S程序原版計(jì)算結(jié)果進(jìn)行對比，對比結(jié)果分別見表7和表8。由對比結(jié)果可知，部分核素存在小于10%的偏差，其可能原因?yàn)榻孛鏀?shù)據(jù)庫替換時(shí)并群處理考慮的中子譜的差異以及截面加工所基于的評價(jià)庫版本的差異。

表6 CPAP數(shù)據(jù)庫的測試算例Table 6 List of test cases for database of CPAPcode

表7和表8中，Ns表示使用ORIGEN-S程序原數(shù)據(jù)庫計(jì)算出的核素的原子密度，cm-3；Nc表示使用CPAP程序中的數(shù)據(jù)庫替ORIGEN-S程序數(shù)據(jù)庫后計(jì)算出的核素的原子數(shù)量，cm-3。delta_c=（Nc-Ns）/Ns。

表7 數(shù)據(jù)庫測試時(shí)AISI304活化例題的計(jì)算結(jié)果Table 7 Calculation result ofAISI304 activation during database test

3.2.3 測試結(jié)論

由表7和表8中測試結(jié)果可知，對于所測試的算例，使用CPAP程序的截面庫和衰變數(shù)據(jù)庫所得絕大部分核素的計(jì)算結(jié)果與ORIGEN-S程序自帶數(shù)據(jù)庫所得計(jì)算結(jié)果之間的偏差在1%以內(nèi)，少量核素存在不超過10%的偏差?？紤]到測試時(shí)截面數(shù)據(jù)庫替換過程中并群處理考慮的中子譜的差異以及截面加工所基于的評價(jià)庫版本的差異，以上偏差可認(rèn)為合理。

4 結(jié)語

本文基于堆外設(shè)備和墻體材料中活化產(chǎn)物的產(chǎn)生和遷移的機(jī)理建立了堆外設(shè)備和墻體材料中活化產(chǎn)物源項(xiàng)的計(jì)算模型，并編制了計(jì)算程序CPAP1.0版?；诘湫筒牧系幕罨}，采用國外同類型軟件對比計(jì)算的方式對CPAP程序的算法和數(shù)據(jù)庫分別進(jìn)行了測試，測試結(jié)果顯示：CPAP程序計(jì)算精度可滿足工程設(shè)計(jì)的需要。CPAP程序具有人機(jī)界面友好以及求解器可選的優(yōu)點(diǎn)，可廣泛應(yīng)用于壓水堆核電站的設(shè)計(jì)、運(yùn)行和退役階段，解決堆外設(shè)備和墻體材料中活化源項(xiàng)計(jì)算的相關(guān)問題。

表8 數(shù)據(jù)庫測試時(shí)混凝土活化例題的計(jì)算結(jié)果Table 8 Calculation result of concrete activation during database test