黒体放射の放射モードのエネルギーを量子化することでプランクの公式を導出し、全周波数にわたる積分にリーマンゼータ関数の値$\zeta(4)$が現れ、シュテファン=ボルツマンの法則が導かれる数理構造を整理します。
19世紀末の物理学において、熱平衡にある空洞放射(黒体放射)のスペクトル分布を古典統計力学と電磁気学から導く試みは、高周波領域でエネルギー密度が無限大に発散する「紫外破綻」という重大な困難に直面しました。プランクによるエネルギー量子化の仮説はこの困難を克服し、量子論の幕開けを告げる契機となりました。同時に、このプランクの放射式から全放射エネルギーを求める積分計算は、純粋数学の対象であるリーマンゼータ関数と深く結びついています。
本記事では、低周波のレイリー=ジーンズの公式、高周波のウィーンの公式、および全周波数をカバーするプランクの公式の3者が、「モード密度×平均エネルギー」という共通構造を持つことを確認します。次に、連続的なボルツマン分布による古典平均エネルギーの積分計算と、離散的なエネルギー準位に基づく級数計算を対比し、エネルギー量子化がどのようにプランクの公式を導くかを具体的に示します。さらに、プランクの公式を全周波数にわたって積分することでリーマンゼータ関数$\zeta(4)$が現れ、シュテファン=ボルツマンの法則が導出される過程を計算します。最後に、パーセヴァルの等式を用いて$\zeta(4)=\pi^4/90$を初等的に導出します。
微積分、広義積分と部分積分、およびフーリエ級数展開の基礎を前提とします。本記事では空洞放射モードの幾何学的計数とエネルギー量子化に伴う代数・積分計算に焦点を当て、熱力学ポテンシャルや量子電磁力学の微視的詳細には立ち入りません。
黒体放射のスペクトルを表す式として、ここでは3つの公式を比べます。現在の立場から見ると、低周波数側はレイリー=ジーンズの公式、高周波数側はウィーンの公式が表します。
$$ u(\nu, T) = \frac{8\pi\nu^2}{c^3} kT $$
$$ u(\nu, T) = \frac{8\pi h\nu^3}{c^3} \frac{1}{e^{h\nu/kT}} $$
ウィーンの公式は、歴史的には経験的な定数を用いて表されました。ここではプランクの公式と比べるため、その定数を$h$と$k$で表した形を用います。
全周波数を1つの式で表すのがプランクの公式です。
$$ u(\nu, T) = \frac{8\pi h\nu^3}{c^3} \frac{1}{e^{h\nu/kT} - 1} $$
ここで$u(\nu, T)$は、周波数$\nu$、絶対温度$T$(ケルビン)における分光エネルギー密度です。これは単位周波数幅当たりの量で、$u(\nu, T)\,d\nu$が周波数区間$[\nu, \nu+d\nu]$に含まれる単位体積当たりのエネルギーを表します。また、$k$はボルツマン定数、$c$は光速度、$h$はプランク定数です。
低周波極限と高周波極限において、プランクの公式から他の公式が得られます。
低周波極限$h\nu \ll kT$では、$e^{h\nu/kT} \approx 1 + h\nu/kT$より
$$ u(\nu, T) = \frac{8\pi h\nu^3}{c^3} \frac{1}{e^{h\nu/kT} - 1} \approx \frac{8\pi h\nu^3}{c^3} \frac{1}{h\nu/kT} = \frac{8\pi\nu^2}{c^3} kT $$
高周波極限$h\nu \gg kT$では、$e^{h\nu/kT} - 1 \approx e^{h\nu/kT}$より
$$ u(\nu, T) = \frac{8\pi h\nu^3}{c^3} \frac{1}{e^{h\nu/kT} - 1} \approx \frac{8\pi h\nu^3}{c^3} \frac{1}{e^{h\nu/kT}} $$
これらの公式は共通の構造を持っています。
$$ u(\nu,T) = \frac{8\pi\nu^2}{c^3} \langle E \rangle $$
この式は2つの部分からなります。
モード密度:$\dfrac{8\pi\nu^2}{c^3}$
周波数$\nu$における単位体積・単位周波数幅当たりのモード(振動パターン)の数を表します。
平均エネルギー:$\langle E \rangle$
周波数$\nu$におけるモードが持つ平均エネルギーです。
モード密度は、3次元空間内に閉じ込められた電磁波の定常波の数を数えることで導出されます。ボルツマン定数$k$との混同を避けるため、本節では波数を$\kappa$で表します。
一辺$L$の完全導体壁で囲まれた立方体空洞を考えます。壁面で電場の接線成分が$0$になるという境界条件により、定常波の各方向の波数は$\pi/L$刻みになります。定常波では波数の符号だけが異なるものは同じモードなので、成分が非負($\kappa_x, \kappa_y, \kappa_z \geq 0$)のものだけを数えます。これは波数空間の第1八分空間(全空間の$1/8$)に当たります。空洞が十分大きいとして、格子点の数を球殻の体積で近似します。境界上のモードの細かな補正は省略します。
波数空間で波数の大きさが$\kappa$と$\kappa+d\kappa$の間にある状態数は、球殻の体積要素
$$ d\left(\frac{4}{3}\pi \kappa^3\right) = 4\pi \kappa^2 d\kappa $$
を状態1つあたりの体積$(\pi/L)^3$と第1象限への制限$8$で割った値になります。電磁波には2つの独立した偏光自由度があるため、モード数は
$$ 4\pi \kappa^2 d\kappa \cdot \frac{1}{(\pi/L)^3} \cdot \frac{1}{8} \cdot 2 =\frac{L^3 \kappa^2}{\pi^2} d\kappa $$
となります。これを体積$L^3$で割って単位体積あたりに直し、関係式$\kappa=2\pi\nu/c$を用いて周波数に変換すれば、$d\nu$の係数としてモード密度が得られます。
$$ \frac{\kappa^2}{\pi^2} d\kappa =\left(\frac{2\pi\nu}{c}\right)^2 \frac{1}{\pi^2} \left(\frac{2\pi}{c}d\nu\right) =\frac{8\pi\nu^2}{c^3} d\nu $$
各公式の違いは、平均エネルギー$\langle E \rangle$にあります。
レイリー=ジーンズの公式:$\langle E \rangle = kT$
ウィーンの公式:$\langle E \rangle = \dfrac{h\nu}{e^{h\nu/kT}}$
プランクの公式:$\langle E \rangle = \dfrac{h\nu}{e^{h\nu/kT}-1}$
前述のように、プランクの公式の高周波極限がウィーンの公式です。
レイリー=ジーンズの公式とプランクの公式の違いは、単に低周波極限というだけに留まらず、取り得るエネルギーの違いに由来します。
レイリー=ジーンズの公式の平均エネルギーは、古典統計力学による振動子の平均エネルギーから求められます。
$$ \langle E \rangle = \frac{\int_0^{\infty} E \cdot e^{-E/kT} dE}{\int_0^{\infty} e^{-E/kT} dE} = kT $$
プランクの公式の平均エネルギーは、エネルギーを離散的な値$E_n=nh\nu$($n = 0,1,2,\ldots$)に制限することにより(エネルギーの量子化)、計算を積分から離散準位についての総和に変えることで求められます。ここでは零点エネルギーを除き、基底状態から測った熱的な励起エネルギーを$E_n$としています。
$$ \langle E \rangle = \frac{\sum_{n=0}^{\infty} E_n \cdot e^{-E_n/kT}}{\sum_{n=0}^{\infty} e^{-E_n/kT}} = \dfrac{h\nu}{e^{h\nu/kT}-1} $$
歴史的には、プランクは1900年にまず実験データを説明できる放射式を提案し、その統計的な導出の中でエネルギー要素$h\nu$を導入しました。これは、光そのものをエネルギー量子として捉える1905年のアインシュタインの光量子仮説とは区別されます。本記事では、各モードのエネルギーを量子化する現代的な立場から導出します。
レイリー=ジーンズの公式の平均エネルギーは、統計力学におけるボルツマン分布によって求められます。
統計力学では、熱平衡状態にある系においてエネルギー$E$を持つ各微視的状態の重みは、ボルツマン因子$e^{-E/kT}$に比例します。これがボルツマン分布の基本的な考え方です。ここで扱う1つの古典的調和振動子では、単位エネルギー幅当たりの状態の数(状態密度)がエネルギーによらず一定なので、エネルギーの確率密度$P(E)$もボルツマン因子に比例します。$P(E)\,dE$が、エネルギーが区間$[E, E+dE]$に入る確率です。
$$ P(E) \propto e^{-E/kT} $$
この表現を確率密度として使うため、全確率が$1$になるように規格化します。
$$ P(E) = \frac{e^{-E/kT}}{\int_0^{\infty} e^{-E/kT} dE} $$
物理量の平均値は、その量に確率密度を掛けて積分することで求められます。これにより平均エネルギー$\langle E \rangle$が求まります。
$$ \langle E \rangle = \int_0^{\infty} E \cdot P(E) dE = \frac{\int_0^{\infty} E \cdot e^{-E/kT} dE}{\int_0^{\infty} e^{-E/kT} dE} $$
変数変換$x = \dfrac{E}{kT}$より、分母の積分を計算します。
$$ \int_0^{\infty} e^{-E/kT} dE = kT \int_0^{\infty} e^{-x} dx = kT [-e^{-x}]_0^{\infty} = kT $$
同じ変数変換で分子の積分も計算します。
$$ \int_0^{\infty} E \cdot e^{-E/kT} dE = \int_0^{\infty} xkT \cdot e^{-x} kT \, dx = (kT)^2 \int_0^{\infty} x \cdot e^{-x} dx $$
積分の部分は、部分積分で計算できます。表面項(中辺第1項)が消える($0$になる)のに注意してください。
$$ \int_0^{\infty} x \cdot e^{-x} dx = [x \cdot (-e^{-x})]_0^{\infty} + \int_0^{\infty} e^{-x} dx = 0+[-e^{-x}]_0^{\infty} = 1 $$
ここまでの結果を$\langle E \rangle$の式に代入します。
$$ \langle E \rangle = \frac{\int_0^{\infty} E \cdot e^{-E/kT} dE}{\int_0^{\infty} e^{-E/kT} dE} = \frac{(kT)^2}{kT}=kT $$
プランクの公式の平均エネルギーを計算します。
$$ \langle E \rangle = \frac{\sum_{n=0}^{\infty} E_n \cdot e^{-E_n/kT}}{\sum_{n=0}^{\infty} e^{-E_n/kT}} = \frac{\sum_{n=0}^{\infty} nh\nu \cdot e^{-nh\nu/kT}}{\sum_{n=0}^{\infty} e^{-nh\nu/kT}} $$
$q = e^{-h\nu/kT}\ (0< q<1)$とおいて、分母について等比数列の和の公式を適用します。
$$ \sum_{n=0}^{\infty} e^{-nh\nu/kT} = \sum_{n=0}^{\infty} q^n = \frac{1}{1-q} $$
分子を整理します。
$$ \sum_{n=0}^{\infty} nh\nu \cdot e^{-nh\nu/kT} = h\nu \sum_{n=0}^{\infty} n q^n $$
総和の部分は、$q\dfrac{d}{dq}q^n=nq^n$より計算できます。べき級数は収束半径の内部($|q|<1$)で項別に微分できることを用います。
$$ \sum_{n=0}^{\infty} n q^n = q \frac{d}{dq}\sum_{n=0}^{\infty} q^n = q \frac{d}{dq}\frac{1}{1-q} = \frac{q}{(1-q)^2} $$
ここまでの結果を$\langle E \rangle$の式に代入します。
$$ \langle E \rangle = \frac{h\nu \cdot \frac{q}{(1-q)^2}}{\frac{1}{1-q}} = h\nu \left( \frac{q}{1-q} \right) = h\nu \left( \frac{1}{\frac{1}{q}-1} \right) = \frac{h\nu}{e^{h\nu/kT}-1} $$
プランクの公式は、周波数$\nu$における単位体積・単位周波数幅当たりのエネルギー密度を表します。
$$ u(\nu, T) = \frac{8\pi h\nu^3}{c^3} \frac{1}{e^{h\nu/kT} - 1} $$
これをすべての周波数にわたって積分することで全エネルギー密度$u(T)$が得られ、ゼータ関数が現れます。
$$ u(T) = \int_0^\infty u(\nu, T) d\nu = \frac{48\pi k^4 T^4}{c^3 h^3} \zeta(4) $$
この積分を計算するための変数変換として$x = \dfrac{h\nu}{kT}$とおけば、$\nu = \dfrac{kT}{h}x$となります。
$$ u(T) = \int_0^\infty \frac{8\pi h}{c^3} \left(\dfrac{kT}{h}x\right)^3 \frac{1}{e^x - 1} d\left(\frac{kT}{h}x\right) = \frac{8\pi k^4 T^4}{c^3 h^3} \int_0^\infty \frac{x^3}{e^x - 1} dx $$
積分の中を計算します。$x>0$では$0< e^{-x}<1$なので、等比級数で展開できます。
$$ \frac{1}{e^x - 1} = \frac{1}{e^x(1 - e^{-x})} = \frac{e^{-x}}{1 - e^{-x}} = e^{-x} \sum_{n=0}^{\infty} e^{-nx} = \sum_{n=1}^{\infty} e^{-nx} $$
これを積分に戻して計算を進めます。各項$x^3e^{-nx}$は非負なので、積分と総和の順序を交換できます(単調収束定理によります。証明には立ち入りません)。
$$ \int_0^\infty \frac{x^3}{e^x - 1} dx = \int_0^\infty x^3 \sum_{n=1}^{\infty} e^{-nx} dx = \sum_{n=1}^{\infty} \int_0^\infty x^3 e^{-nx} dx $$
積分の部分に、部分積分を3回適用します。表面項はすべて消えます。
$$ \begin{alignedat}{2} \int_0^\infty x^3 e^{-nx} dx &= 0 - \int_0^\infty 3x^2 \left(-\frac{1}{n}e^{-nx}\right) dx &&= \frac{3}{n}\int_0^\infty x^2 e^{-nx} dx \\ \int_0^\infty x^2 e^{-nx} dx &= 0 - \int_0^\infty 2x \left(-\frac{1}{n}e^{-nx}\right) dx &&= \frac{2}{n}\int_0^\infty x e^{-nx} dx \\ \int_0^\infty x e^{-nx} dx &= 0 - \int_0^\infty \left(-\frac{1}{n}e^{-nx}\right) dx &&= \frac{1}{n}\int_0^\infty e^{-nx} dx \\ \end{alignedat} $$
最後の積分は直接計算できます。
$$ \int_0^\infty e^{-nx} dx = \left[ -\frac{1}{n}e^{-nx} \right]_0^\infty = \frac{1}{n} $$
この結果から逆にたどります。
$$ \begin{alignedat}{2} \int_0^\infty x e^{-nx} dx &= \frac{1}{n} \cdot \frac{1}{n} &&= \frac{1}{n^2} \\ \int_0^\infty x^2 e^{-nx} dx &= \frac{2}{n} \cdot \frac{1}{n^2} &&= \frac{2}{n^3} \\ \int_0^\infty x^3 e^{-nx} dx &= \frac{3}{n} \cdot \frac{2}{n^3} &&= \frac{6}{n^4} \\ \int_0^\infty \frac{x^3}{e^x - 1} dx &= \sum_{n=1}^{\infty} \frac{6}{n^4} &&= 6\zeta(4) \end{alignedat} $$
ここでリーマンゼータ関数が現れました。
実数$s>1$に対して、次のように定める。
$$
\zeta(s) = \sum_{n=1}^{\infty} \frac{1}{n^s}
$$
$\zeta(4) = \dfrac{\pi^4}{90}$を代入します(この値の導出は次節で行います)。
$$ u(T) = \frac{8\pi k^4 T^4}{c^3 h^3} \int_0^\infty \frac{x^3}{e^x - 1} dx = \frac{8\pi k^4 T^4}{c^3 h^3} \cdot 6 \cdot \frac{\pi^4}{90} = \frac{8\pi^5 k^4 T^4}{15 c^3 h^3} $$
これはシュテファン=ボルツマンの法則に対応します。
$$ u(T) = \dfrac{4\sigma}{c}T^4 $$
$$ \sigma=\dfrac{2\pi^5 k^4}{15 c^2 h^3} $$
シュテファン=ボルツマンの法則は、通常は黒体の表面から単位面積・単位時間当たりに放射されるエネルギー(放射発散度)$M(T)$について$M(T)=\sigma T^4$の形で述べられます。等方的な空洞放射では$M(T)=\dfrac{c}{4}u(T)$の関係があり、上の形と同じ内容になります。この関係の導出には立ち入りません。
$\zeta(4) = \dfrac{\pi^4}{90}$を示すには、フーリエ級数展開によるパーセヴァルの等式を用いる方法が一般的です。
フーリエ係数とパーセヴァルの等式は、次の規約で用います。
区間$[-\pi, \pi]$上の2乗可積分な実数値関数$f$のフーリエ係数を
$$ a_n = \frac{1}{\pi}\int_{-\pi}^{\pi} f(t)\cos(nt)\,dt, \qquad b_n = \frac{1}{\pi}\int_{-\pi}^{\pi} f(t)\sin(nt)\,dt $$
とし、フーリエ級数を$\dfrac{a_0}{2} + \displaystyle\sum_{n=1}^\infty \bigl(a_n\cos(nx) + b_n\sin(nx)\bigr)$とすると、次が成り立つ。
$$ \frac{1}{2\pi}\int_{-\pi}^{\pi} f(x)^2\,dx = \frac{a_0^2}{4} + \frac{1}{2}\sum_{n=1}^\infty \left(a_n^2 + b_n^2\right) $$
右辺の係数$1/4$と$1/2$は、定数関数・余弦・正弦の$L^2$ノルムがそれぞれ$\sqrt{2\pi}$、$\sqrt{\pi}$、$\sqrt{\pi}$であることを織り込んだものです。正規直交基底$\{e_j\}$を用いて$\|f\|_2^2 = \sum_j \langle f, e_j\rangle^2$と書く形と同じ等式です。
区間$[-\pi, \pi]$で定義された関数$x^2$を考えます。これは偶関数なので$b_n = 0$となり、フーリエ級数は余弦項のみで表されます。
$$ \begin{aligned} x^2 &= \frac{1}{2\pi} \int_{-\pi}^\pi t^2 \, dt + \sum_{n=1}^\infty \left(\frac{1}{\pi} \int_{-\pi}^\pi t^2 \cos(nt) \, dt\right) \cos(nx) \\ &= \frac{\pi^2}{3} + \sum_{n=1}^\infty \frac{4 (-1)^n}{n^2} \cos(nx) \end{aligned} $$
すなわち$a_0 = \dfrac{2\pi^2}{3}$、$a_n = \dfrac{4(-1)^n}{n^2}$です。この等式は少なくとも$L^2$の意味で成り立ち、以下の計算にはそれで十分です。なお、この関数では$2\pi$周期の延長が連続で区分的に滑らかなため、各点でも成り立ちます。
パーセヴァルの等式を用いて、フーリエ係数の2乗和と$L^2$ノルムの2乗を結びつけます。
$$ \begin{aligned} \frac{1}{2\pi} \int_{-\pi}^\pi x^4 \, dx &= \frac{1}{4}\left(\frac{2\pi^2}{3}\right)^2 + \frac{1}{2} \sum_{n=1}^\infty \left(\frac{4 (-1)^n}{n^2}\right)^2 \\ &= \frac{\pi^4}{9} + \frac{1}{2} \sum_{n=1}^\infty \frac{16}{n^4} \\ \frac{\pi^4}{5} &= \frac{\pi^4}{9} + 8\zeta(4) \\ 8\zeta(4) &= \frac{\pi^4}{5} - \frac{\pi^4}{9} \\ \zeta(4) &= \frac{\pi^4}{8}\left(\frac{1}{5}-\frac{1}{9}\right)=\frac{\pi^4}{90} \end{aligned} $$
本記事では、黒体放射を表すレイリー=ジーンズの公式・ウィーンの公式・プランクの公式が、いずれもモード密度と平均エネルギーの積という共通の構造を持つことを確認しました。平均エネルギーをボルツマン分布で求めるとき、エネルギーを連続的に扱って積分するとレイリー=ジーンズの公式の$kT$が得られ、エネルギーを$E_n = nh\nu$に量子化して総和を取るとプランクの公式が得られます。
さらに、プランクの公式を全周波数にわたって積分すると、リーマンゼータ関数$\zeta(4)$が現れ、シュテファン=ボルツマンの法則が得られました。$\zeta(4)=\pi^4/90$は、パーセヴァルの等式を用いて導きました。
$$ u(\nu,T) = \frac{8\pi\nu^2}{c^3} \langle E \rangle $$
$$ \langle E \rangle = \frac{\sum_{n=0}^{\infty} nh\nu \cdot e^{-nh\nu/kT}}{\sum_{n=0}^{\infty} e^{-nh\nu/kT}} = \frac{h\nu}{e^{h\nu/kT}-1} $$
$$ u(T) = \int_0^\infty u(\nu, T)\,d\nu = \frac{48\pi k^4 T^4}{c^3 h^3} \zeta(4) = \frac{8\pi^5 k^4 T^4}{15 c^3 h^3} = \frac{4\sigma}{c} T^4 $$