SuperMC可視化方法及其在ITER Clite模型上的驗(yàn)證

曹 佩，甘 佺，曹承龍，郭 里

SuperMC可視化方法及其在ITER Clite模型上的驗(yàn)證

曹 佩1,2，甘 佺2,3,*，曹承龍1,2，郭 里4

（1. 中國科學(xué)技術(shù)大學(xué)，安徽 合肥 230027；2. 中國科學(xué)院合肥物質(zhì)研究院 核能安全技術(shù)研究所，安徽 合肥 230031；3. 中子科學(xué)國際研究院，山東 青島 266199；4. 國家核應(yīng)急技術(shù)支持中心，北京 100080）

隨著反應(yīng)堆模型設(shè)計的日益復(fù)雜，不直觀地物理參數(shù)設(shè)置方式費(fèi)時且容易出錯，大大增加了蒙卡輸運(yùn)計算前期的計算成本。為了增強(qiáng)蒙特卡羅輸運(yùn)軟件的可視化程度，本研究以中子輸運(yùn)設(shè)計與安全評價軟件系統(tǒng)SuperMC為平臺，提出了一種基于數(shù)據(jù)預(yù)處理的權(quán)窗下邊界值分布可視化方法，對生成的權(quán)窗下邊界值的可靠性進(jìn)行實(shí)時的檢查。本文在ITER Clite模型上進(jìn)行可視化方法的驗(yàn)證，結(jié)果表明：借助已有的源信息分布可視化方法，可以直觀地檢查源位置、能量以及方向分布信息設(shè)置的正確性；通過可視化地分析權(quán)窗下邊界值分布的覆蓋趨勢和平滑度，能夠及時終止生成權(quán)窗的迭代過程，輸運(yùn)計算結(jié)果也證明使用可視化方法選擇的權(quán)窗下邊界值可以替代人工經(jīng)驗(yàn)，進(jìn)一步提升蒙卡軟件的易用性。

ITER Clite；源信息分布；權(quán)窗下邊界值分布；可視化方法

隨著聚變工程研究的深入，聚變系統(tǒng)的設(shè)計變得越來越詳細(xì)和復(fù)雜[1,2]。在屏蔽計算中，源信息和權(quán)窗下邊界值比較復(fù)雜并且設(shè)置方式可能比較繁瑣。手動設(shè)置方法費(fèi)時且容易出錯，很難快速得到高質(zhì)量的權(quán)窗下邊界值來提高計算結(jié)果的收斂性。例如，在分析ITER中生物屏蔽插件的中子屏蔽能力時[3-5]，計算內(nèi)層中子分布信息需要花費(fèi)大量時間。源信息設(shè)置的正確性和合理性完全取決于用戶的理解和嘗試性計算。另一方面，在ITER屏蔽計算中，近年來采用的全局權(quán)窗法在加快計算收斂速度上具有突出的優(yōu)勢[6,7]。但是，在其生成的WWINP文件中，權(quán)窗下邊界值多且不可讀。因此，如果源信息和權(quán)窗下邊界值的分布能夠以直觀的形式顯示出來，將為用戶提供一種方便的有效檢查方式。

可視化技術(shù)作為數(shù)值計算中的一種數(shù)據(jù)分析的手段，以其可視性和靈活性的優(yōu)點(diǎn)越來越受到人們的關(guān)注。許多工具（如VisIt、Voxler、Paraview和Amira）目前都能夠可視化大型數(shù)據(jù)集，但大多數(shù)都用在3D建模中。而眾所周知的MCNP程序本身具備一定二維可視化的能力，但需要借助ADVANTG 軟件實(shí)現(xiàn)權(quán)窗下邊界值的三維可視化，以此來檢查權(quán)窗生成器設(shè)置的合理性。應(yīng)用過程中需要反復(fù)在不同程序間進(jìn)行切換，容易出錯，易用性差。

本文研究基于SuperMC平臺，在前期實(shí)現(xiàn)復(fù)雜源信息可視化方法的基礎(chǔ)上[14]，為了解決無法高效檢查生成的權(quán)窗下邊界值是否合理的問題，提出了權(quán)窗下邊界值的三維可視化方法，實(shí)時監(jiān)控權(quán)窗下邊界值的覆蓋范圍和趨勢。并通過ITER Clite例題進(jìn)行了測試。本文的蒙特卡羅計算使用SuperMC進(jìn)行。

1　可視化方法

1.1　自適應(yīng)可視化源信息方法

對于聚變系統(tǒng)模型，源抽樣粒子的分布是不規(guī)則的，目前的聚變模型中包含至少15 000個實(shí)體和29 000個面定義以及超過百萬量級的源粒子。因此，大量的幾何數(shù)據(jù)和粒子數(shù)據(jù)被一一渲染，這將嚴(yán)重影響整個可視化的速度。為了加快源粒子的可視化速度，本研究使用已實(shí)現(xiàn)的基于自適應(yīng)分段樹的源粒子可視化方法。它使用分段結(jié)構(gòu)來組織幾何信息，如圖 1所示。

圖1　分段結(jié)構(gòu)

此外，分段與分段之間具有繼承關(guān)系，通過修改父節(jié)點(diǎn)的分段的屬性，就可以更改其子節(jié)點(diǎn)相應(yīng)的屬性，大大提高了可視化方法的性能。

源抽樣粒子具有位置，能量和方向等特征信息，為了自由切換源抽樣粒子的可視化場景，如圖2所示，展示了分割樹的自適應(yīng)創(chuàng)建過程，通過對葉分段的特征進(jìn)行不斷的提取，形成分段樹結(jié)構(gòu)。

圖2　自適應(yīng)分段樹生成過程

1.2　權(quán)窗下邊界值可視化方法

使用SuperMC全局權(quán)窗生成方法能夠?qū)崿F(xiàn)蒙特卡羅計算的快速收斂，該方法根據(jù)一個粒子進(jìn)入網(wǎng)格后對粒子密度分布均勻性所做的貢獻(xiàn)計算得到每個網(wǎng)格的重要性，再依據(jù)網(wǎng)格重要性計算出權(quán)窗參數(shù)。同時，通過自適應(yīng)迭代計算逐步優(yōu)化權(quán)窗參數(shù)，將粒子均勻地輸運(yùn)到整個模型空間，從而降低全局計算結(jié)果的方差[15]。該方法已經(jīng)在ITER Clite模型上有了很好的驗(yàn)證[7]。然而為了生成全局權(quán)窗下邊界值，計算需要迭代上百次直到覆蓋比率達(dá)到設(shè)定的百分比。但大多數(shù)應(yīng)用場景下，只需要迭代到目標(biāo)計算區(qū)域收斂。目前權(quán)重窗口文件中大約有280 000個下邊界值，數(shù)據(jù)太多，無法直接檢查權(quán)窗下邊界值的可用性。

當(dāng)粒子將要進(jìn)入的第1個網(wǎng)格中，為了確保粒子能最大程度的保留，粒子的權(quán)重最好滿足如下公式（3）：

綜合以上公式（2）和（3），得出公式（4）：

綜上所述，本文相鄰兩個網(wǎng)格之間的權(quán)窗下邊界值比率計算如下：

其中公式（5）中，i,j,k表示三維空間中坐標(biāo)為（，，）的網(wǎng)格的權(quán)窗下邊界值，i,j,k+1表示三維空間中坐標(biāo)為（，，1）的權(quán)窗下邊界值，是沿方向的兩個相鄰網(wǎng)格之間的比值。隨后，在公式2中對比值進(jìn)行比較，最終形成平滑度矩陣Ratio。得到矩陣后，本文將其映射到三維空間去。

原始權(quán)窗下邊界值按一定規(guī)則放置在輸出文件中。為了呈現(xiàn)它，需要對每一個位置的權(quán)窗下邊界值進(jìn)行處理并存儲在合適的結(jié)構(gòu)中。如圖3所示是在可視化時每一個存儲單元WW，其中表示權(quán)窗下邊界值，（，，）表示位置信息，Ratio_r，Ratio_z和Ratio_t分別表示沿半徑方向的比率，沿軸方向的比率以及沿角度方向的比率。圖中顯示的是圓柱坐標(biāo)系下的數(shù)據(jù)存儲結(jié)構(gòu)，直角坐標(biāo)系下與其相似，不再贅述。之后WW結(jié)構(gòu)被存放在一個list里面。

圖3　權(quán)窗參數(shù)的存儲結(jié)構(gòu)

通過對上述存儲結(jié)構(gòu)中的數(shù)據(jù)進(jìn)行可視化，可以直觀地顯示權(quán)窗下邊界值，從而及時監(jiān)測其分布和趨勢。

2　ITER基準(zhǔn)模型

2.1　ITER Clite模型

圖4為本文使用的基本參考模型ITER Clite[16]。圖4所示為等離子體源區(qū)不規(guī)則的幾何結(jié)構(gòu)，瀏覽計算文件，源粒子在幾何結(jié)構(gòu)中分布不均勻。有50個分布區(qū)間，每個區(qū)間又劃分為33個分段，很難檢查設(shè)置的源粒子位置分布是否符合預(yù)期。

圖4　ITER Clite 模型

為了驗(yàn)證基于SuperMC的參數(shù)可視化方法在聚變堆中實(shí)施局部計算的效果，本研究選擇TF和PF線圈計算中子通量。圖5為ITER超導(dǎo)磁體系統(tǒng)，其中PF線圈的主要作用是產(chǎn)生軸向磁場約束等離子體，TF線圈的作用相當(dāng)于控制等離子體位置[17]的平衡、成形和穩(wěn)定。本文在保留TF線圈以及PF線圈結(jié)構(gòu)的同時，對該模型的外圍結(jié)構(gòu)進(jìn)行了簡化，以保證驗(yàn)證的快速進(jìn)行。在計算之前，本文使用第1小節(jié)描述的可視化方法對設(shè)置的源信息和權(quán)窗下邊界值分別進(jìn)行快速檢查和生成監(jiān)測。

圖5　ITER超導(dǎo)磁系統(tǒng)

2.2　源信息可視化

系統(tǒng)中的抽樣粒子起源于設(shè)置的源區(qū)域。蒙特卡羅輸運(yùn)計算之前，源粒子位置，能量，方向等的分布特征進(jìn)行提取。根據(jù)圖4所示，等離子體粒子位置主要分布在①地區(qū)，它是一個不規(guī)則的形狀，粒子能量服從高斯分布和方向是各向同性的。

圖6為等離子體源粒子的三維分布，圖6（a）為源粒子的位置和能量分布。與圖6（b）相比，107個粒子可以更充分地填充等離子體區(qū)域。能量統(tǒng)計結(jié)果如圖7所示，從中可以看出大量的粒子能量在14～14.2 MeV之間。但曲線1的波動幅度比曲線2大，因此有107個粒子的圖更符合標(biāo)準(zhǔn)高斯分布；方向分布如圖8所示，圖8（a）和圖8（b）均為各向同性分布，但抽樣粒子數(shù)為107時圖像呈現(xiàn)出較好的粒子方向分布均勻性。

綜上所述，當(dāng)粒子數(shù)大于 107時，等離子體源的位置、能量和方向分布會更符合輸運(yùn)計算的要求，這說明當(dāng)前源位置、能量以及方向的設(shè)置是合理的，在實(shí)際運(yùn)行中可能根據(jù)不同的實(shí)驗(yàn)要求，將設(shè)置不同的粒子數(shù)。文獻(xiàn)也表明，ITER運(yùn)行時的中子源強(qiáng)最高可達(dá)到1.8×1021n/s[1]。

圖6　ITER Clite源抽樣粒子位置分布

圖7　ITER Clite模型源抽樣粒子能量分布

圖8　ITER Clite模型源抽樣粒子方向分布

2.3　權(quán)窗下邊界值可視化

SuperMC通過迭代生成全局的權(quán)窗下邊界值。理論上，只有當(dāng)權(quán)窗下邊界值覆蓋率達(dá)到設(shè)定值時，才會終止權(quán)窗生成器的計算。經(jīng)驗(yàn)表明，當(dāng)生成的權(quán)窗下邊界值覆蓋目標(biāo)計數(shù)區(qū)域時，進(jìn)行輸運(yùn)計算便可以得到較好的結(jié)果。因此，對生成的權(quán)窗下邊界值進(jìn)行可視化，以決定是否截斷生成過程是非常必要的。

如圖9所示是采用全局權(quán)窗法生成權(quán)窗下邊界值的變化過程。從第10代到第48代，權(quán)窗下邊界值的覆蓋度向外擴(kuò)展，迭代到第47代時，已經(jīng)覆蓋了PF和TF兩個線圈的面積。

圖9　采用全局權(quán)窗法生成權(quán)窗下邊界值

在本例題中，粒子是由內(nèi)部向外部擴(kuò)散，因此主要關(guān)心沿半徑方向的比值。圖10為方向比值分布。如前所述，如果相鄰網(wǎng)格之間的比率可以保持在0.2～5之間，權(quán)窗下邊界值是平滑的。從圖中可以看出，權(quán)窗下邊界值覆蓋的趨勢符合預(yù)期，經(jīng)過47次迭代，目標(biāo)區(qū)域被完全覆蓋，說明這種劃分可以使粒子沿預(yù)期方向飛行并覆蓋整個模型。結(jié)合圖9和圖10可知，第47次迭代的覆蓋面積與第48次迭代基本一致。因此，選擇第47次迭代的結(jié)果作為輸運(yùn)計算的權(quán)窗下邊界值，及時終止迭代程序，減少計算資源的浪費(fèi)。

圖10　相鄰網(wǎng)格之間的比值

為檢驗(yàn)上述可視化方法得到的源信息以及權(quán)窗下邊界值的可靠性，在相同抽樣粒子數(shù)下進(jìn)行了三次輸運(yùn)計算，分別為使用第47代權(quán)窗下邊界值、使用最后一代權(quán)窗下邊界值和未使用窗口的TF線圈和PF線圈的中子通量。TF線圈和PF線圈的中子通量是計算ITER其他運(yùn)行時物理量的重要參數(shù)，快速準(zhǔn)確地計算這些物理量對提高ITER輻射屏蔽設(shè)計效率有重要的作用。表1列出了計算的TF線圈、PF線圈的中子通量和相對統(tǒng)計誤差表，同時本文使用FOM品質(zhì)因子來對收斂效果進(jìn)行定量評估，F(xiàn)OM定義如公式（7）所示：

式中：——輸運(yùn)計算的時間；

σ——第個計數(shù)柵元的統(tǒng)計誤差；

——該次輸運(yùn)計算的計數(shù)柵元總數(shù)。

觀察圖9和圖10所示的可視化結(jié)果，使用全局方法生成的第47代權(quán)窗下邊界值可以提供良好的覆蓋面積和平滑度。同時，由表1可以明顯看出，使用第47代權(quán)窗下邊界值進(jìn)行輸運(yùn)計算，結(jié)果的收斂性足夠好，接近使用最后一代權(quán)窗下邊界值的計算收斂性。

表1　TF線圈和PF線圈的中子通量計算結(jié)果

由表2可以看出，在輸運(yùn)計算中粒子數(shù)相同的情況下，采用第47代權(quán)窗下邊界值的輸運(yùn)計算時間與采用最后一代權(quán)窗下邊界值的輸運(yùn)計算時間相似，收斂效果也接近。雖然沒有使用全局權(quán)窗下邊界值的計算時間短，但品質(zhì)因子FOM的最終結(jié)果與使用權(quán)窗下邊界值的結(jié)果有很大差異。為了達(dá)到與使用權(quán)窗同樣的收斂效果，需要將粒子數(shù)增加至少100倍，對計算資源的要求會加大，同時計算時間也會增加至少100倍。結(jié)合表1和表2的數(shù)據(jù)可以看出，本文提出并實(shí)施的方法在對源采樣粒子合理性的預(yù)判斷和對權(quán)窗下邊界值生成的監(jiān)控方面是可行的。

表2　測試結(jié)果總結(jié)

3　結(jié)論

本研究實(shí)現(xiàn)了權(quán)窗下邊界值分布的可視化方法并將其集成到SuperMC中，并且同時在ITER Clite模型上對復(fù)雜源信息和權(quán)窗下邊界值分布可視化進(jìn)行了測試。使用可視化技術(shù)，直觀地在ITER Clite這種復(fù)雜模型上檢查源分布信息，包括源抽樣粒子的位置、能量和方向。此外，通過監(jiān)測權(quán)窗下邊界值的生成過程中的變化趨勢、覆蓋范圍以及平滑度，可以及時終止迭代程序，減少不必要的迭代次數(shù)，節(jié)省計算資源。結(jié)果進(jìn)一步表明，利用本文提出的可視化方法對源設(shè)置正確性進(jìn)行預(yù)判和對權(quán)窗下邊界值進(jìn)行檢查是可行的。然而，蒙卡輸運(yùn)計算中還有很多需要檢查的復(fù)雜參數(shù)和數(shù)據(jù)分布，未來將對蒙特卡羅軟件的可視化功能進(jìn)一步完善。

［1］ Y.Wu. Fusion Neutronics［M］．Atomic Energy Press，2016.

［2］ R.Pampin，A.Davis，J.Izquierdo，et al.Developments and needs in nuclear analysis of fusion technology［J］．Fusion Eng.Des.，2013，88：454-460.

［3］ Lee，Yi-Kang.Radiation shielding calculate-ons for a 3D ITER benchmark model using TRIPOLI-4? Monte Carlo code［J］．Fusion Eng.Des.，2018，136：602-607.

［4］ H.L.Swami，S.Sharma，A.K.Shaw，et al.Neutron design optimization of ITER TBM port#2 bio-shield plug［J］．Fusion Eng.Des.2018，137：49-55.

［5］ 莊思璇，宋婧，楊琪，等．基于SuperMC的ITER下窗口生物屏蔽插件屏蔽分析［J］．核科學(xué)與工程，2017，37（002）：235-241.

［6］ A Highland Technology V680 TDC（Time to Digital Converter）Driver for EPICS［EB/OL］．http：//www. epics.org/Tdc/Tdc.htm，2005-2-2.

［7］ He，Peng，Zhang，et al.Verification of Su-perMC with ITER C-Lite neutronic model［J］．Fusion Engineering & Design，2016.

［8］ ORNL. ADVANTG［EB/OL］．https：//www.ornl.gov/ onramp/ advantg.

［9］ Mosher，S.W，Johnson，S.R，Bevill，A.M et al. ADVANTG-An Automated Variance Reduction Parameter Generator，Rev.1［R］．2015.

［10］ RSICC.The Visual Editor for MCNP［EB/OL］．http：// mcnpvised.com/visualeditor/visualeditor.html

［11］ Q.Gan，B.Wu，S.Yu，et al.CAD-based hierarchical geometry conversion method for modeling of fission reactor cores［J］．Annals of Nuclear Energy，2016，94：369-375.

［12］WU Yican.Multi-functional Neutrinics Calculation Methodology and Program for Nuclear Design and Radiation Safety Evaluation［J］．Fusion Science and Technology，2018：1-9.

［13］ Y.Wu，J.Song，H.Zheng，et al.CAD-Based Monte Carlo Program for Integrated Simulation of Nuclear System SuperMC［J］．Annals of Nuclear Energy，2015，82：161-168.

［14］曹佩，龍鵬程，甘佺，等．蒙特卡羅粒子輸運(yùn)計算中源可視化分析方法及應(yīng)用［J］．核技術(shù)，2019，42（02）：64-71.

［15］吳宜燦，張澍，宋婧，等．一種基于粒子密度不均勻性的迭代蒙卡全局權(quán)窗參數(shù)生成方法：CN 201610685275［P］．2017-01-04.

［16］ D.Leichtle，B.Colling，M.Fabbri，et al.The ITER t-okomak neutronics reference model C-Model［J］．Fusion Eng.Des，2018，136：742-746.

［17］ N. Mitchell，A.Devred.The ITER magnet system：configuration and construction status［J］．Fusion Eng.Des，2017，123：17-25.

The Visualization Method of SuperMC and Its Validation on ITER Clite Model

CAO Pei1,2，GAN Quan2,3,*，CAO Chenglong1,2，GUO Li4

（1. University of Science and Technology of China，Hefei of Anhui Prov.230027，China；（2. Institute of Nuclear Energy Safety Technology，HFIPS，Chinese Academy of Sciences，Hefei of Anhui Prov.230031，China；3. International Academy of Neutron Sciences，Qingdao of Shandong Prov.266199，China；4. National Nuclear Emergency Response Technical Support Center，Beijing，100080，China）

As reactor model designs become more complex，the unintuitive physical parameter setting method is time-consuming and error-prone，which greatly increases the calculation cost in the early stage of Monte Carlo transport calculation. In order to enhance the visualization of the modeling process for Monte Carlo transport code. This study takes Super Multi-functional Calculation Program for Nuclear Design and Safety Evaluation（SuperMC）as the platform，a visual modeling method based on data pre-processing for weight window bounds distribution was proposed to enhance the visualization of Monte Carlo transport software. Combing with existed approach for visualizing complex source distribution，the ITER Clite model with complex physical parameters was adopted to show the effect of the visualization methods. With the help of visualization method，the correctness of source location，energy and direction distribution information can be visually checked. Running code can be terminated according to visually analyzing the coverage area and smoothness of the weight window bounds distribution. The results of transport calculation also prove that the weight window bounds selected by the visualization method can replace the manual experience and further improve the usability of the Monte Carlo software.

ITER Clite；Source information distribution；Weight window bounds distribution；Visualization methods

TL329+.2

A

0258-0918（2021）05-0899-08

2010-02-22

ITER中子學(xué)分析項(xiàng)目

曹佩（1994—），女，河南南陽人，研究生，現(xiàn)主要從事反應(yīng)堆模擬以及堆芯監(jiān)測技術(shù)方面研究

甘佺，E-mail：quangan@fds.org.cn