EMPbridge软件应用
齿轮感应淬火的多时间尺度电磁-热多物理耦合分析
摘要：感应加热淬火是一项重要的热处理工艺，具有节能、高效的优点。为进一步优化感应加热流程，提高运行效率，基于自主开发的电磁及多物理耦合仿真平台EMPbridge，对电磁-热多物理耦合问题进行了仿真分析。针对电磁变化周期远小于感应加热升温时间（秒级）所导致的时间多尺度问题，采用频域、时域混合的算法，在频域中求解电磁场得到涡流损耗，进而耦合到热场中作为热源对热场进行时域求解，大幅降低了计算量，缩短了计算时间。考虑了磁导率和电导率随温度的变化，并详细阐述了电磁-热双向耦合问题的计算流程，仿真结果展示了感应加热参数与齿轮加热效果的关系，为工程实践和设备优化提供了指导。
关键词：多时间尺度；多物理耦合；感应加热；材料非线性
0. 引言
感应加热淬火是机械加工中一项至关重要的工艺，相比于传统的淬火方式，它具有操作简便、加热速度快、能源利用率高，更加节能环保等优势。为了实现更高的能源利用效率，需要借助于计算机模拟仿真技术对感应加热装置进行设计与优化。通过电磁-热多物理耦合仿真计算得到工件内涡流分布和温度场分布，根据温度场随时间的变化过程可以了解感应加热的效果，以此为依据可以对感应加热频率、激励线圈几何形状、激励线圈与工件的相对位置等参数进行优化设计。数值仿真相比于传统的现场调试可以大幅减少资金与时间的投入，因此开发相关仿真算法与专业软件有助于企业节约人力和时间成本。
本文基于自主开发的电磁及多物理仿真平台EMPbridge，对齿轮感应淬火中涉及的电磁-热多物理耦合过程进行了仿真，并考虑了材料参数随温度的变化，为后续进一步开发感应加热专用仿真和优化软件提供了参考。
1. 感应加热过程的电磁-热多物理耦合分析
1.1 电磁场数学模型
感应加热中涉及涡流损耗的计算，需要求解磁准静态场，其微分方程形式如下：
`\nabla\times\bb{E} = -\frac{\partial \bb{B}}{\partial t}` `(1)`
`\nabla\times\bb{H}=\bb{J}` `(2)`
`\nabla\cdot\bb{D}=\rho_e` `(3)`
`\nabla\cdot\bb{B}=0` `(4)`
式中：E为电场强度（V∙m^(-1)）；H为磁场强度（A∙m^(-1)）；D为电位移矢量（C∙m^(-2)）；B为磁感应强度（Wb∙m^(-2)）；ρ_e为电荷密度（C∙m^(-3)）；J为电流密度（A∙m^(-2)）。
由本构关系和电流密度的定义可以得到材料特性与电磁场的关系：
`\bb{B}=\mu\bb{H}` `(5)`
`\bb{J}=\sigma\bb{E}` `(6)`
其中μ为磁导率，σ为电导率，本文中它们是随温度变化的参数，其变化曲线通常通过测量得到。
由于感应加热的激励电流一般为正弦电流，假定材料参数不随磁场和电场变化，上式可以写为复数形式在频域求解：
`\nabla\times\dot{\bb{E}}=-j\omega\dot{\bb{B}}` `(7)`
`\nabla\times\dot{\bb{H}}=\dot{\bb{J}}` `(8)`
`\nabla\cdot\dot{\bb{D}}=\rho_e` `(9)`
`\nabla\cdot\dot{\bb{B}}=0` `(10)`
式中带上标的`\dot{\bb{E}}`, `\dot{\bb{D}}`, `\dot{\bb{H}}`, `\dot{\bb{B}}`, `\dot{\bb{J}}`分别表示各场量的复数形式，`ω`为激励电流的角频率。
结合相应的边界条件，在频域中求解以上磁准静态场问题，即可得到在感应加热频率下激励电流在齿轮工件中产生的涡流分布，进而得到相应的涡流损耗。
1.2 热场数学模型
感应加热过程中的热场由瞬态热传导方程描述：
`\rho C_p\cdot\frac{\partial T}{\partial t}-\nabla\cdot\lambda\nabla T=Q` `(11)`
式中：`\rho`为密度（`kg\cdot m^{-3}`）；`C_p`为比热容（`J\cdot kg^{-1}\cdot K^{-1}`）；T为温度（`K`）； `\lambda`为导热系数（`W\cdot m^{-1}\cdot K^{-1}`）；`Q`为热源密度（`W{\cdot m}^{-3}`）。热源来自于感应涡流产生的涡流损耗，其表达式为：
`Q_e=\frac{1}{2}\sigma^{-1}\dot{\bb{J}}\cdot\dot{\bb{J}^\ast}` `(12)`
结合边界条件和温度的初始条件，在时域中求解以上热扩散方程，可以得到各个时刻对应的温度场分布。
1.3 电磁-热多物理耦合过程
感应加热中涉及的电磁-热多物理耦合过程为双向耦合，即交变电磁场在工件中感应出涡流，涡流损耗作为热源引起温度场的变化，温度场变化进一步引起材料参数的变化，而材料参数的变化又导致电磁问题中涡流场分布的变化，其耦合示意图如图1所示。
图1 感应加热中电磁-热双向耦合示意图
在仿真计算中，电磁场问题在频域中求解，热场问题在时域中进行求解，其耦合求解过程如图2所示，即首先在求解区域给定一个初始温度分布并根据初始温度时的材料参数计算初始时刻（`t=t_0`）时的频域电磁场，获得涡流分布，然后将涡流损耗作为热源代入，计算时间步长为`\Delta t`的时域热场问题，得到`t_1`时刻的温度场分布。根据此温度场分布，结合材料参数（`\mu`,`\sigma`）随温度变化的曲线，可以得到`t_1`时刻温度变化后的材料参数。由于材料参数的变化，涡流分布也相应变化，因此`t_1`时刻的电磁场需要利用新的材料参数重新计算，并得到新的涡流损耗，以代入热场方程计算下一时刻的温度场分布。之后的过程依此类推，直到到达瞬态热场设置的仿真时长为止。
图2 电磁-热耦合求解过程示意图
值得注意的是，此过程为显式算法，即下一时间步的解可以通过前一时间步的解直接计算得到，因此时间步长的取值不能太大，否则单个时间步内温度变化太大会导致材料参数无法满足线性近似条件。同样，时间步长也不宜选取得过小，否则将增加不必要的计算量和计算时间。可根据某时刻磁导率和电导率的梯度决定时间步长的大小，梯度大（材料参数变化快）时时间步长小，反之，时间步长大。
2. 齿轮淬火数值分析算例
2.1 齿轮淬火仿真模型的建立
如图3所示，利用EMPbridge平台建立齿轮的几何模型。齿轮的几何尺寸采用模数为2.5的齿轮，齿数为26，可以计算得到齿轮的直径约为65 mm。
图3 齿轮几何模型示意图。齿轮的直径：65 mm，模数：2.5，齿数：26。
齿轮的材料是常用的45#钢，其通过测量得到磁导率和电导率参数随温度的变化分别如图4和图5所示。
图4 45#钢相对磁导率随温度的变化
图5 45#钢电导率随温度的变化
关于齿轮的网格剖分，需要考虑到涡流透入深度。涡流透入深度`\delta`受感应加热频率和材料的磁导率、电导率参数的共同影响，其计算公式为：
`\delta=\sqrt{\frac{2}{\omega\mu\sigma}}` `(13)`
仿真中选用的感应加热电流频率为50 kHz，在室温20℃时，根据图3和图4查到的相对磁导率和电导率可以计算得到涡流透入深度约为：0.07 mm。当工件升温至900℃时，由于温度已经高于居里温度，工件快速退磁，相对磁导率下降到1左右，且电导率也随温度下降，此时的涡流透入深度增加到2.4 mm左右。在透入深度范围内涡流幅值成指数形式下降，因此在齿轮表面透入深度范围内需要对网格进行精细剖分。本仿真的网格剖分如图6所示，在齿轮的齿部进行精细剖分，而在远离齿轮边缘的部分采用相对较粗的网格剖分尺寸。
图6 齿轮网格剖分示意图
2.2 涡流与温度场的变化规律
随着感应加热过程的进行，齿轮齿部的温度逐渐上升，而温度的上升将引起材料参数的改变，从而导致涡流场分布产生变化。如图7所示，在初始时刻，温度为室温，涡流透入深度很浅，涡流场集中于齿轮边缘。
图7 初始时刻(t=0s)时电流密度的云图分布
当齿轮齿根部位的温度上升到700℃时，此时的材料的磁导率到达居里温度附近，根据图3的曲线可知，到达居里点后温度的小幅升高将导致磁导率的大幅下降，此时透入深度增加，涡流场分布向齿轮内部扩散，如图8所示。
图8 齿根部位温度到达居里点后涡流场分布的变化
随着感应加热的进一步进行，齿轮的整个齿部都已经升温至居里温度以上，此时涡流分布再次沿齿轮外沿分布，但是透入深度增加到2.4mm左右，如图9所示。
图9 齿部温度均高于居里点后涡流场分布示意
由以上结果可见，涡流分布随温度的上升而变化，并反过来影响温度的分布，因此感应加热过程需要考虑电磁-热的双向耦合过程。
2.3 齿轮材料参数随温度的变化
由于感应加热过程中，齿轮齿部的温度高，而齿轮内部的温度低，因此齿轮不同部位的材料参数并不一致，而是与温度场的分布相关。如图10所示，在加热过程中，温度高的部位，磁导率和电导率下降明显，而齿轮内部由于温度升高不明显，其材料参数与室温时基本相同。对比图10和图8可以进一步印证1.3节中所阐述的电磁-热双向耦合过程。图10中齿根部位的温度升高导致此区域的磁导率和电导率大幅下降，根据公式（10）可知透入深度将升高，观察图8可以发现齿根部位的电流密度升高，涡流透入齿轮齿根部位。
图10 齿轮材料参数随温度的分布
3 结论
本文针对齿轮淬火加热的实际工程问题，对感应加热中涉及的多时间尺度电磁-热多物理耦合过程进行了剖析。在仿真中考虑了电导率和磁导率参数随温度的非线性变化，阐述了利用频域-时域混合的方法对电磁-热双向耦合问题进行仿真计算的方法。通过分析仿真结果，给出了涡流场分布受热场分布影响的规律，为齿轮感应加热的优化设计提供了参考。
附录：电磁-热耦合仿真结果
电流密度 (A/m2)
相对磁导率
温度 (℃)
电导率 (S/m)
参考：郑天宇, 周亚星, 靳剑钊. 齿轮感应淬火的多时间尺度电磁-热多物理耦合分析. 第18届中国CAE工程分析技术年会. 厦门. 2022