基于神經(jīng)網(wǎng)絡(luò)的球床堆在線核材料衡算方法

劉若晴，馬 濤，張立國(guó)

(清華大學(xué) 核能與新能源技術(shù)研究院，先進(jìn)核能技術(shù)協(xié)同創(chuàng)新中心，先進(jìn)反應(yīng)堆工程與安全教育部重點(diǎn)實(shí)驗(yàn)室，北京 100084)

作為第4代先進(jìn)核能系統(tǒng)的典型代表，球床式高溫氣冷堆(簡(jiǎn)稱球床堆)具有國(guó)際公認(rèn)的良好安全特性[1-2]。隨著球床式高溫氣冷堆的商業(yè)開發(fā)，核安保與核不擴(kuò)散問題備受關(guān)注[3-4]。對(duì)球床式高溫氣冷堆的乏燃料進(jìn)行超鈾含量的準(zhǔn)確計(jì)算，是核安全保障的必然要求[5-6]。

10 MW高溫氣冷實(shí)驗(yàn)堆(HTR-10)是由清華大學(xué)核能與新能源技術(shù)研究院研發(fā)的、主要用于實(shí)驗(yàn)的球床式高溫氣冷堆[7-8]。與傳統(tǒng)反應(yīng)堆相比，該實(shí)驗(yàn)堆的核燃料以大量燃料球的形式依靠重力作用在反應(yīng)堆內(nèi)流動(dòng)，具有不停堆換料的特點(diǎn)[9]。燃料球穿越堆底后，由在線燃耗測(cè)量系統(tǒng)對(duì)其燃耗、部分活度容易測(cè)量的核素活度(簡(jiǎn)稱易測(cè)核素活度)進(jìn)行測(cè)量[10]。燃料球的燃耗達(dá)到卸料燃耗整定值將會(huì)被卸出堆芯并成為乏燃料球，否則返回堆芯再次循環(huán)[11]。

對(duì)于輕水反應(yīng)堆，核燃料在反應(yīng)堆內(nèi)的位置可跟蹤標(biāo)記，其進(jìn)堆時(shí)間、出堆時(shí)間可準(zhǔn)確掌握，卸出堆芯成為乏燃料后，結(jié)合反應(yīng)堆運(yùn)行歷史，可用核素生成與耗減計(jì)算軟件算出乏燃料的燃耗和各種核素的活度，即可根據(jù)初態(tài)和過程信息，利用軟件正向計(jì)算末態(tài)信息[12-14]。對(duì)于球床式高溫氣冷堆，針對(duì)特定的乏燃料球，可用在線燃耗測(cè)量系統(tǒng)測(cè)量其末態(tài)的易測(cè)核素活度，而無法掌握其進(jìn)堆時(shí)間、在反應(yīng)堆內(nèi)的流動(dòng)情況和循環(huán)次數(shù)等信息，故無法直接應(yīng)用核素生成與耗減計(jì)算軟件計(jì)算其燃耗和各種核素活度。本工作的難點(diǎn)在于，需根據(jù)一部分已知末態(tài)信息(易測(cè)核素活度)計(jì)算另一部分未知末態(tài)信息(不易測(cè)核素活度)，而根據(jù)已知末態(tài)信息難以反推出初態(tài)和過程信息，故依賴于初態(tài)和過程信息進(jìn)行正向計(jì)算來得到未知末態(tài)信息的路徑不可行。因而直接通過計(jì)算并結(jié)合少量試驗(yàn)分析的方法確定輻照后燃料中超鈾核素含量的方法并不適用。

考慮到核反應(yīng)過程中各種核素的生成和耗減彼此關(guān)聯(lián)，本文使用神經(jīng)網(wǎng)絡(luò)方法，針對(duì)在線燃耗測(cè)量系統(tǒng)中易測(cè)核素活度與不易測(cè)核素活度間的對(duì)應(yīng)關(guān)系進(jìn)行建模，通過易測(cè)核素活度估算輻照后球形燃料元件中超鈾元素含量。

1 基本原理

在反應(yīng)堆的平衡堆芯中，各種核素的生成和耗減滿足以下方程[15]：

(λp+φσp)Xpp=1，…，N

(1)

其中：λq為第q種核素的衰變常量，s-1；σp為第p種核素譜平均中子吸收截面，cm2；lpq為第q種核素衰變成第p種核素的份額；fpr為第r種核素吸收中子后轉(zhuǎn)變?yōu)榈趐種核素的份額；φ為按位置、能量平均的中子注量率，cm-2·s-1；Xq、Xr和Xp分別為第q、r和p種核素的含量；t為時(shí)間。

式(1)中，λq、σp、lpq、fpr與反應(yīng)堆運(yùn)行歷史無關(guān)，φ與反應(yīng)堆運(yùn)行歷史相關(guān)。文獻(xiàn)[16]采用VSOP程序結(jié)合核測(cè)數(shù)據(jù)和校正的反應(yīng)堆功率，追蹤評(píng)價(jià)HTR-10的運(yùn)行過程，給出不同時(shí)段堆芯中子注量率分布。

一方面，對(duì)于HTR-10，堆芯連續(xù)地裝載新元件，同時(shí)卸出乏燃料元件，最終長(zhǎng)期處于穩(wěn)定的平衡堆芯狀態(tài)[17]，不同核素的生成和耗減滿足該方程，具有彼此關(guān)聯(lián)的特性。另一方面，該微分方程的復(fù)雜性和耦合性決定了以乏燃料球中的易測(cè)核素活度為已知量，難以通過直接求逆問題計(jì)算不易測(cè)核素的含量。復(fù)雜運(yùn)行歷史條件下，特別是含有長(zhǎng)時(shí)間停堆時(shí)段，直接針對(duì)式(1)求逆將面臨不適定問題[18]。

由于未知量與已知量的關(guān)聯(lián)性，以及問題直接求解的復(fù)雜性，本文采用深度神經(jīng)網(wǎng)絡(luò)對(duì)易測(cè)核素活度與不易測(cè)核素活度之間的關(guān)系進(jìn)行正向建模[19]。深度神經(jīng)網(wǎng)絡(luò)作為一種多層前饋神經(jīng)網(wǎng)絡(luò)[20]，采用誤差反向傳播算法[21]，能在反復(fù)迭代過程中調(diào)整網(wǎng)絡(luò)權(quán)重和偏置，對(duì)于存在確定內(nèi)部機(jī)制的問題具有強(qiáng)大的擬合能力，在求解復(fù)雜問題上具有明顯優(yōu)勢(shì)。

2 深度神經(jīng)網(wǎng)絡(luò)模型的建立

2.1 數(shù)據(jù)準(zhǔn)備

文獻(xiàn)[22]給出HTR-10在線燃耗測(cè)量系統(tǒng)中較易測(cè)量的核素，共有11種，包括95Zr、106Rh、137Cs、103Ru、140Ba、140La、134Cs、95Nbm、133Xem、125Sb和131Tem，其計(jì)數(shù)率分別為101.6、30.27、21.94、20.83、10.2、10.04、9.2、8.87、8.87、5.31和5.26 s-1。包括超鈾元素同位素在內(nèi)的其余核素活度不易測(cè)量。因此，在深度神經(jīng)網(wǎng)絡(luò)訓(xùn)練過程中無法直接采用不易測(cè)核素活度的實(shí)際測(cè)量值。受此技術(shù)條件所限，本文神經(jīng)網(wǎng)絡(luò)訓(xùn)練數(shù)據(jù)采用軟件模擬值。

文獻(xiàn)[23]應(yīng)用球流模擬軟件LAMMPS模擬了HTR-10大量燃料球在堆芯內(nèi)的流動(dòng)，結(jié)合HTR-10實(shí)際運(yùn)行功率歷史，應(yīng)用軟件KORIGEN模擬了燃料球進(jìn)入堆芯后的核素生成與耗減過程，計(jì)算出乏燃料球的燃耗及75種核素的活度，其中包括11種易測(cè)核素活度。

多種參數(shù)計(jì)算值與實(shí)驗(yàn)值的對(duì)照表明，軟件模擬結(jié)果能較好地反映燃料球在HTR-10堆芯內(nèi)實(shí)際流動(dòng)和燃耗情況，軟件模擬結(jié)果與實(shí)際物理過程的吻合程度較高[23]。因此，本文以LAMMPS和KORIGEN模擬的大量燃料球的燃耗及各種核素活度為數(shù)據(jù)基礎(chǔ)建立神經(jīng)網(wǎng)絡(luò)。其中11種易測(cè)核素活度作為神經(jīng)網(wǎng)絡(luò)的輸入端數(shù)據(jù)，燃耗和其余64種不易測(cè)核素活度作為神經(jīng)網(wǎng)絡(luò)的輸出端數(shù)據(jù)。每個(gè)燃料球的燃耗和75種核素活度為1組數(shù)據(jù)，每組數(shù)據(jù)包含11個(gè)輸入數(shù)據(jù)和65個(gè)輸出數(shù)據(jù)。

在LAMMPS和KORIGEN中，導(dǎo)入HTR-10 2003年初至2007年5月中旬的實(shí)際運(yùn)行功率歷史數(shù)據(jù)以為KORIGEN進(jìn)行核反應(yīng)模擬提供實(shí)時(shí)功率數(shù)據(jù)，令20 000個(gè)模擬的燃料球先后進(jìn)入堆芯，并對(duì)球流和核反應(yīng)過程進(jìn)行模擬，得出了20 000組模擬數(shù)據(jù)。其中19 500組模擬數(shù)據(jù)作為訓(xùn)練組，用于深度神經(jīng)網(wǎng)絡(luò)的訓(xùn)練；其余500組模擬數(shù)據(jù)作為校驗(yàn)組，用于驗(yàn)證該深度神經(jīng)網(wǎng)絡(luò)的性能。

由于不同核素活度數(shù)量級(jí)相差較大，為避免數(shù)值計(jì)算中病態(tài)方程組的發(fā)散問題[24]，用式(2)對(duì)20 000個(gè)燃料球的燃耗和75種核素活度進(jìn)行歸一化處理：

i=1，…，20 000，j=1，…，76

(2)

經(jīng)歸一化處理，燃耗和所有核素的活度變?yōu)榫禐?、方差為1的分布，有利于深度神經(jīng)網(wǎng)絡(luò)初值賦值和最終結(jié)果收斂。

2.2 網(wǎng)絡(luò)結(jié)構(gòu)與參數(shù)設(shè)置

深度神經(jīng)網(wǎng)絡(luò)包括輸入層、隱含層和輸出層[25]。輸入層節(jié)點(diǎn)數(shù)與輸入數(shù)據(jù)個(gè)數(shù)一致，為11；輸出層節(jié)點(diǎn)數(shù)與輸出數(shù)據(jù)個(gè)數(shù)一致，為65。經(jīng)反復(fù)摸索調(diào)整，將隱含層的層數(shù)設(shè)置為8，每層節(jié)點(diǎn)數(shù)分別為100、300、500、500、500、300、100和80。

在梯度下降方法的選擇上，采用批量隨機(jī)梯度下降方法[26]。經(jīng)探索優(yōu)化，每批量的數(shù)據(jù)組數(shù)設(shè)置為20。若每批量的數(shù)據(jù)組數(shù)較大，有利于早期快速收斂，但其缺點(diǎn)在于精度不高，且收斂接近完成時(shí)波動(dòng)較大。

學(xué)習(xí)效率設(shè)置為0.001。研究發(fā)現(xiàn)，學(xué)習(xí)效率設(shè)置過大，將導(dǎo)致深度神經(jīng)網(wǎng)絡(luò)訓(xùn)練時(shí)無法收斂，設(shè)置過小則不利于提高收斂速度。

模型訓(xùn)練中的懲罰函數(shù)[26]采用下式：

(3)

深度神經(jīng)網(wǎng)絡(luò)的權(quán)重和偏置，通過誤差反向傳播算法進(jìn)行更新。深度神經(jīng)網(wǎng)絡(luò)基于目前業(yè)界流行的、支持GPU加速的Tensorflow框架，使用Python編程語言進(jìn)行構(gòu)建[27-29]。

2.3 模型訓(xùn)練

將19 500組模擬數(shù)據(jù)導(dǎo)入深度神經(jīng)網(wǎng)絡(luò)進(jìn)行訓(xùn)練，發(fā)現(xiàn)收斂曲線(圖1)在初期下降十分迅速，后期下降緩慢，這一結(jié)果符合批量隨機(jī)梯度下降方法的特征。

圖1 深度神經(jīng)網(wǎng)絡(luò)的收斂曲線Fig.1 Convergence curve of deep neural network

3 結(jié)果與討論

訓(xùn)練完成后，將校驗(yàn)組數(shù)據(jù)導(dǎo)入所生成的深度神經(jīng)網(wǎng)絡(luò)模型，每組數(shù)據(jù)中的11種易測(cè)核素活度作為模型的輸入數(shù)據(jù)，經(jīng)計(jì)算，模型會(huì)輸出燃耗和不易測(cè)超鈾元素同位素的活度預(yù)測(cè)值。

圖2 239Pu活度的深度神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)值與校驗(yàn)組原始值的對(duì)比Fig.2 239Pu activity comparison between prediction value by deep neural network and original value

將深度神經(jīng)網(wǎng)絡(luò)模型的預(yù)測(cè)值與校驗(yàn)組數(shù)據(jù)中相應(yīng)的原始值進(jìn)行比較，發(fā)現(xiàn)二者吻合程度較好。以239Pu核素為例，如圖2所示，模型預(yù)測(cè)數(shù)據(jù)與其對(duì)應(yīng)的測(cè)試數(shù)據(jù)基本分布在y=x直線附近，模型對(duì)于239Pu的活度預(yù)測(cè)比較準(zhǔn)確。

燃耗預(yù)測(cè)如圖3所示，模型預(yù)測(cè)的燃耗與對(duì)應(yīng)的原始數(shù)據(jù)中燃耗組成的點(diǎn)也基本分布在y=x直線附近。因此，模型對(duì)于燃耗的預(yù)測(cè)也比較準(zhǔn)確。

圖3 燃耗的深度神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)值與校驗(yàn)組原始值的對(duì)比Fig.3 Burnup comparison between prediction value by deep neural network and original value

本文后續(xù)計(jì)算中的誤差均指神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)值與LAMMPS和KORIGEN軟件模擬值之間的誤差，而非神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)值與該數(shù)據(jù)實(shí)測(cè)值間的誤差。

使用下式求各超鈾元素同位素活度的神經(jīng)網(wǎng)絡(luò)平均相對(duì)誤差：

(4)

使用下式求超鈾元素同位素活度的相關(guān)系數(shù)R2：

(5)

(6)

燃料球燃耗和239Pu活度的相對(duì)誤差如圖4所示。

對(duì)燃料球燃耗和超鈾元素活度的部分統(tǒng)計(jì)學(xué)指標(biāo)進(jìn)行計(jì)算，計(jì)算結(jié)果列于表1。其中，相關(guān)系數(shù)R2基本在1.0附近，燃耗、236U活度、239Pu活度相對(duì)誤差在±5%以內(nèi)的燃料球份額超過99%，相對(duì)誤差在±2.5%以內(nèi)的燃料球份額超過81.6%。可看出，深度神經(jīng)網(wǎng)絡(luò)模型對(duì)軟件模擬數(shù)據(jù)的預(yù)測(cè)精度總體上較好。

圖4 燃料球燃耗和239Pu活度的相對(duì)誤差Fig.4 Relative error of burnup and 239Pu activity

表1 超鈾元素同位素活度預(yù)測(cè)準(zhǔn)確性的統(tǒng)計(jì)Table 1 Statistics of model predicted transuranic isotope activity

4 結(jié)論

本文基于深度神經(jīng)網(wǎng)絡(luò)，使用HTR-10 2003年初至2007年5月中旬的運(yùn)行功率歷史數(shù)據(jù)作為KORIGEN核反應(yīng)模擬的功率數(shù)據(jù)進(jìn)行超鈾元素同位素含量計(jì)算。用11種易測(cè)核素的活度推算出的各超鈾元素同位素活度與校驗(yàn)組中對(duì)應(yīng)值符合程度很好，如校驗(yàn)組數(shù)據(jù)中239Pu預(yù)測(cè)結(jié)果相對(duì)誤差在±2.5%以內(nèi)的份額為81.6%，預(yù)測(cè)精度較高。該結(jié)果說明本文提出的方法處理球床堆燃料中核素含量估計(jì)問題時(shí)具有較好的發(fā)展?jié)摿ΑＴ摲椒ㄒ灿芯窒扌浴Ｊ紫龋灰诇y(cè)核素活度缺乏實(shí)測(cè)數(shù)據(jù)，該方法采用軟件模擬的數(shù)據(jù)作為基礎(chǔ)，除深度神經(jīng)網(wǎng)絡(luò)模型對(duì)軟件模擬數(shù)據(jù)的擬合不到位造成的誤差外，還可能疊加軟件對(duì)實(shí)際物理過程模擬不到位造成的誤差；其次，易測(cè)核素中不同核素活度的測(cè)量難易程度和精度存在差別，使用較多的易測(cè)核素活度作為輸入端數(shù)據(jù)，增加了測(cè)量工作量，降低了測(cè)量精度，不利于模型的實(shí)際應(yīng)用。因此，改進(jìn)LAMMPS和KORIGEN等軟件對(duì)實(shí)際物理過程的模擬精度、改進(jìn)神經(jīng)網(wǎng)絡(luò)模型對(duì)軟件模擬數(shù)據(jù)的擬合精度、以盡可能少的易測(cè)量核素活度實(shí)現(xiàn)對(duì)超鈾元素同位素活度的估算，將是下一步的研究目標(biāo)。