做超声检测仿真的人十有八九在第一次跑完COMSOL压电模型后看着满屏波场云图都会发懵说好的表面波呢怎么到处都是乱七八糟的体波和边界反射我从1.5MHz压电片激励、亚克力楔块转换表面波的模型里踩完一轮坑之后最大的感受是——这类仿真最容易翻车的不是软件操作而是物理参数和建模简化之间的匹配。这篇把我搭建COMSOL压电超声表面波裂纹检测模型的完整过程写下来包括1.5MHz频率和亚克力楔块角度怎么推出来的网格和时域步长怎么设才能保住瑞利波以及裂纹反射信号怎么从时域图里判出来。打算复现或者准备改参数做自己方案的可以直接顺着这个思路往下走。1. 为什么用楔块从“贴片”到“表面波”的物理逻辑1.1 表面波在钢板裂纹检测中的不可替代性先说一个反直觉的结论压电片直接贴在钢板上通电后也能让钢板振动但很难激发出可用于裂纹检测的瑞利波。真正有效的做法是把压电片贴在亚克力楔块的斜面上让纵波斜入射到钢表面利用波型转换把能量“交给”表面波。瑞利波是沿固-空界面传播的表面弹性波能量集中在约一个波长深度内。1.5MHz在钢中波长大约2mm也就是说能量基本锁在钢表面2mm左右的薄层里。这个特性对表面裂纹检测非常关键表面波可以沿着表面长距离传播扫查面积大碰到表面开口裂纹时一部分能量会被反射回来一部分透射过去还有一部分散射进钢板内部。反射回波在时域上非常清晰直接对应裂纹的位置。换成纵波或横波做同样的事就没这么舒服。体积波在钢板里会撞到底面再反射回来底面回波、裂纹回波、模式转换波混在一起时域信号一团糟。表面波把能量限定在表面附近等于天然滤掉了底面干扰。这也是为什么金属材料表面裂纹、焊缝表面缺陷、滚动接触疲劳裂纹这类问题常规超声检测方案里都有表面波的身影。1.2 1.5MHz的波长、裂纹尺度与灵敏度三角平衡频率选择背后是一个三角平衡频率高波长短对小裂纹灵敏但衰减大、对表面状态敏感频率低波长长衰减小但裂纹如果太浅反射信号根本出不来。钢中瑞利波速度大约在2910到2980m/s这个范围按2950m/s估算1.5MHz对应的波长约为1.97mm。经验上表面波对深度小于波长的四分之一到八分之一的开口裂纹基本不敏感可靠检测的裂纹深度大约在半个波长的量级。所以1.5MHz下的有效检测窗口在0.5到2mm深度的表面裂纹这正好覆盖了多数金属零部件表面萌生裂纹的尺度范围。频率钢中瑞利波波长可靠检测深度主要限制0.5MHz约6mm1.5mm以上浅裂纹不敏感1.5MHz约2mm0.5-2mm衰减中等常用折中2.5MHz约1.2mm0.3-1mm对表面粗糙度敏感衰减大另外1.5MHz还和压电片厚度有直接关系。PZT压电片工作在厚度伸缩模态时厚度大约是半波长量级1.5MHz对应1.2到1.5mm厚度的片子市面上很容易买到加工和粘接也都方便。你在仿真里如果随便填一个压电片厚度没有和激励频率匹配压电片的振动效率会低得离谱波形也会拖得很奇怪。1.3 亚克力楔块的临界角推演斯涅尔定律的工程实践楔块是整个装置里最容易被低估的零件。很多人觉得它就是把压电片垫斜至于斜多少度随便取个45度或者60度就行。结果是表面波成分很低模型里看到的大半是泄漏进钢板的体波。原理是这样的压电片贴在楔块斜面上发出的纵波在楔块里垂直于斜面方向传播斜入射到楔块与钢的交界面。入射角如果很小反射和折射都以纵波为主随着入射角增大折射波逐渐向掠射方向偏转。当入射角达到临界角时折射角等于90度折射波正好沿着表面方向传播这一方向的能量才可能激发瑞利波。临界角由斯涅尔定律决定sinθc / sin90° c楔块 / c钢表面。这里c楔块取亚克力中的纵波速度约2730m/sc钢表面取钢的瑞利波速度约2950到2980m/s。算下来sinθc 2730 / 2950 ≈ 0.925θc ≈ 67.7度也就是说亚克力楔块的斜面角要接近68度表面波的激发效率才最高。我仿真里习惯取68度起步然后做参数化扫描把65、67.5、70、72.5度各跑一遍观察钢表面法向位移的幅值变化。实测下来角度差个5度表面波幅值能掉两三成效果非常直观。至于为什么用亚克力不是随便选的。亚克力也就是PMMA纵波速度2730m/s正好把临界角推到一个方便加工的斜角范围密度1180kg/m³、声阻抗约3.2MRayl和钢的46MRayl比虽然差很多但楔块本来就不靠声阻抗匹配来传波而是靠入射角把纵波转换成表面波。再加上好加工、便宜、成型精度高自然成了首选。2. 模型搭建从压电片到钢板的几何与材料分配2.1 二维平面应变建模的合理性边界首版模型用二维平面应变来做是最划算的选择。表面波沿钢板长度方向传播数值解关注的也是这个方向的波包运动和裂纹反射和垂直于截面的方向基本无关。只要压电片宽度方向尺寸远大于波长比如片宽10mm而波长2mm平面应变假设就成立二维仿真的计算量比三维低一到两个数量级。参数数值说明钢板长度L_plate150mm端部回波不污染前40μs钢板厚度h_plate10mm厚度远大于波长半空间近似压电片长度L_pzt5mm沿声传播方向的接触长度压电片厚度t_pzt1.2mm1.5MHz厚度伸缩模态楔块高度h_wedge10mm楔块内部混响可控楔块角度theta_wedge68°亚克力/钢瑞利波临界角裂纹深度crack_depth1mm约半个波长裂纹宽度crack_width0.1mm远小于波长裂纹位置x_crack距楔块前沿再往前20mm处便于分离时域信号几何上就三个域压电片一个矩形楔块一个直角梯形钢板一个矩形。压电片的一条边和楔块斜边共线楔块底边和钢板上表面共线。建模时用Form Union把三个域合成一个联合体共享边界自动合并波在界面处自然连续。不需要额外设置接触或者配对这一点后面还会细说。2.2 PZT压电片的极化方向与材料坐标系压电材料用COMSOL材料库里的PZT-5H没问题但一定要检查极化方向。PZT-5H在材料库里的默认极化轴是Z轴而二维模型通常把截面放在X-Y平面。如果你没做任何处理几何里压电片厚度方向是Y材料坐标系里的极化轴却是Z电场激励方向与极化方向完全错位压电片根本不会按预期振动。解决办法是建一个旋转坐标系绕某个轴旋转90度把Z轴转到压电片厚度方向上。然后在压电片域里指定弹性、耦合、介电矩阵都使用这个旋转坐标系的材料。很多人模型跑出来位移量级不对、波形起振很怪多半就是漏了这一步。工作模态方面1.5MHz激励下压电片工作在厚度d33模式电场沿厚度方向应变也沿厚度方向。1.2mm厚度接近PZT-5H在1.5MHz下的半波长共振厚度振动效率高。要注意的是二维平面应变会把平面外也就是压电片宽度方向的应变约束为零这会略微抬高压电片在厚度方向的视在刚度。如果你的目的是算导纳曲线要做实验对标建议改用三维模型如果只是看表面波传播和裂纹回波二维的结果定性完全够用。2.3 亚克力楔块与钢板的声学参数设置楔块和钢板都用各向同性线弹性材料。PMMA的关键参数密度1180kg/m³杨氏模量约3.2GPa泊松比约0.37。钢的参数密度7850kg/m³杨氏模量210GPa泊松比约0.29。材料参数填进去之后可以顺手验算一下波速是否合理。从杨氏模量和泊松比推出来的钢中纵波速度应该在5920m/s左右横波在3220m/s左右瑞利波落在2950m/s上下。如果你填的参数算出来纵波只有4000多那要回头检查单位或者弹性模量是否填错了。这类模型里单位错了波速全部漂移裂纹回波时间对不上排查起来相当痛苦。阻尼建议先不要急着加。我第一版模型加了约等于零的材料阻尼等波速、反射时间窗都验证对了再往钢域里补很小的阻尼调整幅度衰减。如果一上来就加瑞利阻尼波形幅度和相位都会变得奇怪到时候你根本分不清是物理问题还是阻尼参数问题。3. 物理场耦合和激励加载的实操细节3.1 压电效应的耦合方式选择物理场推荐用固体力学加静电再接一个压电效应多物理场耦合。域分配上压电片用压电材料本构楔块和钢用线弹性本构静电接口只作用在压电片域。这个组合最干净也最容易控制每一步求解的内容。COMSOL里也有预置的压电设备接口本质上就是把固体力学、静电和压电耦合打包在一起。新手直接选预置接口也行但对后续想加非线性或者特殊边界会不太方便。我习惯手动搭原因很简单——后处理里我想单独查看哪个物理量就在哪个接口下面找不会因为预置接口的默认设置产生额外约束。这里有一个非常关键的界面设置楔块底面和钢板上表面之间波要能无损传过去。如果用Form Assembly装配需要在边界对上设置连续性条件让两侧的位移和面力连续如果像我前面说的用Form Union几何上就是一个连续体材料不同但位移天然连续。很多人在Assembly模式下忘了配对或者设置成了接触条件结果表面波压根建立不起来。钢板的其余表面按自由边界处理。现实中钢板侧表面暴露在空气里声阻抗差很大超声波在表面反射所以自由边界反而是最接近实测的设定。3.2 突发激励信号设计窗函数、周期数与频带宽度激励表不要直接写正弦函数。如果从t0突然开始加正弦电压压电片会受到一个阶跃冲击激发出一大堆宽频伪波频带远远超出1.5MHz网格根本响应不了那些短波长分量结果就是波形拖尾、边界反射乱成一团。我用的激励是5周期汉宁窗调制的tone burst表达式是这样V V0 * sin(2*pi*fc*t) * 0.5 * (1 - cos(2*pi*fc*t/n))其中fc1.5MHzn5且只让t在0到n/fc范围内生效也就是3.33微秒的突发脉冲。汉宁窗的作用是让幅值从0平滑上升到峰值再平滑回到0避免波包两端出现陡峭的截止边沿。周期数n的选择也有讲究n越小波包越短时域分辨率越好但频带越宽n越大频带越窄网格压力小但波包太长直达波和裂纹回波容易重叠。5周期是两者之间比较稳的折中。在压电片的上表面用终端边界加载电压型激励幅值V0取10V到50V都行线性小信号下波形形状不随幅值变化只影响绝对幅值下表面设置为接地。上表面终端一定要只选压电片的边界不要扩展到楔块边否则几何角点上会出现奇怪的电荷分布。3.3 电极、接地与电荷奇异点的处理电极在二维模型里就是一条完整边界不需要单独画成几何域。真正会出问题的是角点处的电荷奇异压电片边界两端如果刚好是直角尖角终端电压加载后角点区域的电荷密度会发散给导纳计算和个别点的应力带来无物理意义的高峰。处理办法有两个。一是在几何上给压电片角点做0.2mm左右的倒圆角网格在角点处加密一点二是只在后处理里避开尖角的极近邻区域。倒圆角成本很低我建议直接做。另一个实操细节是初始值。瞬态研究默认从零状态开始汉宁窗激励正好从0缓慢上升不会产生阶跃。如果你用了普通正弦或者方波型激励那就需要在激励函数里加一个上升沿不然前几个时间步的加速度冲击会把网格搅出严重噪声。所有设置检查完之后进入研究前先跑一个5微秒的短试算看一眼压电片附近的波场是否对称建立。如果这个阶段就不正常不要急着跑到40微秒省下来的时间够你排查三四遍。4. 网格和求解器超声仿真最容易翻车的环节4.1 每波长网格数的硬约束超声波传播仿真的第一条规矩就是网格尺寸必须能解析最短波长。表面波在钢中的最小相速度约2.9mm/μs1.5MHz对应波长约1.95mmPMMA中纵波速度更低波长约1.8mm。整块模型里要照顾的最短波长就是1.8到1.95mm。二阶单元每个波长至少要6到8个网格一阶单元要求更高要10到12个。按二阶三角形单元来算λ/8就是0.24mm所以钢和楔块区域最大网格取0.25mm。但表面波的能量集中在表面以下一个波长范围也就是最表层0到2mm。我会在钢板表面附近单独定义一层厚约2mm的细化区网格尺寸压到0.15mm裂纹是0.1mm宽的窄槽附近网格要再细化到0.05到0.08mm。区域推荐网格尺寸原因钢板整体0.25mmλ_R/8二阶单元保底钢板表层2mm0.15mm表面波能量集中区裂纹及其附近0.05-0.08mm保证窄槽几何和散射精度楔块0.2-0.3mmPMMA波长约1.8mm网格尺寸的过渡不能跳变。0.15mm的区域不要直接接一个1mm的粗网格区中间留几层逐渐变大的过渡带否则会产生人工反射界面。这个道理跟光在不同介质里折射一样网格尺寸突变就是声学阻抗突变波会被“假裂纹”反射。4.2 时域步长的CFL条件与BDF稳定性时域步长由模型中传播最快的波决定钢中纵波约5920m/s远超瑞利波。显式时间步进要满足CFL条件粗算一下0.2mm的网格对应步长上限大约0.2mm除以5920m/s结果是33.8ns。网格细化到0.15mm之后步长还要更小所以我一般直接取20到30ns。如果用COMSOL默认的BDF隐式求解器理论上可以放宽到50ns左右但每一步迭代求解的成本也不低总账算下来并不比显式便宜。我习惯首版先设Δt25ns直接跑40微秒对应1600步二维模型就算网格有3万单元BDF求解压力也不大。COMSOL 6.4版本里加了弹性波时域显式接口专门针对这种瞬态波传播场景大模型的内存占用和求解效率都更好。如果你的模型后期要转成三维强烈建议从6.4这个接口入手。不过压电激励部分还是要靠固体力学和静电耦合跨接口耦合复杂度会上去二维这个规模用BDF反而省心。求解器的相对容差默认1e-3基本够用。如果关注的是反射波包的高频细节可以放宽到1e-4计算时间增加得不多但绝对不要一上来就调1e-6白白浪费算力。4.3 边界反射抑制PML与低反射边界钢板的左右两个端面如果没有处理表面波传到端面会猛烈反射回来在时域信号里形成和裂纹回波竞争的大幅波包。最简单的处理办法是加长模型尺寸让端面回波晚于关注时间窗。按表面波速度2.95mm/μs算150mm长的钢板端面反射波要走约100微秒才能回到楔块附近而我们的关注时间窗取前40微秒端面回波压根不会进窗。这也是为什么我把钢板长度设为150mm而不是更短的30mm。模型加长带来的网格增加完全可以接受二维求解很快。如果想进一步压缩模型尺寸可以用PML完美匹配层。PML厚度取一个波长左右也就是2到3mm。要注意的是PML在COMSOL里对坐标方向有严格要求二维平面应变下如果坐标没设对PML不但不吸收还会产生虚假反射。低反射边界条件看起来省事但对斜入射波并不是完全吸收体波斜着打到边界上残留误差比较大我不建议作为首选项。5. 裂纹响应的提取与后处理验证5.1 裂纹建模方式的选择表面裂纹最直接的建模方式是开一个矩形窄槽宽度0.1mm深度1mm短边比波长的十分之一还小。实际裂纹通常比这更窄更尖但超声散射行为的主导参数是裂纹深度和开口面积0.1mm和0.05mm宽度算出来的反射系数差异很小而对网格的影响却非常大。这个0.1mm宽的槽在COMSOL里就是一个退化的矩形域网格细化到0.05mm可以完整解析。槽内部不填材料对波来说就是空腔。表面波传到这里垂直位移迅速变化一部分能量反射一部分透射和真实表面开口裂纹的散射规律是一致的。如果你想模拟更理想的零宽度裂纹那就不能靠挖槽得把几何分裂成两个域并用装配让裂纹两侧的网格节点不连续。这个操作在二维里也能做但网格和配对设置麻烦不少。我建议首版先用窄槽方案把时域信号跑出来确认裂纹回波时间窗没问题再考虑要不要做理想裂纹。5.2 时域信号的判读直达波、裂纹反射与底面混响计算完成后在钢板表面、距离楔块入射中心50mm的位置设置一个探针记录表面法向位移分量uy随时间的变化。这个点既是接收点也是你观察裂纹回波的“站台”。根据表面波速度2.95mm/μs接收点距离入射点50mm直达表面波大约在17微秒到达。裂纹开在我设定的楔块前沿再往前20mm处也就是距离接收点20mm反射波要从入射点出发越过接收点打到裂纹再返回接收点总路程是50加20再加20等于90mm所以裂纹反射波大约在30.5微秒到达。这个时间窗和直达波相隔13.5微秒波形完全分得开。拿到时域曲线后先找直达波的最大幅值A0再在裂纹回波时间窗内找峰值Acrack两者比值就是反射系数。1mm深、1.5MHz下的表面波反射系数经验值在0.2到0.4之间如果你仿真算出来是这个量级说明建模基本合理。如果裂纹回波只有直达波的百分之几多半是网格太粗或者楔块角度偏了表面波没有真正建立起来。5.3 用速度场和椭圆极化验证表面波模式光看时间曲线还不够必须确认模型里传播的确实是从楔块耦合出来的瑞利波而不是别的东西。判断标准有三条波速约等于2.9mm/μs而不是5.9mm/μs的纵波表面法向位移幅值明显大于切向位移幅值波场能量集中在一个波长深度内随深度快速衰减。在COMSOL后处理里把表面法向位移云图和位移箭头叠加显示开一个2D动画。如果见到波前从楔块接触区沿着表面扩散表面质点的运动轨迹呈椭圆形状法向分量大于切向分量那就是经典的瑞利波特征。为了量化检查我习惯在距离表面不同深度处画竖向截线看位移幅值随深度的分布。瑞利波在1个波长深度处幅值应该衰减到表面的三分之一以内如果传了很深还有明显幅值说明波型已经转成体波或者板波了。钢板厚度10mm、频率1.5MHz时板内还会出现兰姆波成分特别是时间窗拉长之后板面的上下表面来回反射会让波形越来越复杂。所以判读时只关注前几微秒的主波包和第一道裂纹回波后面的密集震荡不要过度解读。6. 实测调参心得从仿真到与实验对标的几个坑6.1 材料阻尼对波幅衰减的影响仿真里如果材料阻尼设为零表面波传播几十毫米后幅值几乎不衰减这跟真实钢有差距。真实钢材有内部摩擦和晶粒散射波传播过程会有一定损耗。但问题在于很多人在模型里加的瑞利阻尼系数都是拍脑袋填的结果就是把裂纹回波吃得一干二净。我建议的流程是先零阻尼跑通把所有时间窗、波形判读都确定好然后给钢域加等效损耗因子0.02左右的阻尼看波幅衰减比例最后如果手头有实验数据再微调阻尼系数直到仿真回波幅度和实测接近。瑞利阻尼的参数换算比较绕β系数通常在1e-9到1e-8这个量级具体跟频率和材料有关。关键是不要让阻尼参数主导你的物理结论。6.2 耦合剂与被测表面到底用连续性还是接触这个坑我见得太多了。有人在楔块底面和钢板上表面之间设了接触条件结果波形完全出不来或者出来一堆高频振荡然后开始怀疑是网格不够细。其实问题就出在接触模型上。真实检测中楔块和钢之间有一层很薄的耦合剂作用是让超声波穿过去。理想情况下这层耦合剂极薄声学上等效为界面两侧位移连续、应力连续。所以在仿真里这层界面就该用“连续性”边界条件而不是“接触”。接触条件是要模拟两个物体之间可能分开、可能滑移的力学状态它带来的接触刚度和摩擦行为对超声传播完全是噪声。如果你用Form Union把楔块和钢建成联合体界面自动连续根本不涉及这个问题。只有用Form Assembly时才需要手动指定边界配对为连续性。3D模型如果从CAD导入默认走了Assembly这一步特别容易踩雷。6.3 压电片尺寸与楔块接触区域的次生波形干扰最后一个容易被忽略的干扰源是压电片和楔块自身产生的寄生波。压电片宽度如果太大除了厚度方向的d33主振动还会激发横向谐振模式匀波包上叠加高频纹波。在1.5MHz下压电片长度取5mm我做过参数扫描5mm相比3mm和8mm主波包最干净。太短的片子激励能量不够太长的片子次生波形明显变多。楔块内部反射也会产生固定间隔的尾随波包。激励纵波在楔块里发出后一部分能量会从楔块的自由面反射回来再次耦合进钢板表现为主直达波后面的一个滞后波包。如果这个混响波恰好落在裂纹回波时间窗里判读就会出问题。我解决这个问题用了一个很土但很有效的办法先跑一个无裂纹模型把完整时域信号保存为基线再跑有裂纹模型把两条曲线相减。相减之后楔块混响、边界反射、激励泄漏等系统性干扰全部被剥掉剩下的就是纯裂纹散射贡献。这个基线差分法在实验超声里是常规操作放到仿真里同样好用强烈建议做裂纹定量分析的朋友试一试。模型跑通之后我习惯把楔块角度重新做一遍参数化扫描观察不同角度下钢表面法向位移的幅值变化。这个扫描能同时验证物理逻辑和数值实现波形最干净、幅值最大的角度通常正好落在用斯涅尔定律推出的临界角附近。如果仿真和理论推出来的最优角差很多那就要回头检查材料参数或者物理场设置哪里出了问题。做这类超声仿真最重要的一步永远不是软件操作而是搞清楚波到底要在哪里形成、以什么模式传播先把这个物理图景想明白模型只是水到渠成的事。