颗粒包覆过程，即在颗粒表面均匀包覆上一层功能材料，具有颗粒保护、表面性能改善、提高环境相容性等用途[1]。其广泛运用于工业生产、制药行业，同时也是第四代先进核能系统——高温气冷堆(high-temperature gas-cooled reactor, HTGR)包覆燃料颗粒制备的重要环节之一。HTGR所用的燃料元件由TRISO(tristructural-isotropic)包覆燃料颗粒和石墨基体组成，TRISO包覆燃料颗粒的4层包覆层(由内至外分别为疏松热解炭层、内致密热解炭层、碳化硅层、外致密热解炭层)通过流化床-化学气相沉积法(fluidized bed-chemical vapor deposition, FB-CVD)连续制备得到[2-3]，如图 1所示。在FB-CVD颗粒包覆过程中，颗粒在高速气流的作用下处于流态化，而气体反应物通过载带的形式进入流化床，在高温区发生化学反应，产生超细粉末沉积在颗粒表面[4]。
|
| 图 1 用FB-CVD制备高温气冷堆TRISO包覆燃料颗粒示意图 |
FB-CVD颗粒包覆过程复杂，包括一系列复杂的转化，如气体分解产生超细粉末、超细粉末团聚或者沉积在颗粒表面、颗粒表面形成一层包覆层成为包覆颗粒等，具有典型的多尺度特征和多学科交叉的特点[5]，目前国内外尚没有完整的理论揭示其完整过程的机理。如图 2所示，FB-CVD颗粒包覆过程中材料属性、工艺参数、转换和产品属性之间层层关联。不同的材料属性，在不同的工艺参数下会发生不同的转化，从而生成不同的产品。
|
| 图 2 FB-CVD颗粒包覆过程中材料属性、工艺参数、转化和产品属性之间的联系 |
FB-CVD颗粒包覆过程也包括微观尺度、介观尺度、颗粒尺度及宏观尺度等多尺度范围，如图 3所示，这些尺度之间相互联系、互相影响，在本构层次上重叠。设计微尺度的颗粒特征可以改善中尺度的本构行为，如对颗粒进行表面改性可以改善颗粒的流动性和分散性。将设备运行参数与材料属性联系起来的本构模型可以预测产品的转化。研究FB-CVD颗粒包覆过程，相比于实验方法成本昂贵、有时候有辐射和高温的危险，计算模拟方法计算速度快、成本低，而且可以全方位获取颗粒包覆过程信息。将微观和宏观方法相结合的多尺度建模对于实现实际的反应器放大和控制非常重要，有利于理解和操控流化床中颗粒包覆过程。
|
| 图 3 FB-CVD颗粒包覆过程的时间范围和尺寸范围 |
目前有许多研究颗粒包覆过程的模型，如表面更新模型(surface renewal model)[6-7]、Monte Carlo模拟方法(Monte Carlo method)[8-9]、分区群体平衡模型(compartment-based population balance model)[10-11]、DEM-C-PBM模型(discrete element method-compartment-population balance model)[12]、CFD-DEM(computational fluid dynamics-discrete element method)模型[13-14]、耦合球谐函数的CFD-DEM模型(CFD-DEM & spherical harmonic formulation)[15]及离散液滴法(discrete droplet method)[16]等。然而，这些方法存在两大局限性。其一，这些方法都是反应器尺度，即宏观尺度，不能直接用于研究多尺度的流化床-化学气相沉积颗粒包覆过程(这包括化学反应等微观尺度)；其二，这些现象模型基于分区概念，将容器内分为包覆区和非包覆区，由于一个容器内各个位置都可能存在待沉积的组分，因此基于这种分区概念的模型在实际的包覆过程中是非常不精确的[17]。Gregor等[18]基于离散颗粒法建立了用于模拟流化床中化学气相沉积的模型，但模型中假设由化学反应引起的局部前驱体的浓度变化及气体分子数目的变化不会显著影响流体动力学，这是很不精确的，并且会导致一些错误。
根据文[14, 19]以及本文课题组之前的研究[17, 20-22]发现，计算流体力学-离散单元(CFD-DEM)方法可以很好地模拟喷动流化床中气-固两相流行为。然而，流化床-化学气相沉积颗粒包覆过程中还涉及化学反应、组分传输导致的温度场和浓度场的变化。因此，需要建立一个多物理场多尺度的耦合模型，综合考虑流化床中由于颗粒流化导致的速度场的变化，壁面加热和化学反应导致的温度场的变化，以及组成传输与沉积导致浓度场的变化，这样才能更加准确地模拟流化床-化学气相沉积颗粒包覆过程。目前，本文基于流化床中的速度场、温度场、浓度场和化学气相沉积(chemical vapor deposition, CVD)模型建立了一个初步的CFD-DEM-CVD多物理场多尺度耦合模型来模拟颗粒包覆过程，并分别提出了化学反应流(chemical reaction flow, CRF)模型、群体平衡模型(population balance model, PBM)及内嵌式颗粒(inner bonding particle, IBP)模型来一步步优化模型，让模拟结果更加准确。本文首先分析了流化床-化学气相沉积颗粒包覆过程的特点；然后提出了颗粒包覆模型，并对模型进行优化；接着将包覆模型用于实际的案例计算，探究通入气速、通入气体种类、壁面温度等条件对颗粒包覆速率的影响；最后根据模拟计算结果得出颗粒包覆过程的一些影响因素，并指出颗粒包覆模型未来的发展方向。
1 流化床-化学气相沉积颗粒包覆模型对于流化床-化学气相沉积颗粒包覆过程，首先基于流化床中的速度场、温度场、浓度场和CVD模型建立了一个初步的CFD-DEM-CVD多物理场多尺度耦合模型来模拟颗粒包覆过程；然后分别提出了CRF模型、PBM模型及IBP模型来一步步优化模型，让模拟结果更加准确。
1.1 CFD-DEM-CVD模型流化床-化学气相沉积颗粒包覆过程涉及气-固两相流运动带来的速度场的变化，壁面加热及化学反应带来的温度场的变化，以及化学反应和组分传输带来的浓度场的变化，这是一个多物理场多尺度过程。根据研究发现，速度场可以通过CFD-DEM模型求解，温度场可以通过热量传输模型求解，浓度场可以通过组分传输方程求解[17, 21]。
1.1.1 速度场流化床中气-固两相流可以通过CFD-DEM模型求解，并得到广泛验证[19]。在CFD-DEM模型中，气相运动通过连续性方程和Navier-Stokes方程求解，固体颗粒相运动通过Newton第二定律和转动方程求解。
| $ \frac{\partial }{{\partial t}}\left( {{\varepsilon _{\rm{g}}}{\rho _{\rm{g}}}} \right) + \nabla \cdot \left( {{\varepsilon _{\rm{g}}}{\rho _{\rm{g}}}{\mathit{\boldsymbol{u}}_{\rm{g}}}} \right) = 0. $ | (1) |
其中：εg为气体体积分数，ρg为气体密度，ug为气体速度。
对于气相，Navier-Stokes方程为[23]
| $ \begin{array}{c} \frac{\partial }{{\partial t}}\left( {{\varepsilon _{\rm{g}}}{\rho _{\rm{g}}}{\mathit{\boldsymbol{u}}_{\rm{g}}}} \right) + \nabla \cdot \left( {{\varepsilon _{\rm{g}}}{\rho _{\rm{g}}}{\mathit{\boldsymbol{u}}_{\rm{g}}}{\mathit{\boldsymbol{u}}_{\rm{g}}}} \right) = \\ - {\varepsilon _{\rm{g}}}\nabla {P_{\rm{g}}} - {F_{{\rm{pg}}}} + \nabla \cdot \left( {{\varepsilon _{\rm{g}}}{\mathit{\boldsymbol{\tau }}_{\rm{g}}}} \right) + {\varepsilon _{\rm{g}}}{\rho _{\rm{g}}}g. \end{array} $ | (2) |
其中：Pg为气相压力，Fp-g为气-固交换系数，τg为气相应力张量，g为重力加速度。
颗粒的运动可以分为平动和转动，在气-固两相流系统中，颗粒i的平动方程为[24-25]
| $ {m_i}\frac{{{\rm{d}}{\mathit{\boldsymbol{v}}_i}}}{{{\rm{d}}t}} = {\mathit{\boldsymbol{F}}_{{\rm{drag}}, i}} + {m_i}g + \sum\limits_{j = 1}^{{n_i}} {\left( {{\mathit{\boldsymbol{F}}_{c, ij}} + {\mathit{\boldsymbol{F}}_{d, ij}}} \right)} . $ | (3) |
其中：mi和vi分别为颗粒i的质量和速度，ni为与颗粒i接触的颗粒数目，Fdrag, i为曳力，Fc, ij和Fd, ij分别为颗粒j对颗粒i的接触力和阻尼力。
| $ {I_i}\frac{{{\rm{d}}{\mathit{\boldsymbol{\omega }}_i}}}{{{\rm{d}}t}} = \sum\limits_{j = 1}^{{n_i}} {\left( {{\mathit{\boldsymbol{\tau }}_{c, ij}} + {\mathit{\boldsymbol{\tau }}_{r, ij}}} \right)} . $ | (4) |
其中：Ii为颗粒i的转动惯量，ωi为颗粒i的转动角速度，τc, ij和τr, ij分别为颗粒i和颗粒j之间的转矩和滚动摩擦转矩。本文选择Hertz-Mindlin模型描述颗粒之间的碰撞。
1.1.2 温度场热量传递主要包括对流、传导及辐射3种形式。在温度场中，热量的传递主要发生在流体与流体、流体与颗粒、流体与壁面、颗粒与颗粒及颗粒与壁面之间。流体与流体之间的热量传递可以通过能量守恒方程求解[17]，
| $ \begin{array}{c} \frac{\partial }{{\partial t}}\left( {{\rho _{\rm{g}}}E} \right) + \nabla \cdot \left( {{\mathit{\boldsymbol{u}}_{\rm{g}}}\left( {{\rho _{\rm{g}}}E + p} \right)} \right) = \\ \nabla \cdot \left( {{k_{{\rm{eff}}}}\nabla T - \sum\limits_j {{h_j}} {J_j} + \left( {{{\bar {\bar \tau} }_{{\rm{eff}}}} \cdot {\mathit{\boldsymbol{u}}_{\rm{g}}}} \right)} \right) + {S_{\rm{h}}}. \end{array} $ | (5) |
其中: ρg为气体密度，ug为气体速度，keff为有效热导率，hj为组分j的显焓，Jj为组分j的扩散通量。式(5)右边的前3项分别表示由于传导、组分扩散和黏性耗散引起的能量转移，Sh表示由于化学反应产生的能量。
气相和固相之间的流体动力耦合体现在动量守恒方程和能量守恒方程中的相间传输系数之中。对于气-固之间热量传递，只考虑热量的传导，热量传导速率q为[17]
| $ q = {h_{{\rm{gsm}}}}A\left( {{T_{\rm{g}}} - {T_{\rm{p}}}} \right). $ | (6) |
其中: hgsm为气-固相热量传输系数，A为颗粒表面积，Tg为气体温度，Tp为颗粒温度。本文选择Li-Mason模型用来描述气-固相间热量传输系数。
流体-壁面之间的热量传递为[26]
| $ \begin{array}{c} q = {h_{\rm{f}}}\left( {{T_{\rm{w}}} - {T_{\rm{f}}}} \right) + {q_{{\rm{rad}}}} = \\ {h_{{\rm{ext}}}}\left( {{T_{{\rm{ext}}}} - {T_{\rm{w}}}} \right) + {\varepsilon _{{\rm{ext}}}}\sigma \left( {T_\infty ^4 - T_{\rm{w}}^4} \right). \end{array} $ | (7) |
其中：hf、Tw、Tf和qrad分别为流体侧局部传热系数、壁面温度、局部流体温度和辐射热通量，hext为外部传热系数，εext为壁面热辐射系数，Text为外部水槽温度，σ为Stefan-Boltzmann常数，T∞为辐射源温度。式(7)右侧的第一项代表对流传热，第二项代表辐射传热。
同样地，对于颗粒-壁面之间的热量传递可以表示为[27]
| $ q = \frac{{{k_{\rm{p}}}}}{{\Delta n}}\left( {{T_{\rm{w}}} - {T_{\rm{p}}}} \right) + {q_{{\rm{rad}}}}. $ | (8) |
其中：kp为颗粒的热导率，Δn为颗粒与壁面之间的距离，Tp为局部颗粒温度，Tw为壁面温度，qrad为辐射热通量。
颗粒-颗粒之间的热量传递可以采用线性模型表示，
| $ {Q_{ij}} = 2{k_{\rm{S}}}{\left[ {\frac{{3{F_N}{r^*}}}{{4{E^*}}}} \right]^{1/3}}\left( {{T_j} - {T_i}} \right). $ | (9) |
其中：Ti和Tj分别为颗粒i和j的温度，kS、r*和E*分别为颗粒材料的热导率、颗粒几何平均半径及颗粒i和j的有效Young模量。
根据Kirchhoff定律和Stefan-Boltzmann定律，颗粒之间的辐射热传递可以表示为[17]
| $ {Q_{\rm{R}}} = {\varepsilon _{\rm{p}}}{\rm{ \mathsf{ π} }}d_{\rm{p}}^2\left( {{\rm{ \mathsf{ π} }}I - \sigma {n^2}T_{\rm{p}}^4} \right). $ | (10) |
其中：QR为辐射热传递，εp为热辐射系数，dp为颗粒直径，I为颗粒表面的辐射强度，σ为Stefan-Boltzmann常数(约为5.67×10-8 W/(m2·K4))，n为流体的折射率，Tp为均匀温度。
1.1.3 浓度场流化床中化学反应和组分沉积会影响浓度场，一般使用组分传输方程(即对流扩散方程)来描述浓度场的变化，其形式如下[17]：
| $ \frac{\partial }{{\partial t}}\left( {{\rho _{\rm{g}}}{Y_i}} \right) + \nabla \cdot \left( {{\rho _{\rm{g}}}{\mathit{\boldsymbol{u}}_{\rm{g}}}{Y_i}} \right) = - \nabla \cdot {J_i} + {R_i}. $ | (11) |
其中：Yi为组分i局部质量分数，Ri表示由化学反应引起的组分i的净生成率，Ji为由于浓度和温度的梯度而引起的组分i的扩散通量。一般采用Fick定律来描述由于浓度梯度引起的质量扩散，组分i的扩散通量表示为
| $ {J_i} = - {\rho _{\rm{g}}}{D_{i, m}}\nabla {Y_i} - {D_{T, i}}\frac{{\nabla T}}{T}. $ | (12) |
其中：Di, m为混合物中组分i的质量扩散系数，DT, i为热扩散系数。
本文为了简化，一开始只研究流化床中通入乙炔高温裂解生成碳粉体沉积在颗粒表面，并假设该反应只有C(s)和H2产生，
| $ {{\rm{C}}_2}{{\rm{H}}_2}({\rm{g}}) \to 2{\rm{C}}({\rm{s}}) + {{\rm{H}}_2}({\rm{g}}) $ | (13) |
式(13)化学反应采用集总的反应动力学方程和参数，如下所示[28]：
| $ r = k\frac{{c_{{{\rm{C}}_2}{{\rm{H}}_2}}^{1.9}}}{{1 + 18{c_{{{\rm{H}}_2}}}}}. $ | (14) |
其中，反应常数k为
| $ k = A \cdot {{\rm{e}}^{\left( {\frac{{ - {E_A}}}{{RT}}} \right)}}. $ | (15) |
其中：R为摩尔气体常数(8.31 J/(mol·K))，T为热力学温度，A=5.5×106 mol·m3·s，EA=165.0 kJ/mol。
1.1.4 CVD模型颗粒包覆过程涉及晶体的生长机理和颗粒的生长过程，是非常复杂的。由于对沉积过程的细节缺少了解，基础模型中本文采用吹扫理论(即颗粒运动-黏附-沉积理论)来描述该过程，如图 4所示。在沉积过程中，假设颗粒表面的组分沉积是各向同性的，颗粒的包覆层非常薄以至于可以认为包覆颗粒的密度在沉积过程中是常数保持不变。
|
| 图 4 基于颗粒运动-黏附-沉积理论的颗粒包覆示意图 |
借鉴Czok提出的在颗粒表面包覆铝的沉积模型[18]，并结合颗粒运动-黏附-沉积理论，颗粒沉积速率正比于颗粒吹扫表面积、气体和固体相对速度差及周围的有效碳粉体浓度，沉积模型可以表示为
| $ \rho 4{\rm{ \mathsf{ π} }}r_{\rm{p}}^2\frac{{{\rm{d}}r}}{{{\rm{d}}t}} = {k_{\rm{a}}} \cdot {k_{\rm{d}}}{\rm{ \mathsf{ π} }}r_{\rm{p}}^2\left| {{\mathit{\boldsymbol{u}}_{\rm{g}}} - {\mathit{\boldsymbol{v}}_{\rm{p}}}} \right|{C_{{\rm{eff}}}}. $ | (16) |
其中：rp为颗粒半径，ka为黏附系数，kd为沉积系数，ug为气体速度，vp为颗粒速度，Ceff为颗粒周围的有效碳浓度。黏附系数ka和沉积系数kd与颗粒的材料、尺度、温度等因素有关，在初步研究中本文根据实验结果假设ka和kd都是常数，并且ka=1.00，kd=1.00×10-3。颗粒速度vp可以通过DEM框架求解，气体速度ug和颗粒周围的有效碳浓度Ceff可以通过CFD框架求解。
1.1.5 CFD-DEM-CVD多物理场多尺度耦合模型综合考虑速度场、温度场、浓度场和CVD模型，从而可以建立CFD-DEM-CVD多物理场多尺度耦合模型来研究流化床-化学气相沉积颗粒包覆过程，如图 5所示。在CFD模块中，基于质量守恒方程、动量守恒方程及能量守恒方程等可以求解气相运动和化学反应；在DEM模块中，基于曳力模型、接触力模型、Newton第二定律及转动方程等可以求解颗粒的运动；在CVD模块中，基于颗粒运动-黏附-沉积理论及晶体生长机制可以求解组分沉积过程。同时，不同组分的传输方程确保在整个模拟系统中保证质量、动量及能量守恒。
|
| 图 5 CFD-DEM-CVD多物理场多尺度耦合模型示意图 |
1.2 CRF模型 1.2.1 化学反应流机理
在流化床-化学气相沉积颗粒包覆过程中，一般通入乙炔进行碳的包覆。在1.1节提出的CFD-DEM-CVD模型中，假设乙炔在高温分解下只有C(s)和H2产生。显然，这种假设是不够精确的，乙炔在高温裂解中会产生C(s)、H2以及各种碳氢化合物[28-30]，如CH4、C2H4、C4H4及C6H6等，并且这些碳氢化合物会对流化床中的温度场、浓度场及组分传输和沉积产生影响。根据研究发现，乙炔的高温裂解主要涉及9种化学反应和7种产物(C2H2、C(s)、H2、CH4、C2H4、C4H4及C6H6)，如图 6所示。
|
| 图 6 乙炔高温裂解的主要反应路径及产物 |
乙炔裂解具体的每步反应和相关的指前因子Ai及表观活化能EA, i如表 1所示。根据Khan的研究[28]，每步反应的反应速率如式(17)所示。
| $ {\rm{反应速率}}:\left\{ {\begin{array}{*{20}{l}} {{\rm{R}}1:{r_1} = {k_1}\frac{{c_{{{\rm{C}}_2}{{\rm{H}}_2}}^{1.9}}}{{1 + 18{c_{{{\rm{H}}_2}}}}}, }\\ {{\rm{R}}2:{r_2} = {k_2} \cdot {c_{{{\rm{C}}_2}{{\rm{H}}_2}}} \cdot c_{{{\rm{H}}_2}}^{0.36}, }\\ {{\rm{R}}3:{r_3} = {k_3} \cdot c_{{{\rm{C}}_2}{{\rm{H}}_4}}^{0.5}, }\\ {{\rm{R}}4:{r_4} = {k_4} \cdot c_{{{\rm{C}}_2}{{\rm{H}}_2}}^{0.35} \cdot c_{{{\rm{H}}_2}}^{0.22}, }\\ {{\rm{R}}5:{r_5} = {k_5} \cdot c_{{\rm{C}}{{\rm{H}}_4}}^{0.21}, }\\ {{\rm{R}}6:{r_6} = {k_6} \cdot c_{{{\rm{C}}_2}{{\rm{H}}_2}}^{1.6}, }\\ {{\rm{R}}7:{r_7} = {k_7} \cdot c_{{{\rm{C}}_4}{{\rm{H}}_4}}^{0.75}, }\\ {{\rm{R}}8:{r_8} = {k_8} \cdot c_{{{\rm{C}}_2}{{\rm{H}}_2}}^{1.3} \cdot c_{{{\rm{C}}_4}{{\rm{H}}_4}}^{0.6}, }\\ {{\rm{R}}9:{r_9} = {k_9}\frac{{c_{{{\rm{C}}_6}{{\rm{H}}_6}}^{0.75}}}{{1 + 22{c_{{{\rm{H}}_2}}}}}.} \end{array}} \right. $ | (17) |
| 序号 | 反应方程式 | Ai/(s-1或m3·mol-1·s-1) | EA, i /(kJ·mol-1) |
| R1 | C2H2→2C(s)+H2 | 5.5×106 | 165.0 |
| R2 | C2H2+H2→C2H4 | 4.4×103 | 103.0 |
| R3 | C2H4→C2H2+H2 | 3.8×107 | 200.0 |
| R4 | C2H2+3H2→2CH4 | 1.4×105 | 150.0 |
| R5 | 2CH4→C2H2+3H2 | 8.6×106 | 195.0 |
| R6 | C2H2+C2H2→C4H4 | 1.2×105 | 120.7 |
| R7 | C4H4→C2H2+C2H2 | 1.0×1015 | 335.2 |
| R8 | C4H4+C2H2→C6H6 | 1.8×103 | 64.5 |
| R9 | C6H6→6C(s)+3H2 | 1.0×103 | 75.0 |
| 注：一级反应Ai单位为s-1，二级反应Ai单位为m3·mol-1·s-1。 | |||
其中，ki是每步反应的速率常数，可以通过Arrhenius公式得到：
| $ {k_i} = {A_i} \cdot {{\rm{e}}^{\left( {\frac{{ - {E_{A, i}}}}{{RT}}} \right)}}. $ | (18) |
其中：Ai为指前因子，EA, i为表观活化能，R为摩尔气体常数(8.31 J/(mol·K))，T为热力学温度。
根据基元反应动力学和Arrhenius公式，利用MATLAB编程和DETCHEM软件计算可以得到在923~1 323 K范围内乙炔裂解生成各产物的摩尔比。如图 7所示，模拟计算得到的各产物的摩尔比与Khan和Buchholz的实验结果[28, 30]相符合，说明了本文计算的可靠性。可以发现，随着温度升高，C2H2摩尔比下降，C(s)和H2的摩尔比上升，CH4、C2H4、C4H4及C6H6等碳氢化合物的摩尔比基本都小于5%。
|
| 图 7 乙炔高温裂解生成的各产物的摩尔比模拟结果和实验结果对比 |
1.2.2 组分沉积过程机理
如图 8所示，乙炔在高温区裂解产生C(s)、H2以及各种碳氢化合物，然后碳粉体沉积在基体颗粒表面。实际的碳粉体沉积过程非常复杂，包括碳粉体的沉积过程、黏附过程、脱附过程以及各种气体的清除效应，此外，反应产生的氢气会在碳层和沉积层表面形成非常稳固的碳氢络合物阻止碳粉体的沉积，即氢气的抑制效应，如图 9和10所示。
|
| 图 8 乙炔高温裂解和碳粉体沉积示意图 |
|
| 图 9 碳粉体沉积详细过程示意图 |
|
| 图 10 碳原子在基体表面沉积、黏附及脱附过程示意图 |
根据分子热力学理论，当体系温度为绝对零度0 K时，沉积系数应当为0，黏附系数应当为1，脱附系数应当为0；当沉积表面温度为无穷大时，黏附系数应当为0，脱附系数应当为1。结合本文研究，可以得到沉积系数公式[31]如下：
| $ {k_{\rm{d}}} = 0.8579 \times 2{\rm{ \mathsf{ π} }}{r_{\rm{p}}}T. $ | (19) |
其中：kd为沉积系数，rp为颗粒半径，T为沉积基体表面热力学温度。
黏附系数公式[32]如下：
| $ {k_{\rm{a}}} = \frac{1}{{1 + {{\rm{e}}^{\left( {\frac{T}{{1083}} - \frac{{3287.3}}{T}} \right)}}}}. $ | (20) |
其中, ka为黏附系数。
脱附系数公式如下：
| $ {\eta _{{\rm{des}}}} = \frac{r}{{r + {r_{\rm{p}}}}}{{\rm{e}}^{\frac{{ - {E_{\rm{d}}}}}{{4.3531RT}}}}. $ | (21) |
其中：ηdes为脱附系数，r为沉积层厚度，Ed为活化能，R为摩尔气体常数(8.31 J/(mol·K))。
图 11展示了乙炔裂解产物对碳粉体沉积的影响。根据研究发现，氢气抑制效应系数可以表示为
| $ {\eta _{{\rm{inh }}}} = {{\rm{e}}^{ - \sum\limits_{i = 1}^n {{I_i}} }}, \left( {i代表{\rm{C}}{{\rm{H}}_4}, {{\rm{C}}_2}{{\rm{H}}_2}, {{\rm{C}}_2}{{\rm{H}}_4}, {{\rm{C}}_6}{{\rm{H}}_6}} \right). $ | (22) |
|
| 图 11 乙炔裂解产物对碳粉体沉积的影响示意图 |
其中[33]，
| $ {I_{{\rm{C}}{{\rm{H}}_4}}} = \frac{{\frac{{{k_4}}}{{{k_5}}}}}{{\frac{{{k_4}}}{{{k_5}}} + \frac{{{c_{{{\rm{H}}_2}}}}}{{{c_{{\rm{C}}{{\rm{H}}_4}}}}}}}, {I_{{{\rm{C}}_2}{{\rm{H}}_2}}} = \frac{{\frac{{{k_7}}}{{{k_1}}}}}{{\frac{{{k_7}}}{{{k_1}}} + \frac{{{c_{{{\rm{H}}_2}}}}}{{{c_{{{\rm{C}}_2}{{\rm{H}}_2}}}}}}}, $ | (23) |
| $ {I_{{{\rm{C}}_2}{{\rm{H}}_4}}} = \frac{{\frac{{{k_2}}}{{{k_3}}}}}{{\frac{{{k_2}}}{{{k_3}}} + \frac{{{c_{{{\rm{H}}_2}}}}}{{{c_{{{\rm{C}}_2}{{\rm{H}}_4}}}}}}}, {I_{{{\rm{C}}_6}{{\rm{H}}_6}}} = \frac{{\frac{{{k_8}}}{{{k_9}}}}}{{\frac{{{k_8}}}{{{k_9}}} + \frac{{{c_{{{\rm{H}}_2}}}}}{{{c_{{{\rm{C}}_6}{{\rm{H}}_6}}}}}}}. $ | (24) |
其中：ki为速率常数，ci为气体i的浓度。
气体的清除效应系数可以表示为
| $ \begin{array}{c} {\eta _{{\rm{cle }}}} = \frac{r}{{r + {r_{\rm{p}}}}}{{\rm{e}}^{ - \sum\limits_{i = 1}^n {{\rm{C}}{{\rm{l}}_i}} }}, \\ {\rm{ ( }}i代表{{\rm{C}}_2}{{\rm{H}}_2}, {\rm{C}}{{\rm{H}}_4}, {{\rm{C}}_2}{{\rm{H}}_4}, {{\rm{C}}_4}{{\rm{H}}_4}, {{\rm{C}}_6}{{\rm{H}}_6}, {{\rm{H}}_2}, {\rm{Ar )}}{\rm{. }} \end{array} $ | (25) |
其中，
| $ \begin{array}{c} {\rm{C}}{{\rm{l}}_i} = \frac{1}{{1 + \frac{{{{\left| {{\mathit{\boldsymbol{u}}_{{{\rm{g}}_i}}} - {\mathit{\boldsymbol{v}}_{\rm{p}}}} \right|}^2}{c_{{{\rm{g}}_i}}}}}{{\left| {{\mathit{\boldsymbol{u}}_{\rm{g}}}} \right|\left| {{\mathit{\boldsymbol{v}}_{\rm{p}}}} \right|{c_{\rm{g}}}}}}}, \\ {\rm{ (i代表}}{{\rm{C}}_2}{{\rm{H}}_2}, {\rm{C}}{{\rm{H}}_4}, {{\rm{C}}_2}{{\rm{H}}_4}, {{\rm{C}}_4}{{\rm{H}}_4}, {{\rm{C}}_6}{{\rm{H}}_6}, {{\rm{H}}_2}, {\rm{Ar )}}{\rm{. }} \end{array} $ | (26) |
其中：ugi为气体i的速度，vp为颗粒速度，cgi为气体i的浓度，ug为入口气体速度，cg为入口气体浓度。
综合上面碳粉体沉积的详细过程，并结合颗粒运动-黏附-沉积理论，可以得到乙炔裂解-生产碳粉-气相沉积的模型公式如下：
| $ \begin{array}{c} 4{\rm{ \mathsf{ π} }}r_{\rm{p}}^2\rho \frac{{{\rm{d}}r}}{{{\rm{d}}t}} = {k_{\rm{d}}} \cdot {k_{\rm{a}}} \cdot \left( {1 - {\eta _{{\rm{des}}}}} \right) \cdot \left( {1 - {\eta _{{\rm{inh}}}}} \right) \cdot \\ \left( {1 - {\eta _{{\rm{cle}}}}} \right) \cdot 4{\rm{ \mathsf{ π} }}r_{\rm{p}}^2\left| {{\mathit{\boldsymbol{u}}_{{{\rm{g}}_i}}} - {\mathit{\boldsymbol{v}}_{\rm{p}}}} \right|{C_{{\rm{eff}}}}. \end{array} $ | (27) |
其中：rp、ugi、vp和Ceff分别为颗粒半径、气体i的速度、颗粒速度及颗粒周围的有效碳浓度。
综合乙炔裂解的主要反应路径和生成物，并考虑碳粉体沉积的详细过程机理，基于基元反应动力学、分子热力学和颗粒运动-黏附-沉积理论可以建立CRF模型，从而求解物质的分解和沉积过程。
1.3 PBM模型CRF模型中没有考虑乙炔裂解生成碳粉体的尺径大小分布，认为是均匀大小。显然，这在实际中是不准确的，实际生成的碳粉体尺径大小可能彼此不同，从而在基体颗粒上的沉积也不同。根据Doraiswami的研究[34]，纳米颗粒的尺径分布可以用PBM模型来描述，
| $ \begin{array}{c} \frac{{\partial n(L;\mathit{\boldsymbol{x}}, t)}}{{\partial t}} + \nabla \cdot [\mathit{\boldsymbol{u}}n(L;\mathit{\boldsymbol{x}}, t)] - \\ \nabla \cdot \left[ {{\mathit{\Gamma }_{{\rm{ef}}}}\nabla n(L;\mathit{\boldsymbol{x}}, t)} \right] = \\ - \frac{\partial }{{\partial t}}[G(L)n(L;\mathit{\boldsymbol{x}}, t)] + \\ B(L;\mathit{\boldsymbol{x}}, t) + D(L;\mathit{\boldsymbol{x}}, t). \end{array} $ | (28) |
其中：n(L; x, t)为数目密度函数，G(L)为生长速率，B(L; x, t)为由于团聚引起的生长率，D(L; x, t)为由于破损引起的减少率。
当确定成核、生长、团聚和破损模型后，可以用初始条件和边界条件求解该模型[35]，如图 12所示。生长速率G(L)与成核和生长动力学有关，是热解产物过饱和度的函数。为了简化研究，将纳米粒子种子的半径设置为5 nm，并将6个不同的纳米粒子半径分别设置为5、10、20、50、100和200 nm，以表示本模型中的所有纳米粒子。生长速率G(L)可以采用两步模型法[36]计算。在本研究中，采用了Luo模型[37]来模拟纳米颗粒团聚和破损。
|
| 图 12 计算组分颗粒尺径分布示意图 |
1.4 IBP模型
在之前的研究中，由于包覆层很薄，假设在包覆过程中包覆颗粒密度保持不变。但随着包覆层厚度的增加，显然这种假设是不准确的。为此，本文在模拟中构造了一个内嵌式虚拟颗粒(inner bonding particle)，赋予其质量m2和体积V2，并且该虚拟颗粒与流体没有相互作用，如图 13所示。
|
| 图 13 内嵌式颗粒示意图 |
包覆后颗粒质量和密度与内嵌式虚拟颗粒的关系为
| $ \left\{ {\begin{array}{*{20}{l}} {m = {m_1} + {m_2}, }\\ {V = {V_1}.} \end{array}} \right. $ | (29) |
得到包覆后颗粒的密度为
| $ \begin{array}{c} \rho = \frac{m}{V} = \frac{{{m_1} + {m_2}}}{{{V_1}}} = {\rho _1} + \frac{{{m_2}}}{{{V_1}}} = \\ {\rho _1} + \frac{{{\rho _2}{V_2}}}{{{V_1}}} = {\rho _1} + {\rho _2}\frac{{{V_2}}}{{{V_1}}}. \end{array} $ | (30) |
可以看出，在模拟过程中不改变ρ1和ρ2，只要改变V1(通过包覆模型改变)，然后调整V2(虚拟的内嵌式颗粒，和流体没有相互作用)，使得整个颗粒ρ变得和实际一样，从而实现包覆过程中颗粒密度的改变。由此建立了IBP模型，从而求解包覆颗粒密度的改变。
综上所述，在流化床-化学气相沉积颗粒包覆过程的研究中，通过对颗粒包覆过程机理的探究和模型的逐步优化，可以得到一个相对来说比较准确的流化床-化学气相沉积颗粒包覆模型。
2 几何结构和模拟参数 2.1 几何结构和参数模拟计算中本文使用的是单孔矩形槽喷动流化床，如图 14a所示。喷嘴的直径D0为9 mm，底部厚度Di为15 mm，倒锥角γ为60°，床宽Dc为150 mm，床厚L为15 mm。
|
| 图 14 单孔矩形槽喷动流化床结构 |
2.2 模拟参数设置
对于流化床-化学气相沉积颗粒包覆模型的计算，其CFD模块利用FluentTM软件计算，DEM模块利用EDEMTM软件计算，CVD模块、CRF模块、PBM模块及IBP模块通过Fluent软件中的自定义函数(UDF)模块编程实现。Fluent软件中的UDF和EDEM软件中的应用程序接口(API)将这些模块连接起来，从而求解流化床-化学气相沉积颗粒包覆模型。
由于单孔矩形槽喷动流化床底部结构复杂，模拟中构造了非结构性网格，如图 14b和14c所示。首先用单相流对网格独立性进行了测试，网格总数是12 969，使用的CFD网格面积最小为9.70×10-6 m2，最大为4.04×10-5 m2。在气体入口采用速度边界条件，气体出口采用自由压力边界条件，壁面采用无滑移边界条件，具体模拟参数设置如表 2所示。
| 参数 | 数值 | |
| 气相 | 密度(ρg) | 1.62 g/L |
| 黏度(μ) | 2.02×10-5 Pa·s | |
| 气速(vm) | 11 m/s | |
| 颗粒相 | 颗粒数目 | 15 000 |
| 颗粒半径(rp) | 1.0 mm | |
| 颗粒密度(ρp) | 2.38 g/cm3 | |
| 模拟参数 | Poisson比 | 0.25 |
| Young模量 | 108 N/m2 | |
| 恢复系数 | 0.90 | |
| 颗粒碰撞模型 | Hertz-Mindlin模型 | |
| 曳力模型 | Di Felice模型 | |
| 滑动摩擦系数 | 0.5 | |
| 滚动摩擦系数 | 10-5 | |
| 模拟时间 | 4.60 s | |
| CFD时间步长 | 10-5 s | |
| DEM时间步长 | 10-6 s |
在本文研究中，流化床入口通入乙炔、氩气和氢气，3种气体通入的比例如表 3所示。入口气速为11 m/s，温度为室温298 K；流化床壁面温度的变化如表 3所示。其中Ug表示入口气速，Ug, min表示最小入口喷动气速。实验中，当颗粒位于初始状态时，单孔矩形槽喷动流化床的最小表面喷动气速约为0.28 m/s，相当于最小入口喷动气速Ug, min为4.6 m/s。
| 案例 | 入口气速Ug/(m·s-1) | 壁面温度/K | 质量分数比(C2H2∶Ar∶H2) |
| 案例1 | 11.0 | 1 223 | 0.20∶0.80∶0 |
| 案例2 | 11.0 | 1 273 | 0.20∶0.80∶0 |
| 案例3 | 11.0 | 1 323 | 0.20∶0.80∶0 |
| 案例4 | 11.0 | 1 273 | 0.20∶0.79∶0.01 |
| 案例5 | 11.0 | 1 273 | 0.20∶0.40∶0.40 |
| 案例6 | 11.0 | 1 273 | 0.20∶0∶0.80 |
3 模拟结果及讨论
基于本文提出的包覆模型，将其用于流化床-化学气相沉积颗粒包覆过程模拟，并计算了表 3中的6个案例。下文将详细分析颗粒的包覆过程，并比较不同条件下颗粒的包覆速率、颗粒混合度及气体浓度分布等。
3.1 颗粒包覆过程以表 3中的案例2为例来分析颗粒的包覆过程。在模拟开始，未被包覆的颗粒被标记为蓝色，随着包覆的进行和越来越多的碳粉体沉积在颗粒表面，颗粒逐渐变为绿色，然后变为红色。如图 15所示，随着包覆时间的增加，越来越多的颗粒变为绿色和红色，说明越来越多的颗粒都被包覆上一层热解碳层。从图 15可以看出，2.0 s后基本上所有的颗粒都被包覆上一层热解碳层，4.5 s后大部分颗粒都被包上了一层约2.00×10-3 mg的热解碳层。
|
| 图 15 颗粒包覆过程 |
如图 16所示，随着包覆时间的增加，颗粒的包覆层厚度越来越大。在0.50 s时，大部分颗粒的包覆层厚度介于9.92×10-6~1.50×10-5 mm，并且最大的包覆层厚度为4.18×10-3 mm；在1.50 s时, 大部分颗粒的包覆层厚度介于1.98×10-5~3.62×10-5 mm, 并且最大的包覆层厚度为7.98×10-3 mm；在3.00 s时, 大部分颗粒的包覆层厚度介于3.97×10-5~5.53×10-5 mm，并且最大的包覆层厚度为8.91×10-3 mm；在4.50 s时, 大部分颗粒的包覆层厚度介于5.79×10-5~8.10×10-5 mm, 并且最大的包覆层厚度为9.40×10-3 mm。
|
| 图 16 不同时刻颗粒包覆层厚度的颗粒数目分布 |
如图 17所示，所有颗粒的平均包覆层质量和厚度随着包覆时间增加而增加，但不是随着时间线性增加的，因为脱附效应和气体清除效应会随着包覆层厚度增加而增加，导致颗粒包覆的速率会逐渐变缓慢，这与实际是相符合的[38]。在0~0.71 s, 0.71~2.94 s和2.94~4.61 s三段时间范围内，颗粒的平均包覆层质量增速分别为1.0×10-3、4.83×10-4及2.52×10-4 mg/s，颗粒的平均包覆层厚度增速分别为3.33×10-5、1.60×10-5及8.34×10-6 mm/s。
|
| 图 17 所有颗粒的包覆层平均质量和包覆层平均厚度随时间的变化 |
实验中，本文在三维喷动流化床中对半径为1.0 mm的ZrO2颗粒进行包覆(如图 18所示)，入口通入乙炔和氩气，氩气和乙炔的质量分数之比为0.80∶0.20，壁面温度为1 273 K，与模拟中采用同样的条件。在实验中，颗粒包覆150 s增重9.585 27×10-5 g，包覆速率为6.39×10-4 mg/s；在模拟中，颗粒包覆3.00 s增重1.819 3×10-3 mg，包覆速率为6.06×10-4 mg/s。实验和模拟结果的相对误差为
| $ \begin{array}{c} {\eta _{{\rm{d}}\mathit{r}}} = \frac{{\left| {{\rm{\Delta }}{r_{\rm{s}}} - {\rm{\Delta }}{r_{\rm{c}}}} \right|}}{{{\rm{\Delta }}{r_{\rm{c}}}}} \times 100{\rm{\% }} = \\ \frac{{\left| {6.06 \times {{10}^{ - 4}}{\rm{mg}}/{\rm{s}} - 6.39 \times {{10}^{ - 4}}{\rm{mg}}/{\rm{s}}} \right|}}{{6.39 \times {{10}^{ - 4}}{\rm{mg}}/{\rm{s}}}} \times 100{\rm{\% }} = \\ 5.16{\rm{\% }}. \end{array} $ | (31) |
|
| 图 18 流化床-化学气相沉积颗粒包覆实验中颗粒包覆前后对比 |
实验和模拟结果的相对误差为5.16%，说明模拟结果具有可靠性，也验证了本文提出的包覆模型的准确性和可靠性，可以用于流化床-化学气相沉积颗粒包覆过程的研究，同时也可以为实际工业中颗粒包覆过程和流化床设计提供指导，并可以在将来发展成为一种用于研究颗粒包覆过程的基础理论模型。
3.2 沉积组分尺径分布根据模拟结果可以获得乙炔裂解生成的碳粉体纳米颗粒的尺径分布，如图 19所示，给出了流化床出口处纳米颗粒的尺径分布。可以发现，在流化床入口仅形成了5 nm的纳米颗粒核，但随着颗粒生长和团聚，在流化床出口处纳米颗粒尺径有了显著增加。碳纳米颗粒尺径的分布也可以通过电子显微镜分析从实验结果中获得，如图 20所示。通过TEM图像分析可以发现乙炔裂解生成的碳纳米颗粒的尺径大小在5~20 nm，与模拟结果定性上一致。定量和详细的比较需要接下来进一步的参数测试和分析。
|
| 图 19 流化床出口处纳米颗粒的尺径分布 |
3.3 包覆颗粒密度
根据本文的模型计算，可以得到包覆颗粒平均密度随包覆时间的变化情况，如图 21所示。可以发现，包覆颗粒的平均密度随着包覆时间的增加而逐渐减小，由于包覆的热解炭层密度要小于基体颗粒密度，因此随着热解炭层增厚，包覆颗粒的密度逐渐减小。但由于包覆的热解炭层非常薄，因此包覆颗粒的密度的变化也是非常微小，基本可以忽略不计，CFD-DEM-CVD模型一开始的假设，即包覆颗粒密度保持不变，是可以接受的。但随着包覆的热解炭层厚度增加，包覆颗粒的密度变化是需要考虑的。
|
| 图 21 包覆颗粒平均密度随包覆时间的变化 |
3.4 包覆速率
根据本文模型，研究了壁面温度和通入氢气的量对颗粒包覆速率的影响。如图 22a所示，随着壁面温度上升，颗粒包覆速率增加。由于温度越高，乙炔裂解速率越快越充分，产生更多的碳粉体沉积在颗粒表面，因此颗粒包覆速率随着温度上升而增加。当壁面温度从1 223 K上升到1 323 K时，颗粒平均包覆速率从9.276×10-4 mg/s增加为9.955×10-4 mg/s。如图 22b所示，随着通入氢气的质量分数增加，颗粒包覆速率下降。由于通入的氢气一方面会抑制乙炔的裂解，另一方面会在碳层和沉积层表面形成非常稳固的碳氢络合物阻止碳粉体的沉积。因此，颗粒包覆速率随着通入氢气的质量分数增加而下降。当流化床入口通入乙炔的质量分数从0.01变为0.80时，颗粒平均包覆速率从3.109×10-4 mg/s降为1.847×10-4 mg/s。
|
| 图 22 颗粒包覆速率随不同温度和不同氢气质量分数的变化情况 |
3.5 颗粒混合
根据本文模型，也研究了壁面温度和通入氢气的量对流化床中颗粒混合的影响。采用的混合模型如下所示[22]：
| $ M = \frac{{\sum\limits_{i = 1}^n {\left[ {\sum\limits_{j = 1}^m {\left( {{r_{ij}} - d} \right)} } \right]} }}{{\sum\limits_{i = 1}^n {\left[ {\sum\limits_{k = 1}^m {\left( {{r_{ik}} - d} \right)} } \right]} }}. $ | (32) |
其中：M表示混合因子，rij表示颗粒i和与其接触的颗粒j的中心距离(即接触颗粒对ij之间的距离)，m表示颗粒i的配位数，rik表示颗粒i与其随机颗粒k的中心距离(即随机颗粒对ik之间的距离)，n表示流化床中总的颗粒数目。本文模型的物理意义为：用所有颗粒与其一开始选定的接触颗粒之间的距离比上所有颗粒与其一开始选定的随机颗粒之间的距离的比值来衡量体系中颗粒的混合程度。
图 23a显示了流化床中颗粒混合因子随不同温度的变化情况。可以看出在1 223、1 273及1 323 K 3种不同的温度下，颗粒混合因子随时间的变化曲线基本重合，说明当壁面温度变化很小时，流化床中颗粒的混合行为几乎不受壁面温度的影响。图 23b显示了流化床中颗粒混合因子随不同氢气质量分数的变化情况。可以看出，随着通入氢气的质量分数增加，流化床中颗粒混合开始明显变得缓慢。说明通入的氢气可以影响流化床中颗粒的流化、喷动及混合。
|
| 图 23 流化床中颗粒混合因子随温度和氢气质量分数的变化 |
3.6 压力场和速度场
图 24显示了1.90 s时案例2中气体动压力、气体体积分数及颗粒速度的分布。可以发现，在流化床入口处气体动压力最大，随着位置远离入口，气体动压力急剧减小。在流化床中形成了一个大的喷泉区，并且在喷泉区颗粒的包覆速度最快。相比于壁面附近，在入口附近和喷泉区附近的颗粒速度相对较大。
|
| 图 24 某一时刻(1.90 s)时的气体动压力、气体体积分数及颗粒速度分布 |
4 结论及展望
包覆型燃料颗粒由内部燃料核芯与外部包覆层组成，具有固有安全性，可应用于高温气冷堆等新型第四代反应堆，也可用于制备全陶瓷型燃料元件(FCM)等，其由流化床-化学气相沉积法制备而成。本文首先分析了流化床-化学气相沉积颗粒包覆过程的特点，然后综合流化床中的速度场、温度场、浓度场及CVD模型建立了初步的CFD-DEM-CVD多物理场多尺度耦合模型用于颗粒包覆过程的研究，并分别提出了CRF、PBM及IBP模型对颗粒包覆模型进行优化，最后将模型运用于流化床-化学气相沉积颗粒包覆过程的计算，并将模拟结果和实验结果对比。得到的主要结论有：
1) 建立的CRF模型可以用于模拟物质的分解和沉积过程，PBM模型可以求解沉积组分的尺径分布，IBP模型可以求解包覆颗粒密度的改变。
2) 当流化床壁面温度增加，颗粒包覆速率增加；当通入的氢气质量分数增加，颗粒包覆速率下降。
3) 当流化床壁面温度变化不大时，颗粒混合基本不受影响；通入氢气质量分数越大，颗粒流化、喷动、混合越缓慢。
4) 建立的包覆模型可以作为一种基础理论模型用于颗粒包覆过程的研究，并为实际的工业生产提供理论指导。
本文建立的包覆模型可以用于流化床-化学气相沉积颗粒包覆过程的研究，但仍有一些问题亟需解决：
1) 计算时间问题：本文模型可以用于一般简单流化床中颗粒包覆过程的模拟计算，但计算时间较长。为了实现对具有复杂几何结构的真实包覆炉中的大量颗粒的包覆行为进行模拟计算，目前的模型算法及计算平台需要进一步改进才能实现更快的计算。
2) 过程机理问题：目前流化床-化学气相沉积颗粒包覆过程中的一些过程机理尚没有完全探究清楚，为了进一步优化颗粒包覆模型、提高计算精度，需要将实验和模拟相结合，进一步探究颗粒包覆过程的机理，从而建立更加精确的包覆模型。
3) 工程放大问题：目前研究的流化床尺寸是实验室研究尺寸，从实验室研究到工业化生产需要解决逐步放大的问题，特定尺寸的流化床反应器中的多尺度耦合规律在放大尺寸下的迁移及变化规律是放大过程中需要特别重视的问题。
| [1] |
陈猛, 刘马林, 唐亚平, 等. 颗粒包覆过程的数值模拟方法比较研究[J]. 中国粉体技术, 2018, 24(1): 7-15. CHEN M, LIU M L, TANG Y P, et al. Comparative study of numerical simulation method for particle coating process[J]. China Powder Science and Technology, 2018, 24(1): 7-15. (in Chinese) |
| [2] |
LIU M L, LIU R Z, LIU B, et al. Preparation of the coated nuclear fuel particle using the fluidized bed-chemical vapor deposition (FB-CVD) method[J]. Procedia Engineering, 2015, 102: 1890-1895. DOI:10.1016/j.proeng.2015.01.328 |
| [3] |
SAWA K, UETA S. Research and development on HTGR fuel in the HTTR project[J]. Nuclear Engineering and Design, 2004, 233(1-3): 163-172. DOI:10.1016/j.nucengdes.2004.08.006 |
| [4] |
PANNALA S, DAW C S, FINNEY C E A, et al. Simulating the dynamics of spouted-bed nuclear fuel coaters[J]. Chemical Vapor Deposition, 2007, 13(9): 481-490. DOI:10.1002/cvde.200606562 |
| [5] |
刘荣正, 刘马林, 邵友林, 等. 流化床-化学气相沉积技术的应用及研究进展[J]. 化工进展, 2016, 35(5): 1263-1272. LIU R Z, LIU M L, SHAO Y L, et al. Application and research progress of fluidized bed-chemical vapor deposition technology[J]. Chemical Industry and Engineering Progress, 2016, 35(5): 1263-1272. (in Chinese) |
| [6] |
MANN U, RUBINOVITCH M, CROSBY E J. Characterization and analysis of continuous recycle systems[J]. AIChE Journal, 1979, 25(5): 873-882. DOI:10.1002/aic.690250516 |
| [7] |
SHELUKAR S, HO J, ZEGA J, et al. Identification and characterization of factors controlling tablet coating uniformity in a Wurster coating process[J]. Powder Technology, 2000, 110(1-2): 29-36. DOI:10.1016/S0032-5910(99)00265-X |
| [8] |
PANDEY P, KATAKDAUNDE M, TURTON R. Modeling weight variability in a pan coating process using Monte Carlo simulations[J]. AAPS PharmSciTech, 2006, 7(4): E2-E11. DOI:10.1208/pt070483 |
| [9] |
KUSHAARI K, PANDEY P, SONG Y X, et al. Monte Carlo simulations to determine coating uniformity in a Wurster fluidized bed coating process[J]. Powder Technology, 2006, 166(2): 81-90. DOI:10.1016/j.powtec.2006.05.001 |
| [10] |
MARONGA S J, WNUKOWSKI P. Modelling of the three-domain fluidized-bed particulate coating process[J]. Chemical Engineering Science, 1997, 52(17): 2915-2925. DOI:10.1016/S0009-2509(97)00112-7 |
| [11] |
KUMAR R, FREIREICH B, WASSGREN C. DEM-compartment-population balance model for particle coating in a horizontal rotating drum[J]. Chemical Engineering Science, 2015, 125: 144-157. DOI:10.1016/j.ces.2014.06.046 |
| [12] |
LI J F, FREIREICH B, WSAAGREN C, et al. A general compartment-based population balance model for particle coating and layered granulation[J]. AIChE Journal, 2012, 58(5): 1397-1408. DOI:10.1002/aic.12678 |
| [13] |
TOSCHKOFF G, JUST S, FUNKE A, et al. Spray models for discrete element simulations of particle coating processes[J]. Chemical Engineering Science, 2013, 101: 603-614. DOI:10.1016/j.ces.2013.06.051 |
| [14] |
BÖRNER M, BVCK A, TSOTSAS E. DEM-CFD investigation of particle residence time distribution in top-spray fluidised bed granulation[J]. Chemical Engineering Science, 2017, 161: 187-197. DOI:10.1016/j.ces.2016.12.020 |
| [15] |
HILTON J E, YING D Y, CLEARY P W. Modelling spray coating using a combined CFD-DEM and spherical harmonic formulation[J]. Chemical Engineering Science, 2013, 99: 141-160. DOI:10.1016/j.ces.2013.05.051 |
| [16] |
NOROUZI H R. Simulation of pellet coating in Wurster coaters[J]. International Journal of Pharmaceutics, 2020, 590: 119931. DOI:10.1016/j.ijpharm.2020.119931 |
| [17] |
LIU M L, CHEN M, LI T J, et al. CFD-DEM-CVD multi-physical field coupling model for simulating particle coating process in spout bed[J]. Particuology, 2019, 42: 67-78. DOI:10.1016/j.partic.2018.03.011 |
| [18] |
CZOK G, YE M, KUIPERS J A M H, et al. Modeling of chemical vapor deposition in a fluidized bed reactor based on discrete particle simulation[J]. International Journal of Chemical Reactor Engineering, 2005, 3(1): A57. |
| [19] |
FRIES L, ANTONYUK S, HEINRICH S, et al. Collision dynamics in fluidised bed granulators: A DEM-CFD study[J]. Chemical Engineering Science, 2013, 86: 108-123. DOI:10.1016/j.ces.2012.06.026 |
| [20] |
LIU M L, WEN Y Y, LIU R Z, et al. Investigation of fluidization behavior of high density particle in spouted bed using CFD-DEM coupling method[J]. Powder Technology, 2015, 280: 72-82. DOI:10.1016/j.powtec.2015.04.042 |
| [21] |
CHEN M, LIU M L, TANG Y P. Comparison of Euler-Euler and Euler-Lagrange approaches for simulating gas-solid flows in a multiple-spouted bed[J]. International Journal of Chemical Reactor Engineering, 2019, 17(7): 20180254. |
| [22] |
CHEN M, LIU M L, LI T J, et al. A novel mixing index and its application in particle mixing behavior study in multiple-spouted bed[J]. Powder Technology, 2018, 339: 167-181. DOI:10.1016/j.powtec.2018.08.036 |
| [23] |
GRYCZKA O, HEINRICH S, DEEN N G, et al. Characterization and CFD-modeling of the hydrodynamics of a prismatic spouted bed apparatus[J]. Chemical Engineering Science, 2009, 64(14): 3352-3375. DOI:10.1016/j.ces.2009.04.020 |
| [24] |
XU B H, YU A B. Numerical simulation of the gas-solid flow in a fluidized bed by combining discrete particle method with computational fluid dynamics[J]. Chemical Engineering Science, 1997, 52(16): 2785-2809. DOI:10.1016/S0009-2509(97)00081-X |
| [25] |
CHU K W, YU A B. Numerical simulation of complex particle-fluid flows[J]. Powder Technology, 2008, 179(3): 104-114. DOI:10.1016/j.powtec.2007.06.017 |
| [26] |
LI J, MASON D J. A computational investigation of transient heat transfer in pneumatic transport of granular particles[J]. Powder Technology, 2000, 112(3): 273-282. DOI:10.1016/S0032-5910(00)00302-8 |
| [27] |
CHAUDHURI B, MUZZIO F J, TOMASSONE M S. Modeling of heat transfer in granular flow in rotating vessels[J]. Chemical Engineering Science, 2006, 61(19): 6348-6360. DOI:10.1016/j.ces.2006.05.034 |
| [28] |
KHAN R U, BAJOHR S, GRAF F, et al. Modeling of acetylene pyrolysis under steel vacuum carburizing conditions in a tubular flow reactor[J]. Molecules, 2007, 12(3): 290-296. DOI:10.3390/12030290 |
| [29] |
KIEFER J H, VON DRASEK W A. The mechanism of the homogeneous pyrolysis of acetylene[J]. International Journal of Chemical Kinetics, 1990, 22(7): 747-786. DOI:10.1002/kin.550220710 |
| [30] |
BUCHHOLZ D, KHAN R U, BAJOHR S, et al. Computational fluid dynamics modeling of acetylene pyrolysis for vacuum carburizing of steel[J]. Industrial & Engineering Chemistry Research, 2010, 49(3): 1130-1137. |
| [31] |
ZHANG W G, HVTTINGER K J. CVD of SiC from Methyltrichlorosilane. Part Ⅱ: Composition of the gas phase and the deposit[J]. Chemical Vapor Deposition, 2001, 7(4): 173-181. DOI:10.1002/1521-3862(200107)7:4<173::AID-CVDE173>3.0.CO;2-X |
| [32] |
臧卫国. 一个新的黏附系数公式及其应用[C]//中国宇航学会结构强度与环境工程专委会航天空间环境工程信息网学术讨论会. 北京, 中国: 中国宇航学会, 2005: 181-186. ZANG W G. A new formula of adhesion coefficient and its application[C]//Academic Symposium on the Aerospace Environmental Engineering Information Network of the Special Committee of Structural Strength and Environmental Engineering of Chinese Astronautical Society. Beijing, China: CSA, 2005: 181-186. (in Chinese) |
| [33] |
LI A J, DEUTSCHMANN O. Transient modeling of chemical vapor infiltration of methane using multi-step reaction and deposition models[J]. Chemical Engineering Science, 2007, 62(18-20): 4976-4982. DOI:10.1016/j.ces.2007.01.069 |
| [34] |
RAMKRISHNA D. Population balances: Theory and applications to particulate systems in engineering[M]. San Diego, California: Academic Press, 2000.
|
| [35] |
KIM Y, KIM H U, SHIN Y, et al. Modeling of silicon nanoparticle formation in inductively coupled plasma using a modified collision frequency function[J]. Journal of Mechanical Science and Technology, 2014, 28(11): 4693-4703. DOI:10.1007/s12206-014-1036-z |
| [36] |
CHENG J C, YANG C, MAO Z S, et al. CFD modeling of nucleation, growth, aggregation, and breakage in continuous precipitation of barium sulfate in a stirred tank[J]. Industrial & Engineering Chemistry Research, 2009, 48(15): 6992-7003. |
| [37] |
LUO H, SVENDSEN H F. Theoretical model for drop and bubble breakup in turbulent dispersions[J]. AIChE Journal, 1996, 42(5): 1225-1233. DOI:10.1002/aic.690420505 |
| [38] |
FARNOUSH H, MOHANDESI J A, FATMEHSARI D H. Effect of particle size on the electrophoretic deposition of hydroxyapatite coatings: A kinetic study based on a statistical analysis[J]. International Journal of Applied Ceramic Technology, 2013, 10(1): 87-96. DOI:10.1111/j.1744-7402.2012.02818.x |
| [39] |
LIU M L, LIU B, SHAO Y L. The study on pyrolytic carbon powder in the coating process of fuel particle for high-temperature gas-cooled reactor[C]//18th International Conference on Nuclear Engineering. Xi'an, China: ASME, 2010: 545-549.
|