統計力学⑫ デバイ比熱
デバイ模型を用いたデバイ比熱・アインシュタイン比熱との違いは?
今回はデバイ模型を紹介します。アインシュタイン模型ではすべての振動子が同じ振動数で振動しているという点が非現実的なのでした。そこで、振動数についても一様ではない分布があるとして考えます。波数から状態数を考える
周期的境界条件を考えてみます。3次元で、周期が全方向に共通で$L$だとすれば、$$
\left\{
\begin{align*}
k_x&=\dfrac{2n_x\pi}{L}\\
k_y&=\dfrac{2n_y\pi}{L}\\
k_z&=\dfrac{2n_z\pi}{L}
\end{align*}
\right.$$
また、波数の大きさ(ノルム)$k$は以下のようにあらわせます。
\begin{align*}
k=\sqrt{k_x^2+k_y^2+k_z^2}
\end{align*}
デバイ近似を行う
実際は波数と角振動数には固体物理の格子振動などからもっと厳密に求めるべきですが、以下のような関係が成り立つと近似的に仮定しましょう。ただし、$v$は伝搬速度\begin{align*}
\omega\approx vk
\end{align*}
この近似をデバイ近似といいます。さて、調和振動子一個のエネルギーは、$\omega$に比例します。つまり、すべての$\omega$を許したら、無限のエネルギーをもつような場合も考えなければいけないので現実的ではないです。よって、限界値$\omega_D$を持つと仮定しましょう。
状態密度を求める
このとき、$k$も限界値を持つとして、$\omega_D=vk_D$を満たす$k_D$を考えます。さて、この波数$k_D$以下の取りうる波数の数について考えます。\begin{align*}
\sum_{k\leq k_D}1
\end{align*}
波数空間内では波数体積ある一つの軸方向($i=x,y,z$)に関して
\begin{align*}
\Delta k_i=\dfrac{2(n_i+1)\pi}{L}-\dfrac{2n_i\pi}{L}=\dfrac{2\pi}{L}
\end{align*}
の幅があることになります。つまり、
\begin{align*}
\left(\dfrac{L}{2\pi}\right)^3\sum_{k\leq k_D}\Delta k_x \Delta k_y \Delta k_z\to \left(\dfrac{L}{2\pi}\right)^3 \int_0^{k_D} dk_x dk_ydk_z
\end{align*}
極座標に置き換えて、角度依存はないものとすれば(波数に関してどの方向にも対称とすれば)
\begin{align*}
\left(\dfrac{L}{2\pi}\right)^3\int_0^{k_D}4\pi k^2dk
&=\left(\dfrac{L}{2\pi}\right)^3\dfrac{4\pi k^3_D}{3} \\
&=\dfrac{L^3\omega_D^3}{6\pi^2v^3}
\end{align*}
これが角振動数が$\omega_D$以下の状態数です。さて、この式を$\omega_D$で微分すれば状態密度がでて、$\omega_D$を$\omega$で置き換えると、
\begin{align}
D(\omega)
&=\dfrac{L^3\omega^2}{2\pi^2 v^3} \label{eq:1}
\end{align}
また、$\omega_D$以下の状態数が今考えている調和振動子の数$N$と一致する(つまり、エネルギーを低いほうからつめていった基底状態を考える)とすれば、
\begin{align}
\dfrac{L^3\omega_D^3}{6\pi^2 v^3}=N \label{eq:2}
\end{align}
が得られます。さて、\eqref{eq:1}から\eqref{eq:2}をもちいて、$L^3/\pi^2v^3$を消去すると、
\begin{align}
D(\omega)=\dfrac{3N}{\omega_D^3}\omega^2 \label{eq:3}
\end{align}
ときれいな形になりました。
分配関数を計算する
分配関数の性質として、系1,2の分配関数$Z_1,Z_2$について、それぞれの系が独立ならば系1,2を合わせた分配関数は、\begin{align*}
Z=Z_1Z_2
\end{align*}
と表せます。これを用いて解析します。角振動数$\omega$の調和振動子の分配関数$Z_0(\omega)$を考えます。前回の記事の結果を借りると、
\begin{align}
Z_0(\omega)=e^{-\frac{3}{2}\beta\hbar\omega}\left(\dfrac{1}{1-e^{-\beta\hbar\omega}}\right)^3 \label{eq:4}
\end{align}
ここで、角振動数が$\omega$から$\omega+\Delta\omega$の間にある振動子の数は、$D(\omega)\Delta\omega$となるので、全体の分配関数$Z$は、
\begin{align*}
Z=\prod_{0\leq \omega\leq \omega_D} \left\{Z_0(\omega)\right\}^{D(\omega)\Delta\omega}
\end{align*}
いま、$\Delta\omega$が指数部にあり、さらに積の形なので積分として処理ができません。以下のように対処します。
\begin{align*}
\ln{Z}&=\sum_{0\leq\omega\leq \omega_D}D(\omega)\ln{Z_0(\omega)}\Delta\omega \\
&\to \int_0^{\omega_D}D(\omega)\ln{Z_0(\omega)}d\omega
\end{align*}
\eqref{eq:3},\eqref{eq:4}を用いて、
\begin{align*}
\ln{Z}&=\dfrac{3N}{\omega^3_D}\int_0^{\omega_D} d\omega\ \omega^2 \ln{\left\{e^{-\frac{3}{2}\beta\hbar\omega}\left(\dfrac{1}{1-e^{-\beta\hbar\omega}}\right)^3\right\}} \\
&=\dfrac{3N}{\omega_D^3}\int_0^{\omega_D}d\omega \ \omega^2\left\{-\dfrac{3}{2}\beta\hbar\omega-3\ln(1-e^{-\beta\hbar\omega})\right\} \\
&=-\left[\dfrac{9N\beta\hbar\omega_D}{8}+\dfrac{9N}{\omega_D^3}\int_0^{\omega_D}d\omega\ \omega^2\ln{(1-e^{-\beta\hbar\omega})}\right]
\end{align*}
分配関数から内部エネルギーを求める
内部エネルギー$E$を分配関数$Z$から導出しましょう。公式を用いて、\begin{align*}
E
&=-\dfrac{\partial}{\partial \beta}\ln{Z} \\
&=\dfrac{\partial}{\partial \beta}\left[\dfrac{9N\beta\hbar\omega_D}{8}+\dfrac{9N}{\omega_D^3}\int_0^{\omega_D}d\omega\ \omega^2\ln{(1-e^{-\beta\hbar\omega})}\right] \\
&=\dfrac{9N\hbar\omega_D}{8}+\dfrac{9N}{\omega_D^3}\int_0^{\omega_D}d\omega\ \omega^2\dfrac{\hbar\omega e^{-\beta\hbar\omega}}{1-e^{-\beta\hbar\omega}} \\
&=\dfrac{9N\hbar\omega_D}{8}+\dfrac{9N\hbar}{\omega_D^3}\int_0^{\omega_D} \dfrac{\omega^3}{e^{\beta\hbar\omega}-1}d\omega
\end{align*}
というようにもとまりました。
比熱を計算する
比熱を計算します。\begin{align*}
C_V&=\dfrac{\partial E}{\partial T} \\
&=\dfrac{d\beta}{dT}\dfrac{\partial E}{\partial\beta} \\
&=-\dfrac{1}{k_BT^2}\dfrac{\partial}{\partial\beta}\left\{\dfrac{9N\hbar\omega_D}{8}+\dfrac{9N\hbar}{\omega_D^3}\int_0^{\omega_D} \dfrac{\omega^3}{e^{\beta\hbar\omega}-1}d\omega \right\} \\
&=-\dfrac{9N\hbar}{k_BT^2\omega_D^3}\dfrac{\partial}{\partial \beta}\left(\int_0^{\omega_D}\dfrac{\omega^3}{e^{\beta\hbar\omega}-1}d\omega \right) \\
&=-\dfrac{9N\hbar}{k_BT^2\omega_D^3}\left\{\int_0^{\omega_D}\left(-\dfrac{\omega^3\cdot e^{\beta\hbar\omega} \hbar\omega}{(e^{\beta\hbar\omega}-1)^2}\right)d\omega \right\} \\
&=\dfrac{9Nk_B\hbar^2 \beta^2}{\omega_D^3}\left\{ \int_0^{\omega_D}\dfrac{\omega^4\cdot e^{\beta\hbar\omega}}{(e^{\beta\hbar\omega}-1)^2}d\omega \right\}
\end{align*}
ただし、ここでは$\beta$,$\hbar$などの余計な係数が出てきてうっとうしいですね。$x=\beta\hbar\omega$として置換します。$\beta\hbar d\omega=dx$なので、
\begin{align*}
C_V&=\dfrac{9Nk_B}{(\beta\hbar\omega_D)^3}\int_0^{\omega_D}\dfrac{(\beta\hbar\omega)^4e^{\beta\hbar\omega}}{(e^{\beta\hbar\omega}-1)^2}\beta\hbar d\omega \\
&=\dfrac{9Nk_B}{(\beta\hbar\omega_D)^3}\int_0^{\beta\hbar\omega_D}\dfrac{x^4e^x}{(e^x-1)^2}dx
\end{align*}
デバイ温度を導入する
先ほどから出ている$\beta\hbar\omega_D$という量ですが、$x$が指数であり無次元量ということ、\begin{align*}
x=\beta\hbar\omega_D&=\dfrac{\hbar\omega_D}{k_BT}
\end{align*}
とできることから、温度の次元を持つ量
\begin{align*}
T_D=\dfrac{\hbar\omega_D}{k_B}
\end{align*}
を定義できます。これをデバイ温度といいます。$x$$=\beta\hbar\omega_D$$=T_D/T$で、
\begin{align*}
C_V=9Nk_B\left(\dfrac{T}{T_D}\right)^3 \int_0^{\frac{T_D}{T}} \dfrac{x^4e^x}{(e^x-1)^2}dx
\end{align*}
低温極限を考える
低温極限を取ると、$T_D/T\to \infty$となります。公式\begin{align*}
\int_{-\infty}^\infty \dfrac{x^4e^x}{(e^x-1)^2}dx=\dfrac{4\pi^4}{15}
\end{align*}
を用いれば、
\begin{align*}
C_V=9Nk_B\left(\dfrac{T}{T_D}\right)^3\dfrac{4\pi^4}{15}=\dfrac{12\pi^4Nk_B}{5}\left(\dfrac{T}{T_D}\right)^3
\end{align*}
となり、$T^3$に比例する実験結果と一致しました。
高温極限を考える
高温極限、つまり、$T/T_D\to 0$を考えます。$e^x$をマクローリン展開して、最低次数の項まで見ると、\begin{align*}
C_V&=9Nk_B\left(\dfrac{T}{T_D}\right)^3 \int_0^{\frac{T_D}{T}} \dfrac{x^4e^x}{(e^x-1)} dx\\
&\approx 9Nk_B\left(\dfrac{T}{T_D}\right)^3 \int_0^{\frac{T_D}{T}}x^2\ dx \\
&=9Nk_B
\end{align*}
これでデュロン・プティの法則と一致することがわかります。