一维非线性弹性波分岔‒谐波耦合机理与亚谐波激发阈值研究

Bifurcation‒Harmonic Coupling Mechanism and Subharmonic Threshold in One-dimensional Nonlinear Elastic Waves

摘要

本文基于非线性系统分岔演化规律阐释非线性超声检测谐波的产生机理，探讨了高次谐波、亚谐波与非线性系统分岔之间的关联。以一维黏弹性杆中的非线性弹性波为研究对象，采用有限元离散化方法，求解非线性弹性波随外激励幅值变化时的分岔行为，通过构建每个分岔阶段响应信号的相轨图和频谱图，获得弹性波分岔行为与非线性谐波之间的关系。结果表明，系统发生叉式分岔时，非线性弹性波失去对称性，产生偶数次谐波，倍周期分岔与亚谐波存在对应关系。在此基础上，以1/2亚谐波为例，分析其激发阈值条件，发现该阈值在由激发参数“幅值‒频率”构成的平面内呈“V”字形，在亚谐波二倍频率下激发时，产生1/2亚谐波所需的激励幅值最小。对含塑性损伤聚乙烯板进行的实验研究表明，1/2亚谐波的激励阈值条件与仿真分析结果一致，验证了非线性超声检测中谐波分量与系统分岔之间的紧密关联，为阐明谐波的产生条件提供了新思路，也为聚乙烯材料的早期损伤识别提供了参考。

Abstract

Objective Nonlinear harmonics are sensitive indicators of early-stage material degradation and microdamage. However, their generation strongly depends on the excitation frequency, amplitude, and nonlinear dynamic state of the material, leading to empirical selection of excitation parameters and ambiguous physical interpretations of harmonic signals. This study aims to reveal the inherent mechanism of harmonic generation from the perspective of nonlinear system bifurcation evolution in nonlinear ultrasonic testing, investigate the coupling mechanism among high-order harmonics, subharmonics, and system bifurcations, characterize the excitation threshold for the generation of the 1/2 subharmonic, and provide a theoretical basis for optimizing parameters for nonlinear ultrasonic detection. Methods Taking nonlinear elastic waves propagating in a one-dimensional viscoelastic rod as the research object, this study constructs a wave equation incorporating quadratic and cubic nonlinear terms and linear viscoelastic damping. A finite element discretization method is adopted to establish a 10-element numerical model with a fixed left end and sinusoidal excitation applied at the right end. By sweeping the external excitation amplitude, bifurcation diagrams of elastic waves over the full excitation amplitude range are obtained. Time-domain responses, phase portraits, and frequency spectra at each bifurcation stage are analyzed to quantitatively characterize the coupling between elastic wave bifurcation behavior and nonlinear harmonics. Two typical bifurcation modes, period-doubling bifurcation and pitchfork bifurcation, are discussed separately. On this basis, the amplitude‒frequency threshold characteristics of the 1/2 subharmonic are systematically investigated through multi-frequency scanning simulations. Finally, experimental tests are conducted on polyethylene plates with plastic damage to verify the simulated threshold characteristics of the subharmonic. Results and Discussions Numerical results show that pitchfork bifurcation breaks the central symmetry of nonlinear elastic-wave responses, thereby introducing even-order harmonics into the frequency spectra, whereas period-doubling bifurcation has a direct correspondence with the generation of subharmonics. After each period-doubling event, fractional subharmonics (1/2 and 1/4) emerge in the spectrum, while chaotic motion broadens the subharmonic peaks and reduces their identifiability. The amplitude of the 1/2 subharmonic exhibits stable sensitivity over a broader bifurcation range than that of the 1/4 subharmonic, making it more suitable for characterizing variations in material nonlinearity. The effective excitation zone of the 1/2 subharmonic exhibits a typical V-shaped distribution in the amplitude‒frequency phase plane, and the minimum excitation amplitude required to generate the 1/2 subharmonic occurs near twice the natural frequency of the structure. Experimental data obtained from plastically damaged polyethylene specimens are highly consistent with the simulated threshold characteristics of the 1/2 subharmonic, confirming the strong correlation between harmonic components and system bifurcations in nonlinear ultrasonic testing. Conclusions This study elucidates the dynamic correlation mechanism between nonlinear elastic wave bifurcations and harmonics and quantifies the V-shaped amplitude‒frequency excitation threshold of the 1/2 subharmonic. The results demonstrate that even-order harmonics emerge only after the symmetry breaking associated with pitchfork bifurcation, while subharmonics have stricter excitation threshold requirements than high-order harmonics. The good agreement between the experimental and numerical results validates the proposed bifurcation‒harmonic relationship. This work provides a novel analytical perspective for interpreting the conditions governing harmonic generation and offers a theoretical basis for optimizing excitation parameters and identifying early-stage microdamage in nonlinear ultrasonic testing.

Graphical abstract

关键词

非线性黏弹性杆 / 弹性波 / 非线性超声 / 分岔 / 谐波

Key words

nonlinear viscoelastic rod / elastic waves / nonlinear ultrasonic / bifurcation / harmonic

本刊网刊

20世纪50年代，研究人员发现材料中的微观缺陷与弹性波在材料中传播时产生的非线性效应密切相关。基于该非线性效应发展起来的非线性超声无损检测方法得到了迅速发展^[1‒3]。作为一种常用指标，非线性谐波被广泛应用于材料早期力学性能退化和微观损伤的检测中。一般认为，在线性材料中激发特定频率的超声波时，接收信号与激发信号具有相同的频率；但对于非线性材料，接收信号会产生与激发信号不同的额外频率。当额外频率是激发频率的整数倍时，称为超谐波；而当额外频率是激发频率的分数倍时，称为亚谐波。当材料含有早期微观损伤（如早期疲劳、残余应力等），可通过固体材料的二阶弹性常数、三阶弹性常数，甚至高阶弹性常数（塑性损伤则引入塑性常数）等非线性项修正线性化的本构关系来描述^[1]，利用非线性超声谐波特性可实现材料早期的微观损伤识别，该方法已被用于材料疲劳^[4‒6]、高温蠕变^[7]、塑性变形^[8]、残余应力^[9]和应变腐蚀^[10]等引起的材料损伤。例如，Frouin等^[11]利用超声波二次谐波发生技术研究了钛合金的疲劳损伤。研究发现，在90%屈服应力以下，当试样承受40%的寿命周期载荷时，由二次谐波计算的非线性参数的振幅将增加180%。这表明二次谐波对材料的疲劳特性非常敏感。Li等^[12]利用二次谐波对沥青混凝土的热损伤和修复进行评估，获得了优化的热修复效率；通过能量法和非线性超声法研究了沥青混合料的梯度加热愈合性能和最佳愈合效率。Li等^[13]利用非线性兰姆波分析了AISI 316L不锈钢的微观组织演变，得出了马氏体转变、晶粒细化与二次非线性系数之间的关系。Mao等^[14]利用二次谐波评估了检测热处理过程中碳的渗透深度。研究表明，非线性谐波具有广阔的应用前景。

但在很长一段时间内，非线性超声检测的研究多以实验研究为主，目的是找到非线性谐波与损伤之间的关系，对机理研究的缺乏严重阻碍了非线性检测方法的应用。为此，许多学者针对非线性谐波的产生开展了深入的研究。例如，Shan等^[15]认为，二次谐波通常与材料的三阶Bragg弹性常数有关，而三次谐波则与材料的三阶和四阶Bragg弹性常数有关，同时给出了两种三次谐波的产生机制：一方面，由于材料具有四阶弹性常数，三次谐波可以直接由基波产生；另一方面，由于三阶弹性常数先产生二次谐波，随后二次谐波与基波发生混频进而产生三次谐波。Deng等^[16]指出，只有当基波和二次谐波的相速度相互匹配时，能量才会从基波转化为二次谐波。王芝^[17]和Wang^[18]等认为，在实际测量中，放大器、耦合介质和换能器的非线性也可能产生高次谐波，但不会产生亚谐波，其将裂纹视为一个单自由度系统，建立了亚谐波产生的理论模型，得到了亚谐波产生的激发条件。Zhang等^[19]研究了由疲劳裂纹非线性和边界非线性产生的亚谐波阈值区域，发现可激发亚谐波的临界条件在“频率‒幅值”平面上均分布在一个“V”形区域。尽管上述研究为非线性弹性波的产生提供了一定的理论解释，但仍然缺乏对谐波产生的定量分析，这影响了非线性超声波在无损检测中的推广应用。

早在1980年，朱位秋^[20]推导了黏弹性杆中的弹性波，指出在一定条件下，黏弹性杆中的弹性波可化简为杜芬方程。1985年，Parlitz等^[21]通过讨论杜芬方程的共振曲线、分岔图与相轨图，指出分岔存在特定结构的周期性重现现象与系统的非线性共振密切相关。2014年，唐元璋等^[22]应用半数值近似解析方法分析了Duffing 振子系统的倍周期分岔规律，用增量谐波平衡法获得高精度高次谐波半数值近似解析解。Ohara等^[23]将超声亚谐波的特殊激励条件归因于波形失真，导致了倍周期分岔，并产生了激励频率减半的波形。Eftekhari等^[24]指出，对于具有内在非线性特征的系统，采用非线性损伤检测更为适用，并进一步评估了基于分岔边界分析进行损伤检测的可行性。本文聚焦一维黏弹性杆中的非线性弹性波，将非线性超声视为一种非线性动力学系统，尝试从非线性系统分岔理论探讨谐波的产生机理，通过仿真和实验研究验证了理论的有效性，为非线性超声应用提供理论依据。

1 有限元仿真

考虑一维黏弹性杆，其材料属性具有二次、三次非线性，本构关系可表示为：

\sigma =E\epsilon +{\alpha }_{\mathrm{2}}{\epsilon }^{\mathrm{2}}+{\alpha }_{\mathrm{3}}{\epsilon }^{\mathrm{3}}+E^{\text{'}}\frac{\partial \epsilon }{\partial t}（1）

式中：\sigma为应力；\epsilon为应变，有\epsilon =\frac{\partial u}{\partial x}（u为x方向的位移）；E为弹性模量；{\alpha }_{\mathrm{2}}为二次非线性系数；{\alpha }_{\mathrm{3}}为三次非线性系数；E^{\text{'}}为法向黏性系数；E^{\text{'}}\frac{\partial \epsilon }{\partial t}表示线性黏弹性阻尼，其中t为时间。设杆件的几何尺寸和材料参数如表1所示^[25]。

将式（1）代入波动方程中可得：

\begin{array}{l}\rho \frac{{\partial }^{\mathrm{2}}u}{\partial t^{\mathrm{2}}}=E\frac{{\partial }^{\mathrm{2}}u}{\partial x^{\mathrm{2}}}+\mathrm{2}{\alpha }_{\mathrm{2}}\frac{\partial u}{\partial x}·\frac{{\partial }^{\mathrm{2}}u}{\partial x^{\mathrm{2}}}+\\ \mathrm{3}{\alpha }_{\mathrm{3}}{\left(\frac{\partial u}{\partial x}\right)}^{\mathrm{2}}\frac{{\partial }^{\mathrm{2}}u}{\partial x^{\mathrm{2}}}+E^{\text{'}}\frac{{\partial }^{\mathrm{3}}u}{\partial x^{\mathrm{2}}\partial t}\end{array}（2）

为了求解该波动方程，将杆件离散为10个单元，如图1所示，设杆件左端为固定约束，右端受简谐载荷F,F=F_{\mathrm{0}}\mathrm{s}\mathrm{i}\mathrm{n}(\omega t)，F_{\mathrm{0}}表示载荷幅值，\omega表示圆频率。有限元离散化后得出动力学方程为：

\mathit{M}\ddot{\mathit{U}}+\mathit{C}\dot{\mathit{U}}+\mathit{K}\mathit{U}+\mathit{P}\mathit{U}+\mathit{Q}\mathit{U}=\mathit{F}(t)（3）

式中： U 为位移u的列向量；\dot{\mathit{U}}=\frac{\partial \mathit{U}}{\partial t}，为 U 对t的一阶导数；\ddot{\mathit{U}}=\frac{{\partial }^{\mathrm{2}}\mathit{U}}{\partial t^{\mathrm{2}}}，为 U 对t的二阶导数；\mathit{M}为质量矩阵；\mathit{C}为阻尼矩阵；\mathit{K}为整体刚度矩阵；\mathit{P}为非线性二次平方项矩阵；\mathit{Q}为非线性三次平方项矩阵。为了便于求解，以密度、弹性模量、长度为基本物理量，分别记为{\rho }_r、E_r和L_r，对其他物理量做如下无量纲化处理：

\begin{array}{l}{\rho }^{\mathrm{*}}=\frac{\rho }{{\rho }_r},L^{\mathrm{*}}=\frac{L}{L_r},E^{\mathrm{*}}=\frac{E}{E_r},u^{\mathrm{*}}=\frac{u}{L_r},A^{\mathrm{*}}=\frac{A}{L_r{}^{\mathrm{2}}},\\ {\alpha }_n{}^{\mathrm{*}}=\frac{{\alpha }_n}{E_r},F^{\mathrm{*}}=\frac{F}{E_r{L}_r{}^{\mathrm{2}}},t^{\mathrm{*}}=\frac{t}{L_r}\sqrt[]{\frac{E_r}{{\rho }_r}}\mathrm{。}\end{array}

其中，n=2，3。为了便于描述，下文中所涉及到的物理量均为无量纲化后的物理量，因此所有量的上标*均不再标注。经上述无量纲处理之后，各参量为：

\rho =\mathrm{1},E=\mathrm{0.33},{\alpha }_{\mathrm{2}}=-\mathrm{1},{\alpha }_{\mathrm{3}}=\mathrm{0.67},L=\mathrm{1},\omega =\mathrm{1}（4）

图2为系统随驱动力幅值F的分岔图。将无量纲化处理后的参数代入式（3），假定F_{\mathrm{0}}\in[0.000 01, 3.500 00]，且以步长0.000 01变化，初值均为0（所有节点位移和速度均为0），其全局分岔图如图2（a）所示，从图2（a）中可以看出，随着驱动幅值的变化，系统经历了周期、分岔、混沌等多个状态，表明系统具有丰富的动力学特性。为了更清晰地观察该系统的非线性特征，选择图2（a）中区域Ⅰ和区域Ⅱ进行放大，驱动力幅值范围分别为(0.029 0, 0.031 1)和(2.400 0, 3.500 0)，其结果分别如图2（b）和图2（c）所示。依据图2（b）和图2（c）中的分岔点，将驱动力F分为6个区间，其对应的系统状态如表2所示，表示随驱动力变化时系统表现出的多样性。

2 分析与讨论

2.1 倍周期分岔

图3为系统进入倍周期分岔后的响应。图2（b）所示为系统从倍周期分岔到混沌的路径，当0.029 01<F<0.029 12时，不失一般性，选择F = 0.029 06，其响应信号、相位轨迹和频谱分别如图3（a）～（c）所示。响应信号上有一个极值，如图3（a）所示。相位轨迹上显示一个环路，如图3（b）所示，这表明系统处于周期一状态。如图3（c）所示，图中纵坐标S表示频域幅值，可见，频谱图中只有整数倍的高次谐波。当0.029 13< F <0.029 60时，选择F = 0.029 50，其响应信号、相位轨迹和频谱分别如图3（d）～（f）所示。如图3（d）所示，响应信号上有两个不同的极值。相位轨迹上显示两个环路，如图3（e）所示，这表明系统处于周期二状态。可见，非线性系统的周期加倍，则其频率f减半。因此，1/2次谐波的频率及其整数倍出现在频谱中，如图3（f）所示。当0.029 61<F<0.030 22时，系统将再次经历倍周期分岔，处于周期四状态。不失一般性，选择F = 0.030 00，其响应信号、相位轨迹和频谱分别如图3（g）～（i）所示。响应信号上有4个不同的极值，如图3（g）所示。相位轨迹有4个环路，如图3（h）所示。1/4次谐波的频率及其整数倍出现在频谱中，如图3（i）所示。当0.031 02<F<0.031 14时，选择F = 0.031 10，响应信号、相位轨迹和频谱分别如图3（j）～（l）所示。如图3（j）所示，响应的周期性消失。相位轨迹包含许多混沌回路，如图3（k）所示，这些回路堆叠在一起，表明系统在这种情况下已经进入混沌状态。在混沌状态下，亚谐波的频率变得不那么明显，如图3（l）所示，但仍然可以清楚地观察到激励频率及其倍增频率。

2.2 叉式分岔

图4为叉式分岔时系统的响应。如图2（c）所示，F = 2.966 00为叉式分岔点，当2.966 00<F<3.500 00时，不妨令F=3.200 00，并将其代入式（3）进行求解，其时域、相轨和频域结果分别如图4（a）～（c）所示。从图4（a）～（c）可以看出，信号响应幅值关于0点对称，在每个周期结束时，响应仅有一个极大值。相轨迹是一个椭圆形轨道，轨道沿原点中心对称，轨迹中只包含一个环路，说明系统处于周期一状态。对其时域信号进行频谱分析可知，其频谱结果上仅有奇数次谐波，而无偶数次谐波。这种只有奇次谐波的函数被称为奇谐波函数^[26]。只有奇谐波的原因是响应的运动轨迹关于原点中心对称，也就是信号延时半个周期之后与原信号关于原点对称：

z(t)=-z\left(t+\frac{T}{\mathrm{2}}\right)（5）

式中，T为响应信号的周期，z(t)表示响应信号。当响应信号根据z(t)进行傅里叶级数展开时，可以写成：

z(t)=Z_{\mathrm{0}}+{\displaystyle \sum _{k=\mathrm{1}}^{\infty }}{Z}_{k}sin(k\Omega t+{\phi }_k)（6）

式中，Z_{\mathrm{0}}表示响应信号的直流分量，Z_k表示第k阶谐波的幅值，\Omega表示基频圆频率，{\phi }_k表示第k阶谐波的相位。如果延迟z(t)半个周期并取负数，可以得到：

\begin{array}{l}-z\left(t+\frac{T}{\mathrm{2}}\right)=-Z_{\mathrm{0}}-{\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}{Z}_k\mathrm{s}\mathrm{i}\mathrm{n}\left[k\Omega \left(t+\frac{T}{\mathrm{2}}\right)+{\phi }_k\right]=\\ -Z_{\mathrm{0}}-{\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}{Z}_k\mathrm{s}\mathrm{i}\mathrm{n}(k\Omega t+{\phi }_k+k\mathrm{\pi })=\\ -Z_{\mathrm{0}}+{\displaystyle \sum _{k=\mathrm{1}}^{\mathrm{\infty }}}{(-\mathrm{1})}^{k+\mathrm{1}}{Z}_k\mathrm{s}\mathrm{i}\mathrm{n}(k\Omega t+{\phi }_k)\end{array}（7）

显然，为了满足式（5），必须满足Z_{\mathrm{0}}=\mathrm{0}，以及当k是奇数时，Z_k\ne \mathrm{0}；当k是偶数时，Z_k=\mathrm{0}。因此，在傅里叶级数展开中，只保留了奇次项。

当2.696 00<F<2.966 00时，系统已经发生了叉式分岔，不妨设F=2.800 00，并将其代入式（3）进行求解，其时域、相轨和频域结果分别如图4（d）～（f）所示。从时域信号上来看，在每个周期结束时，信号仅有一个极大值，说明系统仍处于周期一的状态，但半波对称特性消失。相应地，相轨迹不再关于原点中心对称，且初始条件取值不同，会出现两种不同的稳定解，即吸引子1和吸引子2，此时的系统敏感依赖于初值的选取，体现了非线性系统的初始条件敏感性。频谱图中出现了包含偶数次谐波在内的完整的整数次谐波，如图4（f）所示，这也与响应半波对称性丧失相对应，使得傅里叶级数展开后出现偶次谐波项。这说明随着驱动力幅值的减小，响应的原点对称性消失，方程从周期一（单吸引子）状态转变为周期一（双吸引子）状态时，结构响应中除了奇次谐波外，还将出现偶次谐波，结构响应的稳定性下降。由此可见，F = 2.966 00是该系统的一个叉式分岔点，在该点附近，若F>2.966 00，系统的响应具有沿原点中心对称的特性，响应信号中将只有奇次谐波，没有偶次谐波，方程的解相对稳定；若F<2.996 00，系统响应将失去原点对称性，响应信号中出现完整的整数次谐波，说明系统已经进入到周期一（双吸引子）状态。一般认为，材料发生非线性损伤会出现二次谐波，从而可以采用二次谐波构造非线性参数^[27]，用以刻画材料损伤积累程度。但从系统分岔结果来看，在叉式分岔前，二次谐波并不会出现。

2.3 亚谐波的激励条件

由以上分析可知，有限元法显示了系统的分岔行为与非线性谐波紧密相关，只有在特定条件下才能激发特定的谐波。亚谐波是表征系统非线性特性的重要分量，其中幅值常用来表征材料的损伤程度。由前面的研究可知，系统只有在发生倍周期分岔后才会产生亚谐波，如周期二时出现1/2亚谐波、周期四状态出现1/4亚谐波，当系统进入混沌态后亚谐波频谱变宽。亚谐波的幅值与系统的分岔特性有紧密联系，如图5所示。当系统处于周期一状态时，1/2亚谐波和1/4亚谐波幅值几乎为0，当系统进入周期二状态时，1/2亚谐波和1/4亚谐波突然增大，1/2亚谐波近似线性增加，但1/4亚谐波表现为近似水平线。当系统进入周期四状态时，1/2亚谐波近似线性减小，随着系统不断发生倍周期分岔，其幅值也会发生相应的突变，当系统进入混沌态后，1/2亚谐波幅值发生波动，系统失去稳定性。同理，当系统进入周期四状态时，1/4亚谐波先快速增加，后变化率变缓，系统再次发生倍周期分岔时，其幅值也会发生相应的突变，当系统进入混沌态后，1/4亚谐波幅值发生波动，系统失去稳定性。

从图5（b）可以看出，两个亚谐波幅值会随着驱动力幅值变化而变化，表明其在一定范围内可用于表征系统参数的变化情况，1/4亚谐波幅值仅在周期四状态内对驱动力幅值敏感，而1/2亚谐波幅值在周期二和周期四状态内均对驱动力幅值敏感，因此1/2亚谐波更适合用于表征系统参数的变化情况，更适合用于材料的微观损伤检测。

从图5（b）还可以看出，亚谐波的出现具有阈值现象，当驱动力幅值控制在一定范围内时，才能激发出亚谐波。由前述分析可知，系统产生亚谐波时说明系统发生了倍周期分岔，因此这里以分岔图的方式分析亚谐波的激发阈值，分别绘制系统关于驱动力幅值的分岔图，通过研究系统倍周期分岔行为确定亚谐波的激励条件。图6所示为系统频率分别取为\omega =\mathrm{0.5}、0.8、1.2、1.5时的倍周期分岔结果，图中箭头所指为系统进入周期二状态时的最小激励幅值，可以看出，随着频率增加，系统进入周期二状态所需要的最小激励幅值呈现出先减小后增大的趋势。

为了获得完整的1/2亚谐波的激励频率‒激励幅值的阈值区间，设置激励频率\omega \in [\mathrm{0.5,1.5}]范围内，以0.1为步长增量，求解共计11个频率下的分岔结果，记录对应频率下出现1/2亚谐波的最小激励幅值，并将其绘制在图7中。从图7可以看到，1/2亚谐波的共振区域在\omega =\mathrm{0.8}处有一个极小值，在\omega =\mathrm{1.6}处有一个极小值。然后，对图1所示的杆结构进行扫频，结果表明\omega =\mathrm{0.8}是结构的一阶固有频率。由此可见，对于1/2亚谐波而言，在结构的一阶固有频率与二阶固有频率附近存在激励幅值的极小值，存在1/2亚谐波的区域在整个激励频率‒幅值区间内呈现“V”字形。

3 实验研究

3.1 试样制备

以聚乙烯为非线性黏弹性材料的代表，通过实验研究，验证非线性黏弹性材料中的谐波特性，探寻亚谐波的激励频率与激励幅值范围。试样为具有黏弹性特性的聚乙烯板材，完整试样长为200 mm，宽为20 mm，厚度为4 mm，试样尺寸如图8所示。

需要指出的是，尽管聚乙烯材料具备一定的非线性特征，但其是一种弱非线性材料，为了便于研究亚谐波的阈值特性，本文通过对聚乙烯材料制造塑性损伤，以增强其非线性特征。具体来说，采用WDW‒300 kN微机控制电子万能试验机对聚乙烯试样进行单轴拉伸，拉伸过程中损伤逐渐累积然后试样断裂。研究发现，当应变达到0.25以上时，材料完全进入塑性阶段，且试样颜色发生明显变化，出现局部发白现象，这是常见的聚合物的银纹现象，表明材料出现塑性损伤^[28]。因此，将损伤试样的拉伸应变定在0.3。然后，对聚乙烯试样进行超声检测试验。在试样中间位置，沿厚度方向上下各粘贴一个压电片，采用Tektronix AFG3235C任意函数发生器输出正弦信号，输出的信号直接作用到一侧压电片上，信号透过聚乙烯板内传播后到达另一侧压电片，再经过MDO4054B数字示波器采集响应数据进行分析。本文的扫频实验以及亚谐波阈值实验均采用此实验方案，示意图如图9所示。

3.2 激发亚谐波的阈值条件

由于亚谐波与材料固有频率有关，对试件进行扫频，结果如图10所示。从图10可以看出：频率小于473 kHz时，随着频率增加，响应信号幅值增大；频率大于473 kHz时，幅值随着频率增大而减小。由此可知，该试样的固有频率为473 kHz。根据共振原理，假定1/2亚谐波频率为473 kHz，当激励频率为亚谐波频率二倍时更容易激励出亚谐波，考虑测量误差，实验中选定5个靠近固有频率的二倍频率作为激励频率：944、945、946、947、948 kHz。激励信号幅值从200 mV开始，每100 mV记录一次数据，直到幅值增加到1 500 mV。以激励频率为946 kHz、激励幅值为800 mV的信号为例进行说明，如图11所示，可以看到时域图中有两个不同的极值，说明系统处于二倍周期状态，相应的频谱如图11（b）所示，出现了1/2亚谐波，这与第2.1节图3的结果一致。图11（b）中，Q为试验测得的电压值，mV。图11（b）还表明，含损伤的聚乙烯可以激励出1/2亚谐波，与利用高次谐波检测损伤类似，因此利用亚谐波响应构建损伤表征指标并进行损伤检测同样具有可行性。进一步根据频谱图获取聚乙烯试样中1/2亚谐波存在的频率‒幅值区间，如图12所示。

从图12可以看出：随着激励频率从943 kHz逐渐增大到949 kHz，出现亚谐波所需要的最小激励幅值先减小后增大，最低点出现在固有频率的二倍频即946 kHz处，这与第2.3节仿真结果图7所示基本一致，亚谐波出现所对应的响应信号幅值‒频率取值呈V字，在激励频率为亚谐波的二倍频率处为最低点，在此处的激励幅值更小，在该频率下更容易激励出1/2亚谐波，这一结论与文献[19,29]的结论一致。这为实际检测中亚谐波的激发条件和参数选取提供了理论基础，也为聚乙烯材料的塑性损伤识别提供了依据。

4 结 论

本文以一维非线性黏弹性杆为研究对象，从非线性动力学角度解释了一维非线性超声出现超谐波和亚谐波的条件，为非线性超声检测的激励参数选取提供了理论依据。主要结论如下：

1） 有限元仿真研究表明，非线性系统的分岔特性与超谐波和亚谐波具有强关联性。具体来说，叉式分岔前仅存在奇数次谐波，叉式分岔后出现偶数次谐波，这说明二次谐波的出现也是有条件的。倍周期分岔中，发生n次倍周期分岔，就对应产生1/n次亚谐波，系统进入混沌态时，系统响应的频谱结果中亚谐波的频带宽度将增加，不利于观察，但超谐波依然可见，这表明亚谐波有着比超谐波更加严格的阈值条件。

2） 有限元仿真研究还发现1/2亚谐波的激发条件幅值‒频率区间呈“V”字形。具体来说，通过对幅值‒频率区间的系统分岔特性进行扫描，得到了1/2亚谐波的激发条件，结果表明，出现亚谐波所需要的最小激励幅值先减小后增大，最低点出现在固有频率的二倍频处，整体呈“V”字形。

3） 实验研究表明，具有塑性损伤的聚乙烯板材作为典型的非线性黏弹性材料，可用于亚谐波阈值条件的研究。此外，在损伤的聚乙烯板材中，1/2亚谐波在频率‒幅值区间内呈“V”字形，在二倍频处激发亚谐波所需的激励幅值最小，这验证了仿真研究的结果。

参考文献

[1] Shui Guoshuang, Wang Yuesheng, Qu Jianmin.Advances in nondestructive test and evaluation of material degradation using nonlinear ultrasound[J].Advances in Mechanics,2005,35(1):52‒68.

[2] 税国双,汪越胜,曲建民.材料力学性能退化的超声无损检测与评价[J].力学进展,2005,35(1):52‒68.

[3] Yost W T, Cantrell J H.Materials characterization using acoustic nonlinearity parameters and harmonic generation:Engineering materials[M]//Review of Progress in Quantitative Nondestructive Evaluation.Boston:Springer,1990:1669‒1676. doi:10.1007/978-1-4684-5772-8_215

[4] Yun Hongguang, Rayhana R, Pant S,et al.Nonlinear ultrasonic testing and data analytics for damage characterization:A review[J].Measurement,2021,186:110155. doi:10.1016/j.measurement.2021.110155

[5] Nagy P B.Fatigue damage assessment by nonlinear ultrasonic materials characterization[J].Ultrasonics,1998,36(1/2/3/4/5):375‒381. doi:10.1016/s0041-624x(97)00040-1

[6] Cantrell J H, Yost W T.Nonlinear ultrasonic characterization of fatigue microstructures[J].International Journal of Fatigue,2001,23:487‒490. doi:10.1016/s0142-1123(01)00162-1

[7] Malfense Fierro G P, Ciampa F, Ginzburg D,et al.Nonlinear ultrasound modelling and validation of fatigue damage[J].Journal of Sound and Vibration,2015,343:121‒130. doi:10.1016/j.jsv.2014.10.008

[8] Wang Xiao, Wang Xue, Niu Xiaoguang,et al.Application of nonlinear ultrasonic technique to characterize the creep damage in ASME T92 steel welded joints[J].NDT&E International,2018,98:8‒16. doi:10.1016/j.ndteint.2018.04.006

[9] Pruell C, Kim J Y, Qu J,et al.A nonlinear-guided wave technique for evaluating plasticity-driven material damage in a metal plate[J].NDT&E International,2009,42(3):199‒203. doi:10.1016/j.ndteint.2008.09.009

[10] Liu Minghe, Kim J Y, Jacobs L,et al.Experimental study of nonlinear Rayleigh wave propagation in shot-peened aluminum plates—Feasibility of measuring residual stress[J].NDT&E International,2011,44(1):67‒74. doi:10.1016/j.ndteint.2010.09.008

[11] Zeitvogel D T, Matlack K H, Kim J Y,et al.Characterization of stress corrosion cracking in carbon steel using nonlinear Rayleigh surface waves[J].NDT&E International,2014,62:144‒152. doi:10.1016/j.ndteint.2013.12.005

[12] Frouin J, Sathish S, Matikas T E,et al.Ultrasonic linear and nonlinear behavior of fatigued Ti‒6Al‒4V[J].Journal of Materials Research,1999,14(4):12. doi:10.1557/jmr.1999.0176

[13] Li Long, Chen Jun.Characterization of optimal healing efficiency based on a new nonlinear ultrasonic method for the process of gradient heating healing in asphalt mixture[J].NDT&E International,2022,130:102674. doi:10.1016/j.ndteint.2022.102674

[14] Li Jianxun, Wang Minghang, Chen Haofeng,et al.Investigation on characteristics of tensile damage and microstructure evolution of steel AISI 316L by nonlinear ultrasonic Lamb wave[J].International Journal of Pressure Vessels and Piping,2022,198:104680. doi:10.1016/j.ijpvp.2022.104680

[15] Mao Hanling, Li Qingzhu, Mao Hanying,et al.Nonlinear ultrasonic characterization of carburized case depth[J].NDT&E International,2020,112:102244. doi:10.1016/j.ndteint.2020.102244

[16] Shan Shengbo, Cheng Li.Mixed third harmonic shear horizontal wave generation:Interaction between primary shear horizontal wave and second harmonic Lamb wave[J].Smart Materials and Structures,2019,28(8):085042. doi:10.1088/1361-665x/ab1fce

[17] Deng Mingxi, Xiang Yanxun, Liu Liangbing.Time-domain analysis and experimental examination of cumulative second-harmonic generation by primary Lamb wave propagation[J].Journal of Applied Physics,2011,109(11):113525. doi:10.1063/1.3592672

[18] Wang Zhi, Wu Weiliang, Xiao Li,et al.Damage detection for closed fatigue crack using subharmonic component[J].Journal of Vibration Engineering,2015,28(4):567‒573.

[19] 王芝,吴维亮,肖黎,等.闭合型疲劳裂纹的亚谐波损伤识别方法[J].振动工程学报,2015,28(4):567‒573.

[20] Wang Zhi, Qu Wenzhong, Xiao Li.Nonlinear ultrasonic damage detection for fatigue crack using subharmonic component[C].EWSHM-7th European Workshop on Structural Health Monitoring.Nantes:Innia Rennes-Bretagne Atlantique,2014.

[21] Zhang Mengyang, Xiao Li, Qu Wenzhong,et al.Damage detection of fatigue cracks under nonlinear boundary condition using subharmonic resonance[J].Ultrasonics,2017,77:152‒159. doi:10.1016/j.ultras.2017.02.001

[22] Zhu Weiqiu.Nonlinear wavesinelastic rods[J].Chinese Journal of Solid Mechanics,1980(2):247‒253.

[23] 朱位秋.弹性杆中的非线性波[J].固体力学学报,1980(2):247‒253.

[24] Parlitz U, Lauterborn W.Superstructure in the bifurcation set of the Duffing equation[J].Physics Letters A,1985,107(8):351‒355. doi:10.1016/0375-9601(85)90687-5

[25] Tang Yuanzhang, Lou Jingjun, Weng Xuetao,et al.Spectrum property of period-doubling bifurcations of Duffing oscillator[J].Journal of Vibration and Shock,2014,33(2):60‒63.

[26] 唐元璋,楼京俊,翁雪涛,等.Duffing振子倍周期分岔谱特性[J].振动与冲击,2014,33(2):60‒63.

[27] Ohara Y, Mihara T, Yamanaka K.Effect of adhesion force between crack planes on subharmonic and DC responses in nonlinear ultrasound[J].Ultrasonics,2006,44(2):194‒199. doi:10.1016/j.ultras.2005.10.006

[28] Eftekhari S A, Bakhtiari‒Nejad F, Dowell E H.Bifurcation boundary analysis as a nonlinear damage detection feature:Does it work?[J].Journal of Fluids and Structures,2011,27(2):297‒310. doi:10.1016/j.jfluidstructs.2010.11.006

[29] Hamilton M F, Blackstock D T.Nonlinear acoustics[M].San Diego:Springer Cham,1998.

[30] Skar S.Odd harmonic periodic solutions of systems with odd nonlinearities[J].IEEE Transactions on Circuits and Systems,1985,32(12):1281‒1285. doi:10.1109/tcs.1985.1085670

[31] Matlack K H, Kim J Y, Jacobs L J,et al.Review of second harmonic generation measurement techniques for material state determination in metals[J].Journal of Nondestructive Evaluation,2014,34(1):273. doi:10.1007/s10921-014-0273-5

[32] Wang Jing.Numerical and experimental study on 1/2 subharnonic propagation in nonlinear viscoelastic materials[D].Taiyuan:Taiyuan University of Science and Technology,2022.

[33] 王晶.粘弹性材料中的1/2亚谐波数值及实验研究[D].太原:太原科技大学,2022.

[34] Sojahrood A J, Earl R, Kolios M C,et al.Investigation of the 1/2 order subharmonic emissions of the period-2 oscillations of an ultrasonically excited bubble[J].Physics Letters A,2020,384(23):126446. doi:10.1016/j.physleta.2020.126446

[35] 【校友简介】. doi:10.1016/j.physleta.2020.126446

基金资助

东莞市社会发展科技项目(20231800940332)