随着国家对环境保护要求的日趋严格和装置设备科学合理设计的要求，从脱硫装置含H2S酸气中回收硫磺越来越受到重视[1-3]。国内天然气净化厂建造了许多硫磺回收装置用于处理含H2S的酸性气体，其中Claus工艺应用最为广泛。其原理是使1/3(φ)的H2S燃烧并生成SO2，再将生成的SO2与剩余的H2S经Claus催化反应生成硫磺和水。2004年，全球以Claus工艺生产的硫磺量已超过4.6×107 t。因此，提高Claus硫磺回收装置的硫回收率对于减少环境污染具有非常重要的意义[4-6]。
面临环保要求日益严格的现状，硫磺回收装置在天然气加工企业中的地位日益提升。在硫磺回收装置的主要设备中，燃烧炉、再热炉、反应器、换热器等的可靠性、稳定性及先进性已成为制约装置长周期安全运行的关键因素[7-8]。在生产过程中，因工况变化等原因，可能造成再热炉局部超温，从而带来安全问题[9-12]。燃料气在燃烧时，若炉内温度分布不均匀，易造成局部超温现象，烧毁炉体。为解决此问题，工程上采用引入保护风的方式降低炉内局部温度，防止再热炉内壁超温，从而保证炉体的安全。因此，有必要对再热炉内的流动燃烧过程开展研究。但由于实验工作量和测试等原因，要进行实验研究较为困难[13]。因此，以生产中的再热炉为实物模型，采用FLUENT软件对其进行了数值模拟研究，研究结果具有重要的工程指导意义。
图 1为天然气生产过程中使用的H-1402(1-6)-H01型再热炉计算模型，其总长度约10 m，炉体前部半径为0.7 m，后部半径为1 m。其工艺过程为：在前部燃烧器部分，空气从圆筒周围均匀流入，经导流板后进入燃烧室；中间为燃料气入口，燃料气进入燃烧器后与空气混合，在燃烧室中燃烧。燃烧室前的小孔为保护风进口，均匀分布在燃烧器的圆筒周围, 其作用是保护燃烧器和附近炉壁面，防止高温的产生。过程气从圆筒周围均匀流入，经炉壁上的开孔流入炉体内部。在炉体中间部分，高温的燃烧气体与过程气发生混合，过程气被加热。
为使研究结果对工程有指导意义，对再热炉实体进行简化后，建立了相应的计算模型。
对于燃烧炉内气体的流动混合问题，采用RNG k-ε双方程模型来模拟，其连续方程、动量方程、能量方程及k-ε方程的通用形式见式(1)。
式中, Φ为自变量；ρ为密度，kg/m3；ГΦ为扩散系数；SΦ为相应的源项。各自变量Φ及ГΦ、SΦ列于表 1。
式中，x，y，z为直角坐标系的3个方向，V为速度，m/s；p为压力，Pa；μeff为有效黏度，Pa·s；H为焓，kJ；λ为导热系数，W/(m·K)；c为比热，kJ/K；λ为容积发热率，kJ/m3；经验常数可取C1=1.5~2.2，C2=0.4~0.5，σε, eff=1.2。
其中，Gk的计算见式(2)。
式中，μeff=μt+μl，μt=Cμρk2/ε；μl为分子动力黏度，Pa·s；Cμ=0.09；ρ为密度，kg/m3；k为紊流动能，kJ；ε为紊流耗散率，kJ/s。
采用涡旋(Electrical Board Description, 以下简称EBD)和概率密度函数(Probability Density Function，以下简称PDF)模型模拟燃烧过程。
辐射过程求解采用P-1辐射模型[14], 对于辐射热流qr可得：
式中，qr为辐射热流，W/m2；a为吸收系数；G为入射辐射能量密度，W/m2；σ为斯蒂芬-波尔兹曼常数，σ=5.67×10-8W/(m2·K)；T为温度，K。
式(3)可直接代入能量方程，从而得到由于辐射所引起的热量源或汇。
几何模型建立后，利用Gambit软件对其进行网格划分。由于模型的几何形状复杂，计算区域较大，为能最大限度地减少数值误差和控制计算区域网格划分单元数，加快数值计算的收敛速度，需要将整个模型分割成若干个子模块。如图 2所示，计算区域分成了燃烧器、炉体前部、过程气流入部分、炉体后部等5个子区域及保护风入口。对模型中大部分子模块使用了六面体单元划分网格，网格尺寸为20 mm。燃烧器和过程气流入部分区域结构较为复杂，使用了四面体单元划分网格。对于燃料气流道较小部分，网格尺寸为2 mm。在进行区域离散计算的时候，采用了容积拆分的方法，将整个计算区域分离成9个模块，分别为燃烧器主体、燃料气细流道、燃烧器与炉体接口、过程气流入前段、过程气流入中段、过程气流入后段、过渡段、燃烧器后部、出口部分，网格尺寸从2~20 mm不等，各模块之间分别用界面进行连接。由于进行再热炉实体模型计算，网格数目较多，计算占用CPU时间非常长。在对计算精度和计算时间进行综合考虑后，通过收敛性计算表明，计算网格约1 000万时具有较好的收敛性和计算精度。通过炉壁散热与燃料燃烧释放热量相比要小得多，可以忽略不计。因此，在模拟计算中，壁面可以按绝热考虑。
模型的边界条件：设置天然气管道进口截面为燃气质量进口，空气管道进口截面为空气质量进口，过程气以及保护气进口截面为质量进口，燃烧室出口截面为质量出口，燃烧室的壁面设置成绝热壁面边界条件。根据现场运行参数，具体计算参数设置见表 2，过程气和燃料气的组成见表 3和表 4。
对于再热炉炉内复杂的三维湍流对流扩散、组分混合及燃烧、辐射的流动传热耦合问题，采用Baliga & Patankar提出的控制容积有限元法(Control Volume based Finite Element Method, 简称CVFEM)进行求解。控制方程的对流项采用Rhie和Chow提出的二阶上风格式进行离散，扩散项采用线性插值多项式进行计算，速度与压力的耦合采用SIMPLE算法处理。由于模拟研究的再热炉内速度场、浓度场以及温度场相互耦合，采用欠松弛因子进行计算。根据需要，速度和压力的欠松弛因子取为0.1~0.2。
为全面了解保护风的作用，针对过程气流量最大的工况，分别模拟了保护风流量为0、2 000 m3/h、4 000 m3/h和6 000 m3/h时的炉内燃烧工况。在不同的保护风量条件下，再热炉炉内的速度场分布是不同的。炉内速度场分布的模拟结果如图 3(无保护风时)和图 4(保护风量为6 000 m3/h时)所示(图 3~图 4中速度单位均为m/s)。
从图 3中可以看出，燃料气在炉内轴向3.5 m范围内燃烧，无保护风作用。图 4表明保护风可以一直流动到壁面过程气入口处。在不同流量下，保护风具有不同的流动范围，这必将影响炉内温度场分布，从而影响壁面温度分布。
图 5为无保护风时，再热炉炉内的温度场分布。图 6为保护风量为6 000 m3/h时，再热炉炉内的温度场分布(图 5、图 6中温度单位均为K)。根据流场分析可知，保护风随着流量的不同，其向前流动的范围不同，直接影响了炉体壁面的温度场分布。明显可以看出有保护风时，炉体内的温度场分布受到了显著的影响。无保护风时，图 5表明从燃烧器喷出口处到炉体约3.5 m处，壁面处温度较高。图 6表明保护风流量为6 000 m3/h时，保护风从壁面圆周处流入，直接降低了流入口附近的壁面温度，且保护风可以保持一直向前流动，与壁面上过程气入口处的过程气混合，全面保护了再热炉前端高温炉体部分。
前述计算假设壁面绝热，在分析炉体内部温度和壁面温度时，可以用靠近壁面的气体温度作为炉体壁面的温度来分析。图 7为不同保护风流量下的壁面温度分布图。
无保护风时，炉墙内壁面温度最高达到1 500 K，入口温度在0.75 m处就达到1 100 K。显然，在实际生产过程中，由于工况变化等原因，炉体内部完全可能出现局部超温现象，造成衬里耐火材料的损坏。当保护风流量为2 000 m3/h时，炉墙内壁面温度最高为1 100 K左右，入口温度在1.5 m处迅速升高到1 100 K。当保护风流量为6 000 m3/h时，炉墙内壁面温度最高为750 K左右，入口温度在2.2 m左右才升高至不到800 K。当保护风流量为6 000 m3/h时，由于流量足够大，可沿着壁面一直“冲”到过程气入口处，对壁面的保护范围更广。分别取过程气流量为最小、中间、正常时所对应的最佳保护风流量，拟合得到保护风流量随过程气流量变化的函数关系：
式中，VB为保护风流量，m3/h；VG为过程气总流量，m3/h。
采用FLUENT软件对Claus硫磺回收过程中再热炉内流动燃烧过程及保护风作用进行了数值模拟研究，得出以下结论：
(1) 针对H-1402(1-6)-H01型再热炉进行的数值模拟研究，对天然气工业其他大型燃烧装置的研究有一定的指导作用。
(2) 随着保护风量的增加，炉墙内壁面温度有较大程度的下降。在模拟研究条件下，过程气流量最大时，保护风流量为6 000 m3/h，炉墙内壁面温度最高从无保护风时的1 500 K下降为低于800 K。保护风在有效保护壁面的同时，对工艺生产没有造成影响。
(3) 通过拟合，得到保护风流量随过程气流量变化的函数关系，对实际生产有一定的指导意义。
2014, Vol. 43 Issue (6): 590-594
2014, Vol. 43 Issue (6): 590-594