トンネル効果の数値解を計算する方法

トンネル効果では、解析解より数値解のほうが簡単に取り扱えます。

微分方程式
Ψ''+(E-U)Ψ=0

について解くだけです。
ここで、数値解を求める際に、微分を差分に直しておきましょう。

1階微分 は
と表現できるため、それをさらにxで微分したものは


と、表現できます。

そうすると、微分方程式は差分方程式

となり
Ψ(n+1)の式

のような数列の漸化式に書き変えることができます。

ここで大事になってくるのが、初期値(境界値)です。

僕が知らないだけかもしれませんが、波動方程式の数値解では、進行波と後退波を分けて計算することができないようです。

そこで、僕の場合は、あらかじめ進行波しかないと解析解でわかっていた透過波のほうから計算してみました。

いうなれば、


ではなく

これを計算したわけです。どこが違うのかおわかりでしょうか。
Ψの添え字n+1とn-1が逆転しているわけです。

ここで透過波は進行波しかないため、自由粒子のように、波動関数の絶対値の2乗が常に一定になることが条件となります。

この条件を満たす初期値というのがありまして
4月19日の日記を参考にします。
Ψ0を一番端、Ψ1を端から2番目の波動関数として

Ψ1/Ψ0=1±i(E-U)dx

一番端のΨはなんでもいいのですが、(E-U)Δx^2に条件が加わり、2番目のΨがなんでもいいわけではなくなるのです。
そうしないと楕円偏光や直線偏光のようになって、波動関数の絶対値の2乗が安定してくれないのです。


実際の計算を見てみましょう

xの行数が多すぎるので、中略して表示しました。
Uが全部ゼロのように見えるのですが
中間のxに、U0の値を逐次入れてます。

一番下のΨ数値はとりあえず1として
下から2番目は
Ψ1={1±i(E-U)dx}Ψ0

を当てはめています。

波動関数の値が複素数に広がっているため、複素数アドインを用いています。
四則演算も関数として記述されるので、最初は慣れないかと思います

COMPLEX(Ψ0,-(E-U)*dx)　(dxとEに絶対参照を忘れずに)

このように書いています。

COMPLEX(Ψ0,(E-U)*dx)

なぜこれではいけないのでしょうか？
(E-U)の符号が反転しています。絶対値が安定するためには、どちらでもいいはずです。

しかし、解析解のほうをexp(ikx)としているので、それに合わせないといけないわけです。
後ろから逆算しているので、螺旋を逆に巻いているのです。

下から3番目以降(以上)では

を用いています。

複素数アドインの関数で書くと

IMSUB(IMPRODUCT(IMSUB(2,(E-U)*(dx^2)),Ψ(n)),Ψ(n+1))
(Eとdxの絶対参照を忘れずに)


計算していくと、波動関数の絶対値|Ψ|は、このようになります。

ちゃんと透過波が安定しているのを確認できますが、透過波だけが安定しすぎていて
全体の大きさが定まっていませんね。

実は、トンネル効果では、解析解でも規格化ができないとされています。
なので仕方なく、トンネルに入る前の進行波を基準にしたりするのですが
先述した通り、数値解では進行波と後退波を区別して計算することができません。
それで透過波を基準に計算を始めたわけですが
本来ならこの図では、透過波が大きくなったり小さくなったりするはずなのに、
まるで透過波を基準に、入射波や反射波などが大きくなったり小さくなったりしています。
次回は、解析解を使って数値解をフィッティングする方法を書きます。

つづく