CN105912504A - 一种机械结构多源冲击载荷识别的压缩感知方法 - Google Patents

Abstract

本发明涉及一种机械结构多源冲击载荷识别的压缩感知方法，用于解决高度欠定系统的多源冲击载荷识别反问题的病态特性。该方法包括以下步骤：1)测量机械结构冲击载荷作用点与机械结构响应点间的频响函数，进而构造感知矩阵；2)采用传感器测量由结构动载荷产生的信号；3)构造多源冲击载荷识别的欠定方程；4)构造基于L1范数的多源冲击载荷识别的压缩感知凸优化模型；5)利用两步迭代阈值算法求解压缩感知优化模型，获得多源冲击载荷的压缩感知解。本发明充分利用冲击载荷的时间和空间的联合稀疏性，适用于识别和定位作用在机械结构的多源冲击载荷，克服了传统的基于L2范数的正则化方法无法求解欠定系统的瓶颈。

Description

一种机械结构多源冲击载荷识别的压缩感知方法

技术领域

本发明属于机械系统状态监测领域，具体涉及一种机械结构多源冲击载荷识别的压缩感知方法。

背景技术

冲击载荷是一类非常重要的动载荷，特别是在复合材料的结构健康监测中。比如用于风力发电的风机叶片在运行和维护过程中，不可避免的遭受到风沙、飞鸟、冰雹、维修工具等外来物的冲击，并且发生频率较高，冲击损伤积累下来会对风机叶片的完整性以及承载能力造成安全隐患；作用于飞机机翼的突变气流、飞鸟、维修工具对机翼的撞击。实时监控这类冲击信号是非常必要的，而直接测量它们是非常困难的。针对冲击载荷这类特殊动载荷，精确地识别其位置甚至比重构时间历程更加重要。传统的冲击载荷位置识别方法比较直观，即采用枚举的策略，在每一个冲击载荷虚拟位置均采用正则化方法识别动载荷，其中响应残差最小的虚拟位置即是冲击载荷作用的真实位置。该方法存在明显缺陷：(1)工作量巨大，每个假定位置都需要应用正则化方法计算一次；(2)无法识别多源冲击载荷；(3)载荷时间历程和位置识别是分离的，不能同时实现重构和定位。

发明内容

基于此，本发明公开了一种机械结构多源冲击载荷识别的压缩感知方法，所述方法包括以下步骤：

S100、测量机械结构冲击载荷作用点j与机械结构响应点i间的频响函数，进而构造感知矩阵；

S200、对机械结构施加冲击载荷并测量冲击载荷的响应y_i；

S300、基于步骤S100和步骤S200构造多源冲击载荷识别的欠定方程；

S400、基于步骤S200中的冲击载荷的响应和步骤S300的欠定方程构造基于L1范数的多源冲击载荷识别的压缩感知凸优化模型；

S500、求解压缩感知凸优化模型，获得多源冲击载荷的压缩感知解。

本发明与现有技术相比具有下列优点：

1.传统的基于L2范数的截断奇异值分解、Tikhonov正则化、函数逼近法等均无法求解欠定系统下的载荷识别问题，本发明充分利用冲击载荷的稀疏性，可获得欠定系统下的唯一解，突破长期以来载荷识别领域欠定系统无法求解的瓶颈。

2.本发明解决了冲击载荷时间历程和位置无法同时识别的难题，彻底改变了传统L2范数框架下的载荷识别思路。

3.采用的两步迭代阈值法具有不涉及传递矩阵求逆运算和不需要明确正则化参数的优点，其迭代过程即是正则化的过程，可高效求解多源冲击载荷的压缩感知模型，且对振动系统初始条件不敏感。

附图说明

图1是本发明一个实施例的方法流程图；

图2是本发明一个实施例中薄壳结构多源冲击载荷识别装置示意图；

图3(a)、图3(b)是本发明一个实施例中力传感器和加速度传感器实测数据，其中，图3(a)为冲击载荷信号，图3(b)为测点的加速度响应；

图4(a)、图4(b)是本发明一个实施例中薄壳结构双源冲击载荷产生的不同噪声水平的加速度响应，其中，图4(a)为第7和8次冲击响应；图4(b)为第13和14次冲击响应；

图5(a)、图5(b)是本发明一个实施例中压缩感知算法TwIST识别的薄壳结构两个真实和七个虚拟载荷，其中，图5(a)为第7和8次连续冲击的压缩感知结果；图5(b)为第13和14次连续冲击的压缩感知结果。

具体实施方式

下面结合附图1-5及具体实施例对本发明作进一步描述，应该强调的是，下述说明仅仅是示例性的，而本发明的应用对象不局限下述示例。

在一个实施例中，本发明公开了一种机械结构多源冲击载荷识别的压缩感知方法，所述方法包括以下步骤：

S100、测量机械结构冲击载荷作用点j与机械结构响应点i间的频响函数，进而构造感知矩阵；

S200、对机械结构施加冲击载荷并测量冲击载荷的响应y_i；

S300、基于步骤S100和步骤S200构造多源冲击载荷识别的欠定方程；

S400、基于步骤S200中的冲击载荷的响应和步骤S300的欠定方程构造基于L1范数的多源冲击载荷识别的压缩感知凸优化模型；

S500、求解压缩感知凸优化模型，获得多源冲击载荷的压缩感知解。

本实施例充分利用冲击载荷的稀疏性，获得欠定系统的唯一解，突破了长期以来载荷识别领域欠定系统无法求解的瓶颈，解决了冲击载荷时间历程和位置无法同时识别的难题，彻底改变了传统L2范数框架下的载荷识别思路。

在一个实施例中，所述步骤S100具体包括以下步骤：

S101、测量机械结构冲击载荷作用点j与机械结构响应点i间的频响函数H_ij(ω)；

S102、对频响函数采用快速傅里叶逆变换得到单位脉冲响应函数h_ij(t)；H_ij(ω)

S102、对单位脉冲响应函数h_ij(t)进行离散化获得传递矩阵H_ij，构造感知矩阵H＝[H_il，…H_ij，…H_iN]；其中，N表示虚拟冲击载荷的数目，ω表示圆频率变量，t表示时间变量。

本实施例中，频响函数的测量方法主要有力锤锤击法和激振器激励法，其中锤击法相对方便，本实施例中优先选择锤击法测量频响函数。

在一个实施例中，采用传感器测量步骤S200中施加于机械结构上的冲击载荷所产生的冲击载荷的响应y_i。

在本实施例中，采用加速度传感器测量作用于机械结构的动载荷产生的相应信号，也可采用速度、位移或应变传感器测量振动响应。

在一个实施例中，所述步骤S300中多源冲击载荷识别的欠定方程为：

$y_i=\left[\begin{array}{ccccc}{H}_{i1},& \mathrm{...}& H_{ij},& \mathrm{...}& H_{iN}\end{array}\right]\left[\begin{array}{c}{f}_1\\ .\\ .\\ .\\ f_j\\ .\\ .\\ .\\ f_N\end{array}\right]$

上式的矩阵-矢量形式为

y＝Hf

其中，f表示多源冲击载荷的压缩感知矢量，所述f仅在1-2个位置有真实载荷而在其他位置元素全部是零，则矢量f是稀疏的。

在本实施例中，由于感知矩阵H的严重病态和高度欠定，且测量响应y中含有噪声，直接对上式采用最小二乘法算法求逆是不可行的，为此，本实施例通过添加稀疏约束条件，构造压缩感知模型，在L1范数下实现上式的求解。

在一个实施例中，所述步骤S400中压缩感知凸优化模型为：

$\underset{f}{\min imize}G\left(f\right)=\frac{1}{2}\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1;$

其中，||·||₂表示向量的L2范数，||·||₁表示向量的L1范数，λ表示正则化参数，H感知矩阵，y为冲击载荷的响应矢量。

在本实施例中，上式第一项具有最小化响应残差的作用，第二项||f||₁具有诱导稀疏解的作用，而正则化参数λ在最小化响应残差和最小化解非零个数||f||₁之间建立了平衡关系。

在一个实施例中，所述步骤S500中利用两步迭代阈值算法求解凸优化模型，具体包括以下步骤：

S501、初始化：

S502、更新当前解f_k+1：

S503、判断当前解f_k+1是否满足如下条件：

G(f_k+1)＞G(f_k)

若满足则继续进入步骤S504，否则进入步骤S505；

S504、单调化处理当前解f_k+1；

S505、判断当前解f_k+1是否满足如下收敛准则：

$\frac{|G\left(f_{k+1}\right)-G\left(f_k\right)|}{G\left(f_{k+1}\right)}\≤\mathrm{\ϵ}$

其中，ε表示终止阈值；

若当前解f_k+1满足上式迭代终止准则，则终止迭代过程，则当前解即为多源冲击载荷的压缩感知解；否则，迭代过程返回步骤S502继续迭代，直到满足上式。

在本实施例中，上式终止准则表示相邻的两次迭代结果f_k和f_k+1的目标函数变化不大时，即可终止迭代过程。其中，终止阈值ε越小，所需的迭代步数越多，耗时也越多，然而过小的终止阈值也是没有必要，在这里ε＝10^-7即可获得高精度的解。

在一个实施例中，所述步骤S502中利用两步迭代阈值算法更新当前解：

f_k+1＝βS_λ(f_k+H^T(y-Hf_k))+(α-β)f_k+(1-α)f_k-1

其中，加权常数α＝1.9874、β＝3.9748。

在本实施例中，将加权常数设置为α＝1.9874和β＝3.9748，是为了平衡相邻迭代解和软阈值滤波结果。本实施例采用的两步迭代阈值法具有不涉及传递矩阵求逆运算和不需要明确正则化参数的优点，其迭代过程即是正则化的过程，可高效求解多源冲击载荷的压缩感知模型，且对振动系统初始条件不敏感。

在一个实施例中，所述步骤S504中利用下式单调化处理当前解f_k+1：

f_k+1＝S_λ(f_k+1)

其中S_λ表示软阈值滤波函数是标准的小波降噪函数。

在一个实施例中，所述步骤S501具体包括：

令非负初始解f₀＝[0，....，0]^T∈Rⁿ、终止阈值ε＝10^-7、正则化参数λ＝0.02||H^Ty||_∞；其中，||·||_∞表示无穷大范数；

初始化k＝1时的解f₁：

f₁＝S_λ(f₀+H^T(y-Hf₀))

其中，S_λ表示软阈值滤波函数是标准的小波降噪函数，定义如下：

$S_{\mathrm{\λ}}\left(f_i\right)=\left\{\begin{array}{cc}\mathrm{sgn}\left(f_i\right)\left(\right|f_i|-\mathrm{\λ}),& \left|f_i\right|>\mathrm{\λ}\\ 0,& otherwise\end{array}\right..$

在一个实施例中，本发明提供一种机械结构多源冲击载荷识别的压缩感知方法，用于解决传统冲击载荷识别方法无法求解载荷识别欠定系统的困境，以高精度和高效地同时识别多源冲击载荷的时间历程和位置。

为达到上述目的，本发明采用的技术方案是，一种机械结构多源冲击载荷识别的压缩感知方法，其特征在于，将当期科学界和工程界广泛关注的压缩感知理论应用到载荷识别领域，采用高效的两步迭代阈值算法求解压缩感知模型，该方法具体包括如下步骤：

1)测量频响函数和构造感知矩阵。采用锤击法测量机械结构冲击载荷作用点j与机械结构响应点i间的频响函数H_ij(ω)，通过快速傅里叶逆变换(Inverse Fast FourierTransform，IFFT)得到单位脉冲响应函数h_ij(t)，进而离散化获得传递矩阵H_ij，构造感知矩阵H＝[H_i1，…H_ij，…H_iN]。其中，N表示虚拟冲击载荷的数目，ω表示圆频率变量，t表示时间变量；

2)施加冲击载荷和测量加速度冲击响应，采用加速度传感器测量由作用于机械结构的冲击载荷产生的响应信号y_i。

3)用一个响应测点数据识别多个位置冲击载荷的，构造多源动载荷识别的欠定方程：

$y_i=\left[\begin{array}{ccccc}{H}_{i1},& \mathrm{...}& H_{ij},& \mathrm{...}& H_{iN}\end{array}\right]\left[\begin{array}{c}{f}_1\\ .\\ .\\ .\\ f_j\\ .\\ .\\ .\\ f_N\end{array}\right]---\left(1\right)$

具紧凑表达形式为：

y＝Hf； (2)

其中：观测矩阵H∈R^n×Nn(N≥2)为严重病态欠定矩阵，y为单一加速度信号。多源载荷矢量f中仅在有限位置有载荷信号而在其他位置元素全部是零，则矢量f是稀疏的，即激励源在时间和空间的分布均是稀疏的。然而，当振动系统是欠定时，传统L2范数框架下的正则化方法的解是不唯一的。

4)针对多源冲击载荷识别，基于冲击载荷时间和空间的联合稀疏性，构造如下的压缩感知模型G(f)：

$\underset{f}{\min imize}G\left(f\right)=\frac{1}{2}\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1;---\left(3\right)$

其中，||·||₂表示向量的L2范数，||·||₁表示向量的L1范数，λ表示正则化参数；上式第一项具有最小化响应残差的作用，第二项||f||₁具有诱导稀疏解的作用，而正则化参数λ在最小化响应残差和最小化解非零个数||f||₁之间建立了平衡关系。

5)正如广泛所熟知的基于L2范数的正则化方法如Tikhonov、TSVD以及函数逼近法均无法求解载荷识别的欠定方程。相反，基于L1范数的正则化方法，在稀疏约束条件下，可以确定欠定系统的唯一稀疏解。两步迭代阈值算法(Two-Step Iterative Shrinkage/Thresholding，TwIST)求解多源冲击载荷压缩感知模型，其具体有如下步骤：

初始化：非负初始解f₀＝[0，....，0]^T∈Rⁿ、终止阈值ε＝10^-7、正则化参数λ＝0.02||H^Ty||_∞。

其中，||·||_∞表示无穷大范数；初始化k＝1时的解f₁：

f₁＝S_λ(f₀+H^T(y-Hf₀))； (4)

其中，软阈值滤波函数是标准的小波降噪函数，定义如下：

$S_{\mathrm{\λ}}\left(f_i\right)=\left\{\begin{array}{cc}\mathrm{sgn}\left(f_i\right)\left(\right|f_i|-\mathrm{\λ}),& \left|f_i\right|>\mathrm{\λ}\\ 0,& othenvise\end{array}\right.;---\left(5\right)$

步骤5.1)：两步迭代阈值算法更新当前解f_k+1：

f_k+1＝βS_λ(f_k+H^T(y-Hf_k))+(α-β)f_k+(1-α)f_k-1 (6)

其中，加权常数设置为α＝1.9874和β＝3.9748，用于平衡相邻迭代解和软阈值滤波结果。

步骤5.2)：判断当前解f_k+1是否满足如下条件

G(f_k+1)＞G(f_k)；(9)

若满足则继续进入步骤5.3)，否则进入步骤5.4)。

步骤5.3)：单调化处理当前解f_k+1：

f_k+1＝S_λ(f_k+1)； (10)

步骤5.4)：判断当前解f_k+1是否满足如下收敛准则：

$\frac{G\left(f_{k+1}\right)-G\left(f_k\right)}{G\left(f_{k+1}\right)}\≤\mathrm{\ϵ};---\left(11\right)$

若当前解f_k+1满足上式迭代终止准则，则终止迭代过程，获得多源冲击载荷的压缩感知解f；否则，迭代过程返回步骤5.1)继续迭代计算，直到满足上式。其中，ε＝10^-7可以保证高精度地重构未知动载荷。

上式终止准则表示相邻的两次迭代结果f_k和f_k+1的目标函数变化不大时，即可终止迭代过程。其中，终止阈值ε越小，所需的迭代步数越多，耗时也越多，然而过小的终止阈值也是没有必要。

在一个实施例中，图1是一种机械结构多源冲击载荷识别的压缩感知方法的流程图，该方法充分利用多源冲击载荷的时间和空间的联合稀疏性，构建多源冲击载荷欠定系统的压缩感知模型，通过两步阈值迭算法对模型求解，实现了冲击载荷时间历程和位置同时识别的目的，具体步骤如下：

1)测量频响函数和构造感知矩阵。采用锤击法测量机械结构冲击载荷作用点j与机械结构响应测点i间的频响函数H_ij(ω)，通过快速傅里叶逆变换(Inverse Fast FourierTransform，IFFT)得到单位脉冲响应函数h_ij(t)，进而离散化获得传递矩阵H_ij，构造感知矩阵H＝[H_i1，…H_ij，…H_iN]。其中，N表示虚拟冲击载荷的数目，ω表示圆频率变量，t表示时间变量；

11)试验示意图如图2所示，薄壳结构基本尺寸参数：圆心角90°、长500mm、半径200mm、厚5mm。薄壳结构表面均匀分布了九个节点，冲击力可作用在九个节点的任意一个或者两个位置。型号为PCB 086C01(灵敏度：12.37mV/N)力锤用于施加冲击载荷，其锤击端内嵌力传感器可保证实时测量冲击力，并将之作为参考信号计算识别载荷的相对误差。一个加速度传感器位于1#节点(型号：PCB 333B32，灵敏度：100mV/g)，实时测量冲击响应；

12)脉冲力锤重复敲击每个虚拟节点五次，同时由LMS SCADASIII数据采集系统同步记录冲击力和加速度信号，获得每个虚拟节点j于响应测点i＝1间的频响函数H_1j(ω)，进而离散化获得传递矩阵H_1j，构造感知矩阵H＝[H₁₁，…H_1j，…H₁₉]；

13)测量系统频响函数时的采样频率为2048Hz，采样时间为1s，数据长度为2050，构造的感知矩阵的维数为2050×18450。可知，该多源冲击载荷识别反问题是严重欠定的；

2)施加冲击载荷和测量加速度冲击响应，采用加速度传感器测量由作用于机械结构的冲击载荷产生的响应信号y_i；

21)来自脉冲力锤的26个连续冲击力交替作用在薄壳结构的5#和8#节点。冲击信号与1#节点的加速度响应信号由LMS SCADASIII数据采集模块同步采样，如图3(a)、图3(b)所示。作用在5#节点的13个冲击力分别被标记为(F1，F3，…，F25)，作用在8#节点的13个冲击力分别被标记为(F2，F4，…，F26)。采样频率同样是2048Hz，连续冲击载荷持续的时间为40s。从图3(a)可知，连续冲击载荷在时域内是非常稀疏的，仅仅在冲击加载区有较大峰值力。对应的1#节点的26个冲击响应如图3(b)所示，每个冲击响应都是快速震荡衰减的，部分相邻冲击响应存在“交叠”现象；

22)本案例考虑两个真实和七个虚拟冲击载荷识别，包含两个相邻的完整的冲击响应的信号(长度为2050)从图3(b)中截取，同时由压缩感知算法TwIST识别作用在九个节点的冲击载荷。以第7、8个相邻冲击和第13、14个相邻冲击为例，对应的加速度冲击响应如图4(a)和(b)所示。从图4(a)可知，两个峰值力分别在t＝0.1567s和t＝0.7856s时达到最大值，且在区间[0，0.1567]s包含有色噪声；从图4(b)可知，两个峰值力分别在t＝0.1924s和t＝0.8096s时达到最大值，且在区间[0，0.1924]s包含有色噪声。可知，两个案例均不满足载荷识别控制方程的零初始条件的假设。

3)用1#响应测点数据识别九个位置冲击载荷的，可构造如下的多源动载荷识别的欠定方程：

$y_1=\left[\begin{array}{ccccc}{H}_{11},& \mathrm{...}& H_{1j},& \mathrm{...}& H_{19}\end{array}\right]\left[\begin{array}{c}{f}_1\\ .\\ .\\ .\\ f_j\\ .\\ .\\ .\\ f_9\end{array}\right];---\left(1\right)$

其紧凑表达形式为：

y＝Hf； (2)

其中：观测矩阵H∈R^2050×18450为严重欠定矩阵，y为单一加速度信号。多源载荷矢量f中仅f₅或者f₈有载荷信号而在其他位置元素全部是零，则矢量f是稀疏的，即激励源在时间和空间的分布均是稀疏的。然而，当振动系统是欠定时，传统L2范数框架下的正则化方法的解是不唯一的。

4)针对多源冲击载荷识别，基于冲击载荷时间和空间的联合稀疏性，构造如下的压缩感知模型G(f)：

$\underset{f}{\min imize}G\left(f\right)=\frac{1}{2}\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1;;---\left(3\right)$

其中，||·||₂表示向量的L2范数，||·||₁表示向量的u范数，λ表示正则化参数；上式第一项具有最小化响应残差的作用，第二项||f||₁具有诱导稀疏解的作用，而正则化参数λ在最小化响应残差和最小化解非零个数||f||₁之间建立了平衡关系。

5)正如广泛所熟知的基于L2范数的正则化方法如Tikhonov、TSVD以及函数逼近法均无法求解载荷识别的欠定方程。相反，基于L1范数的正则化方法，在稀疏约束条件下，可以确定欠定系统的唯一稀疏解。两步迭代阈值算法(Two-Step Iterative Shrinkage/Thresholding，TwIST)求解多源冲击载荷压缩感知模型，其具体有如下步骤：

初始化：非负初始解f₀＝[0，....，0]^T∈Rⁿ、终止容差ε＞0、正则化参数λ＝0.02||H^Ty||_∞。

其中，||·||_∞表示无穷大范数；初始化k＝1时的解f₁：

；f₁＝S_λ(f₀+H^T(y-Hf₀)) (4)

其中，软阈值滤波函数是标准的小波降噪函数，定义如下：

$S_{\mathrm{\λ}}\left(f_i\right)=\left\{\begin{array}{cc}sgn\left(f_i\right)\left(\right|f_i|-\mathrm{\λ}),& \left|f_i\right|>\mathrm{\λ}\\ 0,& otherwise\end{array}\right.;---\left(5\right)$

步骤5.1)：两步迭代阈值算法更新当前解f_k+1：

f_k+1＝βS_λ(f_k+H^T(y-Hf_k))+(α-β)f_k+(1-α)f_k-1 (6)

其中，加权常数设置为α＝1.9874和β＝3.9748，用于平衡相邻迭代解和软阈值滤波结果。

步骤5.2)：判断当前解f_k+1是否满足如下条件

G(f_k+1)＞G(f_k)； (9)

若满足则继续进入步骤5.3)，否则进入步骤5.4)

步骤5.3)：单调化处理当前解f_k+1：

f_k+1＝S_λ(f_k+1)； (10)

步骤5.4)：判断当前解f_k+1是否满足如下收敛准则：

$\frac{|G\left(f_{k+1}\right)-G\left(f_k\right)|}{G\left(f_{k+1}\right)}\≤\mathrm{\ϵ};---\left(11\right)$

若当前解f_k+1满足上式迭代终止准则，则终止迭代过程，获得多源冲击载荷的压缩感知解f；否则，迭代过程返回步骤5.1)继续迭代计算，直到满足上式。其中，ε＝10^-7可以保证高精度地重构未知动载荷。

6)为了定量评价压缩感知算法识别冲击载荷的精度，定义识别的冲击载荷峰值相对误差为：

$\frac{\left|\right|\max \left(f_{measured}\right)-\max \left(f_{identified}\right)|{|}_2}{\left|\right|\max \left(f_{measured}\right)|{|}_2}\×100\%---\left(12\right)$

其中，f_measured和f_identified分别是实测的冲击载荷和应用压缩感知方法重构的冲击载荷。

61)压缩感知算法TwIST应用图4(a)、图4(b)中的响应识别的多源冲击载荷识别结果如图5(a)和图5(b)所示。可知，两种案例中，作用在薄壳结构的第7、8个冲击和第13、14个冲击分别被精准地定位在5#和8#节点，而其他七个位置的虚拟载荷非常小而可忽略不计；

62)至于冲击载荷的峰值力精度，第7、8个冲击载荷实测峰值力分别为66.44N和66.78N，对应压缩感知解的峰值力分别为58.74N和65.33N，则峰值相对误差仅仅为8.85％和2.17％；第13、14个冲击载荷实测峰值力分别为97.31N和105.50N，对应压缩感知解的峰值力分别为91.66N和102.60N，则峰值相对误差仅为5.81％和2.75％。

可知，多源冲击力时间历程被压缩感知算法TwIST从高度不完整不精确的测量中高精度地重构，同时冲击力位置也被压缩感知算法TwIST精准地定位。

以上所述仅为本发明的较佳实施例而己，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种机械结构多源冲击载荷识别的压缩感知方法，其特征在于，所述方法包括以下步骤：

S100、测量机械结构冲击载荷作用点j与机械结构响应点i间的频响函数，进而构造感知矩阵；

S200、对机械结构施加冲击载荷并测量冲击载荷的响应y_i；

S300、基于步骤S100和步骤S200构造多源冲击载荷识别的欠定方程；

S400、基于步骤S200中的冲击载荷的响应和步骤S300的欠定方程构造基于L1范数的多源冲击载荷识别的压缩感知凸优化模型；

S500、求解压缩感知凸优化模型，获得多源冲击载荷的压缩感知解。

2.根据权利要求1所述的方法，其特征在于，所述步骤S100具体包括以下步骤：

优选的，S101、测量机械结构冲击载荷作用点j与机械结构响应点i间的频啊函数H_ij(ω)；

S102、对频响函数H_ij(ω)采用快速傅里叶逆变换得到单位脉冲响应函数h_ij(t)；

S103、对单位脉冲响应函数h_ij(t)进行离散化获得传递矩阵H_ij；

S104、构造感知矩阵H＝[H_i1，…H_ij，…H_iN]；其中，N表示虚拟冲击载荷的数目，ω表示圆频率变量，t表示时间变量。

3.根据权利要求2所述的方法，其特征在于：采用传感器测量步骤S200中施加于机械结构上的冲击载荷所产生的冲击载荷的响应y_i。

4.根据权利要求3所述的方法，其特征在于：所述步骤S300中多源冲击载荷识别的欠定方程为：

$y_i=\[\begin{array}{ccccc}{H}_{i1},& \mathrm{...}& H_{ij},& \mathrm{...}& H_{iN}\end{array}\]\left[\begin{array}{c}{f}_1\\ .\\ .\\ .\\ f_j\\ .\\ .\\ .\\ f_N\end{array}\right]$

上式的矩阵-矢量形式为：

y＝Hf

其中，f表示多源冲击载荷的压缩感知矢量，

5.根据权利要求4所述的方法，其特征在于：所述步骤S400中压缩感知凸优化模型为：

$\underset{f}{\min imize}G\left(f\right)=\frac{1}{2}\left|\right|Hf-y|{|}_2^2+\mathrm{\λ}|\left|f\right|{|}_1;$

其中，||·||₂表示向量的L2范数，||·||₁表示向量的u范数，λ表示正则化参数，H表示感知矩阵，y表示冲击载荷的响应矢量。

6.根据权利要求5所述的方法，其特征在于：所述步骤S500中利用两步迭代阈值算法求解凸优化模型，具体包括以下步骤：

S501、初始化：

S502、更新当前解f_k+1：

S503、判断当前解f_k+1是否满足如下条件：

G(f_k+1)＞G(f_k)

若满足则继续进入步骤S504，否则进入步骤S505；

S504、单调化处理当前解f_k+1；

S505、判断当前解f_k+1是否满足如下收敛准则：

$\frac{|G\left(f_{k+1}\right)-G\left(f_k\right)|}{G\left(f_{k+1}\right)}\≤\mathrm{\ϵ}$

其中，ε表示终止阈值；

若当前解f_k+1满足上式迭代终止准则，则终止迭代过程，则当前解即为多源冲击载荷的压缩感知解；否则，迭代过程返回步骤S502继续迭代，直到满足上式。

7.根据权利要求6所述的方法，其特征在于，所述步骤S501具体包括：

令非负初始解f₀＝[0，....，0]^T∈Rⁿ、终止阈值ε＝10^-7、正则化参数λ＝0.02||H^Ty||_∞；其中，||·||_∞表示无穷大范数；

初始化k＝1时的解f₁：

f₁＝S_λ(f₀+H^T(y-Hf₀))

其中，S_λ表示软阈值滤波函数是标准的小波降噪函数，定义如下：

$S_{\mathrm{\λ}}\left(f_i\right)=\left\{\begin{array}{cc}\mathrm{sgn}\left(f_i\right)\left(\right|f_i|-\mathrm{\λ}),& \left|f_i\right|>\mathrm{\λ}\\ 0,& otherwise\end{array}\right..$

8.根据权利要求7所述的方法，其特征在于：所述步骤S502中利用两步迭代阈值算法更新当前解：

f_k+1＝βS_λ(f_k+H^T(y-Hf_k))+(α-β)f_k+(1-α)f_k-1

其中，加权常数α＝1.9874、β＝3.9748。

9.根据权利要求6所述的方法，其特征在于：所述步骤S504中利用下式单调化处理当前解f_k+1：

f_k+1＝S_λ(f_k+1)

其中S_λ表示软阈值滤波函数。

10.根据权利要求3所述的方法，其特征在于：所述传感器包括：加速度传感器、位移传感器、速度传感器或应变传感器。