基于SLAB的重气泄漏扩散模拟分析

2020-07-14何娟霞周琪勇李昂昂周冬梅

科学技术与工程 2020年17期

何娟霞，刘磊，周琪勇，李昂昂，马野，周冬梅

(1.广西大学资源环境与材料学院，南宁 530004；2.华南理工大学材料科学与工程学院，广州 510641)

汽油、乙醇、液氯等具有燃爆性、毒性、腐蚀性的危险化学品一旦发生泄漏，将会形成重质气云。由于重质气云密度比空气大，在重力作用下，将会沉降到地面，并沿风向缓慢扩散。一旦发生火灾爆炸、中毒窒息事故，将给周边人群生命安全、大气环境及水资源安全带来严重危害[1-3]。

近年来众多学者对重质气云扩散分析采用SLAB模型。姜传胜等[4]通过SLAB模拟与风洞实验的对比，发现结果基本一致，并剖析了重气连续扩散的特点；黄江平等[5]通过改进SLAB模型中风速廓线的计算方法，得到适用于模拟城市重气扩散的SLAB_URBAN模型；朱红亚等[6]通过对SLAB模型预测值与试验观测值的比较，发现数据吻合较好，并在SLAB模型基础上发展了多源重气泄漏扩散模型；刘昕等[7]将SLAB模型用于氨泄漏扩散后下风向浓度分布的模拟；熊发等[8]对SLAB模型大气稳定度及地表粗糙度两个参数进行了优化；Li等[9]考虑气象环境尤其是风向、风速的时变性，提出了SLAB改进模型SLABi。众多学者为SLAB后续研究奠定了基础。

目前基于SLAB模拟下风向三维空间浓度分布的研究较少，在SLAB模型的基础上，结合应急需求，从地平面浓度分布、敏感点浓度随时间变化、浓度曲面随高度变化3个方面对扩散数据进行分析。以便快速、有效预测事故空间影响区域、减少人员伤亡，为事故应急、环境保护提供科学依据。

1 重气云团的形成与判定

1.1 重气云团的形成

重气云团(heavy gas, dense gas)是指分子量大于空气的气体形成的气云或气体本身分子量虽比空气小，但由于从高压、低温环境泄漏等原因形成的比空气重的气云[10-12]。重质气云的形成过程主要有①高压、低温液化气体泄漏到大气中，处于过热状态，迅速气化形成夹杂有液滴的重质气云，如液氯泄漏扩散；②常压液态物质泄漏后，在地面集聚形成液池，进而挥发产生重质气云，如汽油挥发扩散；③气态物质直接在空气中扩散形成重质气云，如硫化氢井喷[13-14]。

1.2 重气云团的判定

重气云团一般用理查森数(Richardson number，Ri)判定[15-16]。理查森数原是研究大气及海洋边界层的重要无量纲参数，后被推广到研究气象学中的对流及大尺度运动中[17]。当初始理查森数Ri0大于临界值R′，则认为是重气云团，用重气扩散模型对其扩散进行模拟；反之，则认为是非重气云团。临界值R′近似取50。Ri0计算公式为[12]

(1)

(2)

式中：g0′=(ρ0-ρa)g/ρa，m/s2；g为重力加速度，m/s2；ρ0为气体初始密度，kg/m3；ρa为大气密度，kg/m3；V0′为气体连续泄漏时所形成的云羽的初始体积通量，m3/s；u为环境风速，m/s；D为泄漏源的特征水平尺寸，m；u*为摩擦速率，m/s；V0为气体瞬时泄漏时所形成云团的初始体积，m3。

2 重气扩散模型

2.1 概述

重气扩散模型主要有唯象模型、箱模型、浅层模型、计算流体力学(computational fluid dynamics, CFD)模型等[18-19]，其中浅层模型结合了CFD模型与一维模型的优点，可以基本保持扩散控制方程的完备性，又能节约计算资源，提高计算效率，模拟结果吻合度较好[13]。

2.2 SLAB模型

SLAB模型由浅层模型发展而来，于1990年美国劳伦斯利佛莫尔国家实验室(Lawrence Livermore National Laboratory，LLNL)开发。SLAB模型考虑到空气卷吸、云团加热、羽流上升、云团曲流等对重气扩散的影响，通过求解质量、能量、动量、组分守恒等平衡方程，得到相应地点、时间的计算浓度。如表1所示，SLAB模型可对4种泄漏类型、3种气云扩散类型进行分析，扩散类型取决于泄漏源持续泄漏时长(continuous source duration，TSD)。

表1 泄漏类型及气云扩散类型Table 1 Types of leakage and gas diffusion

SLAB手册[20](User’s Manual for SLAB)对SLAB模型做了详细的说明，包括数学建模过程、SLAB程序的使用指南、SLAB程序使用案例，如图1所示。数学建模过程包括稳态羽流扩散、瞬时烟团扩散、子模型及时间平均浓度求解。其中子模型包括环境流速剖面模型、卷吸率、热与动量描述、热力学模型、羽流上升模型。图2所示为SLAB模型计算流程，首先确定泄漏源类型，然后选择扩散模型，进而根据扩散方程计算空间平均浓度，最后得到时间平均浓度。

图1 SLAB手册架构Fig.1 Structure of SLAB manual

图2 SLAB模型计算流程Fig.2 Calculation flow of SLAB

3 案例分析

3.1 泄漏场景介绍

以南宁某纸业公司氯气钢瓶整体破裂瞬时泄漏扩散为例。氯(chlorine)具有刺激性、氧化性、毒性，可对大气、水体产生严重危害。美国工业卫生协会(American Industrial Hygiene Association，AIHA)制定的紧急反应计划指南(Emergency Response Planning Guidelines，ERPGs)[21]将氯浓度水平划分为3个等级，如表2所示。本案例取3、9、58、300 mg/m34个浓度标准，分别作为轻危区、中危区、重危区、致死区4个应急区域的划分界限[15, 19]。

该公司周围地形平坦，无高大建筑，厂内员工500人，周边环境风险受体如表3所示。液氯钢瓶容积为715 m3，最高工作压力为1.1 MPa，正常工作温度为-12 ℃，采用单层材质，最大储存量为1 t。

表2 氯气应急浓度标准Table 2 Emergency concentrations of chlorine

表3 泄漏点周边环境风险受体Table 3 Risk receptors near leakage point

3.2 瞬时烟团模型控制方程

氯气钢瓶破裂瞬时泄漏扩散符合SLAB中的瞬时烟团模型。该模型根据质量、动量、能量、组分守恒方程，用空气卷吸的概念衡量湍流对云团扩散的影响，将云团质心的扩散时间作为独立变量。得到t时刻下风向(x,y,z)点的体积浓度值C(x,y,z,t)的计算式为

C(x,y,z,t)=4BxByhC(t)C1(x-

Xc,bx,βx)C1(y,by,βy)C2(z,Zc,σ)

(3)

t时刻的云团体积平均浓度C(t)为

(4)

水平浓度轮廓函数C1(y,by,βy)为

(5)

(6)

附加浓度轮廓函数C1(x-Xc,bx,βx)的对应法则与水平浓度轮廓函数C1(y,by,βy)相同，只需将对应的参数替换即可。垂直浓度轮廓函数C2(z,Zc,σ)为

(7)

(8)

式中：Bx、By分别为云团下风向与侧风向半宽，m；bx、βx为云团下风向半宽参数；by、βy为云团侧风向半宽参数；h为云团高度，m；m(t)为t时刻的云团质量浓度；Ma、Ms分别为干空气、泄漏源物质的分子量；Xc为云团质心坐标，m；Zc为高度参数；erf为误差函数。

3.3 参数确定

广西壮族自治区南宁市属于亚热带季风气候，全年主导风向为东风、南风。表3所述泄漏点北侧人口最为密集，考虑事故后果的严重程度，案例设定事发时风向为南风。SLAB模拟需要输入的参数及其取值情况如表4所示。该环境条件下云团理查森数约为4 794，属重气云团。

3.4 模拟结果及分析

浓度模拟结果都是时间平均值。所谓某一点的时间平均浓度，是指在时间区间(tpk-TAV/2，tpk+TAV/2)内的浓度平均值，如式(9)所示。

(9)

式(9)中：Ctav为时间平均浓度，mg/m3；TAV为浓度平均时间，TAV不应大于云团持续时间TCD，s；tpk为(x,y,z)点浓度达到最大值时的时间，s；C(x,y,z,t)为(x,y,z)点在时间点t的浓度，mg/m3，见式(3)。

3.4.1 浓度分布曲线

根据已有数据点，调用MATLAB中griddata命令的cubic插值法，绘制下风向半边浓度云图，如图3所示。下风向地面处(z=0 m)氯气浓度分布情况如图4所示，图中由外到内4个闭合曲线代表由低到高的4个浓度标准。

结果显示：厂内大部分位于致死区，浓度在300 mg/m3以上，且事故发生后氯气云团随风流扩散较快，厂区最远点(下风向300 m处)浓度到达最大值只需85 s；A村处于事故点下风向位置，约50%面积位于重危区内，受氯气云团影响较大，且居民较多，是本次事故的重点敏感区域；B村距事故地点1 200 m，但由于偏离了下风向位置，只有较少的面积受到影响，承受的风险远小于A村；C、D、E村由于远离下风向位置，几乎不受事故的影响。分析可知，重气云团对下风向区域影响最大，各浓度值影响区域都近似呈现狭长的椭圆形，厂区内人员及A村居民向垂直于下风向的方向逃生，可较快脱离危险区域。

3.4.2 浓度-时间曲线

A村中心点(下风向1 500 m处)的浓度-时间曲线如图5所示。当泄漏发生后，150 s内空气中氯气浓度低于0.000 8 mg/m3；150 ～ 200 s期间，氯气浓度逐渐缓慢升高；在200 ～ 285 s区间，空气中氯气浓度急速上升并达到最大值61.1 mg/m3；此后不断下降，约在450 s左右降至3 mg/m3以下，对人体的伤害逐渐减轻。A村中心点暴露在50～60 mg/m3浓度区间内的时间约为60 s，能刺激人的呼吸道，引起咳嗽，可对健康造成一定的损害，没有严重的生命威胁。

表4 SLAB模拟参数汇总Table 4 Summary of SLAB input parameters

图3 半边浓度云图Fig.3 Half of concentration contour map

图4 下风向地面处氯气浓度平面分布Fig.4 Downwind concentration distribution of chlorine near the ground

图5 A村中心点浓度-时间变化曲线Fig.5 Concentration-time curve in the A village

3.4.3 三维浓度曲面

图6所示为各浓度三维散点轮廓图。图7所示为各浓度曲面侧面轮廓线。可见浓度曲面呈狭长椭圆拱形。近距离处由于云团中氯气的扩散，使浓度曲面的影响高度与影响宽度扩张；远距离处受空气卷吸与湍流影响，导致云团稀释、浓度曲面收缩。各浓度最大影响高度分别为118、82、43、23 m。由于负浮力(重力)效应，加之泄漏源温度较低引起的密度差，使云水平扩散距离远大于纵向爬升高度，浓度曲面呈现扁平型，且近地处浓度曲面随高度变化坡度较大，远地处曲面较平缓。