楊洪新 張競(jìng)宇 李平烈 趙曉澤 譚 琳 肖懿鑫
(華北電力大學(xué)核科學(xué)與工程學(xué)院 北京102206)
在反應(yīng)堆物理的領(lǐng)域中,燃耗計(jì)算主要用來分析重核和裂變產(chǎn)物同位素的原子核密度隨時(shí)間的變化趨勢(shì)以及它們對(duì)輸運(yùn)計(jì)算(如中子通量密度分布和反應(yīng)性)的影響[1]。通過對(duì)燃耗的計(jì)算分析可以得到核素成分與原子核密度隨燃耗的變化以及中子通量密度的密度分布等以供進(jìn)一步研究核燃料的轉(zhuǎn)化及循環(huán)、確定堆芯壽期等,為反應(yīng)堆堆芯燃耗管理、燃料增殖以及燃耗信任制技術(shù)應(yīng)用和乏燃料后處理的研究提供基礎(chǔ)[2]。在反應(yīng)堆源項(xiàng)分析領(lǐng)域,活化計(jì)算主要用來分析材料受中子輻照后產(chǎn)生的放射性核素成分和活度,其計(jì)算模型和方法與燃耗計(jì)算是一致的,僅僅是應(yīng)用領(lǐng)域不同,因此不需特意區(qū)別兩者的差別。對(duì)產(chǎn)生的核素進(jìn)行徑跡追蹤以研究核素的產(chǎn)生路徑及每個(gè)路徑的貢獻(xiàn)統(tǒng)計(jì),有助于了解核素的轉(zhuǎn)化途徑以辨識(shí)重要核素或核反應(yīng),是反應(yīng)堆物理和源項(xiàng)分析中十分重要的內(nèi)容。
現(xiàn)有的反應(yīng)堆燃耗/活化計(jì)算方法一般可以分為兩類,分別是基于單條燃耗鏈求解的線性子鏈方法(Transmutation Trajectory Analysis,TTA)和基于燃耗方程矩陣求解的矩陣指數(shù)方法[3]。其中CINDER[4]點(diǎn)燃耗程序和MCB[5]等蒙特卡羅程序的燃耗計(jì)算模塊使用TTA方法,目前TTA方法已經(jīng)可以實(shí)現(xiàn)以常數(shù)項(xiàng)和基于拉普拉斯變換的多項(xiàng)式項(xiàng)為非齊次項(xiàng)的燃耗方程計(jì)算[6]。ORIGEN[7]點(diǎn)燃耗程序則使用多項(xiàng)式逼近方法。SERPENT[8]蒙特卡羅程序燃耗計(jì)算模塊首先使用了切比雪夫有理近似方法(Chebyshev Rational Approximation Method,CRAM)。歐洲活化程序FISPACT在2007版本中使用了歐拉指數(shù)方法,在后續(xù)版本中使用了GEAR 方法[9]。國內(nèi)計(jì)算程序主要有清華大學(xué)開發(fā)的NUIT程序和華北電力大學(xué)核反應(yīng)堆源項(xiàng)研究組開發(fā)的ABURN程序。ABURN程序基于Fortran語言開發(fā),包含了TTA 以及CRAM 等算法,兼容EAF-2007 和ORIGEN 數(shù)據(jù)庫,程序可實(shí)現(xiàn)活化計(jì)算、衰變計(jì)算、脈沖工況計(jì)算、敏感性及不確定性分析等功能。
表1列出了目前典型的可用于復(fù)雜核素系統(tǒng)計(jì)算的通用點(diǎn)燃耗程序的基本情況。
表1 通用的點(diǎn)燃耗程序Table 1 General ignition consumption program
目前國內(nèi)所開發(fā)的程序大多基于數(shù)值算法,只提供了核素成分和核素存量的計(jì)算結(jié)果,并沒有進(jìn)行核素轉(zhuǎn)變過程的反應(yīng)鏈精細(xì)追蹤,而對(duì)于輻射防護(hù)工作者來說,確定放射性核素的來源并制定相應(yīng)的防護(hù)措施是十分重要的。因此本文基于ABURN程序開發(fā)了徑跡追蹤和貢獻(xiàn)統(tǒng)計(jì)的計(jì)算模塊,可在實(shí)際計(jì)算中考察某核素的產(chǎn)生路徑及貢獻(xiàn)百分比,以進(jìn)一步滿足輻射防護(hù)領(lǐng)域的工程應(yīng)用要求。
TTA的基本原理是將復(fù)雜的燃耗鏈拆解成一條條獨(dú)立的線性鏈,然后再對(duì)每一條燃耗鏈進(jìn)行解析求解,最后綜合各線性鏈中同種核素的計(jì)算結(jié)果,即可得到在復(fù)雜網(wǎng)狀燃耗鏈中的對(duì)應(yīng)某一時(shí)刻的核素存量[10]。TTA方法在理論上可以得到每條鏈的精確解,能夠?qū)μ囟ê怂剡M(jìn)行徑跡追蹤和貢獻(xiàn)統(tǒng)計(jì)研究,所以選取TTA方法對(duì)核素徑跡追蹤進(jìn)行研究,給出指定核素的產(chǎn)生路徑并給出每一條產(chǎn)生該核素的路徑的貢獻(xiàn)統(tǒng)計(jì)。采用TTA 方法求解也可以有效避免矩陣剛性問題引起的數(shù)值不穩(wěn)定性,并且在本質(zhì)上是求出了各條反應(yīng)鏈的解析解,所以結(jié)果也更加準(zhǔn)確。
燃耗計(jì)算即在給定的燃耗步長內(nèi)求解核素原子核密度隨時(shí)間的變化情況,有燃耗方程:
在給定燃耗區(qū),給定燃耗步長內(nèi),式(1)可以簡(jiǎn)化為:
活化計(jì)算主要研究的問題是核反應(yīng)堆的材料受到中子輻照時(shí)的核素成分隨時(shí)間的變化,活化方程的表達(dá)式和燃耗方程比較類似,區(qū)別在于通常活化方程不考慮裂變反應(yīng),其主要關(guān)注核素的衰變和活化反應(yīng),因此有活化方程:
式中:Ni(t)表示第i種核素在t時(shí)刻的存量,cm-3;σi表示第i種核素的中子誘發(fā)反應(yīng)微觀截面,cm2;σij表示第j種核素活化生成第i種核素的中子誘發(fā)微觀截面,cm2;λi表示第i種核素的衰變常量,s-1;λij表示第j種核素活化生成第i種核素的衰變常量,s-1;φ表示中子通量密度,cm-2·s-1。
對(duì)于一條含有n個(gè)核素的核素鏈,假設(shè)只有第一個(gè)核素具有初始核素密度N1(0),利用TTA 方法可以計(jì)算得到第k個(gè)核素經(jīng)過時(shí)間t后的核素存量為:
線性子鏈方法將復(fù)雜的網(wǎng)狀燃耗鏈拆解成一條條簡(jiǎn)單的鏈進(jìn)行計(jì)算,因此實(shí)現(xiàn)徑跡貢獻(xiàn)統(tǒng)計(jì)功能,需要利用回溯法對(duì)計(jì)算中的核素和反應(yīng)信息進(jìn)行搜索,要遍歷整個(gè)解空間搜索所有核素可能的反應(yīng)路徑,并使用截?cái)嗾`差進(jìn)行截?cái)啵?1]。
回溯算法是一種類似枚舉的搜索嘗試過程,按照深度優(yōu)先的準(zhǔn)則進(jìn)行搜索,從一個(gè)頭節(jié)點(diǎn)開始,向核反應(yīng)的縱深方向搜索,直到搜索到某一個(gè)核素不再發(fā)生其他反應(yīng)或滿足預(yù)先設(shè)置的截?cái)鄻?biāo)準(zhǔn)時(shí),回溯到反應(yīng)鏈的上一個(gè)核素檢驗(yàn)有無其他未反應(yīng)的反應(yīng)道,若有則以這個(gè)核素為新節(jié)點(diǎn)繼續(xù)往縱深方向搜索,若無則再回溯到該核素的上一個(gè)核素繼續(xù)進(jìn)行檢驗(yàn),直到回溯到反應(yīng)鏈的頭節(jié)點(diǎn),對(duì)所有初始核素完成相關(guān)反應(yīng)的搜索[12]。
回溯法可以在搜索過程中遍歷整個(gè)解空間,但是真實(shí)的反應(yīng)堆級(jí)別的反應(yīng)網(wǎng)是非常復(fù)雜的,并且如果出現(xiàn)閉環(huán)就會(huì)導(dǎo)致反應(yīng)鏈出現(xiàn)無盡循環(huán),因此需要制定一個(gè)截?cái)鄻?biāo)準(zhǔn),在滿足條件的情況下對(duì)反應(yīng)鏈進(jìn)行截止,剔除不重要的反應(yīng)以提高計(jì)算效率[13]。
定義子鏈中每個(gè)核素的截?cái)啾壤齈k(t),即t時(shí)刻核素鏈中當(dāng)前核素的原子核密度占核素鏈母核初始密度的比例份額[14]:
再根據(jù)線性子鏈的截?cái)鄻?biāo)準(zhǔn)判斷是否截?cái)啵?/p>
假設(shè)存在如下核素A生成核素F的反應(yīng)網(wǎng)絡(luò):
根據(jù)線性鏈方法,可以將復(fù)雜的反應(yīng)網(wǎng)絡(luò)拆解成一條條單獨(dú)的鏈進(jìn)行解析計(jì)算,則核素A生成核素F的反應(yīng)網(wǎng)絡(luò)可以拆解為:
在圖1顯示的反應(yīng)網(wǎng)絡(luò)中,初始核素A共有4條路徑生成核素F,在此反應(yīng)網(wǎng)絡(luò)中,各核素是相互獨(dú)立的,因此共有4 組F核素輸出。在計(jì)算中,核素F的總量可以表示為:
圖1 反應(yīng)網(wǎng)Fig.1 Reaction network
如圖2 所示,即在進(jìn)行活化計(jì)算時(shí)展示的初始核素A生成核素F的路徑,其每條路徑的貢獻(xiàn)百分比ηi即此條鏈生成的F核素量所占F核素總量的份額,即:
圖2 反應(yīng)鏈Fig.2 Reaction chain
ABURN 程序的徑跡追蹤與貢獻(xiàn)統(tǒng)計(jì)功能是基于線性子鏈方法結(jié)合回溯算法得到的,理論上可以對(duì)反應(yīng)中涉及到的所有核素進(jìn)行回溯追蹤并將其路徑及貢獻(xiàn)百分比表示出來,因此具備較高的計(jì)算精度。
本文基于TTA 方法自主開發(fā)了徑跡貢獻(xiàn)統(tǒng)計(jì)計(jì)算模塊,并將其集成到華北電力大學(xué)源項(xiàng)研究組自主開發(fā)的活化程序ABURN 中。為了對(duì)ABURN程序的徑跡追蹤計(jì)算結(jié)果進(jìn)行驗(yàn)證,本文選取IAEA發(fā)布的天然鐵活化基準(zhǔn)題、50Cr 活化基準(zhǔn)題[16]以及UKAEA 發(fā)布的95Zr 衰變基準(zhǔn)題[17],ABURN 程序和歐洲活化計(jì)算程序FISPACT 均使用EAF-2007 數(shù)據(jù)庫進(jìn)行核素徑跡貢獻(xiàn)統(tǒng)計(jì)計(jì)算結(jié)果對(duì)比。
在天然鐵的活化計(jì)算中,初始核素包括鐵的4種同位素,包括54Fe、56Fe、57Fe以及58Fe,能群數(shù)為100群,中子通量密度為3.640 99×1015cm-2·s-1,計(jì)算時(shí)間為1 a,截?cái)鄥?shù)選取1.0×10-29,具體參數(shù)如表2 所示。計(jì)算中能群數(shù)為100,能群分布如表3所示。
表2 天然鐵活化計(jì)算輸入?yún)?shù)Table 2 Input parameters for natural Fe activation calculation
表3 能群分布Table 3 Energy group distribution
續(xù)表
表4 比較了活化反應(yīng)產(chǎn)生的各核素原子核數(shù)目,其中A表示ABURN 程序的計(jì)算結(jié)果,F(xiàn)表示FISPACT 程序的計(jì)算結(jié)果,兩者的差異用相對(duì)誤差δ來表示,相對(duì)誤差δ計(jì)算公式為:
天然鐵活化計(jì)算結(jié)果為:
在計(jì)算結(jié)果中選取52V、51Cr、54Mn、55Fe 以及60Co進(jìn)行徑跡追蹤研究(表5)。
表5 天然鐵徑跡追蹤結(jié)果對(duì)比Table 5 Pathway analysis of Fe activation calculation
對(duì)上述結(jié)果進(jìn)行對(duì)比分析可以發(fā)現(xiàn):對(duì)于核素原子核數(shù)目的計(jì)算結(jié)果,ABURN 程序與FISPACT程序的相對(duì)誤差絕對(duì)值均小于0.5%且大部分小于0.1%,說明ABURN程序的計(jì)算結(jié)果精度較高;對(duì)于目標(biāo)核素的路徑及各路徑貢獻(xiàn)百分比的計(jì)算結(jié)果,ABURN 程序與FISPACT 程序基本一致,說明ABURN徑跡追蹤的計(jì)算結(jié)果較為準(zhǔn)確。
在95Zr 衰變計(jì)算中,初始核素為95Zr,初始核素?cái)?shù)目為1×1020,計(jì)算時(shí)間為5.0×107s,截?cái)鄥?shù)選取1.0×10-29,具體參數(shù)如表6 所示。95Zr 計(jì)算結(jié)果見表7。對(duì)各核素進(jìn)行徑跡追蹤研究,結(jié)果如表8所示。
表6 95Zr衰變計(jì)算輸入?yún)?shù)Table 6 Input parameters of95Zr decay calculation
表7 95Zr衰變計(jì)算結(jié)果對(duì)比Table 7 Comparison of95Zr decay calculation results
表8 95Zr計(jì)算徑跡追蹤結(jié)果對(duì)比Table 8 Pathway analysis of95Zr decay calculation
對(duì)上述結(jié)果進(jìn)行對(duì)比分析,可以發(fā)現(xiàn):對(duì)于核素原子核數(shù)目的計(jì)算結(jié)果,ABURN 程序與FISPACT程序的相對(duì)誤差絕對(duì)值均小于0.01%,說明ABURN程序的計(jì)算結(jié)果精度較高;對(duì)于目標(biāo)核素的路徑及各路徑貢獻(xiàn)百分比的計(jì)算結(jié)果,ABURN 程序與衰變鏈理論計(jì)算結(jié)果基本一致,說明ABURN 徑跡追蹤的計(jì)算結(jié)果較為準(zhǔn)確。
在50Cr 的先活化后衰變基準(zhǔn)題計(jì)算中,初始核素為50Cr,能群數(shù)100,能群分布如表3 所示,中子通量密度為3.640 99×1015cm-2·s-1,計(jì)算時(shí)間為先活化計(jì)算1 a,再衰變計(jì)算0.5 a,截?cái)鄥?shù)選取1.0×10-29,具體參數(shù)如表9 所示。50Cr 先活化后衰變計(jì)算結(jié)果見表10。在計(jì)算結(jié)果中選取39Ar、50V、48Ca以及54Mn進(jìn)行徑跡追蹤研究(表11)。
表9 50Cr先活化后衰變計(jì)算輸入?yún)?shù)Table 9 Input parameters of50Cr activation and decay calculation
表10 50Cr先活化后衰變計(jì)算結(jié)果對(duì)比Table 10 Comparison of50Cr activation and decay calculation results
表11 50Cr先活化后衰變計(jì)算徑跡追蹤結(jié)果對(duì)比Table 11 Pathway analysis of50Cr activation and decay calculation
對(duì)上述結(jié)果進(jìn)行對(duì)比分析,可以發(fā)現(xiàn):對(duì)于核素原子核數(shù)目的計(jì)算結(jié)果,ABURN 程序與FISPACT程序的相對(duì)誤差絕對(duì)值均小于0.1%,說明ABURN程序的計(jì)算結(jié)果精度較高;對(duì)于目標(biāo)核素的路徑及各路徑貢獻(xiàn)百分比的計(jì)算結(jié)果,ABURN 程序與衰變鏈理論計(jì)算結(jié)果基本一致,說明ABURN 徑跡追蹤的計(jì)算結(jié)果較為準(zhǔn)確。
對(duì)活化產(chǎn)生的核素進(jìn)行徑跡追蹤貢獻(xiàn)統(tǒng)計(jì)以研究在實(shí)際反應(yīng)中核素的產(chǎn)生路徑及各路徑的貢獻(xiàn)百分比,是反應(yīng)堆物理和源項(xiàng)分析中十分重要的內(nèi)容。本文采用線性子鏈方法,結(jié)合回溯法,實(shí)現(xiàn)在活化計(jì)算中對(duì)指定核素的產(chǎn)生路徑進(jìn)行追蹤分析,并給出每條路徑的貢獻(xiàn)百分比,并與國外廣泛使用的同類軟件FISPACT進(jìn)行了對(duì)比驗(yàn)證,分別驗(yàn)證了天然鐵活化基準(zhǔn)題,95Zr 衰變基準(zhǔn)題,50Cr 先活化后衰變基準(zhǔn)題等三組基準(zhǔn)題,對(duì)其中部分核素進(jìn)行了徑跡追蹤貢獻(xiàn)統(tǒng)計(jì)研究,核素的路徑一致,各路徑的貢獻(xiàn)百分比也基本一致。
經(jīng)過測(cè)試驗(yàn)證,本文自主編寫的徑跡追蹤程序模塊的計(jì)算結(jié)果與FISPACT程序的基本一致,可以認(rèn)為本文徑跡追蹤方法和程序是可行的,從而為活化計(jì)算中重要核素或核反應(yīng)的篩選提供了有力的仿真工具。
作者貢獻(xiàn)聲明楊洪新:實(shí)施研究,采集數(shù)據(jù),分析解釋數(shù)據(jù),起草文章;張競(jìng)宇:指導(dǎo),技術(shù)材料支持,對(duì)文章的知識(shí)性內(nèi)容進(jìn)行審閱;李平烈:技術(shù)材料支持,統(tǒng)計(jì)分析,內(nèi)容校核;趙曉澤:采集數(shù)據(jù),分析數(shù)據(jù);譚琳:采集數(shù)據(jù),分析數(shù)據(jù);肖懿鑫:采集數(shù)據(jù),分析數(shù)據(jù)。