沈華韻，潘流俊，衷 斌

（北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所，北京 100094）

中子輸運(yùn)外源問題是核能開發(fā)利用中廣泛存在且重要的研究課題，如核裝置的輻射屏蔽研究、加速器驅(qū)動次臨界系統(tǒng)（ADS）和聚變裂變混合堆的研發(fā)、在役核電站裝換料階段的臨界安全分析［1］。與屏蔽問題相比，后兩類問題有如下共同特點(diǎn)：系統(tǒng)中的裂變裝置處于次臨界狀態(tài)，中子注量率直接受源強(qiáng)影響。為此，可將其統(tǒng)稱為外源驅(qū)動次臨界系統(tǒng)。

對外源驅(qū)動次臨界系統(tǒng)的中子物理學(xué)問題的研究，同樣可采用確定論（如SN）方法或蒙特卡羅（MC）方法［2］。但由于外中子源通常只分布在局部空間，使得SN方法的射線效應(yīng)嚴(yán)重影響此類問題的計(jì)算精度。以大亞灣核電站為例，一次源的半徑和高僅0.053cm 和4.00cm，而堆芯的半徑和高分別為161cm 和183cm［1］。另外，SN方法難以精確描述源的復(fù)雜幾何結(jié)構(gòu)，在一定程度上也將影響計(jì)算精度。MC 方法雖不存在上述問題，但計(jì)算效率較低，影響了該方法的應(yīng)用。

雖然以ADS和混合堆為代表的外源驅(qū)動次臨界系統(tǒng)是當(dāng)前的研究熱點(diǎn)，國內(nèi)外均已開展或正在開展大量的研究工作，但這些研究往往關(guān)注于裝置的物理設(shè)計(jì)，而未涉及計(jì)算方法［3-5］。為此，本文針對外源驅(qū)動次臨界系統(tǒng)內(nèi)中子輸運(yùn)過程的特點(diǎn)，提出以中子首次裂變?yōu)轳詈宵c(diǎn)的MC／SN耦合計(jì)算方法。通過充分發(fā)揮MC 方法和SN方法各自的優(yōu)點(diǎn)，實(shí)現(xiàn)外源驅(qū)動問題高效高精度的模擬。在算法研究的基礎(chǔ)上實(shí)現(xiàn)MC／SN耦合輸運(yùn)計(jì)算，并進(jìn)行測試模型的計(jì)算與分析。

1 外源問題

中子輸運(yùn)方程可寫成3種不同的形式：通量型、碰撞密度型和發(fā)射密度型［6］?？紤]到算法討論的便利性，本文將采用發(fā)射密度型積分方程。對于外源問題，該方程可寫為：

其中：P 為相空間r×E×Ω 上任意點(diǎn)；χ（P）為P 點(diǎn)的中子產(chǎn)生率，包括外源中子和中子核反應(yīng)產(chǎn)生的次級中子；S（P）為外源中子產(chǎn)生率；散射核Ks（P′→P）為1個P′的中子經(jīng)過輸運(yùn)并發(fā)生非裂變反應(yīng)產(chǎn)生狀態(tài)為P 的中子數(shù)量的期望值；裂變核Kf（P′→P）為1個P′的中子經(jīng)過輸運(yùn)并發(fā)生裂變反應(yīng)產(chǎn)生狀態(tài)為P 的中子數(shù)量的期望值。式（1）的物理含義非常清晰，即中子發(fā)射密度等于源中子發(fā)射密度、非裂變核反應(yīng)次級中子發(fā)射密度和裂變反應(yīng)中子發(fā)射密度之和。

定義：

其中：

在式（2）兩邊對n 求和，同時利用式（3）和（4），得到：

對比式（1）可知，外源問題的解可表示為級數(shù)形式：

2 MC／SN 耦合算法

2.1 基礎(chǔ)理論

由式（2）和（4）可知，χ′（P）滿足如下方程：

χ″（P）滿足如下方程：

即有

其中，根據(jù)式（3），有：

可見，由式（7）與（8）共同給出的定義完全等價于式（1）。對比式（7）與式（1）可知，χ′（P）與χ（P）對應(yīng)的外源均為S（P），僅在中子裂變反應(yīng)的處理方式上存在差異，計(jì)算χ′（P）時需將裂變處理為吸收。因此，物理上χ′（P）包含兩方面對系統(tǒng)內(nèi)中子分布的貢獻(xiàn)——外源中子的直接貢獻(xiàn)和外源中子通過非裂變反應(yīng)產(chǎn)生次級中子帶來的貢獻(xiàn)。

明確χ′（P）的物理含義后，由式（9）表示的S1（P）的物理含義已非常清晰，具有首次裂變中子源的意義，對應(yīng)的包含兩種來源的裂變反應(yīng)產(chǎn)生的中子——外源中子直接引起的裂變和由外源中子發(fā)生非裂變反應(yīng)產(chǎn)生的次級中子引起的裂變。與實(shí)際外源S（P）相比，S1（P）具有如下特點(diǎn)：運(yùn)動方向上各向同性；分布于整個裂變區(qū)等。

最后，對比式（8）與（1）可知，χ″（P）與χ（P）的差異僅表現(xiàn)在對應(yīng)的源分布上。前者為首次裂變中子源，后者為實(shí)際的中子源。顯然，χ″（P）包含了首次裂變中子通過各種核反應(yīng)引起的對系統(tǒng)內(nèi)中子分布的貢獻(xiàn)。因此，物理上也解釋了式（6）的具體含義。

2.2 耦合算法

由上述分析可知，求解式（1）可轉(zhuǎn)化為耦合求解式（7）和（8），這兩個方程的特點(diǎn)決定了各自適宜采用的方法。對于式（7），由于源S（P）通常具有空間局部分布、幾何結(jié)構(gòu)較復(fù)雜等特點(diǎn)，適合采用MC 方法以保證計(jì)算精度。同時，由于已將裂變處理為吸收，不存在中子增殖過程，使得需要跟蹤的中子歷史較短，不會帶來過多計(jì)算量。對于式（8），由于源S1（P）分布于整個裂變區(qū)，且各向同性，可采用SN方法，在保證計(jì)算精度的同時提高計(jì)算效率。基于此，形成如下以中子首次裂變?yōu)轳詈宵c(diǎn)的MC／SN耦合輸運(yùn)算法。

1）MC 方法求解中子輸運(yùn)方程（式（7））。從實(shí)際源分布S（P）中抽取中子，將裂變處理為吸收，并通過計(jì)數(shù)獲得該中子輸運(yùn)過程對應(yīng)的中子注量率φ′（P）和首次裂變源S1（P）。

2）SN方法求解輸運(yùn)方程（式（8））。此時求解的是對應(yīng)于S1（P）的外源問題，并獲得中子注量率φ″（P）。

3）疊加獲得輸運(yùn)方程（式（1））對應(yīng)的中子注量率φ（P）＝φ′（P）＋φ″（P）。

可見，該耦合輸運(yùn)算法的思想是：將1個難以在兼顧精度和效率條件下實(shí)現(xiàn)求解的外源中子輸運(yùn)問題，轉(zhuǎn)化為另外兩個耦合的外源問題，然后分別采用MC和SN方法實(shí)現(xiàn)耦合問題的求解。最后，需特別指出的是，該算法中的MC與SN耦合是易于實(shí)現(xiàn)且嚴(yán)格精確的，這由如下兩方面保證：1）僅存在MC 到SN的單向單次數(shù)據(jù)傳輸；2）首次裂變中子源S′（P）是各向同性的，可精確實(shí)現(xiàn)角度的離散。

3 算例

為驗(yàn)證耦合算法的有效性，構(gòu)造如下裝置：二維xy 幾何，尺寸為20cm×20cm，x＝0和y＝0處為反射邊界，其他為自由面，材料為富集度10%的濃縮鈾。SN計(jì)算時采用S8，共80個離散方向。另外，假定源中子的能量為14.1 MeV，尺寸為1cm×1cm，考慮不同空間分布的情況。在系列算例中，均以等值線的形式給出14.1 MeV 和2 MeV 附近的中子注量率，并對標(biāo)注作如下約定：MCNP表示僅用MC方法計(jì)算的結(jié)果，作為基準(zhǔn)值；MCS8表示采用耦合算法計(jì)算的結(jié)果；S8表示僅用SN方法計(jì)算的結(jié)果。

3.1 模型1

模型1為僅在中心分布著單一中子源的情況，圖1示出了3種計(jì)算方式的結(jié)果。

從圖1a可見，SN方法的射線效應(yīng)非常明顯。在計(jì)算精度上，無論是源中子所在能量附近還是裂變中子最可幾能量（2 MeV）附近，耦合計(jì)算的結(jié)果均能與MC 符合很好，而SN與MC的結(jié)果有明顯差異。

3.2 模型2

在模型1的基礎(chǔ)上，將外源沿x方向移至x＝10cm 處，得到模型2，圖2示出了3種計(jì)算方式的結(jié)果。

圖2a中同樣可見SN方法的射線效應(yīng)對計(jì)算結(jié)果的影響。在計(jì)算精度上，與模型1一致，耦合計(jì)算的結(jié)果與MC 方法的符合得很好，而SN的結(jié)果與MC方法的有明顯差異。

3.3 模型3

在模型3中，將考慮同時存在兩個中子源的情況，即假定同時存在模型1和2的中子源，圖3示出了3種計(jì)算方式的結(jié)果。

模型3的對比結(jié)果與前兩個一致，耦合計(jì)算的結(jié)果與MC 方法的符合得很好，而SN方法的結(jié)果與MC方法的有較大差異。

4 結(jié)論

圖1 模型1中子注量率等值線對比Fig.1 Comparison of neutron fluence rate contour line of model 1

圖2 模型2中子注量率等值線對比Fig.2 Comparison of neutron fluence rate contour line of model 2

圖3 模型3中子注量率等值線對比Fig.3 Comparison of neutron fluence rate contour line of model 3

根據(jù)外源驅(qū)動次臨界系統(tǒng)中子輸運(yùn)的特點(diǎn)，本文提出了以中子首次裂變?yōu)轳詈宵c(diǎn)的MC／SN耦合算法。在該耦合算法中，MC 方法用于模擬外源中子的早期輸運(yùn)過程，直至發(fā)生裂變；對首次裂變中子的模擬，則采用SN方法。由于首次裂變中子源分布于整個裂變區(qū)，且各向同性，采用SN方法已能保證計(jì)算精度。

不同分布源模型的計(jì)算表明，該耦合算法的計(jì)算結(jié)果與作為基準(zhǔn)的MC 方法結(jié)果具有很好的一致性，而采用單一SN方法計(jì)算的結(jié)果則存在較大的偏差，初步驗(yàn)證了基于中子首次裂變的MC／SN耦合算法在模擬局部源問題時的有效性和正確性。

［1］ 沈華韻，李樹，衷斌，等.堆腔注水改進(jìn)項(xiàng)對堆外探測器的影響分析［R］.北京：北京應(yīng)用物理與計(jì)算數(shù)學(xué)研究所，2012.

［2］ LEWIS E E，MILLER W F.Computational methods of neutron transport［M］.US：American Nuclear Society，Inc.，1993：156-358.

［3］ SOOBY E，BATY A，BENE?O，et al.Candidate molten salt investigation for an accelerator driven subcritical core［J］.Journal of Nuclear Materials，2013，440：298-303.

［4］ LIN Z，CHEN J，GUO W，et al.The conceptual design of electron-accelerator-driven subcritical thorium molten salt system［J］.Energy Procedia，2013，39：267-274.

［5］ SIDDIQUE M T，KIM M H.Conceptual design study of Hyb-WT as fusion-fission hybrid reactor for waste transmutation［J］.Annals of Nuclear Energy，2014，65：299-306.

［6］ 裴鹿成，張孝澤.蒙特卡羅方法及其在粒子輸運(yùn)問題中的應(yīng)用［M］.北京：科學(xué)出版社，1980：213-219.