☰
激光淬火相变模拟:COMSOL多物理场耦合温度场与马氏体预测
2026/10/13 5:11:28 网站建设 项目流程

1. 激光淬火相变模拟的整体设计思路

1.1 为什么激光淬火仿真必须做多物理场耦合

做激光淬火工艺的人都有这种体会:参数窗口大多靠试,试一次切一次样,金相磨完看半天,功率密度调一档再来一轮,工装和工时都耗在试错上。如果能在开机之前把温度场和组织场先算清楚,硬化层深度、马氏体分布甚至残余应力都能提前看到,工艺调试就有了地图,而不是靠猜。这正是 COMSOL 激光淬火相变模拟要解决的核心问题。

激光淬火本质上是一个快速加热-快速自冷的过程:高能激光束扫过金属表面,表层在极短时间内被加热到奥氏体化温度以上,形成奥氏体,随后激光移开,热量迅速向低温基体传导,表层以远大于临界冷却速度的速率冷却,奥氏体转变成马氏体,从而获得高硬度、高耐磨的硬化层。整个过程横跨传热、固态相变和应力应变三个物理场,而且它们之间是强耦合关系:

  • 激光能量输入决定了温度场的分布和演变,温度历史又决定了奥氏体化的程度和马氏体转变的进程。
  • 相变本身会释放潜热,这部分热量会反馈到温度场中,改变局部温度峰值和冷却速度。忽略潜热会导致最高温度偏高、冷速偏快,算出来的硬化层深度和残余应力都不对。
  • 马氏体转变伴随体积膨胀和屈服行为变化,这是后续残余应力分析的源头。

如果只用普通热分析软件做一次纯传热计算,那只能知道哪里超过了奥氏体化温度,得不到马氏体分数分布,也就无法回答“硬化层到底有多深”这个工艺上最关心的问题。COMSOL 的优势在于多物理场耦合是原生机制,传热接口、固体力学接口、域常微分和微分代数方程接口可以绑定在同一套网格上,温度、相变量、应力直接共享节点数据,不需要自己在外部写循环来回倒数据,而且在同一个求解器里处理非线性耦合项,比如相变潜热的反馈项,稳定性比外部耦合好得多。这是我把这类问题放到 COMSOL 里做的主要原因。

1.2 模型简化到什么程度合适,先抓主要矛盾

激光淬火模拟的技术难点不在于“参数多”,而在于“哪些物理必须保留、哪些可以安全地砍掉”。我的经验是,工程上做工艺窗口预测,模型简化的优先级大致是下面这个顺序:

物理机制方面,激光淬火属于固态相变热处理,表面不能熔化。这一点和激光熔覆、激光焊接有本质区别,熔池流动、凝固、元素烧损这类涉及流体和相变界面的复杂问题可以直接排除,大大降低了建模难度。与“熔化-凝固”有关的所有模型都不需要,相变只考虑固态相变即可。

几何与维度方面,旋转体零件,比如轴、套类工件,采用二维轴对称模型就能反映主要热流方向,网格量和计算量只有三维的几十分之一;平面扫描类工件则必须用三维模型,但可以充分利用对称面,只建半个工件,甚至四分之一。扫描路径、工件边界的散热条件越对称,模型越省。

热源模型方面,连续激光淬火最常用的基础模型是高斯面热源,即把激光能量按高斯函数分布在光斑半径内,施加在工件表面。这个模型简单、参数少,工程上对硬化层深度和温度分布的预测精度已经足够。如果想要更细,可以用平顶光束或考虑吸收率随温度变化的修正,但第一步建议先用经典高斯热源把流程跑通,再逐步加复杂度。

相变动力学方面,加热阶段的奥氏体化属于扩散型相变,可以用 JMAK 型方程或简化的松弛方程描述;冷却阶段的马氏体转变属于无扩散型相变,用 K-M 方程描述。这两个模型在 COMSOL 中都可以通过添加状态变量和常微分方程实现,不需要外部程序介入。具体的参数拟合和写法后面会详细展开。

有一个地方需要特别提醒:激光淬火虽然加热速度快,但奥氏体化并不是瞬间完成的,只是奥氏体化的特征时间很短,通常只有零点零几秒到零点几秒。早期有人直接把“超过 Ac3 就认为是全奥氏体”当成判据,这在定性分析时够用,但会导致硬化层深度偏厚,尤其是在光斑边缘温度刚过 Ac3 的区域,那里实际上并没有完成充分的奥氏体化。所以做定量分析时,奥氏体化动力学的引入是必要的。

2. 建模之前必须吃透的热源参数、物性参数和相变方程

2.1 移动热源的数学表达与参数确定

COMSOL 里实现移动激光热源并不复杂,核心思路是在热源表达式中把热源中心坐标写成一个随时间移动的量。以 x 方向扫描为例,面热源公式写成:

Q_laser = eta * P_total / (pi * r0^2) * exp(-2 * ((x - x0 - v * t)^2 + (y - yc)^2) / r0^2)

其中:

  • P_total 是激光器输出功率,单位 W;
  • eta 是材料对激光能量的吸收率,量纲为一;
  • r0 是光斑半径,单位 m;
  • v 是扫描速度,单位 m/s;
  • x0 是热源在 t = 0 时刻的中心 x 坐标;
  • yc 是热源中心的 y 坐标。

这个热源以边界热源的形式施加在工件上表面,单位是 W/m²。注意 COMSOL 的传热接口里,边界热源和域热源的单位不同,边界热源是 W/m²,域热源是 W/m³。如果公式配错了单位,比如把面热源公式直接写进域热源,量级会差出数十倍,温度直接飞到离谱的数值,这是新手最容易踩的坑。

体热源形式则是在深度方向叠加一个指数衰减项:

Q_laser = eta * P_total * beta / (pi * r0^2) * exp(-beta * z) * exp(-2 * ((x - x0 - v * t)^2) / r0^2)

其中 beta 是材料对激光波长的有效吸收系数,金属在红外波段的吸收深度很小,约几十微米量级,所以体热源和面热源的差别在宏观厚度工件上其实不大。我通常先用面热源,简单且稳定。

工程上的参数取值范围大致可以参考下面这张表,具体数值要根据材料和设备标定:

参数典型范围说明
激光功率 P_total1200~3000 W功率太低无法奥氏体化,太高表面会熔化
扫描速度 v5~50 mm/s过慢导致热影响区过大,过快硬化层太浅
光斑半径 r02~5 mm实际光斑能量分布近似高斯,半径按有效直径一半取
吸收率 eta0.3~0.5钢对近红外激光室温吸收率约 0.3,表面涂黑或高温时可达 0.6
初始温度 T0293.15 K车间环境温度
对流换热系数5~15 W/(m²·K)静止空气自然对流,激光淬火主要靠基体自冷,该参数影响很小

吸收率是整个模拟中最不确定的参数。同一台激光器,在抛光表面和喷砂表面作用,实际吸收的能量差一倍都不稀奇。我一般处理方式是:先按典型值 0.4 设一个初始值,算完温度场后,用表面峰值温度或熔池边界与红外测温实测结果对比,反过来修正 eta,这叫反求标定。与其花大力气去查各种文献的吸收率表格,不如自己做一次 5 分钟的温度对比实验来得准。这是工程仿真和学术论文之间的一个重要区别。

2.2 温度相关热物性参数与相变潜热处理

钢铁材料的热导率、比热容、密度都随温度变化,尤其是热导率,从室温到奥氏体区可以下降 30% 以上。如果全程用常数物性,温度场尤其是冷却阶段的冷速会算偏。在 COMSOL 里,处理温度相关物性最方便的做法是在材料节点中用插值函数,把温度作为自变量,物性表做成内插数据,代码块里可以用interp(T, 表名)这样的函数引用。

比热容的处理有个关键点:相变潜热不能漏。奥氏体化和马氏体转变都会释放热量,对钢来说,奥氏体化的潜热大约在 80~100 kJ/kg 量级。这个热量看着不大,但在激光淬火这种快速冷却过程中,相当于给正在冷却的区域额外“补了一脚油”,会明显改变冷却曲线的形状和马氏体的转变量。如果忽略潜热,马氏体转变的起始时间和最终分数都会偏差。

潜热在 COMSOL 里的处理有两种常用方式:

第一种是等效热容法,把潜热折算到比热容上:

cp_eff = cp + L_aus * d(f_aus) / dT

这种做法相当于在相变温度区间给比热容加一个尖峰,物理上好理解,但在 COMSOL 里实现需要计算d(f_aus)/dT,而 f_aus 是一个状态变量,它依赖的是整个温度历史,不是瞬时温度的单值函数,直接用解析导数容易出问题。

第二种是显式潜热源项法,在热方程里加一个体积热源:

Q_latent = -rho * L_aus * d(f_aus, t)

其中d(f_aus, t)是奥氏体体积分数对时间的导数,f_aus 每增加一点,就释放 rho * L_aus * df/dt 的热量。因为 f_aus 在 COMSOL 中本身就是通过常微分方程计算出来的,d(f_aus, t)是现成的变量,所以这个做法实现起来最干净,数值稳定性也好。我推荐用第二种。

热物性方面还需要注意,固态相变带来的组织变化会影响热物性本身,最典型的是奥氏体与马氏体的热导率不同。这一层耦合在激光淬火中影响相对小,除非你研究的钢种对组织极其敏感。工程上先用一套代表平均组织的温度相关物性就够了,真要细化可以后续再加组织加权平均的物性模型。

2.3 奥氏体化动力学与马氏体转变方程的嵌入方法

相变动力学是激光淬火相变模拟的“灵魂”。温度场算得再准,相变模型不对,结果一样废。

先看加热阶段的奥氏体化。奥氏体化是扩散型相变,足够严格的描述是用 JMAK 方程:

f_aus = 1 - exp(-k * t^n)

其中 k 和 n 是从等温转变实验数据拟合出来的参数。但是激光淬火是非等温过程,温度在几秒内从室温飙升到千度以上,直接套用等温 JMAK 公式会带来误差。工程上更实用的做法是引入一个简化的松弛动力学模型:

d(f_aus, t) = (f_eq(T) - f_aus) / tau_aus(T)

当 T 高于 Ac1 时,f_eq 逐步趋近于 1,tau_aus 是奥氏体化的特征时间常数,通常取 0.05~0.2 秒。这个模型虽然不如完整 JMAK 严谨,但胜在参数少、数值稳定,而且对硬化层深度的预测精度已经足够。我一般把 Ac1 和 Ac3 的临界温度设为常数,在某中碳钢上分别取 1000 K 和 1050 K 左右。使用平滑的阶跃函数 flc2hs 来处理相变起始边界,可以避免温度刚跨过临界点时的数值突跳。

再看冷却阶段的马氏体转变。马氏体是无扩散相变,转变量只取决于过冷度,经典 K-M 方程写成:

f_mart = (1 - f_aus) * (1 - exp(-alpha * (Ms - T)))

其中 alpha 是马氏体转变速率常数,钢中经验值大约在 0.011~0.02 1/K 之间;Ms 是马氏体起始温度,常用中碳钢大约 500~600 K。这个方程的物理含义是:当过冷度增大时,马氏体转变量指数地增加,最终趋近于残余奥氏体之外的全部体积分数。

K-M 方程有一个坑:它依赖于当前温度,但马氏体转变是不可逆的。如果温度在 Ms 附近来回振荡,比如模拟中某节点在冷却到 450 K 后又小幅度回升到 480 K,直接用瞬时温度计算 f_mart,马氏体分数会随着温度回升而减少,这完全违反物理。所以必须记录每个材料点经历的历史最低温度 T_min,马氏体转变量只由 T_min 决定,而不是瞬时温度。在 COMSOL 里,T_min 本身也作为一个状态变量,用一阶微分方程更新:

d(T_min, t) = if(T < T_min, (T - T_min) / tau_smooth, 0)

其中 tau_smooth 是一个很小的数值光滑时间,比如 1e-3 秒,作用是避免条件切换时产生刚性。实际实现里还可以用 min 函数,但 COMSOL 的征解器对非光滑的 min 处理不如这种带阻尼的微分方程稳定。

完整的状态变量系统大概是三个:f_aus 记录奥氏体分数,T_min 记录历史最低温度,f_mart 作为代数变量输出马氏体分数。f_aus 和 T_min 在域常微分和微分代数方程接口中定义,f_mart 在变量定义中直接写表达式。这样一层层的依赖关系正好反映了激光淬火的物理过程:先有温度历史,再决定奥氏体化程度,最后在冷却过程中决定马氏体转变量。

3. COMSOL 激光淬火相变模拟实操流程

3.1 几何模型与物理场接口搭建

以一个平板激光淬火为例,几何尺寸取 100mm × 40mm × 20mm,激光沿 x 方向从一端扫描到另一端。这类问题如果工件关于扫描中心线对称,可以只建一半模型,减少一半网格量。几何建模直接画一个长方体即可,不需要任何圆角、倒角,热影响区主要集中在扫描线附近,远离区域对结果影响很小。

物理场选择有两个关键动作:

第一个是传热接口,研究类型选瞬态。初始温度设为 293.15 K,工件底面和四周默认绝热即可,顶面加载激光热源之外,还要加上表面对环境的对流换热,空气自然对流系数取 5 W/(m²·K)。辐射换热在激光淬火的温度范围内影响较小,但在温度接近 1000 K 时可以按等效对流处理,或者用表面辐射条件加上发射率 0.3。

第二个是域常微分和微分代数方程接口,在这个接口中添加奥氏体状态变量 f_aus 和历史最低温度状态变量 T_min。需要注意的是,这两个变量在每个材料点都存在,也就是说它们是空间分布的场变量,而不是单一全局变量。COMSOL 的 DAE 接口天然支持这种“每个节点一个状态”的模式,每次求解时会在每个网格节点上独立更新,这正是我们需要的。

在这个基础上,还要在模型定义节点中创建一组变量和函数:

  • 激光热源 Q_laser 的表达式
  • 潜热源项 Q_latent 的表达式
  • 马氏体分数 f_mart 的代数表达式
  • 奥氏体化速率 d(f_aus, t) 的动力学表达式
  • T_min 的更新微分方程

把这些定义好之后,物理场之间的耦合关系就建立起来了:传热接口提供温度 T,DAE 接口根据 T 更新相变量,相变量又通过潜热源项反馈到传热方程。这是一个双向耦合的温度场-相变场系统。

3.2 关键变量定义与表达式写法

这块是实操里最容易被绕晕的部分,我把 COMSOL 里实际使用的变量和表达式直接列出来,照着抄就能跑通。

在“定义-变量”中创建:

P_total = 1500[W] eta = 0.4 r0 = 3[mm] v = 10[mm/s] x0 = 5[mm] Q_laser = eta * P_total / (pi * r0^2) * exp(-2 * ((x - x0 - v * t)^2 + y^2) / r0^2)

如果模型建了对称一半,y 方向的坐标要根据对称轴调整,热源中心所在的 y 坐标设为对称面位置。

在“定义-变量”中继续定义:

Ac1 = 1000[K] Ac3 = 1050[K] tau_aus = 0.1[s] L_aus = 80[kJ/kg] Ms = 550[K] alpha_km = 0.011[1/K] rho_steel = 7850[kg/m^3]

奥氏体化的松弛动力学,写在 DAE 接口的源项中:

f_aus_eq = if(T > Ac1, min(1, (T - Ac1) / (Ac3 - Ac1)), 0) d(f_aus, t) = (f_aus_eq - f_aus) / tau_aus

这里用(T - Ac1)/(Ac3 - Ac1)做了一个线性过渡,温度在 Ac1 到 Ac3 之间时平衡奥氏体分数从 0 逐渐升到 1,超过 Ac3 后就是 1。这个过渡比硬切割平滑得多,也符合奥氏体化在相变区间渐进完成的物理图像。如果想要更严苛的处理,可以在 f_aus_eq 中再加入温度项,模拟更高温度下奥氏体化速度加快的影响。

潜热源项写在传热接口的域热源中:

Q_latent = -rho_steel * L_aus * d(f_aus, t)

注意负号不是笔误。奥氏体化是一个吸热过程,把珠光体加热到奥氏体需要吸收潜热,所以源项应该是负的,表示这部分热量被相变“吃掉”了。冷却时马氏体转变放热的作用,已经在温度低于 Ms 后的冷却曲线中被“降温放缓”体现出来,采用奥氏体化方向的单一潜热项在工程精度内是合理的近似。

T_min 的更新在 DAE 接口中定义:

d(T_min, t) = if(T < T_min, (T - T_min) / 1e-3[s], 0)

马氏体分数的代数变量:

f_mart = (1 - f_aus) * (1 - exp(-alpha_km * (Ms - T_min))) * if(T_min < Ms, 1, 0)

这里用if(T_min < Ms, 1, 0)确保温度根本没降到 Ms 以下时马氏体为零,非常直观。

最后是在传热接口的边界条件中激活边界热源,把 Q_laser 填进去。另外把 Q_latent 填到域热源中。这样整套耦合就设置完毕。

3.3 网格划分、时间步长与求解器调试

激光淬火模拟的网格策略核心是“局部加密”。光斑半径只有几毫米,而热影响区深度通常只有一到两毫米,如果均匀网格,要么整体网格太密导致计算量爆炸,要么关键区域太粗导致温度尖峰被抹平。我的做法是:

  • 在激光扫描线附近的表面区域做细网格,表面单元边长控制在 0.2 mm 以下,最好 0.1 mm。
  • 深度方向前 1 mm 内划分 5~10 层边界层网格,第一层厚度 0.02 mm,按 1.2 倍递增。这样可以把近表面的温度梯度和相变分布捕捉得比较准。
  • 远离热源的区域使用较粗的扫掠网格,边长 2 mm 即可。

时间步长与空间步长必须匹配。激光以 10 mm/s 的速度扫描,光斑半径 3 mm,光斑在某个固定点上停留的有效时间大约 0.6 秒。如果表面网格边长是 0.2 mm,那么理想时间步长不能超过光斑扫过 0.2 mm 所需的时间,也就是 0.02 秒。可以使用自由时间步长让求解器自适应,但要在求解器设置中把最大时间步长限制在 0.02 秒左右,避免跳过太多细节。

求解器设置建议用 BDF 格式,最大阶设为 2。状态变量 T_min 的微分方程带有 if 条件,是轻度刚性的,一阶或二阶 BDF 比高阶格式更稳定。在求解过程中,如果发现时间步长被压缩到 1e-6 秒的量级还反复不收敛,优先检查相变动力学表达式中是否有未平滑的阶跃切换,比如温度刚过临界点时 f_aus_eq 的不连续会触发求解器的剧烈反应。把硬切换改成 flc2hs 平滑后,计算速度通常能提升一个数量级。

一个非常实用的调试策略是“先热后相变”。第一次跑的时候把相变潜热源项和 DAE 接口都停用,只跑纯热的移动热源模型,确认温度场峰值、云图形态正常,再打开相变耦合。这样如果结果异常,就能迅速定位是热源设置问题还是相变模型问题。我在实际项目里一直都是这个流程,省掉了很多排查时间。

4. 结果解读、硬化层判据与常见问题排查

4.1 从模拟数据中提取硬化层深度

模拟计算完成之后,最直接的问题就是:硬化层在哪里,有多深?

最常用的工程判据是马氏体体积分数阈值。对大多数中碳钢,马氏体体积分数大于 0.5 的区域可以认为达到淬硬要求,硬度值通常在 500 HV 以上;如果要求更严格,可以取 0.8 作为判定线。这个阈值要根据实际钢种的硬度-组织关系标定,但作为初判已经足够。

在 COMSOL 后处理中提取硬化层深度,我一般这样操作:

第一步,创建“二维截面”数据集,在垂直于扫描方向的 yz 平面上切一刀,截面位置取扫描线中心。第二步,在截面上绘制 f_mart 的云图,可以清楚地看到硬化区呈月牙形,靠近表面处马氏体分数最高,沿深度方向逐渐衰减。第三步,创建“一维线图”,横坐标取深度方向,纵坐标取 f_mart,曲线从表面往下延伸。找到 f_mart 等于 0.5 对应的深度,就是硬化层深度。

另一种互补的判据是最大温度判据。可以用状态变量 T_max 记录每个材料点经历的最高温度,然后后处理中绘制 T_max > Ac3 的区域轮廓,这个轮廓代表“理论上发生过完全奥氏体化的区域”。它通常比马氏体分数 0.5 的区域略大,因为奥氏体化不充分时温度可能超过 Ac3 但保温时间不足。两者对比能看出动力学模型带来的修正量,这是个很有价值的诊断信息。

举一个具体的模拟结果示例,参数取值 P_total = 1500 W、v = 10 mm/s、r0 = 3 mm、eta = 0.4:

工况表面峰值温度 (K)硬化层深度 (f_mart > 0.5)硬化区最大宽度 (mm)
工况 115201.16.2
工况 2 P=2000 W17801.77.4
工况 3 v=20 mm/s13800.85.5

从这个表里可以直观看到,功率提升 33%,硬化层深度增加约 55%;速度加倍,硬化层深度反而下降。这些趋势和实际工艺经验完全一致,说明模型物理上是可靠的。有了这两张曲线和云图,工艺工程师可以直接拿来做参数窗口设计。

4.2 常见病态结果与修复方法

我总结激光淬火相变模拟里最常见的四类“翻车现场”:

第一类是温度场出现锯齿状振荡,或者云图上有规律分布的亮点。这是典型的网格与时间步长失配。光斑在网格间移动时,如果时间步长太大,热源中心会从网格节点跳到节点,导致温度场呈脉冲式更新。修复方法简单:加密表面网格,同时把最大时间步长缩小到光斑扫过一个网格所需的时间以内。

第二类是马氏体分数出现负值或者大于 1。这个问题几乎总是出在 T_min 的更新方程或者 f_mart 表达式上。检查 f_mart 的公式里是否用了瞬时温度而不是历史最低温度,检查 T_min 的更新是否在温度回升时错误地往回拉。如果还是有问题,用 min(max(f_mart, 0), 1) 做一层保险,至少物理上不会出现超界值。

第三类是表面峰值温度高到离谱,比如超过 3000 K 甚至更高。原因通常有三个方向:激光热源的功率密度单位不对,边界热源 W/m² 和域热源 W/m³ 搞混;吸收率 eta 设定过高;或者潜热源项的符号反了。先停用相变做纯热模型验证,如果纯热模型温度就飞到几万度,那一定是热源本身的单位或系数有问题。

第四类是求解器时间步长被压缩到极小、计算卡死。多数情况下是相变状态变量的条件切换太突兀。比如if(T > Ac1, 1, 0)这种硬阶跃会让数值雅可比矩阵发生剧烈抖动,BDF 求解器为了满足误差容限只能不断缩小步长。修复方式是使用 flc2hs 平滑函数,把阶跃变成 10K 宽度的平滑过渡。这个修复对计算性能的提升立竿见影。

4.3 与实验对标:模拟不是用来“好看”的

模拟的最终价值是对标实验、指导工艺,所以验证环节不能省。我常用的对标流程有三条线:

第一是温度场对标。在工件近表面埋热电偶,实测几个特征点的温度循环曲线,与模拟曲线对比。重点看升温峰值和冷却段斜率。如果表面峰值温度对不上,优先修正吸收率 eta;如果冷却斜率对不上,检查材料热导率和边界换热设置。温度曲线是唯一能够定量验证传热模型的硬证据,值得多花时间做。

第二是组织对标。淬火后切样、磨抛、腐蚀、金相观察。模拟中 f_mart 等值线轮廓应该与金相试样的白亮硬化层边界吻合。如果模拟预测的硬化层深度偏浅,奥氏体化动力学参数 tau_aus 可能偏大,也就是奥氏体化速率设慢了;如果偏深,则检查 Ac1、Ac3 的取值是否偏低。需要说明的是,模拟预测的组织边界往往比金相看到的更“糊”一些,原因是模型把奥氏体化处理成了渐变过程,而金相观察到的边界人眼会主观锐化,两者对及时时候需要留有判断余地。

第三是硬度对标。显微硬度梯度是最终指标。把马氏体体积分数和硬度建立线性关系是工程上常用的桥接手段,虽然钢的硬度还受含碳量和回火状态影响,但对同一牌号材料、同一工艺窗口,用硬度-马氏体分数拟合曲线来预测硬化层深度已经足够。实际对标时,先测表面硬度和硬化层深度,反向校验 f_mart 阈值;再用这个阈值去预测新工况,如果新工况预测也准,说明整套模型参数已经算标定合格。

常见偏差来源的优先级排序是:激光吸收率不确定性排第一,其次是相变动力学参数拟合偏差,再次是材料物性的温度相关性精度。边界对流换热系数排在最后,激光淬火的自冷主要靠基体导热,表面对流换热影响很小。

5. 进阶扩展:从相变结果到残余应力与工艺优化

5.1 在热-相变耦合基础上加入应力场

如果做完了热-相变的耦合模拟,再往前走一步就是残余应力。激光淬火在表面形成马氏体的同时,伴随着体积膨胀和组织屈服行为的剧烈变化,淬火后表面一般会有压应力或拉应力,这直接影响零件的疲劳寿命和开裂倾向。

在 COMSOL 中把固体力学接口加进来,与传热接口和 DAE 接口联合求解,数学上不复杂,但数值上非常难收敛。原因在于弹塑性本构本身是强非线性,再加上马氏体转变引起的体积应变和相变塑性,方程组的刚性比纯热-相变问题高一个量级。工程上我不建议一步到位直接做全耦合。更稳妥的策略是分两个阶段:

第一阶段,只算热-相变场,并导出每一时刻的温度场和相变场分布;第二阶段,把温度场作为载荷加载到固体力学接口中,再把马氏体转变的体积膨胀处理成等效热膨胀,计算应力和变形。当然,COMSOL 本身支持在一个模型里同时求解,但当模型规模较大时,分阶段操作不管是对调试效率还是对结果可解释性都更好。

马氏体转变的体积膨胀是残余应力模拟的关键输入。钢从奥氏体到马氏体的比体积约增加 1%~2%,可以在固体力学接口中把热膨胀系数设置成依赖于 f_mart 的等效值,比如每增加 1% 马氏体,体积膨胀对应线性应变约 0.003~0.005。相变塑性更加复杂,涉及加载方向与相变交互作用,工程上建议先忽略或者以等效系数合并到体积膨胀项中。

5.2 多道搭接扫描与工艺参数优化思路

实际零件的激光淬火很少是单道扫描,更多是多道搭接。多道扫描时,第二道激光的热循环会部分回火第一道形成的马氏体,产生软化带。如果不考虑回火效应,模拟预测的硬化层会偏乐观。用 COMSOL 做多道扫描时,可以扫描路径设为两段或多段移动热源,让时间轴上热源中心连续切换。计算时间线性增加,网格策略不变。

回火软化目前没有简单的经验公式直接嵌进 COMSOL,但可以通过两阶段近似:先完整模拟多道扫描的温度场和组织场,找到第二道激光对第一道硬化区造成二次加热的最高温度及时间,再根据回火参数-硬度的经验关系修正硬度预测。这是一种近似工程化处理,工程上够用且易于解释。

参数优化方面,COMSOL 自带的参数化扫描可以自动批量计算功率、速度、光斑半径的不同组合。典型做法是固定材料牌号和工件尺寸,对 P、v、r0、eta 做四因素三水平正交参数组,每个工况输出硬化层深度、表面峰值温度和热影响区宽度,然后构建工艺窗口映射。这套流程跑完,原本要 30 次试切实验的活,可能只需要 5 次验证实验就够了。

我在实际项目中,最终的工艺窗口通常是以“硬化层深度 ≥ 1.0 mm,同时表面峰值温度低于液相线”为约束,以加工效率和能耗为优化目标拉出来的。这个窗口可以直接交给车间进行首件验证。仿真不是替代实验,而是把实验次数从“盲试”压缩到“确认”,这才是激光淬火相变模拟在工程里最实际的价值。

关于这套模拟方法,我自己的体会是:最大的坑往往不是相变模型本身,而是基础传热参数没校准就开始跑复杂耦合。先把热源、吸收率、物性这三样搞扎实,相变模型哪怕用简化版,结果也基本靠谱。另一个实用建议是牢牢记住“先纯热、后热-相变、再热-相变-应力”的渐进流程,每一步都确认前一阶段结果正常后再往后走。这样既省时间,也容易在出问题时缩小排查范围。

需要专业的网站建设服务?

联系我们获取免费的网站建设咨询和方案报价,让我们帮助您实现业务目标

立即咨询