基于鸡群算法的危险化学品泄漏扩散事故反演方法

陈国华 胡雪薇 周利兴 门金坤 陈清光

(1.华南理工大学 机械与汽车工程学院，广东 广州 510640；2.广东省安全生产科学技术研究院，广东 广州 510640； 3.广东省安全生产科技协同创新中心，广东 广州 510640)

由化工厂危化品泄漏扩散引发的火灾爆炸事故时有发生，例如：2010年中石油兰州石化公司罐区一球罐出料管开裂，碳四物料大量泄漏，遇明火爆炸，造成6人死亡，6人受伤[1]；2018年河北省张家口市中国化工集团盛华化工有限公司氯乙烯气柜发生泄漏，泄漏的氯乙烯扩散到厂外遇明火爆炸，造成24人死亡，22人受伤[2]；2019年三门峡市义马气化厂空气液化装置倒塌，导致液氧贮槽破裂，大量液氧泄漏，泄漏源周围可燃物在液氧或富氧条件下发生爆炸、燃烧，造成15人死亡，16人受伤，直接经济损失8 170万元[3]。危化品泄漏扩散导致的火灾爆炸事故对人民生命财产安全构成了巨大损失，因此快速准确地定位泄漏源位置及计算泄露源强成为切断“泄漏-火灾-爆炸”连锁事故发生的关键。

在实际事故场景中，危险化学品气体扩散浓度可由传感器监测得到，但并不能根据该数据直接确定泄漏源位置及强度，因此，需要通过反演方法对泄漏气体浓度进行反演推算，进而确定泄漏源强及位置。目前，国内外的泄漏扩散事故反演方法可分为基于统计理论的方法与基于优化理论的方法两类。在基于统计理论的方法方面，Yee[4]使用贝叶斯概率理论推导出了源的数量，并表征了每个源参数的后概率密度函数，然后通过正时拉格朗日随机模型建立了多源分布与检测器阵列测量浓度之间的映射；Keats等[5]利用马尔科夫链蒙特卡罗方法对源参数的后验分布进行采样，验证了该方法在城市环境的高扰动流场反演问题上的可行性；Iovino等[6]利用统计归纳方法对城市中苯的来源进行辨识，验证了统计方法在泄漏源反演问题上的有效性。

在基于优化理论的方法方面，Lushi等[7]使用线性最小二乘间接优化方法逆向估计泄漏源的排放率；Ma等[8- 9]比较了粒子群优化算法、蚁群算法和萤火虫算法对气体发射源的识别能力，提出一种基于浓度分布最优相关匹配的智能优化方法，该方法考虑了大气条件和噪声数据对算法精度的影响，对数据精度的依赖性较小；Qiu等[10]使用人工神经网络进行扩散模拟，采用基于期望最大化的粒子群优化算法进行反演，有效地加速了收敛过程，但超出训练范围的神经网络预测结果不够准确。Cho等[11]和Kim等[12]使用随机森林分类器来追踪可疑泄漏源，但溯源结果准确度最高仅达86.5%。Chen等[13]采用粒子群优化算法对泄漏源项进行反演，算法计算速度、精度及稳定性较强。

基于统计理论的方法充分利用先验信息给出反演参数的不确定度和概率分布[14]，但采用这类方法进行事故反演前，需获得参数的先验分布，并且参数的后验分布抽样过程需较大的时间成本，需要提高算法的计算效率。基于优化理论的方法一般只需提供现场的监测数据，更适用于数据有限、现场情况复杂、要求快速响应的情况。

文中将高斯烟羽模型作为泄漏扩散的计算模型，以气体浓度传感器监测的气体浓度数据与模型计算值的误差最小化为优化目标，运用鸡群(Chicken Swarm Optimization，CSO)算法对泄漏扩散事故进行反演求解。由于初始解集对鸡群算法的性能有着较大影响，采用加权质心法生成优越的初始种群，从而在最初寻优时有效加速种群收敛，以提升算法求解效率。

1 泄漏扩散反演方法

1.1 坐标系建立与传感器布设

以化工厂内的泄漏扩散事故为研究背景，为快速有效地获得气体浓度数据及浓度场分布情况，在厂区内重点罐区及管道区域按矩阵阵列布设传感器，并以化工厂西南角作为原点、正东方向作为x轴正方向、正北方向作为y轴正方向、垂直于水平面且向上作为z轴正方向建立地理坐标系，以此得到各监测点所在地理坐标系中的坐标。气体浓度传感器在地理坐标系下的布设示意图见图1。

由于危化品泄漏扩散的浓度场与风场信息密切相关，气体的扩散会受风速、风向以及湍流等因素影响，因此在厂区内设置风速风向仪，以获取稳态的风速风向相关数据。将地理坐标系逆时针旋转角度θ可得标准风向坐标系，从而将气体浓度传感器地理坐标系位置转化为标准风向坐标系位置，如图2所示。坐标系转化关系可表示为

图1 气体浓度传感器在地理坐标系下的布设

(1)

式中，(xG，yG，zG)为气体浓度传感器在地理坐标系中的位置坐标，(xw，yw，zw)为气体浓度传感器在风向坐标系中的位置坐标。

图2 地理坐标转换示意图

1.2 基于高斯烟羽模型的危化品泄漏扩散反演模型

受泄漏物质危险性、毒害性及实验场地等因素的限制，危险化学品泄露扩散实验难以开展，故文中通过气体扩散模型模拟突发危化品泄漏事故中的气体泄漏扩散过程，获得传感器接收到的气体浓度“监测数据”。模拟气体扩散主要采用高斯模型、箱式模型、拉格朗日/欧拉模型以及计算流体力学模型[15]，其中，高斯模型可用于非重气扩散浓度计算，具有模型简单、计算量小、可快速计算出较好结果等优点，适用于实时应急信息系统。因此，文中将高斯烟羽模型作为泄漏扩散的计算模型[16]，表示为

(2)

σyi(x)、σzi(x)随泄漏源与传感器间的距离x而变化，若仅提供一组气体浓度数据和相应的地理坐标信息，在泄漏源位置、泄漏点距传感器距离均未知的情况下，无法确定σyi(x)、σzi(x)。高斯烟羽模型不具备求解此类复杂组合优化问题的能力，因此，将危化品泄漏扩散事故反演问题转化为一个优化目标最小化的问题展开研究。以气体浓度传感器阵列中n个测量值发生明显变化的气体浓度传感器的实际监测值与理论监测值之间的误差绝对值最小化为目标函数，表示为

(3)

1.3 鸡群算法

由于初始解集的输入对鸡群算法的结果精度有一定影响，而加权质心定位算法可通过权重缩小定位范围，故采用加权质心法生成初始解集，以提高算法的精确度并加速收敛。在标准风向坐标系中，基于高斯烟羽模型的气体浓度分布如图3所示，在此基础上建立加权质心公式。

图3 气体浓度分布图

初始解集规模为pop，任意j∈[1，pop]，其对应解集空间(Q0，xw0，yw0，zw0)通过式(4)、(5)求得：

(4)

(5)

式中，L为监测点间的距离，i为传感器索引，j为初始鸡群索引，ωi为第i个监测点所占权重，xi为第i个监测点的横坐标值，yi为第i个监测点的纵坐标值，Qmax为源强最大值，PRE为气体浓度传感器最小精度。

获得初始解集后，对鸡群个体进行适应度计算，根据适应度值进行排序并确定等级秩序。不同等级的鸡按相应规则更新其位置，位置更新规则[17]如下：当索引为i的个体是公鸡时，在搜寻的过程中拥有更大的搜寻范围，根据式(6)、(7)更新位置：

(6)

(7)

当索引为i的个体是母鸡时，跟随其聚类中的公鸡搜寻食物，其位置更新如式(8)-(10)所示：

(8)

(9)

S2=exp(fr2-fi)

(10)

式中，Rand为0～1之间的随机数，S1为公鸡的索引，S2为小鸡的索引。

当索引为i的个体是小鸡时，小鸡在与其存在母子关系的母鸡周围更新位置，如式(11)所示：

(11)

表1 鸡群算法伪代码

算法伪代码中，f(x)为适应度函数，输出结果为最优解。Rank(pop)是对鸡群根据适应度函数值进行排序，Divide(pop)是将鸡群分成不同组别，并确定小鸡跟随母鸡关系。

2 仿真与分析

泄漏扩散事故反演模型中鸡群算法的相关参数选取如表2所示。M、pop是由下述多组测试仿真所确定的最优值，在进行测试仿真的过程中，每次只改变1个参数，通过泄漏源位置、源强与预设值的误差大小来判断参数取值是否合理。Pr、Ph、Pc、FL、G参考Yu等[17]的建议取值；平均风速u、大气稳定度AT须结合厂区实际环境取值，文中假设事故场景发生在广州市某化工厂；Q0参考Pang等[18]的建议取值。上述参数的设置可以使算法在短时间内找到泄漏源位置并确定泄漏源强，因此，认为上述参数设置足以满足文中仿真需要。

表2 仿真参数[17- 18]

2.1 迭代次数对收敛过程的影响

在算法优化过程中，运算时间与迭代次数成正比，迭代次数直接决定了算法的计算量。算法所得最优解会随着迭代次数的递增而逐渐靠近真实最优解，所以迭代次数的优化选择对最终收敛结果有着决定性的影响。如图4所示，优化计算过程中，目标函数、泄漏源强和泄漏源位置(x，y，z)分别在40、50、70、50、100次左右收敛。相较于其他群智能算法收敛速度慢的情况，该算法计算的5个参数基本都收敛于100次以下，证明了鸡群算法快速收敛的优越性。此外，CSO算法的收敛值都以极小的误差接近于目标值，表明算法未因收敛快而陷入局部最优，反映了其快速准确的特性。

2.2 种群数量对收敛过程的影响

鸡群算法受初始设置的种群数量影响较大[17]，种群数量设置过小会导致无法收敛或陷入局部最优，种群数量设置过大有利于找到全局解，但计算成本也会大大增加。因此，若要在现有算法结构下找到算法的最优计算状态，就需要根据需求选择一个合适的种群数量。本节选择种群数量分别为50、100和200进行计算，结果如图5所示。种群数量为50的曲线在迭代70次左右收敛，种群数量为100和200的曲线在迭代30次左右收敛。相较于种群数量为50的情况，种群数量为100时整个优化计算过程的性能明显提高，而种群数量为100和200时的曲线无明显差异，综合考虑计算效率、成本与必要性后，选择种群数量为100进行优化反演。

2.3 监测点网格对计算精度的影响

如图6所示：监测点网格设置为2×2时，因网格数量过少算法不能得到正确解；监测点网格设置为3×3时，受实际场景中风场等因素的影响，获得有效监测值的监测点数量较少，使得监测浓度与计算浓度差异较小，算法快速达到终止条件，以致搜索结果并不能满足实际要求。

(a)目标函数的收敛过程

(b)泄漏源强的收敛过程

(c)泄漏位置x的收敛过程

(d)泄漏位置y的收敛过程

(e)泄漏位置z的收敛过程

图5 种群数量对收敛过程的影响

(a)监测点网格对源强及位置的影响

而当网格设置为5×5及以上时，监测网络中能够提供足够的监测数据供反演过程利用，保障了数据的充分性，反演得到的源强结果能够很好地趋近于真实源强。相较于3×3网格时的参数计算结果，5×5网格所得的位置误差下降了0.3 m左右，源强误差下降了0.6‰左右，并且计算时间仅为1.3 s。而从5×5网格开始，增加监测点网格后计算位置误差下降小于0.1 m，源强误差下降小于0.1‰，同时计算成本显著增加，11×11网格时的计算时间最高达6.9 s。因此，根据实际情形需求，综合考虑监测成本和计算成本，选择以5×5为合理布置的监测点网格。

2.4 鲁棒性验证

在真实的事故场景中，监测气体浓度传感器可能会由于外界因素干扰出现“异常值”。若输入数据中出现异常值，算法依旧可以实现优化目的，即异常值对算法计算出的结果干扰不明显，说明该算法的鲁棒性较强。本节采用人工输入有限个随机噪声数据的方法来测试鸡群算法的鲁棒性。具体做法是：在25个监测点中随机选取1至5个点进行不同噪声强度条件(编号1至14)下的噪声数据干扰。具体的噪声数据见表3。

表3 噪声数据1)

输入噪声数据后，鸡群算法的精度变化如图7所示。若25个监测点中同时有5个点产生超过自身数据5%的噪声，则源强误差达到11.67%(小于15%)，处于可接受范围内。若同时有3个点产生超过自身数据15%的噪声，对泄漏位置x、y两个参数的精确度影响较大。但若异常值的数量与噪声强度组合超过以上两组数据，误差不可接受。若25个监测点中仅有一个或两个监测点数据异常，即使噪声值达到原值的50%，泄漏源位置x、y和源强Q的计算结果的绝对误差分别不高于4 m、10 m和100 g/s。在实际场景中，多个传感器同时受外部因素影响导致出现异常值的情况很少，所以，一般情况下，鸡群算法的鲁棒性可以满足突发事故现场的泄漏源定位需求。

(a)噪声对泄漏位置x的影响

(b)噪声对泄漏位置y的影响

(c)噪声对源强的影响

2.5 算法性能对比

粒子群优化(Particle Swarm Optimization，PSO)算法和遗传算法(Genetic Algorithm，GA)是泄漏源扩散事故反演中使用较多的两种智能算法，定位较为准确且反演速度快。CSO、PSO及GA算法都是启发式智能优化算法，理论上都可以搜索到全局最优解，但寻优机理不同。3种算法的反演精度如图8所示。采用CSO、PSO、GA算法时，泄漏位置x的反演结果误差分别为0.12、3.00、7.00 m，泄漏位置y的反演结果误差分别为0.01、2.00、8.00 m，泄漏位置z的反演结果误差均超过1.00 m，波动范围分别为4.00～5.50 m、3.50～7.00 m、4.50～7.00 m。由此，可计算出CSO、PSO、GA算法在x方向上的反演结果相对误差分别为0.24%、6.00%、14.0%，在y方向上的反演结果相对误差分别为0.02%、4.00%、16.00%。CSO算法的在泄漏位置x与y方向上的反演精度显著高于PSO、GA算法，且通过(x，y)坐标已能准确确定泄漏源的平面位置。此外，在z方向上，3种算法的误差均超过1.00 m，但CSO算法的反演结果波动范围相较于PSO、GA算法更小，这是因为CSO算法主要依靠“鸡群行为”进行仿生模拟，小鸡在母鸡周围觅食，母鸡跟随公鸡进行觅食，通过“代”的更替扩大搜索空间；而GA算法是依靠“变异”操作跳脱局部极值、扩大搜索空间，以种群适应性为评价准则，所以实现群体收敛需要的时间较长；PSO算法则是依靠粒子与种群的经验反馈来确定搜索速度和方向，隐含并行搜索的特性，故搜索效率较高。3种算法的计算成本(时间)如图8(d)所示。可以发现，CSO算法的计算时间为0.7 s，PSO、GA算法的计算时间分别为0.3与4.5 s，GA算法的计算成本约为CSO算法的6.43倍，说明CSO算法可极大减少计算成本。此外，CSO算法的计算时间相较于PSO算法仅增加0.4 s，可忽略不计。综合考虑反演精度与时间成本，可以认为，CSO算法在危化品泄漏扩散反演中的表现较为突出。

(a)泄漏位置x的反演结果对比

(b)泄漏位置y的反演结果对比

(c)泄漏位置z的反演结果对比

(d)计算时间对比

2.6 优化结果分布

智能优化算法一般是建立在生物智能或物理现象基础上的启发式算法，由于计算结果的随机性是启发式算法的固有属性，文中针对同一案例进行了100次仿真验证，通过比较仿真结果与预设泄漏参数的误差大小，判断反演效果的好坏。泄漏位置x的反演误差小于0.12 m，泄漏源强的反演误差小于1.2 g/s。泄漏位置y的反演效果最好(见图9(b))，正负误差为0.01 m，且大部分数据分布在50 m左右。泄漏位置z的反演效果最差(见图9(c))，正负误差为1.00 m，且大部分数据取分布在4.40 m左右。z方向产生较大误差的主要原因可能是：文中选取的泄漏点高度为5 m，监测点高度均为0.8 m，而非重气危化品气体泄漏后不会大规模向地面下沉，导致部署在0.8 m高度处的气体浓度传感器不能准确得到气体浓度数值，从而使得z方向上的精确度不高。但从实际情况看，事故发生初期仅需确定发生泄漏的装置或储罐，无需确定泄漏的具体坐标点，因此根据优化结果，文中方法满足危化品泄漏扩散事故反演的实际需要。

(a)泄漏位置x的反演结果

(b)泄漏位置y的反演结果

(c)泄漏位置z的反演结果

(d)源强的反演结果

3 结语

文中提出的危险化学品泄漏扩散事故反演方法采用阵列式布设方式布设气体浓度传感器，同时考虑了实际风速、风向等风场信息的影响，并将传统的地理坐标系转化为标准风向坐标系，简化了计算的复杂度。以气体泄漏实际监测数据与高斯烟羽泄漏扩散模型构建未知事故源参数优化模型，基于鸡群算法对模型进行优化求解，并采用加权质心算法创建较高质量的初始解集，进一步提高算法精确度与求解效率。

文中还研究了迭代次数、种群数量对算法收敛过程的影响，通过改变单一变量进行参数测试，综合考虑实际事故场景与反演结果的精度，讨论了监测点网格划分对反演精度的影响，确定5×5网格为最佳网格划分方式；并通过算法鲁棒性验证、与PSO和GA算法进行对比，发现采用鸡群算法时泄漏位置z的最大误差在1.00 m范围内，可满足突发事故现场的应急监测及泄漏源定位需要。

优化结果表明，基于鸡群算法的泄漏扩散事故反演方法能够实现化工厂内部气体泄漏源的快速准确定位以及泄漏源强的准确预估，可满足突发事故现场应急监测及泄漏源定位的需求，为应急决策提供及时有效的技术支持。后续研究中，将通过开展实验研究及事故案例分析等来对CSO算法的有效性进行进一步验证。