Hirano: ページの作成:「== 方針 == 　ポテンシャル等 $f \left( \frac{1}{r} \right)$ 　で定義される関数は、 $r \ll 1$ 　において発散を生じる…」

2026-05-26T04:13:36Z

ページの作成:「== 方針 == 　ポテンシャル等 <math> f \left( \frac{1}{r} \right) </math>　で定義される関数は、<math>r \ll 1 </math>　において発散を生じる…」

新規ページ

== 方針 ==
　ポテンシャル等 <math> f \left( \frac{1}{r} \right) </math>　で定義される関数は、<math>r \ll 1 </math>　において発散を生じる。これを防ぐため、<math>f\left( \frac{1}{r} \right) \to f\left( \frac{1}{r^{mod} } \right)</math> へ置き換えることを考える。<math> \xi=\frac{1}{r^{mod} }</math> とし、<math> f \left( \frac{1}{r^{mod} } \right) </math> の<math> r</math> 微分を考えると、

<math>\begin{align}
\frac{df}{dr} &=\frac{df}{d \xi} \frac{d \xi }{dr} \\

\frac{\partial^2 f}{\partial r^2 } &=\frac{d^2 f}{d \xi^2 } \left( \frac{d \xi} {dr} \right)^2+\frac{df}{d \xi} \frac{d^2 \xi}{dr^2} \\

\frac{\partial^3 f}{\partial r^3 } &=\frac{d^3 f}{d \xi^3 } \left( \frac{d \xi}{dr} \right)^3 + 3\frac{d^2 f}{d \xi^2 } \frac{d \xi}{dr} \frac{d^2 \xi}{dr^2 } + \frac{df}{d \xi} \frac{d^3 \xi}{dr^3 } \\
\end{align}</math>

となっている。つまり、例えば、'''LJ'''ポテンシャルなら、

<math>\begin{align}

U_{LJ} (\xi)&=4 \varepsilon (\sigma^{12} \xi^{12} - \sigma^6 \xi^6 ) \\

\mathbf{F}_{LJ} (\mathbf{r})&=-4 \varepsilon (12 \sigma^{12} \xi^{11} -6 \sigma^6 \xi^5 ) \frac{ d \xi}{dr} \mathbf{e}

\end{align}</math>

相互作用テンソルなら、

<math>\begin{align}

T^{(0)}&=\frac{\xi}{4 \pi \varepsilon_0 } \\

\mathbf{T}^{(1) }&=\frac{1}{4 \pi \varepsilon_0 } \frac{d\xi} {dr} \mathbf{e} \\

\mathbf{T}^{(2) }&=\frac{1}{4 \pi \varepsilon_0 } \left( \frac{d^2 \xi}{dr^2 } \mathbf{ee}+ \frac{d\xi}{dr} \mathbf{e}^{(2) } \right) \\

T_{ijk}^{(3) } &= \frac{1}{4 \pi \varepsilon_0 } \left( \frac{d^3 \xi}{dr^3 } e_i e_j e_k+\frac{\partial^2 \xi}{\partial r^2 } (e_{ij}^{(2) } e_k+e_{ki}^{ (2) } e_j+e_{jk}^{ (2) } e_i )+\frac{ \partial \xi}{\partial r} e_{ijk}^{ (3) }\right)

\end{align}</math>

である。従って、<math>\xi =\frac{1}{r^{mod} } </math>　及びその微分をあらかじめ計算しておくことで、プログラムの変更を最小限にして<math> \frac{1}{r}</math>　の発散を阻止することができる。また、長距離部分のカットオフも同様の方法で組み込むことができる。

== ξの形状 ==
　微分の計算のしやすさから、

<math>
\xi(r)=\frac{1}{r^{mod} } =
\begin{cases}

\sum_{n=0}^8 a_n r^n & r < r_1 \\
\frac{1}{r} & r_1 < r < r_2 \\
\sum_{n=0}^8 b_n (r-r_3 )^n & r_2 < r < r_3 \\
b_0 & r_3 < r

\end{cases}</math>

と置くのが良い。この関数形では、<math> r > r_3</math> 　において ξ が0でないが、エネルギーの基準点を

<math>
U^{mod} (\xi)=U(\xi)-U(b_0 )
</math>

と補正することにより、トラジェクトリー計算には影響がない。（多項式の次数を上げることで<math> b_0=0 </math>　とすることも可能だが、その場合、ξ の高次の微分でスイッチング関数の影響が無視できなくなる。）

　<math>a_n,b_n </math>　は、4回微分が連続という条件を課すことで、

<math>\begin{align}
a_n &=c_n (0,r_1 ) \\
b_n &=c_n (r_3,r_2 ) \\
c_0 (r_1,r_2 )&=-\frac{ -70r_2^4-35r_2^3 r_d-15r_2^2 r_d^2-5r_2 r_d^3-r_d^4}{70r_2^5 } \\
c_1 (r_1,r_2 )&=0 \\
c_2 (r_1,r_2 )&=0 \\
c_3 (r_1,r_2 )&=0 \\
c_4 (r_1,r_2 )&=0 \\
c_5 (r_1,r_2 )&=-\frac{ 35r_2^3+30r_2^2 r_d+15r_2 r_d^2+4r_d^3}{5r_2^5 r_d^4 } \\
c_6 (r_1,r_2 )&=-\frac{ -14r_2^3-13r_2^2 r_d-7r_2 r_d^2-2r_d^3}{r_2^5 r_d^5 } \\
c_7 (r_1,r_2 )&=-\frac{70r_2^3+68r_2^2 r_d+39r_2 r_d^2+12r_d^3}{7r_2^5 r_d^6 } \\
c_8 (r_1,r_2 )&=-\frac{-5r_2^3-5r_2^2 r_d-3r_2 r_d^2-r_d^3}{2r_2^5 r_d^7 } \\

\end{align}</math>

となる。なお、<math>r_d=r_2-r_1</math> である。（導出は[[ファイル:switch.nb]] 参照。関数の形状や微分形もこちらを参照。 3階微分はやや波うちが大きいが、<math>r_d</math>　を大きく取ることで改善される。）

== switch.nb ==
上の <math>1/r</math> および <math>1/r^{mod}</math> の形をプロットする mathematica スクリプトである。スウィッチの境界は <math>r_c - r_s = 5.5</math> A, <math>r_c = 11.5</math> A としている。実行したプロット結果は、Kikkawa et al. ''J. Am. Chem. Soc.'' 2015, 137, 8022-8025 の Supporting Information 中の Figure S3 にも記載されている。

1/r の発散の防止について - 版の履歴

Hirano: ページの作成:「== 方針 == ポテンシャル等 f \left( \frac{1}{r} \right) で定義される関数は、r \ll 1 において発散を生じる…」

Hirano: ページの作成:「== 方針 == 　ポテンシャル等 $f \left( \frac{1}{r} \right)$ 　で定義される関数は、 $r \ll 1$ 　において発散を生じる…」