斯蒂芬·玻尔兹曼定律 — 辐射传热基础与CAE实现

分类: 熱解析 > 輻射 | 综合版 2026-04-12
Stefan-Boltzmann law: blackbody spectral radiance curves at multiple temperatures showing T⁴ dependence of total emissive power
ステファン・ボルツマンの法則 — 黒体放射の全放射エネルギーは絶対温度の4乗に比例する

理论与物理

概述 — 为什么是4次方

🧑‍🎓

斯特藩-玻尔兹曼定律说的是“与温度的4次方成正比”对吧。为什么不是3次方或5次方，而正好是4次方呢？有什么物理上的原因吗？

🎓

问得好。结论是，将普朗克辐射定律在全波长上积分，自然就会得到T⁴。

普朗克的光谱辐射亮度是

$$ B(\lambda, T) = \frac{2hc^2}{\lambda^5}\,\frac{1}{e^{hc/(\lambda k_B T)} - 1} $$

将其在全波长 $\lambda = 0 \sim \infty$ 上积分，求得总辐射能量 $E_b$。使用变量变换 $x = hc/(\lambda k_B T)$，积分变为

$$ E_b = \frac{2\pi k_B^4}{h^3 c^2}\,T^4 \int_0^\infty \frac{x^3}{e^x - 1}\,dx $$

这个定积分是一个确定值 $\Gamma(4)\,\zeta(4) = 6 \times \pi^4/90 = \pi^4/15$。所以得出 $E_b \propto T^4$。4次方来源于“光子在三维空间中的状态密度与 $\nu^2$ 成正比”以及“玻色-爱因斯坦分布的能量平均”，是空间维度和量子统计的结果。

🧑‍🎓

原来如此，不是单纯的实验公式，而是可以从量子力学推导出来的啊。那么，斯特藩-玻尔兹曼常数 $\sigma$ 的值也能从理论上确定吗？

🎓

是的。$\sigma$ 可以用基础物理常数表示：

$$ \sigma = \frac{2\pi^5 k_B^4}{15\,h^3\,c^2} = 5.670374419 \times 10^{-8} \;\text{W/(m}^2\text{K}^4\text{)} $$

2019年SI单位重新定义之后，$k_B$, $h$, $c$ 都成为了定义值，所以 $\sigma$ 也成为了精确确定的值。在CAE教材中，通常使用四舍五入后的 $5.67 \times 10^{-8}$，但了解这一点在精度讨论时会很有用。

控制方程

斯特藩-玻尔兹曼定律的基本形式如下。

黑体的总辐射能量
$$ E_b = \sigma\,T^4 $$

其中 $E_b$ 是单位面积的总辐射能量 [W/m²]，$T$ 是绝对温度 [K]，$\sigma = 5.67 \times 10^{-8}$ W/(m²K⁴) 是斯特藩-玻尔兹曼常数。

实际物体并非黑体，因此引入发射率（emissivity） $\varepsilon$：

$$ E = \varepsilon\,\sigma\,T^4, \qquad 0 \le \varepsilon \le 1 $$

以下列出典型的发射率。

材料温度范围发射率 ε
氧化钢500〜1200 °C0.70〜0.85
抛光铝20〜500 °C0.04〜0.08
耐火砖500〜1000 °C0.75〜0.93
黑体涂料全温度范围0.94〜0.97
玻璃（平板玻璃）20〜300 °C0.90〜0.95
MLI（多层隔热材料）-200〜200 °C0.01〜0.03（有效值）

灰体的辐射交换

🧑‍🎓

在实际的CAE中，物体之间热量是“相互传递”的吧？那种情况怎么计算呢？

🎓

两个表面之间的净辐射热交换，由两个表面的温差决定：

$$ q_{1\to2} = \varepsilon\,\sigma\,A\,(T_1^4 - T_2^4) $$

这是灰体的基本公式。这里重要的是 $T^4$ 的差值，而不是 $(T_1 - T_2)^4$。这种非线性使得CAE变得棘手。

例如，计算铸坯（$T_1 = 1200°C = 1473\,\text{K}$）与周围环境（$T_2 = 30°C = 303\,\text{K}$）的情况。

$$ q = 0.8 \times 5.67 \times 10^{-8} \times (1473^4 - 303^4) \approx 213\,\text{kW/m}^2 $$

另一方面，相同条件下的自然对流，取 $h \approx 10$ W/(m²K)，则 $q_{\text{conv}} = 10 \times 1170 \approx 12$ kW/m²。也就是说辐射大约是对流的18倍。在1200°C的世界里，辐射是绝对主导的。

🧑‍🎓

18倍！那么在像钢铁连铸这样的高温工艺中，如果忽略了辐射，凝固壳的厚度预测就会完全不准了吧。

🎓

没错。在连铸中，铸坯离开结晶器后，其表面立即发生辐射冷却，凝固壳从外侧开始生长。如果不能准确建模 $\sigma\varepsilon T^4$，凝固壳厚度的预测就会出错，进而导致漏钢（Breakout）的风险。这是钢铁厂最害怕的事故之一。

多表面间的辐射交换与角系数

封闭空间（腔体）内的多表面间辐射交换，使用角系数（View Factor） $F_{ij}$ 来描述。

$$ q_i = \sum_{j=1}^{N} A_i\,F_{ij}\,\sigma\,(T_i^4 - T_j^4) \qquad (\text{黑体情况}) $$

对于灰体、漫射表面，则使用辐射度法（Radiosity Method）：

$$ J_i = \varepsilon_i\,\sigma\,T_i^4 + (1-\varepsilon_i)\sum_{j=1}^{N} F_{ij}\,J_j $$

其中 $J_i$ 是表面 $i$ 的辐射度（辐射+反射能量的总和）。这构成了一个N元联立方程组，在CAE求解器中通过矩阵运算求解。角系数具有以下重要性质：

  • 互易性： $A_i F_{ij} = A_j F_{ji}$
  • 完整性： $\sum_{j=1}^{N} F_{ij} = 1$（封闭腔体情况）
  • 凸面： $F_{ii} = 0$（凸面看不到自身）

辐射主导条件 — 何时辐射成为主导

🧑‍🎓

反过来说，温度低的话辐射就可以忽略吗？大概多高的温度是分界线呢？

🎓

一个实用的经验法则是表面温度600°C（约900K）。超过这个温度，辐射的贡献会急剧增加。T⁴的威力是惊人的，温度变为2倍，辐射能量就变为16倍。

但是，在真空中对流为零，所以对于航天器的热设计，即使在常温下，辐射也是唯一的散热手段。国际空间站（ISS）的辐射器在大约280K左右运行，仅靠 $\sigma T^4$ 来决定散热量。正因为温度低，所以需要增大面积。

辐射 vs 对流的支配领域
温度范围辐射的贡献代表性的工程系统
< 200 °C通常很小（5〜15%）电子设备冷却、建筑空调
200〜600 °C不可忽视（20〜50%）排气管、热处理炉的预热区
600〜1000 °C主导（50〜80%）玻璃成型、陶瓷烧成
> 1000 °C压倒性（80〜95%）钢铁制造（连铸、轧制）、火箭喷管
真空环境100%（唯一的散热手段）人造卫星、空间站

数值解法与实现

T⁴项的线性化

🧑‍🎓

T⁴是强非线性的吧。在FEM求解器中是怎么处理的呢？

🎓

有两种方法。最常用的是转换为辐射传热系数。

$$ q = \varepsilon\sigma(T_s^4 - T_\infty^4) = h_{\text{rad}}(T_s - T_\infty) $$

其中辐射传热系数 $h_{\text{rad}}$ 为

$$ h_{\text{rad}} = \varepsilon\sigma(T_s^2 + T_\infty^2)(T_s + T_\infty) $$

$h_{\text{rad}}$ 使用上一次迭代的温度计算，并通过牛顿-拉夫逊法更新。因为它可以与对流的 $h_{\text{conv}}$ 相加，所以很容易纳入现有的传热系数框架。

🧑‍🎓

原来如此，是线性化后迭代求解啊。那另一种方法呢？

🎓

另一种是将T⁴的微分直接纳入雅可比矩阵的方法。

$$ \frac{\partial q}{\partial T} = 4\,\varepsilon\,\sigma\,T^3 $$

将其放入NR法的切线刚度矩阵中。Abaqus的*RADIATION或Ansys的SRDOPT在内部使用了这种方法。收敛速度快，但如果初始温度估计不好，容易发散。

有限元法中的辐射边界条件

将辐射边界条件纳入热传导方程的弱形式。能量方程为：

$$ \rho c_p \frac{\partial T}{\partial t} = \nabla \cdot (k \nabla T) + \dot{Q} $$

辐射面 $\Gamma_r$ 上的边界条件：

$$ -k\,\frac{\partial T}{\partial n}\bigg|_{\Gamma_r} = \varepsilon\,\sigma\,(T^4 - T_\infty^4) + h_{\text{conv}}\,(T - T_\infty) $$

通过伽辽金法纳入弱形式后，在单元级别产生以下非线性残差向量：

$$ \mathbf{R}_e^{\text{rad}} = \int_{\Gamma_r^e} \varepsilon\,\sigma\,(T^4 - T_\infty^4)\,\mathbf{N}^T\,d\Gamma $$

其中 $\mathbf{N}$ 是形函数向量。切线矩阵为：

$$ \mathbf{K}_e^{\text{rad}} = \frac{\partial \mathbf{R}_e^{\text{rad}}}{\partial \mathbf{T}_e} = \int_{\Gamma_r^e} 4\,\varepsilon\,\sigma\,T^3\,\mathbf{N}^T\mathbf{N}\,d\Gamma $$

耦合求解策略

🧑‍🎓

同时处理对流和辐射时，在构建分析方面有什么技巧吗？

🎓

主要有三种策略。

  1. 弱耦合（Sequential）: CFD计算对流系数 → FEM求解热传导+辐射 → 将温度反馈给CFD。实现简单但收敛慢。
  2. 强耦合（Monolithic）: 将流体、固体、辐射汇总到一个方程组中同时求解。Ansys Fluent的S2S（面到面）模型就是这种方式。收敛快但内存消耗大。
  3. DO/MC模型: 当存在参与性介质（如燃烧气体）时，使用离散坐标法或蒙特卡洛法求解辐射传递方程（RTE）。S-B定律作为壁面边界条件使用。

实践中，“先在没有辐射的情况下收敛流场，然后再启用辐射模型”比较稳定。如果一开始就全部启用，有发散的风险。

实现时的注意点：绝对温度

S-B定律中的 $T$ 必须是绝对温度 [K]。如果CAE软件的温度单位仍然是 °C，那么 $T^4$ 的计算从根本上就是错误的。例如，本想输入1000°C，却输入 $T=1000$，那么与正确值 $(1273)^4 = 2.62 \times 10^{12}$ 相比，$(1000)^4 = 1.00 \times 10^{12}$ 会导致辐射能量被低估约2.6倍。

实践指南

分析流程

🧑‍🎓

实际进行辐射分析时，应该按什么步骤推进呢？

🎓

基本流程如下：

  1. 识别辐射面： 找出高温面以及它“能看到”的面。同时确认是否存在遮挡。
  2. 设置发射率： 确认是否有温度依赖性。氧化状态会导致其大幅变化，不要盲目相信文献值。
  3. 计算角系数： 使用求解器的自动计算（如半立方体法），或者如果有解析解，可用于验证。
  4. 网格考虑： 辐射面存在T⁴的非线性，因此需要细化温度梯度陡峭的部分（边缘、角落）。
  5. 设置初始温度： 将接近实际运行温度的值设为初始条件。如果从300K开始，要达到1500K，迭代次数会爆炸式增长。
  6. 分阶段求解： 按对流 → 对流+辐射的顺序启用。不要一开始就全部启用。

案例1: 连铸的辐射冷却

🧑‍🎓

刚才提到的连铸话题，能再具体讲讲吗？实际的CAE分析中会构建什么样的模型呢？

🎓

在连铸中，从结晶器拉出的铸坯会通过二次冷却区。冷却模式分为三个区域：

  1. 喷淋冷却区： 水喷淋引起的强制对流（$h = 500\sim2000$ W/(m²K)）占主导。这里S-B辐射约占整体的10〜20%。
  2. 辐射冷却区： 喷淋结束后，铸坯温度仍有900〜1100°C。在这个区间，$q_{\text{rad}} = \varepsilon\sigma(T_s^4 - T_\infty^4)$ 承担了80%以上的冷却。
  3. 辊子接触区： 与导向辊的接触传导。接触面积小，但在局部很重要。

在CAE模型中，通常采用沿铸造方向的拉格朗日追踪（切片法）。切出铸坯的截面2D模型，按照拉坯速度使其通过各个冷却区。凝固潜热 $L = 270$ kJ/kg 左右通过焓法处理。

🧑‍🎓

在设置发射率时有什么需要注意的吗？钢材的表面状态不同，差别很大吧？

🎓

很敏锐的观察。铸坯表面的氧化皮（氧化皮膜）生长会导致发射率随时间变化。结晶器正下方 $\varepsilon \approx 0.4$（有金属光泽），但在二次冷却区出口，氧化进展后 $\varepsilon \approx 0.8$。在工程实践中，最佳做法是用温度依赖表格来输入这种变化。如果坚持使用固定值 $\varepsilon = 0.8$，会高估上游的冷却速度，导致凝固壳厚度预测出错。

案例2: 航天器的热设计

在太空中，由于对流为零，散热完全依赖于S-B定律的辐射。设计的基本公式是：

$$ Q_{\text{dissipation}} = \varepsilon\,\sigma\,A_{\text{rad}}\,(T_{\text{rad}}^4 - T_{\text{space}}^4) $$

太空的辐射温度 $T_{\text{space}} \approx 3\,\text{K}$（宇宙微波背景辐射）可以忽略不计，因此实际上仅靠 $\varepsilon\sigma A T_{\text{rad}}^4$ 进行设计。国际空间站的辐射器以 $\varepsilon \approx 0.9$，$T_{\text{rad}} \approx 280\,\text{K}$ 运行。$280^4 = 6.15 \times 10^9$ 乘以 $\sigma$ 得到约 $350\,\text{W/m}^2$ — 这就是辐射器每平方米的散热能力。

根据太阳入射热（$\alpha_s \cdot S_0 \approx 0.3 \times 1361 = 408\,\text{W/m}^2$）的热平衡，决定辐射器面积和运行温度。这里 $\alpha_s / \varepsilon$ 比（太阳光吸收率与红外发射率之比）成为热设计的最重要参数。

案例3: 工业炉的炉壁设计

🧑‍🎓

工业炉的设计中也使用S-B定律吗？

🎓

工业炉是辐射传热的宝库。在超过1000°C的炉内，工件（加热对象）获得的热量有70〜90%来自辐射。为了对炉壁、加热器、工件之间的多表面辐射交换进行建模，会使用区域法（Hottel的区域模型）或辐射度法。

特别是在燃气炉中，燃烧气体（CO₂, H₂O）本身参与辐射，仅靠S-B定律是不够的。需要使用加权灰气体和模型（WSGG）计算气体的发射率，并求解辐射传递方程。Ansys Fluent的DO模型或STAR-CCM+的DO/S2S模型支持此功能。

软件比较

主要CAE工具的辐射功能

工具辐射模型角系数计算参与性介质特点
Ansys Mechanical

🔧 相关模拟器

通过交互式参数体验这一理论 → 模拟器集