随着世界经济与科技的迅猛发展,特别是石油和化工原料品种的不断增加,生产、运输和仓储等存在爆炸性气体安全隐患的场所也在不断扩大,其中桶罐区的危险重质气体泄漏导致的事故频繁发生。为保证安全高效的生产,世界各国正致力于主动提高具有防止爆炸的仓储运输车辆和叉车的应用及防爆性能。目前的研究主要集中在以CO₂为代表的重质气体泄漏方面,如肖淑衡等^[1]用CO₂模拟危险性重质气体做空间泄漏扩散实验,周宁等^[2]利用CO₂模拟罐区重质气体泄漏扩散的基本规律,宋黎等^[3]以CO₂作为示踪气体研究室内有害气体浓度场分布。针对一种特定的重质气体进行研究的较少,尤其对CS₂的泄漏扩散研究。CS₂作为一种特殊的重质气体,通常作为添加剂、溶剂和磺化剂应用于石油化工领域,具有高易燃性及高神经毒性,属于危险品^[4]。它的蒸气压为53.32 kPa/28 ℃,熔点为-110.8 ℃,相对分子质量为76.14,相对密度2.64(空气=1),在常温下易挥发,与空气形成爆炸极限为1.3%～50%爆炸性混合物,静电荷、火星、高温、氧化剂等因子都会使CS₂有燃烧危险^[5]。^因此当事故发生后无法准确预测事故影响范围与后果。针对以上情况,本研究采用Fluent软件对CS₂桶罐储存室进行CS₂气体泄漏的数值模拟计算,研究在障碍物存在、不同通风速度、不同泄漏口等条件下CS₂泄漏扩散的规律及其浓度分布情况,为防爆叉车上传感器合理的空间布局,以及制定防火防爆措施提供参考。

从20世纪中后期至今,国内外学者一直在危险气体泄漏扩散的模拟研究方面不断地探索前行,其间提出了许多扩散计算模型,如BM模型^[6]、Gaussian模型^[7]、箱及相似模型^[8]、Sutton模型^[9]和浅层模型^[10]等。国内学者借鉴前人的经验并通过自身研究,提出了板块模型^[11],这些模型都适用于模拟气体泄漏扩散。其中Gaussian模型适用于点源扩散,以统计学方法为基础研究扩散质的浓度分布,它的优点是开发较早,技术较成熟,模型简单,便于计算,应用广泛; 缺点是只适用于轻质气体扩散且误差大。BM模型通过大量简单关系式和图表来描述扩散行为,作为筛选的基准模型适用于中性或重质气体研究。Sutton模型采用湍流扩散统计理论,只适用于中性气体,而且精度较低。箱及相似模型在概念清晰、计算精度等方面优于Gaussian模型,尤其是重质气体的扩散模拟,但模型通常涉及不连续界面,存在较大的不确定性。

上述的理论模型都有其适用范围和局限性,而三维传递现象模型运用CFD方法模拟^[12],通过数值求解室内空气流动的质量、动量、能量和危险气体扩散的控制方程^[13],得到危险气体浓度分布情况,特别适用于重质气体的模拟。因此本研究采用三维传递现象模型对CS₂气体的扩散过程及影响因素进行研究。

连续、动量、能量方程^[14]求解方法如下。

连续性方程(质量守恒方程):

(ρ)/(t)+((ρu_x))/(x)+((ρu_y))/(y)+((ρu_z))/(z)=0。(1)

式(1)中:ρ为密度,kg/m³; u_x、u_y、u_z分别为x、y、z 3个方向的速度分量,m/s; t为时间,s。

可压缩性流体的动量方程:

ρ(du_x)/(dt)=ρf_x+(p)/(x)+/(x){μ[2(u_x)/(x)-2/3((u_x)/(x)+(u_y)/(y)+(u_z)/(z))]}+/(y)[μ((u_x)/(y)+(u_y)/(x))]+/(z)[μ((u_z)/(x)+(u_x)/(z))],(2)

ρ(du_y)/(dt)=ρf_y+(p)/(y)+/(y){μ[2(y)/(y)-2/3((u_x)/(x)+(u_y)/(y)+(u_z)/(z))]}+/(z)[μ((u_y)/(z)+(u_z)/(y))]+/(x)[μ((u_x)/(y)+(u_y)/(x))],(3)

ρ(du_z)/(dt)=ρf_z+(p)/(z)+/(z){μ[2(u_z)/(z)-2/3((u_x)/(x)+(u_y)/(y)+(u_z)/(z))]}+/(x)[μ((u_z)/(x)+(u_x)/(z))]+/(y)[μ((u_y)/(z)+(u_z)/(y))]。(4)

式(2)～(4)中:f_x、f_y、f_z分别为x、y、z 3个方向的摩擦力分量; μ为流体与壁面摩擦系数。

能量方程:

Δ·[k_effΔT-∑_jh_jJ_j+(τ_eff·u)]=((ρE))/(t)+Δ·[u(ρE+p)]。(5)

式(5)中:h_j为组分j的焓,J/kg,h_j=∫^T_TrefC_p,jdT,其中T_ref=298.15 K; k_eff为有效热导系数,W/(m·K); k_t为湍流热传导系数; J_j为组分j的扩散通量; E为流体微团的总能,J/kg,为内能、动能、势能三者之和,即E=h-p/ρ+(u²)/2,其中h为焓,J/kg; u为速度矢量。

标准k-ε模型^[14]:

((ρk))/(t)+((ρku_i))/(x_i)=/(x_j)[(μ+(u_t)/(σ_k))(k)/(x_j)]+G_k+G_b-ρ ε-Y_M+S_k,(6)

((ρε))/(t)+((ρu_i))/(x_i)=C_1εε/k(G_k+C_3ε+G_b)/(x_j)[(μ+(μ_t)/(σ_ε))(ε)/(x_j)]-C_2ερ(ε²)/k+S_ε。(7)

式(6)～(7)中:k为湍流脉动动能; ε为湍流脉动动能的能耗率; σ_k、σ_ε分别为不同普朗特数,Fluent中默认值分别为1、1.3; G_k与平均速度梯度引起的湍动能有关; G_b与浮力影响引起的湍动能有关; Y_M由湍流脉动膨胀对总的耗散率的影响决定; C_1ε、C_2ε、C_3ε为经验值,Fluent中一般取值为1.44、1.92、0.09。

桶罐区CS₂泄漏扩散可认为是多组分气体相互作用的湍流,因此除上述基本方程外还需多组分输运方程来模拟CS₂泄漏扩散过程。

组分传输方程^[15]:

((ρc_a))/(t)+/(x_i)(ρu_ic_a)=/(x_i)(D_a((ρc_a))/(x_i))。(8)

式(8)中:c_a、ρc_a、D_a分别代表组分a的体积分数、质量分数和扩散系数。

马赫数是气体流动速度与音速的比值,记为O,在本研究中,马赫数小于0.4,即低速扩散流动,因此CS₂气体流动时需满足如下流动特性方程^[16]:

(dA)/A=(O ²-1)(dv)/v。(9)

图1 数值计算物理模型
Fig.1 Numerical computational physics model

式(9)中:dA为流管横截面积的微变化值; dv为气流速度的微变化值。

以某CS₂桶罐储存室为研究对象,数值模拟物理模型见图1,尺寸为20 m×12 m×8 m。假定图1中4个矩形框分别为多个CS₂桶罐规则摆放的桶罐区,其尺寸为3 m×4.5 m×1.8 m,通风口大小为3 m×3 m,出风口为3 m×1.8 m,泄漏直径为20 cm。泄漏口坐标分别为A(7.5,8.5,0.9)、B(6.25,7.5,0.9)、C(6.25,9.5,0.9)、D(3,8.5,0.9)、E(14.75,8.5,0.9)、F(6.25,8.5,1.8),除泄漏口所在的桶罐区外,其余桶罐区都视为障碍物存在。

数值模拟的成败首先在于给定模型计算域的大小及相关壁面的边界条件,同时不同的气体存在物理和化学性质上的差异,因此不同气体泄漏扩散的边界条件是不同的。以存储100%质量分数的二硫化碳桶罐室作为研究对象,分别以A、B、C、D、E、F 6个不同的泄漏口作为CS₂泄漏的单源头,对整个室内空间进行CS₂气体泄漏扩散模拟计算。

边界条件设定如下:1)进口边界条件,将泄漏口设置为质量入口,以泄漏速率为0.137 6 kg/s非稳态连续泄漏,方向垂直于入口截面; 2)出口边界条件,压力出口为通风口和出风口,压强为一个标准大气压,设气流方向为水平方向; 3)通风条件,利用大功率风扇强制对流扩散; 4)壁面条件。采用标准壁面函数处理,设置壁面类型为wall,温度为300 K,室内外不存在热传递。

本研究针对泄漏口距地面0.05 m的平面爆炸区域进行对比,图2是在障碍物影响下A处泄漏口经过500 s后的CS₂爆炸区域示意图。在无障碍物时,由于泄漏气体在初速度及重力的作用下,以泄漏口为中心沿四周地表扩散,CS₂桶罐储存室地表空间基本处于危险区,尤其靠近泄漏口一侧较为严重; 而在障碍物存在的情况下,CS₂爆炸区域范围比无障碍物形成的要大,且各区域CS₂的积聚情况更严重。这是由于障碍物的存在对障碍物周围的气流运动产生了影响,即改变了CS₂气体的扩散路径,同时由于CS₂相对分子质量是空气的2倍多,一般沿地面扩散,障碍物的存在,使CS₂气体扩散空间减少,扩散阻力增大,易积聚在障碍物周围。同时根据空气动力学,高浓度区的气体会向低浓度区运动,即CS₂气体会绕过障碍物向前进行扩散运动,因桶罐区内墙壁具有阻碍作用,促使CS₂气体往障碍物后侧方流动。这说明障碍物的存在短期内会增速滞留CS₂气体,同时也增大了CS₂桶罐区潜在的局部易爆区域面积。

图2 A处泄漏经过500 s后有无障碍物爆炸浓度区域
Fig.2 Area of explosive concentration after 500 s at the location of A

when obstacles being present or absent

分别取A、B、C、D、E、F 6个位置作为泄漏源,泄漏速度为0.013 74 kg/s,持续泄漏时间为500 s,入口风速为1 m/s。由图3可知,当CS₂桶罐区分别从A、B、C、D 4个水平方向泄漏时,由于CS₂的重质效应,CS₂泄漏口附近地表快速积聚,易使其浓度达到爆炸极限范围,提高CS₂桶罐储存室爆炸的危险系数,从通风口到出风口整个CS₂桶罐区危险气体的含量先增大后减小。

当E、F处发生泄漏时,CS₂桶罐储存室大部分空间都处在爆炸极限范围内,但CS₂的含量较低,这是由于在E、F处泄漏,泄漏气体本身的初始势能较大,同时受重力作用,将势能转化成动能,泄漏气体的扩散速度相对于水平泄漏口较快,不易在局部地区积聚。

图3 经过500 s后各处泄漏的CS2含量爆炸区域
Fig.3 Explosive area of CS2 concentration at different locations after 500 s

取A作为泄漏源,泄漏方向与通风口方向一致,通风口风速v分别为0、1、3、5 m/s,泄漏速度为0.013 74 kg/s,持续泄漏时间为500 s。由图4可知,随着通风速度的增大,泄漏气体达到爆炸的含量之区域面积明显缩减,且爆炸的CS₂含量区域对应泄漏气体的含量也有明显的降低,主要集中在泄漏源一侧,说明高的通风速度有利于泄漏气体的扩散,减少泄漏气体危险域。这是由于泄漏气体在扩散过程中受到外界风速的干扰,室内空气的流动不断卷吸泄漏气体,提高了稀释泄漏气体的速度,在空气气流的作用下能加快从出风口排除,同时通风口和出风口连线之间无障碍物,风速会加快此处气流运动,阻碍泄漏气体往另一侧扩散,风速越大,阻碍作用越大。

图4 经过500 s后不同速度的CS2含量爆炸区域
Fig.4 Explosive area of CS2 concentration at different velocities after 500 s

本研究运用Fluent软件分别从障碍物、泄漏口位置、通风速度不同角度模拟分析CS₂气体在桶罐室泄漏扩散情况,得到相关规律及建议如下。

1)CS₂气体的扩散方式具有重质气体沿地表扩散、速度慢等基本特征。2)当CS₂桶罐储存室存在障碍物时,对CS₂气体扩散具有阻碍作用,泄漏的CS₂气体易积聚在障碍物周围,增大了爆炸事故发生的概率。3)不同泄漏口形成的爆炸区域不同,CS₂气体沿着桶罐区四周泄漏时,CS₂气体局部积聚较严重,整个CS₂桶罐储存室基本处在爆炸极限范围内; 而沿着垂直方向泄漏时,CS₂气体扩散速度较快,在整个CS₂桶罐储存室分布较均匀且含量较低。4)在一定风速范围内,CS₂气体的扩散速度与风速大小成正相关,随着风速的增加,达到CS₂气体爆炸极限的范围越来越小。5)综合CS₂泄漏扩散的基本情况,可以将CS₂传感器多点位,低平面(窄幅垂直面)合理布局在防爆叉车上。6)利用障碍物易改变流场作用与CS₂的易积聚特性,可以筑堤堵截泄漏二硫化碳或者引流到安全地点。

鉴于CS₂的危险性及试验设施的缺乏,无法进行真实状态下CS₂泄漏扩散实验,仿真结果势必与现实情况存在一定偏差。因此,如有条件在安全措施完备的情况下,可以进行试验,获取真实数据。本文的研究思路与结论可以为其它危爆气体的主动检测提供借鉴。