下料方式和流量分配对D4合成SiO2的影响
张寒1, 黄泽宇2, 张国君1, 欧阳葆华1, 吴学民1, 郑丽丽2    
1. 湖北菲利华石英玻璃股份有限公司, 荆州 434000;
2. 清华大学 航天航空学院, 北京 100084
摘要:燃烧器设计和流量控制对燃烧稳定性、SiO2生成、疏松体品质和关键组分分布有重要影响。该文通过对八甲基环四硅氧烷(OMCTS,D4)制备SiO2疏松体过程的数值研究,分析了不同下料方式,以及D4和氧气流量分配对SiO2生成效率和沉积面上关键组分分布的影响。结果表明:采用多点下料方式可以明显提高SiO2生成率和沉积均匀性;内外D4进口流速接近时,SiO2生成率较高,并且OH浓度有所降低;根据D4质量流量分配对氧气进口流量分配进行调整,能够提高SiO2生成率,但相应的OH浓度也有所增加;但在保持D4质量流量不变的情况下,采用多点下料,导致D4进口速度降低,一定程度上降低了沉积效率。采用多点下料时,可以考虑适当提高D4流量,控制火焰温度,有利于提高SiO2生成和沉积效率。
关键词化学气相沉积    反应流数值模拟    SiO2疏松体    
Effect of feeding mode and flow distribution on SiO2 synthesis from D4
ZHANG Han1, HUANG Zeyu2, ZHANG Guojun1, OUYANG Baohua1, WU Xuemin1, ZHENG Lili2    
1. Hubei Feilihua Quartz Glass Co., Ltd., Jingzhou 434000, China;
2. School of Aerospace Engineering, Tsinghua University, Beijing 100084, China
Abstract: Burner design and flow rate control strongly impact combustion stability, SiO2 generation, soot quality, and key flue gas components. This paper presents a numerical study of SiO2 soot synthesis by the CVD method with octamethylcyclotetra-siloxane (also called 'OMCTS' or 'D4') as the precursor. The study investigates the effects of various feeding modes and D4 and oxygen mass flow rates on the SiO2 generation efficiency and the distribution of key components on the surface. The SiO2 generation rate and deposition uniformity can be significantly improved by using the multi-point feeding mode which increases the SiO2 generation rate and reduces the OH concentration even with similar D4 inlet velocities. Proper control of the oxygen flow rate relative to the D4 mass flow rate improves the SiO2 generation rate, but also increases the OH concentration. The multi-point feeding mode with the same D4 mass flow rate can reduce the D4 inlet velocity with then reduces the deposition efficiency some. Multi-point feeding with the proper increased D4 flow rate and flame temperature control improves the SiO2 generation rate and deposition efficiency.
Key words: chemical vapor deposition    reactive flow numerical simulation    SiO2 soot    

高纯熔石英玻璃因良好的特性在高端光学系统、电子通信、航空航天等领域有广泛应用,如含有一定羟基(—OH)成分的石英玻璃可作为红外探测器材料。在实际工业中,化学气相沉积(chemical vapor deposition, CVD)工艺是制备高纯熔石英玻璃的常用方法,使用无毒环保的有机硅D4替代传统的SiCl4作为原料合成SiO2成为研究重点。根据制备流程,熔石英玻璃制备方法又主要分为直接法和间接法。直接法为通过CVD工艺直接生成液态SiO2液滴,沉积凝固后形成熔石英玻璃;间接法则先制备SiO2疏松体,再进行重熔再生制备玻璃。相比直接法,间接法的优点是可以制备低羟基(—OH)高纯熔石英玻璃[1],并依据需求调制不同OH含量玻璃。在制备疏松体工艺中,燃烧器的设计、D4下料方式及流量分配对燃烧稳定性、燃烧效率有至关重要的影响,并且直接影响SiO2颗粒的沉积均匀性,从而影响下一步合成石英玻璃的重熔再生。

Dobbin等[2]最早提出以硅烷为原料合成熔石英玻璃,用HMDS、OMCTS(D4)取代SiCl4合成熔石英。Briesen等[3]研究了不同原料浓度和成分对合成的SiO2颗粒表面积的影响。聂兰舰等[4]用含硅化合物成功制备无氯、高光谱透过率、高光学均匀性、结构稳定的石英玻璃。这些工作主要基于单点下料技术,而单点下料系统中,反应物形成的焰火宽度小,适用于小口径玻璃合成。为了实现大尺寸(口径)合成,多点下料技术近期获得了重视。张国君等[5]提出一种新型多点下料石英玻璃燃烧器,可以提高流速均匀性和燃烧效率,制备高产量、高质量的SiO2粉料。黄若杰等[6]提出了一种高速合成石英玻璃多点下料燃烧器,能够提高气流稳定性,相比单点下料有效提升产品品质。

上述研究主要针对D4合成石英玻璃或SiO2颗粒形成,以及多点下料燃烧器的设计,但对于合成疏松体工艺过程控制参数影响SiO2颗粒形成和疏松体合成缺乏系统研究,导致无法实现工艺过程的优化控制,难以获得粒径均匀、表面光滑的疏松体。为了实现稳定制备大尺寸SiO2疏松体,有必要针对多点下料技术,研究关键控制参数对工艺过程的影响。本文借助数值模拟方法,研究了以D4为原料制备SiO2疏松体时,下料方式(即单点下料和多点下料)以及流量分配对SiO2生成和关键组分分布的影响,为实际生产提供理论指导。

1 物理模型和数值仿真模型

图 1a为通过CVD方法制备SiO2疏松体的装置示意图。首先通过氢氧火焰预热炉膛,当基底达到预定温度,原料D4由载气氮气送入炉内。D4在氢氧火焰中快速反应,生成SiO2颗粒。SiO2颗粒在火焰射流作用下沉积在基板上,并随着基板旋转向下移动,逐层沉积形成疏松体。反应生成的废气通过抽气口排出。燃烧器设置中心和外环2处原料进口,氢气和氧气进口交替布置,如图 1b所示。

图 1 装置几何结构示意图

系统涉及D4多组分反应的湍流流动,假设:忽略湍流脉动对化学反应热释放的影响;制备疏松体过程近似为准稳态过程;所有组分为气态组分,不考虑颗粒生成的影响。质量、动量、能量和组分输运可以用以下Favre平均的N-S方程一般形式表示如下:

$ \frac{{\partial \bar \rho }}{{\partial t}} + \frac{\partial }{{\partial {x_j}}}\left( {\bar \rho {{\tilde u}_j}} \right) = 0, $ (1)
$ \begin{array}{*{20}{l}} {\frac{{\partial \left( {\bar \rho {{\tilde u}_i}} \right)}}{{\partial t}} + \frac{\partial }{{\partial {x_j}}}\left( {\bar \rho {{\tilde u}_j}{{\tilde u}_i}} \right) = - \frac{\partial }{{\partial {x_i}}}\left( {\bar p + \frac{2}{3}\bar \rho \tilde k{\delta _{ij}}} \right) + }\\ {\frac{\partial }{{\partial {x_j}}}\left( {\left( {\mu + {\mu _t}} \right)\left( {\frac{{\partial {{\tilde u}_i}}}{{\partial {x_j}}} + \frac{{\partial {{\tilde u}_j}}}{{\partial {x_i}}} - \frac{2}{3}\frac{{\partial {{\tilde u}_k}}}{{\partial {x_k}}}{\delta _{ij}}} \right)} \right) + \bar \rho {g_i}, } \end{array} $ (2)
$ \frac{{\partial \left( {\bar \rho {{\tilde Y}_i}} \right)}}{{\partial t}} + \frac{\partial }{{\partial {x_j}}}\left( {\bar \rho {{\tilde u}_j}{{\tilde Y}_i}} \right) = \frac{\partial }{{\partial {x_j}}}\left( {\left( {\bar \rho {D_i} + \frac{{{\mu _t}}}{{S{c_t}}}} \right)\frac{{\partial {{\tilde Y}_i}}}{{\partial {x_j}}}} \right) + {{\dot \omega }_i}, $ (3)
$ \begin{align} & \frac{\partial (\bar{\rho }\tilde{h})}{\partial t}+\frac{\partial }{\partial {{x}_{j}}}\left( \bar{\rho }{{{\tilde{u}}}_{j}}\tilde{h} \right)=\frac{\partial }{\partial {{x}_{j}}}\left( \left( \frac{{{\mu }_{t}}}{P{{r}_{t}}}+\frac{\lambda }{{{C}_{p}}} \right)\frac{\partial \tilde{h}}{\partial {{x}_{j}}} \right)+ \\ & \frac{\partial }{\partial {{x}_{j}}}\left( \sum\limits_{a}{{{h}_{a}}}\left( \bar{\rho }{{D}_{a}}-\frac{\lambda }{{{C}_{p}}} \right)\frac{\partial {{{\tilde{Y}}}_{a}}}{\partial {{x}_{j}}} \right)+\nabla \cdot {{{\bar{\dot{Q}}}}_{\text{rad }}}. \\ \end{align} $ (4)

描述湍流输运采用k-ε模型[7],湍动能k和湍动能耗散率ε由以下输运方程获得:

$ \begin{array}{l} \frac{{\partial (\bar \rho \tilde k)}}{{\partial t}} + \frac{\partial }{{\partial {x_j}}}\left( {\bar \rho {{\tilde u}_j}\tilde k} \right) = \frac{\partial }{{\partial {x_j}}}\left( {\frac{{{\mu _t}}}{{{\sigma _k}}}\frac{{\partial \tilde k}}{{\partial {x_j}}}} \right) + \\ {\mu _t}\left( {\frac{{\partial {{\tilde u}_i}}}{{\partial {x_j}}} + \frac{{\partial {{\tilde u}_j}}}{{\partial {x_i}}}} \right)\frac{{\partial {{\tilde u}_i}}}{{\partial {x_j}}} - {C_D}\bar \rho \tilde k\frac{{\partial {{\tilde u}_j}}}{{\partial {x_j}}} - \bar \rho \tilde \varepsilon , \end{array} $ (5)
$ \begin{array}{*{20}{c}} {\frac{{\partial (\bar \rho \tilde \varepsilon )}}{{\partial t}} + \frac{\partial }{{\partial {x_j}}}\left( {\bar \rho {{\tilde u}_j}\tilde \varepsilon } \right) = \frac{\partial }{{\partial {x_j}}}\left( {\frac{{{\mu _t}}}{{{\sigma _\varepsilon }}}\frac{{\partial \tilde \varepsilon }}{{\partial {x_j}}}} \right) + }\\ {{C_1}\frac{{\tilde \varepsilon }}{{\tilde k}}\left[ {{\mu _t}\left( {\frac{{\partial {{\tilde u}_i}}}{{\partial {x_j}}} + \frac{{\partial {{\tilde u}_j}}}{{\partial {x_i}}}} \right)\frac{{\partial {{\tilde u}_i}}}{{\partial {x_j}}} - {C_D}\bar \rho \tilde k\frac{{\partial {{\tilde u}_j}}}{{\partial {x_j}}}} \right] - {C_2}\bar \rho \frac{{{{\tilde \varepsilon }^2}}}{{\tilde k}}.} \end{array} $ (6)

组分输运和能源方程中的化学反应源项通过有限速率模型计算,化学反应机理采用24组分、39步基元反应的简化反应机理[8],能量方程中的辐射源项通过离散坐标(DO)模型[9]计算。详细参数和模型介绍可见文[10-11]。

为了重点研究燃烧器下料方式和流量分配对火焰射流中组分输运的影响,选取气相反应区作为计算域,将计算模型简化为二维旋转轴对称模型,如图 1所示,疏松体界面采用壁面边界条件。所有壁面的流动边界条件均为无滑移边界,温度边界条件设置为有厚度壁面对流换热,对流换热系数为8 W/(m2·K)。炉子底部和抽气出口设置为压力出口,表压分别为0和-25 Pa。燃烧器进口为质量流量进口。原料进口温度设置为453 K,为D4的蒸发温度,氢氧进口温度设置为303 K,出口环境温度设置为303 K。输运方程使用商用软件Fluent计算。

2 结果与讨论 2.1 D4下料方式和流量分配的影响

为了比较D4单点下料和多点下料,以及多点下料时不同流量分配对火焰和SiO2生成的影响,计算模拟过程中氢氧摩尔比约为2.4,保持D4/载气总质量流量不变,调整原料进口质量流量分配,工况设置见表 1,其中D4占总质量流量的54%,其余为载气氮气。Case 1为单点下料,Case 2—4为多点下料,D4质量流量逐渐由中心进口向外围进口转移。其中,Case 2中内外D4进口的质量流量相同,Case 4中内外D4进口的速度相同。由于本次研究中没有考虑颗粒生成的影响,只考虑气相反应,而计算过程中火焰下游由于氧气不足,存在SiO与SiO2相互转化的现象,因此在分析组分分布时将含Si氧化物作为整体进行分析比较。

表 1 调节D4质量流量工况设置
CaseD4/氮气进口1D4/氮气进口2备注
kg·s-1kg·s-1
Case 10.000 280单点下料
Case 20.000 140.000 14相同质量流量
Case 30.000 060.000 22
Case 40.000 030.000 25相同进口速度

通过火焰射流的速度分布云图(见图 2)可以发现,受湍流作用以及各进口射流相互影响,火焰脉动较大。射流速度、质量流量搭配、进口间距对火焰稳定性有较大影响。如图 3火焰温度云图对比结果所示,相比单点下料,采用多点下料后火焰温度明显升高。这主要是由于采用单点下料时,射流中心D4与载气流量大,射流中心流速大,如图 2a所示,D4分解后的含碳组分与氢氧火焰混合不充分,导致反应不充分。保持D4总质量流量不变,采用多点下料,原料进口流速明显下降,受周围氢氧射流的卷吸混合作用,D4热解后的含碳组分反应更加充分,火焰温度更高。随着中心进口D4质量流量减少,外围进口流量增加,外围局部氧气不足,火焰下游温度略有下降。

图 2 (网络版彩图)不同D4流量分配下火焰射流速度云图

图 3 (网络版彩图)不同D4流量分配下火焰温度云图

OH浓度和均匀性对高纯熔石英玻璃光学特性有重要影响,在疏松体制备过程中也需要关注并控制火焰中和沉积面附近OH浓度。与温度分布类似,采用多点下料后火焰中的OH浓度要明显高于单点下料的情况,如图 4所示。单点下料时火焰下游温度低,OH进一步反应生成H2O;而多点下料时火焰温度高,氢氧火焰中本身的OH浓度增加。另外D4分解后的烃类产物反应更充分,一定程度上增加了OH浓度。但随着D4流量由中心进口向外围进口转移,由于火焰温度有所降低,使得火焰中OH浓度也有所下降。通过沉积面上OH浓度分布(见图 5),也可以看出Case 2,即D4内外进口等质量流率的情况下,沉积面上OH浓度最高,而其对应的火焰温度也最高,如图 3b所示。Case 1中沉积面上OH浓度最低,而Case 3和4中OH浓度比较接近。

图 4 (网络版彩图)不同D4流量分配下火焰中OH浓度分布云图

图 5 (网络版彩图)不同D4流量分配下沉积面上OH浓度分布

通过火焰中的含Si氧化物浓度分布(见图 6)可以发现,随着中心进口D4流量下降,外围D4流量增加,射流中心含Si氧化物浓度明显下降,外围含Si氧化物浓度逐渐增加,含Si氧化物浓度分布逐渐变宽。通过沉积面上含Si氧化物浓度分布(见图 7)可以发现,采用多点下料后,Case 2中沉积面上含Si氧化物浓度要高于Case 1单点下料的情况。随着D4流量向外围转移,沉积面上含Si氧化物浓度分布更加均匀。持续增加外围D4流量(Case 4),沉积面中心含Si氧化物浓度有所下降,但周围含Si氧化物几乎保持不变。这一定程度上是由于外围原料进口横截面积大于中心进口横截面积,因此随着中心原料进口D4流量向外围转移,中心射流速度明显降低,但外围射流速度增加不明显,沉积效率有所下降。并且,由于D4生成SiO2的主要反应是放热反应,而原料进口流速低,D4反应过快,会导致D4核心反应区温度偏高,一定程度上限制了SiO2的生成。

图 6 (网络版彩图)不同D4流量分配下含Si氧化物浓度分布云图

图 7 (网络版彩图)不同D4流量分配下沉积面上含Si氧化物浓度分布

为了比较不同工况对整体SiO2生成效率的影响,定义:

$ {\rm{Si}}{{\rm{O}}_2}生成率 = \frac{{{\rm{Si}}{{\rm{O}}_2}总反应生成速率}}{{{\rm{D}}4总质量流量}} \ \ . $ (7)

图 8为不同D4流量分配下,SiO2生成率的比较。可以发现,采用多点下料,D4反应更充分,SiO2生成率明显高于单点下料的情况。但随着D4质量流量进一步由中心进口向外围进口转移,SiO2生成率上升缓慢,并最终有所下降。这主要是由于本研究中只改变了D4质量流量的分配,但没有改变氧气流量,随着外围D4质量流量增加,局部氧气略显不足,导致外围SiO2生成受到限制。

图 8 不同D4流量分配下SiO2生成率变化

图 9是实际生产试验中采用单点下料和多点下料时疏松体表面实物照片。可以看出,采用多点下料后,沉积面SiO2疏松体颗粒沉积均匀性明显改善。采用单点下料,生成的SiO2颗粒集中在火焰射流中心,容易发生烧结团聚现象,随着射流到达沉积面后,由于粒径变化大,导致沉积面不均匀,气泡增多。疏松体颗粒粒径不均匀和夹杂气泡不利于后续的脱羟和重熔再生,影响纯度,并且导致玻璃体出现气泡等缺陷。采用多点下料后,颗粒粒径均匀性明显提高,夹杂气泡现象也明显减少。

图 9 采用不同下料方式所得到的SiO2疏松体沉积面实物照片

研究结果表明,采用多点下料可以有效提高SiO2生成率和沉积面颗粒均匀性。当内外D4进口速度接近时,即根据进口横截面积分配流量,可以获得较大的SiO2生成率,同时OH浓度也得到控制。但在采用多点下料时,应当注意保证D4进口流速,控制反应核心区温度,同时需要根据D4流量调节氧气的供给。

2.2 氧气流量分配的影响

在调节D4流量分配后,氧气和氢气的流量、流速应该与D4的流量、流速相匹配。这部分研究以Case 4为基础,由于Case 4中内外D4进口速度几乎相同,但外围原料进口D4质量流量较大,为了进一步优化流量配置,提高SiO2生成率,这里保持氧气总质量流量不变,减少氧气进口1(中心D4进口外侧)的流量,增加氧气进口2(外围D4进口外侧)的流量,研究在多点下料条件下,不同氧气流量分配对SiO2生成的影响。相关工况设置见表 2

表 2 调节氧气质量流量工况设置
Case氧气进口1/(kg·s-1)氧气进口2/(kg·s-1)
Case 40.002 650.002 61
Case 50.002 570.002 69
Case 60.002 490.002 77
Case 70.002 410.002 85

通过沉积面上OH浓度分布(见图 10)可以发现,相比Case 4中氧气进口1和2质量流量和速度接近的情况,氧气质量流量向外围转移后,沉积面上OH浓度有不同程度的上升。这是因为氧气流量调整后,反应更加充分。随着氧气流量向外围转移,沉积面上OH浓度呈现先增加再减少的变化趋势。Case 6中沉积面上OH浓度最高,反应相对更加充分。Case 5中沉积面上OH浓度由中心向四周下降更加陡峭,而随着外围氧气浓度增加,Case 7中沉积面上OH浓度变化比较平缓。

图 10 (网络版彩图)不同氧气流量分配下沉积面OH浓度分布

随着氧气流量向外围转移,火焰中含Si氧化物浓度分布变化不大,如图 11所示。Case 6中沉积面上含Si氧化物浓度略低,其余变化不是很明显。在当前D4流量分配下,氧气流量分配变化对沉积面含Si氧化物浓度分布影响不大。

图 11 (网络版彩图)不同氧气流量分配下沉积面上含Si氧化物浓度分布

通过比较整体SiO2生成率(见图 12),可以发现,随着内圈氧气进口流量向外围转移,SiO2生成率有所增加,但继续增加外围氧气流量,射流中心可能出现局部氧气不足,在Case 7中SiO2生成率有所下降。对于SiO2生成率而言,Case 6附近氧气流量分配已接近该D4流量分配下的最优值。但Case 6中沉积面上OH浓度相对较高,含Si氧化物反而略低,原因是此时外围D4反应比较分散,D4分解后的低碳产物反应更充分,使得下游OH浓度相对较高,而SiO2沉积效率则有所下降。

图 12 不同氧气流量分配下SiO2生成率变化

通过研究可知,根据D4质量流量分配适当调整氧气进口流量分配,能够提高SiO2生成率,但相应的OH浓度也有所增加。对SiO2沉积效率、火焰稳定性的进一步优化还需要综合考虑流量、流速以及化学反应热释放的影响。

3 结论

本文通过数值模拟方法,分析了以D4为原料制备SiO2疏松体时,不同下料方式以及D4和氧气流量分配对SiO2生成和关键组分分布的影响。

研究表明,采用多点下料方式可以明显提高SiO2生成率和沉积均匀性。内外D4进口根据进口横截面积分配流量,可以在获得较高SiO2生成率的同时一定程度上控制OH浓度。根据D4质量流量分配适当调整氧气进口流量分配,能够提高SiO2生成率,但相应的OH浓度也有所增加。在保持D4质量流量不变的情况下,采用多点下料会导致D4进口速度降低,降低沉积效率。所以,相比单点下料,采用多点下料时可以考虑适当提高D4流量,控制火焰温度,有利于提高SiO2生成率和沉积效率,并且降低OH浓度。

总之,对D4制备SiO2疏松体多点下料工艺进一步优化,需要综合考虑燃烧器结构如管道直径、间距等对射流稳定性、燃料和氧化剂扩散混合情况的影响,参数控制中流量、流速和化学反应热释放对SiO2生成率和沉积效率以及火焰稳定性的影响。

参考文献
[1]
聂兰舰, 宋学富, 向在奎, 等.一种间接合成石英玻璃的方法及其专用设备以及一种石英玻璃: 102583977A. 2012-07-18.
NIE L J, SONG X F, XIANG Z K, et al. A method of indirect synthesis of silica glass and its special equipment and a kind of silica glass: 102583977A. 2012-07-18. (in Chinese)
[2]
DOBBINS M, MCLAY R E. Method of making fused silica by decomposing siloxanes: US5043002. 1991-08-27.
[3]
BRIESEN H, FUHRMANN A, PRATSINIS S E. The effect of precursor in flame synthesis of SiO2[J]. Chemical Engineering Science, 1998, 53(24): 4105-4112. DOI:10.1016/S0009-2509(98)00219-X
[4]
聂兰舰, 王玉芬, 向在奎, 等.石英玻璃的制备方法: 106277750A. 2017-01-04.
NIE L J, WANG Y F, XIANG Z K, et al. Method of making silica glass: 106277750A. 2017-01-04. (in Chinese)
[5]
张国君, 李得军, 王巍, 等.一种石英玻璃燃烧器: 206799423U. 2017-12-26.
ZHANG G J, LI D J, WANG W, et al. A silica glass burner: 206799423U. 2017-12-26. (in Chinese)
[6]
黄若杰, 王巍, 刘志龙, 等.高速合成石英玻璃生产用多点下料燃烧器: 207760219U. 2018-08-24.
HUANG R J, WANG W, LIU Z L, et al. Multi-point feeding burner for high-speed synthetic silica glass production: 207760219U. 2018-08-24. (in Chinese)
[7]
LAUNDER B E, SPALDING D B. Lectures in mathermatical models of turbulence[M]. London: Academic Press, 1972.
[8]
CHEN J, HUANG Z Y, ZHENG L L. Detailed OMCTS/H2/O2 chemical kinetics mechanism and its application in the fused silica glass synthesis[J]. Journal of Non-Crystalline Solids, 2018, 481: 33-40. DOI:10.1016/j.jnoncrysol.2017.10.014
[9]
CHUI E H, RAITHBY G D. Computation of radiant heat transfer on a non-orthogonal mesh using the finite-volume method[J]. Numerical Heat Transfer, Part B:Fundamentals, 1993, 23(3): 269-288. DOI:10.1080/10407799308914901
[10]
HUANG Z Y, HUANG Y S, ZHENG L L, et al. Effects of precursors on silica particle generation in CVD synthesis for fused silica glass[J]. Journal of Non-Crystalline Solids, 2018, 499: 86-94. DOI:10.1016/j.jnoncrysol.2018.07.020
[11]
黄泽宇, 黄耀松, 陈靖, 等. 入口条件对D4合成石英玻璃的影响[J]. 燃烧科学与技术, 2018, 24(6): 477-486.
HUANG Z Y, HUANG Y S, CHEN J, et al. Effects of entrance conditions on fused silica glass synthesis using D4[J]. Journal of Combustion Science and Technology, 2018, 24(6): 477-486. (in Chinese)