CN109345032A - 基于动态裂纹数目的粒子滤波多裂纹扩展预测方法 - Google Patents

Abstract

本发明公开了一种基于动态裂纹数目的粒子滤波多裂纹扩展预测方法，属于故障预测与健康管理技术领域。本发明采用在线导波结构健康监测方法在线监测结构中裂纹的萌生和扩展，当监测到新的裂纹萌生，更新多裂纹扩展状态方程中裂纹分量和裂纹扩展方程数目，构建多裂纹耦合下的裂纹扩展规律，并更新粒子集。结合导波结构健康监测方法获得的在线观测值，通过正则化粒子滤波在线估计多裂纹扩展状态向量的后验概率密度函数。通过模型参数的后验估计更新裂纹扩展模型，以裂纹长度的后验估计为起点预测各个裂纹在未来时刻的扩展轨迹。本发明可以有效地实现结构中多裂纹的在线监测和预测，在结构多裂纹扩展预测方面具有广泛的应用前景。

Description

基于动态裂纹数目的粒子滤波多裂纹扩展预测方法

技术领域

本发明涉及一种基于动态裂纹数目的粒子滤波多裂纹扩展预测方法，属于故障预测与健康管理技术领域。

背景技术

故障预测与健康管理(Prognostics and Health Management，PHM)技术通过传感器在线获取工程系统相关的信息，基于该信息对系统的健康状态进行实时诊断和预测。PHM技术对于保障结构安全性和可靠性，制定最佳的运行和维护策略具有重要的理论意义和工程应用价值。

工程结构作为工程系统的承载平台，其损伤诊断和预测方法对于PHM技术至关重要。在交变服役载荷作用下，工程结构的主要损伤形式为疲劳裂纹。对于老龄化结构，单个结构单元中常常会出现多个疲劳裂纹并扩展的情况。多个裂纹的扩展和聚合可能使得结构突然失效，造成灾难性事故。然而，由于疲劳裂纹扩展过程是一个包含各种不确定性的随机过程，比如材料微观结构的不确定性、服役载荷的不确定性、以及环境参数不确定性等。这些不确定性使得疲劳裂纹的起裂和扩展存在很大的分散性。同时多个裂纹之间的相互作用，使得多裂纹起裂和扩展问题比结构中存在单个裂纹的情况更加复杂，包含更多的不确定性。近年来，基于结构健康监测方法和粒子滤波的裂纹扩展预测方法得到越来越多的关注，其通过贝叶斯滤波理论融合结构当前的裂纹长度观测值和先验裂纹扩展模型，以消除上述不确定性因素对于疲劳裂纹扩展预测的影响。

在结构健康监测方法中，基于导波的结构健康监测方法具有小损伤敏感，可实现区域监测等优点，得到广泛关注并开始逐渐应用于工程实际中。因此，结合导波结构健康监测和粒子滤波方法实现多裂纹扩展预测非常具有应用前景。但是对于单一结构中存在多个裂纹的问题，每个裂纹的起裂和扩展并不是在同一个时刻，确定裂纹的萌生，以及表征多裂纹的扩展规律对多裂纹的疲劳裂纹扩展预测至关重要。

发明内容

本发明针对结构中的多裂纹扩展可靠预测问题，提出一种基于动态裂纹数目的粒子滤波多裂纹扩展预测方法，旨在结合导波结构健康监测方法对结构疲劳裂纹的萌生和扩展进行监测，同时动态更新多裂纹扩展状态方程中裂纹分量个数，结合正则化粒子滤波理论实现结构多裂纹损伤状态的在线监测和预测。

本发明为解决其技术问题采用如下技术方案：

一种基于动态裂纹数目的粒子滤波多裂纹扩展预测方法，包括如下步骤：

(1)离线状态下，基于实验数据建立结构多裂纹扩展状态和结构损伤程度指数之间的映射关系；从反映结构损伤程度变化的导波信号中提取多种损伤因子，组成多维导波损伤因子向量；建立损伤因子向量和结构损伤程度指数的多维非线性映射关系；

(2)在结构服役的初始时刻k＝0，对算法进行初始化；设定初始裂纹的个数n＝1；根据结构的几何形式和服役条件确定结构中单个裂纹的裂纹扩展状态方程；在结构上布置压电传感网络，通过导波结构健康监测方法在线监测结构疲劳裂纹的萌生和扩展；

(3)在结构服役过程中的当前时刻k，计算每个粒子中各个裂纹分量的形状修正系数；结合定义的多裂纹状态方程，从重要性密度函数中采样得到新的粒子集；

(4)通过压电传感器网络激励和采集导波信号，判断监测区域是否萌生新的裂纹；如果没有监测到裂纹萌生，则直接进行步骤(5)；如果监测到新的裂纹萌生，则动态更新多裂纹扩展状态方程中的裂纹分量数，以及裂纹扩展方程数目；同时更新粒子集，新增裂纹长度分量，其裂纹长度初始化为当前导波结构健康监测方法的裂纹监测阈值；更新裂纹总的裂纹数目n＝n+1，然后进行步骤(5)；

(5)对于每个粒子，计算每个粒子对应的结构损伤程度指数；从导波结构健康监测方法得到的导波信号中提取损伤因子向量，计算结构损伤程度指数的观测似然概率密度函数；通过该观测似然概率密度函数给每个粒子加权，并对所有粒子的权值进行归一化；

(6)通过粒子集和相应的权值近似各个裂纹分量的裂纹长度，以及裂纹扩展模型参数的后验概率密度函数；通过裂纹扩展模型参数的后验估计更新裂纹扩展模型；

(7)以各个裂纹的裂纹长度的后验估计作为起点，采用参数更新后的裂纹扩展模型对多个裂纹在未来时刻的裂纹扩展进行预测；

(8)根据每个粒子的归一化权值，使用系统重采样算法进行重采样，并对粒子进行正则化处理；

(9)进行到下一时刻k＝k+1，并重复步骤(3)，(4)，(5)，(6)，(7)，(8)，(9)。

步骤(1)中所述结构损伤程度指数定义为标量D，其为结构所有裂纹长度a＝[a₁,a₂,...,a_z]的函数，如下所示，

D＝ζ(a)

其中：a₁为第1个裂纹的长度，a₂为第2个裂纹的长度，a_z为第z个裂纹的长度，z取结构中最多存在的裂纹数目；ζ(·)为通过试验数据得到的映射函数，对于共线裂纹，取ζ(·)为所有裂纹的长度的线性求和；对于非共线裂纹，ζ(·)为通过试验数据训练得到的非线性映射模型；

同时在试验中获取特定损伤程度对应的导波结构健康监测数据，从导波信号中提取多种损伤因子，组成导波损伤因子向量y，建立导波向量和结构损伤程度之间的多维非线性映射关系如下式所示

y＝g(D)+ν

其中：y为导波损伤因子向量，ν为观测噪声，g(·)为多维非线性映射。

步骤(3)中从重要性密度函数中采样得到新的粒子集的步骤如下：

(a)根据结构应力强度因子手册计算多裂纹耦合条件下，每个粒子中各个裂纹分量的形状修正系数其中分别为k时刻第i个粒子的第1个，第2个，第n个裂纹分量的形状修正系数；

(b)基于各个裂纹的形状修正系数，计算每个粒子的应力强度因子幅和应力强度因峰值，获得k时刻粒子的裂纹长度分量如下所示，

其中：和分别表示k时刻第i个粒子中第1个，第2个，第n个裂纹的长度，和分别表示k-1时刻第i个粒子中第1个，第2个，第n个裂纹的长度，ω_1,k-1，ω_2,k-1，和ω_n,k-1为服从高斯分布的裂纹扩展噪声，exp(·)为自然指数函数，ΔN循环载荷数步长， 和分别为第i个粒子中第1个，第2个，第n个裂纹分量的裂纹扩展速率，其通过上述NASGRO模型计算得到；

(c)粒子中的参数向量由下式得到，

其中：为k时刻第i个粒子的裂纹扩展模型参数分量，为k-1时刻第i个粒子的裂纹扩展模型参数分量。

步骤(4)中对裂纹扩展状态方程和粒子集的更新方式如下，

如果监测到新的裂纹萌生，首先在已有n个裂纹的基础上，新增第n+1个裂纹的裂纹长度分量和裂纹扩展方程，下式所示为更新后的多裂纹扩展状态方程

其中：k为离散的时刻，n为当前时刻裂纹个数，a_1,k为k时刻第1个裂纹的裂纹长度，a_2,k为k时刻第2个裂纹的裂纹长度，a_n,k为k时刻第n个裂纹的裂纹长度，a_n+1,k为k时刻第n+1个裂纹的裂纹长度；a_1,k-1为k-1时刻第1个裂纹的裂纹长度，a_2,k-1为k-1时刻第2个裂纹的裂纹长度，a_n,k-1为k-1时刻第n个裂纹的裂纹长度，a_n+1,k-1为k时刻第n+1个裂纹的裂纹长度；ω_n+1,k-1为新增的第n+1个裂纹的裂纹扩展噪声，(da₁/dN)|_k-1，(da₂/dN)|_k-1，(da_n/dN)|_k-1和(da_n+1/dN)|_k-1分别为第1个，第2个，第n个，第n+1个裂纹的裂纹扩展速率；

然后对粒子集进行更新，每个粒子中的裂纹长度向量由更新为其中新增第n+1个裂纹的裂纹长度分量该分量的初始值设定为当前导波结构健康监测方法的裂纹监测阈值a_th，并且将变量n设置为n＝n+1。

本发明的有益效果如下：

本发明提出了一种基于动态裂纹数目的粒子滤波多裂纹扩展预测方法，采用在线导波结构健康监测方法在线监测结构中裂纹的萌生和扩展，动态更新裂纹分量和裂纹扩展方程数目，解决了单一结构中新裂纹萌生的监测以及多裂纹扩展规律表征的问题，结合了导波结构健康监测方法获取的在线多裂纹观测值，能有效用于结构的多裂纹扩展的在线预测。

附图说明

图1为本发明提出的多裂纹扩展预测方法流程图。

图2为实施例中的结构尺寸和传感器布置图。

图3为实施例中不同试件的多裂纹扩展轨迹图。

图4为实施例中典型的导波信号直达波图。

图5为实施例中从导波信号直达波中提取的散射信号能量损伤因子图。

图6为实施例中从导波信号直达波中提取的归一化互相关矩损伤因子图。

图7为实施例中从导波信号直达波中提取的能量差损伤因子图。

图8为实施例中从导波信号直达波中提取的散射信号能量损伤因子速率图。

图9为实施例中在N＝30740个载荷循环时的裂纹扩展预测结果图。

图10为实施例中在N＝56800个载荷循环时的多裂纹预测结果图。

具体实施方式

下面结合附图对本发明创造的技术方案进行详细说明。

图1所示为基于动态裂纹数目的粒子滤波多裂纹扩展预测方法流程图，其方法步骤如下，

(1)导波结构健康监测方法提取的损伤因子与结构中存在的所有裂纹长度相关，在离线状态下，定义一个标量D表征结构的损伤程度。进行相关的试验，根据试验数据建立结构中多个裂纹长度a＝[a₁,a₂,...,a_z]和结构损伤程度指数的关系，如下式所示

D＝ζ(a)

其中：a₁为第1个裂纹的长度，a₂为第2个裂纹的长度，a_z为第z个裂纹的长度，z取结构中最多存在的裂纹数目。ζ(·)为通过试验数据得到的映射函数。对于共线裂纹，取ζ(·)为所有裂纹的长度的线性求和。对于非共线裂纹，ζ(·)为通过试验数据训练得到的非线性映射模型，比如高斯过程模型，支持向量机。

同时在试验中获取特定损伤程度对应的导波结构健康监测数据。从导波信号中提取多种损伤因子，组成导波损伤因子向量y，建立导波向量和结构损伤程度之间的多维非线性映射关系如下式所示

y＝g(D)+ν

其中：y为导波损伤因子向量，ν为观测噪声，g(·)为多维非线性映射。

(2)在结构服役的初始时刻k＝0，对算法进行初始化，设定初始裂纹的个数n＝1；根据结构的几何形式和服役条件确定结构中单个裂纹的裂纹扩展状态方程，如下所示

其中：a_1,k为k时刻第1个裂纹的裂纹长度，a_1,k-1为k-1时刻第1个裂纹的裂纹长度，ω_1,k-1为裂纹扩展噪声，N为载荷循环数，ΔN为载荷循环数步长，θ_k-1为k-1时刻裂纹扩展模型参数向量，θ_k为k时刻裂纹扩展模型参数向量，(da₁/dN)|_k-1为第1个裂纹分量的裂纹扩展速率，其通过NASGRO裂纹扩展模型计算，如下所示。

其中：m，p，q为经验材料参数，a为裂纹长度，f是裂纹张开函数，R为应力比，K_max,k-1为k-1时刻的应力强度因子峰值，ΔK_k-1为应力强度因子幅，K_C和ΔK_th分别是断裂韧度和应力强度因子幅阈值，C_k-1为材料特性参数。材料特性参数C_k-1作为待估计的随机变量以表征结构间的裂纹扩展分散性，即θ_k-1＝[C_k-1]，其服从初始分布p(C)。

上式中应力强度因子幅和应力强度因子峰值由下式计算，

其中：a_k-1为k-1时刻的裂纹长度，ΔS为结构承受的应力幅值，S_max为结构承受的应力峰值，F_k-1为k-1时刻该裂纹的裂纹形状修正系数，与结构和裂纹的几何形式相关。

设定粒子个数N_s，对粒子集进行初始化，其中表示初始时刻k＝0的第i个粒子，共N_s个粒子。表示第i个粒子的裂纹长度分量，通过导波结构健康监测方法的裂纹监测阈值a_th初始化。表示第i个粒子的模型参数分量，其从裂纹扩展模型参数的初始分布p(C₀)随机采样得到此外，第i个粒子的初始权值为

(3)在结构服役过程中的k时刻(k＝1,2,3...)，从重要性密度函数中采样得到新的粒子集其中为时刻k-1的第i个粒子，y_k为k时刻提取得到的导波损伤因子向量，x_k为k时刻的状态向量。该采样过程通过如下过程实现：

首先根据结构应力强度因子手册计算多裂纹耦合条件下，每个粒子中各个裂纹分量的形状修正系数 为k时刻第i个粒子的第1个裂纹分量的形状修正系数，为k时刻第i个粒子的第2个裂纹分量的形状修正系数，为k时刻第i个粒子的第n个裂纹分量的形状修正系数。

基于各个裂纹的形状修正系数，计算每个粒子的应力强度因子幅和应力强度因峰值，获得k时刻粒子的裂纹长度分量如下所示，

其中：和分别表示k时刻第i个粒子中第1个，第2个，第n个裂纹的长度，和分别表示k-1时刻第i个粒子中第1个，第2个，第n个裂纹的长度，ω_1,k-1，ω_2,k-1，和ω_n,k-1为服从高斯分布的裂纹扩展噪声，和分别为第i个粒子中第1个，第2个，第n个裂纹分量的裂纹扩展速率，其通过上述NASGRO模型计算得到；

同时粒子中的参数向量由下式得到，

其中：为k-1时刻第i个粒子的裂纹扩展模型参数分量。

该采样过程表征上一时刻粒子的裂纹长度的演化过程。

(4)通过压电传感器网络激励和采集导波信号，判断监测区域的裂纹是否萌生。

如果监测到新的裂纹萌生，首先在已有n个裂纹的基础上，新增第n+1个裂纹的裂纹长度分量和裂纹扩展方程，下式所示为更新后的多裂纹扩展状态方程。

其中：k为离散的时刻，n为当前时刻裂纹个数，a_1,k为k时刻第1个裂纹的裂纹长度，a_2,k为k时刻第2个裂纹的裂纹长度，a_n,k为k时刻第n个裂纹的裂纹长度，a_n+1,k为k时刻第n+1个裂纹的裂纹长度；a_1,k-1为k-1时刻第1个裂纹的裂纹长度，a_2,k-1为k-1时刻第2个裂纹的裂纹长度，a_n,k-1为k-1时刻第n个裂纹的裂纹长度，a_n+1,k-1为k时刻第n+1个裂纹的裂纹长度；ω_n+1,k-1为新增的第n+1个裂纹的裂纹扩展噪声，(da₁/dN)|_k-1，(da₂/dN)|_k-1，(da_n/dN)|_k-1和(da_n+1/dN)|_k-1分别为第1个，第2个，第n个，第n+1个裂纹的裂纹扩展速率。

然后对粒子集进行更新，每个粒子中的裂纹长度向量由更新为其中新增第n+1个裂纹的裂纹长度分量该分量的初始值设定为当前导波结构健康监测方法的裂纹监测阈值a_th。并且将变量n设置为n＝n+1，并进行步骤(5)。

如果没有监测到裂纹，直接进行步骤(5)。

(5)计算每个粒子的损伤程度指数，

其中：为第i个粒子在k时刻的损伤程度指数，为第i个粒子在k时刻的裂纹长度向量。

在基于导波结构健康监测的裂纹萌生和扩展监测过程中，从导波信号中提取当前时刻与结构损伤程度相关的损伤因子向量y_k，根据导波向量和结构损伤程度之间的多维非线性映射得到当前结构损伤程度指数的观测似然函数p(y_k|x_k)，然后采用下式计算每个粒子的权值，

其中：为k时刻第i个粒子的权值，为k-1时刻第i个粒子的权值，为第i个粒子的观测似然值。计算得到所有粒子的权值后，对所有粒子的权值进行归一化处理，如下式所示，

其中：为为k时刻第i个粒子的归一化权值，

则状态向量x_k的后验概率密度函数p(x_k|y_1:k)由粒子集以及对应的归一化权值近似如下，

其中：y_1:k＝{y₁,y₂,...,y_k}，y₁和y₂分别是k＝1和k＝2时刻获取的损伤因子向量，δ为狄利克雷函数，其表达式如下，

其中：x_k为k时刻状态向量，为k时刻第i个粒子。

每个裂纹的裂纹长度后验概率密度函数表示为，

其中a_j,k(j＝1,2,...,n)为结构第j个裂纹的裂纹长度，为第i个粒子在k时刻的第j个裂纹的长度。裂纹扩展模型参数的后验概率密度函数表示为，

其中：θ_k为k时刻裂纹扩展模型参数向量，为k时刻第i个粒子的模型参数分量。

裂纹长度和裂纹扩展模型参数的最大后验概率估计由下式计算，

其中：为结构第j个裂纹的裂纹长度后验估计，为裂纹扩展模型参数的最大后验概率估计。

(6)以裂纹长度的最大后验估计结果为起点，使用裂纹扩展模型参数的最大后验概率估计更新裂纹扩展模型，将每个裂纹的裂纹长度后验估计结果向未来时刻投影得到未来每个裂纹的扩展轨迹预测结果，如下式所示，

a_j,k+p＝a_j,k+p-1+exp(ω_j,k+p-1)(da_j/dN)|_k+p-1ΔN

其中：a_j,k+p为未来k+p时刻的裂纹长度，a_j,k+p-1为未来k+p-1时刻的裂纹长度，ω_j,k+p-1为裂纹扩展噪声；计算裂纹扩展速率(da_j/dN)|_k+p-1时，裂纹扩展模型参数取裂纹扩展预测起点为裂纹长度的最大后验概率估计结果

(7)根据归一化权值对所有的粒子进行重采样，采用系统重采样算法，其步骤如下：

(a)计算所有粒子权值的累积数列M，

其中：为第1个粒子的归一化权值，为第2个粒子的归一化权值，为第i个粒子的归一化权值，为第N_s个粒子的权值。

(b)从均匀分布U(0,1)中随机采样得到随机值u

(c)初始化，设定粒子上标i＝1

(d)计算重采样权值q＝u+(i-1)/N_s

(e)顺序查找数列M，找到M中第一个大于q的元素，记其上标为τ

(f)第i个重采样后的粒子为原始粒子集中的第τ个粒子该重采样后的粒子的权值为

(g)i＝i+1，重复步骤(d)(e)(f)(g)直到i＝N_s

(8)对重采样后的粒子集进行正则化处理。正则化过程中的核密度函数选取为Epanechnikov核函数Φ(x)，其表达式如下，

其中：n_x为状态向量的维数，为单位球的体积，||x||为x的模值，同时核函数的窗宽h选取为，

其中：常数A通过下式计算

正则化步的步骤如下，

(a)计算粒子集合的经验协方差函数S_k。

(b)计算矩阵Z_k，使得Z_kZ_k ^T＝S_k，其中Z_k ^T表示矩阵的转置。

(c)令粒子上标为i＝1。

(d)从标准Epanechnikov核中采样得到随机变量α。

(e)正则化后的第i个粒子

(f)i＝i+1，重复步骤(c)(d)(e)(f)直到i＝N_s。

(g)将正则化处理后的粒子集作为新的粒子集进行下一个时刻的迭代。

(9)重复步骤(6)(7)(8)(9)(10)(11)。

本实施例中以真实的孔边裂纹试件的多裂纹扩展预测为例来具体说明本发明方法的具体实施过程。

如图2所示的AL2024铝制试件，其厚度为3mm。结构中央有直径6mm的通孔。在孔的两边有两个应力集中区域，分别为区域A和区域B。这两个区域可能产生裂纹并扩展。在区域A附近的孔边缘预制长度为2mm的切口，使得区域A的裂纹比区域B的裂纹先开始扩展。结构所承受的疲劳载荷为常幅正弦载荷，应力峰值为S_max＝30Mpa，应力比为R＝0.1。

图3所示为多个试件的疲劳裂纹扩展结果，记为试件T1，T2，T3，T4，T5。同时每个试件的区域A和区域B都萌生裂纹并扩展，记区域A的裂纹扩展为“裂纹A”，其裂纹长度为a₁，区域B的裂纹扩展为“区域B”，其裂纹长度为a₂。由于切口的存在，裂纹A从初始时刻就开始扩展，裂纹B经过一定数目的循环载荷后才开始扩展。

在本实施例中，裂纹A和裂纹B是共线裂纹，结构的损伤程度指数选取为总裂纹长度，D＝a₁+a₂。图4所示为典型的导波信号直达波随着结构损伤程度指数的变化。图5，图6，图7所示分别为从导波信号中提取的散射信号能量损伤因子、归一化互相关矩损伤因子、和能量差损伤因子。记散射信号能量损伤因子为y_1,k，归一化互相关矩损伤因子为y_2,k，能量差损伤因子为y_3,k。

在初始时刻k＝0，初始化裂纹分量数目为n＝1，其裂纹扩展状态方程定义如下，

根据文献查询该结构材料的材料参数，初始化该状态方程的参数为，m＝3.2，p＝0.25，q＝1，f＝0.342， 模型参数C的初始分布p(C)为log(C_k)|_k＝0～N(-11.621,0.1²)，裂纹扩展噪声的标准差为σ_ω＝0.5。裂纹A初始裂纹长度分布为以切口长度为均值的高斯分布N(2,0.1²)，循环载荷数步长为ΔN＝20。

依据应力强度因子手册，采用有限元方法计算不同裂纹长度组合下的裂纹形状修正系数其中为第t个训练样本的裂纹A长度，为第t个训练样本的裂纹B长度，F₁ ^(t)为第t个训练样本的裂纹A的形状修正系数，F₁ ^(t)为第t个训练样本的裂纹B的形状修正系数，S为训练样本个数。通过S个训练样本训练3层前馈人工神经网络得到下式所示的形状修正系数代理模型，

[F_1,k,F_2,k]＝NN(a_1,k,a_2,k)

其中：a_1,k为裂纹A的裂纹长度，a_2,k为裂纹B的裂纹长度，F_1,k为裂纹A的形状修正系数，F_2,k为裂纹B的形状修正系数，NN(·)为训练得到的3层前馈人工神经网络。

在本实施例中，以试件T1作为当前目标结构，试件T2，T3，T4，T5的数据作为提前获取的历史数据以训练观测方程。建立3阶多元多项式的观测模型如下所示，

式中为最小二乘法计算得到的多项式回归系数，t₁，t₂，t₃为求和的中间变量，为散射信号能量损伤因子y_1,k的t₃次方，为归一化互相关矩损伤因子y_2,k的t₂次方，为能量差损伤因子为y_3,k的t₁次方。v为零均值高斯噪声，其标准差同样通过最小二乘方法得到σ_v＝0.402。

设定粒子数为N_s＝200，初始化粒子集其中从裂纹长度的初始分布N(2,0.1²)中随机采样得到，从模型参数的初始分布N(-11.621,0.1²)中随机采样得到。然后将每个粒子的权值初始化为

在试件T1处于服役过程中，通过基于主动导波的结构健康监测方法对结构T1的裂纹萌生和扩展状态进行在线监测。每当获得新的导波信号，计算导波信号的散射信号能量损伤因子，归一化互相关矩损伤因子，以及能量差损伤因子。根据信号的归一化互相关损伤因子变化速率判断新的裂纹是否萌生，如图8所示。当萌生新的裂纹后，损伤因子速率发生突变。同时计算结构损伤程度的观测似然函数，并求解裂纹长度和裂纹扩展模型的后验估计值，在此基础上对每个裂纹进行裂纹扩展预测。图9所示为在N＝30740个载荷循环得到的裂纹A的后验估计和预测结果，此时未监测到裂纹B扩展。如图8所示，在N＝54000个载荷循环监测到新的裂纹出现，在状态方程中新增裂纹分量，并更新粒子集中新裂纹的长度分量，裂纹长度初始化为监测阈值a_th＝2mm，然后进行裂纹长度后验估计和裂纹扩展预测。图10所示为监测到新的裂纹后，在N＝56800个载荷循环得到的多裂纹扩展的后验估计和预测结果。

Claims (4)

1.一种基于动态裂纹数目的粒子滤波多裂纹扩展预测方法，其特征在于，包括如下步骤：

(1)离线状态下，基于实验数据建立结构多裂纹扩展状态和结构损伤程度指数之间的映射关系；从反映结构损伤程度变化的导波信号中提取多种损伤因子，组成多维导波损伤因子向量；建立损伤因子向量和结构损伤程度指数的多维非线性映射关系；

(2)在结构服役的初始时刻k＝0，对算法进行初始化；设定初始裂纹的个数n＝1；根据结构的几何形式和服役条件确定结构中单个裂纹的裂纹扩展状态方程；在结构上布置压电传感网络，通过导波结构健康监测方法在线监测结构疲劳裂纹的萌生和扩展；

(3)在结构服役过程中的当前时刻k，计算每个粒子中各个裂纹分量的形状修正系数；结合定义的多裂纹状态方程，从重要性密度函数中采样得到新的粒子集；

(4)通过压电传感器网络激励和采集导波信号，判断监测区域是否萌生新的裂纹；如果没有监测到裂纹萌生，则直接进行步骤(5)；如果监测到新的裂纹萌生，则动态更新多裂纹扩展状态方程中的裂纹分量数，以及裂纹扩展方程数目；并同时更新粒子集，新增裂纹长度分量，其裂纹长度初始化为当前导波结构健康监测方法的裂纹监测阈值；更新裂纹总的裂纹数目n＝n+1，然后进行步骤(5)；

(5)对于每个粒子，计算其对应的结构损伤程度指数；从导波结构健康监测方法得到的导波信号中提取损伤因子向量，计算结构损伤程度指数的观测似然概率密度函数；通过该观测似然概率密度函数给每个粒子加权，并对所有粒子的权值进行归一化；

(6)通过粒子集和相应的权值近似各个裂纹分量的裂纹长度，以及裂纹扩展模型参数的后验概率密度函数；通过裂纹扩展模型参数的后验估计更新裂纹扩展模型；

(7)以各个裂纹的裂纹长度的后验估计作为起点，采用参数更新后的裂纹扩展模型对多个裂纹在未来时刻的裂纹扩展进行预测；

(8)根据每个粒子的归一化权值，使用系统重采样算法进行重采样，并对粒子进行正则化处理；

(9)进行到下一时刻k＝k+1，并重复步骤(3)，(4)，(5)，(6)，(7)，(8)，(9)。

2.根据权利要求1所述的基于动态裂纹数目的粒子滤波多裂纹扩展预测方法，其特征在于：步骤(1)中所述结构损伤程度指数定义为标量D，其为结构所有裂纹长度a＝[a₁,a₂,...,a_z]的函数，如下所示，

D＝ζ(a)

其中：a₁为第1个裂纹的长度，a₂为第2个裂纹的长度，a_z为第z个裂纹的长度，z取结构中最多存在的裂纹数目；ζ(·)为通过试验数据得到的映射函数，对于共线裂纹，取ζ(·)为所有裂纹的长度的线性求和；对于非共线裂纹，ζ(·)为通过试验数据训练得到的非线性映射模型；

同时在试验中获取特定损伤程度对应的导波结构健康监测数据，从导波信号中提取多种损伤因子，组成导波损伤因子向量y，建立导波向量和结构损伤程度之间的多维非线性映射关系如下式所示

y＝g(D)+ν

其中：y为导波损伤因子向量，ν为观测噪声，g(·)为多维非线性映射。

3.根据权利要求1所述的基于动态裂纹数目的粒子滤波多裂纹扩展预测方法，其特征在于：步骤(3)中从重要性密度函数中采样得到新的粒子集的步骤如下：

(a)根据结构应力强度因子手册计算多裂纹耦合条件下，每个粒子中各个裂纹分量的形状修正系数其中分别为k时刻第i个粒子的第1个，第2个，第n个裂纹分量的形状修正系数；

(b)基于各个裂纹的形状修正系数，计算每个粒子的应力强度因子幅和应力强度因峰值，获得k时刻粒子的裂纹长度分量如下所示，

其中：和分别表示k时刻第i个粒子中第1个，第2个，第n个裂纹的长度，和分别表示k-1时刻第i个粒子中第1个，第2个，第n个裂纹的长度，ω_1,k-1，ω_2,k-1，和ω_n,k-1为服从高斯分布的裂纹扩展噪声，exp(·)为自然指数函数，ΔN循环载荷数步长， 和分别为第i个粒子中第1个，第2个，第n个裂纹分量的裂纹扩展速率，其通过上述NASGRO模型计算得到；

(c)粒子中的参数向量由下式得到，

其中：为k时刻第i个粒子的裂纹扩展模型参数分量，为k-1时刻第i个粒子的裂纹扩展模型参数分量。

4.根据权利要求3所述的基于动态裂纹数目的粒子滤波多裂纹扩展预测方法，其特征在于：步骤(4)中对裂纹扩展状态方程和粒子集的更新方式如下，

如果监测到新的裂纹萌生，首先在已有n个裂纹的基础上，新增第n+1个裂纹的裂纹长度分量和裂纹扩展方程，下式所示为更新后的多裂纹扩展状态方程

其中：k为离散的时刻，n为当前时刻裂纹个数，a_1,k为k时刻第1个裂纹的裂纹长度，a_2,k为k时刻第2个裂纹的裂纹长度，a_n,k为k时刻第n个裂纹的裂纹长度，a_n+1,k为k时刻第n+1个裂纹的裂纹长度；a_1,k-1为k-1时刻第1个裂纹的裂纹长度，a_2,k-1为k-1时刻第2个裂纹的裂纹长度，a_n,k-1为k-1时刻第n个裂纹的裂纹长度，a_n+1,k-1为k时刻第n+1个裂纹的裂纹长度；ω_n+1,k-1为新增的第n+1个裂纹的裂纹扩展噪声，(da₁/dN)|_k-1，(da₂/dN)|_k-1，(da_n/dN)|_k-1和(da_n+1/dN)|_k-1分别为第1个，第2个，第n个，第n+1个裂纹的裂纹扩展速率；

然后对粒子集进行更新，每个粒子中的裂纹长度向量由更新为其中新增第n+1个裂纹的裂纹长度分量该分量的初始值设定为当前导波结构健康监测方法的裂纹监测阈值a_th，并且将变量n设置为n＝n+1。