非线性表面波的二阶微扰解及特性分析*

曾胜洋 贾璐 张书增‡ 李雄兵 王猛

1)(中南大学交通运输工程学院,长沙 410075)

2)(清华大学航天航空学院,应用力学教育部重点实验室,北京 100084)

3)(清华大学,柔性电子技术研究中心,北京 100084)

为解决非线性声表面波的求解难题,本文从二阶非线性各向同性介质的超弹性本构方程出发,采用位移势函数法,建立二维表面波的非线性势函数控制方程;通过微扰法推导非线性表面波的准线性解和绝对非线性系数,讨论表面波二次谐波解的主要组成部分;并建立模拟非线性表面波传播的有限元模型,位移幅值的仿真结果与理论符合良好,验证了本文非线性表面波理论的准确性.根据微扰解的数值结果,探讨了非线性表面波的传播以及非线性系数的特性,结果表明: 表面波二次谐波由累积项及非累积项组成,前者与表面波纵波分量自相互作用相关,但当初始条件和传播距离相同时,该部分谐波幅值比纯纵波的二次谐波幅值大;此外,纵波和表面波的非线性系数存在正比关系,该比例关系由材料的二阶弹性系数确定.本文探究的非线性表面波的传播特性及其绝对非线性系数的定义表达式,对指导非线性表面波的实际应用具有一定意义.

1 引言

近年来,非线性超声理论和技术得到了广泛的研究,并在工业无损检测领域开展了相关应用[1].材料的非线性系数对高阶弹性常数变化和声-微结构作用十分敏感,因此通过测量非线性超声的高次谐波,可有效检测和评价疲劳损伤、热损伤、粘接分层和微裂纹等缺陷信息[2,3].当前非线性超声检测主要采用纵波[4]、表面波[5]、兰姆波[6]以及导波[7]等波形,其中非线性表面波具有传播距离远、能量集中和单边测量的特性,更适用于实际检测应用,因此相关理论和实验研究备受关注.

在非线性表面波的应用方面,国内外学者开展了一系列工作.例如,通过非线性表面波检测玻璃的表面刮擦损伤,体现了其快速评价表面损伤的潜力[8];相关实验表明非线性表面波系数与钛合金的塑性形变程度呈正相关[9];此外,利用空耦探头测量混凝土中的非线性表面波,可有效区分因碱硅反应导致的不同程度的损伤[10].然而,当前非线性表面波的实验研究普遍测量损伤试块的相对非线性系数,并需要相同实验条件下的初始试样测量结果作基准,通过观测该系数的变化做出相对性评价[11].相关研究表明,通过测量材料的绝对非线性系数,不仅可以对材料弹性常数等特性进行评价,还能够对材料的损伤状态等进行定量或定性分析[3].因此利用绝对表面波非线性系数进行材料特性或损伤评价具有更大的优势,而研究非线性表面波的相关理论并提出对应测量方法是保证其应用的基础.

纵波非线性特性已得到深入的研究和广泛认知[12,13],但非线性表面波非线性的研究仍存在不同的认识,受当前相关理论差异的影响,难以实现表面波绝对非线性系数的计算和进一步的实验测量.针对非线性声波问题主要有两类研究方法,一是数值法,以有限差分法或有限元法为代表[14],目前多用于复杂结构的声信号模拟与定性分析[15];另一类是解析法,通常可以给出理想问题的准确解或近似解.前者先将声波方程离散化,输入相关材料参数,再求解出瞬态声场;虽可获得特定复杂问题的数值解,但难以建立各物理参数之间的具体联系,即探究相关问题的物理本质.相比之下,求非线性表面波问题的解析解,对于非线性表面波特性、计算乃至测量表面波绝对非线性系数具有更重要的意义.

针对利用解析法进行非线性表面波求解的难题,有学者通过能流密度[16]等方法进行推导,但重点研究的是表面波的谐波产生和波形演化等特性,研究中包含了多组表面波非线性系数,且并未给出具体的检测方法.Herrmann等[17]基于声表面波的纵波分量推导出表面波非线性系数,并给出该系数的测量表达式,在工程实践中得到一定应用.但Masurkar和Tse[18]指出仅利用表面波纵波分量推导并不能得到真实的非线性系数,他们认为横波分量对表面波非线性存在贡献且不容忽视,并提出一种基于测量谐波幅值的表面波非线性系数,然而该方法与数值法类似,无法解释谐波随距离增加的累积效应,且不能描述表面波的谐波传播特性.

针对当前研究中的不足,本文在准线性理论成立条件下开展非线性表面波传播特性研究,并解决表面波非线性系数定义等问题.为了能够应用解析法对非线性表面波方程进行求解,引入位移势函数描述非线性表面波的波动控制方程,结合表面波的边界条件和微扰法求得表面波基波和二次谐波解,进而研究表面波非线性系数的定义和测量表达式.为验证理论求解的正确性,引入有限元模型仿真对非线性表面波的传播进行模拟.通过数值计算分析表面波、纵波二次谐波的传播特性并探究几种非线性系数的联系,为后续非线性表面波技术的发展和应用提供理论依据.

2 理论

2.1 二阶非线性表面波传播的控制方程

本文主要研究在半无限大、各向同性的固体表面传播的最典型的声表面波,即瑞利波.如图1 的笛卡尔坐标系所示,非线性表面波沿x1方向传播,并在x3方向迅速衰减,忽略声波传播时衍射的影响,即各物理量不随坐标x2变化.当该表面波传播时满足平面应变条件,质点的运动方程表示如下:

图1 表面波传播坐标系Fig.1.The coordinate for the propagation of the surface wave.

式中Pij,j表示∂Pij/∂xj,ui,tt表示∂2ui/∂t2(该 简写标记方式将贯穿本文);Pij为第一类皮奥拉-基尔霍夫应力(Piola-Kirchhoff stress,PK 应力),ui为xi方向的位移,t为时间,ρ0为材料密度;式中i,j满足爱因斯坦求和约定.针对常见金属材料的超弹性本构特性,采用Landau和Liftshitz[19]提出的应变能表达式来描述:

式中W是弹性体每单位质量比应变能;λ和µ是材料的拉梅常数;A,B和C是材料的三阶弹性常数;e是格林应变张量.方程(1)中的PK 应力与应变能之间的关系为

式中F是变形梯度,Fijδij+ui,j.将方程(2)和方程(3)代入到方程(1),仅考虑到二阶非线性项,可得非线性表面波质点位移的运动方程:

式中α=1,β=3 或α=3,β=1.方程(4)等号右侧表示为非线性项,其中包含A,B和C的项体现了材料本构的非线性;包含λ和µ的项体现了应力应变关系的几何非线性.当此项为0 时,方程(4)表示为线性表面波的质点位移运动方程.表面波在自由边界传播时,即在x30 处,其应力满足:

式中P13和P33是应力分量,二者均包括线性项和非线性项,在方程中分别用上标L和NL 表示,具体表达式为

方程(4)—方程(10)构成了二维非线性表面波质点运动位移场问题的控制方程.

2.2 微扰法求解

引入位移势函数和微扰法推导非线性表面波位移场的准线性解,即二次谐波幅值远小于基波幅值,且谐波幅值处于累积范围中.引入纵波势φ和横波势ϕ,他与位移有如下关系:

根据二阶微扰近似法,将方程的解表示为线性解和二阶非线性解部分,即ϕϕ(1)+ϕ(2),φφ(1)+φ(2),且ϕ(2)≪ϕ(1),φ(2)≪φ(1),结合方程(4)—方程(11),可将描述质点位移的运动方程表示为

方程在x30 满足的边界条件为

同理,将方程(11)代入方程(8)和方程(10)可得应力分量非线性项的具体表达式:

方程(12)、方程(13)和方程(16)描述了线性表面波的位移运动方程和边界条件问题,方程的解析解在相关文献[20]中描述如下:

二次势函数谐波分量是基波解的修正,即方程(14)、方程(15)和方程(17)的解.求解之前,首先需获得仅含势函数φ(2)和ϕ(2)的解耦方程,可将线性表面波解方程(22)代入方程(14)和方程(15),接着方程(14)对x1求导,方程(15)对x3求导,二者相加得到方程(25);方程(14)对x3求导,方程(15)对x1求导,二者相减得到方程(26):

其中R1-R4表示如下:

其中,κcl/ct为纵横波声速比,βl-3-(2A+6B+2C)/(λ+2µ)为纵波非线性系数,βt-(λ+2µ+A/2+B)/µ为横波非线性系数[21].方程(25)和方程(26)为非齐次偏微分方程,可使用叠加原理[22]得到特解如下:

至此,方程(22)、方程(31)和方程(32)组成了非线性表面波传播的准线性解,可见仅纵波位移势φ(1)自谐振产生的二次分量(即平方项)具有累积效果,这一结论与非线性表面波超声中可检测的二次谐波来自于纵波分量的自相互作用[23]的认知相符合.

需要注意到前面求解线性问题时曾利用线性边界条件方程(16)确定表面波声速;而表面波在各向同性材料是非频散波[16],且微扰法暗含瑞利波波速不变的假设,故不再对非线性边界条件方程(17)进行研究.

2.3 基于微扰解的表面波非线性系数

忽略方程(4)含下标β的项并令α1 即可得到一维非线性纵波的质点位移运动方程如下:

因此,纵波非线性系数测量表达式为

通过对比纵波非线性系数的定义式及测量表达式,开展对表面波非线性系数的研究.在非线性表面波实验中,测量信号一般是x30 表面的面外位移分量.结合方程(18)、方程(19)和方程(22)—(27),位移u3可采用纵波势函数幅值Aφ表示:

式中,右边第一部分A1代表基波,第二部分A2代表累积效应的二次谐波,第三部分代表无累积效应的二次谐波.其中Dconst表示如下:

在实际应用中,换能器的接收信号中二次谐波相对较弱,且二次以上的高次谐波更难获得,所以求得的摄动解仍具有参考意义.此外,需要说明的是,本文仅考虑在一定的传播距离内进行测量,即可检测随距离增加具有累积效应的二次谐波成分.因此可以参考方程(35)给出非线性表面波的测量表达式:

结合方程(27)、方程(36)和方程(38),可以得到非线性表面波实验中测量的非线性系数(定义式)为

实际上,消去表达式中的波数kr,以及η和ξ中包含的角频率ω和密度ρ0,可以看到方程(39)与纵波非线性系数βl相似,仅由材料的二、三阶弹性常数确定.

3 仿真分析

本节基于COMSOL MultiphysicsR有限元软件验证第2 节中得到的非线性表面波解析解,主要应用该软件中“Murnaghan material”材料模型探究材料非线性与声学非线性问题.该模型需要指定材料二、三阶弹性常数数值以确定材料应力应变关系,与前述理论推导得到非线性表面波仅由材料的二、三阶弹性常数确定的结论相符,且已有研究证实了利用该软件验证非线性超声理论问题的可行性[24].

3.1 表面波有限元模型

忽略表面波传播过程中的衰减和衍射影响,即建立非衰减材料中的二维有限元模型,如图2 所示,选用纯铝作为仿真对象,具体材料参数如表1 所列.为了增强非线性声波的仿真效果,提高信噪比,仿真中将材料三阶弹性常数设置为表1 中材料三阶弹性常数的10倍[25](即βl173.699,βr292.197).根据方程(24)计算得到线性表面波声速为2892.22 m/s,采用位移幅值为10 nm,频率为2 MHz 的8 个周期正弦波信号作为点声源.模型宽度为40 mm,点源设置在左端点10 mm 处,避免将边界网格奇点作为激发源,同时降低反向传播的声波的干扰,且能保证足够的正向传播距离;模型厚度大于表面波传播距离的2 倍,以消除反射的纵/横波的影响.求解区域材料属性设置为超弹性,采用直接求解器中速度最快的Pardiso 求解器求解.最后,在表面点源后以1 mm 间隔放置探针记录x3方向位移,采用FFT 算法将记录的时域位移信号转换为频域信号进行分析.

表1 纯铝和三种铝合金的密度和弹性常数Table 1.Densities and elastic coefficients of pure aluminum and three aluminum alloys.

图2 非线性表面波二维有限元模型示意图Fig.2.Schematic diagram of two-dimensional (2D)finite element model for nonlinear surface wave.

为保证计算精度和可靠性,将近表面黄色区域(厚度大于4 倍波长)最大网格尺寸dmax设置为最大解析频率4 MHz 时表面波波长的1/20,即3.615×10—2mm,最小网格尺寸dmin设置为dmax/2,其余区域划分为预定义极细自由三角形网格.时间步长 Δt由关联最小网格尺寸和时间步长的CFL参数确定:

CFL 一般取值0.2,则 Δt取值为1.25×10—9s.

3.2 有限元模型检验

为了验证本文激励信号波形、材料参数和网格时间步长设置的有效性,开展非线性纵波传播数值模拟的预验证实验.采用线源和周期性边界条件,使纵波能以平面波形式传播.相关传播模型如图3 所示,试块厚度设为1 个纵波波长以节省计算量.与3.1 节相比,仿真材料参数数值设置不变.时间步长不变,网格参数按照纵波波长做相应调整,探针布置在试块中间记录x1方向位移.

图3 非线性纵波二维有限元模型示意图Fig.3.Schematic diagram of 2D finite element model for the nonlinear longitudinal wave.

探针波形信号以及纵波中基波和二次谐波随传播距离的变化情况如图4 所示.基波幅值在传播过程保持不变,拟合幅值是10 nm,与所模拟纵波为平面波的设置相符;二次谐波幅值随距离线性增加,需要说明的是,非0 截距来自数值离散误差.

图4 纵波仿真结果图 (a)传播距离15 mm 处探针的信号图;(b)基波和二次谐波幅值随距离变化图,蓝色点是仿真结果,黄色实线是线性拟合结果Fig.4.Simulation results for the nonlinear longitudinal wave: (a)Typical signal at 15 mm propagation distance;(b)plots of fundamental wave and second harmonic amplitude versus propagation distance,where blue points denote the simulation results,and the yellow solid line denotes the fitting line.

根据方程(34),二次波幅值-传播距离的斜率理论系数应为利用仿真结果拟合值以及βl和kl的理论值可反算出基波幅值,再将其与仿真基波拟合幅值进行对比验证.预测基波幅值为9.98 nm,而仿真结果为10 nm.纵波超声非线性具有较为成熟的理论,仿真结果与理论一致,由此验证了该仿真模型在非线性超声波传播特性模拟中的有效性.

3.3 非线性表面波仿真

本小节对非线性表面波进行仿真研究.仿真中得到的非线性表面波的典型波形信号,以及不同位置处基波和二次波面外分量幅值随传播距离的变化如图5 所示,二次波幅值随传播距离增加近似线性增长.

图5 表面波仿真结果 (a)传播距离15 mm 处探针的信号图;(b)基波和二次谐波幅值随距离变化图,蓝色点是仿真结果,黄色实线是拟合结果Fig.5.Simulation results for the nonlinear surface wave: (a)Typical signal at 15 mm propagation distance;(b)plots of fundamental wave and second harmonic amplitude versus propagation distance,where blue points denote the simulation results,and the yellow solid line denotes the fitting line.

首先,根据方程(41)表示的波形重叠法确定不同距离点信号的传播时间,计算得到非线性表面波声速[29]:

其中,f1和f2是距离为d的两个时域信号序列.计算得到非线性声表面波速度为2892 m/s,验证了2.2 节关于表面波传播声速不变的假设.同时,通过声速可以确定图5(a)所示波包为表面波.假设波包中的二次谐波仅源于纵波分量的纯纵波非线性特性,根据3.2 节该条件下计算所得的二次谐波-距离斜率理论系数为0.121,与图5(b)所示有明显差异,表明该二次谐波应来自非线性表面波传播特性,且该过程与纯纵波非线性效应区别较大.最后,根据方程(38)和方程(39),二次波幅值-传播距离斜率理论系数为据此反算出基波幅值为1.452 nm,与基波拟合幅值相比符合良好.综上所述,有限元仿真结果验证了本文提出的非线性表面波理论的准确性.

4 非线性表面波特性的数值计算与分析

4.1 谐波幅值的数值结果

首先利用微扰解的数值结果分析非线性表面波谐波传播特性,并与非线性纵波谐波对比.根据方程(36)和方程(34),采用表1 中纯铝的材料参数,分别计算不同初始条件和传播距离的表面波面外分量,以及同等条件下纵波传播时二次谐波幅值,结果如表2 所列.

表2 初始条件和传播距离对面外分量和纵波二次谐波幅值的影响Table 2.The effects of initial conditions and propagation distance on the amplitudes of out-of-plane component and the second harmonic longitudinal wave.

由表2 可知,累积二次表面波幅值A2和纵波二次幅值均随基波幅值、频率或传播距离的变化呈正相关,具有相同的传播特性.然而,在相同条件下,由于表面波波数较大,A2幅值会比高一个数量级;在实际中,低衰减和能量集中等特性将使二次表面波可测量性更强.因此,非线性表面波技术相比纵波具有更高的信噪比.非积累二次谐波幅值幅值比A2小两个数量级,且可能包含边界耦合非线性[30]和微扰法误差,因此实际测量中可以忽略.

4.2 非线性系数的数值结果

本小节将根据表1 中的材料参数,计算纵波和表面波非线性系数数值结果并进行对比分析.在此之前,首先列出当前学界应用较多的几种表面波非线性系数定义,并纳入比较范围.Herrmann等[17]基于非线性表面波纵波分量的近似推导得出:

结合方程(38),其定义的表面波非线性系数可表示为

Shull等[31]基于Zabolotskaya[16]的理论体系,以质点速度谐波幅值vn(n1,2,···)作为研究变量,参考流体情况,定义表面波非线性系数β11为

参数S11和ζ的详细表达式见文献[31],vn和有如下关系[16]:

Masurkar和Tse[18]对二阶非线性应力应变关系做进一步近似,提出一种基于谐波实际测量幅值的归一化非线性系数δR.该系数实际是将信号按基波幅值归一化,表示如下:

式中F(kr,η,ξ)是特征函数,具体表达式见文献[18].

将表1 中4 种金属的材料参数代入方程(39)、方程(43)、方程(46)和方程(47)中计算上述相关的非线性系数,结果如表3 所列.

表3 表1 中4 种材料的非线性系数Table 3.Nonlinear parameters of four materials listed in Table 1.

由表3 可知,βr数值比βl大,而的数值最小.这是由于推导过程主要是基于的质点速度慢变幅值,且其推导中存在近似,不适用于实际测量对象是面外位移分量的情况.δR和βr数值相近,但对比方程(47)和方程(39),不考虑衰减和衍射效应时,δR值会随初始条件和传播距离的变化而变化;反之,当衰减或衍射影响成为主导因素时,其数值可能会出现单个或多个局部峰值,该现象还有待深入理论研究.此外,由方程(39)可知,表面波非线性系数与纵波非线性系数存在确定数值关系,与材料的二阶弹性常数相关,随材料变化而改变.综上,结果表明,非线性表面波传播特性与纵波相似,但在初始条件相同时表面波谐波幅值比纵波谐波幅值高,因此在实际测量中信噪比更高.

5 总结

本文首先采用微扰法研究了准线性条件下材料中的非线性表面波传播特性,引入势函数描述非线性表面波的波动控制方程,推导了基波和二次谐波解,并得到了表面波非线性系数定义式.结果显示,由纵波位移势自谐振产生的累积项幅值比纵波的二次谐波幅值大,而且具有能量集中、衰减低等特性,使非线性表面波检测技术在实际应用中比非线性纵波具有更高的信噪比;由横波分量引起的谐波项幅值没有累积效应,因此对表面波二次谐波的影响可以忽略.二维表面波有限元模型位移幅值仿真结果与理论值符合良好,证实了本文理论推导过程的有效性以及所得谐波微扰解的准确性.最后,根据微扰解的数值结果探究了当前几种表面波相关的非线性系数之间的联系,对比计算结果可知表面波非线性系数和纵波非线性系数存在正比关系,并且该系数只与材料二次弹性常数有关.目前仅限于理论推导和数值模拟分析,后续将继续研究实际测量相对于当前理想情况的修正理论,从而研发一套可靠的测量设备及方法,实现利用表面波绝对非线性系数评价材料状态.