, 牛风雷华北电力大学 核科学与工程学院，北京 102206
收稿日期：2022-12-02
基金项目：国家自然科学基金资助项目(12105101, 12027813)；核反应堆系统设计技术重点实验室基金资助项目(KFKT-05-FW-HT-20220010)
作者简介：吴浩(1989—)，男，讲师。E-mail: wuhao@ncepu.edu.cn
摘要：在高温气冷堆(high-temperature gas-cooled reactor，HTGR)堆芯球床中，燃料球间的辐射换热是重要的传热模式，与堆芯固有安全特性密切相关。该文利用机器学习方法提出了球床颗粒间辐射角系数智能预测方案，其中基础计算模型基于角系数显式解析表达式，合理描述了球床热辐射特性随球心距变化规律和周围颗粒球平均阻挡作用，用于快速计算球床堆中辐射角系数的核心主导部分。利用高温气冷堆示范项目(HTR-PM)球床堆积结构和光线追踪方法，建立了高温堆球床高精度角系数大数据集，共包含1.66×107条角系数工况，覆盖了球床各种局部结构。利用大数据训练后的梯度提升决策树模型有效提升了角系数预测精度，综合基础计算模型后角系数回归系数超过0.999。该文成果为高温气冷堆球流传热研究、堆芯优化和热工流体分析提供了高效的辐射传热计算方法。
关键词：高温气冷堆球床热辐射颗粒尺度角系数梯度提升决策树机器学习
Machine learning model of radiation heat transfer in the high-temperature nuclear pebble bed
WU Hao
, NIU FengleiSchool of Nuclear Science and Engineering, North China Electric Power University, Beijing 102206, China
Abstract: [Objective] In high-temperature gas-cooled reactors (HTGR), experimental results show that thermal radiation plays a significant role in the heat transfer process and it is highly related to the inherent safety of the nuclear reactor. The HTGR core is a dense pebble bed of single-sized fuel spheres. It is a fundamental and very challenging task to discuss the particle-scale radiative heat transfer in a pebble bed. In this study, using a machine learning approach, an artificial intelligent regression model is developed and trained on a large dataset to predict the obstructed view factor between fuel spheres in a nuclear pebble bed. [Methods] Comparing the numerical and experimental results, the radiative transfer equation (RTE) method significantly underestimates the radiative heat flux in dense pebble beds. The Lagrangian discrete element method (DEM) is often applied in pebble flow simulations. Thus, particle-particle interactions in the DEM framework are used in this present study. The view factor in the base model is calculated using an explicit analytical expression with the elliptic integral function. The model reasonably describes the effect of the distance between spheres and the average obstruction contribution of the surroundings. The obstruction function is proposed by fitting the numerical results. The analytical base model is used to efficiently obtain the dominant parts of the view factor. Furthermore, for the machine learning model, a large dataset of the view factor is established by the particle-scale DEM packing of the HTR-PM, an HTGR demonstration plant, and the thermal ray tracing method. Following preprocessing, the dataset contains a total of 16.6 million records under various conditions in a pebble bed. The gradient boosting decision tree (GBDT) model is used to learn the rules for view factor regression. The model input is the Cartesian coordinates of the sphere and its surrounding ones. The model output is the difference between the ground truth of the view factor and the analytical base model's prediction. 80% of the dataset is used for training, and 20% is left for validation. The mean square error (MSE) and coefficient of determination are selected to evaluate the machine learning model. The GBDT model was trained using the open-source software LightGBM, and the hyperparameter tuning was performed in the FLAML platform to find the best model parameters. [Results] The MSE of the trained machine learning model decreases gradually as the model complexity increases. The analytical base model provides a generally satisfactory forecast of the view factor in the pebble bed, and the gradient boosting decision tree model trained by big data greatly improved the prediction accuracy. With the base model, the coefficients of determination of the trained machine learning model are greater than 0.999. [Conclusions] This study presents an efficient artificial intelligent model for obstructed view factor prediction for heat transfer research, parameter optimization, and thermal hydraulic analysis of the nuclear pebble bed. The trained machine learning model can also be used in effective thermal conductivity analysis, and it is feasible to be coupled with the CFD-DEM simulations of conduction and heat convection in large-scale nuclear pebble beds.
Key words: high-temperature gas-cooled reactorpebble bedthermal radiationparticle scaleobstructed view factorgradient boosting decision treemachine learning
模块式高温气冷堆(high-temperature gas-cooled reactor，HTGR)是具有固有安全特性的第四代先进球床式核动力反应堆[1-2]。堆芯在严重事故工况下仅需依靠自然的热传递就能将堆芯余热有效载出。近年来，高温球床传热成为国内外****的研究焦点[3]。
德国、南非和中国先后开展了大型球床高温传热实验[4-6]，测量了不同温度下的球床平均有效导热系数，结果指出高温(超过700 ℃)下球床堆辐射换热发挥主导作用。早期颗粒辐射研究主要集中在辐射传输方程上，将研究对象简化为均匀的吸收、发射和散射介质，主要应用在大气、粉尘和气溶胶等稀相颗粒系统。针对稠密的高温球床，Gusarov[7]指出辐射传输方程预测结果显著低于球床辐射等效导热系数的实验值。
在微观层次研究中，离散单元法广泛应用在球床球流运动、堆积和搅混中[2, 8]。由于球床堆堆积结构复杂，颗粒间的辐射换热受到周围球的阻挡，直接计算辐射换热角系数极为烦琐。在计算流体力学-离散单元法(computational fluid dynamics-discrete element method, CFD-DEM) 仿真中常采用近距离截断模型[9-10]，但难以描述远程颗粒对之间的辐射特征。Johnson等[11]提出的径向距离近似方案仍无法反映随机堆积下的球床局部热辐射特征。
机器学习的发展为解决高温球床辐射传热提供了新的研究视角，以梯度提升决策树(gradient boosting decision tree, GBDT)[12]算法为代表的人工智能方法能够高效捕获复杂传热问题的非线性特征[13]。基于以上背景，本研究利用机器学习中的GBDT算法，结合高温球床堆中颗粒间辐射规律，训练开发颗粒间辐射角系数智能回归预测模型。
1 机器学习模型针对典型的10 MW高温气冷实验堆(HTR-10)堆芯[14]，利用离散单元法获取球床随机堆积结构，如图 1所示。在高温堆球床的颗粒尺度辐射模型中，任意2个燃料球间的净辐射换热量Qr，ij可表示为
| $Q_{\mathrm{r}, i j}=\varepsilon \sigma_{\mathrm{b}} A_i X_{i j}\left(T_j^4-T_i^4\right).$ | (1) |
 |
| 图 1 球床随机堆积结构示意图 |
| 图选项 |
其中：ε、σb和Ai分别表示颗粒表面发射率、黑体辐射常数和颗粒i的表面积，高温堆中颗粒发射率取0.8，球形燃料元件设计直径为60 mm；Ti和Tj分别表示颗粒i和j表面的绝对温度；Xij表示受到其他周围颗粒阻挡下的颗粒间角系数。
此外，基于连续介质动力学假设，大型高温球床宏观辐射有效导热系数为
| $k_{\mathrm{r}}=\varepsilon \sigma T^3 \cdot \frac{2}{3} \rho \bar{A} \frac{1}{M} \sum\limits_i \sum\limits_j X_{i j} r_{i j}^2 .$ | (2) |
可以看出，在球床辐射换热模型中无论是微观颗粒尺度计算还是宏观分析，仅Xij需要根据相对堆积结构计算，存在辐射换热关系的颗粒不仅包括周围相互接触的颗粒，还包括距离较远的颗粒。通过计算发现，需要考虑周围3~4倍直径范围内约110个颗粒的辐射作用。颗粒间角系数计算极为复杂，利用传统数值方法在中央处理器(central processing unit，CPU)上计算1×104个颗粒间角系数矩阵约需1.6×107 s，难以用于分析超过42×104个燃料球的大型高温堆球床。此时燃料球间热辐射角系数的机器学习计算模型可表示为
| $X_n=B_{0, n}+\sum\limits_{k=1}^N f_{k, n} .$ | (3) |
基础计算模型主要基于角系数的解析理论。利用角系数的定义，目前仅针对无阻挡的2个等径颗粒球情况存在显式解析解，具体可表示为
| $\begin{gathered}X_0\left(r_{12}\right)=\frac{1}{2}- \\\frac{1}{6 \pi}\left[\left(4-r_{12}^2\right) K\left(\frac{2}{r_{12}}\right)+\left(r_{12}^2+4\right) E\left(\frac{2}{r_{12}}\right)\right] .\end{gathered}$ | (4) |
| $\begin{gathered}X_0\left(r_{12}\right)=\frac{1}{2}- \\\frac{1}{3 \pi r_{12}}\left[\left(r_{12}^2-4\right) K\left(\frac{r_{12}}{2}\right)+\left(r_{12}^2+4\right) E\left(\frac{r_{12}}{2}\right)\right] .\end{gathered}$ | (5) |
针对3个颗粒球情况(n=3)，角系数的主要影响因素包括球心距r12和外部颗粒球的阻挡，如图 2所示。角系数X12的半解析基础计算模型为
| $B_{0, 3}=X_0\left(r_{12}\right) \cdot g\left(s_3\right) \text {. }$ | (6) |
 |
| 图 2 球床中多个颗粒相对位置阻挡示意图 |
| 图选项 |
其中：g(s)为颗粒阻挡函数；s3表示外部第3个颗粒中心到颗粒1和2连线的无量纲距离，即s3= x/d，d表示颗粒直径，x表示距离矢量。同样，s4和s5表示外部第4和第5个颗粒到连线的无量纲距离。当s3≥1时，外部颗粒距离较远，不参与颗粒1和2之间的辐射换热，即g=1；当s3=0时，外部颗粒刚好在颗粒1和2中间，完全阻挡了两者之间的热辐射换热，即g=0。结合光线追踪数值计算结果(见图 3)，阻挡函数可近似拟合为
| $g(s)=5.1358 s^5-12.023 s^4+6.6396 s^3+1.248 s^2 .$ | (7) |
 |
| 图 3 角系数基础计算模型阻挡函数 |
| 图选项 |
该函数描述了外部颗粒对辐射换热角系数的阻挡作用，在计算中无需考虑较远的周围颗粒的影响，具体评估详见第2节。若周围存在多个颗粒，阻挡作用更为复杂，难以用简单表达式准确计算。这里近似取乘积作用，对于n≥4，多球阻挡的情况(见图 2)，角系数X12的基础计算模型调整为
| $B_{0, n}=X_0\left(r_{12}\right) \cdot\left[\prod\limits_{m=3}^n g\left(s_m\right)\right]^{1.4} .$ | (8) |
高温球床辐射换热的基础计算模型可以捕捉角系数的主要成分，包括核心的沿着球心距离衰减的规律和周围颗粒的平均阻挡作用。由于只采用了解析表达式，计算速度快，在精度需求不高的情况下可直接采用基础计算模型获取球床颗粒间角系数的近似值。GBDT模型用于修正角系数基础模型并进一步提高计算精度，其模型结构为
| $F_N=\sum\limits_{k=1}^N f_{k, n}=f_{N, n}+\sum\limits_{k=1}^{N-1} f_{k, n}=f_{N, n}+F_{N-1} .$ | (9) |
| $C_m=X_n\left(\boldsymbol{x}_m\right)-B_{0, n}\left(\boldsymbol{x}_m\right) .$ | (10) |
| $Y_m=\frac{C_m-\bar{C}_m}{\sigma}.$ | (11) |
| $\mathrm{MSE}=\frac{1}{M} \sum\limits_{m=1}^M\left(Y_m-\hat{Y}_m\right)^2$ | (12) |
| $R_f^2=1-\frac{\sum\left(Y_m-\hat{Y}_m\right)^2}{\sum\left(Y_m-\bar{Y}\right)^2}=1-\mathrm{MSE} .$ | (13) |
2 结果与讨论2.1 球床角系数数据集本研究基于模块式高温气冷堆示范项目HTR-PM堆芯[16]球床离散堆积结构(见图 4)获取用于训练机器学习模型的角系数数据集。球床直径3 m，随机填充约42×104个直径为60 mm的燃料球。球床局部燃料球可与大量颗粒球产生辐射换热，任意颗粒对间的角系数均会受到多个周围颗粒阻挡, 故在数据集处理中删除周围不参与颗粒间辐射换热、sm>1的颗粒球，进而建立含有不同颗粒数量的子数据集。
 |
| 图 4 HTR-PM堆芯球床 |
| 图选项 |
对于机器学习模型的输出值，利用光线追踪算法[17]，在英伟达GeForce RTX 3090显卡平台上并行计算每种工况下的热辐射角系数。计算中光线数量均固定为1.0×108，在无阻挡的2个颗粒间角系数数值解与解析解的最大绝对偏差为7.8×10－6，如图 5所示。
 |
| 图 5 HTR-PM球床中2个颗粒间角系数 |
| 图选项 |
由于高温球床中燃料颗粒密集堆积，周围颗粒对角系数阻挡显著。对比数值结果发现，距离大于8倍颗粒半径(r12≥8R)的颗粒对间辐射换热极为微弱(ΔX < 10-4)，角系数近似为0。在颗粒数量超过12时，最外部颗粒对角系数的影响极小(ΔX < 10-4)，故实际计算中可近似忽略n≥12的外部颗粒影响。综合以上规律，结合HTR-PM球床离散堆积特性，分别建立了n=3，4，…，11的多球工况，且r12 < 8R的9个角系数子数据集。每个子数据集的输入维度固定，输出维度均为1，HTR-PM球床共1.66×107条角系数工况，具体信息如表 1所示。
表 1 高温球床角系数数据集
| n | 输入维度 | 数据量 |
| 3 | 9 | 14 156 |
| 4 | 12 | 99 349 |
| 5 | 15 | 499 567 |
| 6 | 18 | 1 388 908 |
| 7 | 21 | 1 691 701 |
| 8 | 24 | 1 777 568 |
| 9 | 27 | 2 555 826 |
| 10 | 30 | 3 755 957 |
| 11 | 33 | 4 864 188 |
表选项
2.2 模型复杂度与回归系数大数据集中的80%数据用于训练模型，剩余20%用于验证。为避免梯度提升决策树模型过拟合，需要合理限制模型的复杂度。梯度提升决策树模型复杂度可定义为
| $C=N \cdot L.$ | (14) |
 |
| 图 6 梯度提升决策树模型的复杂度 |
| 图选项 |
表 2 球床角系数机器学习模型回归系数
| n | RB2 | N | L | Rf2 |
| 3 | 0.999 | 22 | 24 | 0.971 |
| 4 | 0.998 | 103 | 18 | 0.990 |
| 5 | 0.998 | 628 | 15 | 0.992 |
| 6 | 0.997 | 1 197 | 13 | 0.992 |
| 7 | 0.997 | 1 929 | 10 | 0.992 |
| 8 | 0.994 | 4 165 | 8 | 0.994 |
| 9 | 0.981 | 2 425 | 15 | 0.991 |
| 10 | 0.954 | 5 578 | 10 | 0.987 |
| 11 | 0.914 | 8 549 | 7 | 0.978 |
表选项
利用基础计算模型直接计算辐射角系数，角系数大数据集上的回归系数在颗粒数目较少(n < 9)时精度较高，回归系数大于0.994，而在n≥9时精度有所下降。最优的梯度提升决策树模型用于预测角系数真实值与基础模型的差值，其预测回归系数Rf2均在0.978以上。综合基础预测模型和机器学习模型，角系数预测回归系数在n=3，4，…，11时均大于0.999，同等硬件(英特尔酷睿i9-10900K)上实现了平均13.5×104倍加速。
3 结论本文针对高温气冷堆球床辐射换热过程开展机器学习模型研究，基础计算模型和GBDT模型用于高效预测考虑周围阻挡条件下的颗粒间角系数。
基础计算模型基于角系数的显式解析理论，合理描述了角系数随球心距的衰减和周围颗粒的平均阻挡作用，可快速计算球床堆中辐射角系数的核心主导部分。
利用HTR-PM球床随机堆积结构，结合光线追踪技术建立了高温堆球床高精度角系数大数据集，共包含1.66×107条角系数工况，最大绝对误差约7.8×10－6。角系数大数据集很好地覆盖了各种常规局部结构。
GBDT模型用于优化并高效提升角系数预测精度，训练后的模型能够反映颗粒间辐射角系数计算规则。通过限制模型的复杂度避免了机器学习模型过拟合现象。综合机器学习模型角系数回归系数超过0.999，实现了13.5×104倍加速效果，可用于在颗粒尺度上分析大规模高温球床辐射换热规律。
参考文献
| [1] | 张作义, 吴宗鑫, 王大中, 等. 我国高温气冷堆发展战略研究[J]. 中国工程科学, 2019, 21(1): 12-19. ANG Z Y, WU Z X, WANG D Z, et al. Development strategy of high temperature gas cooled reactor in China[J]. Strategic Study of CAE, 2019, 21(1): 12-19. (in Chinese) |
| [2] | 姜胜耀, 桂南, 杨星团, 等. 球床式高温气冷堆球流及球床辐射的理论与实验研究进展[J]. 原子能科学技术, 2019, 53(10): 1918-1929. JIANG S Y, GUI N, YANG X T, et al. Theoretical and experimental research progress on pebble flow and effective thermal conductivity in pebble-type HTGR[J]. Atomic Energy Science and Technology, 2019, 53(10): 1918-1929. (in Chinese) |
| [3] | PENG Z B, DOROODCHI E, MOGHTADERI B. Heat transfer modelling in discrete element method (DEM)-based simulations of thermal processes: Theory and model development[J]. Progress in Energy and Combustion Science, 2020, 79: 100847. DOI:10.1016/j.pecs.2020.100847 |
| [4] | ROUSSEAU P G, DU TOIT C G, VAN ANTWERPEN W, et al. Separate effects tests to determine the effective thermal conductivity in the PBMR HTTU test facility[J]. Nuclear Engineering and Design, 2014, 271: 444-458. DOI:10.1016/j.nucengdes.2013.12.015 |
| [5] | WU Y Y, REN C, YANG X T, et al. Repeatable experimental measurements of effective thermal diffusivity and conductivity of pebble bed under vacuum and helium conditions[J]. International Journal of Heat and Mass Transfer, 2019, 141: 204-216. DOI:10.1016/j.ijheatmasstransfer.2019.06.071 |
| [6] | 任成, 杨星团, 李聪新, 等. 高温堆球床等效导热系数测量实验加热系统设计及可行性验证[J]. 清华大学学报(自然科学版), 2014, 54(8): 1068-1072. REN C, YANG X T, LI C X, et al. Heating system design and validation for the pebble bed effective thermal conductivity experiment in a high temperature gas-cooled reactor[J]. Journal of Tsinghua University (Science and Technology), 2014, 54(8): 1068-1072. (in Chinese) |
| [7] | GUSAROV A V. Statistical approach to radiative transfer in the heterogeneous media of thin-wall morphology-Ⅱ: Applications[J]. Journal of Heat Transfer, 2019, 141(1): 012701. DOI:10.1115/1.4040958 |
| [8] | GUI N, HUANG X L, YANG X T, et al. HTR-PM-based 3D pebble flow simulation on the effects of base angle, recirculation mode and coefficient of friction[J]. Annals of Nuclear Energy, 2020, 143: 107442. DOI:10.1016/j.anucene.2020.107442 |
| [9] | MEHRABIAN R, SHIEHNEJADHESAR A, SCHARLER R, et al. Multi-physics modelling of packed bed biomass combustion[J]. Fuel, 2014, 122: 164-178. DOI:10.1016/j.fuel.2014.01.027 |
| [10] | ZHOU Z Y, YU A B, ZULLI P. Particle scale study of heat transfer in packed and bubbling fluidized beds[J]. AIChE Journal, 2009, 55(4): 868-884. DOI:10.1002/aic.11823 |
| [11] | JOHNSON E F, TARI I, BAKER D. Radiative heat transfer in the discrete element method using distance based approximations[J]. Powder Technology, 2021, 380: 164-182. |
| [12] | KE G L, MENG Q, FINLEY T, et al. LightGBM: A highly efficient gradient boosting decision tree[C]// Proceedings of the 31st International Conference on Neural Information Processing Systems. Long Beach, USA: ACM, 2017: 3149-3157. |
| [13] | HUGHES M T, KINI G, GARIMELLA S. Status, challenges, and potential for machine learning in understanding and applying heat transfer phenomena[J]. Journal of Heat Transfer, 2021, 143(12): 120802. |
| [14] | GAO Z Y, SHI L. Thermal hydraulic calculation of the HTR-10 for the initial and equilibrium core[J]. Nuclear Engineering and Design, 2002, 218(1-3): 51-64. |
| [15] | WANG C, WU Q, WEIMER M, et al. FLAML: A fast and lightweight AutoML library[C]// Proceedings of 4th Machine Learning and Systems Conference. San Jose, USA, 2021: 434-447. |
| [16] | ZHANG Z Y, DONG Y J, LI F, et al. The Shandong Shidao Bay 200 MWe high-temperature gas-cooled reactor pebble-bed module (HTR-PM) demonstration power plant: An engineering and technological innovation[J]. Engineering, 2016, 2(1): 112-118. |
| [17] | HOWELL J R, DAUN K J. The past and future of the Monte Carlo method in thermal radiation transfer[J]. Journal of Heat Transfer, 2021, 143(10): 100801. |
