加速寿命试验设计与评估软件ALT511及其应用(二)

2024-01-08 10:13黄首清刘庆海李芳勇张兆霖
航天器环境工程 2023年6期
关键词:幂律指数分布概率密度函数

黄首清,刘庆海*,李芳勇,杨 勇,遇 今,张兆霖

(1.航天机电产品环境可靠性试验技术北京市重点实验室; 2.可靠性与环境工程技术重点实验室;3.北京卫星环境工程研究所:北京 100094)

0 引言

单应力多个应力水平加速寿命试验是经典的加速寿命试验之一;其可根据所设计的多个应力水平的加速寿命试验失效时间数据,在失效机理不变的前提下外推其他应力水平的寿命。单应力多个应力水平加速寿命试验涉及到应力水平的数量和大小的选取、各应力水平寿命分布拟合、失效机理一致性判别、加速模型的选用和拟合、加速模型与寿命分布模型的关联,以及寿命外推和可靠寿命计算等多个环节,对相关技术和工程经验要求较高,因此,简化应用复杂度以及开展软件设计研究是推广该技术的重要途径和研究方向。

在单应力多个应力水平加速寿命试验领域的相关研究主要聚焦在试验的精细化设计和数据的高精度分析[1-2]。早期的研究集中在应力水平和样本量的选择。例如:Meeker 等[3]认为对于3 个应力水平的情况,从低应力水平到高应力水平的样本量比例应为4:2:1,且中等应力水平应为低应力水平和高应力水平的平均值;而中国国家标准GB 2689.1—1981 [4]规定,在缺少加速因子相关历史数据的情况下,需要设计多个应力水平的加速寿命试验,应力水平一般不少于4 个,每个应力水平下的样本量不少于5 个;Yang[5]提出了如何设计寿命满足威布尔分布的4 个应力水平加速寿命试验;一些学者基于蒙特卡罗方法提出了应力水平、样本量等设计要素的优化设计方法[6-7]。针对加速寿命试验数据,主要采用统计学方法进行分析,例如图估计法[8]、简单线性无偏估计法[9]、最好线性无偏估计法[10]、最小二乘法[11]、极大似然估计法[12]和贝叶斯法[13]。为了进一步减小加速寿命试验的分析误差,一些方法被用来对现有模型进行修正。Mazzuchi 等[14]提出采用动态线性模型与贝叶斯综合方法来降低分析结果误差;Watkins[15]提出一种修正的极大似然估计法,降低了计算的复杂度;殷毅超[16]提出了基于不精确概率理论的加速寿命试验数据分析方法。从工程应用的角度分析,一些传统的方法辅助工程经验可以较好地解决单应力多个应力水平加速寿命试验设计和评估问题,而一些改进和优化方法仍停留在研究层面或者因求解方法过于复杂而难以推广应用。

在软件工具方面,美国的ReliaSoft 软件[17]、JMP Pro 软件[18]、Minitab 软件[19]和中国的ReliaQube软件[20]均具有加速寿命试验模块,可以满足单应力多个应力水平加速寿命试验设计和分析需求,但在适用性和易用性方面存在一些不足。以Minitab 软件为例,它包含的加速模型仅有线性模型、自然对数模型等少数几种,不包括阿伦尼斯模型等常用模型,并且在结果的可视化显示方面存在不足。除上述软件外,我国还有零星的用于研究目的的加速寿命试验软件,但其研发和应用尚处于起步阶段,例如:重庆工业自动化仪表研究所开发的加速寿命试验数据分析软件[21]具备寿命分布拟合和加速模型参数估计等功能;电子科技大学开发了加速寿命试验方案优化设计系统[22],在试验条件设计方面很有特色,并融入了自主提出的改进的随机游走算法。然而,这些软件在可视化显示、工程经验融合、分析结果验证等方面还存在一些不足,尚无实际的工程应用案例。

本文将从软件实现、计算方法和实例应用的角度,对北京卫星环境工程研究所自研的ALT511软件[23]中单应力多个应力水平加速寿命试验设计与评估功能进行系统和全面的阐述,以期对我国加速寿命试验软件的设计和应用提供有益参考和支撑。

1 单应力多个应力水平模块介绍

1.1 “导入数据”选项卡

如图1 所示,“导入数据”选项卡主要用于读入各个应力水平下的失效时间数据;点击下方的“导入数据”按钮会弹出选择数据文件窗口,选定目标数据文件后,数据将以表格形式进行可视化呈现。本软件可处理4 组应力水平数据,并可兼容3 组和2 组应力水平的情况。由于航天器产品通常价格昂贵且样本量有限,所以多个应力水平的加速寿命试验直接应用于单机产品的情况并不多见,它更适用于涂层等材料以及轴承、电机等零件。应力水平建议不少于4 个,每个应力水平下的样本量建议不少于5 个。

图1 “导入数据” 选项卡界面Fig.1 Interface of the tab “data import”

1.2 “每组应力水平寿命分布”选项卡

如图2 所示,“每组应力水平寿命分布”选项卡主要针对指数分布和威布尔分布进行数据拟合并绘制概率密度函数;可针对指数分布给出各个应力水平下的失效率,针对威布尔分布给出各个应力水平下的形状参数和特征寿命。

图2 “每组应力水平寿命分布”选项卡界面Fig.2 Interface of the tab “life distribution under each stress level”

分布类型通过界面左上方的下拉菜单选择;根据选择的分布类型,界面自动更新显示相应的曲线和寿命分布模型参数。建议偶然性失效的产品(例如星载电源产品)选用指数分布;损耗性失效的产品(例如控制力矩陀螺轴承、舱外电动工具按钮、超声电机、太阳帆板驱动机构导电滑环等)选用威布尔分布。

1.3 “加速模型拟合”选项卡

“加速模型拟合”选项卡可选择逆幂律模型和阿伦尼斯模型2 种加速模型;针对指数分布可对不同应力水平的平均故障前工作时间(MTTF)进行拟合建模,针对威布尔分布可对不同应力水平的特征寿命进行拟合建模,最终将拟合结果绘制成曲线;绘图界面中同时显示拟合点和拟合线,同时根据最小二乘法计算并给出拟合相关系数,并显示加速模型公式以及公式中模型常数的计算结果。需要指出,本软件绘制加速模型曲线共4 种情况,即2 种寿命分布和2 种加速模型的排列组合,分别是指数分布逆幂律加速模型、指数分布阿伦尼斯加速模型、威布尔分布逆幂律加速模型和威布尔分布阿伦尼斯加速模型(如图3 所示)。

图3 “加速模型拟合”选项卡界面Fig.3 Interface of the tab “acceleration model fitting”

1.4 “指定应力水平下寿命值”选项卡

“指定应力水平下寿命值”选项卡主要用来计算可靠寿命、失效概率密度函数和可靠度函数。如图4 所示,在界面中输入指定应力水平、可靠度、置信度,然后在下拉菜单中选择是否考虑置信度,即可给出可靠寿命,并自动绘制或更新对应的概率密度函数曲线和可靠度函数曲线。

图4 “指定应力水平下寿命值”选项卡界面Fig.4 Interface of the tab “reliable life value under specified stress level”

2 计算方法

2.1 每组应力水平寿命分布

1)指数分布

对于指数分布,主要进行失效率参数估计,进而得到失效概率密度函数[24]。每组应力水平下失效率可根据极大似然估计法计算,

式中:下标i指第i个应力水平;Ti为累计试验时间;ri为失效数; λˆi为失效率估计值。平均失效前工作时间为

第i个应力水平下失效概率密度函数为

2)威布尔分布

对于威布尔分布,主要进行形状参数和特征寿命两个参数的估计,进而得到失效概率密度函数[24]。对于失效数≥3 的情况,按照极大似然估计法通过求解方程(4)计算每组应力水平下的形状参数。

式中:mi为形状参数;ri为失效数;tij为第j个试验件的试验时间;ni为样本量。

对于失效数<3 的情况,根据工程经验确定形状参数mi。第i个应力水平下特征寿命估计值 ηˆi按照式(5)和式(6)计算:

第i个应力水平下失效概率密度函数为

2.2 加速模型拟合

以下分别针对阿伦尼斯模型和逆幂律模型介绍如何拟合不同应力水平下的寿命。

1)阿伦尼斯模型

用阿伦尼斯模型表达寿命与应力之间的关系:

式中:a和b为模型常数;S为应力;L为应力水平是S时的寿命。这里的寿命L对于指数分布为平均失效前工作时间θ,对于威布尔分布为特征寿命η。对式(8)等号两边取对数进行线性化处理得

式中:Y=lnL;A=lna;B=b;X=1/S。根据最小二乘法,可得常数B和A的估计值分别为:

其中:Xi和Yi分别为第i个应力水平下的X和Y取值;X¯ 和Y¯分别为X和Y的平均值;k为应力水平的数量。根据计算出来的Aˆ 和Bˆ可以反求加速模型的常数a和b。拟合相关系数ρ为

其中,ρ越接近1,代表选择该加速模型与试验数据的匹配程度越高。

2)逆幂律模型

用逆幂律模型表达的寿命与应力之间的关系为

式中:c和α为模型常数。同样,这里的寿命对于指数分布为平均失效前工作时间 θˆi,对于威布尔分布为特征寿命 ηˆi。对式(13)等号两边取对数同样可得到形如式(9)的表达式,其中,Y=lnL,A=lnc,B=-α,X=lnS。类似阿伦尼斯模型拟合,可按照式(10)和式(11)分别计算常数B和A,并反求加速模型的常数c和α。进而可根据式(12)计算拟合相关系数,评估逆幂律模型与试验数据的匹配程度。

2.3 指定应力水平下寿命值

根据2.2 节得到的加速模型,可计算指定应力水平下的寿命。这里的寿命对于指数分布为平均失效前工作时间θ,对于威布尔分布为特征寿命η。以下介绍在不考虑置信度和考虑置信度下,如何得到失效概率密度函数、可靠度函数和可靠寿命。

1)不考虑置信度

在不考虑置信度的情况下,可以得到指数分布的失效概率密度函数f(t)和可靠度函数R(t)表达式,分别见式(14)和式(15),并可绘制曲线。

威布尔分布的失效概率密度函数f(t)和可靠度函数R(t)表达式分别见式(16)和式(17),并可绘制曲线。

指数分布和威布尔分布的可靠寿命tR表达式可分别由式(15)和式(17)推导得出,分别为

2)考虑置信度

在考虑置信度的情况下,对于指数分布,采用式(20)[24]得到平均失效前工作时间的置信下限θL:

图5 置信寿命系数与失效数、置信度、形状参数的关系Fig.5 Impact of failure numbers, confidence level and shape parameter on confidence life factor

考虑置信度的失效概率上限密度函数fu(t)和可靠度下限函数RL(t)表达式分别为式(21)和式(22),可靠寿命下限tRL表达式为式(23)。

对于威布尔分布,采用式(24)得到特征寿命的置信下限ηL:

对比式(24)和式(20),可以发现指数分布是威布尔分布形状参数m¯取1 的特殊情况。将式(24)中ηL/η记为威布尔分布的置信寿命系数,表征考虑置信度下特征寿命的降低程度,从图5(b)~(d)可以看出:置信寿命系数随着失效数的增加而增大、随着置信度的增加而减小,这一规律与指数分布的一致;并且随着m¯增加,置信寿命系数显著降低。

类似指数分布,可以得到考虑置信度的失效概率上限密度函数fu(t)和可靠度下限函数RL(t)表达式,分别见式(25)和式(26);可靠寿命下限tRL表达式见式(27)。

3 软件应用案例

以某航天器超声电机为对象,给出单应力多个应力水平加速寿命试验设计和评估的软件应用案例。

根据失效模式与影响分析(FMEA),以负载为加速应力,设计比正常应力水平0.02 N·m 高且比破坏极限0.22 N·m 低的4 个应力水平开展加速寿命试验并记录试验件的失效时间。如图6 所示,将0.03 、0.07 、0.16 和0.21 N·m 共4 个应力水平下共20 个试验件的失效数据读入软件并显示。

图6 某航天器超声电机在4 个应力水平下的失效时间Fig.6 Failure time of a spacecraft ultrasonic motor under four stress levels

又根据FMEA,该航天器超声电机在负载应力作用下主要失效机理为磨损失效,具有典型的累积损伤失效特征,这样各个应力水平下的寿命应满足威布尔分布。选择分布类型后,软件可自动计算4 个应力水平下的形状参数、特征寿命,并绘制失效概率密度函数曲线(如图7 所示)。以应力水平1(0.03 N·m)为例,寿命满足形状参数为5.1、特征寿命为5941 h 的威布尔分布。此外,观察4 个应力水平的寿命分布,随着负载的增大,特征寿命逐渐降低,说明加大负载能很好起到缩短试验时间的加速效果;而4 个应力水平的形状参数相差不大,提示各个应力水平的失效机理基本一致,各组应力水平下的试验数据有效,可用于建立加速模型。

图7 4 组应力水平下的某航天器超声电机寿命分布Fig.7 Life distribution of a spacecraft ultrasonic motor under four stress levels

基于4 个应力水平下的特征寿命数据进行加速模型拟合,结果如图8 所示,可见:当选择阿伦尼斯模型时,算得模型常数a和b分别为32.22 和0.171,相关系数为0.894 4;而选择逆幂律模型时,算得模型常数c和α分别为0.519 3 和2.776,相关系数为0.971 9。由于逆幂律模型的相关系数更接近1,说明数据拟合的效果更好,所以选择逆幂律模型作为航天器超声电机的加速模型。这一加速模型表达式与图9 所示Excel 软件的拟合结果一致,证明加速模型拟合算法是正确的。但Excel 软件没有针对阿伦尼斯模型的拟合功能,Minitab 软件的加速寿命试验模块也未提供阿伦尼斯加速模型和逆幂律加速模型的拟合功能,这从一个侧面反衬出本软件的专业性。

图8 某航天器超声电机加速模型拟合Fig.8 Acceleration model fitting for a spacecraft ultrasonic motor

图9 Excel 软件拟合的逆幂律加速模型Fig.9 Inverse power law acceleration model fitted by Excel

获得加速模型后,可得到指定应力水平、可靠度和置信度要求下的可靠寿命。如图10 所示,假设0.02 N·m 为正常应力水平,这一负载下,不考虑置信度时可靠度0.95 对应的可靠寿命为15 180 h,这一计算结果与图11 所示的Minitab 软件计算结果几乎一致。

图10 某航天器超声电机可靠寿命计算Fig.10 Calculation of the reliable life of a spacecraft ultrasonic motor

图11 Minitab 软件计算的可靠度0.95 的可靠寿命Fig.11 Reliable life under 0.95 reliability calculated by Minitab

此外,本软件提供了更丰富的概率密度函数曲线和可靠度函数曲线显示功能,例如:设置置信度0.7,考虑置信度条件下,本软件计算的可靠度0.95 对应的可靠寿命为14 760 h。加速寿命试验中各试验件最长试验时长为7100 h、最短试验时长12.5 h,均远小于要验证的寿命指标要求。而且完成多个应力水平的加速寿命试验后,后续针对同类产品可以只做一个短耗时的高应力水平加速寿命试验,即可验证寿命指标。

4 结束语

本文介绍了加速寿命试验设计与评估软件ALT511 的单应力多个应力水平模块,包括:详细说明了“导入数据”、“每组应力水平寿命分布”、“加速模型拟合”和“指定应力水平下寿命值”4 个选项卡;给出了基于多应力水平失效数据的指数分布和威布尔分布拟合方法、阿伦尼斯模型和逆幂律模型的加速模型拟合方法,以及指定应力水平下考虑和不考虑置信度两种情况的可靠寿命、概率密度函数和可靠度函数计算方法。以某航天器超声电机加速寿命试验为例,软件应用结果表明:ALT511 软件的加速模型拟合结果与Excel 软件拟合结果一致,逆幂律模型常数c和α分别为0.519 3 和2.776,相关系数为0.971 9;可靠寿命计算结果与Minitab 软件计算结果几乎一致。

后续还将进一步研究多应力加速寿命试验、加速常数灵敏度分析、加速模型可信度分析,以及这些技术的软件实现。

猜你喜欢
幂律指数分布概率密度函数
幂分布的有效估计*
已知f(x)如何求F(x)
指数分布抽样基本定理及在指数分布参数统计推断中的应用
四川地区降水幂律指数研究
幂律流底泥的质量输移和流场
二元Weinman型指数分布随机变量之和、差、积、商及比率的分布
基于概率密度函数的控制系统性能评价
非高斯随机分布系统自适应控制算法的研究
幂律谱模型原子钟钟差仿真与噪声类型辨识
基于Fibonacci法求幂律模式流变参数最优值