去年年底接了个活儿，要在COMSOL里把一个相控阵探头的三维声场算明白。说实话，相控阵声场模拟这活儿，说难不难，说简单也不简单。难的地方在于三维模型一跑起来计算量蹭蹭涨，内存和收敛问题轮着来；简单的地方在于，核心物理其实就是一个叠加原理，每个阵元发出去的波在空间里互相干涉，算清楚相位差，声压分布就出来了。真正折磨人的反而是那些不起眼的细节：PML（完美匹配层）厚度够不够、网格在波长尺度上划得够不够细、每个阵元的相位激励怎么快速切换、后处理怎么把三维声压分布直观地提取出来。
最近这阵子在COMSOL里折腾三维声压分布，把几个挺有意思的操作技巧摸出来了。这篇文章不打算讲教科书里那套推导，纯粹是从实操角度记录下来，把网格划分、边界条件、相位计算、后处理这些环节里踩过的坑和验证过好用的方法，一并写出来给大伙儿参考。涉及的工具是COMSOL Multiphysics的压力声学模块，模型是水下或者固体耦合场景里常见的相控阵列，核心就两个字：叠加。下面正式开始。
1. 相控阵声场模拟的核心原理与COMSOL实现路径
1.1 相控阵波束成形的物理原理
相控阵这个名字听起来高大上，拆开来看并不复杂。它的本质就是让排成一排或一个面的一堆小阵元，各自发出同频率的声波，但每个阵元的发射时间或者相位不一样。这些波在空间里遇到一起，就会发生干涉。相位对齐的方向，波峰叠加波峰，声压增强；相位对不上的方向，波峰撞波谷，声压被抵消。于是整个阵列的声能量就被“捏”成一束，指向你想要的方向。这就是波束成形。
在仿真里实现这件事，只需在激励条件上做文章。对于连续波（CW）激励，给每个阵元一个不同的初始相位，等价于在频域里给每个阵元入口设置一个复数形式的振动速度或者压力边界条件。比如我要让声束聚焦到空间某个点，那么每个阵元到焦点距离不一样，声波传播过去需要的时间不一样，所以我反过来在每个阵元上施加一个补偿相位，让所有阵元发出的波到达焦点时刚好同相。这个补偿相位的计算公式很简单：
[ \phi_i = \frac{2\pi}{\lambda} \left( r_i - r_0 \right) ]
其中 ( r_i ) 是第 i 个阵元到焦点的距离，( r_0 ) 是参考距离（通常取阵列中心到焦点的距离），( \lambda ) 是介质中的波长。算出来的 ( \phi_i ) 就是第 i 个阵元的激励相位偏移，单位是弧度。把这个相位套进 COMSOL 的边界条件里，求解器解出来的声压场自然就带上了聚焦效果。
这里面有个容易忽略的地方：相位差计算用的是阵元中心到焦点的距离还是阵元表面每个网格点到焦点的距离？实际操作中，大多数情况下用阵元中心到焦点的距离做单值激励就够了，因为单个阵元尺寸通常不超过半个波长，表面各点到焦点距离差带来的相位差很小。但如果你用的是大尺寸阵元或者阵元工作时长比较极端，那就得在边界条件里写空间相关的相位表达式，不能偷懒用常数。
1.2 COMSOL 中实现相控阵的模块选择与整体思路
COMSOL 里做声场模拟，第一个问题是选哪个物理场接口。纯声学问题，比如不考虑结构振动、不关心阵元压电换能器内部的机电耦合，直接用“压力声学-频域”接口就够了。它求解的是亥姆霍兹方程，未知量是复声压，速度和声压的关系通过声学本构方程给出。计算量适中，是做相控阵声场的主力接口。
如果你的模型里除了声场还要考虑阵列面板的结构振动，那就得用“声-结构耦合”或者“压电-声耦合”。压电换能器阵列的完整模型会用到 COMSOL 的压电器件接口，这个接口能同时算压电陶瓷的逆压电效应和周围的声场，计算量比纯声学接口大很多。我做相控阵声场模拟的习惯是：先搞一个纯声学的快速模型，把相位关系、聚焦效果、栅瓣问题验证清楚，再决定要不要上完整压电耦合。
三维模型相比二维模型优势很明显，能真实反映阵列的几何布局和旁瓣分布，但代价是计算资源。二维模型一个网格数量可能就几万，三维模型轻松破百万。COMSOL 多物理场耦合的框架在这里反而帮了忙，你可以先用二维模型快速扫参数，锁定几个关键参数之后再用三维模型精细验证。模型建好之后整体思路按下面几条走：
- 几何：建立阵列阵元、流体域，外圈设置 PML 吸收层。
- 材料：流体域给水的声学属性，密度和声速；阵元结构按实际材料给，但纯声学模型里阵元表面只给边界条件，不需要建实体。
- 边界：阵元表面设置法向加速度或振动速度边界，其余固体边界设硬声场边界，计算域外围设 PML。
- 网格：控制最大单元尺寸与波长之比，通常一个波长内至少 6 个单元，推荐 8 到 10 个。
- 求解：频域直接求解器 + PARDISO，如果做相位扫描就用参数化扫描。
- 后处理：三维切片、等值面、声压级提取、轴向声压分布曲线。
2. 几何建模、材料与边界条件的实操细节
2.1 阵元几何的建立与参数化设计
COMSOL 的几何建模功能在同类软件里算好用的，但做阵列这种大量重复结构时，关键不是画图，而是学会用参数化阵列。先定义全局参数，阵元间距 ( d )，阵元边长 ( a )，阵列行数 ( M )，列数 ( N )。然后用工作平面里的矩形阵列功能，一次性把所有阵元画出来，最后用“形成联合体”做布尔操作，方便后续给每个阵元单独设边界条件。
给每个阵元单独设边界条件这件事，是新手最容易卡住的地方。如果用“形成联合体”，所有阵元表面会变成一个整体边界对象，你没法单独选中某一个阵元施加不同相位。解决办法有两个：一个是建立几何的时候选择“形成装配体”，让每个阵元表面保留独立边界；另一个是直接用“显式选择”按空间坐标框选特定边界。
我的实际经验是：形成联合体时 COMSOL 会自动合并重合边界，适合所有阵元激励相同的简单场景；形成装配体则保留了每个阵元表面的独立边界，适合作相位扫描。如果做的是大规模阵列，后面还要批量设置边界条件，最好用“形成装配体”加“显式选择”。
又在装配体模式下，流体域和阵元之间是接触边界而不是联合边界，需要在多物理场里手动添加“声-结构边界”之类的耦合条件。纯声学模型不需要这个步骤，因为边界条件直接在流体边界上设置即可。不过要注意，装配体模式下网格划分时界面处网格要不要对应一致，如果只关心声场，可以开“非一致网格”选项，界面两侧各自剖分，求解器会自动处理映射关系。
2.2 材料参数与边界条件
材料参数看起来简单，但三维模型里一个微小错误就能让结果偏差巨大。声学模拟里最核心的两个参数是声速 ( c ) 和密度 ( \rho )。水的声速大约 1480 m/s，密度 1000 kg/m³，但注意这是 20°C 左右的数据。如果模拟的是人体组织、工业油液或者其他液体，务必去查对应温度或压力下的准确值。因为波长 ( \lambda = c / f )，声速差 1%，波长就差 1%，网格尺寸和相位补偿全部跟着偏。
边界条件的设置直接影响物理正确性。阵元表面我一般选“法向加速度”边界，给定一个复数幅值来表示相位：
[ a_n = A_i e^{j\phi_i} ]
其中 ( A_i ) 是第 i 个阵元的加速度幅值，( \phi_i ) 是相位。在 COMSOL 的边界条件面板里，直接在“法向加速度”输入框内写表达式，比如A0*exp(j*phi_i)，其中phi_i是全局参数或变量。如果做相位扫描，把phi_i定义成参数扫描的参数名，一次就能扫完整阵列的所有相位组合。
除了阵元表面，计算域最外圈如果不做处理，声波会被硬边界反射回来，干涉条纹就彻底乱了。这一点极其重要。传统的做法是给外圈设“阻抗边界”来近似吸收，但效果远不如 PML 来得干净。PML 必须紧贴计算域外边界，厚度至少取半个波长，推荐取 1 到 2 个波长。PML 内部网格要扫掠或者映射，不能随便自由剖分，不然吸收效果打折扣。
2.3 PML 完美匹配层的使用要点
PML 用得好不好，直接决定三维声场算得准不准。它的工作原理是在边界区域人为构造一个吸收媒质，声波进入后按指数衰减，反射率非常低。COMSOL 里添加 PML 的操作不复杂，但有几个要点必须注意。
第一次跑三维模型的时候，我直接用自由四面体网格把整个域剖了一遍，包括 PML 区域，结果远场看起来总有奇怪的小波纹。后来才发现 PML 层里必须用扫掠网格生成，至少 5 层单元，并且每层厚度要均匀。这相当于给 PML 设置了“梯度折射率”，如果单元大小忽大忽小，吸收效果就不稳定。
PML 厚度的选择也很讲究。厚度太薄，低频段吸收不干净；厚度太厚，网格量猛增。经验值是 PML 厚度取计算域最大频率对应波长的 1 倍，单元层数取 8 到 10 层。模拟水中 500 kHz 的信号时，波长大概 3 mm，PML 厚度 3 mm 到 6 mm 就够了。如果遇到低频且尺寸受限的情况，PML 至少也不能低于半个波长，否则从边界面反射回来的波会让干涉图样面目全非。
PML 外表面默认是硬边界还是软边界，这个也要检查。COMSOL 的 PML 功能通常会把外层设置成硬声场边界，但其实最外层应该是无反射边界。我在实际使用中发现，PML 区域内部已经把波衰减到很低的水平，外层边界啥类型影响不大。不过如果 PML 厚度不够，外层边界条件就会真实地反射出残余波，这就只能靠加厚 PML 来解决。
3. 网格划分与求解设置的工程经验
3.1 网格尺寸与频率的关系
声学模拟的精度几乎被网格尺寸决定。标准的经验准则是：一个波长内至少要 6 个单元，也就是最大单元尺寸不能超过波长的 1/6。我做声场模拟的习惯更保守，最大单元尺寸取波长的 1/8 到 1/10。单元越细，色散误差越小，声束方向和旁瓣位置越准。
举个例子。水中 500 kHz 信号，声速 1480 m/s，波长约 3 mm。1/6 波长就是 0.5 mm，1/10 波长就是 0.3 mm。对于三维模型，0.3 mm 的单元尺寸意味着每立方厘米流体域要剖出几十万到上百万个单元，计算压力不小。这时可以采用区域化网格：阵元附近的近场区域加密到 1/10 波长，远离阵列的传播区域放宽到 1/6 波长，PML 区域用扫掠网格独立剖分。这样既保精度又控规模。
怎么验证网格够不够细？最简单的办法是加密网格跑同一个模型，对比关键位置（比如焦点的声压幅值）的变化。如果两次结果相差小于 1%，说明网格已经收敛了。我经常把网格加密前后的声压幅值变化当成一个算例的自检指标，如果不达标，宁可多跑几小时也要把网格重划，否则后面所有云图都是自欺欺人。
3.2 求解器类型选择与迭代设置
COMSOL 做三维频域声学，默认的求解器通常是直接求解器。三维问题自由度动辄几百万，直接求解器的内存消耗非常大。这里有个取舍：直接求解器（如 PARDISO）鲁棒性强，几乎不用调参数就能收敛；迭代求解器（如 GMRES）内存占用小很多，但需要搭配一个好的预条件子，否则可能不收敛或者收敛很慢。
我个人的习惯是：自由度在 200 万以内直接用 PARDISO，省心；超过 300 万就考虑迭代求解器。但迭代求解器的预条件子设置比较讲究，声学问题里常用“左侧预条件 + 多网格”的组合。实际跑下来，对于纯亥姆霍兹方程，多网格预条件加上适当的平滑器能获得接近线性的收敛速度。
相位参数扫描是相控阵模拟的常见需求，扫描变量就是各阵元的相位组合。COMSOL 的参数化扫描功能会让每一个参数组合都重新组装一遍刚度矩阵，如果扫描步数多，计算时间就是线性叠加的。我一般在扫描之前先关掉“存储解”以外不必要的输出，只保存按需提取的声压数据，大幅减少磁盘占用。还可以利用扫描时“使用上一步解作为初值”的选项，让相邻相位组合的解相互接力，加快迭代收敛。
下面是一个简化的 COMSOL 相位扫描案例代码，展示如何在全局变量里定义阵元相位并在扫描时切换。
// 全局参数定义 // f = 500[kHz], c = 1480[m/s], lambda = c/f // N = 8 (阵元数量) // phi_focus 由扫描参数 slot 控制 // 变量定义部分 // phi_i = 2*pi*(sqrt((x_i)^2 + (y_i)^2 + F^2) - F)/lambda实际建模时，相位表达式要写进每个阵元的法向加速度边界条件里，形如A0*exp(j*phi_i)。在 COMSOL 里可以直接引用全局参数和变量，不需要额外写程序。
3.3 三维模型跑不动的降阶思路
三维模型跑不动是个高频问题，遇到这种情况不要硬扛，先想想能不能降维。很多阵列几何有对称性，时间上可以用 1/2 或 1/4 对称模型，只需要在对称面上加对称边界条件。比如方形阵列四象限对称时，只建 1/4 模型，网格量直接降到原来的 1/4，计算时间大幅缩短。
如果对称性用不上，还可以考虑用二维轴对称模型代替部分三维验证。圆形阵列或环形阵列具有旋转对称性，只要焦点在对称轴上，声场就是轴对称的。用二维轴对称模型直接模拟三维声场在轴上的分布，计算量小几个数量级。不过焦点在偏离轴线位置时，轴对称条件不成立，这时候就老老实实回到三维。
4. 三维声压分布的后处理与可视化
4.1 声压数据提取与声压级计算
求解完成后，COMSOL 里默认得到的变量是复声压 ( p )。要得到声压幅值，需要对复声压取模：abs(p)。但工程上更常用的是声压级（SPL），单位 dB，计算公式是：
[ SPL = 20 \log_{10} \frac{p_{rms}}{p_{ref}} ]
在水中通常取参考声压 ( p_{ref} = 1 \mu Pa )，空气中取 ( 20 \mu Pa )。在 COMSOL 后处理里可以新建一个表达式，直接写成20*log10(abs(p)/1e-6)，然后绘制这个表达式的体切片、等值面或三维云图，出来就是直观的声压级分布。
只画云图还不够，工程上更关注沿某个方向上的声压分布曲线。我从实践中总结了一套标准操作：在模型里定义一个截线或截点集，比如穿过焦点的轴向线段，然后用“一维绘图组”提取该线段上的声压幅值或声压级。这个操作能直接量出焦点位置、旁瓣位置、-6 dB 波束宽度，方便和实验数据对比。
提取“焦点增益”时需要注意，COMSOL 中某个点上的复声压是复数，直接取模即可。但如果模型里使用了 PML，焦点离 PML 太近会造成数据失真，因此计算域设置时保证焦点到 PML 的距离至少有 3 到 5 个波长。
4.2 三维切片、等值面与方向性图
三维声压分布最适合的可视化方式不是等值面，而是切片。沿着 x、y、z 轴各切一刀，两个正交平面和一个水平面组合起来，就能看出波束的形状和聚焦情况。切片颜色用预定义色表，把单位改成 dB，对比度会好很多。
如果想要更漂亮的“声束”效果，可以用等值面显示某个声压级阈值以上的区域。比如只显示 -6 dB 以上的区域，就能看出波束的椭球形状。这个可视化方法非常适合向非声学专业的人展示聚焦效果，一目了然。
方向性图（波束图）是相控阵模拟里必出的一张图。做法是在距阵列中心一定半径的球面上取一个圆周，提取该圆周上的声压幅值，然后以角度为横坐标、dB 值为纵坐标画极坐标图或直角坐标图。从这张图能清楚看到主瓣、栅瓣和旁瓣的位置。
5. 常见问题与排查技巧实录
5.1 迭代不收敛时先查什么
COMSOL 里报“求解器未收敛”，很多人的第一反应是调求解器设置。我的经验是，先别动求解器，先检查网格和 PML。最常见的原因是 PML 内部网格是自由剖分的，而不是扫掠网格，导致吸收效果差，边界把能量反射回来，数值震荡累积最终不收敛。把 PML 网格改成扫掠后，大多情况立刻就好了。
另一个高频原因是材料参数填错导致波长计算偏差。比如水里的声速写成空气中的 343 m/s，那波长就变成实际值的 1/4 左右，网格相对变粗，数值色散增加，长期迭代就会发散。所以不收敛时，先复核材料参数，再检查网格尺寸是否满足每波长至少 6 个单元的要求，最后才是动求解器。
如果上面都没问题，那可能是几何里存在非常小的特征尺寸，导致网格剖分时生成了劣质单元。这时在“大小”设置里把“狭窄区域分辨率”调高，或者手动指定最小单元尺寸，就能解决。
5.2 几何布尔操作报错与装配体模式
COMSOL 做阵列几何时经常报“转换为 CAD 内核时不支持的拓扑”，这通常是阵列单元之间的重复面导致联合体操作失败。解决办法很简单：改用装配体模式，或者先把阵列生成组件，再用“布置”功能摆放。装配体模式虽然接触边界多，但胜在稳定，后续加边界条件也灵活。
再说一个细节：装配体模式下阵元表面和流体域是独立网格，边界条件直接设在流体边界上。如果要模拟实际阵元振动，还要在边界上添加结构耦合；但纯声学模拟就完全不需要，给法向加速度即可。
5.3 内存不足的应对策略
三维声场模拟吃内存吃得厉害，尤其是用直接求解器。我遇到过几次内存溢出，后来总结了三个实用的减负方案：
第一，用对称模型缩减计算域。对称面设对称边界条件，求解自由度成倍下降。第二，用迭代求解器替代直接求解器，同时开多网格预条件，内存占用能降到一半以下。第三，把多余的几何细节删掉。比如阵列单元背后的支撑结构、倒角、圆角，这些几何细节对声场影响极小，但会让网格量暴涨，纯粹浪费资源。
还有一个冷门技巧：把模型的存储格式改成“稀疏矩阵”和“只有非零元素”，能省不少内存。COMSOL 默认会存很多中间变量，关掉“保留临时解”选项也可以腾出大量空间。
5.4 多阵元激励设置效率太低怎么办
阵列一多，每个阵元都要单独设一条边界条件和相位表达式，手动操作非常繁琐。COMSOL 支持右键点击“边界条件”下的“自动”功能，也可以借助“显式选择”按编号或坐标批量选中边界。还可以利用“全局定义”里的“变量”功能，一次性定义所有阵元的相位变量，然后每条边界只用引用对应的变量。
批量设置边界条件还有一个更高效的办法：用 COMSOL 的 Java 或 MATLAB API 脚本控制。比如对 64 阵元的阵列，写个循环给每个阵元设置法向加速度边界。这个操作对不熟悉编程的人有点门槛，但在阵元数量大、相位组合多的场景下，脚本一劳永逸。
6. 个人实操建议
最后再分享几个我踩过坑之后形成的习惯，不一定写进文档里，但对提高成功率特别有帮助。
第一，跑全尺寸三维模型之前，务必先用二维或二维轴对称模型把相位关系和焦点位置验证一遍。原因很简单，三维模型调试周期太长，一个边界条件写错可能要一天之后才能发现。二维模型几分钟就能出结果，先拿它把物理机制摸清楚，三维模型一次跑对的概率会高很多。
第二，相位表达式里用到的坐标参考点要固定。我遇到过一个很隐蔽的问题，使用全局坐标和局部坐标混用，导致阵列左侧和右侧阵元的相位符号相反，聚焦位置完全错位。后来我规定所有阵元的相位全部基于全局坐标系计算，就再也没出过类似问题。
第三，保存结果时尽量保存复声压而不是只存幅值。因为后处理阶段你可能还要提取相位信息，或者计算声强、方向性图。如果只存了幅值，这些后处理全没戏，需要重新跑一遍模型，白白浪费时间。
我在实际做相控阵声场模拟的过程中，感受最深的一点就是：这个项目几乎不会栽在物理公式上，反而容易栽在网格、边界和几何设置这些看起来不起眼的工程细节上。按照上面这套流程走一遍，多数情况下能顺利拿到干净准确的三维声压分布。希望这份实操笔记对正在折腾 COMSOL 相控阵模拟的你有帮助，也欢迎交流各自踩过的坑。
