基于survival signature和子集模拟的复杂系统可靠性灵敏度分析

吴建新，李禹雄，杨国栋，黄贤振

(1.铁正检测科技有限公司，济南250014；2.东北大学机械工程与自动化学院，沈阳110819；3.中车青岛四方车辆研究所有限公司，青岛266031)

工业系统的规模和复杂程度随着科技的进步不断上升，对系统可靠性的要求也日益提高[1-3]。可靠性灵敏度定义为系统的可靠度对分布均值或标准差的偏导数，能直观地反映出系统中某元件的相关参数对整个系统可靠性的影响程度，有助于认清系统或机构的薄弱点[4，5]。由于现实中的系统连接结构复杂繁琐，元件数量多，传统的可靠性分析工具存在一定困难。Coolen团队在signature的基础上，针对多种分布元件的复杂系统提出了survival signature的概念，省略了推导系统结构函数的复杂过程，大大提高可靠性分析的效率[6，7]。当前主流可靠性灵敏度分析主要分为解析法和数字模拟，较为常用的如一次二阶矩法[8]、蒙特卡罗法[9，10]等应用在复杂系统中会受到一定的局限。一次二阶矩法的计算效果会受到系统功能函数形式的限制，对于非线性程度较高的功能函数往往精确度不够理想；传统的蒙特卡罗法虽然能够通过简易的编程获得灵敏度，但却需要生成大量的样本才能达到所需的精度，计算效率低，在工程实际应用中受限。Au等[11，12]提出的子集模拟法通过将小概率事件转化为多个条件概率的连乘，有效地减少采样次数，显著提高蒙特卡罗仿真的效率。本文将结合survival signature和子集模拟法，提出一种针对复杂系统的元件可靠性灵敏度分析方法，大大简化系统连接结构的分析，显著减少采样次数，提高模拟速度。文中的数值算例验证了该方法的计算效率和精度。

1 基于survival signature的系统可靠性分析

已知一个系统中包含m 个元件，向量x ＝(x1，x2，…，xm) 表示系统中元件某一时刻的工作状态，其中xi＝1或0分别代表第i个元件正常工作或失效。定义结构函数φ x( ) 值为1或0分别代表系统正常或失效。现假设某系统由K(K ≥2)种元件组成，mk为第k 种元件的总数。同种元件服从独立同分布，且假设各个元件的失效是相互独立的。定义向量x ＝(x1，x2，…，xK) 为系统所有元件的工作状态，其中，xk＝代表第k 种元件的工作状态。当第k 种元件中有lk个正常工作时(k＝1，2，…，K )，survival signature可写作[6]

其中，S l k代表已知各类元件正常工作数为l1，l2，…，l K的所有可能的状态向量集合Φ (l1，l2，…，l K)本质是条件概率，代表系统中各类元件分别有l1，l2，…，l K个保持正常工作状态时，系统能够正常工作的概率。因此，根据全概率公式，系统某一时刻的可靠度可表示为[6]

其中，P{C k(t)＝l k} 为t时刻下，第k种元件中有l k个元件正常工作这一事件的概率。若在t时刻第k种元件的可靠度为R k(t)，得

2 子集模拟法

子集模拟法通过引入中间事件，将目标小概率事件转化为多概率乘积

子集模拟的原理如图1所示，通过引入一系列临界值确定每一层中间失效域F，通过对中间失效域的模拟最终确定目标事件的失效概率。由于每一层的条件事件概率都远大于目标事件概率，因此模拟每一层的条件概率时所需的样本点数大大降低，从而显著提高模拟效率。生成中间条件样本时采用了Metropolis-Hastings算法[13]:选取特定的建议分布生成样本的备选状态并按照一定的概率接受或拒绝该状态，之后按照是否落入失效域用备选状态取代原状态，生成服从平稳的目标分布的马尔可夫链，得到所需样本。

中间概率事件的个数和大小影响着子集模拟法的效率。中间概率较大、中间事件较多时，模拟总层数会随之增加；反之，总模拟层数减少，但每层需要更多的样本点进行模拟[12]。实际操作中，通常将中间条件概率设定为0.1～0.3的固定值，根据中间概率值选择合适的每层失效域的阙值和总模拟层数。

3 基于子集模拟和survival signature的可靠性灵敏度模拟

3.1 基于子集模拟求解各元件可靠度

为了分析系统的可靠度和各元件的可靠性灵敏度，首先需要得到各元件的可靠度。元件在投入使用的初期，其失效概率较小，因此采用子集模拟法对各元件失效率进行分析来提高计算效率。通过子集模拟法模拟t时刻某一元件的可靠度流程如下:

(1)确定中间条件概率p0、最大模拟层数、功能函数g (T )＝T－t和失效事件F ＝{T:g (T ) ＜0} ；

(2)在第一层模拟中，生成N 个服从元件寿命分布的随机样本，按功能函数值升序排列。取前Np0个样本落入失效域，并根据第Np0－1个样本点确定第一层失效事件F1＝{T:g (T ) ＜b1} 的阙值b1；

(3)对第i层进行模拟(i≥2)，以落入上一个失效域的Np0个样本点为种子，根据M-H 算法，生成服从本层条件概率分布的样本，使样本总数再次达到N ；按功能函数值升序排列，取前Np0个样本作为落入失效域的样本，并根据第Np0－1个样本点确定第i层失效事件Fi＝{T:g (T ) ＜bi} 的阙值bi；

(4)重复步骤(3)，直到第m 层时的阙值bm＜0，此时令bm＝0。 统计落入失效域的样本点个数Nm，此时目标事件的失效概率即为:

3.2 系统可靠度

模拟仿真过程中，服从同种寿命分布的一类元件在t 时刻的可靠度模拟值也不会完全相同，因此，不能直接使用式(3)求解P {Ckt( )＝lk}。 由于上文已通过子集模拟法求得元件在t 时刻的失效率Fi，从而可知元件可靠度Ri＝1－Fi，可采用排列组合的方法就可求得某一类元件中有若干元件正常工作的概率。举例来说，假设第k 种元件共5 个，k 种元件中有3 个正常刚工作的概率可以表示为:P {Ckt( )＝3}＝R1R2R3(1－R4) (1 －R5)＋…＋(1 －R1) (1－R2)R3R4R5，将结果带入式(2)中即可得到系统在t时刻的可靠度。

3.3 元件可靠性灵敏度

元件的可靠性灵敏度定义为系统的可靠度对随机参数均值或标准差的偏导数

其中，RS(t)和Ri(t) 分别代表在t时刻系统和第i个元件的可靠度；θi代表第i个元件的寿命分布参数，即分布均值μ 或标准差σ。 偏导数ƏRS(t)/ƏRi(t) 表示第i个元件的Birnbaum 重要度(Birnbaum importance)[14]。根据全概率公式，系统在t时刻的可靠度可以写作

其中，h (1i，R (t)) 和h (0i，R (t)) 分别为第i个元件正常工作和第i个元件完全失效条件下原系统的可靠度。对等式两边求偏导可得到重要度的计算公式

根据式(7)求解每一个元件的重要度时都需要计算一次h (1i，R (t)) 和h (0i，R (t)) ，较为繁琐，为简化计算，可以将式(6)带入式(7)中，得到元件重要度的另外一种表达方式

其中，RS(t)和Ri(t) 已经求得，只需额外求出h (0i，R (t)) 就能得到第i 个元件的重要度。找到原系统第i个元件失效条件下的等效系统，求解该系统的可靠度即为h (0i，R (t)) 。 例如在图2中，在包含5个元件的桥式系统中求解h (03，R (t)) ，只需将元件3视为断路，将原系统转化为4元件的串并联系统并求解该系统的可靠度即可。求解等效系统可靠度同样可采用survival signature。根据前文可知，若系统由K 种元件组成(K ≥2)，mk为第k 种元件的总数，且假设所分析的元件i属于第h 种元件，则

其中，Φ＇l(1，l2，…，lK) 代表元件i失效时的等价系统的survival signature。将上式带入式(8)中即可得到元件重要度。偏导数表示第i个元件可靠度对分布均值或标准差的偏导。根据元件分布寿命表达式，可以将偏导数转化为通过可靠度R i(t)表达的形式。以Weibull分布为例来说明转化过程。服从Weibull分布的元件在t时刻的可靠度为:，其中a为比例参数，b为形状参数。分别对两个参数求偏导，得

将上式进行一定变化，得

由于已知Weibull分布的均值和标准差为

其中，Γ(x) 为伽马函数，即欧拉第二积分，是阶乘函数在实数与复数上扩展的一类函数，得

可求出均值μ和标准差σ分别对于比例参数a和形状数b的偏导并整合为矩阵。根据Jacobian矩阵的性质，对该矩阵求逆即可得到a、b分别对μ和σ的偏导[15]

接着通过下式即可得到基于元件可靠度对均值μ和标准差σ的偏导

至此，元件的可靠性灵敏度已完全转化为元件可靠度的表达式。

4 算例分析

某系统共两类元件，其中第1类元件(元件1、元件2、元件5)服从a＝3、b＝4二参数Weibull分布；第2类元件(元件3、元件4、元件6)服从a＝4、b＝2的二参数Weibull分布。通过式(2)计算该系统的survival signature，结果如表1。系统连接结构如图3所示。

表1 算例系统的survival signature

首先通过子集模拟法得到各元件的可靠度，通过3.2所述方法得到系统基于时间的可靠度曲线如图4所示，且图像显示可靠度的计算结果与实际解析值基本拟合。分别计算出各元件失效时等效系统的survival signature如表2～表7所示。

表2 元件1失效时的等效系统

表3 元件2失效时的等效系统

表4 元件3失效时的等效系统

表5 元件4失效时的等效系统

表6 元件5失效时的等效系统

表7 元件6失效时的等效系统

将各等效系统的survival signature和各元件的可靠度带入式(9)，可以得到各元件的重要度(Birnbaum importance)曲线如图5 所示。元件1的重要度要高于其他元件，说明其在系统中具有关键作用。但需注意重要度只能反映元件所处位置对系统可靠度的影响，并没有考虑元件自身寿命分布带来的影响，因此须进一步求解灵敏度数值。由于已知元件可靠度服从Weibull分布，因此将式(10)带入式(5)，即可得到各元件的可靠性灵敏度曲线如图6和7所示。

由图6和7中的曲线可知，相比其他元件，元件3和6的寿命分布均值和方差的灵敏度绝对值更大，因此其对于系统可靠性的影响程度是最大的，在进行可靠性优化等工作时应对其进行着重考虑。选取t＝0.5时各个元件的失效率模拟效果来说明此方法的效率。具体数据见表8。

表8 子集模拟法的精度和采样效率

由表8可知，子集模拟法的采样效率要远高于直接蒙特卡罗，其模拟精确度仍然较为理想，即使适当提高样本点数量或进行重复试验，其采样总数仍然少于直接蒙特卡罗。直接蒙特卡罗所需的样本点会随着目标概率的缩小而呈指数增长，因此子集模拟法在处理小量级的事件概率时，计算效率的优势会更加明显。

5 结论

本文结合了survival signature和子集模拟法，针对复杂系统提出了一种新的系统可靠性和元件灵敏度模拟方法。survival signature简化了系统的可靠性分析，避免了传统方法中求解系统结构函数的复杂运算，对于结构复杂元件种类繁多的系统尤为有效；而子集模拟法则大大减少了所需的随机样本点数量，在保障足够精度的前提下提高了模拟的效率。通过包含多元件的复杂系统算例计算，验证了该方法的精度和效率，体现了良好的适用性。