CN109117500B - 一种基于薄层离散计算层状板中Lamb波频散曲线的方法 - Google Patents

一种基于薄层离散计算层状板中Lamb波频散曲线的方法 Download PDF

Info

Publication number
CN109117500B
CN109117500B CN201810720164.XA CN201810720164A CN109117500B CN 109117500 B CN109117500 B CN 109117500B CN 201810720164 A CN201810720164 A CN 201810720164A CN 109117500 B CN109117500 B CN 109117500B
Authority
CN
China
Prior art keywords
thin layer
matrix
layer
calculating
wave
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201810720164.XA
Other languages
English (en)
Other versions
CN109117500A (zh
Inventor
柴华友
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Wuhan Institute of Technology
Original Assignee
Wuhan Institute of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Wuhan Institute of Technology filed Critical Wuhan Institute of Technology
Priority to CN201810720164.XA priority Critical patent/CN109117500B/zh
Publication of CN109117500A publication Critical patent/CN109117500A/zh
Application granted granted Critical
Publication of CN109117500B publication Critical patent/CN109117500B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Investigating Or Analyzing Materials By The Use Of Ultrasonic Waves (AREA)

Abstract

本发明涉及一种基于薄层离散计算层状板中Lamb波频散曲线的方法，包括以下步骤，在层状板的下方构建一虚拟介质层，并在虚拟介质层底面设置刚性基；对层状板和虚拟介质层进行水平向离散，得到多个薄层；分别计算每一薄层的刚度矩阵Kn；将各薄层的刚度矩阵Kn集成为总刚度矩阵K；由各薄层的交界面位移及外力分量分别构筑薄层的交界面的总位移矢量U及总外力矢量F，则有KU＝F；在自由状态下，总外力矢量F＝0，KU＝F可写为KU＝0，求解KU＝0的特征值，并对特征值进行筛选，由筛选出的特征值计算出Lamb波的相速度。本发明将层状板中Lamb波频散曲线传统矩阵行列式根搜索法转化成代数矩阵特征值求解方法，可以有效克服根搜索算法存在搜索范围不收敛、根遗漏问题。

Description

一种基于薄层离散计算层状板中Lamb波频散曲线的方法
技术领域
本发明可用于复合材料板、混凝土板、金属材料板中Lamb传播特性分析，Lamb波传播特性可广泛应用于板介质力学参数分析、板耦合介质选择、层耦合效果分析以及板缺陷探测。
背景技术
在单层板或层状板的表面施加动荷载，板的动力响应与板中Lamb波有关，研究板中Lamb波频散特性有助于对板中波场分析。Lamb波频散特性与层状板的分层结构及层材料力学参数有关，利用Lamb波频散特性反过来可以对板层结构及层材料力学参数分析，当板中存在裂纹、空洞或异质体等异常，Lamb波会在这些异常处反射，不同模态Lamb波反射特征不同，由Lamb波频散特性可以分析反射波特征，并由此可确定板中异常位置、程度及范围。对均匀单层板，已有全频率域Lamb波对称及反对称模态频散方程，由此可得到Lamb波对称及反对称模态频散曲线(即相速度随频率变化)。对两层及以上层状板，由系数矩阵、传递矩阵法或刚度矩阵法，利用自由边界条件得到矩阵行列式，对行列式求解可以得到层状板中Lamb波各模态频散曲线。然而，矩阵行列式一些矩阵元素是波数及厚度超越函数，行列式求解需要根搜索方法，即：设定目标函数，通过优化方法，不断收敛搜索范围，当前、后步目标函数差值小于设定值，得到行列式根。由于存在多个根，根搜索非常耗时，当搜索起始值设置不当，甚至会出现遗漏根和搜索范围不收敛的问题。
发明内容
本发明所要解决的技术问题是提供一种基于薄层离散计算层状板中Lamb波频散曲线的方法，可以克服现有层状板中Lamb波频散曲线计算矩阵方法行列式求根时搜索范围不收敛和根遗漏的问题，从而可以实现高效快速得到Lamb波频散曲线的目的。
本发明解决上述技术问题的技术方案如下：一种基于薄层离散计算层状板中Lamb波频散曲线的方法，包括以下步骤，
S1，在总厚度为H的层状板的下方构建一虚拟介质层，并在虚拟介质层底面设置刚性基；
S2，对层状板和虚拟介质层构成的整体结构进行水平向离散，得到多个薄层；
S3，分别计算每一薄层的刚度矩阵Kn；
S4，根据各薄层的交界面位移连续，将各薄层的刚度矩阵Kn集成为总刚度矩阵K；
S5，由各薄层的交界面位移及外力分量分别构筑薄层的交界面的总位移矢量U及总外力矢量F，则有
KU＝F (1)；
S6，在自由状态下，总外力矢量F＝0，式(1)可写为
KU＝0 (1-1)，
用矩阵分解方法求解式(1-1)的特征值，并对特征值进行筛选，筛选出与Lamb波特性对应的特征值，由筛选出的特征值计算出Lamb波的相速度。
其中：
在所述S1中，所述虚拟介质层的剪切波速cbs＝0.1m/s，密度ρbs＝1.25kg/m3，泊松比vbs＝0.45，即以虚拟介质模拟空气，且所述虚拟介质层的厚度与层状板的厚度相同。
设Lamb波频散曲线分析频率范围用符号[fmin，fmax]表示，其中，fmin及fmax分别代表频率范围的下限及上限，
所述S2具体为，
S21，由层状板中各层板材料的泊松比νs及剪切波速cs，按关系式
Figure GDA0003890626610000031
计算层状板中各层板材料的纵波波速cp，且层状板中各层板材料的纵波速cp中的最小值用符号cp，min表示；
S22，由层状板中各层板材料的纵波速cp中的最小值cp，min，以及频率范围[fmin，fmax]中最大频率fmax，计算纵波最小波长
λmin＝cp，min/fmax (3)；
S23，按照厚度h≤λmin/20对层状板和虚拟介质层构成的整体结构进行水平向离散。
由第n薄层的上、下交界面位移分量及应力分量分别构建该薄层的广义位移向量un及广义应力向量pn，其中
un＝[ux,n,ux,n+1,iuz,n,iuz,n+1]T (4)，
pn＝[-σxz,n,σxz,n+1,-iσz,n,iσz,n+1]T (5)，
其中，ux,n和uz,n分别表示第n薄层的上交界面水平向位移分量及竖直向位移分量，ux,n+1和uz,n+1分别表示第n薄层的下交界面水平向位移分量及竖直向位移分量，σxz,n和σz,n分别表示第n薄层的上交界面剪切应力及法向应力，σxz,n+1和σz,n+1分别表示第n薄层的下交界面剪切应力及法向应力，虚数
Figure GDA0003890626610000032
所述S3具体为，薄层内位移由该薄层的上、下交界面的位移线性插值计算，利用虚功原理得到薄层刚度矩阵，第n薄层的刚度矩阵Kn可表示为，
Kn＝Ank2+Bnk+Cn (6)，
式中k是波数，An、Bn和Cn分别为，
Figure GDA0003890626610000041
Figure GDA0003890626610000042
Figure GDA0003890626610000043
其中，ω是角频率，hn为第n薄层的厚度，μn和λn均为第n薄层的拉梅常数，ρn为第n薄层的密度,矩阵Bxzn T是矩阵Bxzn的转折。
所述S4具体为，根据薄层的交界面位移连续条件，将由式(6)表示的各薄层刚度矩阵集合成总刚度矩阵K，总刚度矩阵K可表示为
Figure GDA0003890626610000044
在式(7)中，矩阵
Figure GDA0003890626610000045
是式(6-1)中各薄层子矩阵Axn的集成；矩阵
Figure GDA0003890626610000046
是式(6-1)中各薄层子矩阵Azn的集成；矩阵
Figure GDA0003890626610000047
是式(6-2)中各薄层子矩阵Bxzn的集成；矩阵
Figure GDA0003890626610000048
是式(6-3)中各薄层子矩阵Cxn的集成；矩阵
Figure GDA0003890626610000049
是式(6-3)中各薄层子矩阵Czn的集成。
在所述S5中，薄层的交界面总位移矢量U用矢量
Figure GDA00038906266100000410
代替，
Figure GDA00038906266100000411
如下式
Figure GDA00038906266100000412
其中，矢量Φx和Φz分别为
Φx＝[ux,1,ux,2,…,ux,N-1,u,x,N]T (8-1)，
Φz＝i[uz,1,uz,2,…,uz,N-1,u,z,N]T (8-2)，
上标符号T表示转折，N为薄层的交界面总数，ux,n和uz,n分别表示第n薄层的上交界面的水平向和竖直向位移分量，n＝1,2,3...，N；
利用刚性基位移为零条件，式(1-1)可改写为
Figure GDA0003890626610000051
其中
Figure GDA0003890626610000052
为零矢量。
在所述S6中，具体用矩阵分解法求取式(1-1-1)的特征值；
所述S6具体为，
S61，给定频率后，由式(6-1)～(6-3)计算各薄层子矩阵Axn,Azn,Bxzn,Cxn,Czn，将各薄层子矩阵分别集成得到矩阵
Figure GDA0003890626610000053
并利用数学工具软件Matlab中的函数eig计算出式(1-1-1)的特征值1/k2，其中，特征值1/k2有2N个，2N个特征值1/k2对应有4N个特征波数k，设km为4N个特征波数中的任一个，且m＝1，2，3，...，4N；
S62，筛选出实部大于零、虚部小于零且实部是虚部绝对值十倍以上的特征波数，设筛选出来的特征波数有M个，同时设kj为筛选出来的M个特征波数中的任一个，且j＝1,2,3，...，M；
S63，根据筛选出来的特征波数kj的实部按下式计算Lamb波的相速度cj
cj＝ω/Re(kj) (9)，
这里，符号Re表示实部；
S64，将Lamb波的相速度cj由小到大排序，并重新编号为c1<c2<...<cM，这些值依次对应于该频率Lamb波不同阶次模态相速度；
S65，改变频率，重复以上步骤S61～S64，得到在[fmin,fmax]范围内不同频率下Lamb波各模态相速度，进而得到Lamb波各模态频散曲线。
在利用Matlab中的函数eig计算式(1-1-1)的特征值1/k2的同时，还得到特征向量
Figure GDA0003890626610000061
根据筛选出来的特征波数所对应的特征向量
Figure GDA0003890626610000062
中的位移分量矢量，得到Lamb波各模态竖直向位移沿厚度方向变化形状。
在单层或以板中心平面为对称面的板中，由竖直向模态位移形状判断板中Lamb波模态为对称还是反对称类型。
本发明的有益效果是：本发明一种基于薄层离散计算层状板中Lamb波频散曲线的方法就是在层状板下方构建虚拟介质层并在虚拟介质层底面施加刚性基，将层状板及虚拟介质层离散成厚度相对波长很小的多个薄层，这样，可将薄层刚度矩阵简化成波数的代数矩阵，给定频率及材料特性参数后，通过矩阵分解方法得到特征值，由特征值可计算相速度，由特征向量可得模态位移沿厚度方向变化形状。对计算相速度筛选、排序可得到Lamb波各模态相速度，对单层板或以板中心面为对称的板，由模态位移形状可确定模态是对称还是反对称类型；本发明将现有矩阵行列式根搜索法转化成用矩阵分解法求解特征值及特征向量方法，可以有效克服现有层状板中Lamb波频散曲线计算根搜索算法存在搜索范围不收敛、根遗漏问题，从而实现高效快速得到Lamb波频散曲线的目的。
附图说明
图1为本发明的流程图；
图2为本发明中由层状板到构建虚拟介质层再到薄层离散的状态示意图；
图3为本发明中Lamb波各模态相速度随频率变化曲线图；
图4为本发明中150Hz处Lamb波各模态竖直向位移沿厚度方向变化形状图；
图5为本发明中Lamb波各模态归一化相速度随无量纲厚度变化计算曲线图；
图6为本发明中Lamb波各模态归一化相速度随无量纲厚度变化理论曲线图。
具体实施方式
以下结合附图对本发明的原理和特征进行描述，所举实例只用于解释本发明，并非用于限定本发明的范围。如图1所示，一种基于薄层离散计算层状板中Lamb波频散曲线的方法，包括以下步骤，
S1，在总厚度为H的层状板的下方构建一虚拟介质层，并在虚拟介质层底面设置刚性基；
S2，对层状板和虚拟介质层构成的整体结构进行水平向离散，得到多个薄层；
S3，分别计算每一薄层的刚度矩阵Kn；
S4，根据各薄层的交界面位移连续，将各薄层的刚度矩阵Kn集成为总刚度矩阵K；
S5，由各薄层的交界面位移及外力分量分别构筑薄层的交界面的总位移矢量U及总外力矢量F，则有
KU＝F (1)；
S6，在自由状态下，总外力矢量F＝0，式(1)可写为
KU＝0 (1-1)，
用矩阵分解方法求解式(1-1)的特征值，并对特征值进行筛选，筛选出与Lamb波特性对应的特征值，由筛选出的特征值计算Lamb波的相速度。
以下内容为对本发明方法中的各个步骤进行详细的解释说明。
在所述S1中，所述虚拟介质层的剪切波速cbs＝0.1m/s，密度ρbs＝1.25kg/m3，泊松比vbs＝0.45，即以虚拟介质模拟空气，且所述虚拟介质层的厚度与层状板的厚度相同。
设Lamb波频散曲线计算频率范围用符号[fmin，fmax]表示，其中，fmin及fmax分别代表频率范围的下限及上限，
所述S2具体为，
S21，由层状板中各层板材料的泊松比νs及剪切波速cs，按关系式
Figure GDA0003890626610000081
计算层状板中各层板材料的纵波波速cp，且层状板中各层板材料的纵波速cp中的最小值用符号cp，min表示；
S22，由层状板中各层板材料的纵波速cp中的最小值cp，min，以及频率范围[fmin，fmax]中最大频率fmax，计算纵波最小波长
λmin＝cp，min/fmax (3)；
S23，按照厚度h≤λmin/20对层状板和虚拟介质层构成的整体结构进行水平向离散。
在S23中，对整体结构(整体结构为层状板+虚拟介质层)进行水平向离散；这里的离散就是分层的意思，按厚度h离散，就是分层后每层的厚度为h。图2给出了由层状板到构建虚拟介质层再到薄层离散的状态示意图。
在所述S3中：
由第n薄层的上、下交界面位移分量及应力分量分别构建该薄层的广义位移向量un及广义应力向量pn，其中
un＝[ux,n,ux,n+1,iuz,n,iuz,n+1]T (4)，
pn＝[-σxz,n,σxz,n+1,-iσz,n,iσz,n+1]T (5)，
其中，ux,n和uz,n分别表示第n薄层的上交界面水平向位移分量及竖直向位移分量，ux,n+1和uz,n+1分别表示第n薄层的下交界面水平向位移分量及竖直向位移分量，σxz,n和σz,n分别表示第n薄层的上交界面剪切应力及法向应力，σxz,n+1和σz,n+1分别表示第n薄层的下交界面剪切应力及法向应力，虚数
Figure GDA0003890626610000095
(图2给出薄层编号及薄层交界面编号，其中，图2左边的为薄层编号，右边的为薄层交界面编号)
所述S3具体为，薄层内位移由该薄层的上、下交界面的位移线性插值计算，利用虚功原理得到薄层刚度矩阵，第n薄层的刚度矩阵Kn可表示为，
Kn＝Ank2+Bnk+Cn (6)，
式中k是波数，矩阵An、Bn和Cn分别为，
Figure GDA0003890626610000091
Figure GDA0003890626610000092
Figure GDA0003890626610000093
其中，ω是角频率，hn为第n薄层的厚度，μn和λn均为第n薄层的拉梅常数，ρn为第n薄层的密度,矩阵Bxzn T是矩阵Bxzn的转折。
在上述公式(6-1)、(6-2)和(6-3)中，子矩阵Axn、Azn、Bxzn、Cxn、Czn和零矩阵0分别为
Figure GDA0003890626610000094
Figure GDA0003890626610000101
Figure GDA0003890626610000102
Figure GDA0003890626610000103
Figure GDA0003890626610000104
Figure GDA0003890626610000105
在所述S3中，将薄层内位移用层上、下交界面位移线性插值的过程为：假设第n薄层的上、下交界面位移矢量分别为
Figure GDA0003890626610000106
Figure GDA0003890626610000107
则该薄层内任一位置位移矢量
Figure GDA0003890626610000108
为
Figure GDA0003890626610000109
在式(6-12)中，
Figure GDA00038906266100001010
Figure GDA00038906266100001011
薄层厚度hn＝zn+1-zn。
在所述S4中：
所述S4具体为，根据薄层的交界面位移连续条件，将由式(6)表示的各薄层刚度矩阵集合成总刚度矩阵K，
且总刚度矩阵K的计算公式为
Figure GDA00038906266100001012
在式(7)中，矩阵
Figure GDA0003890626610000111
是式(6-1)中各薄层子矩阵Axn的集成，各薄层子矩阵Axn的表达式如式(6-4)；矩阵
Figure GDA00038906266100001112
是式(6-1)中各薄层子矩阵Azn的集成，各薄层子矩阵Azn的表达式如式(6-5)；矩阵
Figure GDA0003890626610000112
是式(6-2)中各薄子矩阵Bxzn的集成，各薄层子矩阵Bxzn的表达式如式(6-6)；矩阵
Figure GDA0003890626610000113
是式(6-3)中各薄层子矩阵Cxn的集成，各薄层子矩阵Cxn的表达式如式(6-7)；矩阵
Figure GDA0003890626610000114
是式(6-3)中各薄层子矩阵Czn的集成，各薄层子矩阵Czn的表达式如式(6-8)。在所述S5中：
薄层的交界面总位移矢量U用矢量
Figure GDA0003890626610000115
代替，
Figure GDA0003890626610000116
如下式
Figure GDA0003890626610000117
其中，矢量Φx和Φz分别为
Φx＝[ux,1,ux,2,…,ux,N-1,u,x,N]T (8-1)，
Φz＝i[uz,1,uz,2,…,uz,N-1,u,z,N]T (8-2)，
上标符号T表示转折，N为薄层的交界面总数，ux,n和uz,n分别表示第n薄层的上交界面的水平向和竖直向位移分量，n＝1,2,3...，N；
利用刚性基位移为零条件，式(1-1)可改写为
Figure GDA0003890626610000118
其中，
Figure GDA0003890626610000119
为零矢量。
在所述S6中，具体用矩阵分解法求取式(1-1-1)的特征值；
所述S6具体为，
S61，给定预设的频率后，由式(6-1)～(6-3)计算各薄层子矩阵Axn，Azn，Bxzn，Cxn，Czn，将各薄层子矩阵分别集成得到矩阵
Figure GDA00038906266100001110
Figure GDA00038906266100001113
利用数学工具软件Matlab中的函数eig计算出式(1-1-1)的特征值1/k2，其中，由于N个薄层共有2N个自由度，所以特征值1/k2有2N个，2N个特征值1/k2对应有4N个特征波数k，设km为4N个特征波数中的任一个，且m＝1，2，3，...，4N；
S62，在4N个特征波数中，其中有一半(2N)特征波数的实部为负值，实部为负值的特征波数对应由远处向中心传播的波，由于只考虑由中心向外传播的波，这些实部为负值的特征波数应舍去；舍去这些实部为负值的特征波数后，余下一半的特征波数的实部大于零且虚部小于零，为了确保波可以传播较远，再从余下一半中筛选出实部是虚部绝对值十倍以上的特征波数，即Re(kj)＞-10Im(kj)，这里符号Re和Im分别表示实部和虚部，用符号kj(j＝1,2,…M)表示筛选出来的特征波数，M是筛选出来的特征波数的数量；
S63，根据筛选出来的特征波数kj的实部得到Lamb波的相速度cj，即，
cj＝ω/Re(kj) (9)，
另外，由筛选出来的特征波数kj的虚部可以得到衰减系数αj，其中αj＝-Im(kj)；
S64，将Lamb波的相速度cj由小到大排序，并重新编号为c1<c2<...<cM，这些值依次对应于该频率Lamb波不同阶次模态相速度；
S65，改变频率，重复以上步骤S61～S64，得到在[fmin,fmax]范围内不同频率下对应Lamb波各模态相速度，进而得到Lamb波不同模态频散曲线。例如，将频率范围[fmin,fmax]等间隔离散成L个点，第j个点对应频率fj＝fmin+(j-1)Δf，(j＝1,2，…L)，这里Δf＝(fmax-fmin)/(L-1)为相邻点频率间隔，按上式步骤得到频率fj处相速度，重复以上步骤S61～S64进而得到在频率范围[fmin，fmax]内，Lamb波各模态频散曲线。
在利用Matlab中的函数eig计算式(1-1-1)的特征值1/k2同时，还得到特征向量
Figure GDA0003890626610000121
根据保留下来的特征波数所对应的特征向量
Figure GDA0003890626610000122
中的位移分量矢量，得到Lamb波各模态竖直向位移沿厚度方向变化形状。
在单层或以板中心平面为对称面的板中，由模态竖直向位移形状判断板中Lamb波模态为对称还是反对称类型。
本发明通过矩阵分解方法得到总刚度矩阵K的特征值及特征向量，特征值对应于特征波数平方的倒数，特征向量对应模态位移形状。对特征波数筛选、排序可计算Lamb波各模态相速度，对单层板或以板中心面为对称的板，由模态位移形状可确定模态是对称还是反对称类型。
在本具体实施例中：
单层板是层状板的一种特例，单层板中的对称及反对称Lamb波有解析解，可以计算对称及反对称Lamb波频散曲线，以下通过本方法计算结果与理论曲线比较验证本发明方法。
A:一均匀板厚度1m,材料剪切波速cs＝130m/s,密度ρs＝1800kg/m3,泊松比νs＝0.3。
B:虚拟介质层的剪切波速cbs＝0.1m/s，密度ρbs＝1.25kg/m3，泊松比νbs＝0.45。
C:取虚拟介质厚度1m，在其底面施加刚性基。
D:板纵波速按式(2)计算得
Figure GDA0003890626610000131
分析频率范围取为[5,500]Hz，按式(3)计算纵波最小波长λmin＝243/500≈0.486m；
E:薄层厚度取0.02m，满足薄层厚度h<λmin/20条件，对板及虚拟介质层水平向离散。
F:按式(1-1-1)求Lamb波的特征值并进行筛选。
G:按式(9)计算Lamb波的相速度随频率变化曲线如图3，在频率150Hz处，Lamb波各模态竖直向位移沿厚度方向变化形状如图4，由此可判断模态为对称还是反对称类型，图4中A0、S0分别表示零阶反对称、对称类型，A1、S1分别表示一阶反对称、对称类型。
H:按瑞利波速cR与剪切波速cs及泊松比νs回归关系式
Figure GDA0003890626610000132
由板剪切波速cs和泊松比νs计算板瑞利波速cR，不同频率瑞利波波长λR(f)＝cR/f。将板厚度对不同频率波长λR(f)归一化，将Lamb波各模态相速度对剪切波速cs＝130m/s进行归一化，归一化Lamb波各模态计算频散曲线见图5，理论曲线如图6所示，两者计算结果一致，本发明方法有效性得到验证。
在实际工程中，振源能量一般介于5～10000Hz有限频率范围内，这样，只需关注激发频率范围而不需关注全频率范围(低频趋于0，高频趋于无穷)内的Lamb波频散特性。本发明就是在于解决在工程应用频率范围内层状板中Lamb波频散曲线计算。本发明目的在于克服现有矩阵行列式求解层状板中Lamb波频散曲线困难，为层状板振动信号分析、解释以及利用Lamb波进行板结构探伤提供一种高效快速计算Lamb波频散曲线方法。
本发明就是在层状板下方构建虚拟介质层及在虚拟介质层底面施加刚性基，将层状板及虚拟介质层离散成厚度相对波长很小薄层，这样，可将薄层刚度矩阵简化成波数的代数矩阵，给定频率及材料参数后，通过矩阵分解方法得到总刚度矩阵的特征值及特征向量，特征值对应于波数平方的倒数，特征向量对应模态位移形状。对波数筛选、排序可计算Lamb波各模态相速度，对以板中心面为对称的板，由模态位移形状可确定模态是对称还是反对称类型。本发明将Lamb波频散曲线计算现有的行列式根搜索法转化成矩阵分解方法，通过矩阵分解方法计算代数矩阵的特征值及特征向量，从而避免了现有矩阵行列式根搜索方法存在各种问题。
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (9)

1.一种基于薄层离散计算层状板中Lamb波频散曲线的方法，其特征在于：包括以下步骤，
S1，在总厚度为H的层状板的下方构建一虚拟介质层，并在虚拟介质层底面设置刚性基；
S2，对层状板和虚拟介质层构成的整体结构进行水平向离散，得到多个薄层；
S3，分别计算每一薄层的刚度矩阵Kn；
S4，根据各薄层的交界面位移连续，将各薄层的刚度矩阵Kn集成为总刚度矩阵K；
S5，由各薄层的交界面位移及外力分量分别构筑薄层的交界面的总位移矢量U及总外力矢量F，则有
KU＝F (1)；
S6，在自由状态下，总外力矢量F＝0，式(1)可写为
KU＝0 (1-1)，
用矩阵分解方法求解式(1-1)的特征值，并对特征值进行筛选，筛选出与Lamb波特性对应的特征值，由筛选出的特征值计算出Lamb波的相速度。
2.根据权利要求1所述的一种基于薄层离散计算层状板中Lamb波频散曲线的方法，其特征在于：在所述S1中，所述虚拟介质层的剪切波速cbs＝0.1m/s，密度ρbs＝1.25kg/m3，泊松比vbs＝0.45，即以虚拟介质模拟空气，且所述虚拟介质层的厚度与层状板的厚度相同。
3.根据权利要求1或2所述的一种基于薄层离散计算层状板中Lamb波频散曲线的方法，其特征在于：设Lamb波频散曲线计算频率范围用符号[fmin，fmax]表示，其中，fmin及fmax分别代表频率范围的下限及上限，
所述S2具体为，
S21，由层状板中各层板材料的泊松比νs及剪切波速cs，按关系式
Figure FDA0003890626600000021
计算层状板中各层板材料的纵波波速cp，且层状板中各层板材料的纵波速cp中的最小值用符号cp，min表示；
S22，由层状板中各层板材料的纵波速cp中的最小值cp，min，以及频率范围[fmin，fmax]中最大频率fmax，计算纵波最小波长
λmin＝cp，min/fmax (3)；
S23，按照厚度h≤λmin/20对层状板和虚拟介质层构成的整体结构进行水平向离散。
4.根据权利要求1或2所述的一种基于薄层离散计算层状板中Lamb波频散曲线的方法，其特征在于：
由第n薄层的上、下交界面位移分量及应力分量分别构建该薄层的广义位移向量un及广义应力向量pn，其中
un＝[ux,n,ux,n+1,iuz,n,iuz,n+1]T (4)，
pn＝[-σxz,n,σxz,n+1,-iσz,n,iσz,n+1]T (5)，
其中，ux,n和uz,n分别表示第n薄层的上交界面水平向位移分量及竖直向位移分量，ux,n+1和uz,n+1分别表示第n薄层的下交界面水平向位移分量及竖直向位移分量，σxz,n和σz,n分别表示第n薄层的上交界面剪切应力及法向应力，σxz,n+1和σz,n+1分别表示第n薄层的下交界面剪切应力及法向应力，虚数
Figure FDA0003890626600000022
所述S3具体为，薄层内位移由该薄层的上、下交界面的位移线性插值计算，利用虚功原理得到薄层刚度矩阵，第n薄层的刚度矩阵Kn可表示为，Kn＝Ank2+Bnk+Cn (6)，
式中k是波数，矩阵An、Bn和Cn分别为，
Figure FDA0003890626600000031
Figure FDA0003890626600000032
Figure FDA0003890626600000033
其中，ω是角频率，hn为第n薄层的厚度，μn和λn均为第n薄层的拉梅常数，ρn为第n薄层的密度,矩阵Bxzn T是矩阵Bxzn的转折。
5.根据权利要求4所述的一种基于薄层离散计算层状板中Lamb波频散曲线的方法，其特征在于：
所述S4具体为，根据薄层交界面位移连续条件，将由式(6)表示的各薄层刚度矩阵集合成总刚度矩阵K，总刚度矩阵K可表示为
Figure FDA0003890626600000034
在式(7)中，矩阵
Figure FDA0003890626600000035
是式(6-1)中各薄层子矩阵Axn的集成；矩阵
Figure FDA0003890626600000036
是式(6-1)中各薄层子矩阵Azn的集成；矩阵
Figure FDA0003890626600000037
是式(6-2)中各薄层子矩阵Bxzn的集成；矩阵
Figure FDA0003890626600000038
是式(6-3)中各薄层子矩阵Cxn的集成；矩阵
Figure FDA0003890626600000039
是式(6-3)中各薄层子矩阵Czn的集成。
6.根据权利要求5所述的一种基于薄层离散计算层状板中Lamb波频散曲线的方法，其特征在于：
在所述S5中，薄层的交界面总位移矢量U用矢量
Figure FDA0003890626600000041
代替，
Figure FDA0003890626600000042
如下式
Figure FDA0003890626600000043
其中，矢量Φx和Φz分别为
Φx＝[ux,1,ux,2,…,ux,N-1,u,x,N]T (8-1)，
Φz＝i[uz,1,uz,2,…,uz,N-1,u,z,N]T (8-2)，
上标符号T表示转折，N为薄层的交界面总数，ux,n和uz,n分别表示第n薄层的上交界面的水平向和竖直向位移分量，n＝1,2,3...，N；
利用刚性基位移为零条件，式(1-1)可改写为
Figure FDA0003890626600000044
其中
Figure FDA0003890626600000045
为零矢量。
7.根据权利要求6所述的一种基于薄层离散计算层状板中Lamb波频散曲线的方法，其特征在于：在所述S6中，具体用矩阵分解法求取式(1-1-1)的特征值；
所述S6具体为，
S61，给定频率后，由式(6-1)～(6-3)计算各薄层子矩阵Axn，Azn，Bxzn，Cxn，Czn，将各薄层子矩阵分别集成得到矩阵
Figure FDA0003890626600000046
并利用数学工具软件Matlab中的函数eig计算出式(1-1-1)的特征值1/k2，其中，特征值1/k2有2N个，2N个特征值1/k2对应有4N个特征波数k，设km为4N个特征波数中的任一个，且m＝1，2，3，...，4N；
S62，筛选出实部大于零、虚部小于零且实部是虚部绝对值十倍以上的特征波数，设筛选出来的特征波数有M个，同时设kj为筛选出来的M个特征波数中的任一个，且j＝1,2,3，...，M；
S63，根据筛选出来的特征波数kj的实部按下式计算Lamb波的相速度cj
cj＝ω/Re(kj) (9),
这里，符号Re表示实部；
S64，将Lamb波的相速度cj由小到大排序，并重新编号为c1<c2<...<cM，这些值依次对应于该频率Lamb波不同阶次模态相速度；
S65，改变频率，重复以上步骤S61～S64，得到在[fmin,fmax]范围内不同频率Lamb波各模态相速度，进而得到Lamb波各模态频散曲线。
8.根据权利要求7所述的一种基于薄层离散计算层状板中Lamb波频散曲线的方法，其特征在于：在利用数学工具软件Matlab中的函数eig计算式(1-1-1)的特征值1/k2的同时，还得到特征向量
Figure FDA0003890626600000051
根据筛选出来的特征波数所对应的特征向量
Figure FDA0003890626600000052
中的位移分量矢量，得到Lamb波各模态竖直向位移沿厚度方向变化形状。
9.根据权利要求8所述的一种基于薄层离散计算层状板中Lamb波频散曲线的方法，其特征在于：在单层或以板中心平面为对称面的板中，由竖直向模态位移形状判断板中Lamb波模态为对称还是反对称类型。
CN201810720164.XA 2018-07-03 2018-07-03 一种基于薄层离散计算层状板中Lamb波频散曲线的方法 Active CN109117500B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810720164.XA CN109117500B (zh) 2018-07-03 2018-07-03 一种基于薄层离散计算层状板中Lamb波频散曲线的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810720164.XA CN109117500B (zh) 2018-07-03 2018-07-03 一种基于薄层离散计算层状板中Lamb波频散曲线的方法

Publications (2)

Publication Number Publication Date
CN109117500A CN109117500A (zh) 2019-01-01
CN109117500B true CN109117500B (zh) 2022-11-22

Family

ID=64822154

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810720164.XA Active CN109117500B (zh) 2018-07-03 2018-07-03 一种基于薄层离散计算层状板中Lamb波频散曲线的方法

Country Status (1)

Country Link
CN (1) CN109117500B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109752262B (zh) * 2019-01-18 2020-10-27 中国水利水电科学研究院 一种基于原位相对密度确定覆盖层土体动剪模量参数的方法
CN110795820B (zh) * 2019-09-24 2021-07-16 武汉大学 工程结构裂纹问题求解方法以及装置
CN111307940B (zh) * 2020-04-01 2023-04-07 东北电力大学 一种金属管道周向导波激励频率区间的确定方法
CN114088818B (zh) * 2021-11-16 2024-03-22 南京工业大学 一种识别全域刚度的超声导波方法及系统
CN115982522B (zh) * 2023-01-10 2023-09-26 西南交通大学 一种弯曲波在周期性加劲板结构中的传播带隙计算方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8355893B2 (en) * 2008-12-12 2013-01-15 Wisconsin Alumni Research Foundation Method and system for analysis and shape optimization of physical structures using a computerized algebraic dual representation implicit dimensional reduction
CN106680375B (zh) * 2016-11-25 2019-08-02 中国商用飞机有限责任公司 用于确定材料的弹性模量的空气耦合超声检测方法
CN107818209B (zh) * 2017-10-26 2021-02-02 哈尔滨工程大学 一种弹性板结构的振动分析方法

Also Published As

Publication number Publication date
CN109117500A (zh) 2019-01-01

Similar Documents

Publication Publication Date Title
CN109117500B (zh) 一种基于薄层离散计算层状板中Lamb波频散曲线的方法
Pandey et al. Analysis of functionally graded sandwich plates using a higher-order layerwise theory
Chandra et al. Vibro-acoustic response and sound transmission loss analysis of functionally graded plates
Ferreira et al. Analysis of composite plates using higher-order shear deformation theory and a finite point formulation based on the multiquadric radial basis function method
Awrejcewicz et al. Chaotic dynamics of flexible Euler-Bernoulli beams
Ferreira et al. Buckling analysis of isotropic and laminated plates by radial basis functions according to a higher-order shear deformation theory
Hajlaoui et al. Nonlinear dynamics analysis of FGM shell structures with a higher order shear strain enhanced solid-shell element
Cunfu et al. The propagation of coupled Lamb waves in multilayered arbitrary anisotropic composite laminates
Ma et al. A symplectic analytical wave based method for the wave propagation and steady state forced vibration of rectangular thin plates
Narwariya et al. Harmonic analysis of moderately thick symmetric cross-ply laminated composite plate using FEM
Vu et al. A refined sin hyperbolic shear deformation theory for sandwich FG plates by enhanced meshfree with new correlation function
Swaminathan et al. Vibration and stability characteristics of functionally graded sandwich plates with/without porosity subjected to localized edge loadings
Ding et al. Wave propagation in a periodic elastic-piezoelectric axial-bending coupled beam
Javed et al. Free vibration of cross-ply laminated plates based on higher-order shear deformation theory
Kapuria et al. AC 1‐continuous time domain spectral finite element for wave propagation analysis of Euler–Bernoulli beams
Pichler et al. A mode-based meta-model for the frequency response functions of uncertain structural systems
Balcı et al. Free vibration analysis of a laminated composite beam with various boundary conditions
Burlayenko A continuum shell element in layerwise models for free vibration analysis of FGM sandwich panels
Van et al. Static bending and free vibration analysis of functionally graded porous plates laid on elastic foundation using the meshless method
Shahid et al. Dynamic response of flexible viscoelastic kerf structures of freeform shapes
Arenas Matrix method for estimating the sound power radiated from a vibrating plate for noise control engineering applications
Onyeka et al. Study on stability analysis of rectangular plates section using a three-dimensional plate theory with polynomial function
Parrinello et al. Modal density of rectangular structures in a wide frequency range
Cao et al. A VAM-based equivalent model for random vibration of composite sandwich plate with arrowhead-on cores
Ozer et al. Development of an equivalent shell finite element for modelling damped multi-layered composite structures

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant