3. 瞬态热传导中的时间步长选择与验证
问题设定与解析基准
考虑一根长 L = 0.1 m 的均质一维导热杆，初始温度均匀为 T₀ = 20 °C，两端突然施加 T_b = 100 °C 的恒温边界。材料取铜的典型教学参数：密度 ρ = 8960 kg/m³、比热 c_p = 385 J/(kg·K)、导热系数 k = 400 W/(m·K)，则热扩散率 α = k/(ρc_p) ≈ 1.16×10⁻⁴ m²/s。
该问题的解析解为无穷级数：
T(x,t) = T_b + (T₀ - T_b) · Σ_{n=1,3,5,...} (4/(nπ)) · sin(nπx/L) · exp(-n²π²αt/L²)特征时间 τ = L²/(π²α) ≈ 8.7 s，可作为时间步长选择的参考尺度。下面用 COMSOL 的“传热传质 > 固体传热”接口建立一维模型，并对比不同时间步长下的数值解。
建模要点与时间步长设置
- 几何：一维线段，长度 0.1 m；网格用“极细化”映射网格，约 200 个单元，确保空间离散误差远小于时间离散误差。
- 边界：左端
x=0与右端x=L均设为温度边界T = 100 °C。 - 初始值：
T = 20 °C。 - 研究类型：瞬态，时间范围
0 ~ 30 s。 - 时间步长策略：分别采用
Δt = 2 s（粗）、Δt = 0.5 s（中）、Δt = 0.1 s（细），并开启“严格时间步进”以排除自适应步长的干扰。
不同时间步长的温度演化对比
取杆中点 x = L/2 作为监测点，记录其温度随时间的变化。下图给出三种步长下的数值曲线与解析解的对比关系，以及误差随步长缩小的收敛趋势。
误差分析与步长选择依据
从图中可以读出三个关键现象：
Δt = 2 s时，早期温度上升段明显滞后于解析解，最大误差约 8.5 °C，占总温升的 10% 以上，不可接受。Δt = 0.5 s时，曲线已接近解析解，最大误差约 1.2 °C，约占总温升的 1.5%。Δt = 0.1 s时，数值解与解析解几乎重合，最大误差约 0.05 °C，远小于工程允许的 0.5 °C。
误差随步长缩小的比例约为 7 倍（步长缩小 4 倍），接近二阶收敛，说明 COMSOL 在该问题中使用了 BDF2 格式。若以“最大误差 ≤ 0.5 °C”为验收标准，Δt = 0.5 s 已满足要求；若需更严格的 0.1 °C 标准，则应选 Δt = 0.1 s。
步长选择经验规则：对于扩散型瞬态问题，建议初始时间步长取 Δt ≤ τ/20，其中 τ = L²/(π²α) 为特征时间。本例 τ ≈ 8.7 s，则 Δt ≤ 0.44 s，与 Δt = 0.5 s 的结论一致。若关注早期快速变化阶段，可进一步取 Δt ≤ τ/50。
验证方法与操作检查表
| 检查项 | 操作 | 通过判据 |
|---|---|---|
| 空间离散独立性 | 将网格从 200 单元加密到 400 单元，重算 Δt = 0.1 s | 中点温度变化 < 0.01 °C |
| 时间步长收敛性 | 依次用 2 s、0.5 s、0.1 s 计算，记录最大误差 | 误差随步长缩小而单调下降 |
| 能量守恒检查 | 在结果中计算域内总热量随时间的变化 | 与边界流入热量之差 < 1% |
| 解析解对比 | 用上述级数公式取前 20 项，在相同时间点采样 | 细步长解与解析解偏差 < 0.1 °C |
| 求解器日志 | 查看“时间步进”日志中的步长与迭代次数 | 无步长被反复截断的警告 |
效率与精度的权衡结论
在本例中，Δt = 0.5 s 是兼顾精度与效率的推荐值：计算 30 s 物理时间仅需 60 个时间步，单次求解在普通笔记本上约数秒完成，而最大误差已控制在 1.2 °C 以内。若后续需要耦合其他物理场或进行参数扫描，可先用 Δt = 0.5 s 快速筛选，再对关键工况用 Δt = 0.1 s 复核。这样既避免了 Δt = 2 s 的严重失真，也避免了全程使用 Δt = 0.1 s 带来的不必要计算量。