HOME > natural science Laboratory > コンピュータ・シミュレーション講座 > 計算物理学

【量子力学再入門9】
シュレディンガー方程式の数値解法(スペクトル法)

文責:遠藤 理平 (2018年5月11日) カテゴリ:仮想物理実験室(286)計算物理学(143)

量子力学とは原子や電子といったナノスケールにおける粒子の運動を記述する物理学の一分野です。 メートルスケールで成り立つ「ニュートンの運動方程式」に対して、量子力学では「シュレディンガー方程式」と呼ばれる基礎方程式が成り立つと考えられています。これまでは解析的に得られている固有関数の重ねあわせで波束の運動をシミュレーションしました。 しかしながら、シュレディンガー方程式で固有関数が解析的に得られるのは非常に稀な系だけで、多くの系では数値的に計算する必要があります。 そこで本稿ではスペクトル法と呼ばれるシュレディンガー方程式の数値解法を導入します。

スペクトル法による解法

シュレディンガー方程式は波動方程式なので解は一般の有限差分法による数値計算を行うことができます。本稿では有限差分法よりも高精度な計算が可能なスペクトル法というアルゴリズムを用います。スペクトル法というのは、求めたい波動関数を正規直交系で展開した展開係数の時間発展を計算するというアルゴリズムで、最も簡単なのはフーリエ級数展開を用いるものです。 シュレディンガー方程式から展開係数の時間発展を表す常微分方程式を導きます。波動関数とポテンシャル項を

\psi(x,t) =\sum\limits_{n} a_n(t) \, e^{ik_nx}
V(x,t) =\sum\limits_{l} v_l(t) \, e^{ik_lx}

と波数空間で展開して、時間依存性はその係数部分に含めます(-L/2\leq x\leq L/2)。これらを元のシュレディンガー方程式に代入します。

i\hbar\sum\limits_{n} \frac{da_n(t)}{dt} \,e^{ik_nx} = \sum\limits_{n}\frac{\hbar^2k_n^2}{2m}\, a_n(t)e^{ik_nx} +\sum\limits_{n}\sum\limits_{l} v_l(t)a_n(t)\, e^{ik_nx}e^{ik_lx}

この関係式からx依存性を排除するために、両辺に\exp^{-ik_{n'}x}を掛けてxで積分を行います。直交関係から

i\hbar\sum\limits_{n} \frac{da_n(t)}{dt}\, \delta_{nn'} = \sum\limits_{n}\frac{\hbar^2k_n^2}{2m}\, a_n(t)\delta_{nn'} +\sum\limits_{n}\sum\limits_{l} v_l(t)a_n(t)\, \, \delta_{n+l,n'}

となるので、最終的に次の表式が導かれます。

i\hbar\,\frac{da_n(t)}{dt} = \frac{\hbar^2k_n^2}{2m}\, a_n(t) +\sum\limits_{l} v_l(t)a_{n-l}(t)

これはn=0~N-1までのN個の連立常微分方程式となっています。本稿ではこれをルンゲ・クッタ法を用いて解きます。


シミュレーション結果

中心エネルギーE0 = 0 の波束の時間発展と解析解との差である計算誤差を示します。

波動関数の空間分布

20180511.gif

計算誤差

20180511.gif

計算パラメータ

空間スケール 1E-11[m];//横軸の値
時間スケール 1E-14[s];//1コマの時間
数値計算の時間間隔1E-16[s]


考察と次の課題

・ポテンシャルが存在しない系で動作は確認できた。
・次は箱型ポテンシャル障壁に波束を衝突させる系で誤差チェックを行う。

プログラムソース(C++)

http://www.natural-science.or.jp/files/physics/QuantumPhysics_OpenMP.zip
※VisualC++、GCC(MinGW)で動作確認しています。



▲このページのトップNPO法人 natural science トップ

関連記事

仮想物理実験室







計算物理学







▲このページのトップNPO法人 natural science トップ