不确定条件下航空不安全事件灵敏度分析的M onte-Carlo方法

陈浩然，崔利杰，2，*，任博，2，张贾奎

（1.空军工程大学 装备管理与无人机工程学院，西安710051； 2.光电控制技术重点实验室，洛阳471000）

灵敏度分析（Sensitivity Analysis，AS）是研究一个系统的状态或输出变化率对系统参数或周围条件变化的敏感程度的方法［1-2］。航空器作为一个组件构成多、逻辑交互复杂、不确定因素广泛的巨系统，开展灵敏度分析可以有效指导安全性设计与维护工作，提高航空安全水平。

国内外学者针对航空领域的灵敏度分析方法进行了大量探索，Li等［3］通过多变量概率积分变换分析输入不确定性对多变量输出整个分布的影响，得到飞机机翼旋转轴模型的灵敏度指数；Cao等［4-5］通过建立飞机发动机高保真热力学模型，得到了发动机参数的全局灵敏度；Zentner等［6］在不确定条件下分析变量数据来源，得到了随机变量和主观变量对核涡轮机振动的灵敏程度；金燕和刘少军［7］结合神经网络和一次二阶矩法对航空滚动轴承的影响因素开展灵敏度分析；权凌霄等［8］基于ANASYS软件平台和多目标遗传算法分析了航空液压管路支架参数对振动响应的灵敏度；张马兰等［9］结合区间数学和贝叶斯神经网络方法，对航空公司安全管理体系的指标变量进行灵敏度分析；陈超［10］、锁斌［11］等通过模糊理论、概率包络对认知不确定性进行描述，提出了针对航空复杂系统的灵敏度计算方法。上述研究主要针对航空组件本身开展灵敏度分析，少有针对事件级进行深入研究，无法直观地表征组件及其分布参数对航空安全态势的影响程度。

由此，本文充分考虑随机事件和主观事件的不确定性，基于航空不安全事件Bow-tie模型，提出一种适合多输出条件的航空安全性函数与灵敏度指标，并采用Monte-Carlo方法开展灵敏度测算，运用实例验证所提指标的合理性和准确性，为提高航空安全水平提供参考。

1 典型航空不安全事件建模

航空安全分析的模型有很多，如事故树分析（Fault Tree Analysis，FTA）和事件树分析（Event Tree Analysis，ETA）是2个借助图形工具危险源和事故后果定量分析的有效工具，但两模型均无法直观建立危险源与事故后果的关系，而Bow-tie模型以顶事件为媒介连接危险源和事故后果，打通了事故原因和事故后果之间的联系，解决了传统模型条块分割严重、直观性不强的问题，是一种新型的安全分析工具。

1.1 Bow-tie模型构建

Bow-tie模型可分为两部分，左侧事故树是一个用逻辑符和事件符连接的由原因到结果的树形图，连结基本事件和顶事件，通过分析基本事件导致顶事件发生的逻辑关系有针对地制定防范措施；右侧事件树是一个由顶事件到后果事故的过程分析图，通过分析顶事件发生后不同的控制事件，演绎推断出所有事故后果。因此，Bow-tie模型将基本事件、顶事件、后果事件和控制事件统一起来，构成了一幅展示事故发生前因后果的可视化领结图，原理如图1所示。

图1 Bow-tie模型原理Fig.1 Principle of Bow-tie model

以往许多学者基于Bow-tie模型开展安全分析过程中，新的模型构建方法不断涌现［12-14］，总结其构造原则主要有以下3点：①模型中只有一个顶事件，是事故树的终点，也是事件树的起点；②模型中最左端是导致顶事件发生的危险源，最右端是不同控制措施下顶事件可能导致的后果事件；③模型中事故树的所有分支均向顶事件汇集，事件树由顶事件拓展至所有事故后果。

1.2 Bow-tie模型求解

由于Bow-tie模型基于事故树和事件树模型而来，其求解方法可在2类树形分析工具计算基础上进行，以图2简化的Bow-tie模型展示其求解过程。

式中：

因此，后果事件OEi的发生概率可表示为n个基本事件与m个控制事件发生概率的函数，即

图2 Bow-tie模型示意图Fig.2 Sketch map of Bow-tie model

2 安全性灵敏度分析方法

由于航空系统不确定性特征显著，传统量化指标往往难以测算，本文将基于Bow-tie模型提出新的航空安全性指标和灵敏度测度。

2.1 航空安全性指标

通常情况下，人们往往更加关注航空领域超过预期严酷度的事故，如人员伤亡、财产损失或环境污染等，因此，航空安全性指标可定义为：在预期坏境下，执行预定航空任务时导致危险后果的可能性低于人们预期值的概率［12］，即

航空活动中，后果事件发生概率阈值已被提前规定，如MIL-HDBK-882D［15］、ARP4754A［16］等标准已经规定了一系列航空事故后果事件严酷度阈值，在给定严酷度下的安全性功能函数为

因此，航空安全性指标可进一步转化安全功能函数的表示形式为

若基本事件之间相互独立，航空安全性指标的概率表达式为

为失效域指示函数；E［·］为数学期望算子；fX（x）表示基本事件和随机变量的联合概率密度函数。

若采用Monte-Carlo方法计算上述指标，则安全性指标可转换样本的估计值为

式中：N为仿真抽样次数。

在单一航空事故后果事件计算的基础上，考虑多个后果发生的可能性，得到多模式下航空器安全性指标计算公式为

综合式（7）、式（8），得到多模式下采用Monte-Carlo方法计算航空安全性指标的样本估计值为

2.2 航空器安全性灵敏度指标

传统的灵敏度测度需要获得事故发生概率的解析函数，逻辑推导困难且计算量大而繁琐，给灵敏度求解过程带来巨大困难。因此，寻求一种适合不确定条件下的灵敏度分析方法，对航空不安全事件量化分析工作具有重要意义。

1）全局灵敏度

全局灵敏度是从平均的角度来衡量输入变量的不确定性对输出的贡献，也被称为输入变量的重要性测度。参考航空领域重要性测度的概念，结合式（9）航空安全性指标，得到第i个基本事件的全局灵敏度Sx（xi）为

Monte-Carlo方法下全局灵敏度的样本估计值可参考式（9），只是第i个基本事件的发生概率有所变化，这里不再赘述。

2）局部灵敏度

为了便于运用Monte-Carlo方法开展抽样估计，将局部灵敏度表达式转换成均值形式为

式（12）可进一步转换成样本均值估计形式为

式中：xj为概率密度函数fX（x）的第j次抽样值。

若n个基本事件之间相互独立，由于基本事件和控制事件相互独立，则联合概率密度fX（x）为n个基本事件概率密度函数与m个控制事件概率密度函数的乘积。由于分布参数θik只与第i个概率密度函数fi（xi）有关，对式（13）进一步变换，得到局部灵敏度表达式为

大量故障统计信息显示［17］，若基本事件为机械件失效，其失效概率服从对数正态分布；若基本事件为电子件失效，则其失效概率服从指数分布。对于服从对数正态分布的失效件，显然有

将式（15）、式（16）代入式（14），得到安全性指标对第i个基本事件的均值μi和标准差σi的局部灵敏度分别为

同样思路，可以得到基本事件失效率服从指数分布时，安全性函数对第i个基本事件均值μi的局部灵敏度为

由于控制事件涉及大量主观因素，将在3.1节给出其概率测算方法，不再讨论其灵敏度测度。

3 航空不安全事件灵敏度求解

航空活动中危险变量的不确定特征显著，若要求解上述提出的安全性指标和灵敏度测度，须提出科学准确的变量分布特征描述方法。

3.1 变量描述

对数据库［17］中大量失效信息统计分析后发现，对于多数机械类和电子类基本事件，可采用概率方法计算其失效概率。而对于数据库中没有的基本事件，可根据其组件类型和故障机理，判断失效概率的分布类型，统计其分布参数求解基本事件发生概率。

而多数电子类产品的MTBF服从指数分布，假定故障率为λ，则分布参数μMTBF=1/λ，因此飞行时间为T时，电子类产品的概率密度函数为

3.2 M on te-Car lo方法仿真流程

根据所提各项指标，结合相关变量描述即可求解Bow-tie模型各类输出。Monte-Carlo方法是一种常用的数字模拟抽样方法，本文以该方法为例开展仿真分析，具体流程见图3。

图3 Monte-Carlo方法的航空安全灵敏度仿真流程Fig.3 Simulation flowchart of aviation safety sensitivity using Monte-Carlo method

由图3可知，由于基本事件失效多是与飞行时间相关，因此求解的航空安全性指标和灵敏度指标随着飞行时间动态变化；考虑到不同基本事件失效概率的差异性，航空器安全性变化率对不同基本事件敏感程度不同。另外，其他高效的数值仿真方法也适用于本文所提各个指标。

4 典型案例

轮胎是航空器起落系统的重要构件，对飞机安全起降有着至关重要的作用。飞机轮胎一旦爆破，往往将造成重大的财产损失或人员伤亡。根据轮胎爆破历史事故统计结果［18］，人员、机械、环境等因素是导致事故发生的主要原因，应对事故的不同控制措施也往往导致不同严酷程度的事故后果。很多学者已经对轮胎爆破事件进行了深入研究，这里基于Cui等［12］构建的Bow-tie模型，验证文中所提各类指标的可行性、合理性，模型中基本事件分布参数及其不确定描述情况见表1，控制事件服从区间分布的不确定性描述见表2。

表1 轮胎爆破基本事件分布参数及其类型描述Tab le 1 Distribu tion param eters and type descrip tion of basic even ts for tire burst acciden t

表2 轮胎爆破事故树控制事件不确定性描述Tab le 2 Uncertain ty descrip tion of control events for event tree of tire burst accident

在模型中，轮胎爆破是唯一的顶事件，共有16个基本事件可能致事故发生；轮胎一旦爆破，共有4种控制事件分别定义如下：SE1为启动应急刹车；SE2为避让飞机和建筑物；SE3为增设隔离网；SE4为启动应急消防措施。不同的控制事件组合由会导致4种不同严酷程度的后果事件：OE1（轻微的）表示停留在跑道上；OE2（较大的）表示机体轻度损伤；OE3（危险的）表示机体严重损伤；OE4（灾难的）表示飞机起火或人员伤亡，根据航空安全领域对于风险事件的通用界定方法［19］，其可接 受 概 率 阈 值 分 别 是10-3、10-5、10-7、10-9。利用表1和表2所提供数据，结合上述分析，计算得到航空安全性指标随时间变化情况见图4。

由图4可知，轮胎爆破事件航空安全性指标随时间不断减小，且在500～600 h发生突变，安全度急剧下降，此时是基本事件失效所致，因此需要特别关注该时间段各基本事件失效导致事故发生的可能性。考虑到灵敏度和安全水平的关系，这里选取450～650 h时间区间，计算得到航空安全性指标对基本事件的全局灵敏度（见图5），局部灵敏度见图6～图8。

图4 不确定条件下轮胎爆破事件航空安全性指标Fig.4 Aviation safety index of tire burst accident under uncertain conditions

图5 基本事件全局灵敏度Fig.5 Global sensitivity of basic events

图6 机械类基本事件分布参数μ局部灵敏度Fig.6 Local sensitivity for distribution parameter μofmechanical events

图7 电子类基本事件分布参数λ局部灵敏度Fig.7 Local sensitivity for distribution parameter λof electronic events

图8 机械类基本事件分布参数σ局部灵敏度Fig.8 Local sensitivity for distribution parameter σofmechanical events

图5显示，轮胎爆破事件的安全性指标对不同基本事件的全局灵敏度存在显著差异，且各个基本事件灵敏度指标随着时间发生变化显著，但各基本事件灵敏度的重要性排序保持不变。图中所有基本事件全局灵敏度均为负值，说明减小基本事件发生概率，能够有效提高航空器安全指数；BE6和BE4由于MTBF较小，是导致危险发生的主要因素，灵敏度指标甚至超过0.8；灵敏度较大的5种基本事件排序为：BE6＞BE4＞BE7＞BE3＞BE5，因此机务维护人员要保证刹车盘、热熔塞、充气嘴、轮胎气压和机轮磨损符合飞行要求。另外，还有一些基本事件对航空器安全性指标敏感程度较小，如全部的电子类事件和BE13～BE16等机械类事件。

图6和图7表明，航空安全性指标对机械类基本事件分布参数μ和电子类基本事件分布参数λ的局部灵敏度均为随时间动态变化的正值且差异巨大：为正值说明航空安全性指标与分布参数μ和λ正相关，增大基本事件分布参数均值能够有效提高航空安全水平；差异巨大是由于电子件和机械件事件MTBF相差太大，电子产品在450～650 h区间内几乎不会发生故障，导致电子类事件的均值灵敏度与机械类事件相比可以忽略不计。机械类基本事件的灵敏度重要性排序为：

BE6＞BE14＞BE4＞BE15＞BE13＞BE7＞BE16＞BE3＞BE5＞BE2＞BE1。另外，600 h后，机械类事件灵敏度整体都有所下降，此时航空安全性指标趋于稳定，而机械类事件的概率密度函数对均值μ的偏导开始下降，导致灵敏度有所下降。

图8显示，机械类事件分布参数σ对航空安全性指标的敏感程度排序与分布参数σ相同，这也恰恰验证了两者的一致性、合理性。另外，安全性指标对分布参数σ的敏感度为负值，说明要提高安全性指标要努力降低机械类基本事件的分布参数标准差σ。

5 结 论

1）本文从典型航空不安全事件的Bow-tie模型出发，提出了多模式航空安全性指标求解方法，进而推导出不确定条件下，针对不安全事件的全局灵敏度和局部灵敏度计算公式。

2）运用概率方法和区间理论对不确定变量进行描述，又采用Monte-Carlo方法对所提指标开展仿真计算，解决了航空器这类复杂系统难以进行灵敏度测算的困难。

3）结合典型航空事故——轮胎爆破事件开展仿真验证分析，通过与实际机务工作数据对比，证明了所提灵敏度指标的合理性、准确性和可操作性，对提高航空器的安全设计水平具有一定的指导意义。

本文主要针对单一顶事件开展安全性和灵敏度分析，而对于航空器这类复杂系统，往往包含大量顶事件，如何使用本文方法开展安全性评估和灵敏度分析有待进一步研究。