CN112069697A - 多应力多失效模式相依竞争条件下加速退化试验统计分析方法 - Google Patents

Abstract

本发明提供了多应力多失效模式相依竞争条件下加速退化试验统计分析方法，包括：获取产品多应力加速退化试验数据；选择合适的分布模型并估计分布参数；建立退化量分布参数与时间之间的模型并估计模型参数；构造多应力加速模型；多失效模式相依竞争可靠度建模；求解模型未知参数，得到产品可靠度分布。本发明在考虑了多应力对寿命特征的耦合影响，退化过程与随机冲击过程的相依性，建立了N应力M失效模式相依竞争条件下的可靠度函数，填补了相关研究领域的空白；基于产品退化量分布和多应力多失效模式相依竞争可靠度模型，首次系统建立N应力M失效模式相依竞争加速退化试验的分析方法，更加符合工程实际，提高了产品寿命和可靠性评估的准确性。

Description

多应力多失效模式相依竞争条件下加速退化试验统计分析 方法

技术领域

本发明属于可靠性工程技术领域，具体涉及多应力多失效模式相依竞争条件下加速退化试验统计分析方法。

背景技术

在时间和成本约束下，基于加速退化试验的寿命预测技术已成为分析解决高可靠长寿命产品可靠度的必然要求。然而，产品在实际服役时通常会受到多种应力的作用，包括工作应力(如机械载荷、电流、电压等)和环境应力(如温度、振动、冲击等)。故施加单一应力的传统加速退化试验不能够真实体现产品实际的应力状态，研究多应力综合作用下产品的加速退化试验统计分析方法，可以更加全面体现产品的服役工况，使寿命预测更加准确。

产品在多应力下的可靠度建模是进行多应力加速退化试验统计分析方法的关键。但是，在多应力综合作用下，通常伴随有多个失效模式(如磨损、侵蚀和疲劳等)，且其失效是多失效模式相依竞争的结果。现有的方法主要集中于单应力作用下多失效模式相依竞争的可靠度建模，而极少有研究多应力作用下多失效模式相依竞争的可靠度建模方法。

竞争失效的一个典型工程实例是电液伺服阀受到磨损(或软失效)和随机冲击载荷(或硬失效故障)过程，这也被称为退化阈值冲击模型。目前，许多学者研究了随机冲击-退化之间的单侧依赖性，即随机冲击过程会导致增加退化量或增加退化率。然而，在工程上退化与冲击过程不仅仅是单侧依赖，而是相互依赖，即每一个随机冲击载荷都会导致退化量或退化率的增加，而冲击强度也在退化的累积下逐渐变大。

因此，亟需一种多应力作用下多失效模式相依竞争条件下的加速退化统计试验分析的新型技术。

发明内容

本发明针对当前多应力加速退化试验统计分析的需要，提出多应力多失效模式相依竞争条件下加速退化试验统计分析方法，所述方法包括：

获取产品多应力加速退化试验数据；所述多应力加速退化试验数据为多个产品在多个应力组合多个退化失效模式测试时间下进行加速退化实验得到的；

根据所述退化试验数据，确定所述退化试验数据对应的分布模型，以及获取所述分布模型的参数向量{ρ_i1m,ρ_i2m,...,ρ_iKm}；其中，ρ_iKm表示第i个应力组合下第m个失效模式的第K个测试时间点的模型参数向量；

根据所述参数向量，确定产品的寿命特征β_im，其中，β_im表示产品在第i个应力组合下第m个失效模式的寿命特征；

根据所述应力组合中包含的冲击应力，构建时间尺度协变量因子，根据所述时间尺度协变量因子和所述寿命特征，构建产品多个应力组合下的多应力加速模型；

在加速退化试验的随机冲击过程中，根据作用于产品的所述冲击应力的强度函数，确定产品第m个失效模式下的总退化量，以及获取冲击应力的冲击幅值；

根据所述冲击幅值达到预先设置的硬失效阈值时确定产品失效的第一失效关系，以及根据所述总退化量达到预先设置的软失效阈值时确定产品失效的第二失效关系；

根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型；

根据所述退化试验数据采用最小二乘算法求出多应力加速模型中的未知参数；

将所述未知参数输入到所述可靠度模型中，输出产品在多应力多失效模式相依竞争条件下的可靠度分布。

进一步的，根据所述应力组合中包含的冲击应力，构建时间尺度协变量因子，根据所述时间尺度协变量因子和所述寿命特征，构建产品在多个应力组合下的多应力加速模型，包括：

产品在多个应力组合下的多应力加速模型表示为：

其中，α₀,α₁,...α_j,...α_u,α_u+1,α_u+2为多应力加速模型的参数；

ξ_i＝(ξ₁,ξ₂,...ξ_N-1)为多应力水平组合S_i的标准化环境应力向量，ξ_j为某一标准化环境应力(1≤j≤N-1)，如标准化温度应力、湿度应力、冲击应力等；Y_j表示第j次冲击产生的退化量；L(t)表示t时间之内的冲击个数。

进一步的，在加速退化试验的随机冲击过程中，根据作用于产品的所述冲击应力的强度函数，确定产品第m个失效模式下的总退化量，以及获取冲击应力的冲击幅值，包括：

产品第m个失效模式下的总退化量X_Sm(t)表示为

其中，X_m表示产品第m个失效模式下的理论退化量；

表示产品第m个失效模式下的初始退化量；θ_m表示产品第m个失效模式下的退化速率，L(t)表示t时间之内的冲击个数。

进一步的，根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型中，

所述第一失效关系和第二失效关系的相依竞争为：产品在多应力作用下发生性能退化，会增加随机冲击的强度；而在随机冲击载荷作用下，既会增大产品的退化特征量，也会同时增大产品的退化速率。

进一步的，根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型，包括：

所述可靠度函数表示为：

其中，D为硬失效阈值，F_W(D)为W_j的累积分布函数，W_j为第j次随机冲击的幅值，j＝1,2,…L；θ_m表示第m个失效模式下的退化速率；λ_j(t)表示时间t内已承受j个随机冲击的强度函数；η为强度因子；λ₀(t)是时间t内已承受0个随机冲击的强度函数；λ₀为初始冲击的初始强度；χ为相依因子；f_θm表示第m个退化模式下的概率密度函数；μ_Y,表示随机冲击载荷所产生退化损伤的均值，σ_Y表示随机冲击载荷所产生退化损伤的标准差。

进一步的，根据所述退化试验数据采用最小二乘算法求出多应力加速模型中的未知参数，包括：

设定产品在试验过程中所承受冲击载荷的个数为L，时间分别为t₁,t₂,…t_L；

(a)、基于0冲击载荷过程退化数据(τ_ik,y_ijkm)，其中τ_ik<t₁；此阶段不包含冲击过程，因此，基于多个应力水平下的分布参数变化速率与多应力水平，采用最小二乘算法求取未知参数α₀,α₁,α₂,…,α_u；

(b)、将步骤(a)所求取的α₁,α₂,…,α_u带入所述多应力加速模型，基于多个冲击载荷过程退化数据(τ_ik,y_ijkm)，其中t₁<τ_ik<t_s，采用最小二乘算法求取未知参数α_u+1,α_u+2。

进一步的，所述退化试验数据服从正态分布。

本发明还提供了多应力多失效模式相依竞争条件下加速退化试验分析装置，所述装置包括：

数据获取模块，用于获取产品多应力加速退化试验数据；所述多应力加速退化试验数据为多个产品在多个应力组合多个退化失效模式测试时间下进行加速退化实验得到的；

数据处理模块，用于根据所述退化试验数据，确定所述退化试验数据对应的分布模型，以及获取所述分布模型的参数向量{ρ_i1m,ρ_i2m,...,ρ_iKm}；其中，ρ_iKm表示第i个应力组合下第m个失效模式的第K个测试时间点的模型参数向量；根据所述参数向量，确定产品的寿命特征β_im，其中，β_im表示产品在第i个应力组合下第m个失效模式的寿命特征；

多应力加速模型构建模块，用于根据所述应力组合中包含的冲击应力，构建时间尺度协变量因子，根据所述时间尺度协变量因子和所述寿命特征，构建产品多个应力组合下的多应力加速模型；

可靠度模型构建模块，用于在加速退化试验的随机冲击过程中，根据作用于产品的所述冲击应力的强度函数，确定产品第m个失效模式下的总退化量，以及获取冲击应力的冲击幅值；根据所述冲击幅值达到预先设置的硬失效阈值时确定产品失效的第一失效关系，以及根据所述总退化量达到预先设置的软失效阈值时确定产品失效的第二失效关系；根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型；

输出模块，用于根据所述退化试验数据采用最小二乘算法求出多应力加速模型中的未知参数；将所述未知参数输入到所述可靠度模型中，输出产品在多应力多失效模式相依竞争条件下的可靠度分布。

本发明还提供了一种计算机设备，包括存储器和处理器，所述存储器存储有计算机程序，其特征在于，所述处理器执行所述计算机程序时实现上述方法的步骤。

本发明还提供了一种计算机可读存储介质，其上存储有计算机程序，其特征在于，所述计算机程序被处理器执行时实现上述的方法的步骤。

与现有加速退化试验统计分析方法相比，本方法具有以下显著优点：

(1)本发明在考虑了多应力对寿命特征的耦合影响，退化过程与随机冲击过程的相依性，建立了N应力M失效模式相依竞争条件下的可靠度函数，填补了相关研究领域的空白。

(2)基于产品退化量分布和多应力多失效模式相依竞争可靠度模型，本发明首次系统建立N应力M失效模式相依竞争加速退化试验统计分析方法，因此更加符合工程实际，提高了产品寿命和可靠性评估的准确性。

(3)多应力加速模型参数的多个未知参数不能采用最小二乘算法直接估计，本发明提出二步分析法以解决模型包含多个未知的问题。

附图说明

图1为多应力多失效模式相依竞争下加速退化试验统计分析方法流程图；

图2为随机冲击和退化过程的相依竞争失效图；

图3为多应力组合下产品的磨损量与时间的关系曲线图；

图4为多应力组合下不同时刻磨损量的分布模型参数μ_iK估计值与时间的关系曲线图；

图5为多应力组合下退化量分布参数变化速率β_i估计值与冲击个数的关系曲线图；

图6(a)为磨损量均值为常数时的磨损量分布图；

图6(b)为不同温度应力下的磨损量分布图；

图6(c)为随机冲击载荷作用下的磨损量分布图；

图6(d)为不同湿度应力下的磨损量分布图；

图6(e)为不同温湿度应力下的磨损量分布图；

图6(f)为不同温湿度应力、冲击载荷下的磨损量分布图；

图7为不同应力水平条件下可靠度曲线图；

图8(a)为不同应力水平下的概率密度函数曲线图；

图8(b)为不同应力水平下的失效率函数曲线图。

具体实施方式

为了使本申请的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本申请进行进一步详细说明。应当理解，此处描述的具体实施例仅仅用以解释本申请，并不用于限定本申请。

多应力多失效模式相依竞争条件下加速退化试验统计分析方法，包括以下步骤：

步骤1：获取产品多应力加速退化试验数据，所述多应力加速退化试验数据为多个产品在多个应力组合多个退化失效模式测试时间下进行加速退化实验得到的。

具体的，设计多应力加速退化试验方案，选取适当的加速应力水平及样本量，开展多应力加速退化试验。选择q个加速应力组合S₁,S₂,...,S_q，在每种应力水平组合下随机抽取一定数量的样品进行一定时间的恒定多应力加速退化试验，如在多应力水平组合S_i＝(S₁,S₂,...,S_N)下(i＝1,2,...,q)，随机选取P个样品投入试验。在试验过程中，按照一定间隔τ_i1,τ_i2,...τ_iK对产品退化参数进行K次测试。记S_i应力组合下τ_ik时刻第p个试验样本，第m个退化失效模式的退化数据为(τ_ik,y_ipkm)(p＝1,2,...,P；k＝1,2,K；m＝1,2,...,M-1)。

步骤2：根据所述退化试验数据，确定所述退化试验数据对应的分布模型，以及获取所述分布模型的参数向量。

具体的，在各个测试时间点{τ_i1,τ_i2,...τ_iK}，基于退化量数据选择合适的分布模型F(y；ρ_iKm)，并估计分布模型参数向量{ρ_i1m,ρ_i2m,...,ρ_iKm}。ρ_iKm表示第i个应力组合下第m个失效模式的第K个测试时间点的模型参数向量，如当退化数据服从正态分布时ρ_iKm＝(μ_iKm,σ_iKm)，其中μ_iKm表示均值，σ_iKm表示标准差。

步骤3：根据所述参数向量，确定产品的寿命特征。

具体的，选择合适的模型描述产品各退化量分布参数{ρ_i1m,ρ_i2m,...,ρ_iKm}与测试时间之间{τ_i1,τ_i2,...τ_iK}的关系并估计模型参数β_im，β_im表示第i个应力组合下第m个失效模式的寿命特征，与多应力水平组合S_i相关，可用多应力加速模型表示，如当退化数据服从正态分布，且分布参数与测试时间满足线性关系时，μ_iKm＝β_im·τ_iK。

步骤4：根据所述应力组合中包含的冲击应力，构建时间尺度协变量因子，根据所述时间尺度协变量因子和所述寿命特征，构建产品多个应力组合下的多应力加速模型。

具体的，产品在第i个应力水平组合下第m个退化失效模式的寿命特征β_im，可用多应力加速模型描述，考虑应力耦合的多应力加速模型可表示为

其中，α₀,α₁,...α_u为多应力加速模型参数；

ξ_i＝(ξ₁,ξ₂,...ξ_n)为多应力水平组合S_i的标准化环境应力向量，ξ_i为某一标准化环境应力，如标准化温度应力、湿度应力、冲击应力等。

然而，由于冲击载荷的时间和幅值均具有随机性，多应力加速模型不能表达随机冲击载荷对寿命特征β_im的影响。学者提出时间尺度协变量因子，描述了随机冲击对寿命特征的影响，即

其中，γ₁,γ₂≥0表示多应力加速模型的未知参数，Y_j为第j次随机冲击对产品退化量的影响值，L(t)为时间t内产品承受的非致命性损伤个数，冲击幅值表示为W_j(1≤j≤L)～N(μ_W,σ² _W)。随机冲击对产品的损伤表示为Y_j(1≤j≤L)～N(μ_Y,σ² _Y)，W_j与Y_j是独立的。上式中第一项γ₁L(t)表示冲击载荷的个数对寿命特征的影响，第二项

表示冲击载荷的累积对寿命特征的影响。

结合多应力加速模型和时间尺度协变量因子，可得在N应力综合作用下(包含冲击应力)，产品多个应力组合下的多应力加速模型，表示为

其中，α₀,α₁,...α_j,...α_u,α_u+1,α_u+2为多应力加速模型的参数；

ξ_i＝(ξ₁,ξ₂,...ξ_N-1)为多应力水平组合S_i的标准化环境应力向量，ξ_j为某一标准化环境应力(1≤j≤N-1)，如标准化温度应力、湿度应力、冲击应力等；Y_j表示第j次冲击产生的退化量；L(t)表示t时间之内的冲击个数。

另外，α_u+1＝γ₁，α_u+2＝γ₂。

步骤5：在加速退化试验的随机冲击过程中，根据作用于产品的所述冲击应力的强度函数，确定产品第m个失效模式下的总退化量，以及获取冲击应力的冲击幅值；根据所述冲击幅值达到预先设置的硬失效阈值时确定产品失效的第一失效关系，以及根据所述总退化量达到预先设置的软失效阈值时确定产品失效的第二失效关系；根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型。

具体的，第一失效关系和第二失效关系的相依竞争，表现为：a、产品在多应力作用下发生性能退化，会增加随机冲击的强度；b、在随机冲击载荷作用下，既会增大产品的退化特征量，也会同时增大产品的退化速率，如图2所示。第一失效关系可以为突发失效，如轴和轴套在突然的随机冲击载荷下导致轴或轴套的断裂。第二失效关系可以为退化失效，如轴和轴套之间的磨损达到失效阈值。

假设产品M-1个退化失效的理论退化模型为一般混合效应模型，表示为

其中，X＝(X₁,X₂,…X_M-1)，

表示产品初始退化量，为固定效应参数向量，θ＝(θ₁,θ₂,…θ_M-1)表示产品的退化速率，为随机效应参数向量。该模型在t时刻，X与θ同分布，ε表示随机测量误差，独立同分布，服从正态分布ε～N(0,σ_ε ²)。

在随机冲击载荷作用下，产品的性能退化被分为多个阶段，较为广泛使用的维纳过程模型不再适用。此外，由于Poisson过程中的每个点随机独立于其他点，而在相依竞争条件下，随机冲击过程会受到退化过程的影响，使得每次的随机冲击不再独立，故泊松过程不适宜描述相依条件下的随机冲击过程，Facilitation模型作为一个特殊的马尔可夫过程，可用于描述随机冲击过程，时间t内已承受j个随机冲击的强度函数λ_i(t)为

λ_j(t)＝(1+ηj)λ₀(t)

其中，η为强度因子，λ₀(t)是时间t内已承受0个随机冲击的强度函数，表达式为

λ₀(t)＝λ₀+χX_Sm(t)

其中，λ₀为初始冲击的初始强度，χ为相依因子，表示X_Sm(t)对强度函数的影响。X_Sm(t)表示产品第m个失效模式下的总退化量(m＝1,2,…,M-1)，由理论退化量和随机冲击导致的突然退化量组成，可表示为：

其中，X_m表示产品第m个失效模式下的理论退化量；

表示产品第m个失效模式下的初始退化量；θ_m表示产品第m个失效模式下的退化速率，L(t)表示t时间之内的冲击个数。

P_j表示L(0)＝0状态至L(t)＝j状态的转换概率，即时间t之前产品已经承受的冲击数为j的概率，P_j(t)服从负二项分布，

其中，Λ₀表达式为

产品在第j次随机冲击前未发生突发失效的概率为：

P(W_j＜D)＝F_W(D)

其中，W_j为第j次随机冲击的幅值(j＝1,2,…L)，D为硬失效阈值，F_W(D)为W_j的累积分布函数。

产品的失效是M-1个退化过程达到阈值H_m(m＝1,2,…,M-1)，或者随机冲击载荷达到阈值D之间的相依竞争失效，可靠度函数可表示为

R(t)＝P(W₁＜D,W₂＜D,...,W_n＜D,X_S1(t)＜H₁,...,X_Sm(t)＜H_m)

由于M-1个退化过程相互独立，故式可以表示为

R(t)＝P(W₁＜D,W₂＜D,...,W_n＜D,X_S1(t)＜H₁)

·P(W₁＜D,W₂＜D,...,W_n＜D,X_S2(t)＜H₂)

,...,·P(W₁＜D,W₂＜D,...,W_n＜D,X_Sm(t)＜H_m)

根据全概率公式可得，相依竞争条件下的可靠度函数可表示为

将第m个失效模式下的总退化量X_Sm(t)代入上式，可靠度函数可表示为

为简化表达，令事件

依据条件独立的定义，当退化模型的随机效应参数向量θ、冲击个数L(t)及冲击路径Y_j确定的条件下，令{θ_m＝x_m}＝B_m,(m＝1,2,…,M-1)，

可得可靠度函数

进一步可化简为

即

其中，P[E₁|A,B_m,C](m＝1,2,…,M-1)表示在冲击个数L、随机效应参数向量θ及累积损伤u确定下，产品并未发生硬失效的概率，只与冲击个数L有关，可表示为

其中，Ф是标准正态分布的概率分布函数。

P[E₂|A,B₁,C]表示在L个冲击、参数θ₁及累积损伤u确定下，产品在第一个退化模式下并未发生软失效的概率可表示为

P[E_m+1|A,B_m,C]表示在L个冲击、参数θ_m及累积损伤u确定下，产品在第m个退化模式下并未发生软失效的概率可表示为

P[A|B_m,C](m＝1,2,…,M-1)表示参数θ_m及累积损伤u确定下，在t时间之前产生L个冲击的概率，可表示为

是n个Y_j变量之和

的概率密度函数(j＝0,1,2，n)，假设冲击载荷幅值服从正态分布

因此

则有

将P[E₁|A,B_m,C]、P[E₂|A,B₁,C]、P[E_m+1|A,B_m,C]及P[A|B_m,C]表达式带入可靠度函数，可得M-1个退化失效与突发失效的相依竞争可靠度函数模型的解析表达式表示为

其中，D为硬失效阈值，F_W(D)为W_j的累积分布函数，W_j为第j次随机冲击的幅值，j＝1,2,…L；θ_m表示第m个失效模式下的退化速率；λ_j(t)表示时间t内已承受j个随机冲击的强度函数；η为强度因子；λ₀(t)是时间t内已承受0个随机冲击的强度函数；λ₀为初始冲击的初始强度；χ为相依因子；f_θm表示第m个退化模式下的概率密度函数；μ_Y,表示随机冲击载荷所产生退化损伤的均值，σ_Y表示随机冲击载荷所产生退化损伤的标准差。

步骤6：根据所述退化试验数据采用最小二乘算法求出多应力加速模型中的未知参数；将所述未知参数输入到所述可靠度模型中，输出产品在多应力多失效模式相依竞争条件下的可靠度分布。

具体的，基于退化数据采用最小二乘算法(LSE)求出寿命特征β_im中的未知参数，带入解析表达式即可求到产品在正常应力组合下的可靠度函数。之所以选择LSE估计未知参数而非极大似然估计法(MLE)的原因是：LSE可以用于加速寿命试验，也可以用于基于伪失效和基于退化量分布的加速退化试验；MLE只能用于加速寿命试验和基于伪失效寿命的加速退化试验，不可以用于基于退化量分布的加速退化试验，其原因为在基于退化量分布的加速退化试验分析方法中，退化量分布参数{ρ_i1m,ρ_i2m,...,ρ_iKm}通常是时间的一次函数，斜率与多应力构成加速模型，在多个斜率与加速应力已知的条件下，可以通过最小二乘算法计算出未知参数，而MLE的基础是计算出某时刻的可靠度函数和概率密度函数，与不同时刻所求出的分布参数无关。

需要说明的是，传统LSE是拟合分布参数与时间的一次函数，斜率即为多应力的函数，采用加速模型来描述，此种方法假设了应力与时间的无关性。然而，实际服役环境下，应力通常也是时间的函数，比如冲击载荷，是时间的随机离散函数。本发明所假设的模型参数与多应力的关系，考虑了应力与时间对分布参数的耦合影响。即冲击载荷随时间发生变化，与其他应力的耦合影响着模型参数的变化。因此，加速模型参数的未知参数不能采用最小二乘算法直接估计，本发明提出二步分析法以解决模型包含多个未知的问题，具体为：

假设产品在试验过程中所承受冲击载荷的个数为L，时间分别为t₁,t₂,…t_L。

(a)基于0冲击载荷过程退化数据(τ_ik,y_ijkm)，其中τ_ik<t₁。此阶段不包含冲击过程，因此，可以基于多个应力水平下的分布参数变化速率与多应力水平，采用最小二乘算法求取未知参数α₀,α₁,α₂,…,α_u。

(b)将步骤(a)所估计的α₀,α₁,α₂,…,α_u带入多应力加速模型中，基于多个冲击载荷过程退化数据(τ_ik,y_ijkm)，其中t₁<τ_ik<t_s，采用最小二乘算法求取未知参数α_u+1,α_u+2。

本发明在考虑了多应力对寿命特征的耦合影响，退化过程与随机冲击过程的相依性，建立了N应力M失效模式相依竞争条件下的可靠度函数，填补了相关研究领域的空白。基于产品退化量分布和多应力多失效模式相依竞争可靠度模型，本发明首次系统建立N应力M失效模式相依竞争加速退化试验统计分析方法，因此更加符合工程实际，提高了产品寿命和可靠性评估的准确性。

另外，多应力加速模型参数的多个未知参数不能采用最小二乘算法直接估计，本发明提出二步分析法以解决模型包含多个未知的问题。

本发明研究产品在多应力综合作用下发生多失效模式相依竞争的可靠度建模，包括多应力对退化速率的影响、退化过程与随机冲击过程的相互依赖；然后研究多应力综合作用下的加速退化试验统计分析方法，主要解决可靠度模型多个未知参数的估计问题。

本发明还提供了多应力多失效模式相依竞争条件下加速退化试验分析装置，所述装置包括：

数据获取模块，用于获取产品多应力加速退化试验数据；所述多应力加速退化试验数据为多个产品在多个应力组合多个退化失效模式测试时间下进行加速退化实验得到的；

数据处理模块，用于根据所述退化试验数据，确定所述退化试验数据对应的分布模型，以及获取所述分布模型的参数向量{ρ_i1m,ρ_i2m,...,ρ_iKm}；其中，ρ_iKm表示第i个应力组合下第m个失效模式的第K个测试时间点的模型参数向量；根据所述参数向量，确定产品的寿命特征β_im，其中，β_im表示产品在第i个应力组合下第m个失效模式的寿命特征；

多应力加速模型构建模块，用于根据所述应力组合中包含的冲击应力，构建时间尺度协变量因子，根据所述时间尺度协变量因子和所述寿命特征，构建产品多个应力组合下的多应力加速模型；

可靠度模型构建模块，用于在加速退化试验的随机冲击过程中，根据作用于产品的所述冲击应力的强度函数，确定产品第m个失效模式下的总退化量，以及获取冲击应力的冲击幅值；根据所述冲击幅值达到预先设置的硬失效阈值时确定产品失效的第一失效关系，以及根据所述总退化量达到预先设置的软失效阈值时确定产品失效的第二失效关系；根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型；

输出模块，用于根据所述退化试验数据采用最小二乘算法求出多应力加速模型中的未知参数；将所述未知参数输入到所述可靠度模型中，输出产品在多应力多失效模式相依竞争条件下的可靠度分布。

本发明还提供了一种计算机设备，包括存储器和处理器，所述存储器存储有计算机程序，其特征在于，所述处理器执行所述计算机程序时实现上述方法的步骤。

本发明还提供了一种计算机可读存储介质，其上存储有计算机程序，其特征在于，所述计算机程序被处理器执行时实现上述的方法的步骤。

下面结合典型竞争失效案例——电液伺服阀两个失效模式之间的相依失效，按照图1所示流程说明本方法的实施方式。

步骤1：获取电液伺服阀多应力加速退化试验数据

温度和湿度应力水平的选择依据均匀正交设计，选取四种应力组合，各个应力组合下具有50个样本量，如表1所示。试验时间共480h，测试时间间隔为24h，试验共测试20次。产品在各个应力组合下的仿真数据如图3所示。

表1试验应力水平设置

步骤2：选择合适的分布模型并估计分布参数

采用KS检验不同应力组合下各个时刻的退化数据(磨损量)是否服从正态分布模型F(y；ρ)＝F(y_i；μ_i,σ_i)，(i＝1,2,3,4)表示4个应力组合，检验结果如表2所示。估计不同时刻{τ₁,τ₂,...,τ_K}分布模型参数{ρ_i1,ρ_i2,...,ρ_iK}＝{μ_i1,μ_i2,...,μ_iK}，K＝20，如图4所示。

表2K-S检验结果(5％的置信水平)

步骤3：建立退化量分布参数与时间之间的模型并估计模型参数

由图4可知，分布参数与时间呈多阶段线性关系，因此拟采用多阶段线性模型来描述产品各退化量分布参数{ρ_i1,ρ_i2,...,ρ_iK}＝{μ_i1,μ_i2,...,μ_iK}与测试时间{τ₁,τ₂,...,τ_K}之间的关系，模型显著性检验如表3所示，故可用μ_iK＝β_i·τ_K进行拟合，其中t₁,t₂,...,t₇为冲击载荷的随机到达时间。不同应力组合下退化量分布参数变化速率β_i的估计值，如图5所示。

表3模型显著性检验结果(5％的置信水平)

步骤4：构造多应力加速模型

依据式，考虑三应力综合作用下(含随机冲击应力)的加速模型表达式为

其中，α₁₀,α₁₁,α₁₂,α₁₃表示环境应力对退化速率的影响，α₁₄＝γ₁，α₁₅＝γ₂。ξ₁、ξ₂分别为标准化后的温度、湿度应力。

温度应力选取阿伦尼斯模型，湿度应力选取工业界广泛使用的幂率型湿度模型，故上式中ξ₁、ξ₂的表达式分别为

ξ₁＝ξ(T)＝[log(1/T)-log(1/T₀)]/[log(1/T_H)-log(1/T₀)]

ξ₂＝ξ(RH)＝(logRH-logRH₀)/(logRH_H-logRH₀)

其中，X₀为正常使用下应力水平，X_H为加速应力水平。

步骤5：多失效模式相依竞争可靠度建模

电液伺服阀有两个失效模式，其一是内部轴和套筒之间相互磨损而发生软失效；其二是在随机冲击载荷作用下油液中突然产生较大的磨粒导致污染锁死而发生硬失效。内部轴和套筒之间相互磨损可以通过线性模型描述，污染锁死可以通过Facilitation模型来描述。两个失效模式的相依性表现为：内部轴和套筒之间相互磨损可能产生较大磨粒会导致污染锁死的发生；污染锁死会产生额外的磨粒，增强磨损过程。

由一般混合效应模型可知，轴和套筒的磨损过程可表示为

其中，

表示初始磨损幅值，θ_f表示电液伺服阀的磨损率。

退化失效与突发失效的相依竞争可靠度函数的解析表达式可表示为

步骤6：求解模型未知参数

依据退化数据及多应力加速模型估计模型未知参数。基于0冲击载荷过程退化数据(τ_ik,y_ijkm)，其中τ_ik<t₁，采用最小二乘算法估计出前4个未知参数α₁₀,α₁₁,α₁₂,α₁₃。再将得到的参数带入多应力加速模型，基于多个冲击载荷过程退化数据(τ_ik,y_ijkm)，其中t₁<τ_ik<t_s，采用最小二乘算法求取未知参数α₁₄,α₁₅。

采用均方误差MSE来描述参数的估计精度，均方误差

是真值θ的估计，是指参数估计值与参数真值之差平方的期望值，是衡量“平均误差”的一种较方便的方法，可以评价数据的变化程度，MSE的值越小，说明预测结果具有更好的精确度。通过Monte Carlo仿真方法进行100次仿真，计算待估参数的MSE，如表4所示。

表4模型参数估计值及相对误差

将估计参数的均值带入可靠度函数的解析表达式即可得产品可靠度函数

其中，

此外，基于{μ_i1,μ_i2,...,μ_iK}与测试时间{τ₁,τ₂,...,τ_K}之间的线性关系，可以仿真得到初始磨损量在不同应力状态下的分布。如图6(a)所示为磨损量不随应力发生变化的分布。图6(b)、图6(c)、图6(d)为不同应力下的磨损量分布，可以看出，温度应力、湿度应力及冲击载荷越大，磨损量的分布参数μ_f越大。由图6(e)和图6(f)可以看出，在双应力或多应力组合下，磨损量分布参数相对于单应力有更快的变大趋势。

基于多应力多失效模式相依竞争条件下的可靠度模型，结合不同应力水平条件下的磨损量分布参数μ_f，可以仿真得到产品在不同应力水平下的可靠度曲线，如图7所示。可以看出，温度应力、湿度应力越大，同一时刻的可靠度越低，这是由于应力的增大使得磨损量分布参数μ_f增大，进而使得可靠度降低。概率密度函数能够反映可靠度随时间的变化速率，根据上节所提出的数值分析方法可以得到不同应力水平下的概率密度函数，如图8(a)所示，应力水平越大，概率密度函数在较早的时刻达到顶点，表明可靠度变化速率随着应力水平的增大而增大。同时，由λ(t)＝f(t)/R(t)可以得到失效率函数随时间的变化关系，如图8(b)所示，在同一时刻，应力水平越大，失效率越高，进一步说明了本发明所提模型的正确性。

以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.多应力多失效模式相依竞争条件下加速退化试验统计分析方法，其特征在于，所述方法包括：

获取产品多应力加速退化试验数据；所述多应力加速退化试验数据为多个产品在多个应力组合多个退化失效模式测试时间下进行加速退化实验得到的；

根据所述退化试验数据，确定所述退化试验数据对应的分布模型，以及获取所述分布模型的参数向量{ρ_i1m,ρ_i2m,...,ρ_iKm}；其中，ρ_iKm表示第i个应力组合下第m个失效模式的第K个测试时间点的模型参数向量；

根据所述参数向量，确定产品的寿命特征β_im，其中，β_im表示产品在第i个应力组合下第m个失效模式的寿命特征；

根据所述应力组合中包含的冲击应力，构建时间尺度协变量因子，根据所述时间尺度协变量因子和所述寿命特征，构建产品多个应力组合下的多应力加速模型；

在加速退化试验的随机冲击过程中，根据作用于产品的所述冲击应力的强度函数，确定产品第m个失效模式下的总退化量，以及获取冲击应力的冲击幅值；

根据所述冲击幅值达到预先设置的硬失效阈值时确定产品失效的第一失效关系，以及根据所述总退化量达到预先设置的软失效阈值时确定产品失效的第二失效关系；

根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型；

根据所述退化试验数据采用最小二乘算法求出多应力加速模型中的未知参数；

将所述未知参数输入到所述可靠度模型中，输出产品在多应力多失效模式相依竞争条件下的可靠度分布。

2.根据权利要求1所述的方法，其特征在于，根据所述应力组合中包含的冲击应力，构建时间尺度协变量因子，根据所述时间尺度协变量因子和所述寿命特征，构建产品在多个应力组合下的多应力加速模型，包括：

产品在多个应力组合下的多应力加速模型表示为：

其中，α₀,α₁,...α_j,...α_u,α_u+1,α_u+2为多应力加速模型的参数；

ξ_i＝(ξ₁,ξ₂,...ξ_N-1)为多应力水平组合S_i的标准化环境应力向量，ξ_j为某一标准化环境应力，其中1≤j≤N-1；Y_j表示第j次冲击产生的退化量；L(t)表示t时间之内的冲击个数。

3.根据权利要求1所述的方法，其特征在于，在加速退化试验的随机冲击过程中，根据作用于产品的所述冲击应力的强度函数，确定产品第m个失效模式下的总退化量，以及获取冲击应力的冲击幅值，包括：

产品第m个失效模式下的总退化量X_Sm(t)表示为

其中，X_m表示产品第m个失效模式下的的理论退化量；

表示产品第m个失效模式下的初始退化量；θ_m表示产品第m个失效模式下的退化速率，L(t)表示t时间之内的冲击个数。

4.根据权利要求1所述的方法，其特征在于，根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型中，

所述第一失效关系和第二失效关系的相依竞争表现为：产品在多应力作用下发生性能退化，会增加随机冲击的强度；而在随机冲击载荷作用下，既会增大产品的退化特征量，也会同时增大产品的退化速率。

5.根据权利要求4所述的方法，其特征在于，根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型，包括：

所述可靠度模型表示为：

其中，D为硬失效阈值，F_W(D)为W_j的累积分布函数，W_j为第j次随机冲击的幅值，j＝1,2,…L；θ_m表示第m个失效模式下的退化速率；λ_j(t)表示时间t内已承受j个随机冲击的强度函数；η为强度因子；λ₀(t)是时间t内已承受0个随机冲击的强度函数；λ₀为初始冲击的初始强度；χ为相依因子；f_θm表示第m个退化模式下的概率密度函数；μ_Y,表示随机冲击载荷所产生退化损伤的均值，σ_Y表示随机冲击载荷所产生退化损伤的标准差。

6.根据权利要求1所述的方法，其特征在于，根据所述退化试验数据采用最小二乘算法求出多应力加速模型中的未知参数，包括：

设定产品在试验过程中所承受冲击载荷的个数为L，时间分别为t₁,t₂,…t_L；

(a)、基于0冲击载荷过程退化数据(τ_ik,y_ijkm)，其中τ_ik<t₁；此阶段不包含冲击过程，因此，基于多个应力水平下的分布参数变化速率与多应力水平，采用最小二乘算法求取未知参数α₀,α₁,α₂,…,α_u；

(b)、将步骤(a)所求取的α₀,α₁,α₂,…,α_u带入所述多应力加速模型，基于多个冲击载荷过程退化数据(τ_ik,y_ijkm)，其中t₁<τ_ik<t_s，采用最小二乘算法求取未知参数α_u+1,α_u+2。

7.根据权利要求1所述的方法，其特征在于，所述退化试验数据服从正态分布。

8.多应力多失效模式相依竞争条件下加速退化试验分析装置，其特征在于，所述装置包括：

数据获取模块，用于获取产品多应力加速退化试验数据；所述多应力加速退化试验数据为多个产品在多个应力组合多个退化失效模式测试时间下进行加速退化实验得到的；

数据处理模块，用于根据所述退化试验数据，确定所述退化试验数据对应的分布模型，以及获取所述分布模型的参数向量{ρ_i1m,ρ_i2m,...,ρ_iKm}；其中，ρ_iKm表示第i个应力组合下第m个失效模式的第K个测试时间点的模型参数向量；根据所述参数向量，确定产品的寿命特征β_im，其中，β_im表示产品在第i个应力组合下第m个失效模式的寿命特征；

多应力加速模型构建模块，用于根据所述应力组合中包含的冲击应力，构建时间尺度协变量因子，根据所述时间尺度协变量因子和所述寿命特征，构建产品多个应力组合下的多应力加速模型；

可靠度模型构建模块，用于在加速退化试验的随机冲击过程中，根据作用于产品的所述冲击应力的强度函数，确定产品第m个失效模式下的总退化量，以及获取冲击应力的冲击幅值；根据所述冲击幅值达到预先设置的硬失效阈值时确定产品失效的第一失效关系，以及根据所述总退化量达到预先设置的软失效阈值时确定产品失效的第二失效关系；根据第一失效关系和第二失效关系的相依竞争，构建可靠度模型；

输出模块，用于根据所述退化试验数据采用最小二乘算法求出多应力加速模型中的未知参数；将所述未知参数输入到所述可靠度模型中，输出产品在多应力多失效模式相依竞争条件下的可靠度分布。

9.一种计算机设备，包括存储器和处理器，所述存储器存储有计算机程序，其特征在于，所述处理器执行所述计算机程序时实现权利要求1至7中任一项所述方法的步骤。

10.一种计算机可读存储介质，其上存储有计算机程序，其特征在于，所述计算机程序被处理器执行时实现权利要求1至7中任一项所述的方法的步骤。

Images