大型板结构全场高空间分辨率模态参数辨识 与损伤识别研究

念海波，曹善成，郭 宁，徐 超; NIAN Hai-bo; CAO Shan-cheng; GUO Ning; XU Chao

使用Chrome浏览器效果最佳，继续浏览，你可能不会看到最佳的展示效果，

确定继续浏览么?

复制成功，请在其他浏览器进行阅读

大型板结构全场高空间分辨率模态参数辨识与损伤识别研究 PDF

- ORCID：
念海波 ¹
✉
- ORCID：
曹善成 ^1,2
✉
- ORCID：
郭宁 ¹
- ORCID：
徐超 ^1,2

1. 西北工业大学航天学院，陕西西安 710072； 2. 西北工业大学长三角研究院，江苏太仓215400

中图分类号： O327

最近更新：2024-12-31

DOI：10.16385/j.cnki.issn.1004-4523.2024.12.012

摘要

基于数字图像相关的振动位移测量方法因其环境适应性高、非接触全场测量等优势而被广泛应用于航空航天结构的测试与监测中。然而，基于全场振动位移的损伤识别面临两方面难题：现有高空间分辨率模态分析方法计算效率低；基于模态振型的无基线损伤定位方法损伤特征提取困难且抗噪性能差。为解决上述难题，本文分别提出了一种基于核函数分解与联合主成分分析的频域模态分析方法和一种基于虚拟激励的大型板类结构损伤定位方法。对全场振动测量数据进行奇异值分解，得到含有局部损伤特征的核函数及其系数，在此基础上使用联合主成分分析改进的频域模态分析，获得高空间分辨率和高抗噪性的模态振型；将结构损伤引起的局部动态平衡方程扰动等效为一个虚拟激励力，以实现对损伤位置的识别；根据损伤特征的局部邻近性和噪声特征的空间分布稀疏性，引入层次聚类方法以优化损伤定位的精度；提出一种融合多模态信息的鲁棒损伤定位指标。数值和试验算例结果验证了本文所提方法的准确性和有效性。

关键词

模态分析; 损伤定位; 虚拟激励法; 大型板结构; 层次聚类

引言

随着制造工艺的不断发展，大型空间板结构特别是复合材料层合板因其具有较高比强度、比刚度以及良好的热稳定性和抗腐蚀性能，在航空航天领域得到广泛应用，如运载火箭和在轨航天器外壳与防护层^［

1］。但在服役过程中，复合材料层合板结构承受循环、冲击和热等多种复杂载荷，不可避免地产生损伤、力学性能受到破坏。为此，迫切需要对复合材料层合板结构进行模态参数辨识和损伤监测研究，以确保其在服役周期中的安全性和可靠性。

基于振动的结构损伤识别方法，因测量方法成熟且振动信号包含结构全局模态与局部损伤信息，被广泛地应用于土木和航空航天领域^［

2］。传统基于有限测点的动力学损伤识别研究，可通过频率特性的改变识别裂纹是否存在，但难以实现损伤的定位。近年来，非接触全场振动测量因其可提供高空间分辨率的模态振型信息^{［参考文献 3‑4}3‑4］，非常适用于基于结构特征挠度形状（characteristic deflection shapes，CDS）的损伤定位^{［参考文献 5

百度学术}5］。基于CDS的损伤定位方法通过检测损伤引起的CDS或其空间导数的局部形状变化确定损伤位置。此外，由于模态振型比其他类型的CDS^{［参考文献 6‑7}6‑7］具有更高的信噪比，更适用于损伤定位。

基于高速摄像的全场测量技术可以快速获取工程结构表面的振动位移^［

8］。但由于其空间分辨率高、采样频率高和数据采集时间长，会产生大量的测量数据^{［参考文献 3

百度学术}3，9‑10］。传统工作模态分析（operational modal analysis，OMA）方法如频域分解（frequency domain decomposition，FDD）方法，对于数千个测点的动态位移处理效率较低^{［参考文献 11‑13}11‑13］。为解决这一问题，WANG等^{［参考文献 14

百度学术}14］提出使用图像分解来压缩大量的像素化位移，同时保留重要的模态信息。此类方法通常采用切比雪夫核函数、三角函数或小波等核函数对图像进行分解，进而用少量形状描述因子来表征原始大量的全场振动位移信息。基于形状描述因子，采用传统的OMA方法可以快速高效地辨识模态参数^{［参考文献 15

百度学术}15］。然而，采用传统多项式核函数，如切比雪夫核函数，难以表征损伤引起的局部形状畸变。因此，基于传统多项式形状描述因子辨识的模态振型一般不包含损伤信息，无法用于损伤位置识别。为解决这一问题，本文提出一种自适应核函数构造方法。该方法首先对全场位移测量数据进行奇异值分解（singular value decomposition，SVD），进而采用若干较大的奇异值对应的奇异向量（核函数）来分解与表征全场位移数据，可有效保留损伤位置的畸变信息与模态参数信息。

在基于振动的损伤识别方法中，损伤位置可以通过检测振型或振型空间导数的局部奇异点来确定^［

16‑18］。但是，在缺少健康状态的参考模态振型时，局部损伤的空间特征难以被进一步量化和融合。针对这一问题，本文将虚拟激励（pseudo‑excitation，PE）方法与自适应降噪技术相结合，提出一种无基线损伤定位方法。该方法基于模态振型计算局部动态平衡方程，并通过评估损伤引起的PE来识别损伤位置^{［参考文献 19

百度学术}19］。但该方法需要计算位移场的四阶空间导数，容易受到测量噪声的影响^{［参考文献 20

百度学术}20］。针对这一问题，通常采用小波变换或空间滤波等降噪方法来减少测量噪声的影响^{［参考文献 21‑23}21‑23］，但如何选择降噪算法的正则化参数以确保损伤定位结果更准确仍然是一个关键问题。为此，本文提出一种将层次聚类和局部最小二乘拟合平滑相结合的自适应降噪技术，以提高损伤定位精度。在此基础上，可结合视觉拼接技术，实现工程中大型结构的快速损伤定位。

基于上述背景，本文针对高空间分辨率结构的离面振动数据，首先通过奇异值分解方法构造分解核函数，使用基于形状描述因子的、改进的频域分解法提取结构的模态振型。以此为基础，将结构损伤造成模态振型的局部畸变等效为一个虚拟激励力，以实现损伤定位。为进一步降低噪声影响，提高损伤定位精度，提出通过局部多项式拟合和层次聚类方法优化损伤定位结果。

1 改进的工作模态分析方法

1.1 自适应核函数构造方法

传统的核函数通常具有全局光滑的特性，且独立于获取的全场振动位移数据。此外，位移场通过传统核函数表征时容易丢失局部损伤特征。为解决这些问题，本文提出对全场振动位移数据集进行奇异值分解来构造核函数。SVD作为一种正交矩阵分解方法，实现了从原始矩阵到对角矩阵的变换，通过少量奇异值和奇异向量（主成分）即可反映原矩阵的主要特征。当使用SVD方法处理由位移场列向量组装的矩阵时，所求解的左奇异向量保留了由损伤引起的位移场局部不连续变化的特征，进而可基于奇异向量或基于奇异向量辨识的模态振型实现损伤定位。

基于高速摄像与数字图像相关（digital image correlation，DIC）方法测量的离面位移场如图1所示。通常，每个时刻的位移场可由矩阵 $Y (t) \in ℝ^{l_{1} \times l_{2}}$ 表示，其中 $l_{1}$ 和 $l_{2}$ 为测量平面x和y方向上的测量点数。

图1 基于数字图像相关方法测量的离面位移场

Fig.1 Out‑of‑plane displacement field measured by DIC method

为更高效地处理采集的位移场数据，将位移矩阵 $Y (t)$ 重新排列成一个向量，进而将所有时刻的位移数据集合成一个矩阵 $Z \in ℝ^{m \times N}$ ：

Z = \{z_{1} z_{2} \dots z_{i} \dots z_{N}\}

（1）

式中 $z_{i} \in ℝ^{m \times 1}$ 为一个长度为 $m (m = l_{1} \times l_{2})$ 的列向量。为了消除均值的影响，将矩阵的每一列减去相应列向量的均值：

\tilde{Z} = Z - \bar{Z}

（2）

其中，矩阵 $\bar{Z}$ 每一列的元素相同，其值为 $Z$ 矩阵每一列的均值。然后，对 $\bar{Z}$ 矩阵进行奇异值分解：

\tilde{Z} = U Σ V^{T}

（3）

式中 $U \in ℝ^{m \times m}$ 和 $V \in ℝ^{N \times N}$ 为两个酉矩阵； $Σ$ 为奇异值矩阵，是一个对角矩阵。为提高 $U$ 计算结果的噪声鲁棒性，式（3）可改写为：

\tilde{Z} {\tilde{Z}}^{T} = U Λ_{1} U^{T}

（4）

式中 $Λ_{1} \in ℝ^{m \times m} = Σ Σ^{T}$ ， $\tilde{Z}$ 可以通过前几个较大的奇异值与其对应的奇异向量近似表示：

\tilde{Z} \approx U_{k} Σ_{k} V_{k}^{T}

（5）

式中 $U_{k} \in ℝ^{m \times k}$ 由核函数向量组成，即 $U_{k}$ 每一列代表一个核函数向量；核函数的数量k可以通过重构的位移场与原始测量位移场的相关度确定。通过SVD方法构造的核函数 $U_{k}$ 既保留了占主导地位的空间特征挠度形状，也保留了局部损伤诱导的损伤特征，这使得局部化损伤检测成为可能。

基于振动位移场，使用构造出来的 $U_{k}$ 计算其核函数系数：

P_{k} = U_{k}^{T} \tilde{Z} = Σ_{k} V_{k}^{T}

（6）

式中 $P_{k} \in ℝ^{k \times N}$ 表示核函数系数，可将其称作形状描述因子，通过形状描述因子可以完成模态参数的求解。由于 $k ≪ m$ ，大大减少了工作模态分析方法中所需处理的数据量。

与切比雪夫多项式和泽尼克多项式等经典正交多项式核函数相比，利用奇异值分解构造的正交核函数简单且不受位移场几何形状限制。其次， $U_{k}$ 的确定是通过各奇异向量对构造原始位移场的贡献来衡量的，而传统正交核函数构造方法难以选择最优核函数类型以及核函数数量。此外，传统核函数通常是全局平滑的函数，同时滤除了噪声和局部形状特征，然而局部形状特征是损伤定位的必要信息。

1.2 基于形状描述因子的改进频域分解法

频域分解法可以很容易地确定形状描述因子 $P_{k} \in ℝ^{k \times N}$ 的模态参数。

p (t_{i}) = Φ_{S D} q (t_{i}) + n (t_{i}) \in

（7）

式中 $p (t_{i}) \in ℝ^{k \times 1} (1 \leq i \leq N)$ 为矩阵 $P_{k}$ 的一列； $Φ_{S D} \in ℝ^{k \times N}$ 表示形状描述因子的模态振型矩阵； $q (t_{i})$ 为模态坐标系的动力学响应； $n (t_{i}) \in ℝ^{k \times 1}$ 表示噪声的影响。

当 $Φ_{S D} q (t_{i})$ 与 $n (t_{i})$ 不相关时， $p (t_{i})$ 的协方差矩阵为：

R_{p p} (τ) = Φ_{S D} R_{q q} (τ) Φ_{S D}^{T} + R_{n n} (τ)

（8）

式中 $τ$ 为延迟时间。通过快速傅里叶变换得到功率谱密度（Power Spectral Density，PSD）矩阵：

S_{p p} (ω) = Φ_{S D} S_{q q} (ω) Φ_{S D}^{H} + S_{n n} (ω)

（9）

对 $S_{p p} (ω)$ 矩阵进行奇异值分解，以最大奇异值在频域坐标下的峰值点对应的频率作为辨识的共振频率：

S_{p p} (ω_{r}) = D_{r} Λ_{2} D_{r}^{H}

（10）

式中 $ω_{r}$ 为第r阶共振频率； $D_{r} = [d_{r 1}, d_{r 2}, \dots, d_{r k}]$ $\in ℝ^{k \times k}$ 为一个酉矩阵；向量 $d_{r 1} \in ℝ^{k \times 1}$ 为对角矩阵 $Λ_{2}$ 中的最大奇异值对应的奇异向量，是基于形状描述因子求解的模态振型 $φ_{r S D}$ 。

然而，只通过单个PSD矩阵的 $d_{r 1}$ 估计模态振型不能保证其准确性和鲁棒性。为解决这一问题，采用联合主成分分析方法对多个PSD矩阵同时进行联合对角分解以求解模态振型。

\begin{array}{l} S_{p p} (ω_{r + l}) = Φ_{S D} S_{q q} (ω_{r + l}) Φ_{S D}^{H} + S_{n n} (ω_{r + l}) = \\ D_{r} Λ_{2} (ω_{r + l}) D_{r}^{H} + E_{r + l} \end{array}

（11）

式中 l（l=-L，-L+1，…，L；L为正整数）表示 $ω_{r}$ 邻近频率。对于任何非零L，式（11）中的分解是超定的。且对于任意 $ω_{r + l}$ ，误差项 $E_{r + l}$ 的值不总为0。为求解该问题，将式（11）转化为一个极小化问题：

\begin{array}{l} J (D_{r}, Λ_{2} (ω_{r + l})) = \\ \sum_{l = - L}^{L} ‖S_{p p} (ω_{r + l}) - D_{r} Λ_{2} (ω_{r + l}) D_{r}^{H}‖ \end{array}

（12）

式中 $D_{r}$ 为联合求解得到的通用的酉矩阵，具有统计噪声鲁棒性。此外，利用扩展Jacobi方法^［

24］可有效求解式（12）。在这种情况下，

D_{r}

的第一列向量

d_{r 1}

为模态振型

φ_{r S D}

的抗噪估计。

进而在空间坐标下，频率 $ω_{r}$ 处的模态振型为：

φ_{r} = U_{k} d_{r 1}

（13）

式中 $φ_{r} \in ℝ^{m \times 1}$ 为结构模态振型向量，可以被重新排列为模态振型矩阵 $Ψ_{r} \in ℝ^{l_{1} \times l_{2}}$ 。其整个模态振型求解过程如图2所示。

图2 基于自适应形状描述因子的模态振型求解流程

Fig.2 Solving process of modal shape based on adaptive shape description factor

2 基于虚拟激励的损伤识别方法

2.1 损伤定位指标

以厚度为h的对称正交各向异性层合板为例，频率 $ω$ 处的简谐运动方程如下式所示：

D_{11} \frac{\partial^{4} w}{\partial x^{4}} + 2 (D_{12} + 2 D_{66}) \frac{\partial^{4} w}{\partial x^{2} \partial y^{2}} + D_{22} \frac{\partial^{4} w}{\partial y^{4}} + j C ω w (x, y, ω) - ρ h ω^{2} w (x, y, ω) = f (x, y, ω)

（14）

式中 $w (x, y, ω)$ 和 $f (x, y, ω)$ 为傅里叶变换后频率 $ω$ 处的值； $D_{i j}$ 表示板结构的刚度参数； $ρ$ 为质量密度； $C$ 为黏性阻尼系数。

当外部激励 $f (x, y, ω)$ 不存在时，对于一个完整的板结构，式（14）的左边等于零。如果损伤存在，由于局部刚度、阻尼等参数的变化，式（14）等于零将不成立，可以利用这一特性进行损伤识别。

\begin{array}{l} D I (x, y, ω) = (D_{11} - Δ_{D_{11}}) \frac{\partial^{4} w}{\partial x^{4}} + \\ 2 (D_{12} + 2 D_{66} - Δ_{D_{12}} - 2 Δ_{D_{66}}) \frac{\partial^{4} w}{\partial x^{2} \partial y^{2}} + \\ (D_{22} - Δ_{D_{22}}) \frac{\partial^{4} w}{\partial y^{4}} + j (C - Δ_{C}) ω w (x, y, ω) - \\ (ρ h - Δ_{ρ h}) ω^{2} w (x, y, ω) + j Δ_{C} ω w (x, y, ω) - \\ ∆_{ρ h} ω^{2} w (x, y, ω) + Δ_{D_{11}} \frac{\partial^{4} w}{\partial x^{4}} + \\ 2 (Δ_{D_{12}} + 2 Δ_{D_{66}}) \frac{\partial^{4} w}{\partial x^{2} \partial y^{2}} + Δ_{D_{22}} \frac{\partial^{4} w}{\partial y^{4}} \end{array}

（15）

式中 $Δ_{D_{11}}$ ， $Δ_{D_{12}}$ ， $Δ_{D_{22}}$ ， $Δ_{D_{66}}$ ， $Δ_{C}$ 和 $Δ_{ρ h}$ 分别表示损伤引起的刚度、阻尼和质量等变化。由于局部的运动学方程依旧处于平衡状态，式（15）可简化为：

\begin{array}{l} D I (x, y, ω) = Δ_{D_{11}} \frac{\partial^{4} w}{\partial x^{4}} + 2 (Δ_{D_{12}} + 2 Δ_{D_{66}}) \frac{\partial^{4} w}{\partial x^{2} \partial y^{2}} + \\ Δ_{D_{22}} \frac{\partial^{4} w}{\partial y^{4}} + j Δ_{C} ω w (x, y, ω) - Δ_{ρ h} ω^{2} w (x, y, ω) \end{array}

（16）

式中 $D I (x, y, ω)$ 适用于非共振频率和共振频率，但因共振频率处信噪比较高，通常选择共振频率处进行损伤识别。1.2节中的 $Ψ_{r} \in ℝ^{l_{1} \times l_{2}}$ 可替代式（16）中的 $w (x, y, ω)$ 用于损伤定位指标DI的计算。

由式（16）可以看出，基于虚拟激励的损伤定位方法将损伤区域刚度、阻尼和质量等物理量的变化等效为运动平衡方程中施加在损伤区域的一个虚拟力，并将其作为损伤识别指标。因此，此指标的有效性与损伤区域的大小、形状和位置无关。

对于式（16）提出的损伤定位指标 $D I (x, y, ω)$ ，仍然存在两个需要解决的基本问题： $D_{i j}$ ， $C$ 和 $ρ$ 的值未知，以及其对于噪声的鲁棒性。

2.2 结构参数未知时的损伤定位指标

根据式（16）计算损伤指数时，通常无法获得 $D_{i j}$ ， $C$ 和 $ρ$ 等结构参数值。为了解决这一问题，本文提出一种组合参数估计策略，将式（16）改写为：

\begin{array}{l} D I (x, y, ω) = c_{1} \frac{\partial^{4} w}{\partial x^{4}} + c_{2} \frac{\partial^{4} w}{\partial x^{2} \partial y^{2}} + \\ c_{3} \frac{\partial^{4} w}{\partial y^{4}} + w (x, y, ω) \end{array}

（17）

其中：

\begin{array}{l} c_{1} = \frac{Δ_{D_{11}}}{j Δ_{C} ω - Δ_{ρ h} ω^{2}}, c_{2} = \frac{2 (Δ_{D_{12}} + 2 Δ_{D_{66}})}{j Δ_{C} ω - Δ_{ρ h} ω^{2}}, \\ c_{3} = \frac{∆_{D_{22}}}{j Δ_{C} ω - Δ_{ρ h} ω^{2}} \end{array}

（18）

式中 $c_{i} (i = 1,2, 3$ ）为待定系数。

通常情况下，损伤区域远小于整个检测区域，因此 $D I (x, y, ω) = 0$ 适用于结构的绝大部分位置，可以通过最小二乘法对 $c_{i}$ 进行估计。在确定 $c_{i}$ 后，可以计算出式（17）中定义的损伤定位指标 $D I (x, y, ω)$ 完成损伤定位。但由于测量噪声的影响，对于未损伤区域， $D I$ 的值也不严格等于零，一般在零值附近波动。

为减小测量噪声对 $D I (x, y, ω)$ 的影响，本文提出采用基于局部区域最小二乘拟合的方法对式（17）模态振型高阶导数进行求解。若直接进行四阶导数求解，其拟合的形函数需要包括 $x^{4}, x^{3} y, x^{2} y^{2}, x y^{3}, y^{4}$ 等多个高阶项，拟合的系数为15个，在数据不足的情况下会产生不确定项，影响高阶导数计算精度。因此，提出通过计算两次二阶导数来避免拟合系数过多的现象，二阶导数计算的形函数如下式所示：

w (x, y) = a_{0} + a_{1} x^{2} + a_{2} y^{2} + a_{3} x y + a_{4} x + a_{5} y

（19）

二阶导数计算结果与形函数系数的关系为：

\{\begin{matrix} \frac{\partial^{2} w}{\partial x^{2}} = 2 a_{1} \\ \frac{\partial^{2} w}{\partial y^{2}} = 2 a_{2} \\ \frac{\partial^{2} w}{\partial x \partial y} = a_{3} \end{matrix}

（20）

通过对 $\frac{\partial^{2} w}{\partial x^{2}}, \frac{\partial^{2} w}{\partial y^{2}}$ 再次进行二阶导数计算，可以得到 $\frac{\partial^{4} w}{\partial x^{4}}, \frac{\partial^{4} w}{\partial x^{2} \partial y^{2}}, \frac{\partial^{4} w}{\partial y^{4}}$ 。四阶导数的计算精度以及对噪声的平滑程度与所选取的局部区域大小有关，通过局部区域大小参数 $s$ 来描述局部区域的大小，如图3所示。

图3 导数计算区域大小参数表示

Fig.3 Parameter representation of derivative calculation area size

通过调节参数 $s$ 的值可以降低测量噪声对DI的影响。但如何确定合适的 $s$ 值以获得更准确的DI是一个重要的问题。

此外，采用 $3 σ$ 准则去除 $D I (x, y, ω; s)$ 中的小值，如下式所示：

D I = \{\begin{array}{l} |D I|, |D I| \geq \bar{D I} + 3 σ_{D I} \\ 0, |D I| < \bar{D I} + 3 σ_{D I} \end{array}

（21）

式中 $\bar{D I}$ 和 $σ_{D I}$ 分别为 $D I (x, y, ω; s)$ 的平均值和标准差。

在式（21）中，对于整个结构， $D I (x, y, ω; s)$ 中的元素通常位于 $[\bar{D I} - 3 σ_{D I}, \bar{D I} + 3 σ_{D I}]$ 范围内。但由于测量噪声影响，一些超出此范围的值随机分散在二维表面上，如图4所示。在结构损伤处，损伤引起更多的超出此范围的值出现，如图4所示，这与测量噪声相关的离群值有显著的空间分布区别。

图4 损伤区域与测量噪声离群值分布

Fig.4 Distribution of damage area and outlier values of measurement noise

2.3 自适应降噪方法

$D I (x, y, ω)$ 的损伤定位性能更多地取决于局部区域大小参数 $s$ 。因此，将损伤定位问题转化为优化问题，以提高损伤定位精度为目标，确定最优的s值。在这种情况下，需要定义一个损伤定位精度指标来衡量损伤定位性能的精度。众所周知，局部损伤是空间连续的，如图4所示。因此，本文利用损伤的空间邻近特性，通过层次聚类算法评估损伤定位精度。

该方法采用单关联聚类算法计算 $D I (x, y, ω)$ 任意两个离群值簇之间的最小距离：

d (C_{1}, C_{2}) = \underset{c_{1} \in C_{1}, c_{2} \in C_{2}}{m i n} d (c_{1}, c_{2})

（22）

式中 $C_{1}$ 和 $C_{2}$ 表示两个不同的离群值簇； $d (C_{1}, C_{2})$ 表示 $C_{1}$ 和 $C_{2}$ 两个离群值簇之间的欧氏距离。另外，任意两个离群值簇之间的距离可以用树状图表示，如图5所示。

图5 基于树状图的损伤区域与测量噪声离群值分布

Fig.5 Distribution of damage area and outlier values of measurement noise based on dendrogram

图5中，高度表示 $D I (x, y, ω)$ 的离群值之间的距离。可以看出，大型簇包含大量空间上接近的离群值，这通常是由损伤引起的。相比之下，测量噪声引起的离群值往往在空间上分散，相应的簇仅由少数离群值组成。因此，层次聚类可以用于评估损伤识别精度。为量化损伤定位的精度，提出以噪声诱导离群值簇的数量来指示损伤定位的精度，数量越少精度越高。此外，本文认为簇中离群值个数 $N_{C} \leq 3$ 的簇是由测量噪声引起的。

随着参数 $s$ 的增大，噪声引起的簇的数量逐渐减少，因为当 $s$ 较大时，测量噪声引起的离群值得到了平滑。因此， $s$ 的选取可以通过 $N_{C} \leq 3$ 的离群值聚类的个数来确定。为了避免损伤诱发的局部特征因 $s$ 取值过大而被平滑掉，可以自动选择最小 $s$ 值，使得 $N_{C} \leq 3$ 离群值聚类个数小于给定数值 $N_{o}$ 即完成损伤定位。本文取 $N_{o} = 3$ ，即当最小的 $s$ 值使得 $N_{o} \leq 3$ 时达到损伤定位精度的优化要求。

2.4 多模态损伤指标融合

单一模态振型无法包含全部的损伤信息，因此通过多阶模态损伤定位指标叠加可以更加精确和鲁棒地实现损伤定位。

D I (x, y) = \frac{1}{n} \sum_{i = 1}^{n} D I (x, y, ω_{i})

（23）

式中 $n$ 表示使用模态振型的数量。

综上，基于虚拟激励方法的损伤定位计算流程如图6所示。

图6 基于虚拟激励方法的损伤识别流程图

Fig.6 Damage identification flow chart based on PE method

3 数值仿真验证

本节首先验证基于形状描述因子的模态振型提取方法的可行性，其次对所提出的损伤定位方法的有效性进行验证。

建立一个含分层损伤的四层碳/环氧树脂复合材料板模型，并进行有限元仿真，如图7所示，其左侧固支。板的尺寸为 $(0.5 \times 0.5 \times 0.008)$ $m^{3}$ ，每层铺层厚度相等，均为0.002 m，铺层顺序为［0°，90°，90°，0°］，密度为1500 kg/m³，所用材料各个方向上的工程参数为： $E_{1} = 109.34 G P a$ ， $E_{2} = 8.82 G P a$ ， $E_{3} = 8.82 G P a$ ， $G_{12} = 4.32 G P a$ ， $G_{13} = 4.32 G P a$ ， $G_{23} = 3.20 G P a$ ， $ν_{12} = 0.342$ ， $ν_{13} = 0.342$ ， $ν_{23} = 0.520$ 。方形分层损伤区域位于板结构对称面单侧方向上0°与90°两层之间，损伤区域长度和宽度分别为0.05 m和0.05 m，位于板结构的正中心。在ABAQUS中采用2478个SC8R单元进行网格划分，得到了较好的收敛解。将振幅为5 N的脉冲激励力 $F$ 加在点（0.4013 m，0.4209 m）处。在采样率为2000 Hz的情况下，采集了上表面1307个单元节点2 s的振动位移数据用于工作模态分析。

图7 复合材料分层损伤示意图

Fig.7 Schematic diagram of layered damage of composite material

为验证形状描述因子方法提取的模态振型和基于虚拟激励的损伤定位方法具有噪声鲁棒性，对仿真得到的振动响应添加噪声水平 $n_{l e v e l} = 0.1 %$ （信噪比60 dB）的高斯白噪声，使用第1节提出的改进工作模态振型计算方法计算其前6阶弯曲模态振型，如图8所示。此外，图8中所示模态振型为通过工作模态分析方法获得的工作模态振型，并非有限元等数值方法计算得出的质量归一化的模态振型。

图8 基于自适应形状描述因子的模态振型识别结果

Fig.8 Identification results of modal shapes based on adaptive shape description factor

以第3阶模态振型为例，如果不使用局部区域数据进行最小二乘拟合区域求解4阶导数，而直接采用4阶中心差分进行计算，得到的损伤识别结果如图9所示。

图9 不使用降噪处理的损伤定位结果

Fig.9 Damage localization results using PE method without denoising technique

由图9可以看出，损伤定位结果受到测量噪声的严重污染，无法得到有效的结果。为了减少噪声影响，采用2.3节提出的降噪方法，通过对局部区域大小参数 $s$ 进行寻优计算，图10（a）为参数 $s$ 优化计算后（ $s = 6$ ）的损伤定位结果。图10（b）为参数 $s$ 与 $N_{C} \leq 3$ 的群落个数之间的关系。局部区域大小参数 $s$ 的值可以显著减少测量噪声引起的随机聚类。

图10 损伤定位优化

Fig.10 Damage localization optimization

使用前6阶弯曲模态振型分别进行损伤定位，得到结果如图11所示，其中第1，3，9，10阶模态振型可以较为准确地定位损伤，第4，6阶模态振型因对该区域损伤的敏感程度有限，无法完成损伤区域定位。因此单一模态振型无法提供鲁棒的损伤定位，需要融合多阶模态振型的损伤特征，实现鲁棒的损伤定位。

图11 损伤定位结果

Fig.11 Damage localization results

为实现鲁棒的损伤定位，这里将前6阶弯曲模态振型中的损伤信息融合起来，形成一个集合多阶模态振型的损伤指数，多模态融合损伤识别结果如图12所示。结果表明，融合多阶模态振型损伤指数能够更准确和有效地完成损伤定位。

图12 多阶模态振型融合损伤定位结果

Fig.12 Integrated damage localization results based onmultiple order modal shapes

最后，为了证明传统形函数作为形状描述因子在损伤定位中难以应用，使用二维切比雪夫核函数对全场位移图进行分解。选取前20个二维切比雪夫核函数用以重建振动位移场，基于所提出的方法识别第3阶模态振型及其损伤定位结果如图13所示，通过与图10（a）的损伤识别结果对比表明，将传统多项式作为形状描述因子会使得模态振型丢失局部损伤信息，所得到的模态振型难以用于损伤识别。

图13 基于切比雪夫核函数的第3阶模态振型及其损伤定位结果

Fig.13 The 3rd order modal shape and its damage localization results based on Chebyshev kernel function

4 试验验证

为通过试验验证本文所提基于核函数分解与联合主成分分析的频域模态分析方法和基于虚拟激励的损伤定位方法的有效性与可行性，使用简单金属板结构验证本文提出的“振动位移‑模态振型‑损伤定位”整个识别流程的可行性；为节省时间与成本，使用公开数据集验证所提损伤定位方法在复合材料板上的适用性。

4.1 铝板结构损伤定位试验

为验证所提基于形状描述因子的模态分析方法和基于模态振型自适应降噪的损伤定位方法的有效性，首先开展金属板结构损伤定位试验，试验设计如图14（a）所示。

图14 基于高速摄像的板结构振动试验

Fig.14 Vibration experiment of a cantilever plate structure based on stereo‑cameras

板结构试验采用单边固支的铝板作为试验试件，整体尺寸为410 mm×220 mm×3 mm，不计固支部分的尺寸为330 mm×220 mm×3 mm，在距板结构右侧自由端220 mm位置处设置裂纹损伤（位于测量面背面），裂纹损伤长220 mm，宽2 mm，深2 mm，如图14（b）所示。在测量平面制作散斑图案来增强特征，用以提高振动测量精度，本文散斑点直径大小为2.5 mm，如图15所示。

图15 测量区域空间位置和网格划分

Fig.15 The spatial position and meshing of measurement area

双目相机系统的左相机与右相机均采用FASTCAM Mini AX200高速相机，采用外部TTL信号同步触发，采样频率为4000 Hz，单帧图片分辨率为1024 pixel×1024 pixel。对板结构以锤击的方式施加脉冲激励，采用双目相机系统同步拍摄图像，采样时长为1.8 s。

在计算过程中，需要将测量区域划分为若干子区，通过不同图像子区的匹配完成振动位移的测量。因此，网格划分的密度决定工作模态振型的空间测量密度。过密的网格点不仅会增加计算量，而且会增大测量噪声对高阶位移导数的影响（可通过降噪方法进行改进），过于稀疏的网格密度又存在无法捕捉到局部损伤信息的缺点。因此，在计算效率满足要求的前提下，建议采取较密的测点，以实现损伤的高精度定位。本次试验测量区域大小为285 mm×168 mm，网格密度为60×36，测量区域的空间位置和网格划分如图15所示。

通过双目测量系统完成板结构的振动测量，然后基于1.2节中的改进模态分析方法完成结构模态振型计算，其前4阶弯曲模态振型如图16所示。

图16 前4阶弯曲模态振型计算结果

Fig.16 Calculation results of the first four order bending modal shapes

基于第2节提出的虚拟激励方法，使用前4阶弯曲模态振型完成损伤定位，结果如图17所示，图17（a）识别结果再次表明单一模态振型有时无法有效地定位损伤位置。因此，将前4阶弯曲模态振型损伤信息组合起来，多模态融合的鲁棒损伤识别结果如图18所示。

图17 损伤定位结果

Fig.17 Damage localization results

图18 多阶模态振型融合损伤定位结果

Fig.18 Integrated damage localization results based on multiple order modal shapes

图18表明了本文所提多阶模态振型融合损伤定位方法的有效性。此外，本文所提方法可应用于多种损伤形状、分布位置的损伤检测。

4.2 复合材料层合板损伤识别算例

为了验证基于高空间分辨率模态振型的虚拟激励损伤定位方法在复合材料层合板损伤定位的可行性和有效性，本文使用KATUNIN^［

25］所给出的试验模态测量数据进行测试。

试验试件为尺寸大小为（300×300×2.5） mm³的12层玻璃纤维增强环氧基复合材料层合板，铺层顺序为［0°，60°，-60°，60°，0°］，如图19所示。数据集中包含如下损伤：（1）贯穿平面的裂纹损伤，损伤宽度为1 mm，损伤深度为0.5 mm，裂纹距离板结构边缘最近处为100 mm；（2）正方形区面损伤，损伤深度为0.5 mm；（3）多种损伤，包括一个线形域损伤及两个面损伤，损伤深度均为0.5 mm。试件均采用四边固支，试验中实际测量区域大小为（250×250） mm²，即上下左右各有25 mm区域未进行测量，更为详细的试验信息可参考文献［

25］。

图19 复合材料损伤试件^［

25］（单位：mm）

Fig.19 Damaged specimen of composite material^［

25］（Unit： mm）

将第1~5阶模态振型中的损伤信息组合起来，多模态融合损伤识别结果如图20所示。图20中为实际测量区域尺寸，相较于图19上下左右各有25 mm的减少量。

图20 基于虚拟激励方法的损伤定位结果

Fig.20 Damage localization results based on PE method

图20表明基于前5阶模态振型损伤信息融合的损伤指标可有效识别三个损伤案例中的损伤区域，且适用于多类型与多损伤区域的识别，进而验证本文所提基于虚拟激励的损伤识别方法在复合材料中损伤定位的可行性与有效性。

5 结论

基于高空间分辨率的振动位移场数据，本文研究了无基准结构参数损伤定位的两个关键问题：鲁棒模态振型估计和自适应无基准结构参数的损伤定位方法。对于鲁棒模态振型估计，提出了自适应核函数的构造方法和基于形状描述因子的改进频域分解法，在保留损伤特征的前提下有效提取高空间分辨率模态振型。对于损伤定位，在没有基准结构参数的情况下，采用虚拟激励方法识别损伤诱发局部特征以进行损伤定位。此外，将损伤定位问题转化为优化问题，利用层次聚类和局部最小二乘拟合平滑方法有效降低了测量噪声的影响。其他结论如下：

（1）通过奇异值分解法获得的核函数可有效保留结构局部损伤信息，基于形状描述因子的改进模态方法可以有效地处理数千个测点的振动响应，实现含损伤特征的高分辨率模态振型的辨识。

（2）单阶模态振型可能对某些损伤部位不敏感，无法鲁棒地实现损伤定位，需融合多阶模态振型中的损伤信息，才能保证鲁棒的损伤定位结果。

（3）降噪技术对损伤定位的鲁棒性和准确性起着至关重要的作用。

（4）基于虚拟激励的损伤定位方法仅与损伤区域的物理参数退化有关，适用于板结构裂纹、分层和质量损失等多种损伤形式的损伤定位。

参考文献

王惠芬，杨碧琦，刘刚. 航天器结构材料的应用现状与未来展望［J］. 材料导报， 2018， 32（增刊1）： 395-399. [百度学术]

WANG Huifen， YANG Biqi， LIU Gang. Application status and future prospect of materials for spacecraft structures［J］. Materials Reports， 2018， 32（Sup1）： 395-399. [百度学术]

赵建刚，张玉祥，陈家照，等. 结构损伤识别的广义均布载荷面法［J］. 振动工程学报， 2021， 34（5）： 987-994. [百度学术]

ZHAO Jiangang， ZHANG Yuxiang， CHEN Jiazhao， et al. Method of generalized uniform load surface for structure damage identification［J］. Journal of Vibration Engineering， 2021， 34（5）： 987-994. [百度学术]

BEBERNISS T J， EHRHARDT D A. High-speed 3D digital image correlation vibration measurement： recent advancements and noted limitations［J］. Mechanical Systems and Signal Processing， 2017， 86： 35-48. [百度学术]

NERI P， PAOLI A， RAZIONALE A V， et al. Low-speed cameras system for 3D-dic vibration measurements in the kHz range［J］. Mechanical Systems and Signal Processing， 2022， 162： 108040. [百度学术]

MOLINA-VIEDMA A J， PIECZONKA L， MENDROK K， et al. Damage identification in frame structures using high-speed digital image correlation and local modal filtration［J］. Structural Control & Health Monitoring， 2020， 27（9）： e2586. [百度学术]

HU Y J， GUO W G， ZHU W D， et al. Local damage detection of membranes based on Bayesian operational modal analysis and three-dimensional digital image correlation［J］. Mechanical Systems and Signal Processing， 2019， 131： 633-648. [百度学术]

XU Y F， CHEN D M， ZHU W D， et al. Delamination identification of laminated composite plates using measured mode shapes［J］. Smart Structures and Systems， 2019， 23（2）： 195-205. [百度学术]

陈亚军，孙胜洁，季春明. 三维数字图像相关技术（3D DIC）在材料形变研究中的应用进展［J］. 航空材料学报， 2017， 37（4）： 90-100. [百度学术]

CHEN Yajun， SUN Shengjie， JI Chunming. Development and application of 3D digital image correlation （3D DIC） in deformation measurement of materials［J］. Journal of Aeronautical Materials， 2017， 37（4）： 90-100. [百度学术]

LI M， FENG G， DENG R， et al. Structural vibration mode identification from high-speed camera footages using an adaptive spatial filtering approach［J］. Mechanical Systems and Signal Processing， 2022， 166： 108422. [百度学术]

李建宇，任朝，张丽丽，等. 基于变形场测量数据主元压缩的模型参量反求方法［J］. 计算力学学报， 2020， 37（2）： 226-232. [百度学术]

LI Jianyu， REN Zhao， ZHANG Lili， et al. Inverse method of model parameters based on deformation measurement data using principal component compression［J］. Chinese Journal of Computational Mechanics， 2020， 37（2）： 226-232. [百度学术]

GORJUP D， SLAVIČ J， BOLTEŽAR M. Frequency domain triangulation for full-field 3D operating-deflection-shape identification［J］. Mechanical Systems and Signal Processing， 2019， 133： 106287. [百度学术]

王世东，任伟新. 应变功率谱密度传递比与工作应变模态参数识别［J］. 振动工程学报， 2022， 35（4）： 806-813. [百度学术]

WANG Shidong， REN Weixin. Strain power spectrum density transmissibility and operational strain modal parameter identification［J］. Journal of Vibration Engineering， 2022， 35（4）： 806-813. [百度学术]

陈永高，钟振宇. 桥梁结构模态参数识别中系统阶次的自动辨识算法［J］. 振动工程学报， 2021， 34（4）： 680-689. [百度学术]

CHEN Yonggao， ZHONG Zhenyu. Automatic identification algorithm of system order in modal parameter identification for bridge structures［J］. Journal of Vibration Engineering， 2021， 34（4）： 680-689. [百度学术]

WANG W Z， MOTTERSHEAD J E. Adaptive moment descriptors for full-field strain and displacement measurements［J］. The Journal of Strain Analysis for Engineering Design， 2013， 48（1）： 16-35. [百度学术]

CHANG Y H， WANG W Z， CHANG J Y， et al. Compressed sensing for OMA using full-field vibration images［J］. Mechanical Systems and Signal Processing， 2019， 129： 394-406. [百度学术]

PANDEY A K， BISWAS M， SAMMAN M M. Damage detection from changes in curvature mode shapes［J］. Journal of Sound and Vibration， 1991， 145（2）： 321-332. [百度学术]

DAHAK M， TOUAT N， KHAROUBI M. Damage detection in beam through change in measured frequency and undamaged curvature mode shape［J］. Inverse Problems in Science and Engineering， 2019， 27（1-4）： 89-114. [百度学术]

贺文宇，何健，任伟新. 基于间接法识别的桥梁振型的损伤定位方法［J］. 振动与冲击， 2018， 37（24）： 13-17. [百度学术]

HE Wenyu， HE Jian， REN Weixin. A damage localization method based on indirectly identified mode shapes［J］. Journal of Vibration and Shock， 2018， 37（24）： 13-17. [百度学术]

XU H， ZHOU Q， CAO M S， et al. A dynamic equilibrium-based damage identification method free of structural baseline parameters： experimental validation in a two-dimensional plane structure［J］. Journal of Aerospace Engineering， 2018， 31（6）： 4018081. [百度学术]

ROY K， RAY-CHAUDHURI S. Fundamental mode shape and its derivatives in structural damage localization［J］. Journal of Sound and Vibration， 2013， 332（21）： 5584-5593. [百度学术]

SOLÍS M， ALGABA M， GALVÍN P. Continuous wavelet analysis of mode shapes differences for damage detection［J］. Mechanical Systems and Signal Processing， 2013， 40（2）： 645-666. [百度学术]

CAO M S， QIAO P Z. Integrated wavelet transform and its application to vibration mode shapes for the damage detection of beam-type structures［J］. Smart Materials and Structures， 2008， 17（5）： 055014. [百度学术]

YANG Z B， RADZIENSKI M， KUDELA P， et al. Fourier spectral-based modal curvature analysis and its application to damage detection in beams［J］. Mechanical Systems and Signal Processing， 2017， 84： 763-781. [百度学术]

CAO S C， OUYANG H J. Robust structural damage detection and localization based on joint approximate diagonalization technique in frequency domain［J］. Smart Materials and Structures， 2016， 26（1）： 015005. [百度学术]

KATUNIN A. Nondestructive damage assessment of composite structures based on wavelet analysis of modal curvatures： state-of-the-art review and description of wavelet-base damage assessment benchmark［J］. Shock and Vibration， 2015： 735219. [百度学术]