在激光加工仿真圈子里,COMSOL 几乎是讨论度最高的工具。激光烧蚀、激光熔覆、选区激光熔化(SLM)这三类工艺,文献里能找到大量经典算例,但真正能把这些算例从论文里的参数表“复现”成可运行模型的人,并不多。这篇内容我把三条线放在一起讲,因为它们共享同一套多物理场耦合思维:热源怎么加载、材料状态怎么变化、边界怎么移动、网格怎么处理。理清这套框架之后,你会在 COMSOL 里搭激光加工模型这件事上有一个完整的全局认识,也能直接照着搭出能跑的简化版本。适合正在啃仿真文献的研究生,也适合刚接触这类仿真的工艺工程师。
1. 三大工艺的仿真共性与总体建模思路
1.1 三种工艺的物理主线其实不一样
很多人一开始就把激光加热当成一个“传热问题”来做,这没错,但对烧蚀、熔覆、SLM 来说,只做传热远远不够。三者宏观上的区别非常明显:激光烧蚀是材料在消失,激光熔覆是材料在增加,选区激光熔化是材料在粉末和实体之间切换并逐层累积。
这个区别直接决定了你在 COMSOL 里选择哪些物理场接口。烧蚀必须处理边界的后退运动,对应的是移动网格加传热;熔覆需要处理熔覆层的生长和熔池流动,对应的是变形几何、层流和传热的组合;SLM 除了移动热源之外,还要处理粉末向实体的状态转变,以及逐层铺粉带来的几何激活问题。如果一开始没把这条主线想清楚,后面搭模型很容易陷入“不知道怎么给边界条件”的困境。
| 工艺 | 核心物理过程 | 宏观表现 | COMSOL 必选接口 |
|---|---|---|---|
| 激光烧蚀 | 固相升温、汽化/熔化喷射 | 材料表面后退、形成凹坑 | 传热 + 移动网格(ALE) |
| 激光熔覆 | 基体熔池形成、粉末熔化沉积 | 熔覆层生长、稀释率变化 | 传热 + 变形几何 + 层流(可选) |
| 选区激光熔化 | 粉末选择性熔化、凝固、逐层堆积 | 熔道成形、致密化、缺陷形成 | 传热 + 域激活/材料切换 + 层流(可选) |
1.2 热源几何统一到高斯分布
经典文献里绝大多数激光热源都可以写成高斯分布,这是复现时的第一块基石。面热源公式为:
q(r) = q0 * exp(-2 * r² / r_beam²),其中 q0 = 2P / (π * r_beam²)
P 是激光功率,r_beam 是光斑半径,r 是当前点到光斑中心的距离。需要体热源时,再沿深度方向加指数衰减因子:Q(z) = q_surface * exp(-α * z),α 是材料的吸收系数。
COMSOL 里最常用的做法是在传热模块中加“热通量(Heat Flux)”边界条件,直接在表达式里写入上面的高斯分布。这里有一个极其容易踩的坑:r 必须动态追踪光斑中心。如果激光在某个方向上移动,光斑中心坐标要写成时间的函数,比如 r = sqrt((x - x_beam(t))² + (y - y_beam(t))²)。很多人第一次运行时发现等温线定在那里不动,基本都是因为把 r 写成了到某个固定点的距离。
1.3 多物理场之间怎么串联
温度场是绝对的主变量,流场只在涉及熔池时必要,固体力学场则要看你是否需要残余应力。第一次复现经典算例时,我的建议是先只开传热和热源移动,把温度场跑通、跑稳定,然后再逐步加流场、加移动网格、加相变。上来就把所有物理场全部打开,结果往往是模型还没算几步就因为不收敛或网格扭曲而中止,浪费大量时间去排查一个其实并不复杂的耦合错误。
耦合顺序有个朴素的经验法则:先算一个“纯传热 + 固定网格”的版本,确认热源、材料参数、边界条件给的量级与文献一致;再引入边界动网格或变形几何,观察几何变化;最后才加入流体流动和表面张力效应。每加一层物理场之前,都要先保证上一层结果是合理的。
2. 激光烧蚀复现:从单脉冲到动态形貌演化
2.1 先用二维轴对称搭几何,别直接上三维
单脉冲激光烧蚀在经典算例中非常常见,几何上通常具有较好的对称性,没必要一上来就建三维。二维轴对称模型能把计算量降低一个数量级,而且物理上完全足够。基体建一个矩形域,厚度取烧蚀深度的 3~5 倍,半径或半宽取光斑半径的 5~10 倍,这几点能避免边界反射影响温度场。
材料参数需要准备密度、热导率、比热容、熔化温度、蒸发温度、汽化潜热。纳秒脉冲激光烧蚀金属时,脉冲宽度只有十几到几百纳秒,热影响区小得多,所以激光光斑下的网格要加密到 0.1~0.2 μm 量级,远离光斑的区域逐渐粗化。注意:第一次跑透之前,不要急着用自适应网格,手动画一个合适的渐变网格更能控制变量。
2.2 烧蚀面后退用“指定网格速度”实现
材料表面因为汽化和熔体喷射而后退,这是烧蚀建模的核心。COMSOL 的移动网格(ALE)接口里,可以把边界速度设成由热输入决定的表达式。简化能量平衡下,表面后退速度可以写为:
v_n = (q_abs - q_loss) / [ρ * (c_p * ΔT + L_v)]
其中 q_abs 是表面吸收的热流密度,q_loss 是向体内热传导等损失项,ΔT 是材料从初始温度升到蒸发温度的温度差,L_v 是汽化潜热。分母本质上是单位质量材料从室温变成蒸气所需的总能量。
这个公式直观地反映了物理过程:热流密度足够大时,表面快速后退;热流密度低到不足以汽化材料时,速度自然会趋向于零,也就没有烧蚀发生。实际操作中,方向取边界内法向,并且建议在移动网格设置里打开“自动重新剖分(Automatic Remeshing)”,否则边界位移超过网格尺寸后单元会严重扭曲,求解器很可能直接报错。
2.3 脉冲参数换算与时间步控制
连续激光直接把功率填入热源公式即可,但脉冲激光需要先算能量密度。单脉冲能量密度公式:
F = E_p / A_spot = E_p / (π * r_beam²)
举个例子,单脉冲能量 0.5 mJ,光斑直径 50 μm,则 F = 0.5e-3 / (π * (25e-6)²) ≈ 25.5 J/cm²,这个量级明显高于常见金属的烧蚀阈值,表面会进入烧蚀状态。如果文献给的是平均功率、重复频率和脉冲宽度,还需要先用平均功率除以重复频率得到单脉冲能量,再换算成能量密度。
时间步长也要格外小心。脉冲宽度的 1/10 是最低要求,否则峰值热流会被时间平均抹平,算出来的烧蚀深度会严重偏低。用 COMSOL 时可以把热源写成高斯脉冲形式,例如 q(t) = q0 * exp(-((t - t0)/σ)²),σ 取脉宽的四分之一左右,让热流在脉冲时间内平滑上升和下降。这样既符合物理,也对数值求解更友好。
2.4 烧蚀深度、热影响区与文献对比
复现经典算例时,最常用的验证指标有三个:单脉冲烧蚀深度、烧蚀坑直径、热影响区厚度。烧蚀深度可以直接在后处理中读取边界位移量;热影响区则可以画出温度场中高于再结晶温度的等温线范围,和文献截面图作定性对比。
如果算出来的深度偏深或者偏浅,先检查两处:第一,能量密度有没有算对;第二,表面吸收率是否设置正确。金属对 1064 nm 波长激光的吸收率通常在 0.3~0.5,而不是 1。经典文献很多时候会把吸收率写在“光学参数”或者“材料参数”小节里,复现时很容易漏掉。我实测中最常见的偏差就是这里,因为大家一开始图省事直接给 100% 吸收,结果烧蚀深度和实验对不上。
3. 激光熔覆复现:移动热源与熔池生长
3.1 从二维截面模型入手,先看熔深稀释率
激光熔覆的完整三维熔池模拟非常吃算力,经典复现一般从二维截面模型或带厚度方向对称的三维窄条模型入手。二维模型能看到熔池形状、稀释率、热循环曲线,三维窄条模型(比如宽度取实际熔宽的一半)则可以额外看到表面形貌的走势。
推荐路线是先做二维:基体矩形加熔覆层区域,热源用移动高斯热流,覆层区域用变形几何持续添加,整个模型的输出核心指标是熔深、熔宽、稀释率。稀释率的定义是基体熔深除以熔深与熔覆层高度之和:
D = d_melt / (d_melt + h_clad)
稀释率太大会稀释合金元素,太小又说明冶金结合不充分,文献里通常要求控制在 5%~15% 范围。这个量是最容易和实验对比的,也是调工艺参数的天然靶子。
3.2 熔覆层生长的两种建模方式
第一种是表面生长法,也是最贴近物理的方式。用移动网格接口让熔覆层表面随时间抬高,高度增量速度可以近似写成:
h_dot = m_dot / (ρ * w * v)
其中 m_dot 是送粉速率,w 是熔覆道宽度,v 是扫描速度,ρ 是粉末材料密度。这种方式能刻画出熔覆层的自由表面形貌,但网格变形控制难度较高,需要配合自动重新剖分。
第二种是源项法,在熔覆层区域分别施加质量源和能量源。质量源对应送粉量,能量源对应粉末吸收的激光能量。源项法不需要动网格,数值稳定性明显更好,适合快速扫描参数看趋势。缺点是无法刻画真实的自由表面。我的建议是:方案评估阶段用源项法,最终复现和对比形貌时用表面生长法。
还需要特别注意激光能量的分配。送粉熔覆时,粉末在进入熔池之前就提前吸收了一部分激光能量,所以基体吸收的功率并不是总功率。常用近似分配:基体吸收 P*(1 - η_powder),粉末吸收 P*η_powder,η_powder 一般取 0.2~0.4。这个分配比例对熔深影响很大,可以当作调参的核心旋钮。
3.3 熔池对流做到什么程度
如果只关心温度场和稀释率的大致范围,完全可以不开流场。但经典复现如果涉及熔池中的温度梯度和成分混合分析,就需要把层流、传热、移动网格耦合起来。熔池对流的驱动力主要是表面张力的温度依赖性,也就是 Marangoni 效应。COMSOL 里可以在自由表面添加切向应力条件,或者直接用两相流接口追踪液-气界面。
开了流场之后,时间步长和网格质量要求会显著提升,计算成本成倍增长。比较踏实的顺序是:先跑纯热模型确定合理参数范围,再打开层流做最终验证。不建议从一开始就全物理场全开。
3.4 工艺参数的工程级手算校验
我强烈建议在跑模型之前,先用笔算一遍主要量级。单位长度能量输入 E = P / v,激光熔覆一般取 20~60 J/mm。举个例子,P = 2000 W,v = 10 mm/s,则 E = 200 J/mm,这个值偏高,容易出现基体过度熔化,需要降低功率或提高扫描速度。再看送粉量,熔覆层截面面积约为 A_clad = m_dot / (ρ*v)。取 m_dot = 0.2 g/s,ρ = 7800 kg/m³,v = 10 mm/s,则 A ≈ 2.56 mm²。如果希望得到的熔覆层宽 5 mm、高 0.5 mm,截面面积刚好是 2.5 mm²,与送粉量吻合,说明工艺参数搭配合理。
这种手算核算在复现文献时非常有用。很多时候你发现模型怎么调都和实验对不上,回头一核对,发现是文献里功率或送粉率的单位没看清楚,直接抄错了一个数量级。手算是最高效的防错机制。
4. 选区激光熔化逐层复现:多物理场与逐层激活
4.1 单道模型是经典复现的最佳起点
SLM 的工艺参数非常多,铺粉厚度、激光功率、扫描速度、扫描间距、扫描策略,都会影响成形质量。经典复现几乎无一例外地从单道扫描开始:一层粉末铺在基体上,激光沿一条直线扫过,观察熔池宽度、深度和球化倾向。
单道模型的几何非常简单,但物理过程是完整的:热源移动、粉末熔化、熔池流动、凝固收缩全都有。几何上就是粉末层加基体。这里有一个关键细节:粉末层的有效热导率远低于实体,这个差异必须体现在模型中,否则温度场会完全错误,熔池形状也会差之千里。
4.2 粉末到实体的材料状态切换
SLM 数值模拟最核心的问题是:同一个位置,扫描前是粉末,扫描后是实体,材料属性随温度动态变化。通常定义一个固相分数变量或相指示函数 φ,初始时粉末区域 φ = 0;当温度超过液相线,该单元内的 φ 趋向 1。然后有效热导率写成:
k_eff = k_powder * (1 - φ) + k_solid * φ
密度和比热容也做类似插值。这种方式能较好地模拟粉末烧结和实体熔化共存的状态。
粉末熔化后体积会收缩,粉末堆积密度通常在 58%~60%,凝固后变成致密实体,高度大约降为原来的 40%~50%。如果要用变形几何描述这个收缩,需要额外小心,网格很容易扭曲。稳妥做法是:第一版模型里先只用材料属性切换,不做几何收缩,等结果和文献趋势一致后再加入收缩。
4.3 Marangoni 对流对熔池形貌的影响
SLM 熔池的温度梯度非常大,表面张力随温度变化产生了很强的切向应力,驱动熔池流体从中心高温区流向边缘低温区,最终摊展出更宽的熔池。如果不考虑这个效应,熔池宽度通常会被低估。金属材料的表面张力温度系数 dγ/dT 一般是负值,约为 -0.2e-3 到 -0.5e-3 N/(m·K)。
在 COMSOL 中做这一步,需要对自由表面施加切向边界应力,方向和幅值都要仔细核对。这里易踩的坑是坐标系方向,比如表面法向与切向搞反,计算结果会明显不合理,熔池形貌也会变得很奇怪。建议先用非常简单的等温平板测试 Marangoni 边界是否产生预期的环流,确认无误后再放入完整激光模型。
4.4 从单道到多层多道的算力策略
完整的 SLM 工艺仿真计算量非常大,一个 5 层 5 道的简化模型动辄需要运行几十小时。常用策略是分步走:
先做单道,验证材料参数与热源模型;再做单层多道,观察搭接率是否合理;最后做多层少道,研究层间重熔和残余应力。网格在激光路径附近加密到 5~10 μm,远处必须逐渐粗化。热源方面,高斯面热源优先使用;体热源物理上更准确,但需要额外标定吸收深度,否则容易引入新偏差。
如果网格规模已经很大,还可以考虑使用对称面或周期性扫描路径复用,把计算区域减小一半。这里有个容易让人忽略的点:逐层激活几何的“事件”接口虽然看似方便,但层层之间的状态传递很容易引起数值不稳定。我个人经验是,先用参数化扫描一层一层手动算,跑通之后再考虑自动化。
5. 经典复现路上的高频坑与排查实录
5.1 热源位置不对:坐标系选择错误
这是出现频率最高的错误。症状是等温线始终在初始位置不动,或者热源固定在某个坐标附近移动不了。原因绝大多数是表达式里使用了全局坐标,而不是随时间变化的动态光斑中心。移动热源的正确写法类似于:
q = 2P/(πr_beam²) * exp(-2((x - x0(t))² + (y - y0(t))²)/r_beam²)
其中 x0(t)、y0(t) 是光斑中心的动态位置坐标。为了平滑启动,可以直接把 x0(t) 写成 x_start + v*t 的形式。如果担心起步瞬间热源跳变,再用一个平滑的 step 函数做过渡。测试时可以先给非常慢的扫描速度看等温线是否随热源一起推进,确认正确之后再恢复真实速度。
5.2 不收敛:相变潜热的数值处理
在 COMSOL 中处理相变潜热,常用方法是把潜热折算成表观比热容,在相变温度区间内给一个很大的比热峰值。这个峰值如果设置得过窄,比如金属熔化区间只给 1~2 K,很容易引起迭代震荡。解决办法是先展宽相变区间到 20 K 以上,跑通之后再逐步收窄;或者先用固定温度场跑热源加热,确认整体稳定后再加入相变。一个隐藏技巧是,在流场或压力场出现压力尖峰时,可以先把初始步长压到非常小,让前几步平稳通过相变点。
5.3 移动网格崩溃:变形超限
熔覆和烧蚀模型中最常见的中断原因是网格单元被压溃,尤其是熔覆层持续生长、边界位移很大时。排查顺序是:先检查是否已开启自动重新剖分(Automatic Remeshing);再检查最大单元变形限制是否设置得太宽松;最后确认边界速度的方向是否确实沿法向。如果这些都没问题,可以把移动网格的“平滑”选项切换到各向异性,往往能大幅改善单元质量。
5.4 计算时间爆炸:物理场全开并不划算
如果你同时打开传热、层流、移动网格、相变、Marangoni 对流,一个 0.1 ms 的物理时间可能要跑一整个晚上。实际经验告诉我,凡是一开始就全物理场全开的模型,大概率是在浪费 CPU。正确路径永远是先热模型、后流动模型,分阶段验证。计算过程中可以用“求解器配置”限制 BDF 阶数到 2,这通常比默认的全阶数更快且更稳定。先做一次网格无关性测试确定单元数,不要一上来就用最细网格,那只会让排查问题的时间成倍增加。
5.5 常见问题速查表
| 症状 | 可能原因 | 排查顺序 | 解决办法 |
|---|---|---|---|
| 等温线不动 | 热源坐标写死 | 检查光斑中心表达式 | 改成随时间变化的坐标变量 |
| 烧蚀深度偏小 | 时间步长太大 | 检查能量密度与时步 | 时间步长小于脉冲宽度1/10 |
| 烧蚀深度偏大 | 吸收率设为1 | 核对吸收率参数 | 金属 1064 nm 吸收率取0.3~0.5 |
| 熔覆层形状失真 | 网格变形过度 | 查看自动重新剖分设置 | 开启自动重新剖分并调小单元变形限值 |
| SLM熔池宽度偏窄 | 没开Marangoni对流 | 确认切向应力边界 | 加入表面张力温度系数 |
| 整体计算过慢 | 物理场全部无差别开启 | 检查求解器配置 | 去掉非必要接口,限制BDF阶数 |
6. 最后的实操体会
我个人在实际操作中体会最深的一点是,经典复现真正有价值的环节并不是“把模型跑出一个漂亮的彩色云图”,而是“把文献参数表里每一个数值填进模型时被迫思考它的物理含义”。能量密度为什么取这个量级,稀释率为什么控制在这个范围,Marangoni 数大概有多大,这些问题会逼着你把一个模糊的物理直觉变成清晰的数值表达。
还有一个很实用的小建议:所有热源表达式、材料函数、边界速度,都单独定义成变量或函数,不要直接写数字。后期换参数、做参数化扫描、出图对比都会非常方便,尤其在连续调参的时候,这个习惯能帮你省下大量时间。
最后分享一个经验判断,也是做这类仿真最容易忽略的:文献里给的是总功率还是吸收后的功率,单位是 J/cm² 还是 W/mm²,吸收率有没有单独列出来,这些单位细节决定了最终结果是否对得上。把这几项在建模前先统一换算清楚,比任何高超的数值技巧都管用。