The Stefan–Boltzmann law: the power a black body radiates¶
The power a black body radiates from each square metre of its surface into the half-space above it, the Stefan-Boltzmann law, j = sigma T^4; sigma is lusine.constants', 2 pi^5 k^4 / (15 h^3 c^2), exact.
It is pi times Planck's radiance (lusine.radiation.planck) integrated over every wavelength. A grey body's is its emissivity times this, the caller's to apply. Source: Stefan (1879), Boltzmann (1884); Rybicki & Lightman, Radiative Processes in Astrophysics, 1.4 and 1.5 (eq. 1.43). Checked against a reference in numpy, which has its own tests.
- Source
- Stefan (1879), Boltzmann (1884); Rybicki & Lightman, Radiative Processes in Astrophysics, 1.4 and 1.5 (eq. 1.43).
- Proof
- reference, invariant.
radiated_power.trait
namespace lusine.radiation;
float<W/m2> radiated_power(float<K> temperature)
{
let t2 = temperature * temperature;
return lusine.constants.sigma * t2 * t2;
}
The power a black body radiates from each square metre of its surface, j = σT⁴, with σ read from lusine.constants (2π⁵k⁴/(15h³c²), exact). Reading the constant makes it a binding's kernel, not a standalone function: it runs in a program, as this page runs it.
import os, sys
root = os.getcwd()
while not os.path.isfile(os.path.join(root, 'docsite.json')):
root = os.path.dirname(root)
sys.path.insert(0, os.path.join(root, 'tools'))
import numpy as np
import matplotlib.pyplot as plt
from lusine import trait_run
from lusine.library_sources import parse_rows
EPS = float(np.finfo(np.float32).eps) # Trait's float is binary32
import re
def constant(name):
"""A constant of lusine.constants as the library prints it."""
return float(re.search(r' ' + name + r' = ([-0-9.e]+) ', trait_run.source('constants.trait')).group(1))
The kernel¶
The law as the library states it.
print(trait_run.source('radiation/radiated_power.trait'))
// lusine.radiation.radiated_power: the power a black body radiates from each square metre of its surface into the
// half-space above it, the Stefan-Boltzmann law, j = sigma T^4; sigma is lusine.constants', 2 pi^5 k^4 / (15 h^3 c^2),
// exact. It is pi times Planck's radiance (lusine.radiation.planck) integrated over every wavelength. A grey body's is
// its emissivity times this, the caller's to apply.
// Source: Stefan (1879), Boltzmann (1884); Rybicki & Lightman, Radiative Processes in Astrophysics, 1.4 and 1.5
// (eq. 1.43).
// Checked against a reference in numpy, which has its own tests.
namespace lusine.radiation;
float<W/m2> radiated_power(float<K> temperature)
{
let t2 = temperature * temperature;
return lusine.constants.sigma * t2 * t2;
}
The program¶
Each section runs the kernel over a lattice of items, one a case.
print(trait_run.law_program('radiation/radiated_power.trait', 'radiated_power'))
namespace probe;
lusine.version version = { { number = 1 } };
struct item {
host float<K> temperature;
float<W/m2> result_0;
};
lattice<item, 1> items(spacing = 1 m, boundary = wall);
void main() { items.result_0 = lusine.radiation.radiated_power(items.temperature); }
Reference: a room and the Sun¶
A black body at 300 K radiates 459.3 W/m². The Sun's nominal effective temperature, 5772 K, over its nominal radius, 6.957 × 10⁸ m, radiates its nominal luminosity, 3.828 × 10²⁶ W (IAU 2015 Resolution B3, which defines the three to agree), within 0.1%.
j = trait_run.law('radiation/radiated_power.trait', 'radiated_power', temperature=[300.0, 5772.0])
lum = float(j[1]) * 4 * np.pi * 6.957e8 ** 2
assert abs(j[0] - 459.3) < 0.1 and abs(lum / 3.828e26 - 1) < 1e-3, (j, lum)
print(f'300 K: {j[0]:.2f} W/m2; the Sun: {j[1]:.4e} W/m2, {lum:.4e} W')
300 K: 459.30 W/m2; the Sun: 6.2939e+07 W/m2, 3.8280e+26 W
Invariant: π times Planck's radiance over every wavelength, and T⁴¶
The emitted lusine.radiation.planck at 1500 K, integrated over ln λ from 50 nm to 1 cm, times π, is this kernel's power within 10⁻⁴; twice the temperature radiates sixteen times the power, bit for bit (σ(2T)²(2T)² scales by powers of two without rounding), and j/(σT⁴) is 1 to binary32's rounding with lusine.constants.sigma as printed.
lam = np.geomspace(50e-9, 1e-2, 4001)
b = trait_run.law('radiation/planck.trait', 'planck', wavelength=lam, temperature=1500.0)
y = b * lam
total = np.pi * float(np.sum(0.5 * (y[1:] + y[:-1]) * np.diff(np.log(lam))))
j = float(trait_run.law('radiation/radiated_power.trait', 'radiated_power', temperature=1500.0)[0])
assert abs(total / j - 1) < 1e-4, total / j - 1
ts = np.float32(np.geomspace(1.0, 1e6, 32)).astype(float)
j1 = trait_run.law('radiation/radiated_power.trait', 'radiated_power', temperature=ts)
j2 = trait_run.law('radiation/radiated_power.trait', 'radiated_power', temperature=2 * ts)
assert np.array_equal(j2, 16 * j1)
assert np.all(np.abs(j1 / (constant('sigma') * ts ** 4) - 1) <= 4 * EPS)
print(f'pi * integral / j - 1 = {total / j - 1:+.1e}')
fig, ax = plt.subplots(figsize=(7, 3.2))
ax.loglog(ts, j1)
ax.set_xlabel('T (K)')
ax.set_ylabel('j (W/m2)')
plt.show()
pi * integral / j - 1 = +2.8e-07
