本發(fā)明涉及裂變產(chǎn)物核的燃耗計(jì)算方法，具體涉及一種強(qiáng)中子場裂變產(chǎn)物核燃耗的計(jì)算方法。

背景技術(shù)：

在高密度的強(qiáng)中子場中，裂變產(chǎn)物核可能發(fā)生可觀的二次反應(yīng)，使得產(chǎn)物核的產(chǎn)額(或積存量)發(fā)生改變，因此需要通過燃耗計(jì)算，來確定產(chǎn)額的改變量。

目前除了國內(nèi)錢晶的博士論文編寫了FIRENEQ計(jì)算程序，尚無發(fā)現(xiàn)有公開資料介紹強(qiáng)中子場燃耗的計(jì)算方法。與本發(fā)明接近的是反應(yīng)堆的燃耗計(jì)算程序。反應(yīng)堆燃耗計(jì)算方法主要有兩大類：基于燃耗方程組矩陣求解的數(shù)值計(jì)算方法和基于單燃耗鏈的解析方法。如洛斯阿拉莫斯實(shí)驗(yàn)室在上世紀(jì)80年代開發(fā)的Origen2程序，是基于泰勒展開的解析方法。基于數(shù)值計(jì)算的芬蘭VTT中心開發(fā)的Serpent系統(tǒng)中的燃耗計(jì)算模塊、法國開發(fā)的燃耗程序MENDEL、歐洲開發(fā)的活化計(jì)算系統(tǒng)中的核素存量計(jì)算程序FISPACT，關(guān)于這方面的介紹，可以參看吳明宇的博士論文《基于STEP1.0和MCMG-II的輸運(yùn)-燃耗耦合計(jì)算系統(tǒng)的開發(fā)與研究》。下面介紹的是與本發(fā)明最為接近的基于解析法的CINDER90程序和基于數(shù)值方法的FIRENEQ計(jì)算程序。

CINDER90程序通過對核反應(yīng)網(wǎng)絡(luò)進(jìn)行拆分，得到單個(gè)核素鏈，然后通過解析方法得到每個(gè)步長內(nèi)的核素密度的變化。其關(guān)鍵問題之一是反應(yīng)網(wǎng)絡(luò)的拆分，當(dāng)出現(xiàn)循環(huán)以及循環(huán)嵌套時(shí)，拆分比較困難；關(guān)鍵問題二是確保在步長時(shí)間內(nèi)中子通量的近似不變。通過試算，CINDER90的步長不能小于0.1ns，而在這個(gè)步長內(nèi)，強(qiáng)中子場的中子通量變化很大，所以CINDER90不適合于強(qiáng)中子場的燃耗計(jì)算，計(jì)算結(jié)果會出錯(cuò)。FIRENEQ程序采用數(shù)值的GEAR算法，對整個(gè)核反應(yīng)網(wǎng)絡(luò)和局部網(wǎng)絡(luò)進(jìn)行計(jì)算。裂變產(chǎn)額生的產(chǎn)物核有1300多個(gè)，產(chǎn)物核的參與的反應(yīng)有(n，γ)、(n，2n)、(n，3n)、(n，p)、(n，d)、(n，t)、(n，3He)和(n，α)等，是一個(gè)1300x1300的矩陣數(shù)值求解，計(jì)算一次需要幾個(gè)小時(shí)，甚至需要1天或幾天。

技術(shù)實(shí)現(xiàn)要素：

本發(fā)明的目的是提供一種強(qiáng)中子場裂變產(chǎn)物核燃耗的計(jì)算方法，用于解決強(qiáng)中子場裂變產(chǎn)物核的燃耗修正，并有效地提高計(jì)算效率。

本發(fā)明的技術(shù)方案如下：一種強(qiáng)中子場裂變產(chǎn)物核燃耗的計(jì)算方法，包括如下步驟：

(1)將裂變產(chǎn)生的產(chǎn)物核的核反應(yīng)網(wǎng)絡(luò)拆分為電荷數(shù)Z＝22-72的51個(gè)子網(wǎng)；

(2)建立裂變產(chǎn)物核的子網(wǎng)數(shù)據(jù)文件；

(3)根據(jù)裂變產(chǎn)物核的子網(wǎng)數(shù)據(jù)文件，求解各裂變產(chǎn)物核的核密度隨時(shí)間的變化函數(shù)，裂變產(chǎn)物核i核密度變化微分方程如下：

其中，

表示為第i個(gè)產(chǎn)物核素密度隨時(shí)間變化；

右邊第一項(xiàng)為單位時(shí)間裂變產(chǎn)生i核的個(gè)數(shù)，f(t)為裂變率，y_i為裂變產(chǎn)額；

右邊第二項(xiàng)為i核的消耗項(xiàng)，n(t)為中子注量率，σ_out為譜平均消失截面，包括該核的(n，γ)、(n，2n)、(n，3n)反應(yīng)截面；

右邊第三項(xiàng)為i核的產(chǎn)生項(xiàng)，分別為其它核j通過(n，γ)、(n，2n)、(n，3n)反應(yīng)產(chǎn)生i核的貢獻(xiàn)，相應(yīng)的N_j為其它產(chǎn)物核j核密度，σ_in為相應(yīng)的平均截面。

進(jìn)一步，如上所述的強(qiáng)中子場裂變產(chǎn)物核燃耗的計(jì)算方法，步驟(2)中所述的子網(wǎng)數(shù)據(jù)文件包括相同電荷的裂變產(chǎn)物核的ID和核數(shù)據(jù)，其中，對于電荷數(shù)Z，質(zhì)量數(shù)A，同核異能態(tài)I，其ID＝A×1000+Z×10+I；核數(shù)據(jù)包括(n，γ)、(n，2n)、(n，3n)的反應(yīng)截面、獨(dú)立產(chǎn)額Y。

進(jìn)一步，如上所述的強(qiáng)中子場裂變產(chǎn)物核燃耗的計(jì)算方法，步驟(3)中通過隱式龍格-庫塔方法求解各裂變產(chǎn)物核的核密度隨時(shí)間的變化函數(shù)，所述的隱式龍格-庫塔方法采用PERL代碼軟件包MATH-ODE，用PERL語言編寫CABOR程序?qū)崿F(xiàn)。

本發(fā)明的有益效果如下：采用本發(fā)明的計(jì)算方法，建立相應(yīng)的模型，可以快速地進(jìn)行計(jì)算，計(jì)算一次只需要幾分鐘。速度的提高，得益于整個(gè)網(wǎng)絡(luò)的分割。裂變產(chǎn)物大概有1300多個(gè)，求解整個(gè)網(wǎng)絡(luò)，需要求解矩系數(shù)矩陣為1300*1300的微分方程組；分割后成51個(gè)子網(wǎng)后，每個(gè)子網(wǎng)的產(chǎn)物核只有十幾個(gè)，以Z＝58的子網(wǎng)為例，其核素只有11個(gè)，其矩陣為11*11，約為整個(gè)網(wǎng)絡(luò)的14000分之一。循環(huán)所有51個(gè)子網(wǎng)約為整網(wǎng)的280分之一。根據(jù)這樣粗略的估計(jì)，時(shí)間可以縮短為整個(gè)網(wǎng)絡(luò)事件的1/280。如果求解整個(gè)網(wǎng)絡(luò)需要24小時(shí)，那么單個(gè)子網(wǎng)的時(shí)間約為6秒，所有51個(gè)子網(wǎng)的時(shí)間累加約為5分鐘?？梢?，本發(fā)明可以顯著的提高強(qiáng)中子場裂變產(chǎn)物核燃耗的計(jì)算效率。

附圖說明

圖1為任意裂變產(chǎn)物核在強(qiáng)中子場下生成與消失的可能反應(yīng)路徑圖；

圖2為由裂變產(chǎn)物核、核反應(yīng)、核衰變組成的核反應(yīng)網(wǎng)絡(luò)示意圖，節(jié)點(diǎn)表示產(chǎn)物核，箭頭表示可能的核反應(yīng)(參考圖1)，圖中給出的是質(zhì)量數(shù)范圍A＝118-162的網(wǎng)絡(luò)，虛線帶箭頭的是Z＝58的子網(wǎng)拆分線；

圖3為Ce(Z＝58，A＝142-153)同位素反應(yīng)子網(wǎng)絡(luò)示意圖，中間向右表示(n，γ)反應(yīng)，該子網(wǎng)包含許多循環(huán)，適用于數(shù)值求解；

圖4為本發(fā)明的方法用PERL語言編寫CABOR程序?qū)崿F(xiàn)的流程圖；

圖5為本發(fā)明的具體實(shí)施方案的流程圖。

具體實(shí)施方式

下面結(jié)合附圖和實(shí)施例對本發(fā)明進(jìn)行詳細(xì)的描述。

在強(qiáng)中子場中，i產(chǎn)物核密度變化微分方程為：

其中，表示為第i個(gè)產(chǎn)物核素密度[單位：/b-cm]隨時(shí)間變化，右邊第一項(xiàng)為單位時(shí)間裂變產(chǎn)生i核的個(gè)數(shù)，f(t)為裂變率[1/s]，y_i為裂變產(chǎn)額；右邊第二項(xiàng)為i核的消耗項(xiàng)，n(t)為中子注量率[1/cm²-s]，σ_out為譜平均消失截面，主要包括該核的(n，γ)、(n，2n)、(n，3n)反應(yīng)截面；右邊第三項(xiàng)為i核的產(chǎn)生項(xiàng)，分別為其它核(j)通過(n，γ)、(n，2n)、(n，3n)等反應(yīng)產(chǎn)生i核的貢獻(xiàn)，相應(yīng)的N_j為其它產(chǎn)物核核密度，σ_in為相應(yīng)的平均截面。

每個(gè)產(chǎn)物核可以通過其它核的衰變或中子反應(yīng)產(chǎn)生，也可以通過本身的衰變和核反應(yīng)消失。這些衰變和核反應(yīng)包括β^-衰變、(n，γ)、(n，2n)、(n，3n)、(n，p)、(n，d)、(n，t)、(n，³He)、(n，a)反應(yīng)等，如圖1所示。裂變產(chǎn)生1300多個(gè)產(chǎn)物核，其核素和反應(yīng)組成的核反應(yīng)網(wǎng)絡(luò)如圖2所示，質(zhì)量數(shù)的范圍為A＝66-172，電荷數(shù)Z＝22-72。

強(qiáng)中子場的持續(xù)時(shí)間較短，遠(yuǎn)遠(yuǎn)小于衰變過程，衰變過程可以忽略，微分方程(1)式去掉了衰變過程，其它中子反應(yīng)，如(n，p)、(n，d)、(n，t)、(n，³He)、(n，a)等反應(yīng)截面遠(yuǎn)遠(yuǎn)小于(n，γ)、(n，2n)、(n，3n)反應(yīng)截面之和，因此也可以忽略。按照上述假設(shè)，裂變產(chǎn)物核參加的核過程，電荷數(shù)不變，相應(yīng)的核反應(yīng)網(wǎng)絡(luò)方程可以拆分成相同電荷數(shù)的子網(wǎng)，整網(wǎng)可以拆分為Z＝22-72的51個(gè)子網(wǎng)。拆分前的網(wǎng)絡(luò)范圍為A＝72-172，Z＝22-72，圖2(示意圖)給出了Z＝54-62(Xe-Sm)的部分網(wǎng)絡(luò)。拆分的方法是在Z＝22.5，23.5，...，71.5的地方橫向切開，(圖2中帶箭頭的虛線為Z＝58的子網(wǎng)拆分線)，得到Z＝22，23，...，72的子網(wǎng)，圖3為拆分后的Z＝58的子網(wǎng)，并標(biāo)出了(n，γ)、(n，2n)、(n，3n)的核反應(yīng)；如圖3是拆分后Ce(Z＝58)同位素組成的子網(wǎng)。Ce同位素的反應(yīng)子網(wǎng)，存在諸多循環(huán)的復(fù)雜網(wǎng)絡(luò)，因此很難把該子網(wǎng)拆分成無循環(huán)的單線網(wǎng)絡(luò)，無法用解析的方法給出核密度的隨時(shí)間的變化函數(shù)N(t)。但是，由于該子網(wǎng)只有11個(gè)產(chǎn)物核，所以可以用數(shù)值解法進(jìn)行求解，核素少，可以避免速度慢和剛性問題。常用的隱式的龍格-庫塔方法可以滿足要求。隱式龍格-庫塔方法有現(xiàn)成的代碼，本發(fā)明采用公開免費(fèi)的PERL代碼軟件包MATH-ODE(http：//search.cpan.org/dist/MATH-ODE-0.07/)。

生成核反應(yīng)子網(wǎng)，包括數(shù)據(jù)文件和代碼，本發(fā)明的技術(shù)方案如下(如圖4所示)：

1)從數(shù)據(jù)文件library讀入相同電荷的裂變產(chǎn)物核的ID和核數(shù)據(jù)(包括截面和產(chǎn)額)；

a)對于某一電荷Z，質(zhì)量數(shù)A，同核異能態(tài)I，其身份代碼為ID＝A*1000+Z*10+I；

b)核數(shù)據(jù)包括(n，g)，(n，2n)和(n，3n)截面、獨(dú)立產(chǎn)額Y；

2)根據(jù)產(chǎn)物核ID生成子網(wǎng)嵌入式代碼Pu239_Z.code，Z表示電荷數(shù)，以Z＝58為例，生成的代碼如表1行B1-B27所示，信息包括ID，產(chǎn)額，微分代碼；

3)主代碼中利用require(“Pu239_$Z.code”)，引入子網(wǎng)嵌入式代碼，并進(jìn)行求解；

4)得到計(jì)算結(jié)果。

基于以上原理，用PERL語言編寫了CABOR程序，流程圖如圖5所示。主要步驟如下：

A.參考CINDER90的數(shù)據(jù)格式建立輸入數(shù)據(jù)表，文件名為library；

B.根據(jù)library文件，根據(jù)相同電荷的產(chǎn)物核建立子網(wǎng)，生成子網(wǎng)數(shù)據(jù)文件OD-Z.dat，Z＝51-72.該文件包括反應(yīng)截面和產(chǎn)額；

C.生成子網(wǎng)嵌入式代碼Pu239_Z.code，Z＝51-72，該代碼采用PERL語言，可以嵌入主程序；

D.執(zhí)行主程序Cabor.pl，進(jìn)行計(jì)算，代碼見表1；

E.計(jì)算生成輸出文件sum.out，見表2；

F.結(jié)束。

關(guān)鍵語句如表1所示。以求解Z＝58、Ce133-Ce162核素密度為例，表1列出了核心語句和模塊。行A5插入子網(wǎng)代碼(B1-B27)。行A7-A16初始化微分方程組。行A18-29是求解微分方程組的循環(huán)，每個(gè)循環(huán)完成一個(gè)燃耗步長，燃耗時(shí)間到達(dá)給定的時(shí)間即停止。這樣就可以得到不同時(shí)刻的核密度，然后換算成產(chǎn)額，即可得到產(chǎn)額的變化數(shù)據(jù)，表2為輸出的例子。

通過以上方法，實(shí)現(xiàn)了網(wǎng)絡(luò)分割和子網(wǎng)的數(shù)值求解，最后得到產(chǎn)額的變化數(shù)據(jù)。

表1 Cabor.pl主要代碼(A)和PU239_Z.code嵌入式代碼(B)

表2 Ce的計(jì)算結(jié)結(jié)果文件sum.out

顯然，本領(lǐng)域的技術(shù)人員可以對本發(fā)明進(jìn)行各種改動(dòng)和變型而不脫離本發(fā)明的精神和范圍。這樣，倘若對本發(fā)明的這些修改和變型屬于本發(fā)明權(quán)利要求及其同等技術(shù)的范圍之內(nèi)，則本發(fā)明也意圖包含這些改動(dòng)和變型在內(nèi)。