第六十一篇：热辐射基本定律

摘要

热辐射是热量传递的三种基本方式之一，在航天、能源、冶金、建筑等领域具有重要应用。本教程系统介绍热辐射的基本概念和定律，包括热辐射的本质、电磁波谱特性、黑体辐射模型以及辐射换热的理论基础。通过Python实现黑体辐射光谱计算、辐射能量分布可视化以及不同温度下辐射特性的对比分析。本教程涵盖热辐射的量子理论基础、辐射强度与辐射力的关系、以及辐射换热的基本计算方法。通过本教程的学习，读者将掌握热辐射的基本理论，能够进行辐射换热的定量计算，为后续深入学习辐射传热奠定基础。

关键词

热辐射，黑体辐射，电磁波谱，辐射强度，辐射力，普朗克定律，斯特藩-玻尔兹曼定律，维恩定律

---

1. 热辐射基础理论

1.1 热辐射的本质

热辐射是物体由于温度而发射电磁波的能量传递方式。与导热和对流不同，热辐射不需要介质，可以在真空中传播。

热辐射的特点：

1. 无需介质：通过电磁波传播，可在真空中进行
2. 能量转换：热能 ↔ 辐射能 ↔ 热能
3. 双向传递：高温物体向低温物体辐射，同时低温物体也向高温物体辐射
4. 光谱分布：辐射能量按波长分布，与温度密切相关

电磁波谱：
热辐射主要位于红外波段，波长范围通常为0.1-100 μm：

• 紫外线：0.01-0.4 μm
• 可见光：0.4-0.76 μm
• 红外线：0.76-100 μm
• 远红外：4-100 μm

1.2 黑体辐射模型

黑体是理想化的辐射体，具有以下特性：

• 吸收所有入射辐射（吸收比 α = 1）
• 在相同温度下发射最大辐射能量
• 辐射特性仅与温度有关，与材料无关

黑体辐射的物理基础：
黑体辐射是量子力学的起源之一。1900年，普朗克提出能量量子化假设，成功解释了黑体辐射光谱：

$$E=h\nu =\frac{hc}{\lambda }$$

其中 $h=6.626\times 10^{-34}$ J·s 是普朗克常数。

1.3 辐射强度与辐射力

辐射强度（Radiation Intensity）：
辐射强度 $I$ 定义为单位面积、单位立体角、单位波长的辐射功率：

$$I(\theta ,\varphi ,\lambda )=\frac{dQ}{dA\cdot d\Omega \cdot d\lambda }$$

辐射力（Emissive Power）：
辐射力 $E$ 是单位面积发射的总辐射功率：

$$E={\int }_0^{2\pi }{\int }_0^{\pi /2}I(\theta ,\varphi )\cos \theta \sin \theta d\theta d\varphi$$

对于漫射表面（各向同性辐射）：

$$E=\pi I$$

---

2. 黑体辐射定律

2.1 普朗克定律

普朗克定律描述了黑体辐射的光谱分布：

$$E_{b\lambda }(\lambda ,T)=\frac{2\pi h{c}^2}{{\lambda }^5}\cdot \frac{1}{\exp \left(\frac{hc}{\lambda k_{B}T}\right)-1}$$

或写成：

$$E_{b\lambda }(\lambda ,T)=\frac{C_1}{{\lambda }^5}\cdot \frac{1}{\exp \left(\frac{C_2}{\lambda T}\right)-1}$$

其中：

• $C_1=2\pi h{c}^2=3.742\times 10^8$ W·μm⁴/m²
• $C_2=hc/k_B=1.439\times 10^4$ μm·K

2.2 维恩位移定律

维恩位移定律给出了峰值波长与温度的关系：

KaTeX parse error: Undefined control sequence: \cdotp at position 1: \̲c̲d̲o̲t̲p̲

这表明温度越高，峰值波长越短（向短波方向移动）。

2.3 斯特藩-玻尔兹曼定律

斯特藩-玻尔兹曼定律给出了黑体的总辐射力：

$$E_b=\sigma T^4$$

其中 $\sigma =5.67\times 10^{-8}$ W/(m²·K⁴) 是斯特藩-玻尔兹曼常数。

2.4 兰贝特定律

兰贝特定律描述了漫射表面的方向辐射特性：

$$I(\theta )=I_n=\text{常数}$$

即漫射表面在各个方向的辐射强度相等。

---

3. Python仿真实现

3.1 黑体辐射光谱计算

"""
主题061：热辐射基本定律
黑体辐射光谱计算与可视化
"""
import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate
import os

# 创建输出目录
output_dir = r'd:\文档\500仿真领域\工程仿真\传热学仿真\主题061'
os.makedirs(output_dir, exist_ok=True)

# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei', 'DejaVu Sans']
plt.rcParams['axes.unicode_minus'] = False

# 物理常数
h = 6.626e-34      # 普朗克常数 (J·s)
c = 2.998e8        # 光速 (m/s)
k_B = 1.381e-23    # 玻尔兹曼常数 (J/K)
sigma = 5.67e-8    # 斯特藩-玻尔兹曼常数 (W/m²·K⁴)

# 第一和第二辐射常数
C1 = 2 * np.pi * h * c**2 * 1e24  # W·μm⁴/m²
C2 = h * c / k_B * 1e6            # μm·K

print(f"第一辐射常数 C1 = {C1:.3e} W·μm⁴/m²")
print(f"第二辐射常数 C2 = {C2:.3e} μm·K")

def planck_law(lambda_um, T):
    """
  普朗克黑体辐射定律
    lambda_um: 波长 (μm)
    T: 温度 (K)
    返回: 光谱辐射力 (W/m²·μm)
    """
    lambda_m = lambda_um * 1e-6  # 转换为米
    
    # 普朗克公式
    term1 = (2 * np.pi * h * c**2) / (lambda_m**5)
    term2 = 1 / (np.exp(h*c/(lambda_m*k_B*T)) - 1)
    
    E_lambda = term1 * term2 * 1e-6  # 转换为 W/m²·μm
    return E_lambda

def wien_law(T):
    """维恩位移定律，返回峰值波长 (μm)"""
    return 2898 / T

def stefan_boltzmann(T):
    """斯特藩-玻尔兹曼定律，返回总辐射力 (W/m²)"""
    return sigma * T**4

# 波长范围
lambda_range = np.linspace(0.1, 50, 1000)  # 0.1-50 μm

# 不同温度下的黑体辐射光谱
T_list = [300, 500, 800, 1200, 2000]  # K
colors = ['blue', 'green', 'orange', 'red', 'purple']

fig, axes = plt.subplots(2, 2, figsize=(14, 12))

# 1. 不同温度下的辐射光谱
for T, color in zip(T_list, colors):
    E_spectrum = [planck_law(lam, T) for lam in lambda_range]
    axes[0, 0].plot(lambda_range, E_spectrum, color=color, linewidth=2, label=f'T = {T} K')
    
    # 标记峰值波长
    lambda_max = wien_law(T)
    axes[0, 0].axvline(x=lambda_max, color=color, linestyle='--', alpha=0.5)

axes[0, 0].set_xlabel('Wavelength (μm)', fontsize=11)
axes[0, 0].set_ylabel('Spectral Emissive Power (W/m²·μm)', fontsize=11)
axes[0, 0].set_title('Blackbody Radiation Spectrum', fontsize=12)
axes[0, 0].legend()
axes[0, 0].grid(True, alpha=0.3)
axes[0, 0].set_xlim([0, 20])

# 2. 维恩位移定律验证
T_range = np.linspace(300, 3000, 100)
lambda_max_range = [wien_law(T) for T in T_range]

axes[0, 1].plot(T_range, lambda_max_range, 'b-', linewidth=2)
axes[0, 1].set_xlabel('Temperature (K)', fontsize=11)
axes[0, 1].set_ylabel('Peak Wavelength (μm)', fontsize=11)
axes[0, 1].set_title('Wien Displacement Law', fontsize=12)
axes[0, 1].grid(True, alpha=0.3)

# 添加一些标记点
for T in [300, 5800, 1000]:
    lam_max = wien_law(T)
    axes[0, 1].plot(T, lam_max, 'ro', markersize=8)
    if T == 5800:
        axes[0, 1].annotate(f'Sun (5800K)\nλ_max={lam_max:.2f}μm', 
                           xy=(T, lam_max), xytext=(T-1000, lam_max+2),
                           arrowprops=dict(arrowstyle='->', color='red'))

# 3. 斯特藩-玻尔兹曼定律
T_range_sb = np.linspace(200, 2500, 100)
E_total = [stefan_boltzmann(T) for T in T_range_sb]

axes[1, 0].semilogy(T_range_sb, E_total, 'b-', linewidth=2)
axes[1, 0].set_xlabel('Temperature (K)', fontsize=11)
axes[1, 0].set_ylabel('Total Emissive Power (W/m²)', fontsize=11)
axes[1, 0].set_title('Stefan-Boltzmann Law', fontsize=12)
axes[1, 0].grid(True, alpha=0.3)

# 4. 可见光范围内的辐射能量比例
T_vis_range = np.linspace(500, 6000, 100)
visible_fraction = []

for T in T_vis_range:
    # 可见光范围: 0.4-0.76 μm
    E_vis, _ = integrate.quad(lambda lam: planck_law(lam, T), 0.4, 0.76)
    E_total_T = stefan_boltzmann(T)
    visible_fraction.append(E_vis / E_total_T * 100)

axes[1, 1].plot(T_vis_range, visible_fraction, 'b-', linewidth=2)
axes[1, 1].axvline(x=5800, color='r', linestyle='--', label='Sun (5800K)')
axes[1, 1].set_xlabel('Temperature (K)', fontsize=11)
axes[1, 1].set_ylabel('Visible Light Fraction (%)', fontsize=11)
axes[1, 1].set_title('Visible Light Energy Fraction', fontsize=12)
axes[1, 1].legend()
axes[1, 1].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig(f'{output_dir}/blackbody_radiation.png', dpi=150, bbox_inches='tight')
plt.close()
print("黑体辐射光谱图已保存")

3.2 辐射强度与方向分布

"""
辐射强度方向分布与兰贝特定律验证
"""

# 方向角
theta = np.linspace(0, np.pi/2, 100)

# 漫射表面（兰贝特表面）
I_diffuse = np.ones_like(theta)  # 各方向强度相等

# 镜面反射表面
I_specular = np.zeros_like(theta)
I_specular[0] = 1.0  # 只在反射角方向

# 实际表面（混合）
n_rough = 2  # 粗糙度指数
I_real = np.cos(theta)**n_rough  # 余弦分布

fig, axes = plt.subplots(1, 2, figsize=(14, 6))

# 极坐标图
ax1 = plt.subplot(121, projection='polar')
ax1.plot(theta, I_diffuse, 'b-', linewidth=2, label='Diffuse (Lambertian)')
ax1.plot(theta, I_real, 'r--', linewidth=2, label=f'Real surface (n={n_rough})')
ax1.set_theta_zero_location('N')
ax1.set_theta_direction(-1)
ax1.set_thetamax(90)
ax1.set_title('Directional Radiation Intensity', fontsize=12, pad=20)
ax1.legend(loc='upper right', bbox_to_anchor=(1.3, 1.0))
ax1.grid(True)

# 辐射力计算验证
ax2 = plt.subplot(122)

# 不同表面类型的辐射力
E_diffuse = np.pi * np.ones_like(theta)  # 漫射表面
E_real = np.pi * I_real  # 实际表面

ax2.plot(theta * 180/np.pi, I_diffuse, 'b-', linewidth=2, label='Intensity (Diffuse)')
ax2.plot(theta * 180/np.pi, I_diffuse * np.cos(theta), 'b--', linewidth=2, label='Projected Intensity')
ax2.plot(theta * 180/np.pi, I_real, 'r-', linewidth=2, label=f'Intensity (Real, n={n_rough})')
ax2.plot(theta * 180/np.pi, I_real * np.cos(theta), 'r--', linewidth=2, label='Projected Intensity')

ax2.set_xlabel('Polar Angle (degrees)', fontsize=11)
ax2.set_ylabel('Normalized Intensity', fontsize=11)
ax2.set_title('Radiation Intensity vs Angle', fontsize=12)
ax2.legend()
ax2.grid(True, alpha=0.3)
ax2.set_xlim([0, 90])

plt.tight_layout()
plt.savefig(f'{output_dir}/radiation_intensity.png', dpi=150, bbox_inches='tight')
plt.close()
print("辐射强度方向分布图已保存")

3.3 实际物体辐射特性

"""
实际物体与黑体辐射的对比
"""

# 定义几种实际材料的发射率（随波长变化）
def emissivity_metal(lambda_um):
    """金属的发射率模型（简化）"""
    return 0.1 + 0.05 * np.sin(lambda_um / 2)

def emissivity_ceramic(lambda_um):
    """陶瓷的发射率模型（简化）"""
    return 0.8 + 0.1 * np.cos(lambda_um / 5)

def emissivity_coating(lambda_um):
    """选择性涂层的发射率"""
    # 太阳波段高吸收，红外波段低发射
    if lambda_um < 2.5:
        return 0.9
    else:
        return 0.1

T = 800  # K

# 计算实际物体的辐射光谱
lambda_range = np.linspace(0.5, 20, 500)

E_blackbody = [planck_law(lam, T) for lam in lambda_range]
E_metal = [planck_law(lam, T) * emissivity_metal(lam) for lam in lambda_range]
E_ceramic = [planck_law(lam, T) * emissivity_ceramic(lam) for lam in lambda_range]
E_coating = [planck_law(lam, T) * emissivity_coating(lam) for lam in lambda_range]

fig, axes = plt.subplots(2, 2, figsize=(14, 12))

# 1. 不同材料的辐射光谱对比
axes[0, 0].plot(lambda_range, E_blackbody, 'k-', linewidth=2, label='Blackbody')
axes[0, 0].plot(lambda_range, E_metal, 'b-', linewidth=2, label='Metal')
axes[0, 0].plot(lambda_range, E_ceramic, 'r-', linewidth=2, label='Ceramic')
axes[0, 0].plot(lambda_range, E_coating, 'g-', linewidth=2, label='Selective Coating')
axes[0, 0].set_xlabel('Wavelength (μm)', fontsize=11)
axes[0, 0].set_ylabel('Spectral Emissive Power (W/m²·μm)', fontsize=11)
axes[0, 0].set_title(f'Radiation Spectrum at T = {T} K', fontsize=12)
axes[0, 0].legend()
axes[0, 0].grid(True, alpha=0.3)
axes[0, 0].set_xlim([0, 15])

# 2. 发射率随波长变化
eps_metal = [emissivity_metal(lam) for lam in lambda_range]
eps_ceramic = [emissivity_ceramic(lam) for lam in lambda_range]
eps_coating = [emissivity_coating(lam) for lam in lambda_range]

axes[0, 1].plot(lambda_range, eps_metal, 'b-', linewidth=2, label='Metal')
axes[0, 1].plot(lambda_range, eps_ceramic, 'r-', linewidth=2, label='Ceramic')
axes[0, 1].plot(lambda_range, eps_coating, 'g-', linewidth=2, label='Selective Coating')
axes[0, 1].axhline(y=1.0, color='k', linestyle='--', label='Blackbody')
axes[0, 1].set_xlabel('Wavelength (μm)', fontsize=11)
axes[0, 1].set_ylabel('Emissivity', fontsize=11)
axes[0, 1].set_title('Spectral Emissivity', fontsize=12)
axes[0, 1].legend()
axes[0, 1].grid(True, alpha=0.3)
axes[0, 1].set_ylim([0, 1.2])

# 3. 温度对总发射率的影响
T_range = np.linspace(300, 1500, 100)

# 假设发射率随温度变化
emissivity_avg_metal = 0.15 + 0.0001 * (T_range - 300)
emissivity_avg_ceramic = 0.85 - 0.0001 * (T_range - 300)

axes[1, 0].plot(T_range, emissivity_avg_metal, 'b-', linewidth=2, label='Metal')
axes[1, 0].plot(T_range, emissivity_avg_ceramic, 'r-', linewidth=2, label='Ceramic')
axes[1, 0].axhline(y=1.0, color='k', linestyle='--', label='Blackbody')
axes[1, 0].set_xlabel('Temperature (K)', fontsize=11)
axes[1, 0].set_ylabel('Total Emissivity', fontsize=11)
axes[1, 0].set_title('Temperature Dependence of Emissivity', fontsize=12)
axes[1, 0].legend()
axes[1, 0].grid(True, alpha=0.3)
axes[1, 0].set_ylim([0, 1.2])

# 4. 辐射能量对比（饼图）
E_bb_total = stefan_boltzmann(T)
E_metal_total = E_bb_total * np.mean(emissivity_avg_metal)
E_ceramic_total = E_bb_total * np.mean(emissivity_avg_ceramic)

labels = ['Blackbody', 'Metal', 'Ceramic']
sizes = [E_bb_total, E_metal_total, E_ceramic_total]
colors_pie = ['gold', 'lightblue', 'lightcoral']

axes[1, 1].pie(sizes, labels=labels, colors=colors_pie, autopct='%1.1f%%', startangle=90)
axes[1, 1].set_title(f'Radiation Energy Comparison at T = {T} K', fontsize=12)

plt.tight_layout()
plt.savefig(f'{output_dir}/real_surface_radiation.png', dpi=150, bbox_inches='tight')
plt.close()
print("实际物体辐射特性图已保存")

3.4 辐射换热基础计算

"""
简单几何形状间的辐射换热计算
"""

# 两个平行平板间的辐射换热
def parallel_plates_radiation(T1, T2, eps1, eps2, A=1.0):
    """
    两个无限大平行平板间的辐射换热
    T1, T2: 两板温度 (K)
    eps1, eps2: 两板发射率
    A: 面积 (m²)
    返回: 换热量 (W)
    """
    sigma = 5.67e-8
    
    # 等效辐射热阻
    R_total = (1/eps1 + 1/eps2 - 1) / (A * sigma)
    
    # 换热量
    Q = (T1**4 - T2**4) / R_total
    
    return Q

# 两个同心圆球间的辐射换热
def concentric_spheres_radiation(T1, T2, eps1, eps2, r1, r2):
    """
    两个同心圆球间的辐射换热
    r1: 内球半径 (m)
    r2: 外球半径 (m)
    返回: 换热量 (W)
    """
    sigma = 5.67e-8
    
    A1 = 4 * np.pi * r1**2
    A2 = 4 * np.pi * r2**2
    
    # 等效辐射热阻
    R1 = (1 - eps1) / (eps1 * A1)
    R12 = 1 / A1  # 空间热阻
    R2 = (1 - eps2) / (eps2 * A2)
    
    R_total = R1 + R12 + R2
    
    Q = sigma * (T1**4 - T2**4) / R_total
    
    return Q

# 计算示例
T_hot = 800  # K
T_cold = 300  # K

# 不同发射率组合下的换热量
eps_range = np.linspace(0.1, 1.0, 50)

Q_parallel_same = [parallel_plates_radiation(T_hot, T_cold, eps, eps) for eps in eps_range]
Q_parallel_diff = [parallel_plates_radiation(T_hot, T_cold, 0.9, eps) for eps in eps_range]

fig, axes = plt.subplots(1, 2, figsize=(14, 6))

# 1. 平行平板换热量随发射率变化
axes[0].plot(eps_range, np.array(Q_parallel_same)/1000, 'b-', linewidth=2, label='ε₁ = ε₂')
axes[0].plot(eps_range, np.array(Q_parallel_diff)/1000, 'r--', linewidth=2, label='ε₁ = 0.9, ε₂ = variable')
axes[0].set_xlabel('Emissivity', fontsize=11)
axes[0].set_ylabel('Heat Transfer Rate (kW)', fontsize=11)
axes[0].set_title(f'Radiation Heat Transfer\n(T₁={T_hot}K, T₂={T_cold}K)', fontsize=12)
axes[0].legend()
axes[0].grid(True, alpha=0.3)

# 2. 温度差对换热量的影响
T_hot_range = np.linspace(400, 1000, 100)
Q_temp_diff = [parallel_plates_radiation(T, T_cold, 0.8, 0.8) for T in T_hot_range]

axes[1].plot(T_hot_range - T_cold, np.array(Q_temp_diff)/1000, 'b-', linewidth=2)
axes[1].set_xlabel('Temperature Difference (K)', fontsize=11)
axes[1].set_ylabel('Heat Transfer Rate (kW)', fontsize=11)
axes[1].set_title('Heat Transfer vs Temperature Difference\n(ε₁ = ε₂ = 0.8)', fontsize=12)
axes[1].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig(f'{output_dir}/radiation_heat_transfer.png', dpi=150, bbox_inches='tight')
plt.close()
print("辐射换热计算图已保存")

# 打印一些计算结果
print("\n辐射换热计算结果:")
print(f"黑体辐射力 (800K): {stefan_boltzmann(800):.2f} W/m²")
print(f"黑体辐射力 (300K): {stefan_boltzmann(300):.2f} W/m²")
print(f"平行平板换热量 (ε=0.8): {parallel_plates_radiation(800, 300, 0.8, 0.8):.2f} W")
print(f"峰值波长 (800K): {wien_law(800):.2f} μm")
print(f"峰值波长 (300K): {wien_law(300):.2f} μm")

print("\n" + "="*60)
print("仿真5：辐射传热系数计算")
print("="*60)

# 计算不同几何结构的辐射传热系数

# 1. 无限大平行平板
T1, T2 = 800, 300  # K
epsilon1, epsilon2 = 0.8, 0.8

# 辐射传热系数
h_rad_parallel = 4 * 5.67e-8 * (T1**3 + T1**2*T2 + T1*T2**2 + T2**3) / (1/epsilon1 + 1/epsilon2 - 1)

# 2. 同心圆柱
r1, r2 = 0.1, 0.15  # m
h_rad_cylinder = 4 * 5.67e-8 * (T1**3 + T1**2*T2 + T1*T2**2 + T2**3) / \
                 (1/epsilon1 + r1/r2*(1/epsilon2 - 1))

# 3. 同心球
h_rad_sphere = 4 * 5.67e-8 * (T1**3 + T1**2*T2 + T1*T2**2 + T2**3) / \
               (1/epsilon1 + (r1/r2)**2*(1/epsilon2 - 1))

print(f"辐射传热系数对比 (T1={T1}K, T2={T2}K):")
print(f"  平行平板: {h_rad_parallel:.2f} W/(m²·K)")
print(f"  同心圆柱: {h_rad_cylinder:.2f} W/(m²·K)")
print(f"  同心球: {h_rad_sphere:.2f} W/(m²·K)")

# 温度对辐射传热系数的影响
T1_range = np.linspace(300, 1000, 50)
h_rad_range = 4 * 5.67e-8 * (T1_range**3 + T1_range**2*T2 + T1_range*T2**2 + T2**3) / \
              (1/epsilon1 + 1/epsilon2 - 1)

# 可视化
fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# 辐射传热系数vs温度
ax1 = axes[0]
ax1.plot(T1_range, h_rad_range, 'b-', linewidth=2)
ax1.fill_between(T1_range, h_rad_range, alpha=0.3)
ax1.set_xlabel('Surface Temperature T₁ (K)', fontsize=11)
ax1.set_ylabel('Radiation Heat Transfer Coefficient (W/(m²·K))', fontsize=11)
ax1.set_title('Radiation Heat Transfer Coefficient vs Temperature', fontsize=12, fontweight='bold')
ax1.grid(True, alpha=0.3)

# 几何结构对比
ax2 = axes[1]
geometries = ['Parallel\nPlates', 'Concentric\nCylinders', 'Concentric\nSpheres']
h_values = [h_rad_parallel, h_rad_cylinder, h_rad_sphere]
colors = ['steelblue', 'forestgreen', 'orangered']

bars = ax2.bar(geometries, h_values, color=colors, alpha=0.8)
ax2.set_ylabel('Heat Transfer Coefficient (W/(m²·K))', fontsize=11)
ax2.set_title('Comparison of Radiation Heat Transfer Coefficients', fontsize=12, fontweight='bold')
ax2.grid(True, alpha=0.3, axis='y')

for bar, val in zip(bars, h_values):
    ax2.text(bar.get_x() + bar.get_width()/2., bar.get_height() + 0.5,
             f'{val:.1f}', ha='center', va='bottom', fontsize=11, fontweight='bold')

plt.tight_layout()
plt.savefig(f'{output_dir}/radiation_htc.png', dpi=150, bbox_inches='tight')
plt.close()

print("图5：辐射传热系数已保存")

print("\n" + "="*60)
print("仿真6：辐射屏蔽效应")
print("="*60)

# 辐射屏蔽层效果分析

# 无屏蔽
Q_no_shield = 5.67e-8 * (T1**4 - T2**4) / (1/epsilon1 + 1/epsilon2 - 1)

# 有n层屏蔽
n_shields = np.arange(0, 6)
Q_shielded = []

for n in n_shields:
    if n == 0:
        Q_shielded.append(Q_no_shield)
    else:
        # 假设所有表面发射率相同
        Q_shielded.append(Q_no_shield / (n + 1))

# 屏蔽效率
shield_efficiency = [(1 - Q/Q_no_shield)*100 for Q in Q_shielded]

print(f"无屏蔽换热量: {Q_no_shield:.2f} W/m²")
for i, (n, Q, eff) in enumerate(zip(n_shields[1:], Q_shielded[1:], shield_efficiency[1:])):
    print(f"  {n}层屏蔽: {Q:.2f} W/m², 效率: {eff:.1f}%")

# 可视化
fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# 换热量vs屏蔽层数
ax1 = axes[0]
ax1.plot(n_shields, Q_shielded, 'bo-', linewidth=2, markersize=8)
ax1.fill_between(n_shields, Q_shielded, alpha=0.3)
ax1.set_xlabel('Number of Radiation Shields', fontsize=11)
ax1.set_ylabel('Heat Transfer Rate (W/m²)', fontsize=11)
ax1.set_title('Effect of Radiation Shields', fontsize=12, fontweight='bold')
ax1.grid(True, alpha=0.3)

# 屏蔽效率
ax2 = axes[1]
ax2.bar(n_shields[1:], shield_efficiency[1:], color='forestgreen', alpha=0.8)
ax2.set_xlabel('Number of Radiation Shields', fontsize=11)
ax2.set_ylabel('Shielding Efficiency (%)', fontsize=11)
ax2.set_title('Radiation Shielding Efficiency', fontsize=12, fontweight='bold')
ax2.grid(True, alpha=0.3, axis='y')

for i, eff in enumerate(shield_efficiency[1:]):
    ax2.text(i, eff + 1, f'{eff:.1f}%', ha='center', va='bottom', fontsize=10, fontweight='bold')

plt.tight_layout()
plt.savefig(f'{output_dir}/radiation_shield.png', dpi=150, bbox_inches='tight')
plt.close()

print("图6：辐射屏蔽效应已保存")
print("\n所有仿真完成！")

---

4. 结果分析与讨论

4.1 黑体辐射光谱特性

从仿真结果可以看出：

1. 温度对光谱分布的影响：

   • 温度升高，辐射峰值向短波方向移动（维恩位移定律）
   • 300K时峰值波长约9.7 μm（远红外）
   • 5800K（太阳表面）时峰值波长约0.5 μm（可见光）

2. 总辐射力与温度的关系：

   • 遵循T⁴关系，温度升高导致辐射力急剧增加
   • 800K表面的辐射力是300K表面的约50倍

3. 可见光能量比例：

   • 低温物体（<1000K）主要发射红外辐射
   • 太阳（5800K）约40%的能量在可见光范围

4.2 实际物体与黑体的差异

实际物体的辐射特性与黑体存在显著差异：

1. 发射率小于1：

   • 金属发射率通常较低（0.1-0.3）
   • 非金属发射率较高（0.8-0.95）
   • 发射率随波长和温度变化

2. 选择性表面：

   • 太阳能选择性吸收涂层在可见光区高吸收，红外区低发射
   • 辐射制冷材料在红外大气窗口高发射

4.3 辐射换热的特点

1. 非线性特性：

   • 换热量与温度的四次方差成正比
   • 高温时辐射换热占主导地位

2. 表面特性敏感：

   • 发射率对换热量影响显著
   • 表面粗糙度、氧化层等都会改变辐射特性

---

5. 工程应用

5.1 航天器热控

航天器在太空面临极端热环境：

• 太阳直接辐射：约1367 W/m²
• 地球反照：约300 W/m²
• 地球红外辐射：约200 W/m²
• 太空背景：接近0 K

通过选择合适的表面涂层（发射率和吸收率），可以控制航天器温度。

5.2 太阳能利用

1. 太阳能集热器：

   • 选择性吸收涂层：太阳吸收率>0.95，红外发射率<0.1
   • 提高集热效率

2. 光伏电池：

   • 理解光谱响应特性
   • 优化带隙设计

5.3 高温工业炉

• 辐射是主要传热方式
• 炉壁发射率影响加热效率
• 需要精确计算辐射换热

5.4 建筑节能

• 窗户的辐射特性影响能耗
• Low-E玻璃控制红外辐射
• 外墙材料选择考虑辐射特性

---

6. 总结

本教程系统介绍了热辐射的基本定律，主要内容包括：

1. 理论基础：介绍了热辐射的本质、电磁波谱特性、黑体辐射模型等基本概念

2. 基本定律：详细阐述了普朗克定律、维恩位移定律、斯特藩-玻尔兹曼定律和兰贝特定律

3. 数值仿真：通过Python实现了黑体辐射光谱计算、方向辐射特性分析和辐射换热计算

4. 工程应用：讨论了热辐射理论在航天、太阳能、工业炉和建筑等领域的应用

热辐射是传热学的重要组成部分，理解其基本定律对于解决实际工程问题具有重要意义。通过本教程的学习，读者应能够掌握热辐射的基本计算方法，并能够应用于实际问题的分析。

---