简介:这份资源面向从事红外仿真、物理模拟与科学计算的工程师及科研人员,提供基于MATLAB实现普朗克黑体辐射公式的完整代码集合,帮助解决不同温度与波长组合下辐射出射度的数值计算与可视化问题。压缩包内共5个m文件,整体约3KB,涵盖普朗克公式主计算、辐射度求解、曲线拟合及映射函数等模块,分别承担核心公式运算、结果拟合与数据映射等用途,结构紧凑便于直接调用与二次修改。目前已有4854人学习下载,说明其在相关领域具有较高的参考价值。读者可借助这些脚本快速搭建温度与波长参数扫描流程,得到辐射强度矩阵并绘制三维分布图,同时通过曲线拟合与映射函数模块理解黑体辐射随温度变化的规律,为红外仿真建模、光谱分析及教学演示提供可复用的计算基础。
1. 普朗克公式到底在算什么:从“紫外灾难”说起
如果你把一块铁放进火里烧,它会先变暗红,再变橙黄,最后亮到发白。这个“颜色随温度变化”的现象,背后藏着一个让物理学家头疼了几十年的问题:一个物体在不同温度下,到底会辐射出多少能量、分布在哪些波长上?
普朗克公式就是回答这个问题的数学表达式。它描述的是黑体辐射的光谱辐射度——通俗地说,就是“在某个温度下,某个波长的光,每单位面积、单位波长、单位立体角里携带多少能量”。单位通常是 W/(m²·sr·m) 或 W/(m²·sr·μm)。
在普朗克之前,经典物理给出的公式在短波长区域会算出“能量趋于无穷大”的荒谬结果,被戏称为“紫外灾难”。普朗克在 1900 年提出,能量不是连续的,而是一份一份的,每份能量等于 hν(h 是普朗克常数,ν 是频率)。这个“能量量子化”的假设,不仅修好了公式,还顺手开启了量子力学的大门。
用 MATLAB 来计算普朗克公式,本质上就是把这个物理规律翻译成代码,让你能快速画出不同温度下的辐射曲线,或者反推某个波长对应的辐射度。搞懂它,你就能理解为什么红外测温枪能测体温、为什么太阳看起来是黄的、为什么白炽灯发光效率那么低。
2. 为什么普朗克公式值得用代码算:手算的坑与代码的爽
2.1 手算普朗克公式的三个“翻车点”
普朗克公式长这样:
B(λ, T) = (2hc² / λ⁵) · 1 / (exp(hc / (λkT)) - 1)
其中 h 是普朗克常数,c 是光速,k 是玻尔兹曼常数,λ 是波长,T 是温度。
如果你拿计算器手算,第一个坑是数量级跨度太大。λ 在可见光范围大约是 5×10⁻⁷ 米,λ⁵ 就是 3×10⁻³¹ 量级;而 exp(hc/(λkT)) 里的指数部分,在常温下可能大到几十甚至上百。手算时,一个指数溢出就能让你得到“Inf”或者“0”,完全看不出物理意义。
第二个坑是单位不统一。物理常数用的是米和开尔文,但工程上经常用微米和摄氏度。如果你把 λ=0.5(微米)直接代进公式,而 h、c、k 用的是国际单位制,结果会差 10²⁰ 倍以上。这种错误在血泪经验里太常见了——算出来的曲线形状对,但数值完全没法看。
第三个坑是需要批量计算。你想看 300K、500K、1000K、3000K 四条曲线放在一起对比,手算就得重复几十次,每次还要处理指数和幂次。一旦某个温度下波长范围没选对,曲线就画不全。
2.2 MATLAB 把这三个坑一次填平
MATLAB 的优势在于向量化运算和内置常数。你可以把波长定义成一个数组,比如从 0.1 微米到 20 微米取 1000 个点,然后一行代码算出所有波长对应的辐射度,不需要写循环。
% 定义物理常数(国际单位制) h = 6.62607015e-34; % 普朗克常数,单位 J·s c = 2.99792458e8; % 光速,单位 m/s k = 1.380649e-23; % 玻尔兹曼常数,单位 J/K % 定义波长范围:0.1 到 20 微米,转成米 lambda = linspace(0.1e-6, 20e-6, 1000); % 定义温度 T = 300; % 开尔文 % 普朗克公式(向量化计算) B = (2*h*c^2) ./ (lambda.^5) ./ (exp(h*c ./ (lambda*k*T)) - 1);这段代码里,./和.^是逐元素运算,exp会自动对数组每个元素求指数。算完之后,B就是一个 1000 个元素的数组,直接plot(lambda, B)就能出图。
参数说明:lambda的单位必须是米,因为 h、c、k 都是国际单位制。如果你习惯用微米,记得在代入公式前乘以1e-6。温度T必须是开尔文,摄氏度要加 273.15。
2.3 从“能算”到“算得明白”:几个必须知道的物理量
用代码算只是第一步,看懂结果才是目的。普朗克公式里藏着几个关键推论:
维恩位移定律:辐射度峰值对应的波长 λ_max 满足 λ_max · T ≈ 2898 μm·K。也就是说,温度越高,峰值波长越短。太阳表面约 5800K,峰值在 0.5 微米左右,正好是绿光——所以太阳看起来是白的,但大气散射让天空变蓝。
斯特藩-玻尔兹曼定律:对所有波长积分,总辐射度正比于 T⁴。温度翻一倍,总辐射能量变成 16 倍。这就是为什么高温物体散热极快。
瑞利-金斯极限:当波长很长时,普朗克公式退化成经典公式,和实验吻合。短波长时,量子效应才显现出来。
在 MATLAB 里,你可以用trapz(lambda, B)做数值积分验证 T⁴ 规律,也可以用[~, idx] = max(B)找到峰值位置,再反推 λ_max · T 是否接近 2898。
3. 普朗克公式和相近概念的区别:别把辐射度、辐照度、亮度混成一锅
3.1 辐射度、辐照度、辐射强度:三个“度”到底差在哪
初学热辐射时,最容易搞混的就是这几个词。它们都描述“能量”,但限定条件不同。
辐射度(Radiance):普朗克公式算的就是这个。它描述的是“从单位面积、向单位立体角、在单位波长间隔内发出的能量”。关键词是“单位立体角”和“单位面积”。单位是 W/(m²·sr·m)。
辐照度(Irradiance):描述的是“到达单位面积上的总能量”,不区分方向。单位是 W/m²。如果你把普朗克公式对所有立体角积分,再对所有波长积分,得到的就是黑体表面的总辐照度(也就是斯特藩-玻尔兹曼定律的结果)。
辐射强度(Radiant Intensity):描述的是“点光源向单位立体角发出的能量”,不除以面积。单位是 W/sr。
打个比方:辐射度像“某个具体方向的探照灯亮度”,辐照度像“整个探照灯照到墙上的总光量”,辐射强度像“探照灯本身有多亮”。在 MATLAB 里,如果你只算普朗克公式,得到的是辐射度;要得到辐照度,需要乘以 π(对半球积分)。
3.2 常见误区:普朗克公式不是“黑体总辐射公式”
很多人第一次接触普朗克公式,会以为它算的是“黑体辐射的总能量”。其实不是。它算的是光谱分布——也就是“每个波长贡献多少”。总能量需要积分。
另一个误区是把频率形式和波长形式混用。普朗克公式有两种写法:按频率 ν 的 B_ν 和按波长 λ 的 B_λ。两者之间差一个 λ²/c 的因子,因为 dν = -c/λ² dλ。如果你在 MATLAB 里用频率算了一遍,又用波长算了一遍,发现数值对不上,大概率是忘了这个转换。
还有一个坑是把“黑体”当成“黑色物体”。黑体是指“对所有波长都完全吸收、完全发射”的理想物体,不是颜色黑。太阳、白炽灯灯丝、人体皮肤,在特定波段都可以近似成黑体。
3.3 用 MATLAB 验证区别:一个小实验
你可以写一段代码,分别算辐射度和辐照度,看看差多少。
% 算辐射度(普朗克公式) B = (2*h*c^2) ./ (lambda.^5) ./ (exp(h*c ./ (lambda*k*T)) - 1); % 算辐照度:对半球立体角积分,乘以 pi E = pi * B; % 算总辐照度:对波长积分 E_total = trapz(lambda, E); % 验证斯特藩-玻尔兹曼定律:sigma * T^4 sigma = 5.670374419e-8; E_theory = sigma * T^4; % 输出对比 fprintf('数值积分总辐照度: %.4e W/m²\n', E_total); fprintf('理论值 sigma*T^4: %.4e W/m²\n', E_theory);运行后你会发现两者非常接近(误差来自积分范围和离散化)。这个实验能帮你把“辐射度”和“辐照度”的区别刻进脑子里。
4. 普朗克公式在 MATLAB 里的典型用法:从画曲线到反推温度
4.1 画不同温度下的辐射曲线:一眼看出维恩位移
最直观的用法就是画图。把 300K、1000K、3000K、5800K 四条曲线画在一起,你能立刻看到峰值往短波方向移动。
T_list = [300, 1000, 3000, 5800]; lambda = linspace(0.1e-6, 20e-6, 2000); figure; hold on; for T = T_list B = (2*h*c^2) ./ (lambda.^5) ./ (exp(h*c ./ (lambda*k*T)) - 1); plot(lambda*1e6, B, 'LineWidth', 1.5); end xlabel('波长 (\mum)'); ylabel('光谱辐射度 W/(m^2·sr·m)'); legend('300K', '1000K', '3000K', '5800K'); title('不同温度下的普朗克辐射曲线'); grid on;注意lambda*1e6把米转成微米,方便看横轴。hold on让多条曲线叠加。legend里的温度标签要手动对应。
这张图能直接回答“为什么高温物体会发白光”——5800K 的峰值在可见光中间,而 300K 的峰值在 10 微米左右,属于红外,人眼看不见。
4.2 反推温度:从辐射度算温度
有时候你已知某个波长下的辐射度,想反推温度。普朗克公式对 T 没有解析解,但可以用数值方法。
% 已知:波长 0.65 微米,辐射度 1e12 W/(m^2·sr·m) lambda_target = 0.65e-6; B_target = 1e12; % 定义误差函数 err = @(T) abs((2*h*c^2 / lambda_target^5) / (exp(h*c/(lambda_target*k*T)) - 1) - B_target); % 用 fminbnd 找温度 T_est = fminbnd(err, 100, 10000); fprintf('反推温度: %.2f K\n', T_est);fminbnd在 100 到 10000K 之间找误差最小的温度。这个方法在红外测温里有实际对应——测温枪就是测某个波段的辐射度,再反推温度。
4.3 参数怎么改:波长范围、温度、精度
波长范围:如果你只关心可见光,用linspace(0.38e-6, 0.78e-6, 500)。如果关心红外,用linspace(1e-6, 20e-6, 1000)。范围太宽会导致短波区域点数太少,曲线不平滑。
温度:常温物体用 300K,火焰用 1000-1500K,白炽灯丝用 2500-3000K,太阳用 5800K。温度越高,峰值越靠短波,波长上限可以适当减小。
精度:linspace的点数决定曲线平滑度。1000 个点通常够用,如果要做积分,建议 5000 个点以上。trapz的精度和点数直接相关。
失败时看什么:如果曲线出现NaN,检查exp的参数是否溢出。当h*c/(lambda*k*T)大于 700 左右时,exp会返回Inf,导致1/Inf = 0。这不是错误,是物理上短波辐射度趋近于零。如果出现Inf,检查lambda是否太小或T是否太高。
5. 什么时候该用普朗克公式,什么时候别用
5.1 适用场景:黑体或近似黑体的辐射计算
普朗克公式的前提是黑体——完全吸收、完全发射的理想物体。实际物体用发射率 ε修正:实际辐射度 = ε · B(λ, T)。ε 在 0 到 1 之间,黑体 ε=1。
适合用普朗克公式的场景:
- 白炽灯、卤素灯的灯丝辐射(ε 约 0.3-0.5)
- 太阳辐射(近似 5800K 黑体)
- 红外测温、热成像(人体 ε 约 0.98,接近黑体)
- 高温炉、燃烧产物的辐射分析
- 天文观测中的恒星光谱
不适合的场景:
- 荧光灯、LED——它们不是热辐射,是电致发光,光谱由能带结构决定,普朗克公式不适用
- 激光——单色、相干,完全不同的机制
- 大气辐射——需要考虑气体吸收谱线,不是连续黑体谱
5.2 一个具体例子:用 MATLAB 估算白炽灯效率
白炽灯灯丝温度约 2800K,只有约 5% 的能量落在可见光波段(0.38-0.78 微米),其余都是红外。你可以用 MATLAB 算一下。
T = 2800; lambda = linspace(0.1e-6, 20e-6, 5000); B = (2*h*c^2) ./ (lambda.^5) ./ (exp(h*c ./ (lambda*k*T)) - 1); % 总辐射度(对波长积分) B_total = trapz(lambda, B); % 可见光波段辐射度 visible = lambda >= 0.38e-6 & lambda <= 0.78e-6; B_visible = trapz(lambda(visible), B(visible)); % 效率 eff = B_visible / B_total * 100; fprintf('可见光效率: %.2f%%\n', eff);运行后你会得到 5% 左右。这个数字解释了为什么白炽灯被淘汰——大部分电能变成了热,而不是光。
5.3 我踩过的坑:单位、积分范围、发射率
刚用 MATLAB 算普朗克公式时,我犯过三个错。第一次是波长用了微米但没转米,结果辐射度大了 10²⁰ 倍,曲线形状对但数值完全离谱。第二次是积分范围只到 5 微米,漏掉了长波红外,总能量算少了 30%。第三次是忘了乘发射率,把灯丝当成理想黑体,效率算出来偏高。
从那以后,我每次用普朗克公式都强制走一遍检查:波长单位是不是米?温度是不是开尔文?积分范围覆盖了峰值和长波尾巴吗?发射率要不要乘?这四个问题过一遍,基本不会翻车。
希望帮到你。
本文还有配套的精品资源,点击获取