Weighted Histogram Analysis Method

1.Multiple Histogram Reweighting

カノニカルアンサンブルにおいてある座標構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle Q(\mathbf{r}^N)} に関してバイアスポテンシャル$W_{i}(Q)$ を加えることを 考える。ここでi はトラジェクトリーの番号であり，それぞれのトラジェクトリーで違うバイアスを掛けて いる。この時，それぞれのトラジェクトリーから直接計算されるQ についてのバイアス付きの分布関数は，

$p_{i}(Q)={\frac {1}{Z_{i}}}\int d\mathbf {r} ^{N}exp\left[-\beta \left\{U(\mathbf {r} ^{N})+W_{i}(Q(\mathbf {r} ^{N}))\right\}\right]\delta (Q-Q(\mathbf {r} ^{N}))$

(1)

$={\frac {1}{Z_{i}}}exp\left[-\beta W_{i}(Q)\right]\int d\mathbf {r} ^{N}exp[-\beta U(\mathbf {r} ^{N})]\delta (Q-Q(\mathbf {r} ^{N}))$

(2)

$Z_{i}=\int d\mathbf {r} ^{N}exp[-\beta {U(\mathbf {r} ^{N})+W_{i}(Q(\mathbf {r} ^{N}))}]$

(3)

と表わせる。一方，バイアスなしの場合の分布関数は,

$p_{0}(Q)={\frac {1}{Z_{0}}}\int d\mathbf {r} ^{N}exp[-\beta U(\mathbf {r} ^{N})]\delta (Q-Q(\mathbf {r} ^{N}))$

(4)

$Z_{0}=\int d\mathbf {r} ^{N}exp[-\beta U(\mathbf {r} ^{N})]$

(5)

であるから，バイアス付きの分布関数とは

${\begin{aligned}p_{0}(Q)&={\frac {Z_{i}}{Z_{0}}}exp[\beta W_{i}(Q)]{\frac {1}{Z_{i}}}exp[-\beta W_{i}(Q)]\int d\mathbf {r} ^{N}exp[-\beta U(\mathbf {r} ^{N})]\delta (Q-Q(\mathbf {r} ^{N}))\\&={\frac {Z_{i}}{Z_{0}}}exp[\beta W_{i}(Q)]p_{i}(Q)\end{aligned}}$

(6)

の関係をもつ。

　式6 は理想的には$p_{i}(Q)$ からでも同じ$p_{0}(Q)$ を計算できることを意味している。勿論，有限のトラジェ クトリーでは$p_{i}(Q)$ がほぼ0 となる領域については正しくサンプリングできないため，実際にはそうはいか ない。そこでそれぞれのトラジェクトリーから実際に得られる分布関数

$p_{i}^{est}(Q)\delta Q={\frac {H_{i}(Q)}{M_{i}}}$

(7)

（$H_{i}(Q)$ はQ に関するヒストグラム，$\delta Q$はヒストグラムの幅，$M_{i}$ は全サンプル数）の線形結合

$p_{0}^{est}(Q)=\sum _{i}w_{i}(Q){\frac {Z_{i}}{Z_{0}}}exp[\beta W_{i}(Q)]p_{i}^{est}(Q)$

(8)

$\sum _{i}w_{i}(Q)=1$

(9)

を取り，最も真の分布p0(Q) に近くなるよう重み関数$w_{i}(Q)$ を決定することにする。

　式8 が最も真の分布に近くなるのは、 $p_{0}^{est}(Q)$の分散が最小になる時である。$p_{0}^{est}(Q)$の分散は

$\sigma ^{2}[p_{0}^{est}(Q)]=<p_{0}^{est}(Q)^{2}>-<p_{0}^{est}(Q)>^{2}$

(10)

$=\sum _{i}\sum _{j}w_{i}(Q){\frac {Z_{i}}{Z_{0}}}exp[\beta W_{i}(Q)]w_{j}(Q){\frac {Z_{j}}{Z_{0}}}exp[\beta W_{j}(Q)]\quad \times \quad \left(<p_{i}^{est}(Q)p_{j}^{est}(Q)>-<p_{i}^{est}(Q)><p_{j}^{est}(Q)>\right)$

(11)

$=\sum _{i}\left[w_{i}(Q){\frac {Z_{i}}{Z_{0}}}exp[\beta W_{i}(Q)]\right]^{2}\left(<p_{i}^{est}(Q)^{2}>-<p_{i}^{est}(Q)>^{2}\right)$

(12)

$=\sum _{i}\left[w_{i}(Q){\frac {Z_{i}}{Z_{0}}}exp[\beta W_{i}(Q)]\right]^{2}{\frac {<H_{i}(Q)^{2}>-<H_{i}(Q)>^{2}}{(M_{i}\delta Q)^{2}}}$

(13)

と表わされる。ヒストグラムを作るために使用したサンプルの間に相関がない場合，あるサンプルが$[Q,Q+\delta Q)$ に入る確率，つまり，$H_{i}(Q)$ としてカウントされる確率は$p_{i}(Q)\delta Q$となるから、$H_{i}(Q)=k$ と なる確率は二項分布

$_{M_{i}}C_{k}[p_{i}(Q)\delta Q]^{k}[1-p_{i}(Q)\delta Q]^{M_{i}-k}$

(14)

に従う。従ってその分散は

$<H_{i}(Q)^{2}>-<H_{i}(Q)>^{2}=M_{i}p_{i}(Q)\delta Q[1-p_{i}(Q)\delta Q]\fallingdotseq M_{i}p_{i}(Q)\delta Q$

(15)

となる。式15 を式13 に代入すると，

$\sigma ^{2}[p_{0}^{est}(Q)]=\sum _{i}\left[w_{i}(Q){\frac {Z_{i}}{Z_{0}}}exp[\beta W_{i}(Q)]\right]^{2}{\frac {p_{i}(Q)}{M_{i}\delta Q}}$

(16)

$=p_{0}(Q)\sum _{i}w_{i}^{2}(Q){\frac {Z_{i}}{Z_{0}}}exp[\beta W_{i}(Q)]{\frac {1}{M_{i}\delta Q}}$

(17)

が得られる。後はLagrange の未定乗数法により式9 の条件のもと，式17 を最小化する$w_{i}(Q)$ を求め ればよい。Lagrange の未定乗数法より

${\frac {\partial }{\partial w_{i}(Q)}}\left[\sigma ^{2}[p_{0}^{est}(Q)]-2\alpha \left\{\sum _{i}w_{i}(Q)-1\right\}\right]$

(18)

$=p_{0}(Q)2w_{i}(Q){\frac {Z_{i}}{Z_{0}}}exp[\beta W_{i}(Q)]{\frac {1}{M_{i}\delta Q}}-2\alpha =0$

(19)

従って，

$w_{i}(Q)=\alpha {\frac {1}{p_{0}(Q)}}{\frac {Z_{0}}{Z_{i}}}M_{i}\delta Qexp[-\beta W_{i}(Q)]$

(20)

式9 からα は規格化定数なので

$w_{i}(Q)={\frac {M_{i}exp[-\beta W_{i}(Q)]/Z_{i}}{\sum _{j}M_{j}exp[-\beta W_{j}(Q)]/Zj}}$

(21)

よって，求めるべき分布関数は式8 から

$p_{0}^{est}(Q)={\frac {1}{\delta Q}}{\frac {\sum _{i}H_{i}(Q)}{\sum _{j}M_{j}exp[-\beta W_{j}(Q)]Z_{0}/Zj}}$

(22)

となる。これで未知変数は${\frac {Z_{i}}{Z_{0}}}$ のみとなった。式3,4 から

${\frac {Z_{i}}{Z_{0}}}={\frac {1}{Z_{0}}}\int d\mathbf {r} ^{N}exp[-\beta \{U(\mathbf {r} ^{N})+W_{i}(Q(\mathbf {r} ^{N}))\}]$

(23)

$={\frac {1}{Z_{0}}}\int d\mathbf {r} ^{N}\int dQexp[-\beta U(\mathbf {r} ^{N})]\delta (Q-Q(\mathbf {r} ^{N}))exp[-\beta W_{i}(Q)]$

(24)

$=\int dQp_{0}(Q)exp[-\beta W_{i}(Q)]$

(25)

となるので，$p_{0}(Q)=p_{0}^{est}(Q)$として式22 と式25 をiterative に解くことで$p_{0}^{est}(Q)$を計算できる。自由エネルギー曲線は，

$F(Q)=-k_{B}Tlnp_{0}^{est}(Q)$

(26)

により得られる。

2. 鏡像電荷によるポテンシャル

　誘電体の境界近傍の電荷は鏡像電荷によるポテンシャルを感じる。本章では，そのポテンシャルの大きさ について議論する。

　誘電体$\varepsilon _{0}$,構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \varepsilon_1} の境界を構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle x = 0} にとり，構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \varepsilon_0} 側構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle x = a} の位置に電荷q が存在する系を考える。鏡像電荷の方法によると，観測者にはそれぞれの誘電体構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \varepsilon_0} ,構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \varepsilon_1} 中で図のような電荷q, q', q" が存在しているように見える。

　そこでまずはq', q" をq を用いて表すことを目指す。

　構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle x > 0 } 　の任意の点構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \mathbf{r}} における電場は左図より，

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \mathbf{E}_0(\mathbf{r}) = \frac{1}{4 \pi \varepsilon_0} \left[ \frac{q}{(\mathbf{r} - a \hat{x})^2} \frac{\mathbf{r} - a\hat{x} }{|\mathbf{r} - a\hat{x} |} +\frac{q'}{(\mathbf{r} + a\hat{x})^2 } \frac{\mathbf{r} + a\hat{x} } {|\mathbf{r} + a\hat{x}| } \right] }

(27)

と表わされる。一方構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle x < 0 } 　における電場は右図より，

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \mathbf{E}_1(\mathbf{r}) = \frac{1}{4 \pi \varepsilon_1} \frac{q''}{(\mathbf{r} - a \hat{x})^2} \frac{\mathbf{r} - a\hat{x} }{|\mathbf{r} - a\hat{x} |} }

(28)

となる。

　電場の接続条件（構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \nabla \times E = 0 } 　より）から，境界面において水平方向の電場は等しくなる（構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle E_{0,y}(0, y, z) = E_{1,y}(0, y, z)} ）ので

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \frac{1}{4 \pi \varepsilon_0} \frac{(q + q')y }{(a^2 + y^2 + z^2)^{\frac{3}{2} } } =\frac{1}{4 \pi \varepsilon_1} \frac{q''y }{(a^2 + y^2 + z^2)^{\frac{3}{2} } } }

(29)

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \frac{(q + q')} {\varepsilon_0} = \frac{q''} {\varepsilon_1} }

(30)

が得られる。また，電束密度の接続条件(構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \nabla \cdot D = 0 } 　より）から，境界面において垂直方向の電束密度は等し くなる（構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle D_{0,x}(0, y, z) =D_{1,x}(0, y, z)} ）ので，

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle \frac{1}{4 \pi} \frac{(-q + q')a }{(a^2 + y^2 + z^2)^{\frac{3}{2} } } =\frac{1}{4 \pi} \frac{-q''a }{(a^2 + y^2 + z^2)^{\frac{3}{2} } } }

(31)

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle -q + q' = -q'' }

(32)

が得られる。式30 ，32 から求める電荷は，

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle q'=\frac{ \varepsilon_0 - \varepsilon_1 } {\varepsilon_0 + \varepsilon_1}q }

(33)

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle q''=\frac{ 2 \varepsilon_1 } {\varepsilon_0 + \varepsilon_1}q }

(34)

となる。

　q', q" が求まったので電荷q にかかるポテンシャルを計算しよう。電荷q はq' から

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle F_x(a) =-\frac{1}{4 \pi \varepsilon_0} \frac{ \varepsilon_0 - \varepsilon_1 } {\varepsilon_0 + \varepsilon_1} \frac{q^2}{4a^2} }

(35)

の力を受ける。従ってq を無限遠から構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle x = a } まで持っていくために必要な仕事（＝ポテンシャルエネル ギー）は，

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle U_{mirror}(a) = \int_{\infty}^a F_x(a')da' =\frac{1}{4 \pi \varepsilon_0} \frac{ \varepsilon_0 - \varepsilon_1 } {\varepsilon_0 + \varepsilon_1} \frac{q^2}{4a^2} }

(36)

となる。鏡像電荷の移動には仕事を伴わないので，電荷q とq' が距離2a だけ離れた電荷間に働くポテン シャルエネルギー

構文解析に失敗 (MathML、ただし動作しない場合はSVGかPNGで代替（最新ブラウザーや補助ツールに推奨）: サーバー「https://wikimedia.org/api/rest_v1/」から無効な応答 ("Math extension cannot connect to Restbase."):): {\displaystyle U(a) =\frac{1}{4 \pi \varepsilon_0} \frac{ \varepsilon_0 - \varepsilon_1 } {\varepsilon_0 + \varepsilon_1} \frac{q^2}{2a^2} }

(37)

の半分になることに注意せよ。