クラウジウス–モソッティ則導出にローレンツ空洞は不要 — Ewald 和による導出
Clausius–Mossotti (CM) 則
\varepsilon = 1 + \frac{4\pi n \alpha}{1 - \dfrac{4\pi}{3} n \alpha}
\quad\Leftrightarrow\quad
\frac{\varepsilon - 1}{\varepsilon + 2} = \frac{4\pi}{3} n \alpha \qquad (1)
は、左の形で書けば意味は明白である: 裸の応答 $4\pi n\alpha$ が、局所場のフィードバックで $1/(1 - \frac{4\pi}{3}n\alpha)$ 倍にエンハンスされる、という式である(分母がゼロに近づくのが強誘電不安定 = polarization catastrophe)。教科書は右の同値な形で書くことが多いが、本稿では左の形を主とする。
さて、この式の教科書導出(Kittel、Ashcroft–Mermin、Jackson)は、分子の位置に仮想の球形空洞を掘り、空洞表面のローレンツ場 $(4\pi/3)P$ を足す、というやり方をする。これは昔から気持ちが悪い。物理として必要なのは「自分の作る場は自分に働かない」($j \neq i$)だけのはずで、球はどこから来たのか。
答え: 球は要らない。$j \neq i$ の双極子格子和を Ewald 法で計算すれば CM 則はそのまま出る。空洞構成は、エバルト法以前の時代の方法である。が、なぜか延々と教科書に載せられ教えられている。結晶の無限周期性と表面、長距離力というところに注意を払う必要があり、とても難解である。僕の考えではもっと簡単に以下の方法で教えるべき。これは「クーロン力という長距離力の場合の平均場近似」というべきものであり、一般の平均場近似と並べて理解しておくと良い。そしてその上で「固体中ではダイポールって定義できるの?」という疑問もでてくる。
本稿では、真面目に Ewald 法を遂行することで CM 則が出ることを見る。シナリオをまずストレートに示す:
- 周期的な電荷密度に対するポアソン方程式の解には、一様電場 1 本分の不定性がある。セル平均場 $\mathbf E_{\rm macro}$ は方程式が決めない独立変数である(§1)。
- 格子点に分極率 $\alpha$ の点双極子を置き、自己無撞着条件 $\mathbf p = \alpha(\mathbf E_{\rm macro} + S \mathbf p)$ を立てる。物理は「自分の場は自分に働かない」($j \neq i$)だけ(§2)。
- 格子和 $S$ を Ewald 法で計算する。これは周期性を仮定した $\mathbf E_{\rm periodic}$ の厳密解であり、$\mathbf E_{\rm macro}$ は級数の外に独立変数として立つ(§3)。
- 立方格子で $S$ をまじめに評価すると $S_{ab} = (4\pi/3v)\delta_{ab}$ となり、式 (8) から教科書の CM 則が出る。数値確認込み(§4)。
- 静電容量 $C = \varepsilon C_0$ を直接計算し、後づけで定義した $\varepsilon$ が「電極に溜まる電荷の倍率」という実測量に対応することを見る(§5)。
- おまけとして一つの問いに答える: 同じクーロン力の RPA なのに、なぜ分極は局所場で増強され、電荷は遮蔽で弱められるのか(§6)。答えは $G=0$ の輸出の帰結として一行で出る。
ただしその前に、マクロ電場についての理解の確認(手順 1、§1)が必要である。難しい部分(条件収束の分解、対称性による導出、fixed-$D$ の代数、電子系への拡張)は付録に分離した。
1. 周期系のポアソン方程式 — 静的電場の不定性
格子和に入る前に、微分方程式の側で成り立つ基本的な事実を先に固定しておく。これが後で出てくる「$G = 0$ = 境界条件」の正体である。
単位胞で中性($\int_{\rm cell}\rho , d^3r = 0$)の周期的な電荷密度 $\rho(\mathbf r)$ を与え、静電場を求める:
\nabla\times\mathbf E = 0, \qquad \nabla\cdot\mathbf E = 4\pi\rho \qquad (2)
この解は一意ではない。$\mathbf E_{\rm periodic}$ をセル平均ゼロ($\langle\mathbf E_{\rm periodic}\rangle_{\rm cell} = 0$)の周期解とすると、有界な解の全体は
\mathbf E(\mathbf r) = \mathbf E_{\rm periodic}(\mathbf r) + \mathbf E_0 \qquad (3)
($\mathbf E_0$ は任意の一様ベクトル)。証明は二行で済む: 二つの解の差 $\Delta\mathbf E$ は div も curl もゼロだから各成分が調和関数で、全空間で有界な調和関数は定数に限る(Liouville)。つまり周期的な $\rho$ は $\langle\mathbf E\rangle$ を決定できない。Maxwell 方程式のバルクの情報だけでは、セル平均場が一様ベクトル 1 本分だけ未定のまま残る。
Fourier で見ると同じことがこう見える。$G^2 \phi_{\mathbf G} = 4\pi\rho_{\mathbf G}$ は $\mathbf G \neq 0$ の成分を全て決めるが、$\mathbf G = 0$ では中性条件 $\rho_0 = 0$ により $0 = 0$ となって何も決まらない。さらに強い言い方をすると、一様場のポテンシャル $-\mathbf E_0 \cdot \mathbf r$ は周期関数ではないので Fourier 級数の表現能力の外にある。「$\phi$ を Fourier 級数で書く」ことは「$\mathbf E_0 = 0$ と宣言する」ことと同値であり、周期境界条件の計算コードが暗黙にやっているのはこの宣言である。
ここで場とポテンシャルの見え方の違いを明確にしておく。場の配位としては $\mathbf E_0$ は周期的である: 定数場は自明に周期的だから $\mathbf E(\mathbf r + \mathbf R) = \mathbf E(\mathbf r)$ が成り立ち、$\mathbf E$ は $G = 0$ 係数 $\langle \mathbf E \rangle_{\rm cell} = \mathbf E_0$ を含む Fourier 級数を持つ正当な周期ベクトル場になる。不定性は「周期関数の $G = 0$ 係数を方程式が決めない」という形で級数の中に現れる。一方ポテンシャルは $\phi(\mathbf r + \mathbf R) = \phi(\mathbf r) - \mathbf E_0 \cdot \mathbf R$ で周期的でない。傾いたスロープと同じで、勾配はどのセルでも同一だが、高さは周期的でない。第一原理計算(周期境界条件のコード)では周期的ポテンシャルの勾配を電場とするので、$\mathbf E_0 \neq 0$ のこの $\mathbf E$ は表現の枠内に存在せず、扱えないことになる。
では $\mathbf E_0$ を物理的に決めているのは何か。無限周期系の中だけで考えていると答えは出ない。$\mathbf E_0$ はバルクの $\rho$ に含まれない情報であり、境界の物理的状況を指定して初めて決まる独立変数である。そこで、試料の表面に金属電極を蒸着し、電圧 $V$ をかけた状況を考える。このとき一様場成分は外部回路が制御する量になる。これを $\mathbf E_{\rm macro}$ と書き、式 (3) の $\mathbf E_0$ を置き換えた形をわざわざ書いておく:
\mathbf E(\mathbf r) = \mathbf E_{\rm periodic}(\mathbf r) + \mathbf E_{\rm macro} \qquad (4)
周期性から $V$ と $\mathbf E_{\rm macro}$ の関係が直ちに出る。式 (4) を電極間で線積分すると、$\mathbf E_{\rm periodic}$ の寄与はセル 1 個分の並進で完全に消える: 任意の始点 $\mathbf r$ と格子ベクトル $\mathbf R$ に対し $\int_{\mathbf r}^{\mathbf r + \mathbf R} \mathbf E_{\rm periodic} \cdot d\mathbf l = \phi_{\rm periodic}(\mathbf r) - \phi_{\rm periodic}(\mathbf r + \mathbf R) = 0$($\phi_{\rm periodic}$ が周期関数だから)。したがって試料の厚さ $\mathbf d$ にわたる線積分で残るのは
V = \mathbf E_{\rm macro} \cdot \mathbf d \qquad (5)
だけで、$\mathbf E_{\rm macro}$ は「試料両端の電位差 ÷ 厚さ」という測定可能量そのものである。$V$ は外部回路が直接制御する量だから、この意味で $\mathbf E_{\rm macro}$ は外部からかけた電場である。試料表面には束縛電荷が生じるが、理想導体電極の自由電荷がこれを過不足なく補償するので(界面のガウス則 $\sigma_f = D/4\pi$; 供給電荷の増加が静電容量 $\varepsilon$ 倍の正体)、内部の場は $V/d$ のままである。
§1 のまとめ:
- 周期化した計算で得られるのは $\mathbf E(\mathbf r) = \mathbf E_{\rm periodic}(\mathbf r)$、すなわち $\mathbf E_{\rm macro} = 0$ の解だけである。
- 強誘電体や誘電体においては $\mathbf E_{\rm macro} \neq 0$ の状況を扱う必要がある。マクロ分極 $\mathbf P$ の定義は、むしろ後づけで出てくることになる量であり、本質的な物理量でない。
2. 設定 — 自己無撞着方程式
先に場のクラスを限定しておく。本稿で扱う電場は静電場、すなわち $\mathbf E = -\nabla\phi$ と勾配で書けるもの($\nabla\times\mathbf E = 0$)に限る。誘導起源の横成分は扱わない。この制限の下で §1 の結果を使うと、周期系に許される場は「$G \neq 0$ で決まる周期部分 + 一様場 $\mathbf E_0$」で尽きる。一様場も勾配場だが、そのポテンシャル $-\mathbf E_0 \cdot \mathbf r$ だけが非周期になる、というのが §1 で見た構造である。以下に現れる場もすべて勾配場で、双極子場テンソル $T = \nabla\nabla(1/r)$ 自体が勾配の形をしている($T\mathbf p = -\nabla\phi_{\rm dip}$、$\phi_{\rm dip} = \mathbf p \cdot \mathbf r/r^3$)。
体積 $v$ の単位胞、格子点 $\mathbf{R}_j$ に分極率 $\alpha$ の点双極子を置く。双極子場テンソルは、超関数として正確に書くと
T_{ab}(\mathbf r) = \partial_a \partial_b \frac{1}{r}
= \frac{3 r_a r_b - r^2 \delta_{ab}}{r^5}
- \frac{4\pi}{3} \delta_{ab} \delta^3(\mathbf r) \qquad (6)
第一項は $\mathbf r \neq 0$ の正則部(原点は主値として正則化)、第二項がデルタ部である。デルタ部の係数はトレース勘定で一意に決まる: 正則部は $r \neq 0$ で traceless、一方 $\mathrm{Tr}, T = \nabla^2 (1/r) = -4\pi \delta^3(\mathbf r)$ なので、トレースはデルタ部が全部担うしかなく、等方性から $-(4\pi/3)\delta_{ab}\delta^3(\mathbf r)$ と定まる。(デルタ関数が気持ち悪い場合は電荷分布を微小だが有限な分布とかんがえて(ガウス関数などが便利)、極限であると考えれば良い。実際計算もできるはず。)
物理的には、双極子 $\mathbf p$ を含む任意の球で電場を体積積分すると
\int_V \mathbf E \, dV = -\frac{4\pi}{3} \mathbf p \qquad (7)
となることを主張する項で、一様分極した媒質の微視的平均場と $\mathbf P$ を結ぶ関係($-4\pi P/3$ 系の恒等式)はここから出る。Fourier 側では、正則部の変換が $-4\pi(\hat k_a \hat k_b - \delta_{ab}/3)$、デルタ部が $-(4\pi/3)\delta_{ab}$、和が付録 A に出る非解析的な $-4\pi \hat k_a \hat k_b$ — つまり $k$ 空間の等方成分 $1/3$ を担っているのがデルタ部である。ちなみに磁気双極子では対応する項が $+(8\pi/3)\delta_{ab}\delta^3(\mathbf r)$(Fermi の接触項)で係数も符号も異なり、超微細相互作用として実測に掛かる。
$j \neq i$ の和には $\mathbf R_j \neq 0$ だからデルタ部は入らない。ただし消えてなくなるのではなく、$G = 0$(巨視場・境界条件、付録 A)の側に現れる。自己無撞着条件は
\mathbf p = \alpha \left( \mathbf E_{\rm macro} + \sum_{j \neq 0} T(\mathbf R_j)\, \mathbf p \right)
\equiv \alpha \left( \mathbf E_{\rm macro} + S\, \mathbf p \right) \qquad (8)
これで物理は終わりである。これを真面目に計算すればCM則が出る。空洞も球も出てこない。僕の考えではCM則はこれで教えるべき。
なお式 (8) は、ただの分子場近似である。$\alpha$ を固有分極率(proper / 既約分極率)と見たときの RPA の形 — 応答の Dyson 方程式 $\chi = \tilde\chi(1 - v\tilde\chi)^{-1}$ を点双極子基底で書いたもの — に他ならない。$j \neq i$ が既約性の定義(自分の誘起場への応答は $v$ 鎖の反復側に押し出す)、格子和 $S$ が $v$ の双極子成分に対応する。鎖の構造自体は $\alpha$ が厳密なら厳密で、近似は (i) $\alpha$ の計算水準と (ii) 単位間の交換相関核・重なりの無視、にのみ入る(付録 C)。
3. Ewald 分割 — 周期系を仮定した E_periodic の厳密な解
\frac{1}{r} = \frac{\mathrm{erfc}(\eta r)}{r} + \frac{\mathrm{erf}(\eta r)}{r} \qquad (9)
と分け、erfc 部は実空間で(短距離)、erf 部は Poisson 和で逆格子側へ移す:
\sum_{j\,(\text{全格子点})} \frac{\mathrm{erf}(\eta|\mathbf r-\mathbf R_j|)}{|\mathbf r-\mathbf R_j|}
= \frac{4\pi}{v} \sum_{\mathbf G} \frac{e^{-G^2/4\eta^2}}{G^2}\, e^{i\mathbf G\cdot\mathbf r} \qquad (10)
erf 部は全格子点の和(自己項込み)になっているので、$j \neq 0$ に直すには自己項の $r \to 0$ 極限を引く:
\lim_{r\to 0} \partial_a\partial_b \frac{\mathrm{erf}(\eta r)}{r}
= -\frac{4\eta^3}{3\sqrt{\pi}}\, \delta_{ab} \qquad (11)
($\mathrm{erf}(\eta r)/r = (2\eta/\sqrt\pi)(1 - \eta^2 r^2/3 + \cdots)$ の二階微分)。これが「自分の場は自分に働かない」の Ewald 表現である — 滑らかな周期解を作ってから、その原点での自己寄与を引く。まとめると($\partial_a\partial_b, e^{i\mathbf G\cdot\mathbf r} = -G_a G_b, e^{i\mathbf G\cdot\mathbf r}$ に注意):
S_{ab} = \sum_{\mathbf R \neq 0} H_{ab}(\eta;\mathbf R)
- \frac{4\pi}{v} \sum_{\mathbf G \neq 0} \frac{G_a G_b}{G^2}\, e^{-G^2/4\eta^2}
+ \frac{4\eta^3}{3\sqrt{\pi}}\, \delta_{ab} \qquad (12)
実空間部は $g(r) \equiv \mathrm{erfc}(\eta r)/r$ として
H_{ab} = g''(r)\, \hat r_a \hat r_b + \frac{g'(r)}{r}\, (\delta_{ab} - \hat r_a \hat r_b) \qquad (13)
g' = -\frac{\mathrm{erfc}(\eta r)}{r^2} - \frac{2\eta}{\sqrt{\pi}}\frac{e^{-\eta^2 r^2}}{r},
\qquad
g'' = \frac{2\, \mathrm{erfc}(\eta r)}{r^3} + \frac{4\eta}{\sqrt{\pi}}\frac{e^{-\eta^2 r^2}}{r^2}
+ \frac{4\eta^3}{\sqrt{\pi}}\, e^{-\eta^2 r^2} \qquad (14)
式 (12) は $\eta$ に依存しない(分割パラメータは物理に入らない)。$G = 0$ は含めない — Poisson 和で $G$ 空間の級数表示に移した時点で、ポテンシャルの周期性($\langle\mathbf E\rangle_{\rm cell} = 0$、§1)を仮定していることに注意。規約としての位置づけは付録 A で述べる。エバルド和は、G=0をゼロにする、という条件のもとできっちりクーロン力1/rの和を取る方法である。
少し横道になるが、分子動力学法(タンパク質)などでもクーロンをきっちり計算することは非常に重要である。一般には原子はクーロン力の長距離性を打ち消すように配置しようとする(正電荷があればその周りに負電荷がよってきて打ち消す、さらには残ったダイポールに対しても、それを打ち消すようなダイポールの形に原子が並ぼうとする。高次のモーメントも同様)ーなので最終的配置には長距離クーロンが残らないことが多いのだが、その打ち消すような原子配置を求めるにはちゃんとクーロンの総和を取ってやる必要がある。
4. 立方格子 — 教科書の CM 則
立方格子で、まじめに $S$ を式 (12) から計算すると
S_{ab} = \frac{4\pi}{3v}\, \delta_{ab} \qquad (15)
となる。数値評価の結果($v = 1$ に規格化、実空間 $6^3$ 殻・逆格子 $8^3$ 殻):
| 格子 | $S$ の対角成分 | 期待値 |
|---|---|---|
| sc ($\eta$ = 1, 2, 3) | 4.1887902 ($\eta$ 非依存、非対角 < $10^{-16}$) | $4\pi/3$ = 4.1887902 |
| fcc | 4.1887902 | $4\pi/3$ |
| bcc | 4.1887902 | $4\pi/3$ |
$\eta$ 非依存性込みで $4\pi/3$ が確認できる。実は式 (15) は数値に頼らず、立方対称性とトレース勘定だけで導ける(付録 B)。検証スクリプトも付録 B。
これを式 (8) に入れて $\mathbf p$ について解くと
\mathbf p = \alpha \left( \mathbf E_{\rm macro} + \frac{4\pi}{3v} \mathbf p \right)
\quad\Rightarrow\quad
\mathbf p = \frac{\alpha}{1 - \frac{4\pi}{3} n \alpha}\, \mathbf E_{\rm macro} \qquad (16)
($n = 1/v$)。この分母の逆数 $1/(1 - \frac{4\pi}{3} n \alpha)$ が RPA のエンハンスメント因子である(式 (8) の $S \mathbf p$ 項 — $v$ 鎖 — の幾何級数の和、§2)。
最後に、マクロな分極を $\mathbf P = \mathbf p/v$ で定義する。そして $\varepsilon = 1 + 4\pi P/E_{\rm macro}$ — これも定義である(§1 のまとめで述べた「後づけ」がここにあたる)。代入して
\varepsilon = 1 + \frac{4\pi n \alpha}{1 - \frac{4\pi}{3} n \alpha}
\quad\Leftrightarrow\quad
\frac{\varepsilon - 1}{\varepsilon + 2} = \frac{4\pi}{3} n \alpha \qquad (17)
冒頭の式 (1) に到達した。導出のどこにも空洞は現れない。使ったのは (i) $j \neq i$、(ii) Ewald 正則化、(iii) $G = 0$ = $E_{\rm macro}$ 固定という境界条件、(iv) 立方格子での格子和の評価 $S = (4\pi/3v), \mathbf 1$、の四つである。
5. 静電容量 — 系の物理的応答
$\mathbf P$ も $\varepsilon$ も後づけの定義だった。実測に掛かる応答は、電圧 $V$ を入力、電極に流れ込む自由電荷 $Q$ を出力とする微分静電容量 $C = dQ/dV$ である。ここで $Q$ は、電極を囲む閉曲面 $S$ 上の微視的な全電場(解いた配置が作る場そのもの)のフラックス $Q = (1/4\pi)\oint_S \mathbf E \cdot d\mathbf A$ として定義する。$S$ は点双極子を切らない限り任意でよい: ガウス則よりフラックスは囲んだ全電荷を数えるが、各点双極子は中性なので、$S$ が誘電体側へどれだけ食い込んでも囲まれる正味電荷は電極上の自由電荷だけであり、値は $S$ の取り方に依らない。本稿の線形模型では $Q \propto V$ なので $dQ/dV = Q/V$ に一致する(強誘電体では短絡 $V = 0$ でも $Q \neq 0$ となるため、微分で定義しておく必要がある)。ただし、ここからは周期系の数学の枠を離れて、物理的な考察に移る — 「カット(終端)」は前節までの枠内で数学的に定義された操作ではなく、半無限の格子 + 電極という別の設定に対する物理的な議論である。
問うべきことは一つ: $\mathbf p$ が並んでいるものをカット(終端)すると、どういう電場ができるか。 電極に垂直な方向に間隔 $a$ で並ぶ双極子の柱を電荷対 $\pm q$($q = p/a$)で書き直すと、内部では隣り合う $+q$ と $-q$ が打ち消し合い、カットした端に $\pm q$ — 単位面積あたり $\pm p/v$ — だけが残る。つまり終端の効果は $p/v$ の不連続が作る表面電荷、それだけである(教科書の $\sigma_b = P \cdot \hat{\mathbf n}$ はその言い換え)。正極側では $-p/v$。金属内 $\mathbf E = 0$ のガウス箱から界面の全表面電荷は $E_{\rm macro}/4\pi$ だから、電極の自由電荷は
\sigma_f = \frac{E_{\rm macro}}{4\pi} + \frac{p}{v} \qquad (18)
したがって
C = \frac{\sigma_f A}{V} = \frac{A}{4\pi d} \left( 1 + 4\pi \frac{p/v}{E_{\rm macro}} \right) = \varepsilon\, C_0
\qquad C_0 = \frac{A}{4\pi d} \qquad (19)
括弧の中身は §4 で定義した $\varepsilon$ そのものである。つまり後づけで定義した $\varepsilon$ の操作的な中身は $C/C_0$ — 「同じ電圧を保つのに電極へ余計に流れ込む電荷の倍率」という実測量である。系の物理的応答は $C$ であり、$\varepsilon$ はその梱包にすぎない。なお表面近傍の数層はバルクと分極がずれるので、厳密には $C(d) = \varepsilon A/(4\pi d) + (d$ に依らない表面項$)$ であり、バルクの $\varepsilon$ は容量の厚さスケーリング $\varepsilon = \lim_{d \to \infty} 4\pi d, C(d)/A$ の係数として定義される。上の計算はその主要項を与えている。
6. なぜ分極は「増強」され、電荷は「遮蔽」されるのか
式 (16) の分母 $1 - \frac{4\pi}{3}n\alpha$ は 1 より小さい: 局所場は応答を増強し、極限では発散する(polarization catastrophe = 強誘電不安定の素朴模型)。ところが電荷(モノポール)に対する同じ RPA 構造 $\varepsilon = 1 - v\tilde\chi$ では、密度応答の安定性 $\tilde\chi \le 0$ により分母は 1 より大きく、応答は遮蔽されて弱まる。同じクーロン力、同じ幾何級数なのに、フィードバックの符号が逆になるのはなぜか。
答え: 双極子相互作用にも遮蔽(負のフィードバック)は入っている。それは $k \to 0$ の縦成分 $-4\pi\hat k\hat k$(付録 A)であり、本稿の定式化はまさにその部分を核から外して独立変数 $\mathbf E_{\rm macro}$ に移した。増強とは、負の部分を端子(電極)に輸出した後に残る $+4\pi/3$ の話である。
| 設定 | 結合するチャネル | 核 × 応答 | 帰結 |
|---|---|---|---|
| 電荷 | 縦のみ($\nabla\cdot\mathbf E = 4\pi\rho$) | $v\tilde\chi < 0$ | 遮蔽 |
| 分極、$E_{\rm macro}$ 固定(短絡) | 縦成分は電極が処理 | $+\frac{4\pi}{3}n\alpha > 0$ | 増強、catastrophe |
| 分極、$D$ 固定(開回路) | 縦成分が内部に残る | $-\frac{8\pi}{3}n\alpha < 0$ | 遮蔽(付録 D) |
| スピン(Stoner) | 長距離部が最初から結合しない(スピン密度は電荷中性) | $+I N(E_F) > 0$ | 増強、Stoner 条件 |
- 電荷には選択肢がない。 電荷密度は $\nabla\cdot\mathbf E = 4\pi\rho$ を通じて縦場にしか結合できないから、必ず $-4\pi\hat k\hat k$ チャネルに落ち、必ず遮蔽される。
- 分極には二通りの埋め込みがある。 $\mathbf P \parallel \mathbf k$(縦変調)なら $\nabla\cdot\mathbf P \neq 0$ で束縛電荷を作り、電荷と同じ縦チャネルで遮蔽される。$\mathbf P \perp \mathbf k$(横)なら divergence-free で遮蔽すべき電荷がどこにも生じない。残るのは近接場の協力的整列 — これが増強。モノポールは自分の電荷で自分の首を絞めるが、横分極は電荷を出さないので絞める手段がない。
- 増強は端子条件込みの言明である。 短絡電極が像電荷で反電場(輸出された負部分)を食ってくれるから成り立つ。開回路(fixed-$D$)では同じ結晶が遮蔽側に回り、catastrophe の点でも応答は有限に留まる — 強誘電不安定は「短絡条件の現象」である(代数は付録 D。$\varepsilon$ 自体は同じ CM 値になる)。
両チャネルが一つの系で同時に見えるのが LO–TO 分裂である: LO(縦)フォノンは $+4\pi$ の復元力で硬くなり(遮蔽側)、TO(横)は局所場で軟化する(増強側)。強誘電転移 = TO ソフトモード = polarization catastrophe であり、Lyddane–Sachs–Teller 関係 $\varepsilon(0)/\varepsilon_\infty = \omega_{\rm LO}^2/\omega_{\rm TO}^2$ は文字通り二つのチャネルの比になっている。
ストーナー因子との比較. 磁性側の対応物を並べると、この構造はさらにはっきりする。遍歴磁性のスピン帯磁率
\chi = \frac{\chi_{\rm Pauli}}{1 - I N(E_F)}
は式 (1)・(16) と同じ幾何級数で、分母ゼロが Stoner 条件(強磁性不安定)— polarization catastrophe の完全な対応物である(局在スピン側では Weiss 分子場と Curie–Weiss 則)。ただし増強の成立のしかたが違う。スピン密度 $n_\uparrow - n_\downarrow$ は電荷を運ばないので、クーロン核の長距離部 $4\pi/q^2$ はスピンチャネルに最初から結合しない。残るのは短距離の交換 $I > 0$ だけで、増強は輸出の手続きなしに無料で成立する。誘電体は逆で、核の全体が長距離双極子相互作用だから、負の縦成分を電極(端子条件)に輸出して初めて $+4\pi/3$ が残る。磁性にも長距離双極子相互作用は存在する(反磁場 $-4\pi N\cdot\mathbf M$ — まさに「輸出されるべき負部分」の磁性版)が、交換(Fe で $\sim$ eV)に対して桁違いに小さく($\sim 0.1$ meV/原子)、Stoner 条件にはほぼ効かない。それは磁区形成・形状異方性という別の階層に現れる。強磁性がありふれている一方、真の強誘電は反電場との戦い(電極、ドメイン、遮蔽)を要する — この非対称性の起源は、増強因子が「無料」か「輸出後」かの違いにある。
7. 式 (8) の一般形 — p が定義できない場合
そもそも一般の系では $\mathbf p$ が定義できるわけではない。応答する自由度は、電子密度 $\delta\rho(\mathbf r)$ と原子位置 ${\mathbf u_\kappa}$ である。$\mathbf E_{\rm macro}$ を固定した設定で式 (8) に代わる一般形は、「$\mathbf E_{\rm macro}$ が各自由度に及ぼす力」を種とし、「自由度同士が及ぼし合う場を $G = 0$ を除いたクーロン核で足し上げる」自己無撞着方程式になる。自由度をまとめて $x = (\delta\rho, \mathbf u)$ と書けば
x = \tilde\chi \left( f[\mathbf E_{\rm macro}] + K' x \right) \qquad (20)
- $\tilde\chi$: 各自由度の固有応答。電子側は既約分極関数 $\tilde\chi(\mathbf r, \mathbf r')$($\delta\rho = \tilde\chi, \delta V$)、イオン側は力定数行列の逆(変位応答)。式 (8) の $\alpha$ にあたる。
- $K'$: $G = 0$ 成分を除いた相互作用核。電子–電子は周期クーロン核 $v'$($G \neq 0$)、イオン–イオン・イオン–電子も同じ核の Ewald 表示。式 (8) の $S$ にあたり、§3 の Ewald 和はこの核の点双極子版である。
- $f[\mathbf E_{\rm macro}]$: $\mathbf E_{\rm macro}$ が各自由度に及ぼす力。電子には長波長ポテンシャル(下記の有限 $q \to 0$)、イオンには電荷 × $\mathbf E_{\rm macro}$。
解いた後、誘起双極子(電子分 + イオン変位分)を後づけで束ねれば $\varepsilon$ になり、静的な $\varepsilon$ にはイオン変位の寄与(Born 有効電荷・フォノン)が電子分に加算される。点双極子模型は、$\tilde\chi$ をサイト対角の $\alpha$ に、$K'$ を双極子格子和 $S$ に退化させた特別な場合である。電子側の $\tilde\chi$ の中身(既約性とフェルミ統計)については付録 C。
実際の計算をどうするか。§1 の通り、周期計算は $\mathbf E_{\rm macro} = 0$ の解しか作れないので、$f[\mathbf E_{\rm macro}]$ を直接かけることはできない。やり方は実質二つになる。
- 有限 $q$ 変調($q \to 0$ 極限): 外部電荷(ポテンシャル)を波長 $2\pi/q$ の長周期で変調した摂動としてかけ、周期計算の枠内で応答を解いて、$q \to 0$ の極限を取る。一様場はこの極限としてしか周期計算に乗らないので、原理的にはこれでやるしかない。極限が方向依存になる(付録 A の非解析性; §6 の縦/横の違いはこの方向依存の物理的中身である)点に注意。現代分極理論(Berry 位相)は、この $q \to 0$ の極限を波動関数の言葉で整理したものにあたる。
- スラブモデル: 真空層を挟んだスラブに実際に一様電場をかける(鋸歯状ポテンシャル)。ただし電場下で表面の配置を緩和させる — 分極が最大になるような表面配置に落とす — 必要があり、表面の処理が面倒になる。とはいえ、有限 $q$ 変調の方法が正しいとするなら、緩和の仕方の詳細によらないかもしれない。
8. ローレンツ空洞とは + 参考
ローレンツ空洞の手法は Ewald 法(1921)以前のものであり、自己分極を差し引く方法である。初学者には難解であるし、Ewald 法があるのになぜこんなややこしい方法を行う必要があるのかよくわからない、ということになる。それなのに、なぜかよく電磁気の教科書に取り上げられている。
逆に、本記事の議論 — Ewald 法による CM 則の導出 — は標準教科書には載っていない(Kittel、Ashcroft–Mermin、Jackson はいずれもローレンツ空洞で導く)。ただし部品は古くからあり、専門文献では「球は不要」はほぼ確立した認識である:
- Born & Huang, Dynamical Theory of Crystal Lattices (1954), §30 — 双極子格子和を Ewald 法で厳密に扱った古典。ただし文脈がフォノン(長波長光学振動)なので、CM 則の導出として読まれることは少ない。
- A. Redlack & J. Grindlay, J. Phys. Chem. Solids 36, 73 (1975) — 格子和の条件収束と形状依存の系統的解析。§1・付録 A の主張の直接の典拠。
- S. W. de Leeuw, J. W. Perram & E. R. Smith, Proc. R. Soc. A 373, 27 (1980) — 周期境界下の双極子系の条件収束と境界項(tin-foil)の決定版。
- D. E. Aspnes, Am. J. Phys. 50, 704 (1982) — 微視的立場からの局所場の教育的再導出。
- P. B. Allen, J. Chem. Phys. 120, 2951 (2004) — 格子和ベースで CM 則と有限系の分極を再構成した、本記事に最も近い先行。
したがって本記事の新しさは「球なし」自体ではなく、$\mathbf E_{\rm macro}$ を「方程式が決めない自由度の名前」として定義するところから始め(§1)、Ewald 正則化(§3–4)、静電容量による操作的検証(§5)、増強/遮蔽のチャネル分解(§6)までを一本の筋として書いた構成にある。
付録 A. 条件収束の分解 — 正則化と境界条件は別物
条件収束の原因は Fourier 側で見ると明確になる。双極子相互作用の $k \to 0$ 極限は
\tilde T_{ab}(\mathbf k) \to -4\pi\, \hat k_a \hat k_b \qquad (k \to 0) \qquad (21)
と方向依存で非解析的。この $G = 0$ 成分が巨視場・脱分極場・試料形状の情報の置き場所であり、和の「取り方依存性」の全体がここに局在している。したがって処方は:
- $G \neq 0$: 絶対収束するように正則化して足す(Ewald)。格子の内部構造のみ。
- $G = 0$: 和に入れず、境界条件として別に指定する。本稿では「$E_{\rm macro}$ を固定」(周期ポテンシャル解、分子シミュレーション文献で言う tin-foil / 導体境界)を採る。式 (8) の $\mathbf E_{\rm macro}$ がその項である。
空洞構成はこの分離を球面上で暗黙にやっている: 球内離散和が正則化、空洞表面場 $(4\pi/3)P$ が $G = 0$ の側に相当する。分割面が球であることの任意性は、条件収束級数の正則化として辻褄が合っているから最後に消えるだけで、論理は Ewald 側の方が一本道である。
付録 B. 対称性による導出と検証スクリプト
対称性とトレース勘定による式 (15) の導出. $S$ は($\langle\mathbf E\rangle = 0$ の規約の下で)格子構造だけで決まる対称テンソルである。立方格子では立方群で不変だから $S_{ab} = s,\delta_{ab}$、残るのはスカラー $s$ の決定だけで、これはトレースで取れる。式 (12) のトレースを項ごとに:
- 実空間部: $\mathrm{Tr}, \partial_a\partial_b (\mathrm{erfc}(\eta r)/r) = \nabla^2(\mathrm{erfc}(\eta r)/r) = 4\pi\gamma_\eta(r)$($r \neq 0$)。ここで $\gamma_\eta(r) = (\eta/\sqrt\pi)^3 e^{-\eta^2 r^2}$ は規格化ガウス関数($\mathrm{erf}(\eta r)/r$ が単位ガウス電荷のポテンシャルであることの言い換え)。寄与は $4\pi\sum_{\mathbf R\neq 0}\gamma_\eta(R)$。
- 逆格子部: $-\frac{4\pi}{v}\sum_{\mathbf G\neq 0} e^{-G^2/4\eta^2}$。
- 自己項: $3 \cdot \frac{4\eta^3}{3\sqrt\pi} = 4\pi\gamma_\eta(0)$。
Poisson 和 $\sum_{\mathbf R}\gamma_\eta(R) = \frac{1}{v}\sum_{\mathbf G} e^{-G^2/4\eta^2}$ を実空間部に使うと、$\gamma_\eta(0)$ と $G$ 和がすべて相殺して
\mathrm{Tr}\, S = \frac{4\pi}{v}
\qquad\Longrightarrow\qquad
s = \frac{4\pi}{3v}
数表($4\pi/3 = 4.1887902$、sc/fcc/bcc、$\eta$ 非依存)はこの二行の検算である。$\mathrm{Tr}, S = 4\pi/v$ の起源は、式 (6) のトレース勘定($\mathrm{Tr}, T = -4\pi\delta^3$)と $G = 0$ を落としたことの差額であり、要するに $1/3$ の正体は幾何学的な球ではなく、縦射影子のトレース $\mathrm{Tr}(\hat k\hat k) = 1$ の等方三等分である。
contact 恒等式(物理的読み). 独立な整合性チェックとして、セル平均の側から同じ数が出る。$\langle\mathbf E\rangle = 0$ の規約では並進対称性からどのセルも同一で、セル平均の全電場はゼロ: $\langle\mathbf E_{\rm self}\rangle + \langle\mathbf E_{\rm others}\rangle = 0$。自己場のセル平均は、正則部が立方セルで厳密にゼロ($\int_{\rm cell}(3z^2 - r^2)/r^5, d^3r = 0$、$x,y,z$ 軸の等価性から)、デルタ部が式 (6)・(7) より $-(4\pi/3v)\mathbf p$。よって
\langle\mathbf E_{\rm others}\rangle_{\rm cell} = +\frac{4\pi}{3v}\,\mathbf p
これは §4 の $S\mathbf p$(サイト値)と同じ数である: **ローレンツ補正 = 自己場のセル平均の符号反転(contact 項)**という一行の言明が得られる。なお「サイト値 = セル平均」という一致自体は上の二つの独立な計算の帰結であって、調和関数の平均値定理の類推で直接示せるほど自明ではない(調和関数のセル平均は一般に次数 4 以上の成分で中心値からずれる)ことを注意しておく。
検証スクリプト. そのまま実行できる。§4 の表の再現。
import numpy as np
from scipy.special import erfc
def S_tensor(lat, vol, eta, nreal=6, nrec=8):
S = np.zeros((3, 3))
rng = range(-nreal, nreal + 1)
for i in rng:
for j in rng:
for k in rng:
if i == j == k == 0: continue
R = i*lat[0] + j*lat[1] + k*lat[2]
r = np.linalg.norm(R); rh = R / r
e = np.exp(-eta**2 * r**2)
g1 = -erfc(eta*r)/r**2 - (2*eta/np.sqrt(np.pi))*e/r
g2 = (2*erfc(eta*r)/r**3 + (4*eta/np.sqrt(np.pi))*e/r**2
+ (4*eta**3/np.sqrt(np.pi))*e)
S += g2*np.outer(rh, rh) + (g1/r)*(np.eye(3) - np.outer(rh, rh))
B = 2*np.pi*np.linalg.inv(lat).T
for i in range(-nrec, nrec + 1):
for j in range(-nrec, nrec + 1):
for k in range(-nrec, nrec + 1):
if i == j == k == 0: continue
G = i*B[0] + j*B[1] + k*B[2]; G2 = G @ G
S += -(4*np.pi/vol)*np.outer(G, G)/G2*np.exp(-G2/(4*eta**2))
return S + (4*eta**3/(3*np.sqrt(np.pi)))*np.eye(3)
# sc: 4pi/3 = 4.1887902...
print(np.diag(S_tensor(np.eye(3), 1.0, 2.0)), 4*np.pi/3)
付録 C. 電子密度の応答としての再定式化
CM 則の前提は「分極率 $\alpha$ を持つ分離可能な単位(分子)」である。共有結合結晶(Si など)や電子が非局在な系では、格子点に $\alpha$ を割り当てる分割自体が定義できない。このとき何が残り、何が置き換わるか、あらましだけ書いておく。
何が変わるかと言えば、考える対象である。各セルの双極子 $\mathbf p$ ではなく、電子密度に対する応答 — ポテンシャルの変化に対する誘起密度 $\delta\rho(\mathbf r)$ — を考えることになる。その点が違う。その上で構造はそのまま生き残る。式 (8) の三つ組($\alpha$、$S$、$\mathbf E_{\rm macro}$)は、微視的応答の言葉に置き換わる:
- $\alpha$ → 固有(既約)分極関数 $\tilde\chi(\mathbf r, \mathbf r')$。「自分の誘起場への応答を数えない」という $j \neq i$ は、$\tilde\chi$ の既約性($v$ 鎖で繋がる図を含まない)として定義に昇格する。
- 格子和 $S$ → クーロン核 $v$ の $G \neq 0$ 成分。単位ごとの空間分割の代わりに、分割は $G$ 成分で行う: 巨視場($G = 0$)対 微視的局所場($G \neq 0$)。Dyson 方程式 $\chi = \tilde\chi(1 - v\tilde\chi)^{-1}$ の $v$ 鎖が、点双極子模型で $S$ が担っていた局所場の足し上げにあたる(§2)。
- $\mathbf E_{\rm macro}$ 固定 → $v$ の $G = 0$ 成分を落とし、誘電率は逆誘電行列の head から $\varepsilon_M = 1/[\varepsilon^{-1}]_{00}$ で読む(Adler–Wiser)。CM の局所場補正($\varepsilon + 2$ の分母)に相当するものは、$G \neq 0$ の局所場効果として自動的に含まれる。
ただし、既約性だけでは $\tilde\chi$ の定義として足りない。古典の点双極子は可弁別だから $j \neq i$ というラベルの除外で自己場を除けたが、電子は不可弁別なのでラベルで「自分」を除くことはできない。自己相互作用の除去はフェルミ統計を満たすこと — exchange pair diagram を正しくカウントすること — によって初めて実現される(クーロン力を点相互作用にシュリンクしたときに同一ダイヤグラムになるグループ。符号が違うのでキャンセルする。ハバードモデルだとハミルトニアンレベルで落としてしまうことができる)。RPA のバブル鎖は既約ではあるが交換対を落としており、一電子系ですら自己相関が残る。つまり量子系では、局所場の扱い(既約性)とフェルミ統計(交換対の完全性)が独立な要求として並ぶ。
第一原理計算はこの形で $\varepsilon_\infty$ を出している。あるいは §7 の通り、長周期の $\delta\rho$ で有限 $q$ の場をかけ、$q \to 0$ の極限を取る。分子性結晶の極限では $\tilde\chi$ が分子ブロック対角に近づき、本稿の CM 則に戻る。局所的な「分子の分極」をどうしても再導入したければ、局在軌道(Wannier 化)単位の分極として後づけに定義できるが、それは §1 のまとめ 2 で述べた「後づけ」の別の顔であって、理論の入力ではない。
付録 D. 増強/遮蔽の細部
fixed-$D$ の代数. 開回路では $\langle\mathbf D\rangle = \mathbf D$ を固定し、$\mathbf E_{\rm macro} = \mathbf D - 4\pi\mathbf P$ が従属変数になる。式 (8) に代入すると($x \equiv \frac{4\pi}{3}n\alpha$)
\mathbf p = \alpha\left(\mathbf D - 4\pi\mathbf P + \frac{4\pi}{3}\mathbf P\right)
= \alpha\left(\mathbf D - \frac{8\pi}{3}\mathbf P\right)
\quad\Rightarrow\quad
\mathbf P = \frac{n\alpha}{1 + 2x}\,\mathbf D
分母 $1 + 2x > 1$: 遮蔽。同じ格子・同じ $\alpha$ で、フィードバックの符号が端子条件だけで反転する。誘電率は
\varepsilon = \frac{D}{E_{\rm macro}} = \frac{1 + 2x}{1 - x}
で、$(\varepsilon-1)/(\varepsilon+2) = x$ — CM 則 (17) と同一。物質定数は端子条件によらず、増強/遮蔽は「どの変数を固定して測るか」の言明である。catastrophe($x \to 1$)では fixed-$E$ の応答は発散するが fixed-$D$ は有限に留まる: 開回路では反電場が系を安定化する、という強誘電体の標準的事実の模型版。
縦チャネルの一般定理. 静的な縦誘電関数には $\varepsilon(\mathbf q) \ge 1$ または $\varepsilon(\mathbf q) \le 0$ のみが許され、$0 < \varepsilon < 1$(反遮蔽)は安定性から禁止される(Dolgov–Kirzhnits–Maksimov)。「縦チャネルは増強に連続的には行けない」— 増強が許されるのは $G=0$ を輸出した後の横/fixed-$E$ 側だけ、という §6 の非対称性は定理のレベルでも成り立っている。
注意: 双極子が本質的に協力的なわけではない. 単純立方格子の点双極子の真の基底状態は強誘電配列ではなく反強誘電的配列である(Luttinger–Tisza)。§6 の増強は $k = 0$・fixed-$E$ チャネルの応答についての言明であって、双極子格子のエネルギー最小配置の言明ではない。$k \neq 0$ のチャネルまで含めた不安定性解析には格子和を $S(\mathbf k)$ に一般化する必要があり、その $\mathbf k \to 0$ 非解析項がフォノンの LO–TO 分裂の非解析項と同一物になる。