中子導管輸運特性的數(shù)值計算方法

王 燕 李新喜 黃朝強 陳 波

(中國工程物理研究院核物理與化學研究所，綿陽 621900)

中子已廣泛地應用于材料科學以及基礎(chǔ)科學研究中，中子導管則廣泛地應用于反應堆上，以傳輸和過濾中子[1–3]。中子的輸運特性，是中子光學設(shè)備理論設(shè)計和優(yōu)化的一個重要方面。常用的中子輸運模擬計算程序有Mcstas[4]與Vitess[5]，系Mont Carlo法的中子輸運軟件，可從網(wǎng)上下載獲得。Pusenkov等[1]介紹了他們的數(shù)值計算軟件的優(yōu)越性，相對于Mcstas和Vitess軟件，數(shù)值計算方法能精確計算中子軌跡，尤宜應用于有缺陷的中子光學輸運設(shè)備。而許多中子光學設(shè)備并非十分理想，Mcstas和Vitess軟件難以避免產(chǎn)生輸入誤差。我們基于C++控制臺應用程序編寫了數(shù)值計算程序NCMP，本文利用該程序計算了兩個中子導管的中子輸運特性，通過分析和比較，為獲得程序的計算誤差及可能的進一步開發(fā)提供參考依據(jù)。

1 基本原理

中子導管數(shù)值模擬計算中，導管出口注量率可表示為：

式中，T(l)是中子導管的透射系數(shù)(Mean accepted divergence)[1]，對于矩形中子導管，

式中，H(l)和V(l)分別是橫向和縱向的透射系數(shù)，

式(3)中，θmin(x) ≤ θ ≤ θmax(x), 0 ≤ x ≤ w; 式(4)中，θmin(y) ≤ θ ≤ θmax(y), 0 ≤ y ≤ h。每個中子的波長lm以及其在導管入口坐標xj和飛行角度qi決定其到達出口的可能性，即

式中，R(4παn/lm)為相應波長的反射率，N為中子在導管內(nèi)的反射次數(shù)，an為第n次反射時的掠入射角。

2 結(jié)果與討論

2.1 理想多道彎導管模擬計算

為得到更好的幾何近似，采用了彎曲部分的切分近似，中子導管被分成很多短的、直的部分。短者的高和寬在出入口有所不同，第k段的示意描述見圖1。圖中，是入口和出口參數(shù)，Φwk和Δwk為角度與尺寸的偏差，玻璃基底的厚度為δk。

計算采用Maxwell波長分布無窮大面源FM(λ)=2F0(lT)4l–5e–lT/l，積分通量 F0=1.4×1013n·cm–2·s–1，特征波長lT= 0.168 nm。多道彎導管超鏡因子m=2，全反射率為 1(理想情況)。導管在橫向上(x方向)彎曲，多道也在橫向上排列。彎曲中子導管相關(guān)導管參數(shù)見表1。

表1 計算采用的主要輸入?yún)?shù)Table 1 The main input parameters for calculation

圖2給出了NCMP計算的結(jié)果(與文獻[1]的結(jié)果一致)。對短波長的中子(<0.1 nm)，橫向透射系數(shù)基本為 0，即彎導管對短波長中子的過濾很好。對于長波長中子，波長越大，透射系數(shù)也越大，因為長波長中子的全反射臨界角大。圖2表明，波長大于0.2nm的中子的透射系數(shù)與波長基本成線性關(guān)系，這是因為計算時把中子在導管中的反射小于臨界角的反射率都假設(shè)成 1，而大于臨界角反射率為0。

圖2 NCMP計算彎導管縱向和橫向透射系數(shù)Fig.2 The mean accepted divergence of guide C3 calculated by NCMP.

由式(4)、(5)，長波長中子縱向透射系數(shù)與波長線性關(guān)系可近似為:

式中的系數(shù)kni= 0.0173 rad/nm為天然鎳的中子全反射臨界角與波長比[6]，m為超鏡因子。

對長波長中子的橫向透射系數(shù)，由于多道的影響，須在式(6)中乘一個因子G，G為導管入口中子通過截面積與去掉多道后的總截面積之比，即

式中參數(shù)見表1。則長波長橫向透射系數(shù)可近似為：

式(6)、(7)、(8)是在反射率設(shè)為1的理想情況下，實際的導管反射率并非如此，故需用數(shù)值模擬。

2.2 冷中子導管C3的模擬計算

上文中，將理想導管全反射率假設(shè)為 1，但實際的導管不會如此理想。C3導管的反射率模型為：

式中q為動量轉(zhuǎn)移，m為導管超鏡因子，rmax、rmin分別為導管反射率最大值以及截斷處反射率。

從慢化劑出來的中子源可看作是面源，面源尺寸90 mm(W)×200 mm(H)。其理論中子能譜由俄羅斯圣彼得堡核物理研究所的Pusenkov教授提供，如圖3所示。

圖3 中子注量率計算中使用的液氫冷中子源初始譜Fig.3 The energy spectrum of liquid hydrogen moderator in calculation.

中子導管長1.6 m的直導管入口與慢化劑的距離為2.295 m，接長度為13.45 m、曲率半徑為556m的彎中子導管，再接長度為31.45m的中子直導管。整個導管系統(tǒng)超鏡因子均為 m=1.5，導管光學截面均為 30mm(W)×200mm(H)。表 2是計算程序所用到的C3導管系統(tǒng)參數(shù)，其中導管表面彎曲度以及單段導管間隔厚度精度是本數(shù)值模擬程序為了更準確模擬導管所用到的特色參數(shù)[1]。

表2 C3導管的輸入?yún)?shù)表Table 2 The main input parameters of C3 neutron guide for calculation

圖4給出了對C3冷中子導管出口注量率的計算比較結(jié)果，從圖中可以看出，利用NCMP計算得到的結(jié)果與 Pusenkov教授提供的計算結(jié)果基本上一致，偏差的引起可能是彎曲導管的分段數(shù)與Pusenkov所用分段數(shù)不一致以及精度步長不同，導致對彎曲部分導管的計算精度不一致。從中子注量率的結(jié)果來看，差距非常小。

圖4 出口注量率與俄羅斯Pusenkov教授提供的計算結(jié)果比較Fig.4 Comparion between the calculated flux by Professor V.M.Pusenkov and NCMP.

C3冷中子導管的透射系數(shù)與中子波長的函數(shù)關(guān)系見圖5，可見中子導管對波長大于0.25nm中子具有良好的中子傳輸性能，能較好地完成冷中子的輸運，對短波長中子則能很好地過濾。

圖5 NCMP 計算的C3導管透射系數(shù)與波長關(guān)系曲線Fig 5 The mean accepted divergence in XY plane and index calculated by NCMP.

中子在中子導管內(nèi)的輸運性能，與中子在 XY平面內(nèi)的發(fā)散密切有關(guān)。XY平面內(nèi)的發(fā)散影響中子通量的大小，對于帶有一定準直功能的多道中子導管，導管輸出 XY平面內(nèi)的發(fā)散情況也反映了中子導管的中子輸運特性。C3冷中子導管XY平面內(nèi)的發(fā)散情況計算結(jié)果見圖6。

圖6 NCMP計算的橫向和縱向發(fā)散情況（X，橫向；Y，縱向）Fig.6 The divergence in X-Y plane calculated by NCMP.

圖7給出的是俄羅斯Pusenkov教授的計算結(jié)果，圖6的發(fā)散情況與之基本相同，僅在圖形的光滑度上略有差別，這可能是由于計算采用的精度步長不同所致。

圖7 Pusenkov教授提供的橫向和縱向發(fā)散情況（X，橫向；Y，縱向）Fig.7 The calculation result of divergences in X-Y plane supplied by Professor V.M.Pusenkov.

綜上所述，我們用自行編寫的數(shù)值計算程序NCMP計算了兩個中子導管模型、理想的多道彎導管以及俄羅斯Pusenkov教授提供的冷中子導管C3。

對于理想的多道彎導管，計算結(jié)果與文獻[1]中提供的透射系數(shù)以及理論分析完全相符。

對于冷中子導管 C3，本文計算結(jié)果與文獻[1]基本上吻合(圖4)，僅由于模擬計算所采用的精確度不同而略有差異。對于源與C3導管組成的整個系統(tǒng)的透射系數(shù)，Pusenkov教授未提供相應的計算輸入?yún)?shù)，無法將NCMP程序計算結(jié)果(圖5)作對比；但由圖4中的注量率計算結(jié)果，用NCMP軟件計算的XY平面內(nèi)可接受的發(fā)散及透射系數(shù)是準確的，可供導管優(yōu)化或?qū)Ч芎笞V儀優(yōu)化的輸入特性計算作參考。

3 結(jié)論

本文介紹了 NCMP數(shù)值計算方法模擬程序基本原理和相關(guān)計算，并與文獻[1]以及 Pusenkov教授提供的計算結(jié)果作比較，結(jié)果表明該程序可以進行準確的中子輸運特性模擬計算。由于數(shù)值計算方法的靈活性和準確性，針對一些有缺陷的中子光學輸運設(shè)備，通過模型上的處理和修正，能夠獲得更為精確的計算結(jié)果。因此，該程序可為中子光學輸運設(shè)備的理論模擬計算提供更多的選擇，為相應的中子光學設(shè)備理論設(shè)計或優(yōu)化提供參考。

1 Pusenkov V M, Schebetov A, Gibcus H P M, et al.Numerical calculation of neutron fluxes at the exit of a complex neutron-guide system at IRI, Delft [J].Nucl Instr Meth, 2002, A492: 105–116

2 van Well A A, de Haan V O.The average number of reflections in a curved neutron guide [J], Nucl Instr Meth,1991, A309: 284–286

3 Mildner D F R, Chen-Mayer H H, Lamaze G P, et al.Characterization of a cold neutron beam from a curved guide [J].Nucl Instr Meth, 1998, A413: 341–351

4 http://www.mcstas.org.2008 [OL]

5 http://www.hmi.de/projects/vitess, 2008 [OL]

6 丁大釗, 葉春堂, 趙志祥, 等.中子物理學[M].北京:原子能出版社, 2005 DING Dazhao, YE Chuntang, ZHAO Zhixiang, et al.Neutron Physics [M].Beijing: Atomic Energy Press, 2005