1、文章编号:1009-6094(2023)07-2159-08地震诱发储罐火灾多米诺效应的动态定量风险评估方法研究*陈海翔，唐珑枰，方伟(中国科学技术大学火灾科学国家重点实验室，合肥 230026)摘要:为评估地震对化工储罐区可能造成的火灾多米诺效应的风险，提出了一种基于蒙特卡洛模拟和矩阵计算的动态风险评估方法。该方法考虑了时空演化过程的不确定性和协同效应，可用于多个储罐着火的主要事件场景和高级别多米诺传播的复杂情况，能够得到可能发生的主要事件场景及其发生概率、储罐群全过程失效概率、特定主要事件场景演化细节，进而指导制定更具针对性的应急管理措施。通过案例分析了不同地震和储罐条件，以及多米诺级别对
2、储罐群多米诺失效概率的影响，验证了该方法的有效性，该方法可为多米诺效应的风险应急管理提供理论依据。关键词:安全工程;化学工业安全;Natech 事件;多米诺效应;不确定性;风险分析;蒙特卡洛模拟中图分类号:X937文献标志码:ADOI:10.13637/j issn 1009-6094.2022.0168*收稿日期:2022 02 23作者简介:陈海翔，副研究员，从事火灾及灾害风险评估研究，hxchen ustc edu cn。基金项目:国家重点研发计划课题(2016YFC0800604)0引言由于规模效应需求，化工业倾向于集群化，设施密集1 2，一旦发生极端自然灾害易导致严重的技术事故(Na
3、tech 事故3)，对密集工业设施造成巨大破坏，引发多米诺效应，导致巨灾4。典型事故如2011 年，日本东北部大地震及其引发的海啸冲击福岛，对储罐区造成了毁灭性的破坏5;2017 年，飓风“哈维”引发了风暴潮和降雨，冲击美国休斯敦地区，造成了大量石油化工设施损坏6。近年来Natech 事件导致的多米诺效应越来越受到关注，成为研究热点7。其中地震因其巨大破坏力、诱发的Natech 事件后果严重，引起了重点关注。2008 年，Campedel 等8 对地震引发的 78 个 Natech 事件进行了统计，发现常压储罐在所有受损设备中占比最高(80.3%)。此 外，据 统 计9，19512012 年，
4、约41.8%的多米诺事故发生在储存区域。2011 年，Abdolhamidzadeh 等10 分析了 19102008 年发生的224 起工业事故，发现 43%的多米诺事故由火灾引发。因此，对化工储罐区的常压储罐群因地震诱发火灾多米诺效应进行定量风险评估具有重要意义。当前针对地震引起 Natech 事件诱发多米诺效应的综合研究较为缺乏，更多是把 Natech 事件和多米诺效应分别进行研究。一方面，先前研究已经将地震纳入了 Natech 事件量化风险分析(Quantitativeisk Assessment，QA)的具体框架11 13。如 2003年，Salzano 等11 首先通过数据统计，利
5、用 Probit 模型表征储罐致损概率，简化了地震触发储罐事故的量化风险分析计算;2007 年，Antonioni 等12 最先开发了一个可分析地震触发事故风险的 QA 框架。但上述研究没考虑 Natech 事件可能影响周围设备并引发多米诺效应造成事故升级。另一方面，前人提出并验证了多种多米诺效应的风险量化分析方法14 18。如 2005 年，Cozzani 等14 首次通过使用后果、损坏和升级评估的 Probit 模型，提出了一种针对多米诺效应的 QA 方法;2019 年，华敏等15 率先综合运用动态贝叶斯方法与图示评审技术网络对多米诺效应进行建模。然而，上述研究未考虑 Natech 事件诱
6、发主要事件的不确定性。因此，发展一种综合考虑地震烈度不确定性与 Natech 事件导致多米诺效应的定量评估方法有重要意义。此外，化工储罐区往往涉及大量设施，地震可能造成多个储罐同时破坏的主要事件场景，导致更高级别的多米诺效应。目前研究大多只考虑了一级多米诺效应14 16，且多数基于场景分析的方法无法很好地适应复杂场景15，17 18。事故升级传播过程会涉及多个事故储罐的物理效应同时作用于一个未发生事故储罐的情景(协同效应)，且此未发生事故储罐所受能量可能动态变化。因此，需要考虑过程的协同效应，发展一种能适应复杂场景的动态化定量风险评估方法。本文旨在发展一种以风险为导向，针对地震诱发常压储罐群的
7、火灾多米诺效应的动态化定量风险评估方法。该方法基于矩阵计算和蒙特卡洛模拟，相较于以往研究，考虑了地震诱发主要事件场景的不确定性和过程的协同效应，改进了多米诺升级模型，增加了过程随机性，能便捷地获得储罐群在主要事件场景的多米诺效应细节和全过程的多米诺失效概率，从而更好地指导风险管理活动。1研究方法如图 1 所示，针对化工储罐区地震造成 Natech事件从而诱发火灾多米诺事故，本文提出的动态化9512第 23 卷第 7 期2023 年 7 月安全 与 环 境 学 报Journal of Safety and EnvironmentVol 23No 7Jul，2023定量风险评估方法主要包括 4 个
8、版块:1)初始数据收集;2)地震造成的主要事件场景评估;3)地震造成的火灾多米诺升级评估;4)全过程风险评估。该方法的核心是地震造成的火灾多米诺升级评估板块中基于矩阵计算和蒙特卡洛模拟的动态化计算流程。图 1地震诱发火灾多米诺效应研究方法流程图Fig 1Flowchart of the methodology of fire-induceddomino effect scenarios triggered by earthquakes1.1初始数据收集初始数据的收集包括:1)储罐区目标储罐辨识;2)地震危害描述;3)应急救援行动量化描述。储罐区目标储罐收集的数据主要包括储罐类型、体积、几何特征
9、、储存量、布局及储存物质特性等。地震危害收集的数据包括当地地质参数和气温、风向、湿度等环境因素，以求出特定重现期的地震烈度概率分布和辨识可能造成的 Natech 主要事件场景。应急救援行动方面主要是收集当地的应急计划和缓解行动的历史数据，从而帮助专家确定灾害控制时间的对数正态分布19。对地震危害的定量描述，峰值地面加速度 A 是一个 常 用 烈 度 物 理 量，而 概 率 地 震 危 害 分 析(Probabilistic Seismic Hazard Analysis，PSHA)是一种广泛应用的评估方法，可估计具有一定烈度的地震发生概率20。PSHA 给出了地震危害曲线的泊松超越概率(Pe)
10、。对工业场所的风险评估，50 年(储罐等主要设备的典型使用寿命)重现期较为常用。参考 Campedel 等21 的工作，只考虑50 年重现期，使用下面模型描述地震危害。Pe A =aebA(1)式中 是某给定的地震烈度;a 和 b 是两个依赖于局部地震活动性的常数。对应急救援行动的量化描述，主要是确定应急控制火灾所需时间 tttc。考虑到其存在不确定性，一般用对数正态分布描述19。ln(tttc)N(ttc，2ttc)(2)式中ttc和 2ttc分别是函数 ln(tttc)的均值和方差，与当地救援力量有关，min。1.2地震造成的主要事件场景评估此板块的目的是计算出地震可能造成的主要事件场景及
11、其发生频率。首先计算常压储罐在地震作用下的破坏概率。考虑到计算资源和时间，使用简单快捷的 Probit 模型来计算储罐在地震下的易损性。Probit 单位 Ys计算如下11。Ys=k1+k2ln(100A)(3)式中A 单位是 g;k1和 k2是 Probit 系数。对地震情景下的常压储罐，Probit 系数见表 111。计算出 Ys后，基于累积标准正态分布密度函数 可以得到储罐的破坏概率 Pd。Pd=(Ys 5)(4)储罐被地震破坏后，会根据自身和储存物质的特点，触发不同的主要事件场景，造成不同的后果。在高烈度地震作用下，生命线系统和安全屏障可能受到严重破坏，泄漏的易燃物质着火概率极高22。
12、根据前人研究，通过事件树方法，选择其中引起池火的路径，可得到某次泄漏事件后储罐处于池火状态的概率23。图 2 为地震作用下易燃易爆液体物质的事件树，本文只考虑用灰线标注的池火情况。通过蒙特卡洛模拟，可以得到所有可能出现的主要事表 1地震作用下常压储罐损坏的 Probit 模型系数11 Table 1Probit coefficients for atmospheric steel tankdamage due to earthquake11 储罐状态储存量Probit 系数k1k2锚固的接近满的2.431.54锚固的50%2.421.25非锚固的接近满的0.8331.25非锚固的50%0.83
13、1.250612Vol 23No 7安全 与 环 境 学 报第 23 卷第 7 期图 2易燃易爆液体物质发生泄漏后的事件树23 Fig 2Event tree for flammable and volatile liquidmaterials after loss of containment23 件场景 Sk及对应场景的期望频率 Pk。假设储罐区有 n 个储罐，则 Sk是 1 n 维的矩阵。Sk=(f1，f2，fi，fn)(5)式中fi(i=1，2，n)取值为1 和0，为1 时表示第 i 个储罐处于池火状态，为 0 时则表示不处于池火状态。1.3地震造成的火灾多米诺升级评估此版块针对图 1
14、 第三部分给定主要事件场景的多米诺效应进行概率计算，动态考虑过程协同效应和不确定性，适用于大量储罐中高级别多米诺传播的场景。本文考虑大型常压储罐，假设其中充装的可燃物质在储罐失效着火后，在应急救援控制之前都足够燃烧，不会出现燃料燃尽的情况。计算流程如图 3 所示，重点关注的是多米诺效应升级的判据与每次升级过程更新的关键物理量。如图 3 所示，针对 1.2 节辨识的每个主要事件场景，综合热辐射阈值与储罐状态升级 Probit模型(式(6)(9)、应急控制火灾时间与未着火储罐最小失效时间之间的关系作为升级判据。考虑到应急控制火灾时间和储罐受热辐射失效时间的不确定性，将事故频率的计算简化为确定性的蒙
15、特卡洛抽样模拟计算。此外，使用矩阵计算方法考虑升级过程的协同效应，能简单高效地进行计算24。在给定主要事件场景之后，输入储罐参数、环境参数和应急控制火灾的时间。M 为蒙特卡洛模拟总次数，m 为当前模拟次数，初值 m=1，L 为截止的多米诺传播级别，l 为当前多米诺传播级别。利用升级判据对未着火储罐判断多米诺级别是否升级，如果升级会导致储罐本身的状态变化及其他未着火储罐受到总热辐射量增大，从而更新储罐群的状态 S 和图 3基于蒙特卡洛模拟的动态多米诺概率计算流程Fig 3Flowchart of probability calculation of dynamicdomino effect ba
16、sed on Monte Carlo simulation受到的总热辐射量。可通过 S 与目标罐从着火罐接受到的热辐射量矩阵 0的关系求出24。0通过 ALOHA 软件计算得出。考虑到多米诺升级过程中，可能发生变化，需要对过程进行动态化评估。目前常用储罐状态升级Probit 模型并不能很好地适用于这种情况25。本文使用递归模型来适应储罐所受总热辐射量变化的情况，并对这种情况下 l 级传播储罐 j 剩余失效时间tl，jrttf进行修正，如下式所示。16122023 年 7 月陈海翔，等:地震诱发储罐火灾多米诺效应的动态定量风险评估方法研究Jul，2023tl，jrttf=e1.13ln(l，j)2.67103Vj+9.88/60tjd=0l，jl1，()j1.13(tl1，jrttf tl1)tjd0(6)为考虑随机性，再令 tl，jrttf=N(tl，jrttf，2rttf)，使该变量符合标准差为rttf的正态分布，避免相同条件储罐同时失效升级的情形，使其更符合实际。式(6)中，l，j表示第 j 个储罐在 l 级多米诺传播过程接受到的总热辐射量，Vj表示第 j 个储罐的体积，tl1为
