Lamb波作为一类在薄板结构中传播的弹性导波,常用于检测大型板状结构的缺陷以及评价结构健康状况,在无损检测与结构健康监测领域具有重要的理论研究价值与工程应用前景[1-2]。近半个世纪以来,基于Lamb波的缺陷诊断技术已成为多学科交叉研究的热点领域。相较于体波,导波在传播距离、检测灵敏度及覆盖范围等方面展现出显著优势[3-4],因而在板壳类结构的缺陷表征中更具适用性。
在超声导波结构健康监测技术体系中,弹性导波由换能器激励产生,并以被测结构作为波导介质进行传播[5]。当导波与结构缺陷相互作用时,会引发散射效应、能量衰减及模式转换等物理现象,此类缺陷信号可通过分布式传感器阵列采集。相较于无缺陷状态下的直达波信号,缺陷散射信号通常呈现幅值微弱(通常低于噪声基底10~20 dB)、信噪比低等特征。此外,导波固有的频散特性(相速度与群速度的频率依赖性)及多模态耦合效应(如S0与A0型的共存)进一步加剧了信号解析的复杂性。缺陷特征参量主要包括飞行时间、波包能量衰减率及模态转换系数等,其中飞行时间与传播路径几何参数(群速度、传播距离)的关联性为缺陷定位算法提供了理论依据。
缺陷成像作为超声导波结构健康监测技术体系的核心环节,旨在通过多源信号融合重构缺陷的空间分布信息[6]。延迟求和(delay-and-sum, DAS)方法是由澳大利亚国防科学与技术组织下属飞行器部的Wang以及美国斯坦福大学的Chang等学者于2004年提出[7],现已成为导波缺陷成像中的经典方法之一[8]。美国佐治亚理工学院的Michaels等[9-10]用DAS成像法对铝板结构中的通孔缺陷进行了缺陷成像研究。伦敦帝国理工学院的Khodaei 和Aliabadi[11]对DAS成像方法在铝板结构及复合材料结构中的缺陷检测进行了数值和实验研究。针对复合材料层合板中群速度各向异性的问题,中国飞机强度研究所的王莉等[12]对传统DAS成像法中的速度项进行了修正,提出了群速度修正的DAS成像法,即考虑了群速度的各向异性,实现了更高精度的缺陷成像。然而,由于上述DAS方法均针对等厚度板结构,对变厚度板并不适用。出于承载能力及轻量化等需求,变厚度板已广泛应用于飞机机翼、船舶纵梁等结构,其缺陷检测对保障服役安全至关重要,但目前对于变厚度板的Lamb波缺陷检测研究仍较为缺乏。
对于Lamb波在变厚板中的传播问题,相关学者已经进行了一些有益探索。Fu 和Marical等研究了Lamb波在高斯截面变化的弹性板中的行为,在变截面内观察了不同模式Lamb波的透射和反射特性[13-14]。Moll等[15-16]分析了线性变厚度各向异性波导中的反对称模式导波,发现其符合绝热波特征且其波速随板厚变化而变化,从而导致明显的缺陷定位误差。Deng等[17]研究了Lamb波在具有不同力学参数的锥形黏接层的复合结构中的传播特性。Predoi等[18]通过有限元方法研究了厚度线性变化的铝板中Lamb波的传播特性,重点分析了绝热模式假设的适用性与坡度角对波数变化的影响:当坡度角小于3°时,S0模式的局部波数k(x)与等厚度板的频散曲线高度吻合,支持绝热假设(即波数仅依赖局部厚度);坡度角超过3°时,波数出现显著峰值,表明S0模式向高阶模式转换,且转换点与模态截止频率无关;坡度角增大至5.14°时,波数波动加剧,显示模态耦合效应增强。de Marchi等[19-20]提出了锥形波导中Lamb波群延迟预测的方法,并在已知Lamb波传播路径的条件下实现了变厚度板中Lamb波的频散补偿。Pagneux等[21]提出了一种变厚度板中Lamb波传播模拟方法,在此基础上Moreau等[22]提出了一种测量变厚度波导中超声导波频散曲线的方法。Zhang等[23]观察并计算了Lamb波在线性变化厚度的铝板中的声阴影区。Kang等[24-25]基于扫描激光多普勒测振法获得了Lamb波数据,利用声波波数谱法重构了波导结构的厚度分布图像。文献[26-28]报道了几种不同几何结构变厚度铝板中的Lamb波传播行为,并对其中Lamb波的飞行时间进行了预测,其研究结果表明阶梯式板的厚度分布也可以根据可变厚度板中的Lamb波传播模型来确定。虽然已有研究对变厚度板中Lamb波的传播特性等进行了有益探索,但目前针对变厚度板的Lamb波缺陷检测,尤其是缺陷成像的研究仍鲜见报道,其主要原因在于变厚度板中的Lamb波在传播过程中受波导结构厚度的影响,使得其波速等传播特性不断改变,故常规等厚度波导中的Lamb波成像方法在变厚度板中不再适用。
因此,针对实际应用中变厚度板结构的缺陷检测需求,本文构建了新的缓变厚度板中Lamb波传播模型,通过反向传播相位补偿算法,消除了变厚度板中Lamb波频散效应导致的信号畸变问题,在此基础上,结合延时求和成像算法实现缓变厚度板的缺陷成像。
1 缓变厚度板中Lamb波缺陷成像
1.1 缓变厚度板中Lamb波特性
1.2 缓变厚度板中Lamb波传播模型
$u(t, x)=\frac{A(x, \omega)}{2 \pi} \int_{-\infty}^{+\infty} U(\omega, 0) \mathrm{e}^{\mathrm{i}[\omega t-k(\omega) x]} \mathrm{d} \omega 。$
式中:u(t,x)为接收处时域响应信号;U(ω,0)为频域激励信号;i为虚单位;ω为角频率;x为Lamb波传播距离;A(x,ω)为幅值因子,它与角频率和传播距离相关联;k(ω)为角频率ω处的波数,其值既取决于板的厚度,也取决于频率。
在式(1)中,U(ω,0)为激励信号u(t,0)的傅里叶变换。如果不考虑Lamb波传播过程中的幅度衰减,则上述方程可以表示为
U(ω,x)=U(ω,0)e-ik(ω)x。
式中,U(ω,x)是响应信号u(t,x)的傅里叶变换。从式(2)可以看出,只要知道传播距离x、波数k(ω)和激励位置的信号,就可以预测响应信号波形。此外,基于式(2),当已知x1位置响应信号时,就可以预测x2处的响应信号波形
U(ω,x2)=U(ω,x1)e-ik(ω)(x2-x1)。
Lamb波在等厚度板中的传播模型如上所述,厚度沿传播方向不发生变化时模型的构建较为简单,若厚度沿传播方向发生变化时,上述传播模型不再适用,现考虑Lamb波在变厚度板中的传播情况。
单模态Lamb波在变厚度板中沿x方向传播,z方向为板厚方向,变厚度板沿y方向的尺寸为无限大。对于厚度沿传播方向发生变化的板,波导可以近似为一系列不同厚度的微小等宽微元。对于某一特定微元,其沿传播方向保持恒定的厚度,长度为Δx,如图3所示。 在不考虑模式转换的情况下, 单模态Lamb波在这种变厚度板中的传播可以看作是波沿着一系列不同厚度的微小等宽微元的传播。沿着波的传播方向,可以把前一微元输出的响应信号看作是后一微元的输入激励信号。当Lamb波在厚度发生变化的区域传播时,会同时产生反射波和透射波。Lamb波在相邻的两个不同厚度的微元之间可能会经历多次反射和透射,反射和透射系数都是未知的。因此,为简单起见,只考虑第一次透射波的传播。在下文提出的缺陷成像方法中,将响应信号中包含的反射波和其他透射波视为噪声。由于波在每个微元中的传播可以用式(3)来预测,因此只要x足够小,就可以预测任意微元的响应波形。这样,可以将第N个微元的输出响应信号[30]表示为
式中:k(ω,dj)是厚度为dj的第j个微元内Lamb波的波数;αi是表示透射波能量比例的系数; $\stackrel{-}{k}$ (ω)为广义波数,定义为
当k(ω,dj)不随微元变化时,就意味着板的厚度沿着波的传播方向没有发生变化,广义波数 $\stackrel{-}{k}$ (ω)等于Nk(ω),公式(4)会被简化为公式(2)。当每个微元对应的波数已知,则可以预测任意位置的输出响应信号。因为每个Lamb波模式的频散曲线都是频率和厚度积的函数,所以任意厚度下对应的未知波数曲线可以利用对波数曲线进行一维插值得到
k(ω,dj)=finterp1[ω·d,k(ω,d),ω·dj]。
式中:finterp1是一维插值函数;(ω,k(ω,d))是板厚等于d时频散曲线上的点;(ω,k(ω,dj))是板厚为dj时频散曲线上的点[31]。
1.3 缓变厚度板缺陷成像方法
1.3.1 缓变厚度板中缺陷散射信号传播模型
基于Lamb波传播模型,当已知波的传播路径,就可以预测变厚板中缺陷引起的散射信号。基于缺陷信号预测模型,通过从接收信号中提取的散射信号并进行频散补偿,可构建缺陷成像方法[30]。
考虑具有缺陷的变厚度板,将其视为点状散射源(如图4所示),基于pitch-catch模式对缺陷引起的散射信号进行建模,图中黑色虚线为两点间连线方向,蓝色实线是基于Snell定律计算所得到的Lamb波真实传播路径[32],基于该传播路径上板的厚度信息,再利用式(4)就可以实现Lamb波响应信号的预测。假设变厚度板上有一个激励点和一个接收点,考虑到它们的尺寸比波的传播路径小得多,它们都可以被看作是点状换能器。这里不考虑直达波反射和边界反射,只考虑从激励位置到缺陷点再到接收位置的散射波。与等厚板中散射信号建模的情况类似,缺陷引起的散射信号的过程可以分为2部分:第1部分是从激励位置向缺陷传播的波,第2部分是从缺陷位置向接收位置传播的波。在不考虑缺陷与Lamb波相互作用过程中的模式转换和非线性效应的情况下,变厚板中缺陷引起的散射信号可以建模为
$U_{\mathrm{S}-\mathrm{D}-\mathrm{R}}(\omega)=\alpha U(\omega, 0) \mathrm{e}^{-\mathrm{i} \Delta x\left[\bar{k}_{\mathrm{S}-\mathrm{D}}(\omega)+\bar{k}_{\mathrm{D}-\mathrm{R}}(\omega)\right]} 。$
式中:US-D-R(ω)为频域散射信号,下标S、D、R分别表示激励、缺陷、接收;α表征缺陷散射系数,α参数越大,缺陷越严重; ${\stackrel{-}{k}}_{S-D}$ (ω)和 ${\stackrel{-}{k}}_{D-R}$ (ω)分别是Lamb波从激励位置传播到缺陷位置和从缺陷位置传播到接收位置的广义波数。当波的传播路径已知时, ${\stackrel{-}{k}}_{S-D}$ (ω)和 ${\stackrel{-}{k}}_{D-R}$ (ω)可以由式(5)计算得到。在式(7)中,传播路径进行离散化时,使用相同的单元长度Δx。
由于式(4)中角频率ω与广义波数 $\stackrel{-}{k}$ (ω)之间的非线性关系,因此响应信号中的相位与频率和距离有关。与Lamb波在等厚板中传播的情况类似,群延迟是频率的非线性函数。在等厚度板中,可以使用频散补偿方法来消除非线性关系,从而使多个散射信号在各自飞行时间的移位后处于同相。如前所述,在变厚度板中记录的信号中的相位畸变无法通过常规频散补偿方法消除。如果传播路径已知,则可以将散射信号从接收处反向传播到缺陷点,再反向传播到激励位置。这个反向传播过程可以表示为
$U_{\mathrm{b}}(\omega)=U_{\mathrm{S}-\mathrm{D}-\mathrm{R}}(\omega) \mathrm{e}^{\mathrm{i} \Delta x\left[\bar{k}_{\mathrm{R}-\mathrm{D}}(\omega)+\bar{k}_{\mathrm{D}-\mathrm{S}}(\omega)\right]}$。
式中, ${\stackrel{-}{k}}_{R-D}$ (ω)和 ${\stackrel{-}{k}}_{D-S}$ (ω)分别是Lamb波从接收位置反向传播到缺陷位置和从缺陷位置反向传播到激励位置的总波数。
根据式(5)的定义,广义波数满足可逆性原则, ${\stackrel{-}{k}}_{S-D}$ (ω)= ${\stackrel{-}{k}}_{D-S}$ (ω), ${\stackrel{-}{k}}_{T-D}$ (ω)= ${\stackrel{-}{k}}_{D-T}$ (ω)。因此,将式(7)代入式(8),反向传播后的信号可简化为
Ub(ω)=αU(ω,0)。
式(9)表明反向传播信号的波形与正向的波形相似,两者只相差系数α,通过反向传播信号可以消除频散,使得传感器阵列在反向传播后采集到的所有信号都是同相的。由于求和可以放大缺陷引起的信号成分,同时抑制干扰,所以两个信号可以直接加和来指示缺陷情况。
需要注意,式(8)中的反向传播过程需要接收到的散射信号的传播路径。但由于缺陷位置未知,散射信号的传播路径无法确定。为解决此问题,将感兴趣的成像区域内的每个接收点作为潜在的缺陷源,代入式(8)。如果输入的接收点与缺陷位置重合,则得到的反向传播信号与式(9)中的Ub(ω)相同,如果输入的定位点与实际缺陷位置不重合,则得到的反向传播信号为
$U_{\mathrm{b}}(\omega)=\alpha U(\omega, 0) \mathrm{e}^{-\mathrm{i} \Delta x\left[\bar{k}_{\mathrm{S}-\mathrm{D}}(\omega)-\bar{k}_{\mathrm{P}-\mathrm{S}}(\omega)+\bar{k}_{\mathrm{D}-\mathrm{R}}(\omega)-\bar{k}_{\mathrm{R}-\mathrm{P}}(\omega)\right]} \text { 。 }$
式中, ${\stackrel{-}{k}}_{P-S}$ (ω)和 ${\stackrel{-}{k}}_{R-P}$ (ω)分别为Lamb波从潜在缺陷位置传播到激励位置和从接收位置传播到潜在缺陷位置的广义波数。由于潜在缺陷位置与实际缺陷位置不一致,则 ${\stackrel{-}{k}}_{S-D}$ (ω)≠ ${\stackrel{-}{k}}_{P-S}$ (ω), ${\stackrel{-}{k}}_{D-R}$ (ω)≠ ${\stackrel{-}{k}}_{R-P}$ (ω)。在这种情况下,式(10)中反向传播信号的相位Ub(ω)并不总是与激励的相位相同。由于同相信号之和的幅值很小,所以可以用反向传播的异相信号来表征缺陷。
1.3.2 缓变厚度板中缺陷定位
基于缓变厚度板中Lamb波传播模型以及缺陷信号预测模型,可以进行缺陷定位。该方法的核心在于将缓变厚度板中的Lamb波传播路径离散化为一系列不同厚度的微小等宽微元,利用频散曲线数据,通过式(6)进行线性插值计算每个微元的局部波数,并沿传播路径对总波数积分以构建理论散射信号。
首先构建传感器阵列来激励和接收Lamb波。每个采集到的信号对应唯一的激励点和接收点,称为激励-接收对。为了生成指示缺陷的图像,需要将感兴趣的成像区域离散为M个成像点。每个成像点都会被视为潜在的缺陷位置。假设一共收集到N个信号,对第m个成像点的反向传播信号可以表示为
$U_{\mathrm{b}, n}(\omega, m)=U_{n}(\omega) \mathrm{e}^{-\mathrm{i} \Delta x\left[\bar{k}_{\mathrm{R}-m, n}(\omega)+\bar{k}_{m-\mathrm{S}, n}(\omega)\right]} 。$
式中:Ub,n(ω,m)为第m个成像点的第n个反向传播信号;Un(ω)为频域内第n个采集信号; ${\stackrel{-}{k}}_{R-m,n}$ (ω)和 ${\stackrel{-}{k}}_{m-S,n}$ (ω)分别为第n个收发对接收点到第m个成像点的Lamb波广义波数和从第m个成像点到第n个收发对的激励位置的总波数。因为传播路径已知,所以可以根据公式(5)和公式(6)进行计算。若与实际缺陷位置重合,反向传播信号的相位与激励信号一致,形成相干叠加;否则,信号因相位错位而抵消。
为提取缺陷散射成分,对反向传播信号进行希尔伯特包络分析
UD(t,m)=H[ub,n(t,m)]。
式中,ub,n(t,m)是Ub,n(ω,m)的时域信号对于第i个传感器对的包络信号。最终缺陷指标Em通过时间对齐与能量累加求得,计算出M个成像点的所有缺陷指标后,就可以生成一幅图像来指示缺陷位置:
Em= ${\left[\sum _{i=1}^{M}{y}_{i}\left({T}_{i,m}\right)\right]}^{2}$ 。
式中,yi(Ti,m)表示在像素点m处对应飞行时间的包络信号幅值。 本文提出的变厚度板缺陷成像方法流程图如图5所示。
2 数值仿真
用通孔(孔的直径为8 mm)来模拟缺陷。铝板沿厚度方向对称,其厚度在y轴方向上从2 mm线性增加到4 mm,板厚沿z轴方向不发生变化。试件内所有单元类型皆设置为C3D8R。设置网格单元长度为1 mm,厚度方向划分为4层。在仿真中,采用Abaqus显式求解器进行动态分析。模拟中的时间增量设置为2.5×10-8 s,总时间长度为2×10-4 s。所有接收到的信号都以10 MHz的新采样频率下采样,测点位置与编号及网格缺陷如图7所示。
选择xy平面中-200<x<200和-200<y<200的区域作为成像区域将其划分为一系列尺寸为1 mm × 1 mm的网格。计算每条传播路径中Lamb波的总波数时,将每条传播路径划分为一系列不同厚度的微小等宽微元,沿路径方向的微元长度设为Δx=1 mm。由于本文方法不涉及非线性,因此互易性在模拟中有效,即传感器对i-j的响应信号与传感器对j-i的响应信号相同。因此,为了降低计算复杂度,考虑到互易性,缺陷成像只使用15个唯一的响应信号。将激励Lamb波更改为S0型激励,并保持其他条件相同,所得成像结果如图9c所示。从图9c可见,该结果与A0型激励所得成像结果(图9a)并无较大差别,由此可见,激励模式对本文成像方法并无影响。
不同成像方式的数值模拟结果如图9所示,成像结果使用最大值进行归一化,并进行对数处理。图9中白色圆圈代表通孔缺陷的实际位置。成像幅值的动态范围选择为-15 dB到0 dB。其动态范围的解释为:以所有成像点中的最大幅值作为0 dB,其余成像点幅值若在-15 dB到0 dB范围之内,则不变;若小于-15 dB,则显示为-15 dB。从图9a中可以看出,虽然存在一定干扰,但本文方法可以正确定位缺陷位置,定位误差约为5.8 mm。将变厚度板分别近似为2 mm、3 mm和4 mm等厚度板,成像结果出现较大误差。由表1中5种方式的成像指标数据可以看出:本文成像方法获得的缺陷位置最明显(像素值最大),同时定位误差最小。
表1 不同方法的成像指标Tab.1 Imaging metrics for the four modalities |
| 方法 | 实际缺陷中心对应 的像素值/dB | 定位误差/mm |
|---|---|---|
| A0模式下本文方法 | 19.45 | 5.83 |
| S0模式下本文方法 | 20.38 | 5.66 |
| 等效为2 mm等厚度板 | 11.69 | 24.33 |
| 等效为3 mm等厚度板 | 15.28 | 10.63 |
| 等效为4 mm等厚度板 | 14.96 | 15.52 |
3 实验验证
试样为铝板,其力学参数与数值模拟中所用铝板力学参数相同。板在xy平面上的尺寸为1 000 mm×1 000 mm,厚度在z方向上从2 mm线性增加到8 mm,在板的上表面共粘接9个压电陶瓷(PZT)作为传感器阵列来激发和接收Lamb波,传感器阵列排布和缺陷位置如图10所示。PZT为圆柱形,尺寸为直径8 mm,厚度0.5 mm,在厚度方向的谐振频率为4 MHz,在径向方向的谐振频率为250 kHz。粘贴PZT可以改变局部声阻抗,当入射波与其相互作用时将产生散射信号,所以可以通过粘贴尺寸相同的PZT来模拟缺陷。
实验系统主要包括:Tektronix AFG31022信号发生器,功率放大器ATA-2042,直径为8 mm、厚度为0.5 mm的压电陶瓷圆片以及Tektronix DPO3014数字示波器等。Tektronix AFG31022信号发生器能使用线性调频信号激励。功率放大器ATA-2042能够提供最大60倍的电压增益,实现激励信号的放大功能。
实验流程如下:将自定义激励信号导入信号发生器后,在功率放大器中实现电压放大,传输到激励阵元,信号在板中传播到接收阵元处,接收信号通过示波器采集,采集256次后,将接收信号的均值保存在计算机中。信号发生器产生的激励是线性调频信号(如图11所示),起始频率为0 Hz,带宽为600 kHz,持续时间为100 μs。信号采集以中心位置(图10中1号位置)作为激励位置,其余的位置作为接收位置。在构建的系统下,总共接收了8个信号,信号由采样频率为10 MHz的数字示波器记录。基于线性时不变系统的假设[33-34]恢复了不同中心频率的对应脉冲响应,发现以80 kHz为中心的脉冲响应以A0型为主。因此,对中心频率为80 kHz的5周期汉宁窗调制脉冲的窄带响应可以从相应的线性调频响应中恢复。
在原点处激励时,对比了2号传感器处接收到的实际信号和预测响应信号(见图12)。可以看出,接收到的信号不仅包含损伤散射信号和直达波,还包含模式转换波和边界反射信号。
4 结论
本文提出了一种基于Lamb波的缓变厚板缺陷定位成像方法。首先介绍了2种Lamb波传播模型:等厚度板和变厚度板的Lamb波传播模型。在变厚度板中Lamb波的传播可以被视为在一系列不同的等厚板中传播的Lamb波。通过将波导结构沿传播路径离散化,形成一系列具有不同厚度的等宽微元,在已知激励信号和传播路径时,就可以实现响应信号的预测及其频散补偿。在此基础上提出的延迟求和成像方法,可以实现对各向同性缓变厚度板的缺陷成像,成像中心定位误差小于10 mm。与现有基于等厚度板的延迟求和成像方法相比,中心定位误差精度得到了明显提升。