2025年度 熱統計力学(A科) 第7章
2/3に試験の予定。
ここまでの簡易的なまとめ
これだけではないので注意(特にカルノーサイクル、理想気体のエントロピーなど)。
カノニカルアンサンブル
以前、エントロピー\(S\)は同一のマクロな量に対応する場合の数\(\Omega\)と、以下のような関係があることを示した。
\[\begin{equation}
S = k\log\Omega
\end{equation}\] 孤立系において、\(\Omega\)は時間と共に単調増加し平衡状態で最大をとる場合の数である。系1(\(\Omega_1\))と系2(\(\Omega_2\))をあわせた系での場合の数は\(\Omega_1\Omega_2\)となる。一方、エントロピーは同様の性質をもち、あわせた系では\(S_1 +
S_2\)となる。この２つをつなぐ式が上記のものである。
この式はミクロな変数\(\Omega\)とマクロな変数\(S\)をつなぐ重要な式であった。ここでは、ミクロな状態に注目し、マクロな状態との対応関係を考えていく。
一定の温度をもつ系の位相空間上の確率分布を考える。ここでの位相空間とは、１つの粒子\(i\)の状態を位置\(q_i\), 運動量\(p_i\)により表すことを意味する(例えば理想気体ならば、運動エネルギー\(\frac{1}{2}mv^2\)は、運動量\(p = mv\)を用いて、\(\frac{p^2}{2m}\)と表すことができる。位置と運動量があれば、その粒子がどの位置にあり、どのようなエネルギーをもち、さらに速度の方向かつ振幅も表すことができる)。両者はハミルトアニアンというエネルギー関数を構成する変数であり、位相空間上の点それぞれにいわゆるエネルギーを割り当てることができる。
教科書の図7.1を参照に、系を\(S\)、熱浴を\(R\)とする。ここでは、系と熱浴を合わせた閉じた系を考える。すなわち、\(S\)と\(R\)を合わせた系のエネルギー\(E\)は、 \[\begin{equation}
E = E_{\rm R} + E_{\rm S} = \rm{const.}
\end{equation}\]
であり、閉じた系であるためこれは一定値である。熱浴とは系を等温に保つために熱のやり取りをする巨大な系であり、
\[\begin{equation}
\frac{E_{\rm S}}{E} = \frac{E_{\rm S}}{E_{\rm S} + E_{\rm R}} =
\frac{E_{\rm S}/E_{\rm R}}{1 + E_{\rm S}/E_{\rm R}} \ll 1
\end{equation}\]
とする。ここで、熱浴と接するため温度は一定であるが、粒子は動き続けており、系のエネルギー\(E_{\rm
S}\)は一定ではない。ここで、理想気体であれば内部エネルギーは温度と粒子数が一定ならば決まることと反するように一見みえる。しかし、内部エネルギーは運動エネルギーを粒子間平均したものであり[さらにいうと時間平均したものでもあり]、とある時間スナップショットにおける複数の粒子の運動エネルギーはつねに揺らぐ。そのため、この系のエネルギー\(E_S\)は揺らぐと想定しても矛盾はない。したがって系のエネルギーの実現値は、何かしらの確率分布に従うものと想定できる。系\(S\)のエネルギーが実現値\(E_i\)となる確率\(p_i\)を求める。
\(p_i\)は、エネルギー\(E_i\)[これまでの内部エネルギー\(U\)と異なる。\(E_i\)の期待値が\(U\)]と両立しうる場合の数\(\Omega_{\rm
S}(E_i)\)に比例する(例えば、目が1,1,1,2,3,4のサイコロであれば、1が確率は3/6であるように)。すなわち、
\[\begin{equation}
p_i \propto \Omega_{\rm S}(E_i)
\end{equation}\] 系のエネルギーの実現値が\(E_i\)となるならば，熱浴のエネルギーは\(E -
E_i\)と決まる。すなわち、熱浴におけるミクロな状態の場合の数を考えると、
\[\begin{equation}
p_i \propto \Omega_{\rm R}(E - E_i)
\end{equation}\]
となる。便宜的にエントロピーを考えると、系のエネルギーの実現値が\(E_i\)となるとき、熱浴のエントロピーは、\(E_i\)が\(E\)に比べて十分に小さいと仮定すると、 \[\begin{equation}
S_{\rm R} = k\ln \Omega_{\rm R}(E - E_i) \simeq k\ln \Omega_{\rm R}(E) -
k \frac{\partial}{\partial E}(\ln\Omega_{\rm R}(E)) \times E_i + ... =
k\ln \Omega_{\rm R}(E) - \frac{\partial S_{\rm R}}{\partial E}E_i = k\ln
\Omega_{\rm R}(E) - \frac{E_i}{T}
\end{equation}\] という近似式を得る(\(E_i \ll
1\)を暗に仮定している[ように思える])。
したがって、 \[\begin{equation}
p_i \propto \Omega_{\rm R}(E - E_i) = \Omega_{\rm R}(E)\exp \left( -
\frac{E_i}{kT}\right)
\end{equation}\] を得る。\(E\)は一定であることを想定しているため、\(\Omega_{\rm R}(E)\)も一定。このことから、
\[\begin{equation}
p_i \propto \exp \left( - \frac{E_i}{kT}\right)
\end{equation}\] を得る。さらに、\(p_i\)を、\(\sum_i
p_i = 1\)となるように正規化すると(確率を考えているため)、 \[\begin{equation}
p_i = \frac{\exp \left( - \frac{E_i}{kT}\right)}{\sum_i \exp \left( -
\frac{E_i}{kT}\right)}
\end{equation}\]
を得る。この確率分布を、カノニカル分布と呼ぶ。
多くの場合、\(\beta =
1/(kT)\)として逆温度を定義し、 \[\begin{equation}
p_i = \frac{\exp \left( - \beta E_i\right)}{\sum_i \exp \left( - \beta
E_i\right)}
\end{equation}\]
と書くことが多い。また、カノニカル分布の分母を\(Z\)、すなわち \[\begin{equation}
Z = \sum_i \exp \left( - \beta E_i\right)
\end{equation}\] として、\(Z\)を分配関数と呼ぶ。また、\(E_i\)が連続値である場合は、 \[\begin{equation}
p_i = \frac{\exp \left( - \beta E_i\right)}{\int \exp \left( - \beta
E_i\right) dE_i}
\end{equation}\] \[\begin{equation}
Z = \int \exp \left( - \beta E_i\right) dE_i
\end{equation}\]
として積分を考えれば良い。要は、確率分布の性質を満たすように分配関数を設定すればカノニカル分布として成立する(確率分布の性質は、非負、定義域内で確率変数に関して積分すると1になる、ということ)。
少し脱線
とあるデータ\(\boldsymbol{y}\in\mathbf{R}^{N}\)を計測したとする。このデータには加法的にバイアスのないノイズが加わっており、ノイズの加わっていない観測できない真のデータ\(\boldsymbol{x}\in\mathbf{R}^{N}\)から、\(\boldsymbol{y} = \boldsymbol{x} +
\boldsymbol{\xi}\)として生成されているという確率モデルを考える。このとき、尤度関数は\(p(\boldsymbol{y} |
\boldsymbol{x})\)、事前分布は\(p(\boldsymbol{x})\)、事後分布は \[\begin{equation}
p(\boldsymbol{x} | \boldsymbol{y}) = \frac{p(\boldsymbol{y} |
\boldsymbol{x})p(\boldsymbol{x})}{\int p(\boldsymbol{y} |
\boldsymbol{x})p(\boldsymbol{x}) d \boldsymbol{x}}
\end{equation}\] いま、無情報量事前分布を考えるとすると、 \[\begin{equation}
p(\boldsymbol{x} | \boldsymbol{y}) = \frac{p(\boldsymbol{y} |
\boldsymbol{x})}{\int p(\boldsymbol{y} | \boldsymbol{x}) d
\boldsymbol{x}}
\end{equation}\] これを書き換えると、 \[\begin{equation}
p(\boldsymbol{x} | \boldsymbol{y}) = \frac{\exp(- [-\log
p(\boldsymbol{y} | \boldsymbol{x})])}{\int \exp(- [-\log
p(\boldsymbol{y} | \boldsymbol{x})]) d \boldsymbol{x}}
\end{equation}\] となり、\(E = -\log
p(\boldsymbol{y} | \boldsymbol{x})\)かつ逆温度を\(\beta =
1\)としたときのカノニカル分布と一致する。このことから、統計力学と機械学習(特にベイズ統計)は同一の枠組みとみなすことができ、互いの分野を行き来するような解析が可能となることが期待できる(事実、多くの機械学習手法にたいして統計力学的解析などが行われている)。
(2026. 01. 09ここから):
等温状態にある熱力学系において、エネルギー\(E\)が取りうる状態は以下のカノニカル分布に従う。
\[\begin{equation}
p(E = E_i) = \frac{\exp \left( - \beta E_i\right)}{Z}
\end{equation}\] ただし、\(\beta =
1/(kT)\)は逆温度、\(Z = \sum_i \exp
\left( - \beta
E_i\right)\)は分配関数である。取りうるエネルギーの状態が連続的に変化するならば、
\[\begin{equation}
p(E = E_i) = \frac{\exp \left( - \beta E_i\right)}{Z}, \ Z = \int \exp
\left( - \beta E\right)dE
\end{equation}\]
となる。エネルギーとして抽象的なものを扱ってきたが、以下ではハミルトニアンを用いて速度分布を導出する(ハミルトニアンは位置と運動量の関数であるが、質点かつユークリッド幾何なら運動量でも速度でも大差はない)。
カノニカル分布に基づく、粒子が従う速度分布の導出
カノニカル分布を利用して、理想気体の各粒子が従う速度の分布を求めてみる。1粒子のみを考え、さらに1次元の筒内に動作が拘束されているものと仮定する。このとき、この系のハミルトニアン(ここではエネルギーと同じもの)は、
\[\begin{equation}
H(q, p) = \frac{p^2}{2m}
\end{equation}\]
となる。したがって、カノニカル分布は、粒子1の速度を\(v_1\)、運動量を\(p_1\)とすると、 \[\begin{equation}
p(H(q, p_1) = H(q, p)) = p(v_1 = v) = \frac{1}{Z}\exp \left( - \beta
\frac{1}{2}mv^2\right) = \frac{1}{Z}\exp \left( -\frac{mv^2}{2kT}\right)
\end{equation}\]
ここで、確率変数\(x\)が実現値\(z\)をとる確率が平均\(\mu\)、分散\(\sigma^2\)のガウス分布に従うとき、その確率分布は
\[\begin{equation}
p(x = z) = \frac{1}{\sqrt{2\pi\sigma^2}}\exp \left(
-\frac{1}{2\sigma^2}(z - \mu)^2 \right)
\end{equation}\]
となる。これと照らし合わせると、粒子の速度は、平均0、分散\(kT/m\)のガウス分布であることがわかる。そして、分配関数\(Z\)は、 \[\begin{equation}
Z = \sqrt{2\pi kT/m}
\end{equation}\]
であることがわかる。この分布は1次元に拘束されている場合の話だが、3次元を考えるならば、各次元は独立しているため、
\[\begin{equation}
p(v_{1x} = v_x, v_{1y} = v_y, v_{1z} = v_z) = \left(\frac{m}{2\pi
kT}\right)^{3/2} \exp \left( -\frac{m(v_x^2 + v_y^2 +
v_z^2)}{2kT}\right) = \left(\frac{m}{2\pi kT}\right)^{3/2} \exp \left(
-\frac{m |\boldsymbol{v}|^2 }{2kT}\right)
\end{equation}\]
を得る。この速度分布はマクスウェルの速度分布と呼ばれる(p. 12の例1.2)。0からこれを導くとかなり難しいものの、カノニカル分布を利用すると一瞬で導出することができる。
カノニカル分布に基づく、理想気体の内部エネルギー
マクスウェルの速度分布を利用すると、運動エネルギーの期待値は、 \[\begin{equation}
E[K] = \frac{1}{2}mE[v^2] = \frac{1}{2}m \int v^2 p(v)dv
\end{equation}\] となる。ここで、分散の定義は\(\mathrm{Var}(v) = E[v^2] -
(E[V])^2\)より、\(E[v^2] =
\mathrm{Var}(v) + (E[v])^2\)。上記の\(p(v)\)は平均が0(\(E[v]\)=0)、分散が\(\mathrm{Var}(v) =
kT/m\)のガウス分布であるため、\(E[v^2]
= kT/m\)。
したがって、1粒子1次元の場合は運動エネルギーは、 \[\begin{equation}
E[K] = \frac{1}{2}mE[v^2] = \frac{1}{2}kT
\end{equation}\] となる。これが3次元の場合ならば、\(v_x\), \(v_y\), \(v_z\)は独立なため、 \[\begin{equation}
E[K] = \frac{1}{2}m(E[v_x^2]+E[v_y^2]+E[v_z^2]) = \frac{3}{2}kT
\end{equation}\] となる。\(N\)粒子ならば、各々の粒子が独立なため、
\[\begin{equation}
E[K] = N\times \left( \frac{1}{2}m(E[v_x^2]+E[v_y^2]+E[v_z^2]) \right) =
\frac{3}{2}NkT
\end{equation}\]
となる。理想気体ならば、粒子間の距離は十分に長くポテンシャルエネルギーが0の状況を想定している。したがって、理想気体の内部エネルギー\(U\)は運動エネルギーのみであり、 \[\begin{equation}
U = \frac{3}{2}NkT
\end{equation}\]
となる。これまたかなり難しい導出が必要だったが、カノニカル分布を利用すると容易に導出できる。
分配関数と(ヘルムホルツの)自由エネルギー
カノニカル分布を利用することで、系の平均エネルギー\(U\)を、 \[\begin{equation}
U = \frac{1}{Z}\sum_i E_i \exp \left( - \beta E_i\right)
\end{equation}\]
として計算することができる。また、天下り的ではあるが、 \[\begin{equation}
Z = \sum_i \exp \left( - \beta E_i\right)
\end{equation}\] より、 \[\begin{equation}
\frac{\partial}{\partial \beta}Z = -\sum_i E_i \exp \left( - \beta
E_i\right)
\end{equation}\] となる。すなわち、 \[\begin{equation}
-\frac{\partial}{\partial \beta} \ln Z = \frac{1}{Z}\sum_i E_i \exp
\left( - \beta E_i\right) = U
\end{equation}\] となる。
一方、内部エネルギーの\(S\)への依存性を\(T\)への依存性へとルジャンドル変換することで、以下のヘルムホルツの自由エネルギーを得る。
\[\begin{equation}
F = U - TS = U - \frac{S}{k\beta}
\end{equation}\] ここから、 \[\begin{equation}
\beta F = \beta U - \frac{S}{k}
\end{equation}\] \[\begin{equation}
\frac{\partial}{\partial\beta}(\beta F) = \frac{\partial}{\partial\beta}
(\beta U) = U
\end{equation}\] となり、上記の式を結びつけることで、 \[\begin{equation}
-\frac{\partial}{\partial \beta} \ln Z = U =
\frac{\partial}{\partial\beta}(\beta F)
\end{equation}\]
となり、すなわち以下の超重要な式を得ることができる。 \[\begin{equation}
\ln Z = - \beta F
\end{equation}\] 以下の形で表すことが多い印象である。 \[\begin{equation}
F = U - TS = -kT \ln Z
\end{equation}\]
また、ここからエントロピー\(S\)について求めていく。その前に少しだけ情報理論に関しておさらい。確率分布\(p_i\)が与えられたとき、その事情を観測することで得る情報量は
\[\begin{equation}
I = -\log p_i
\end{equation}\] である(あえて\(\log\)と書いているのは、情報理論では底が自然対数ではない場合も多いため)。エントロピーは
\[\begin{equation}
H = -\sum_ip_i\log p_i
\end{equation}\]
である。単純な例として、ベルヌーイ分布を考える(\(x\)は0か1のバイナリ)。 \[\begin{equation}
p(x | q) = q^{x}(1-q)^(1-x)
\end{equation}\] この情報量は \[\begin{equation}
I = -p = -x\log q - (1-x)\log (1-q)
\end{equation}\] 事象\(x=1\)のときは\(-\log q\)、事象\(x=0\)のときは\(-\log
(1-q)\)という情報量を得る。また、エントロピーは \[\begin{equation}
H = -\sum_{x}q^{x}(1-q)^(1-x) (x\log q + (1-x)\log (1-q)) = -q\log q -
(1-q)\log(1-q)
\end{equation}\] であり、エントロピーを最大化するのは\(q=0.5\)であることが簡単な計算からわかる。すなわち、コイントスならば表も裏も同確率のとき、もっとも予測が困難なときにエントロピーは最大となり、ここからもエントロピーが不確かさの尺度と呼ばれる理由がわかる。
以下は本題であり、エントロピーを求めていく。\(p_i = p(E_i)\)とすると、 \[\begin{equation}
S = k\beta(U - F) = k\beta (\sum_i E_ip_i + \frac{1}{\beta}\ln Z)
\end{equation}\] を得る。以下、\(\sum_i
x_ip_i = \langle x \rangle\)とする。ただし、\(\langle y \rangle =
y\)。分配関数は正規化の定数のため、以下のように書いても変化はない。
\[\begin{equation}
S = k (\beta \langle E \rangle+ \ln Z) = -k \langle \ln\exp(-\beta E) -
\ln Z \rangle
\end{equation}\] \(p_i\)はカノニカル分布であることを踏まえて、
\[\begin{equation}
S = -k \langle \ln p \rangle = -k \sum_i p_i \ln p_i
\end{equation}\]
これより、カノニカル分布においては、熱力学でのエントロピーと情報理論でのエントロピーが一致することがわかる。すなわち、カノニカル分布に従う系においては、エントロピーは不確かさを表す。
カノニカル分布に基づくエントロピーの解釈(時間があれば)
以前考えた1粒子かつ1次元の動きに制約された系に注目する。分配関数は、
\[\begin{equation}
p(H) = \frac{1}{Z}\exp(-\frac{mv^2}{2kT})
\end{equation}\] であることから、ガウス分布との対比より \[\begin{equation}
Z = \sqrt{2\pi kT/m}
\end{equation}\] である。自由エネルギーは、 \[\begin{equation}
F = -kT\log Z = -\frac{kT}{2} (\log (2\pi kT) - \log m)
\end{equation}\] となる。自由エネルギーは、\(F = U =TS\)より \[\begin{equation}
\frac{\partial F}{\partial T} = -S
\end{equation}\] を満たすことから、考えている系のエントロピーは、
\[\begin{equation}
S = \frac{k}{2}\left( \ln (2\pi kT) + 1 - \ln m\right)
\end{equation}\]
となる。例えば、温度を絶対零度に近づけていくとき(\(T \to
0\))、速度分布は分散を0に近づけるガウス分布、すなわち値が0になる確率が有限の値、その他の値になる確率が0となるデルタ関数となる。すなわち、絶対零度への極限では、粒子は運動しない。このとき、エントロピーは最小の値をとる。一方、温度を大きくしていくと速度分布の分散が大きくなり、粒子の速度は様々な値をとるようになり、エントロピーは上昇していくことがわかる。このことから、この系におけるエントロピーとは、「粒子がもつ速度を予測できるか否か」もしくは「粒子が取りうる速度の乱雑さ」などと解釈することができる。
レポート課題: 以下の問題を解け。
(1) 平均\(\mu\)、分散\(\sigma^2\)のガウス分布のエントロピーを求めよ。
(2) ハミルトニアンが\(H(q,p) = \frac{1}{2}aq^2 + \frac{1}{2}bq^2\)で与えられる1次元の系を考える(\(a(> 0)\)、\(b(>0)\)は定数、\(q \in (-\infty, \infty)\), \(p \in (-\infty, \infty)\))。この系の分配関数ならびに自由エネルギーを求めよ。
(2026. 01. 16ここから):
等温状態にある熱力学系において、エネルギー\(E\)が取りうる状態は以下のカノニカル分布に従う。
\[\begin{equation}
p(E = E_i) = \frac{\exp \left( - \beta E_i\right)}{Z}
\end{equation}\] ただし、\(\beta =
1/(kT)\)は逆温度、\(Z = \sum_i \exp
\left( - \beta
E_i\right)\)は分配関数である。取りうるエネルギーの状態が連続的に変化するならば、
\[\begin{equation}
p(E = E_i) = \frac{\exp \left( - \beta E_i\right)}{Z}, \ Z = \int \exp
\left( - \beta E\right)dE
\end{equation}\] となる。先週は、自由エネルギー\(F = U - TS\)は、 \[\begin{equation}
F = -kT\ln Z
\end{equation}\]
であることを学んだ。さらに、カノニカル分布においてはエントロピーは情報理論でのエントロピー\(E[-\log
p]\)と一致すること、すなわちエントロピーとは不確かさの尺度であることを学んだ。
ガウス分布のおさらい
理想気体(+調和振動子)においてカノニカル分布を扱うさい、ガウス分布の知識は必須である。念の為、もう一度おさらいしておく。確率変数\(X\)の実現値が\((x+x+\Delta x)\)となる確率を\(p(X=x)\Delta x = p(x)\Delta
x\)と定義する(\(\Delta \ll
1\))。ガウス分布は以下の確率密度関数で定義される。 \[\begin{equation}
p(x|\mu, \sigma^2) =
\frac{1}{\sqrt{2\pi\sigma^2}}\exp\left(-\frac{1}{2\sigma^2}(x-\mu)^2\right)
\end{equation}\] この関数系からわかるように、\(\mu\)を中心として滑らかに減衰する偶関数である。平均は、
\[\begin{equation}
E[x] =
\int_{-\infty}^{\infty}\frac{1}{\sqrt{2\pi\sigma^2}}x\exp\left(-\frac{1}{2\sigma^2}(x-\mu)^2\right)dx
\end{equation}\] であり、\(x' =
\frac{x-\mu}{\sigma}\)と定義すると、 \[\begin{equation}
E[x] =
\int_{-\infty}^{\infty}\frac{1}{\sqrt{2\pi\sigma^2}}x'\exp\left(-\frac{1}{2}[x']^2\right)\sigma
dx + \mu
\end{equation}\]
偶関数より、右辺第一項は0。したがって、平均は\(\mu\)。分散は、 \[\begin{equation}
\mathrm{Var}(x) =
\int_{-\infty}^{\infty}\frac{1}{\sqrt{2\pi\sigma^2}}(x-\mu)^2\exp\left(-\frac{1}{2\sigma^2}(x-\mu)^2\right)dx
\end{equation}\] 同様に\(x'\)を用いると、 \[\begin{equation}
\mathrm{Var}(x) =
\int_{-\infty}^{\infty}\frac{1}{\sqrt{2\pi\sigma^2}}\sigma^2[x']^2\exp\left(-\frac{1}{2\sigma^2}(x-\mu)^2\right)\sigma
dx' = \sigma^2
\end{equation}\]
となる。この平均と分散はどの分野にせよ必須事項である。
次に規格化定数が\(\sqrt{2\pi\sigma^2}\)であることを示す。平均は簡単のため\(\mu=0\)とする。 \[\begin{equation}
p(x|\sigma^2) =
\frac{1}{\sqrt{2\pi\sigma^2}}\exp\left(-\frac{1}{2\sigma^2}x^2\right)
\end{equation}\] 上記と同様に\(x'\)を定義し、計算を進めると、以下を示すことと等価であることがわかる。
\[\begin{equation}
\int\frac{1}{\sqrt{2\pi}}\exp\left(-\frac{1}{2}x^2\right)dx = 1
\end{equation}\] つまり、\(I =
\int\exp\left(-\frac{1}{2}x^2\right)dx =
\sqrt{2\pi}\)を示せば良い。
\[\begin{equation}
I^2 = \int\int\exp\left(-\frac{1}{2}(x^2+y^2)\right)dxdy
\end{equation}\] より、\(x^2 + y^2 =
r^2\), \(x = r\cos\theta\),
\(y =
r\sin\theta\)とおく。変数変換に留意して(ヤコビアンの行列式がかかる)
\[\begin{equation}
I^2 =
\int_0^{\infty}rdr\int_{0}^{2\pi}d\theta\exp\left(-\frac{1}{2}r^2\right)
= 2\pi\times 1
\end{equation}\] ただし、 \[\begin{equation}
\int_0^{\infty}rdr\exp\left(-\frac{1}{2}r^2\right) = -\int_0^{\infty}dr
\left(\frac{d}{dr}\exp\left(-\frac{1}{2}r^2\right)\right) = 1
\end{equation}\] である。したがって、求めたい\(I\)は、\(I^2 =
2\pi\)より、 \[\begin{equation}
\int\exp\left(-\frac{1}{2}x^2\right)dx = \sqrt{2\pi}
\end{equation}\] 適宜変数変換すれば、 \[\begin{equation}
\int\exp\left(-\frac{1}{2\sigma^2}x^2\right)dx = \sqrt{2\pi\sigma^2}
\end{equation}\]
がわかる。この規格化定数はカノニカル分布では分配関数に相当するため、計算では多用する。
カノニカル分布と状態方程式
再び1粒子の系を考える。体積\(V\)の容器に封入されている状況を考え、真面目に分配関数を計算してみる。カノニカル分布は、
\[\begin{equation}
p(H(\boldsymbol{q}_1, \boldsymbol{p}_1) = H(\boldsymbol{q},
\boldsymbol{p})) = \frac{1}{Z}\exp \left( - \frac{1}{2mkT}(p_x^2 + p_y^2
+ p_z^3) \right)
\end{equation}\] であり、分配関数は本来、位置\(q\)、運動量\(p\)に対して積分をする必要がある。 \[\begin{equation}
Z = \int dq_x \int dq_y \int dq_z \int dp_x \int dp_y \int dp_z \ \exp
\left( - \frac{1}{2mkT}(p_x^2 + p_y^2 + p_z^3) \right)
\end{equation}\] となる。ここで、\(\int
dq_x \int dq_y \int dq_z\)は容器の体積\(V\)となる。
そして、 \[\begin{equation}
\int dp_x \int dp_y \int dp_z \ \exp \left( - \frac{1}{2mkT}(p_x^2 +
p_y^2 + p_z^3) \right) = (2\pi mkT)^{3/2}
\end{equation}\]
となる(この計算は難しいので、ガウス分布が積分して1になることから理解すればよい)。すなわち、
\[\begin{equation}
\int dp_x \int dp_y \int dp_z \ \frac{1}{(2\pi mkT)^{3/2}} \exp \left(
- \frac{1}{2mkT}(p_x^2 + p_y^2 + p_z^3) \right) = 1
\end{equation}\] から計算すれば良い。したがって、最終的に \[\begin{equation}
Z = V (2\pi mkT)^{3/2}
\end{equation}\] を得る。
続いて、\(N\)個の粒子から構成される理想気体を考える。\(N\)個の粒子の間には相互作用がなく、互いに独立しているとみなせる。したがって、カノニカル分布は(単純のため実現値のみを確率分布の括弧内に記す)、
\[\begin{equation}
p(H(\boldsymbol{q}, \boldsymbol{p}) = p(H((\boldsymbol{q}_1,
\boldsymbol{p}_1), H((\boldsymbol{q}_2, \boldsymbol{p}_2), ...,
H((\boldsymbol{q}_N, \boldsymbol{p}_N)) = p(H((\boldsymbol{q}_1,
\boldsymbol{p}_1())p(H((\boldsymbol{q}_2,
\boldsymbol{p}_2))...p(H((\boldsymbol{q}_N, \boldsymbol{p}_N))
\end{equation}\]
である。そして、上記と同じ計算より、この系全体のカノニカル分布に対する分配関数は、
\[\begin{equation}
Z = V^N (2\pi mkT)^{\frac{3N}{2}}
\end{equation}\] となる。したがって、自由エネルギーは \[\begin{equation}
F = -kT\log Z = -kT (N\log V + \frac{3N}{2} \log (2\pi mkT)) = -NkT
\left( \log V + \frac{3N}{2} \log (2\pi mkT) \right)
\end{equation}\] となる。
この系の圧力\(p\)は、\(p = -\frac{\partial F}{\partial
V}\)から計算することができる。なぜならば、 \[\begin{equation}
F = U - TS
\end{equation}\] より、 \[\begin{equation}
\frac{\partial F}{\partial V} = \frac{\partial U}{\partial V} = -p
\end{equation}\] を得る(\(dU = -pdV +
TdS\)より)。したがって、 \[\begin{equation}
p = -\frac{\partial F}{\partial V}
\end{equation}\] となる。
\[\begin{equation}
F = -kT\log Z = -NkT \left( \log V + \frac{3N}{2} \log (2\pi mkT)
\right)
\end{equation}\] より、 \[\begin{equation}
p = -\frac{\partial F}{\partial V} = \frac{NkT}{V}
\end{equation}\] すなわち、理想気体の状態方程式、 \[\begin{equation}
pV = NkT
\end{equation}\]
を得る。これはカノニカル分布により、ミクロな系である1個1個の粒子がもつ位置と運動量から、マクロな系である理想気体が満たす方程式が導けたことを意味する(ちなみに、熱力学自体はある程度確立された後に統計力学が進展したため、これは一致するように枠組みが作られているという考えもある)。
調和振動子
\(N\)個の互いに異なる各周波数\(\omega_i\)で振動する調和振動子がカノニカル分布に従うと考える\(i = 1, ...
,N\)(固体内部で発生する格子振動の単純なモデル)。ハミルトニアンは、
\[\begin{equation}
H(q, p) = \sum_i \left( \frac{p_i^2}{2m} + \frac{1}{2}m\omega_i^2q_i^2
\right)
\end{equation}\]
として与えられる。分配関数は理想気体のときと同様に計算可能であり、自由エネルギーやエントロピーを計算可能である(先週のレポート課題参照)。
ここまでの本講義のまとめ
レポート課題: 以下、期末テストの過去問より。
体積 V の容器に封入されている粒子数 N の理想気体を考える。粒子 \(i\) の位置と運動量を各々 \(\boldsymbol{q}_i = (q_{i,x}, q_{i,y}, q_{i,z})\), \(\boldsymbol{p}_i = (p_{i,x}, p_{i,y}, p_{i,z})\)とする(\(i = 1, ..., N\))。以下の問題を解け。
(1) N=1 のとき、カノニカル分布の分配関数を計算せよ。
(2) N \(\neq\) 1 のとき、カノニカル分布の分配関数を計算せよ。
(3) 上記の問い(2)の状況において、ヘルムホルツの自由エネルギーを求めよ。
(4) 上記の問い(2)の状況において、圧力を求め、理想気体の状態方程式を導け。
イジングモデル
\end{equation}
磁性体のモデルとして、上向きと下向きの2状態をとるスピンを利用した模型を考える。この模型をイジングモデルと呼び、ハミルトニアンは
\[\begin{equation}
H = -\sum_{i,j}J_{i,j}\sigma_i\sigma_j - \sum_ih_i\sigma_i
\end{equation}\] で与えられる(\(i = 1,
..., N\))。ただし、\(\sigma_i\)はスピン変数であり-1か1を取る離散変数である。\(J_{i,j}\)はスピン\(i\)とスピン\(j\)の相互作用を表し、\(h_i\)は外部磁場を表す。特に\(J_{i,j}\)の変化により多彩なモデルを作ることができ、一次元格子上に横並びさせ、隣同士に相互作用するものを1次元イジングモデル、2次元格子状に並べ、近接スピン同士で相互作用するものを2次元イジングモデル、…として定義されることが多く、3次元はまだ解かれていないと言われている印象。また、\(J_{i,j}\)にランダム性をもたせるとスピングラス(ガラスの模型)、パターンを埋め込むと2024年にノーベル物理学賞を受賞したJ.J.Hopfield(+甘利俊一先生)が提案した連想記憶モデルになる。計算方法の多様性も踏まえ(平均場近似、ベーテ近似、動的平均場近似、レプリカ法、ラプラス近似、など)、統計力学の花形的模型である。
簡単のため、\(J_{i,j} = 0, h_{i} =
h\)とする一様な外部磁場のみが与えられるイジングモデルを解析してみよう。ハミルトニアンは、
\[\begin{equation}
H = - h\sum_i\sigma_i
\end{equation}\]
となり、スピン同士に相互作用がない状況となる。分配関数は、 \[\begin{equation}
Z = \sum_{\sigma_i}\exp(\beta h \sum_i\sigma_i) =
\sum_{\sigma_i}\prod_i\exp(\beta h \sigma_i)
\end{equation}\] すべてのスピン変数の値に関して足し算すると、
\[\begin{equation}
Z = (\exp(\beta h)+\exp(-\beta h))^N = 2^N\cosh(\beta h)^N
\end{equation}\] となる。自由エネルギーは、 \[\begin{equation}
F = -kT \ln Z = - \frac{1}{\beta} N\ln \cosh (\beta h) + \mathrm{const.}
\end{equation}\] となる。
多くの場合、統計力学では秩序変数(order
parameter)を考える。代表的なものとして、\(m =
\langle\sigma_i
\rangle\)があり(要はスピン変数の期待値)、上記のモデルにおいては、\(J_{ij} =
J\)において(一様な相互作用)ハミルトニアンを下記のように考えることができる(\(o\)はランダウのオー)。 \[\begin{equation}
H = -J\sum_{i,j}(m+\delta\sigma_i)(m+\delta\sigma_j) -h\sum_i\sigma_i =
-J\sum_{i,j}(m^2+m\delta\sigma_i+m\delta\sigma_j) - h\sum_i\sigma_i +
o(\delta\sigma^2)
\end{equation}\] ここで、\(\sum_{i,j}1\)はスピン同士がどのように繋がっているかという模型の制約に依存する。例えば1次元鎖なら\(N\)なので、以下では\(\sum_{i,j}1=N\)を考える。\(\delta\sigma_i\)は\(\delta\sigma_i = \sigma_i -
m\)として再び考えてみると、ハミルトニアンは、 \[\begin{equation}
H = -JNm^2-2Jm\sum_i(\sigma_i-m) - h\sum_i\sigma_i = JNm^2 - (h +
2Jm)\sum_i\sigma_i
\end{equation}\] となる。つまり、見かけ上の外部磁場\((h + 2Jm)\)が加わるような形となる。
カノニカル分布は、 \[\begin{equation}
p(\boldsymbol{\sigma}) = \frac{1}{Z}\exp(-\beta H) =
\frac{1}{Z}\exp(-\beta JNm^2)\exp(-\beta [h+2Jm]\sum_i\sigma_i)
\end{equation}\] であるため、分配関数は、\(\exp(-\beta
JNm^2)\)の項と見かけ上の外部磁場\((h +
2Jm)\)に依存する部分とになるため、さきほどの相互作用なしの場合を参考にすると、
\[\begin{equation}
Z = \exp(-\beta JNm^2)2^N\cosh(\beta [h + 2Jm])^N
\end{equation}\] となる。
したがって、自由エネルギーは、 \[\begin{equation}
F = -kT\ln Z = - \frac{1}{\beta}\ln Z
\end{equation}\] より、 \[\begin{equation}
F = -\frac{1}{\beta}[-\beta JNm^2+N\ln 2 + N\ln\cosh(\beta [h + 2Jm])] =
JNm^2 -NkT\ln\cosh(\frac{h+2Jm}{kT}) + \mathrm{const.}
\end{equation}\]
自由エネルギーは内部エネルギーやポテンシャルと同様、最小になる値が平衡解となる。したがって、平衡状態での\(m\)は \[\begin{equation}
\frac{\partial F}{\partial m} = 0
\end{equation}\] を満たす。それを解くと、 \[\begin{equation}
\frac{\partial F}{\partial m} = 2JNm -NkT
\frac{2J}{kT}\tanh(\beta[h+2Jm]) = 0
\end{equation}\] したがって、平衡状態ではスピン変数の平均値は、
\[\begin{equation}
m = \tanh(\beta[h+2Jm])
\end{equation}\]
を満たす。この先は厳密解を求めることが困難である。平衡状態を決める方程式の両辺に\(m\)があらわれていることから、この方程式をself-consistent
equation(自己無撞着方程式)と呼ぶ。あとは以下のように数値解などで求めるしかない。
以下では\(h = 0\)を想定している。まずは\(2\beta J = 1\)のときから。self-consistent equationの解として、\(m=0\)しかないことがわかる。すなわち、相互作用の強さが十分でない場合、スピンはばらばらの方向を向いていることを意味し、イジングモデルを磁性体のモデルとしてみなすときには、磁性が帯びない状態であると解釈する。
つづいて\(2\beta J = 5\)のとき。self-consistent equationの解として、\(m=-1, 0, 1\)を取りうることがわかる。すなわち、相互作用\(J\)を強くすると、互いのスピンが同じ値をとり\(m=1\)もしくは\(m=-1\)となることがわかる。数値解を求める際には\(m\)の初期値を決めるが、\(h=0\)のときは\(m=0\)から始めると\(m=0\)も解となる。ただし、\(m=0\)以外から始めると\(m\neq0\)で収束するため、\(m=0\)は不安定解であり、\(m = -1, 1\)が安定解となる。
このようにパラメータを変化させることで挙動が変わることを相転移と呼び、相転移の(比較的)簡単なモデルとしてイジングモデルもしくは磁性体のモデルがある。実験的な相転移は固体、液体、気体などの相が変わることを意味し、これまた例えば温度を変化させることで生じる比熱の不連続な変化により相転移を判定する。
上記では\(J_{ij}\)を想定したが、例えば2次元格子上に文字の\(A\)を描くようなパターンを埋め込み(黒が1、白が-1になるように)、この正解パターン\(\sigma_i = \xi_i\)を\(J_{ij} \propto
\xi_i\xi_j\)とすると、連想記憶モデルとなる。また、\(J_{ij} = J + \zeta_{ij}\)として、\(\zeta_{ij}\)をガウス分布から生成したランダムな変数とするモデルをスピングラスと呼ぶ。さらに、\(\sigma_i\sigma_j\)の部分を\(\cos(\theta_i-\theta_j)\)とすればXYスピンモデルになり、さらには非対称な結合荷重も
考えることが可能である(例えばnon-reciprocal Ising modelは2025年のPRLに,
non-reciprocal XY-spin
modelは2023年のPRLに出版されていたりする)。非対称な結合荷重では一般的に平衡状態が存在せず、動的な変化の解析が必要となり、解析の難しさが1段階上がる。このように、イジングモデルは多岐に渡り、かつ今このときでもまだまだ発展し続けるモデルである。
レポート課題:
上記のイジングモデルにおいて、\(J=0\)かつ\(h\neq0\)の場合を考える。self-consistent equationにおける\(J=0\)の場合を考え、以下のように外部磁場\(h\)が\(h>0\)のときには\(m>0\)、\(h<0\)のときには\(m<0\)となることを示せ。以下の図は\(\beta h = 0.25\)の場合。
以下の図は\(\beta h = -0.25\)の場合。
