一种考虑随机载荷和强度退化的疲劳可靠性分析方法与流程

本发明属于可靠性分析技术领域，尤其针对机电产品疲劳可靠性分析方法。

背景技术：

在载荷的反复作用下，各种结构及其加工而成的产品都会产生疲劳。据统计，疲劳破坏在机械结构破坏中所占比例高达50％-90％，疲劳破坏往往十分突然，并且造成灾难性的事故，引起巨大的经济损失。特别是近40年来，在机械向高温、高速以及大型化方向发展，机械零件所承受的应力越来越高，服役工况越来越复杂，导致疲劳失效事故越来越多。因此，疲劳可靠性分析方法的研究，对预防疲劳破坏事故的发生和经济发展意义重大。

疲劳可靠性的研究目的主要在于从经济性和维修性出发，在规定的工况下，规定的使用时间内，以及完成规定功能的条件下，使得产品因疲劳而失效的可能性降到最低程度。目前结构元件的疲劳分析方法主要有三类：一是疲劳累积损伤模型，该模型基于传统的机构疲劳寿命模型，结合疲劳累积损伤理论和应力统计分布，得到疲劳寿命的相关统计参数，同时该模型比较简单，只能反映疲劳破坏的基本趋势；二是疲劳寿命模型，该模型以疲劳寿命作为参数，结合疲劳累积损伤进行疲劳可靠性分析，但是在随机载荷作用下疲劳寿命的瞬时分布的解析表达式很难给出；三是剩余强度模型，该模型通过对元件剩余强度在疲劳载荷下随时间变化规律的研究，从而给出其寿命参数。目前的研究工作主要关注在以剩余强度为指标的疲劳可靠性分析方法，而直接以安全余量为指标的疲劳可靠性研究尚未开展。

技术实现要素：

本发明为解决上述技术问题，提出了一种考虑随机载荷和强度退化的疲劳可靠性分析方法，以安全余量为指标，通过随机载荷仿真，构建安全余量向量和对应的载荷作用次数向量，从而建立疲劳可靠性分析模型。

本发明采用的技术方案是：一种考虑随机载荷和强度退化的疲劳可靠性分析方法，包括：

S1、根据建立的产品随机载荷数学模型和剩余强度数学模型，构建产品安全余量数学模型；

S2、给定随机载荷作用次数，通过仿真得到安全余量向量和相对应的载荷作用次数向量；

S3、结合步骤S2所得的安全余量向量和相对应的载荷作用次数向量，利用核函数估计建立安全余量向量和相对应的载荷作用次数向量的联合概率密度函数，最终建立考虑随机载荷和强度退化的产品疲劳可靠性分析模型。

进一步地，所述步骤S1包括：

S11、通过分析产品在实际工况中的随机载荷作用情况，建立随机载荷数学模型；

S12、根据步骤S11得到的随机载荷数学模型，建立剩余强度模型；

S13、根据步骤S11得到的随机载荷数学模型，以及步骤S12得到的剩余强度模型，建立产品安全余量模型。

进一步地，所述步骤S2包括：

S21、给定随机载荷作用次数，通过仿真得到安全余量最小值和相对应的载荷作用次数；

S22、若步骤S21所得的安全余量最小值小于零，则用安全余量小于零的第一时间点对应的安全余量值和载荷作用次数作为最小值点对应的安全余量最小值和载荷作用次数；否则保持原来的安全余量最小值和对应的载荷作用次数不变；

S23、重复N次步骤S21、步骤S22的操作，得到N维的安全余量向量，和N维的载荷作用次数向量。

更进一步地，所述步骤S21具体包括以下分步骤：

S211、设定随机载荷作用次数为L₁，则仿真一次得到包含L₁个元素的随机载荷向量；

S212、根据步骤S211得到的随机载荷向量得到随机载荷谱，并根据该随机载荷谱计算得到剩余强度向量；

S213、结合步骤S211得到的随机载荷向量和S212得到的剩余强度向量，根据安全余量数学模型的计算公式，得到安全余量向量；所述的安全余量数学模型的计算公式为：

M(ω)＝R(ω)-s(ω)；

其中，M(ω)表示安全余量，R(ω)表示剩余强度，s(ω)表示随机载荷，ω表示载荷作用次数；

S214、对步骤S213所得的安全余量向量进行比较，选择安全余量的最小值和安全余量最小值对应的载荷作用次数。

更进一步地，步骤S3所述核函数K(v)的表达式为：

其中：v表示核函数参数。

进一步地，所诉步骤S3还包括：

A1、计算安全余量向量的均方差，得到核密度估计中的安全余量向量最佳带宽；

A2、计算载荷作用次数向量的均方差，得到核密度估计中的载荷作用次数向量最佳带宽；

A3、结合步骤A1得到的安全余量向量最佳带宽，步骤A2得到的载荷作用次数向量最佳带宽，以及核函数公式得到安全余量向量和载荷作用次数的联合概率密度函数；

A4、根据步骤A3得到的联合概率密度函数，建立考虑随机载荷和强度退化的产品疲劳可靠性分析模型。

更进一步地，所述步骤A4具体为：

根据安全余量向量和载荷作用次数向量的联合概率密度函数，得到疲劳失效概率；

根据疲劳失效概率，得到考虑随机载荷和强度退化的产品疲劳可靠性分析模型。

更进一步地，所述疲劳失效概率表达式为：

其中，w表示随机载荷作用次数，x为取最小值点时对应的载荷作用次数；y为最小值点对应的安全余量数值。

本发明的有益效果：本发明通过建立产品的安全余量方程，并对一定次数随机载荷进行仿真，计算安全余量最小值和对应的载荷作用次数，利用核函数估计分别建立安全余量向量的概率密度函数、安全余量向量和对应载荷作用次数向量的联合概率密度函数，最终建立考虑随机载荷和强度退化的产品疲劳可靠性分析模型，相比于现有的剩余强度模型，本申请不需要考虑剩余强度和随机载荷概率密度函数类型以及二者之间的关系，直接以安全余量方程为指标进行求解，本申请对全面提高产品的设计水平和制定合适的维修策略具有较高的理论支撑和工程实践意义。

附图说明

图1是本发明疲劳可靠性分析方法的流程图。

图2是本发明实施例中滤波减速器固定齿轮的三维模型图。

图3是本发明实施例中滤波减速器固定齿轮安全余量向量的直方图。

图4是本发明实施例中滤波减速器固定齿轮安全余量向量对应的载荷作用次数向量的直方图。

图5是本发明实施例中可靠度随载荷作用次数的变化图。

具体实施方式

为便于本领域技术人员理解本发明的技术内容，下面结合附图对本发明内容进一步阐释。

以滤波减速器固定齿轮为例进行说明，在以本发明技术方案为前提下进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。

如图1所示为本发明的方案流程图，本发明采用的技术方案是：一种考虑随机载荷和强度退化的疲劳可靠性分析方法，针对如图2所示滤波减速器固定齿轮，包括以下步骤：

S1、建立产品随机载荷和剩余强度数学模型，结合随机载荷模型和剩余强度模型得到产品的安全余量数学模型；具体包括：

S11、分析根据谐波减速器固定齿轮在实际工况中和双联齿轮的啮合情况，固定齿轮所受的随机载荷s(ω)的数学模型如公式(1)所示：

公式(1)中，Z_H为区域系数，本实施例中取值为2.5；Z_E为材料的弹性影响系数，本实施例中取值为189.8；Z_ε为重合度系数，本实施例中取值为0.9642；K为载荷系数，本实施例中取值为0.3635；T(ω)为输出转矩，为随机过程，在给定ω时，假定T(ω)服从正态分布N(50933，4244.1²)；b为固定齿轮齿宽，服从正态分布N(5，0.033²)；d为固定齿轮分度圆直径，服从正态分布N(79.5，0.06²)；μ是滤波减速器固定齿轮与双联齿轮的齿数比，本实施例中取值为0.9056。

S12、结合步骤S11得到的随机载荷作用数学模型，得到固定齿轮的剩余强度R(ω)模型如公式(2)所示：

公式(2)中，δ_s表示固定齿轮的初始强度，本实施例中取值为1092MPa；s(ω)为步骤S11求得的随机载荷；m为S-N曲线(S-N曲线，也称应力－寿命曲线，具体为：以材料标准试件疲劳强度为纵坐标，以疲劳寿命的对数值lg N为横坐标，表示一定循环特征下标准试件的疲劳强度与疲劳寿命之间关系的曲线。)的斜率参数，本实施例中取值为11.4223；C为S-N曲线的材料参数，本实施例中取值为10^39.7032；ω为载荷作用次数，取值范围为[0，10⁶]；c为材料相关的常数，本实施例中取值为2。

S13、结合步骤S11得到的随机载荷数学模型和步骤S12得到的剩余强度模型，固定齿轮的安全余量模型如公式(3)所示：

M(ω)＝R(ω)-s(ω) 公式(3)

S2、给定随机载荷作用次数，通过仿真得到安全余量最小值和相对应的载荷作用次数，若安全余量最小值小于零，则用安全余量小于零的第一时间点对应的安全余量值和载荷作用次数替代最小值点对应的安全余量最小值和载荷作用次数，重复上述仿真N次，得到最小安全余量向量和相对应的载荷作用次数向量；包括：

S21、本实例中，给定随机载荷作用的次数ω为10⁶，通过仿真得到安全余量最小值和相对应的载荷作用次数；具体包括如下分步骤：

S211、当载荷作用次数为10⁶时，仿真一次得到随机载荷向量，包含10⁶个元素；

S212、根据步骤S211得到的随机载荷谱，计算得到剩余强度向量；

S213、结合步骤S211得到的随机载荷向量和S212得到的剩余强度向量，根据安全余量模型计算公式(3)，得到安全余量向量；

S214、对步骤S213得到的安全余量向量进行比较，选择安全余量的最小值和对应的载荷作用次数，本实施例中通过一次仿真后可以得到最小的安全余量为42.75MPa，对应的载荷作用次数为984367；

S22、若步骤S21得到的安全余量最小值小于零，则用安全余量小于零的第一时间点对应的安全余量值和载荷作用次数替代最小值点对应的安全余量最小值和载荷作用次数；否则保持原来的安全余量最小值和对应的载荷作用次数不变；具体包括如下分步骤：

S221、对步骤S21得到的安全余量最小值进行判断；

S222、对于步骤S221的判断结果，如果安全余量最小值大于零，则保持原来的安全余量最小值和对应的载荷作用次数不变；

S223、对于步骤S221的判断结果，如果安全余量小于零，则用安全余量小于零的第一时间点对应的安全余量值和载荷作用次数替代最小值点对应的安全余量最小值和载荷作用次数。

S23、对步骤S21和S22仿真N次，N可根据工程实际中的效率和精度要求进行取值，本实施例中N＝1000，得到最小安全余量向量和相对应的载荷作用次数向量；具体包括如下分步骤：

S231、通过仿真产生N个随机载荷向量和剩余强度向量，得到N个安全余量向量；

S232、对于每一个安全余量向量，重复步骤S22，得到一个N维的安全余量向量M＝[m₁,m₂,···,m_N]和一个相对应的N维载荷作用次数向量ω＝[ω₁,ω₂,···,ω_N]，安全余量向量的直方图如图3所示，对应的载荷作用次数向量的直方图如图4所示。

S3、结合步骤S2得到的安全余量向量和相对应的载荷作用次数向量，利用核函数估计建立考虑随机载荷和强度退化的产品疲劳可靠性分析模型；具体包括：

S31、对于步骤S2得到的N维的安全余量向量和相对应的N维载荷作用次数向量，用核密度函数估计得到二者的联合概率密度函数；具体包括如下分步骤：

S311、用来估计安全余量向量和载荷作用次数的联合概率密度函数的高斯核函数K(ν)如公式(4)所示：

S312、安全余量向量核密度估计的最佳带宽h_M如公式(5)所示：

公式(5)中，N为安全余量向量的样本个数；σ_M表示安全余量向量的样本均方差，本实施例中取值为32.2699；因此，本实施例中取值为安全余量向量核密度估计的最佳带宽h_M＝8.5859。

S313、载荷作用次数向量核密度估计的最佳带宽h_ω如公式(6)所示：

公式(6)中，N为载荷作用次数向量的样本个数；σ_ω表示安全余量向量的样本均方差，本实施例中取值为63950.4709；因此，本实施例中取值为安全余量向量核密度估计的最佳带宽h_ω＝17015。

S314、结合步骤S311、S312和S313，安全余量向量和载荷作用次数的联合概率密度函数f_Mω(x，y)如公式(7)所示：

公式(7)中，为N安全余量向量的样本个数；h_M为步骤S312所得的安全余量向量最佳带宽，本实施例中取值为8.5859；K(·)为步骤S311中的高斯核函数，即表示步骤S311中的高斯核函数中的时的表达形式，同理表示时的表达形式；m_i(i＝1,2,...，N)为安全余量向量M中的第i个安全余量值；h_ω为步骤S313所得的载荷作用次数向量最佳带宽，本实施例中取值为17015；ω_i(i＝1,2,...，N)为载荷作用次数向量ω中的第i个载荷作用次数。

S32、根据步骤S31所得的概率密度函数，建立考虑随机载荷和强度退化的产品疲劳可靠性分析模型；具体包括如下分步骤：

S321、根据安全余量向量和载荷作用次数向量的联合概率密度函数，疲劳失效概率P(w)如公式(8)所示；

在公式(8)中，ω表示载荷作用次数；Pr{·}表示概率。当ω为一定值时，P(w)为一确定值；而对于一变量t∈[0,w]，则P(t)为一随t变化的曲线，P(t)即为失效概率曲线。

S322、根据步骤S321所得疲劳失效概率曲线，疲劳可靠度曲线如图5所示。

疲劳可靠度曲线R(t)如公式(9)所示：

R(t)＝1-P(t) 公式(9)

本领域的普通技术人员将会意识到，这里所述的实施例是为了帮助读者理解本发明的原理，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的权利要求范围之内。