一种超声导波频散曲线计算方法

摘要

本发明提供了一种超声导波频散曲线计算方法，属于超声波测量技术领域。随着待测板结构的厚度增加、结构复杂度增加，在构建待测板结构的有限元模型时，若采用传统有限元分析方法，则需要采用高阶多项式插值函数，而且需要划分较多数量的单元，导致有限元网格的数量增加，计算效率降低。本发明通过采用区间B样条小波单元的形函数作为插值函数进行求解，能够保证计算精度，而且无需划分较多数量的单元，有限元网格的数量减少，能够提高计算效率。

说明书

技术领域

本发明涉及一种超声导波频散曲线计算方法，属于超声波测量技术领域。

背景技术

超声导波检测是在传统超声检测技术发展而来的一种新型无损检测技术。与传统超声体波相比，导波只需要单点激励便可完成对整个被测介质的检测，且导波具有传播速度快、能量衰减小和检测效率高等优点，因此被广泛应用于板结构的快速检测。导波在不同波导介质中的传播特性包含了超声导波的一些基本信息，主要通过建立相关模型对超声导波传播的相速度、群速度频散曲线进行求解，从而为实际检测中导波频率的选取、探头的斜入射角度的选取等提供理论支持。

目前，板结构的频散曲线求解方法主要是通过求解特征方程和有限元法得到频率和波数的关系，进而绘制出被测板结构的相速度和群速度频散曲线。而采用有限元方法时，通常需要对有限元网格进行划分，并采用高阶多项式基函数进行求解。由于有限元网格的尺寸划分一般要控制在λ/20～λ/10(λ为最小波长)，随着板结构的厚度增加，有限元网格数量增加，导致求解过程的计算效率低。

发明内容

本发明的目的在于提供一种超声导波频散曲线计算方法，用于解决板厚度增加导致求解频散曲线的计算效率低的问题。

为了实现上述目的，本发明提供了一种超声导波频散曲线计算方法，包括如下步骤：

S1、在待测板结构的厚度方向上建立有限元模型，并对有限元模型进行单元划分；

S2、利用区间B样条小波单元构造的形函数，计算有限元模型中各单元的单元位移；

S3、根据步骤S2计算出的单元位移，计算对应单元的单元应变；

S4、根据各单元的单元位移和单元应变，计算对应单元的刚度矩阵和质量矩阵；

S5、根据各单元的刚度矩阵和质量矩阵，建立波动方程，并对波动方程进行求解，得到超声导波的波数与频率之间的关系，进而绘制出超声导波在待测板结构中传输时的频散曲线。

随着待测板结构的厚度增加，在构建待测板结构的有限元模型时，若采用传统有限元分析方法，则需要采用高阶多项式插值函数，而且需要划分较多数量的单元，导致有限元网格的数量增加，计算效率降低。本发明通过采用区间B样条小波单元的形函数作为插值函数进行求解，能够保证计算精度，而且无需划分较多数量的单元，有限元网格的数量减少，能够提高计算效率。

进一步地，在上述方法中，步骤S2中所述的区间B样条小波单元构造的形函数通过如下公式得到：

式中，

提供一种具体的方法来获取区间B样条小波单元的形函数，便于实现。

进一步地，在上述方法中，步骤S2中，通过如下公式计算有限元内第e个单元的单元位移：

式中，u

提供一种具体的方法来计算单元位移，便于实现。

进一步地，在上述方法中，步骤S3中，通过如下公式计算各单元的单元应变：

式中，ε

提供一种具体的方法来计算单元应变，便于实现。

进一步地，在上述方法中，步骤S4中，通过将单元位移和单元应变代入哈密顿公式的离散形式中，求出各单元的刚度矩阵和质量矩阵；所述哈密顿公式的离散形式为：

式中，t

提供一种具体的计算方法来求解刚度矩阵和质量矩阵，便于本发明的实施。

进一步地，在上述方法中，将单元位移和单元应变代入哈密顿公式的离散形式后得到如下公式，对该公式进行求解得到对应单元的刚度矩阵和质量矩阵；

式中，δq

进一步地，在上述方法中，步骤S5中，建立波动方程时，引入特定算子符号Γ

式中，T

通过引入特定算子对刚度矩阵和质量矩阵进行化简，降低后续计算过程的计算量，提高计算效率。

进一步地，在上述方法中，通过如下公式将刚度矩阵和质量矩阵整合为整体刚度矩阵和整体质量矩阵，进而通过整体刚度矩阵和整体质量矩阵建立波动方程；

式中，K

进一步地，在上述方法中，所述波动方程为：

[K

式中，κ表示超声导波的波数，i表示虚数符号，K

进一步地，在上述方法中，在求解波动方程时，引入对角矩阵T

式中，下标O为总自由度数。

在求解波动方程时，引入对角矩阵T

附图说明

图1为本发明方法实施例中超声导波频散曲线计算方法的流程框图；

图2为本发明方法实施例中4阶3尺度的区间B样条小波尺度函数示意图；

图3为本发明方法实施例中单元内节点数量与节点坐标的对应关系示意图；

图4为本范发明方法实施例中板结构的半解析小波有限元模型；

图5为本发明方法实施例中采用本发明与采用常规有限元法求解各向均质板结构频散曲线的对比图；

图6为本发明方法实施例中采用本发明与采用常规有限元法求解各向均值板结构频散曲线的时间对比图；

图7为本发明方法实施例中复合材料层合板的模型示意图；

图8为本发明方法实施例中采用本发明与采用常规有限元法求解复合材料层合板相速度频散曲线的对比图；

图9为本发明方法实施例中采用本发明与采用常规有限元法求解复合材料层合板群速度频散曲线的对比图；

图10为本发明方法实施例中采用本发明与采用常规有限元法求解复合材料层合板频散曲线的时间对比图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明了，以下结合附图及实施例，对本发明进行进一步详细说明。

方法实施例：

如图1所示，本发明的超声导波频散曲线计算方法包括如下步骤：

1、利用半解析法，建立待测板结构的有限元模型，有限元模型与待测板结构在厚度方向的尺寸相同。同时根据求解精度的需求，选用不同尺度的区间B样条小波函数，并得到该区间B样条小波单元的尺度函数，以及一维区间B样条小波单元的单元转换矩阵。利用区间B样条小波单元的尺度函数和一维区间B样条小波单元的单元转换矩阵，构建一维区间B样条小波单元的形函数。

如图2所示为待测板结构的半解析小波有限元模型，X方向为超声导波的波传播方向，y-z为待测板结构的截面，其中y方向为无限宽。超声导波的波数为κ，频率为ω。

本实施例中，根据求解精度的需求，选取4阶3尺度的区间B样条小波函数。4阶3尺度区间B样条小波函数的尺度函数共有11个，分别为：0边界尺度函数

区间B样条小波(BSWI)函数是通过给定简单的节点序列并在节点之间使用分段多项式函数来将它们连接在一起，从而获得一定阶的整体平滑度来构造的。对于处在求解区间[a,b]上的任意一维函数f(x)，可以通过线性映射关系ξ＝(x-a)/(b-a)转换到标准求解区间[0,1]，因此只需要在求解区间[0,1]上构造m阶B样条小波函数。为了保证在区间[0,1]内至少有一个内小波，必须满足以下条件：

2

其中，m和j分别为区间B样条小波函数的阶数和尺度参数。通过m阶0尺度的区间B样条小波尺度函数，通过如下递推公式，可以得到任意j尺度的区间B样条小波尺度函数

因此，区间[0,1]上区间B样条小波函数的尺度函数的行向量可以表示为：

其中ξ的取值区间为[0,1]。如图3所示，即为4阶3尺度区间B样条小波函数的尺度函数。

当采用m阶j尺度的区间B样条小波函数作为插值函数构造单元时，单元上的节点数量和每个节点相应的坐标如图4所示。长度为l

定义坐标转换公式为：ξ＝(x-x

将实际坐标值代入坐标转换公式中，可以得到每一个节点坐标x

ξ

当采用m阶j尺度的区间B样条小波函数作为插值函数时，未知场函数u(ξ)可表示为：

式中，

假设物理自由度的列向量为

式中，R

将u

式中，

2、在待测板结构的有限元模型上，利用小波有限元方法进行离散，将有限元模型在厚度方向上的截面Ω划分为有限个单元Ω

待测板结构的各单元中每个点的简谐位移、应力和应变分量可表示为：

u＝[u

应力与应变的本构关系为σ＝Cε，其中，C为待测板结构的弹性常数矩阵。应变和位移的关系式为：

假设超声导波沿传播方向x的位移场为简谐波，通过空间函数表示单元内任意一点的位移，表示为：

当采用m阶j尺度的区间B样条小波函数的尺度函数作为插值函数构造单元时，每个单元内有n+1个节点。利用区间B样条小波单元形函数，表示单元内任意一点的单元位移。单元位移可表示为：

式中，

3、根据单元位移与单元应变的关系，求得每个单元的单元应变。

将单元位移代入应变与位移的关系式中，得出单元的单元应变为：

4、将单元位移和单元应变代入哈密顿公式的离散形式中，分别求得每个单元的刚度矩阵和质量矩阵。

记截面Ω离散化的得到的有限个单元Ω

将单元位移和单元应变代入上式中，可得：

式中，t

式中，

为方便计算，引入特定算子符号Γ

利用特定算子，对刚度矩阵和质量矩阵进行化简。以

Γ

5、将每个单元的刚度矩阵和质量矩阵整合为整体刚度矩阵，代入波动方程中，并对波动方程进行求解，得到相速度和群速度的频散曲线。

通过如下公式将刚度矩阵和质量矩阵整合为整体刚度矩阵和整体质量矩阵。

式中，K

将整体刚度矩阵和整体质量矩阵代入波动方程：

[K

式中，κ表示超声导波的波数，i表示虚数符号，K

由于公式中的虚数i会为求解带来很大难度，因此引入对角矩阵T

式中，下标O为总自由度数。令

通过对波动方程进行求解，给定每一个频率可以求得对应波数k的2M个特征值。所有的特征值中，实数特征值的个数即为对应的波数。利用波数和频率的关系，经过一定的数值变换可绘制出相速度和群速度频散曲线。

以各相同性均质铝板和各相异性复合材料板为例，对采用本发明的效果进行展示。

各相同性均质铝板的材料参数为：厚度1mm、弹性模量68.9Gpa、泊松比0.33、密度2690kg/m

如图7所示为各向异性复合材料板的模型，为了更进一步体现此发明在保证计算精度的前提下计算效率高的优点。采用常规有限元法进行分析时，在厚度方向上划分40个单元，而采用本发明后，在厚度方向上构建8个单元，节点总数相同。从图8和图9可以看出，与采用常规有限元法得出的频散曲线相比，采用本发明后得出的群速度和相速度频散曲线具有同样的精度，而从图10可以看出，采用本发明后，计算时间从476s降低至226s，计算效率大幅增加。