Fourier Grid Hamiltonian法 #
Fourier Grid Hamiltonian method1 (以降FGH法、または省略してFGH) は、離散フーリエ変換のアイデアを 元に設計された時間依存しないシュレーディンガー方程式を高精度に解く数値計算手法です。 位置空間について対角化された基底の線形結合で解を表現することで、少ない基底関数の数で高精度な結果が得られます。
- 文献1に従ってハミルトニアン行列要素を求める。
- 変分原理の直接法の考えでハミルトニアン行列要素を求める。
1. まとめ #
こちらにFGHを実施するwebツールを作成しました。
Solving the 1D Time-Independent Schrödinger Equation
1.1. 問題 #
次の1次元時間依存しないシュレーディンガー方程式を区間$x=[0,L]$, 周期境界条件$\psi(x)=\psi(x+L)$の下で解くことを考えます。
$$ \begin{align} \hat{H}\psi(x)=E\psi(x) \end{align} $$
$$ \begin{align} \hat{H} = -\frac{\hbar^2}{2m}\frac{d^2}{dx^2}+V(x) \end{align} $$
1.2. 解 #
$N$が偶数の場合に限り、FGH法のハミルトニアン行列要素は次になります。$N$が奇数の場合は本文に記載しています。
-
ハミルトニアン行列要素
$$ \begin{align} H_{n,n’}=\left\{ \begin{aligned} &\frac{\hbar^2}{2m}\frac{(-1)^{n-n’}}{2}\frac{1}{\sin^2(\pi (n-n’)/N)} \left(\frac{2\pi}{L}\right)^2, &&(n\ne n’)\\ &\frac{\hbar^2}{2m}\frac{N^2+2}{12} \left(\frac{2\pi}{L}\right)^2 + V(x_n), && (n= n’) \end{aligned} \right. \end{align} $$
$$ \begin{align} H_{n,n’}=H_{n’,n}^* \end{align} $$
-
波動関数
$$ \begin{align} \psi(x) = \sum_{n=0}^{N-1} c_n \varphi_n(x) \end{align} $$
ここで$c_n$は固有ベクトルです。
-
基底関数
$$ \begin{align} \varphi_n(x)=\frac{1}{\sqrt{NL}}\sum_{l=0}^{N-1} e^{i 2\pi l(\frac{x}{L} - \frac{n}{N})} \end{align} $$
2. 問題設定 #
次のハミルトニアンを持つ1次元時間依存しないシュレーディンガー方程式を解くことを考えます。
$$ \begin{align} \hat{H} &= -\frac{\hbar^2}{2m}\frac{d^2}{dx^2}+V(x) \label{e1}\\ &\to \frac{\hbar^2\hat{k}^2}{2m} + V(\hat{x}) \label{e2} \end{align} $$
ここで、$\hbar$はプランク定数、$m$は量子の質量、$\hat{p}=\hbar \hat{k}$であり、 波数$\hat{k}$を位置空間で表現する場合、$\hat{k}=-i\frac{d}{dx}$となります。
3. FGH1の導出 #
3.1. 位置基底における行列要素の導出(連続) #
ブラケット記法において、位置表示・波数表示を次の通りに決めます。
位置 | 波数 | |
---|---|---|
演算子 | $\hat{x}$ | $\hat{k}$ |
固有状態 | $\hat{x}|x\rangle = x|x\rangle$ | $\hat{k}|k\rangle = k|k\rangle$ |
状態の直交性 | $\langle x | x’\rangle = \delta(x-x’)$ | $\langle k | k’\rangle = 2\pi\delta(k-k’)$ |
波動関数 | $\psi(x)=\langle x | \psi\rangle$ | $\psi(k)=\langle k | \psi\rangle$ |
恒等演算子 | $\displaystyle \int dx |x\rangle\langle x |$ | $\displaystyle \int \frac{dk}{2\pi} |k\rangle\langle k | $ |
また、位置表示⇔波数表示の変換の核は下記の通りに書かれます。
$$ \begin{align} \label{e3} \langle x|k \rangle = e^{ikx} \end{align} $$
これらを利用すると、位置の固有状態に対するハミルトニアンの行列要素は次の通りに計算できます。
$$ \begin{align} \langle x|\hat{H}|x’ \rangle &= \langle x|\frac{\hbar^2\hat{k}^2}{2m} + V(\hat{x})|x’ \rangle \\ &= \frac{\hbar^2}{2m}\int \frac{dk}{2\pi}k^2 e^{ik(x-x’)} + V(x)\delta(x-x’) \label{e4} \end{align} $$
3.2. 離散化 #
コンピュータでは連続な量を扱うことは出来ませんので、離散化して考える必要があります。 続いて位置・波数を離散化し、式\eqref{e4}を表現します。
3.2.1. 位置・波数空間の離散化 #
次のように位置・波数を$N$個の点を用いて離散化します。
位置 | 波数 | |
---|---|---|
範囲 | $x =[0,L]$, $L=N\Delta x$ | $k = [-\frac{K}{2},\frac{K}{2}]$, $K=N\Delta k$ |
分点 | $x_n = n\Delta x$ | $k_l = l\Delta k$ |
インデックス | $n=0,1,\cdots,N-1$ | 式\eqref{e6}参照 |
ここで、
$$ \begin{align}\label{e5} \Delta x\Delta k=\frac{2\pi}{N} \end{align} $$
の関係があります。また、波数は下記のようにインデックスを決めます。
$$ \begin{align} l=\left\{ \begin{aligned} &-\frac{N}{2},-\frac{N}{2}+1,\cdots, \frac{N}{2}-1, &&(\text{if $N$ is even})\\ &-\frac{N-1}{2},-\frac{N-1}{2}+1,\cdots, \frac{N-1}{2}, &&(\text{if $N$ is odd}) \\ \end{aligned} \right. \label{e6} \end{align} $$
離散化したことで、周期境界条件が生まれていることに注意してください。つまり、$x=0$と$x=L$は同じ値になります。
3.2.2. 行列要素 #
では式\eqref{e4}を離散化していきます。 クロネッカーのデルタを利用すると、ディラックのデルタが$\delta(x)\to \frac{\delta_{n,n’}}{\Delta x}$として離散化されることも利用して次のように書けます。
$$ \begin{align} \langle x|\hat{H}|x’ \rangle &= \frac{\hbar^2}{2m}\int \frac{dk}{2\pi}k^2 e^{ik(x-x’)} + V(x)\delta(x-x’) \\ &\to \frac{\hbar^2}{2m}\sum_l \frac{\Delta k}{2\pi} k_l^2 e^{i k_l (x_n-x_{n’})} + V(x_n)\frac{\delta_{n,n’}}{\Delta x} \label{e7a}\\ &=\frac{1}{\Delta x}\left[\frac{\hbar^2}{2m} \sum_l l^2 e^{i2\pi l(n-n’)/N} \cdot \frac{\Delta k^2}{N}+ V(x_n)\delta_{n,n’}\right] \\ &=\frac{1}{\Delta x}\left[\frac{\hbar^2}{2m} \frac{1}{N} T^{(2)}_{n-n’,N} \left(\frac{K}{N}\right)^2 + V(x_n)\delta_{n,n’}\right] \label{e7} \\ &=\frac{1}{\Delta x}\left[\frac{\hbar^2}{2m} \frac{1}{N} T^{(2)}_{n-n’,N} \left(\frac{2\pi}{L}\right)^2 + V(x_n)\delta_{n,n’}\right] \label{e7x} \end{align} $$
ここで、$T^{(k)}_{n,N}, (0\leq n<N)$は下記を表します。この和の計算は解析的に実行できますが、紙面を多く使用するので付録に計算結果を示しています。
$$ \begin{align} T^{(k)}_{n,N}= \left\{ \begin{aligned} S_n^{(k)}&\equiv \sum_{l=-{(N-1)/2}}^{(N-1)/2} l^k e^{i2\pi nl/N}, &&(\text{if $N$ is even})\\ R_n^{(k)}&\equiv \sum_{l=-{N/2}}^{(N/2)-1} l^k e^{i2\pi nl/N}, &&(\text{if $N$ is odd}) \\ \end{aligned} \right. \label{e8} \end{align} $$
離散化した表現におけるハミルトニアン行列要素\eqref{e7}が得られたので、これを対角化すれば固有エネルギーが得られます。
4. 変分法(直接法)による導出 #
さて、ハミルトニアン行列要素\eqref{e7}が得られたわけですが、この行列要素を与える離散的な基底関数系で得た結果だ、と考えることもできます。
つまり、本章の目的は波動関数$\psi(x)$を基底関数$\varphi_n(x)$の線形結合で書き、 その基底関数の行列要素が式\eqref{e7}もしくは比例する形になるような$\varphi(x)$の具体的な形を見つけることです。
具体的に、
$$ \begin{align} \psi(x)=\sum_n c_n \varphi_n(x) \end{align} $$
と展開して、
$$ \begin{align} \int \varphi_n^\ast(x) \hat{H} \varphi_{n’}(x) dx \hspace{0.5em}\propto \hspace{0.5em}\text{Eq. \eqref{e7}} \end{align} $$
になるような$\varphi_n(x)$とは何か、という問題です。
なぜこのような表現を見つけるかというと、元々の離散化後のFGHでは分点と分点の間の値が与えられていないため、なにかと都合が悪い場合があるのです。
実際の計算では離散化した分点上の波動関数値だけが得られれば十分な場合がほとんどですが、そうではない状況が起こる場合に対応するため、基底関数による表現を求めます。
4.1. 基底関数 #
式\eqref{e7a}から類推して、式\eqref{e7}と同じ形を与えるような正規直交基底関数を探してみましょう。
4.1.1. 具体的な表現 #
導出過程は類推だけなので、基底関数$\varphi_n(x)$の結果を与えてしまいます。$\varphi_n(x)$は次の形をしています。
$$ \begin{align} \varphi_n(x)&=\frac{1}{\sqrt{NL}}\sum_l e^{i k_l (x-x_n)} \label{e9a}\\ &=\frac{1}{\sqrt{NL}}\sum_l e^{i 2\pi l(\frac{x}{L} - \frac{n}{N})} \label{e9b}\\ \end{align} $$
ここで、$l, x_n, k_l, L, N$は項3.2.1と同じ定義です。
一点、注意しておきますが$N$が偶数の場合、基底関数は複素関数となります。これは$l=-N/2$の対となる$l=N/2$がいないためです。
FGHを使用する場合$x=x_n$の選点上のみに注目することがほとんどであり、選点上では実部しか存在しないので(式\eqref{e12}より)気にする必要は余り多くないでしょう。
4.1.2. 図示 #
基底関数\eqref{e9b}は、次のような形を持ちます。黒い点はFGHに使用する点$x=x_n$です。
実部 | 虚部 | |
---|---|---|
$N=10$ | ||
$N=11$ |
4.1.3. 特徴 #
この$\varphi_n(x)$は次の特徴を持っています。
-
正規直交性
$$ \begin{align} \label{e10} \int_0^L \varphi^{\ast}_n(x)\varphi_m(x) dx = \delta_{n,m} \end{align} $$
-
選点直交性 $$ \begin{align}\label{e11} \varphi_n(x_j)=\sqrt{\frac{N}{L}}\delta_{n,j} \end{align} $$
波動関数が基底関数$\varphi_n$の線形結合で展開される場合、
-
$x$の定義域内の任意の位置$x$における波動関数の値
$$ \begin{align} \psi(x) = \sum_{n=0}^{N-1} c_n\varphi_n(x) \end{align} $$
-
分点の位置$x=x_j$における波動関数の値(※基底関数の選点直交性の利用)
$$ \begin{align} \psi(x_j) = c_j \sqrt{\frac{N}{L}},\hspace{1em}(j=0,1,\cdots,N-1) \label{e12} \end{align} $$
4.2. ハミルトニアン行列要素 #
具体的にハミルトニアン行列要素を計算すると、次のようになります。
$$ \begin{align} H_{n,n’}&\equiv \int_0^L \varphi_n^\ast(x) \hat{H} \varphi_{n’}(x) dx \\ &=\frac{\hbar^2}{2m}\frac{1}{N}T^{(2)}_{n-n’,N}\left(\frac{K}{N}\right)^2 + V(x_n)\delta_{n,n’} \label{e13} \\ &=\frac{\hbar^2}{2m}\frac{1}{N}T^{(2)}_{n-n’,N}\left(\frac{2\pi}{L}\right)^2 + V(x_n)\delta_{n,n’} \end{align} $$
式\eqref{e7}と比較して、$(1/\Delta x)$倍だけ異なりますが、定数倍はシュレーディンガー方程式に関係ないですので同一の物が得られました。
4.3. 規格化 #
$\varphi_n(x)$の線形結合で表現された波動関数の内積は、次の通りに計算できます。
$$ \begin{align} \int_0^L |\psi(x)|^2 dx &= \sum_n c_n^\ast \sum_{n’}c_{n’}\int dx \varphi_n^\ast(x)\varphi_{n’}(x) \\ &= \sum_n |c_n|^2 \end{align} $$
また十分に滑らかな関数$f(x)$に対する内積は次の式で近似できます。
$$ \begin{align} \int_0^L f(x)|\psi(x)|^2 dx &= \sum_n c_n^\ast \sum_{n’}c_{n’}\int dx \varphi_n^\ast(x)f(x)\varphi_{n’}(x) \\ &\approx \sum_n c_n^\ast \sum_{n’}c_{n’} f(x_n)\delta_{n,n’} \\ &= \sum_n |c_n|^2 f(x_n) \end{align} $$
5. 実際の性能 #
調和振動子に対して、実際にFGHで計算した固有値の計算結果を示します。 FGHでハミルトニアン行列要素を生成し、対角化して求めた量子数$n$に属する固有エネルギーを$E_n^{(FGH)}$と書きます。
下記の図は、横軸に基底の数$N$を取り、縦軸に数値計算結果と解析解の差を示しています。
FGH法 (周期境界条件) | $|E_n^{(FGH)}-E_n|\hspace{0.5em} \text{vs} \hspace{0.5em}N$ |
---|---|
$N=10\sim 200$ |
基底の数が増えるに従って、精度(解析解と数値解との差)が指数的に向上することが分かります。 つまり、FGH法で得られる固有値の次数$N$の依存性は
$$ \begin{align} |E^{(FGH)}_n-E_n|\propto c^{-N} \end{align} $$
となります。ここで実数$c>0$は適当な定数です。
比較のために、有限差分法 (Finite-Difference method, 以降FD法) で同じ問題を同じように計算した結果を示します2。
有限差分法で求めた量子数$n$に属する固有エネルギーを$E_n^{(FD)}$と書くと下記の図の結果を得ます。
有限差分法 (固定端条件) | $|E_n^{(FD)}-E_n| \hspace{0.5em}\text{vs}\hspace{0.5em} N$ |
---|---|
$N=10\sim 200$ (FGHと同スケール) |
|
$N=10\sim 10000$ |
有限差分法では次数を10倍にすると精度が1/100倍になるので、精度は$N$のべき乗で向上します。 つまり、有限差分法で得られる固有値の次数$N$の依存性は
$$ \begin{align} |E^{(FGH)}_n-E_n|\propto N^{-2} \end{align} $$
となります。
FGHは$c^{-N}$, FDは${N}^{-c}$であるため、FGHは有限差分法よりも圧倒的に少ない$N$で同じ精度に到達できることが分かります。
6. 数値計算ツール #
6.1. webツール #
こちらにwebツールを作成しました。ブラウザ上でアクセスして、各種パラメータを入力すればFGH法と有限差分法で固有値問題を解き、エネルギー固有値とその固有値に属する固有関数を得ることができます。
Solving the 1D Time-Independent Schrödinger Equation
皆様のブラウザ上で計算・結果出力をするので、入力データをサーバーに送ったり等はしません。 ただし、JavaScript単体で全てを行っているので、計算時間の問題から N=150~200 位が限界かと思います。
6.2. python #
Pythonで書いたプログラムも置いておきます。
※上記プログラムでは、scipy.linalg.eig
を使用していますが、これを対称行列専用のルーチンscipy.linalg.eigh
に変えることで計算速度が2倍程度向上します。
下記の通り、これを
eigenvalues, eigenvectors = sp.linalg.eig(matH, left=False, right=True)
次のように変更してください。
eigenvalues, eigenvectors = sp.linalg.eigh(matH)
7. 応用 #
元々、FGH法を開発する動機としては高速フーリエ変換で使用するグリッド上の点で、位置空間で対角化された基底を探したい、というものがあったそうです1。
推測するに、多くの場合で高速フーリエ変換で用いる分点の位置は、別の方法で用いる分点の位置と一致していないため、両者を結ぶ線形変換が一度挟まってしまい、それが精度や計算速度のボトルネックになってしまったためだと思われます。 それを解消する手段としてFGHが有効であるとなったのでしょう。
波動関数がフーリエ変換出来て何が嬉しいかというと、位置空間と波数空間の行き来が精度の損失なくできるという点です。フーリエ変換を介することで例えば位置空間の微分演算子$d^n/dx^n$は波数空間では$k^n$を掛けるだけになるので、波動関数の微分を考える必要が無くなります。
そのため、例えば運動エネルギー項に含まれる微分演算子の計算が非常に簡潔になり、時間発展など扱いやすくなります。
8. 付録 #
8.1. $\langle x|A(\hat{k})B(\hat{x})|x’\rangle$の計算 #
$$ \begin{align} \langle x|A(\hat{k})B(\hat{x})|x’\rangle &= \langle x|A(\hat{k})\left(\int \frac{dk}{2\pi} |k\rangle\langle k|\right)|x’\rangle B(x’)\\ &= \int \frac{dk}{2\pi} A(k) \langle x|k\rangle\langle k|x’\rangle B(x’)\\ &= \int \frac{dk}{2\pi} A(k) e^{ik(x-x’)}\cdot B(x’) \\ \end{align} $$
- $A(\hat{k})=\hat{I}$の場合
$$ \begin{align} \langle x|B(\hat{x})|x’\rangle &= \int \frac{dk}{2\pi} e^{ik(x-x’)}\cdot B(x’) \\ &= B(x’)\delta(x-x’) \end{align} $$
- $B(\hat{x})=\hat{I}$の場合
$$ \begin{align} \langle x|A(\hat{k})|x’\rangle &= \int \frac{dk}{2\pi} A(k) e^{ik(x-x’)}\\ \end{align} $$
です。
8.2. 指数関数の和の計算 #
8.2.1. Nが偶数の場合 #
$$ \begin{align} S_n^{(k)} &\equiv \sum_{l=-{N/2}}^{(N/2)-1} l^k e^{i2\pi nl/N} \end{align} $$
$$ \begin{align} S^{(0)}_{n}&= \left\{ \begin{aligned} &0, &&(n\ne 0)\\ &N, &&(n=0) \end{aligned} \right. \\ S^{(1)}_{n}&= \left\{ \begin{aligned} &(-1)^n \frac{N}{2}\frac{e^{-i\pi n/N}}{i\sin(\pi n/N)}, &&(n\ne 0)\\ &-\frac{N}{2}, &&(n=0) \end{aligned} \right. \\ S^{(2)}_{n}&= \left\{ \begin{aligned} &(-1)^n \frac{N}{2}\frac{1}{\sin^2(\pi n/N)}, &&(n\ne 0)\\ &\frac{N}{12}(N^2+2), &&(n=0) \end{aligned} \right. \end{align} $$
8.2.2. Nが奇数の場合 #
$$ \begin{align} R_n^{(k)} &\equiv \sum_{l=-{(N-1)/2}}^{(N-1)/2} l^k e^{i2\pi nl/N} \end{align} $$
$$ \begin{align} R^{(0)}_{n}&= \left\{ \begin{aligned} &0, &&(n\ne 0)\\ &N, &&(n=0) \end{aligned} \right. \\ R^{(1)}_{n}&= \left\{ \begin{aligned} &(-1)^n \frac{N}{2}\frac{1}{i\sin(\pi n/N)}, &&(n\ne 0)\\ &0, &&(n=0) \end{aligned} \right. \\ R^{(2)}_{n}&= \left\{ \begin{aligned} &(-1)^n \frac{N}{2}\frac{\cos(\pi n/N)}{\sin^2(\pi n/N)}, &&(n\ne 0)\\ &\frac{N}{12}(N^2-1), &&(n=0) \end{aligned} \right. \end{align} $$
8.2.3. 導出 #
まず下記の関数$S^{(k)}(\alpha)$を計算した後、$\alpha=2\pi n/N$として書き変えることで$S_n^{(k)}$を求めます。 まずは$N$が偶数の時を計算し、その後$N$が奇数の場合を計算します。
$$ \begin{align}\label{ap1} S^{(k)}(\alpha) = \sum_{l=-{N/2}}^{(N/2)-1} l^k e^{i\alpha l} \end{align} $$
両辺を微分することで、次の関係式を得るため、これを利用していきます。
$$ \begin{align} \label{ap2} S^{(k+1)}(\alpha) = -i \frac{\partial}{\partial \alpha} S^{(k)}(\alpha) \end{align} $$
8.2.3.1. $N$が偶数の場合 #
-
$k=0$
等比級数の考えを用いて計算します。これは初項$e^{-iN\alpha/2}$, 公比$e^{i\alpha}$の等比級数です。そのため、次のように計算できます。
$$ \begin{align} S^{(0)}(\alpha) &= \left\{ \begin{aligned} &e^{-i\alpha/2}\frac{\sin(\alpha N/2)}{\sin(\alpha/2)},&&(\alpha\ne 2m\pi)\\ &N,&&(\alpha=2m\pi) \end{aligned} \right. \label{ap3} \end{align} $$
ここで、$m$は任意の整数です。
-
$k=1$
関係式\eqref{ap2}を利用して、次の結果を得ます。
$$ \begin{align} S^{(1)}(\alpha) &= \left\{ \begin{aligned} &-\frac{2ie^{i\alpha}\sin(N\alpha/2) + N(1-e^{i\alpha})\cos(N\alpha/2)}{(1-e^{i\alpha})^2},&&(\alpha\ne 2m\pi)\\ &-\frac{N}{2},&&(\alpha=2m\pi) \end{aligned} \right.\label{ap4} \end{align} $$
適宜、必要に応じで下記の関係式を利用して整理してください。
$$ \begin{align} \left\{ \begin{aligned} 1-e^{i2x}&=-2ie^{ix}\sin{x} \\ 1+e^{i2x}&=2e^{ix}\cos{x} \end{aligned} \right. \end{align} $$
-
$k=2$
関係式\eqref{ap2}を利用して、次の結果を得ます。
$$ \small \begin{align} &S^{(2)}(\alpha) \\ &=\left\{ \begin{aligned} &-\frac{4Ne^{i\alpha}(1-e^{i\alpha})\cos(N\alpha/2)+i[-2(N^2-2)e^{i\alpha}+(N^2+4)e^{2i\alpha}+N^2]\sin(N\alpha/2)}{2(1-e^{i\alpha})^3},&&(\alpha\ne 2m\pi)\\ &\frac{N}{12}(N^2+2),&&(\alpha=2m\pi) \end{aligned} \right. \end{align} $$
8.2.3.2. $N$が奇数の場合 #
- $k=0$
$N$が偶数の場合と同じように、等比級数で考えていくと次の結果を得ることが分かります。
$$ \begin{align} R^{(0)}(\alpha)= \left\{ \begin{aligned} &e^{i\alpha/2}S^{(0)}(\alpha),&&(\alpha\ne 2m\pi)\\ &N,&&(\alpha= 2m\pi)\\ \end{aligned} \right. \end{align} $$
- $k=1$
$$ \begin{align} R^{(1)}(\alpha)= \left\{ \begin{aligned} &e^{i\alpha/2}\left[S^{(1)}(\alpha)+\frac{1}{2}S^{(0)}(\alpha)\right],&&(\alpha\ne 2m\pi)\\ &0,&&(\alpha= 2m\pi)\\ \end{aligned} \right. \end{align} $$
- $k=2$
$$ \begin{align} R^{(2)}(\alpha)= \left\{ \begin{aligned} &e^{i\alpha/2}\left[S^{(2)}(\alpha)+S^{(1)}(\alpha)+\frac{1}{4}S^{(0)}(\alpha)\right],&&(\alpha\ne 2m\pi)\\ &\frac{N}{12}(N^2-1),&&(\alpha= 2m\pi)\\ \end{aligned} \right. \end{align} $$
8.3. 各種の値 #
$$ \begin{align} \int_0^L dx \varphi_n^*(x) \frac{d}{dx} \varphi_{n’}(x) =i\frac{1}{N}T_{n-n’,N}^{(1)}\frac{2\pi}{L} \end{align} $$
$$ \begin{align} \int_0^L dx \varphi_n^*(x) \frac{d^2}{dx^2} \varphi_{n’}(x) = -\frac{1}{N}T_{n-n’,N}^{(2)}\left(\frac{2\pi}{L}\right)^2 \end{align} $$
8.3.1. 一般化 #
$$ \begin{align} \hat{H}=V_2(x)\frac{d^2}{dx^2}+V_1(x)\frac{d}{dx}+V_0(x) \end{align} $$
の時、
$$ \begin{align} H_{n,n’}&\equiv \int_0^L \varphi_n^\ast(x) \hat{H} \varphi_{n’}(x) dx \\ &\approx -V_2(x_{n})\frac{1}{N}T^{(2)}_{n-n’,N}\left(\frac{K}{N}\right)^2+iV_1(x_{n})\frac{1}{N}T^{(1)}_{n-n’,N}\left(\frac{K}{N}\right)+V_0(x_n)\delta_{n,n’} \end{align} $$
となります。更に具体的に書けば、次のように書き表すことができます。
- Nが偶数の時
$$ \small \begin{align} H_{n,n’}=\left\{ \begin{aligned} &-V_2(x_{n})\frac{(-1)^{n-n’}}{2}\frac{1}{\sin^2(\pi (n-n’)/N)} \left(\frac{2\pi}{L}\right)^2 + V_1(x_{n})\frac{(-1)^{n-n’}}{2}\frac{e^{-i\pi (n-n’)/N}}{\sin(\pi (n-n’)/N)}\left(\frac{2\pi}{L}\right), &&(n\ne n’) \\ &-V_2(x_{n})\frac{N^2+2}{12} \left(\frac{2\pi}{L}\right)^2 - \frac{i}{2}V_1(x_{n}) \left(\frac{2\pi}{L}\right) + V_0(x_n), && (n= n’) \end{aligned} \right. \end{align} $$
- Nが奇数の時
$$ \small \begin{align} H_{n,n’}=\left\{ \begin{aligned} &-V_2(x_{n})\frac{(-1)^{n-n’}}{2}\frac{\cos(\pi(n-n’)/N)}{\sin^2(\pi (n-n’)/N)} \left(\frac{2\pi}{L}\right)^2 + V_1(x_{n})\frac{(-1)^{n-n’}}{2}\frac{1}{\sin(\pi (n-n’)/N)}\left(\frac{2\pi}{L}\right), &&(n\ne n’) \\ &-V_2(x_{n})\frac{N^2-1}{12} \left(\frac{2\pi}{L}\right)^2 + V_0(x_n), && (n= n’) \end{aligned} \right. \end{align} $$
$H_{n,n’}$は対称ではなくなり、更に一階微分の係数が複素数になるので、ここまで一般化したものはあまり使わないかもしれません。
9. 参考文献 #
-
C. C. Marston and G. G. Balint‐Kurti, “The Fourier grid Hamiltonian method for bound state eigenvalues and eigenfunctions”, J. Chem. Phys., Vol. 91, No.6 (1989) https://doi.org/10.1063/1.456888" ↩︎ ↩︎ ↩︎ ↩︎
-
FGHとFDで周期境界条件が異なりますが、エネルギー固有値とポテンシャルが一致する転回点$x_\text{tp}$は$n=0$で$\pm 1.0$, $n=4$で$\pm 3.0$に対し、$x=\pm 10$の計算区間を取っているので、境界条件はそこまで重要ではありません。 ↩︎