CN107391822B - 一种基于自适应时间步长的瞬变电磁-温度场耦合计算方法 - Google Patents

Abstract

本发明提出一种基于自适应时间步长的瞬变电磁‑温度场耦合计算方法，采用指数平滑法预测电磁‑温度场耦合时间节点。并在两耦合时间节点间，通过预测‑校正法和响应特征值计算电磁场和温度场最佳离散步长。与传统等步长耦合方法对比，电磁、温度场均使用最佳的离散步长，避免了温度场过频的计算，减小计算时间。最后，以通电铜导环为例，采用自适应步长耦合计算铜导体在交流电下0.1s内温升，比传统等步长耦合方法时间减小20%，证明该方法的有效性。

Description

一种基于自适应时间步长的瞬变电磁-温度场耦合计算方法

技术领域

本发明一种基于自适应时间步长的瞬变电磁-温度场耦合计算方法，涉及电磁温度计算领域。

背景技术

在线圈发射装置、继电器、电机等设备中，运行时产生过高的温升将影响线圈和接头导电性、材料的绝缘性能，甚至会产生破坏性的热膨胀，对于设备运行性能及安全性产生一定的影响。在上述问题中，瞬变的大电流使得导体的温度在较短的时间内迅速增加，因而在计算过程不仅考虑电磁场对温度场的影响，还需考虑温度变化对电磁场材料导电性能的影响，为此需要建立瞬变电磁-温度场顺序强耦合模型。目前采用有限元法求解瞬变电磁-温度场间接耦合时，电磁场和温度场往往采用相同的时间步长。然而在求解过程中，温度场响应较电磁场慢，若采用相同的时间步长计算，使得温度场求解过频，造成了求解时间的增加。

对于两物理场瞬态间接耦合时，物理场响应时间不同的问题。有方法根据现有软件，提出代码耦合概念，对各物理场采用不同的时间离散策略，并在时间的耦合点上进行载荷的传递。也有方法针对在流-固耦合时，对流体区域和固体区域采用不同的时间步长的DFMT-SCSS算法，选取固体求解区域时间步长为流体求解区域的倍数计算。然而上述方法在计算两物理场耦合时，两物理场时间步长只是存在简单的倍数关系，而物理场都采用恒定步长计算，存在两物理场未获得最佳离散策略的情况。为此提出自适应步长耦合概念，即各物理场采用自适应步长算法获取最佳离散时间策略，并在预测的耦合时间节点上自动耦合。可以使各物理场获得最佳的时间离散策略，也避免了响应较慢的物理场过频的计算。

发明内容

为解决上述技术问题，针对瞬态电磁-温度场间接耦合问题，本发明提出一种基于自适应时间步长的瞬变电磁-温度场耦合计算方法。在T_n时刻，如图1所示，首先通过指数平滑法预测下一耦合时间点t_n+1。然后在两耦合时间节点间，通过预测-校正法和响应特征值计算电磁场和温度场最佳离散步长。得到t_n至t_n+1时间段内电磁场和温度场最佳时间离散步长为△t_n ^E和△t_n ^T。电磁场、温度场分别采用△t_n ^E、△t_n ^T计算至耦合时间节点。使各物理场在两耦合时间节点上，以最佳的时间步长计算，避免了响应较慢的物理场过频的计算。

本发明所采用的技术方案是：

一种基于自适应时间步长的瞬变电磁-温度场耦合计算方法，包括以下步骤：

1)、电磁场与温度场的耦合时间点确定：采用温度场触发热量来判断电磁-温度场耦合时间节点，首先电磁-温度场步长耦合三次，用来获取温度场触发热量Q_pre的预测数据，然后采用指数平滑法预测温度场触发热量，通过触发热量判断耦合时间节点；

2)、电磁场自适应步长计算：通过载荷离散误差和响应特征值确定电磁场载荷最佳离散的步长△t_n ^E，当电磁场累积热量达到温度场触发热量时，暂停电磁场计算，计算温度场；

3)、温度场自适应步长计算：通过载荷离散误差和响应特征值确定电磁场载荷最佳离散的步长△t_n ^T，将电磁场计算平均热功率作为载荷加至温度场，当温度场计算时间与电磁场同步时，暂停温度场计算，更新电磁场节点温度，并进行下一时间步的温度场触发热量计算，反复迭代计算至最终时间。

电磁场与温度场的耦合时间点确定包括以下步骤：

步骤1)：采用温度场触发热量来判断电磁场与温度场进行耦合计算时的时间节点，首先电磁-温度场步长耦合三次，用来获取温度场触发热量Q_pre的预测数据；

步骤2)：由式(1)、(2)、(3)计算电磁场计算过程中允许材料变化最大温度。其中线性材料的热功率P与电流密度J的函数关系如式(1)、(2)：

P_n＝∫_V J_n ²ε(T)dV (1)

ε＝aT+ε₀ (2)

式中：P_n为t_n时刻发热功率；V为单元体积；J_n为T_n时刻电流密度；ε为电阻率，a为电阻随温度变化率，ε(T)为电阻率与温度函数，ε₀为0℃时电阻率；∫_V dV表示体积分；

当输入热功率为P_n，由温度变化引起功率计算误差应满足：

|P_n(T_n)-P_n(T_n+△T)|≤γP_n(T_n) (3)

式中：P_n(T_n)为功率与温度函数，如式(1)；γ为由温度变化引起的功率误差，取1％；

式(3)取等号时，可得出t_n时刻温度允许变化最大值△T_max；

步骤3)：温度场触发热量计算，根据t_n时刻最大温度变化△T_max。由t_n，t_n-1，t_n-2时刻输入功率及温度变化，预测t_n+1温度场触发热量

计算t_n，t_n-1，t_n-2时刻温度随热量变化率k_n-1，k_n-2，k_n-3，如式(4)：

式中：Q_n，Q_n-1，Q_n-2分别为t_n，t_n-1，t_n-2时刻发热量；T_n，T_n-1，T_n-2，T_n-3为t_n，t_n-1，t_n-2，t_n-3时刻温度；

最后计算温度场触发热量计算，采用如式(5)的指数平滑方法预测t_n+1时刻温度场变化率K_n+1：

K_n+1＝αK_n+(1-α)K_n-1+(1-α)²K_n-2 (5)

式中：α为指数平滑系数，α＝0.75；k_n-1，k_n-2，k_n-3如式(5)；

温度场触发热量为：

Q^pre＝△T_maxK_n+1 (6)

式中：Q^pre为温度场触发热量；△T_max为允许变化最大温度，见式(4)；K_n+1计算见式(5)；

电磁场自适应步长计算包括以下步骤：

步骤4)：根据载荷离散误差确定电磁场时间步长

t∈(t_n-1,t_n)时，当载荷矩阵P采用线性插值时，等效载荷采用斜坡加载，由载荷离散产生的误差为如式(7)，

式中：

为第n个时间步长的载荷离散误差；f_s(t)为电磁场中电流或电压载荷；△t为离散步长；

对式(7)采用梯形公式积分，得离散误差近似如式(8)；

式中：f_s(ξ)为电磁场中电流或电压载荷；ξ为t_n-1至t_n间常数；∝表示正比于，其余参数如式(7)；

根据式(8)，载荷离散误差近似正比于(△t)²，可将下一步长计算可分为如下两步：

(a)步长预测：采用式(9)根据第n步计算误差，预测第n+1步长△t¹ _n+1；

式中：

为第n+1步预测步长；

为安全系数，

e_tolerance为允许最大误差；

为第n步产生的离散误差；

(b)步长校正：判断当第n+1时间步长△t¹ _n+1所产生误差是否满足e^k+1<e_tolerance，如不满足采用(8)进行修正迭代计算，直至满足e^k<e_tolerance；

式中：△t¹ _n+1为第n+1个时间步长；

为安全系数，

e_tolerance为允许最大误差；e^k+1为第k+1次迭代产生的离散误差；

步骤5)：根据响应特征值确定电磁场时间步长

采用响应特征λ_r值确定稳定时间步长△t_n+1。定义△t_n+1λ为震荡限制条件，当△t_n+1λ>>1时系统处于震荡状态，为保证计算稳定性可取最大步长需满足：

式中：f<1，f为稳定系数；λ_r为响应特征值；u为场量，电磁场中

A为矢量磁位，V为电压；u_n为t_n-1到t_n时间段场量u的变化；K为一阶有限元方程中传导矩阵；C为一阶有限元中容性矩阵；

步骤6)：电磁场在t_n时刻离散时间步长为：

式中：

为电磁场t_n时刻离散时间步长；

为式(10)中载荷离散最大步长；

为特征值确定离散步长；

步骤7)：耦合时间判断：当t_n+1时刻电磁场累积热量满足以下条件之一时，t_n+1为电磁-温度场耦合时间点，暂停电磁场计算，启动温度场计算；

(a)：单元最大热量变化达到步骤5)中单元预测热量阈值时自动启动温度场计算：

式中：Qⁱ为第i个单元的累积热量；

为第i个单元的预测触发热量，由式(6)计算；

(b)：总体热量变化达到步骤5)中单元预测热量阈值时自动启动温度场计算：

式中：Q^pre为温度场总体触发热量，由式(6)计算；Qⁱ为第i个单元的累积热量；β为触发安全系数，β＝0.9；

电磁场自适应步长计算包括以下步骤：

步骤8)：计算温度场自适应时间步长

温度场自适应时间步长计算同步骤4)、步骤5)、步骤6)；

步骤9)：当温度场计算时间到达步骤7)中耦合时间T_n+1时，停止温度场计算，并更新电磁场节点温度；

步骤10)：反复迭代步骤2)～步骤8)，直至计算总时间T_total。

本发明一种基于自适应时间步长的瞬变电磁-温度场耦合计算方法，优点在于：

1)、由于顺序耦合时，在计算电磁场的单个时间步内，并未考虑温度场的变化，因而存在累积的非平衡误差。采用温度场触发热量来判断耦合时间节点，可以有效控制非平衡误差。

2)、与传统等步长耦合方法相比，避免了由于温度场响应时间比电磁场长，而采用统一计算时间步长导致温度场求解次数过多而增加计算时间的问题，减小了计算时间。

附图说明

图1是电磁-温度场自适应耦合示意图。

图2是电磁-温度场自适应步长耦合流程图。

图3是载荷离散误差图。

图4(a)是铜导环结构示意图；

图4(b)是有限元模型图。

图5(a)等步长温度计算温度分布图；

图5(b)自适应步长温度分布图。

图6(a)是1号单元计算结果图；

图6(b)是140号单元计算结果图。

具体实施方式

一种基于自适应时间步长的瞬变电磁-温度场耦合计算方法，包括以下步骤：

1)、初始化。首先建立分析对象的电磁场和温度计算有限元模型，并给定电磁场计算中的磁导率、电阻率和温度场计算的比热容、导热率等材料参数。最后施加相应的初始条件、边界及载荷；

2)、耦合时间点确定。采用温度场触发热量来判断电磁温度场耦合时间节点。首先电磁-温度场步长耦合三次，用来获取温度场触发热量Q_pre的预测数据。然后采用指数平滑法预测温度场触发热量，通过触发热量判断耦合时间节点；

3)、电磁场自适应步长计算。通过载荷离散误差和响应特征值确定电磁场载荷最佳离散的步长△t_n ^E。启动电磁场，采用最佳离散步长计算。当电磁场累积热量达到温度场触发热量时，暂停电磁场计算；

4)、温度场自适应步长计算。同3)中，通过载荷离散误差和响应特征值确定电磁场载荷最佳离散的步长△t_n ^T。将电磁场计算平均热功率作为载荷加至温度场，当温度场计算时间与电磁场同步时，暂停温度场计算。更新电磁场节点温度，并进行下一时间步的温度场触发热量计算。反复迭代计算至最终时间。

具体步骤如下：

步骤1)：建立电磁场计算有限元模型，并给定磁导率、电阻率材料参数，施加载荷。由似稳条件忽略位移电流，时变电磁场方程可写为如下形式的矢量磁位方程：

式中：t为时间(s)，Ω₁为导体区域，Ω₂为非导体区域，A为矢量磁位(Wb/m)；V为电位(V)；σ为电导率(S/m)；μ为磁导率(H/m)；J_s为源电流密度(A/m²)，v为导体的运动速度，

为旋度运算，

为梯度运算，

为散度运算；

采用伽辽金法将上式写成有限元格式：

式中：

R为电磁场阻尼矩阵；S为电磁场系数矩阵；J为磁场载荷矩阵。

步骤2)：建立温度场场计算有限元模型，并给定比热容、导热率材料参数。在只考虑热传导和对流条件下，温度场导热微分方程可写成如下形式：

式中：T为温度(℃)，t为时间(s)，ρ为密度(Kg/m)，C_p为比热容(J/(Kg·K))，k为导热率(W/(m·K))，Q为发热功率W。

为梯度运算，

为散度运算，

为温度对时间的偏导。

初始条件及边界条件为：

式中：G(x,y,z)表示初始温度分布；F(x,y,z)表示恒定温度边界条件；Γ_T表示恒温边界；Γ_q表示散热量边界条件，q为边界发热功率，h为边界对流换热系数，n为边界法线向量。

采用伽辽金法将式(3)，(4)形成有限元格式如下：

式中：C为温度场比热矩阵，K为温度场导热矩阵，P为载荷矩阵。T为有限元节点温度矩阵，

为有限元节点温度梯度矩阵；

步骤3)：采用温度场触发热量来判断电磁温度场耦合时间节点。电磁-温度场采用较小步长耦合三次，用来获取温度场触发热量Q_pre的预测数据。

步骤4)：获取电磁场计算过程中允许材料变化最大温度。线性材料热功率P与电流密度为J关系：

P_n＝∫_V J_n ²ε(T)dV (6)

ε＝aT+ε₀ (7)

式中：P_n为t_n时刻发热功率；V为单元体积；J_n为T_n时刻电流密度；ε为电阻率，a为电阻随温度变化率，ε(T)为电阻率与温度函数，T表示温度，ε₀为0℃时电阻率；∫_V dV表示体积分；

当输入热功率为P_n，由温度变化引起功率计算误差应满足：

|P_n(T_n)-P_n(T_n+△T)|≤γP_n(T_n) (8)

γ为由温度变化引起的功率误差，取1％；式(8)取等号时，可得出T_n时刻温度允许变化最大值△T_max。

步骤5)：温度场触发热量计算。根据T_n时刻最大温度变化△T_max。由T_n，T_n-1，T_n-2时刻输入功率及温度变化，预测T_n+1温度场触发热量

首先计算各时刻输入热量及温度变化，如表1。

表1各时刻热量、温度

再计算各时刻温度随热量变化率，如表2。

表2时刻温度变化率

最后计算温度场触发热量计算。采用指数平滑方法预测t_n+1时刻温度场变化率K_n+1：

K_n+1＝αK_n+(1-α)K_n-1+(1-α)²K_n-2 (9)

α为指数平滑系数，α＝0.75；

温度场触发热量为：

Q^pre＝△T_maxK_n+1 (10)

式中：Q^pre为温度场触发热量；△T_max为允许变化最大温度。

步骤6)：启动电磁场，根据载荷离散误差确定电磁场时间步长

t∈(t_n-1,t_n)时，当载荷矩阵P采用线性插值时，等效载荷采用斜坡加载，如图3所示。由载荷离散产生的误差为如式(13)。

式中：

为第n个时间步长的载荷离散误差；f_s(t)为电磁场中电流或电压载荷；△t为离散步长；

对式(13)采用梯形公式积分，得离散误差近似如式(14)。

式中：f_s(ξ)为电磁场中电流或电压载荷；ξ为t_n-1至t_n间常数；∝表示正比于；

根据式(14)，载荷离散误差近似正比于(△t)²，可将下一步长计算可分为如下两步：

(a)步长预测。根据第n步计算误差，预测第n+1步长△t¹ _n+1。

式中：

为第n+1步预测步长，

为安全系数，

e_tolerance为允许最大误差；

为第n时间步内产生的离散误差。

(b)步长校正。判断当第n+1时间步长△t¹ _n+1所产生误差是否满足e^k+1<e_tolerance。如不满足采用(8)进行修正迭代计算，直至满足e^k<e_tolerance。

式中：△t¹ _n+1为第n+1个时间步长；

为安全系数，

e_tolerance为允许最大误差；e^k+1为第k+1次迭代产生的离散误差；

步骤7)：根据响应特征值确定电磁场时间步长

Hughes提出根据响应特征λ_r值确定计算稳定时间步长△t_n+1的方法，定义△t_n+1λ为震荡限制条件。当△t_n+1λ>>1时系统处于震荡状态。为保证计算稳定性可取最大步长需满足：

式中：f<1，f为稳定系数；λ_r为响应特征值；u为场量，电磁场中

A为矢量磁位，V为电压；u_n为t_n-1到t_n时间段场量u的变化；K为一阶有限元方程中传导矩阵；C为一阶有限元中容性矩阵。

步骤8)：电磁场在t_n时刻离散时间步长为：

式中：

为电磁场t_n时刻离散时间步长；

为式(10)中载荷离散最大步长；

为特征值确定离散步长。

步骤9)：耦合时间判断。当T_n+1时刻电磁场累积热量满足以下条件之一时，T_n+1为电磁-温度场耦合时间点，启动温度场计算。

(a)单元最大热量变化达到步骤5)中单元预测热量阈值时自动启动温度场计算：

式中：Qⁱ为第i个单元的累积热量；

为第i个单元的预测触发热量；

(b)总体热量变化达到步骤5)中单元预测热量阈值时自动启动温度场计算：

式中：Q^pre为温度场总体触发热量，由式(6)计算；Qⁱ为第i个单元的累积热量；β为触发安全系数，β＝0.9。

步骤10)：启动温度场计算，计算温度场自适应时间步长

温度场自适应时间步长计算同步骤6)、步骤7)、步骤8)；

步骤11)：当温度场计算时间到达步骤9)中耦合时间T_n+1时，停止温度场计算，更新电磁场节点温度。

步骤12)：反复迭代步骤4)～步骤10)，直至计算总时间T_total。

实施例，以通电铜导环瞬态温升为例:

本发明以电铜导环热分析为例说明自适应时间步长在电磁-热耦合中的应用，该模型广泛存在于继电器、线圈发射等装置中。将一长10mm，厚2mm铜导环置于空气中，如图4(a)所示)，铜导环的上下及内外表面对流换热系数均为5W/m²。对铜导环施加电流i＝10⁴sin(50πt)，持续时间0.1s，分析铜环的温度变化。

建立轴对称模型，如图5(b)所示，电磁场计算区域包含铜环、空气区域，采用三角划分网格，共988个节点，2061个单元。温度场计算区域为铜环区域，也采用三角形划分网格，共121个节点，200个单元。

当式(15)中电磁场载荷离散误差取为e_tolerance＝0.5％，式(8)、式(21)中取γ＝1％，β＝1。电磁-温度场前三个时间步内采用1ms耦合，随后采用自适应时间步长耦合。最后将自适应时间步长计算结果与定步长计算结果对比，如图5(a)、图5(b)所示，分别为铜导环施加交流载荷0.1s时，采用等步长和自适应时间步长计算所得的温度分布云图。

由图6(a)、图6(b)可见，采用自适应时间步长所得的铜环温度分布与定步长温度分布基本一致。其中1号单元温度最高而相对误差最小为0.35％，140号单元温度最低而相对误差最大为0.53％。为此选取1号单元和140号单元，分析其在整个加热过程中温度及误差的变化规律。

1号和140号在0.1s内自适应时间步长与定步长计算所得单元温度对比如图6(a)、图6(b)所示，铜导块温度近似以0.02s为周期上升。以0.06s～0.08s为例，在区域a和c中，电流幅值较小，温度场升较慢，对应的温度场计算步长较大。而区域b，电流幅值较大，温度上升较快，对应温度场计算时间步长较小，说明电磁-温度场自动步长耦合算法的有效性。在整个计算过程中，自动步长计算与定步长计算误差在在0.7％以内，明了该方法的准确性。

电磁-温度场自适应时间步长与定步长耦合计算性能对比如表3所示。在0.1s内，采用定步长耦合，选取时间步长为1ms时，电磁场和温度场均需计算100次，整体计算时间为235s。而采用自适应步长求解时，电磁场计算93次，温度场计算29次，计算时间为184s，计算时间减小21％。

表3计算性能对比

Claims (1)

1.一种基于自适应时间步长的瞬变电磁-温度场耦合计算方法，其特征在于包括以下步骤：

步骤一、电磁场与温度场的耦合时间点确定：采用温度场触发热量来判断电磁-温度场耦合时间节点，首先电磁-温度场步长耦合三次，用来获取温度场触发热量Q_pre的预测数据，然后采用指数平滑法预测温度场触发热量，通过触发热量判断耦合时间节点；

步骤二、电磁场自适应步长计算：通过载荷离散误差和响应特征值确定电磁场载荷最佳离散的步长△t_n ^E，当电磁场累积热量达到温度场触发热量时，暂停电磁场计算，计算温度场；

步骤三、温度场自适应步长计算：通过载荷离散误差和响应特征值确定电磁场载荷最佳离散的步长△t_n ^T，将电磁场计算平均热功率作为载荷加至温度场，当温度场计算时间与电磁场同步时，暂停温度场计算，更新电磁场节点温度，并进行下一时间步的温度场触发热量计算，反复迭代计算至最终时间；

电磁场与温度场的耦合时间点确定包括以下步骤：

步骤1)：采用温度场触发热量来判断电磁场与温度场进行耦合计算时的时间节点，首先电磁-温度场步长耦合三次，用来获取温度场触发热量Q_pre的预测数据；

步骤2)：由式(1)、(2)、(3)计算电磁场计算过程中允许材料变化最大温度；其中线性材料的热功率P与电流密度J的函数关系如式(1)、(2)：

P_n＝∫_VJ_n ²ε(T)dV (1)

ε＝aT+ε₀ (2)

式中：P_n为t_n时刻发热功率；V为单元体积；J_n为T_n时刻电流密度；ε为电阻率，a为电阻随温度变化率，ε(T)为电阻率与温度函数，ε₀为0℃时电阻率；∫_V dV表示体积分；

当输入热功率为P_n，由温度变化引起功率计算误差应满足：

|P_n(T_n)-P_n(T_n+△T)|≤γP_n(T_n) (3)

式中：P_n(T_n)为功率与温度函数，如式(1)；γ为由温度变化引起的功率误差，取1％；

式(3)取等号时，可得出t_n时刻温度允许变化最大值△T_max；

步骤3)：温度场触发热量计算，根据t_n时刻最大温度变化△T_max，由t_n，t_n-1，t_n-2时刻输入功率及温度变化，预测t_n+1温度场触发热量

计算t_n，t_n-1，t_n-2时刻温度随热量变化率k_n-1，k_n-2，k_n-3，如式(4)：

式中：Q_n，Q_n-1，Q_n-2分别为t_n，t_n-1，t_n-2时刻发热量；T_n，T_n-1，T_n-2，T_n-3为t_n，t_n-1，t_n-2，t_n-3时刻温度；

最后计算温度场触发热量计算，采用如式(5)的指数平滑方法预测t_n+1时刻温度场变化率K_n+1：

K_n+1＝αK_n+(1-α)K_n-1+(1-α)²K_n-2 (5)

式中：α为指数平滑系数，α＝0.75；k_n-1，k_n-2，k_n-3如式(5)；

温度场触发热量为：

Q^pre＝△T_maxK_n+1 (6)

式中：Q^pre为温度场触发热量；△T_max为允许变化最大温度，见式(4)；K_n+1计算见式(5)；

电磁场自适应步长计算包括以下步骤：

步骤4)：根据载荷离散误差确定电磁场时间步长

t∈(t_n-1,t_n)时，当载荷矩阵P采用线性插值时，等效载荷采用斜坡加载，由载荷离散产生的误差为如式(7)，

式中：

为第n个时间步长的载荷离散误差；f_s(t)为电磁场中电流或电压载荷；△t为离散步长；

对式(7)采用梯形公式积分，得离散误差近似如式(8)；

式中：f_s(ξ)为电磁场中电流或电压载荷；ξ为t_n-1至t_n间常数；∝表示正比于，其余参数如式(7)；

根据式(8)，载荷离散误差近似正比于(△t)²，可将下一步长计算可分为如下两步：

(a)步长预测：采用式(9)根据第n步计算误差，预测第n+1步长△t^m _n+1；

式中：

为第n+1步预测步长；

为安全系数，

e_tolerance为允许最大误差；

为第n步产生的离散误差；

(b)步长校正：判断当第n+1时间步长△t¹ _n+1所产生误差是否满足e^k+1<e_tolerance，如不满足采用(8)进行修正迭代计算，直至满足e^k<e_tolerance；

式中：△t¹ _n+1为第n+1个时间步长；

为安全系数，

e_tolerance为允许最大误差；e^k+1为第k+1次迭代产生的离散误差；

步骤5)：根据响应特征值确定电磁场时间步长

采用响应特征λ_r值确定稳定时间步长△t_n+1；定义△t_n+1λ为震荡限制条件，当△t_n+1λ>>1时系统处于震荡状态，为保证计算稳定性可取最大步长需满足：

式中：f<1，f为稳定系数；λ_r为响应特征值；u为场量，电磁场中

A为矢量磁位，V为电压；u_n为t_n-1到t_n时间段场量u的变化；K为一阶有限元方程中传导矩阵；C为一阶有限元中容性矩阵；

步骤6)：电磁场在t_n时刻离散时间步长为：

式中：

为电磁场t_n时刻离散时间步长；

为式(10)中载荷离散最大步长；

为特征值确定离散步长；

步骤7)：耦合时间判断：当t_n+1时刻电磁场累积热量满足以下条件之一时，t_n+1为电磁-温度场耦合时间点，暂停电磁场计算，启动温度场计算；

(a)：单元最大热量变化达到步骤5)中单元预测热量阈值时自动启动温度场计算：

式中：Qⁱ为第i个单元的累积热量；

为第i个单元的预测触发热量，由式(6)计算；

(b)：总体热量变化达到步骤5)中单元预测热量阈值时自动启动温度场计算：

式中：Q^pre为温度场总体触发热量，由式(6)计算；Qⁱ为第i个单元的累积热量；β为触发安全系数，β＝0.9；

温度场自适应步长计算包括以下步骤：

步骤8)：计算温度场自适应时间步长

温度场自适应时间步长计算同步骤4)、步骤5)、步骤6)；

步骤9)：当温度场计算时间到达步骤7)中耦合时间T_n+1时，停止温度场计算，并更新电磁场节点温度；

步骤10)：反复迭代步骤2)～步骤8)，直至计算总时间T_total。

Images