基于埃特金插值的蒸气云爆炸规律研究

郭晓伟，张园园

（1.永城永安矿山安全技术工程有限公司，河南永城 476600；2.辽宁石油化工大学环境与安全工程学院，辽宁抚顺 113001）

当石油化工企业含有CO 或者H2等易燃易爆气体泄漏到空气中，没有立刻遇到点火源，CO 或者H2就会在自由空间与空气混合，形成爆炸性混合物。这时，系统中一旦出现点火源就会发生无约束的蒸气云爆炸。爆炸导致体系温度和压力急剧升高，给石油化工装置造成严重的破坏。国内外学者从不同的角度，主要通过实验分析、模拟分析、事故后果理论分析3 种方法对蒸气云爆炸的危害进行了研究。Thomas[1]进行了无约束的蒸气云爆炸实验，Lind和Whiston[2]对Thomas做的实验进行了验证，证明其实验结果具有可行性。毕明树[3]对一定浓度下球状蒸气云的爆炸后果进行了实验研究。罗正鸿等[4]对一定浓度下圆柱形蒸气云的爆炸后果进行了实验研究。孙博等[5]对蒸气云爆炸后果进行了研究。吕涛[6]进行了内部储罐失效研究。蒸气云爆炸数值模拟主要是以CFD[7-8]为基础。使用的软件主要有CFX，FLACS，Autoreagas，EXSIM等。蒸气云爆炸事故后果的理论分析主要是TNT当量法[9-10]，TNO 模型法[11-12]。综上所述，国内外学者对蒸气云爆炸进行了较为深入和广泛的研究，但是对不同浓度下含有CO或者H2的易燃易爆气体蒸气云爆炸温度和压力变化规律研究的较少。所以，文章以埃特金插值数学原理为基础，采用软件模拟研究不同的体积分数下CO和H2蒸气云爆炸温度和压力变化规律。

1 理论分析

基本理论依据易燃易爆物质爆炸后温度和压力计算的热力学原理及埃特金插值的基本原理和基本公式。

1.1 爆炸温度与压力计算依据

易燃易爆的介质与空气混合后形成爆炸性混合物，爆炸后的温度和压力可以依据热力学的基本原理进行计算。首先计算出反应前所有反应物的热力学能，根据能量守恒，推导出爆炸后的温度；根据理想气体状态方程，推导出爆炸后的压力。需要说明的是，由于爆炸性混合物发生爆炸的时间非常短暂，所以假设反应系统为绝热系统。具体计算公式为：

式中：∑Uc—系统发生化学反应后释放的总热力学能，kJ；∑Uf—反应物反应前在一定的温度下的总热力学能，kJ；n—参与反应的易燃易爆物质的摩尔数，mol；Qr—易燃易爆物质的燃烧热，kJ /mol。

1.2 数值模拟基本原理

采用模拟软件进行数值模拟，依据的数学原理就是埃特金插值。埃特金插值有两个典型的特点，一是逐步线性插值，二是误差事后估计。埃特金插值的优点显而易见，首先是可以通过逐步提高插值的阶来降低插值的误差；其次是具有线性特点，计算比较方便。下面介绍具体的计算原理。

设Pk-1，k（x）、Pk-1，k+1（x）分别为k次插值多项式，该表达式的含义是，他们共同通过点x0，x1，…，xk-1，Pk-1，k（x），还通过点xk，Pk-1，k+1（x）和点xk+1，（xk≠xk+1）。这两个k次插值多项式由于通过的点不相同，其余项也不相同，具体的余项为：

其中，f（k+1）（ξ1）与f（k+1）（ξ2）不相等，但是因为前面两个插值多项式通过的插值点几乎相同，所以这两个导数差别不是很大，若ξ1与ξ2相差很小，则f（k+1）（ξ1）与f（k+1）（ξ2）可认为近似相等，即

根据前面的分析及推导结果，可以推出两个插值多项式余项的比值是：

从上式中解出f（x）得：

即：

通过上式进一步分析，如果将这个误差加上，则可以得到更精确的k＋1次插值多项式，即：

由上式可以看出，两个k次插值多项式Pk-1，k（x）和Pk-1，k+1（x）经线性组合后，便得到k＋1次插值多项式Pk，k+1（x）。这种方法称为埃特金逐步线性插值。

2 数值模拟及分析

2.1 CO 蒸气云爆炸模拟

根据埃特金插值数学原理以及热力学基本原理，编制MATLAB插值程序，进行数值模拟。对体积分数为0.125～0.740的CO形成的无约束的蒸气云爆炸温度峰值和压力峰值进行统计，绘制变化曲线见图1。

图1 不同体积分数的CO蒸气云爆炸温度峰值、压力峰值曲线

由图1 可知，温度峰值和压力峰值随着体积分数的增加而先增大后减小。当体积分数为0.296时，可燃物和助燃物正好达到完全反应浓度配比的要求，所有反应物除惰性气体外都发生化学反应转化成生成物，所以爆炸后温度峰值和压力峰值最大，温度峰值最大值达到3 227 K，压力峰值最大值达到1.25 MPa。当体积分数为0.740时，可燃物和助燃物达到化学反应极限配比浓度，反应不充分导致燃烧释放的热量过少，所以爆炸后温度峰值和压力峰值最小，温度峰值最小值为1 563 K，压力峰值最小值为0.499 MPa。

温度峰值最大值是最小值的2.06倍，压力峰值最大值是最小值的2.5 倍。在体积分数为0.300 时，温度峰值有所下降，但下降的幅度不是很明显，由3 227 K下降到3 211 K，下降了0.49%；压力峰值下降比较明显，由1.250 MPa下降到0.925 MPa，下降了26%。压力的明显下降是由于系统中参与反应物质和剩余物质发生变化导致。

2.2 H2 蒸气云爆炸模拟

根据埃特金插值数学原理以及热力学基本原理，编制插值程序，进行数值模拟。对体积分数为0.04 ～0.75 的H2形成的无约束蒸气云爆炸温度峰值和压力峰值数据进行统计，结果见图2。

图2 不同体积分数的H2蒸气云爆炸温度峰值、压力峰值

由图2 可知，温度峰值和压力峰值随着体积分数的增加先增大后减小。当体积分数为0.30时，可燃物和助燃物正好达到完全反应浓度配比的要求，所有反应物除惰性气体外都发生化学反应转化成生成物，所以爆炸后温度峰值和压力峰值最大，温度峰值最大达到3 508 K，压力峰值最大达到1.360 MPa。当体积分数为0.04时，可燃物和助燃物达到化学反应极限配比浓度，反应不充分导致燃烧释放的热量过少，所以爆炸后温度峰值和压力峰值最小，温度峰值最小为823 K，压力峰值最小为0.283 MPa。

温度峰值最大值是最小值的4.26倍，压力峰值最大值是最小值的4.81 倍。在体积分数为0.32 时，温度峰值有所下降，但下降的幅度不是很明显，由3 512 K下降到3 417 K，下降了2.7%。在体积分数为0.32时，压力峰值下降比较明显，由1.361 MPa下降到0.989 MPa，下降了27.33%，压力的明显下降是由于系统中参与反应物质和剩余物质发生变化导致。

3 结论

通过采用模拟软件建立拟合公式，模拟了CO、H2蒸气云爆炸的温度峰值和压力峰值变化规律，发现CO、H2体积分数直接影响了蒸气云爆炸的温度峰值和压力峰值，具体结论如下：

1）CO体积分数自0.125向0.740增大的过程中，蒸气云爆炸的温度峰值和压力峰值先增大后减小。当体积分数为0.296 时，蒸气云爆炸的温度峰值和压力峰值达到最大值；当体积分数为0.740 时，蒸气云爆炸的温度峰值和压力峰值达到最小值；温度峰值最大值是最小值的2.06倍，压力峰值最大值是最小值的2.5倍。

2）H2体积分数自0.04 向0.75 增大的过程中，蒸气云爆炸的温度峰值和压力峰值先增大后减小。当体积分数为0.30时，蒸气云爆炸的温度峰值和压力峰值达到最大值；当体积分数为0.04时，蒸气云爆炸的温度峰值和压力峰值达到最小值；温度峰值最大值是最小值的4.26倍，压力峰值最大值是最小值的4.81倍。

3）化工企业可根据管道气体性质，参照此方法模拟蒸气云爆炸的温度峰值和压力峰值，选用相应等级的抗高温和抗高压装置设备，保障化工企业的整体安全性。