transfer_matrix.py

プログラムの動作

本プログラムは、伝達行列法（Transfer Matrix Method）を用いて、多層構造（多重量子井戸など）における電子の1次元波動関数および透過確率を計算するPythonスクリプトです。計算対象となるポテンシャルプロファイルは、プログラム内部で生成する多重量子井戸（MQW）モデル、または外部のExcelファイルから読み込んだ任意の構造を指定できます。計算された波動関数の分布やエネルギーに対する透過確率はグラフとして可視化され、各種データはExcelファイルに保存されます。特定のエネルギーに対する波動関数の振る舞いを調べるモードと、広いエネルギー範囲での透過確率のスペクトルを計算するモードの2つを備えています。

原理

本プログラムでは、1定常シュレディンガー方程式を解くために空間を \(N\) 個の微小領域に分割し、各領域における波動関数 \(\Psi\) を前進波と後退波の重ね合わせとして扱います。各領域 \(i\) における電子の波数 \(k_{z, i}\) は、入射エネルギー \(E_z\) 、ポテンシャルエネルギー \(U_i\) 、および相対有効質量 \(m^*_i\) から以下の式で計算されます。

\[k_{z, i} = \frac{\sqrt{2 m^*_i m_e (E_z - U_i) e}}{\hbar}\]

ここで、\(m_e\) は電子の静止質量、 \(e\) は素電荷、 \(\hbar\) はディラック定数です。各領域での波動関数 \(\Psi_i(z)\) は、振幅 \(A_i\) と \(B_i\) を用いて以下のように表されます。

\[\Psi_i(z) = A_i \exp(i k_{z, i} z) + B_i \exp(-i k_{z, i} z)\]

隣接する領域間の境界条件（波動関数とその微分の連続性）から、振幅 \(A_i\) と \(B_i\) を結びつける伝達行列の要素が導出されます。プログラム内では以下の漸化式を用いて隣接領域の振幅を計算します。まず、有効質量と波数の比 \(m_k\) を求めます。

\[m_k = \frac{m^*_i}{m^*_{i-1}} \frac{k_{z, i-1}}{k_{z, i}}\]

これを用いて、係数 \(a_+\) と \(a_-\) を定義します。

\[a_+ = \frac{1 + m_k}{2}, \quad a_- = \frac{1 - m_k}{2}\]

位相因子 \(P\) と \(Q\) を以下のようにおきます。

\[P = \exp(i (k_{z, i-1} - k_{z, i}) z_i)\]

\[Q = \exp(i (k_{z, i-1} + k_{z, i}) z_i)\]

これらにより、隣接する領域の振幅は次のように更新されます。

\[A_i = a_+ P A_{i-1} + a_- Q^{-1} B_{i-1}\]

\[B_i = a_- Q A_{i-1} + a_+ P^{-1} B_{i-1}\]

終端領域に到達した後、系全体の透過確率 \(T\) は、始点と終点の有効質量、波数、および前進波の振幅の比から次のように求められます。

\[T = \left| \frac{m^*_{N-1}}{m^*_0} \frac{k_{z, 0}}{k_{z, N-1}} \frac{A_0}{A_{N-1}} \right|^2\]

必要な非標準ライブラリとインストール方法

本プログラムを実行するためには、数値計算およびデータ処理、グラフ描画のためのライブラリと、Excelファイルを読み書きするためのライブラリが必要です。以下のコマンドを使用して、必要なパッケージをインストールしてください。

pip install numpy pandas matplotlib openpyxl

必要な入力ファイル

引数で外部ポテンシャルファイルを指定する場合、Excel形式（.xlsx）の入力ファイルが必要です。プログラムは指定されたファイルから以下のデータ構造（列の順番）を読み込みます。ヘッダーの有無にかかわらず、列の位置でデータを取得します。

1列目: z座標（単位：オングストローム）
2列目: ポテンシャルエネルギー \(U\) （単位：eV）
3列目: 有効質量 \(m^*\) （電子質量 \(m_e\) に対する比率）

なお、引数に mqw を指定した場合は、プログラム内部で多重量子井戸のポテンシャルデータが自動的に生成されるため、入力ファイルは不要です。

生成される出力ファイル

透過確率計算モード（tr）を実行した場合、カレントディレクトリに以下のExcelファイルが出力されます。

potential.xlsx : 計算に使用された座標、ポテンシャルエネルギー、および有効質量のプロファイルが保存されます。列名は z (A), U (eV), m* (me) です。
transmission.xlsx : 対数スケールで等間隔に設定された各エネルギーに対する透過確率のデータが保存されます。列名は E (eV), T です。

コマンドラインでの使用例 (Usage)

コマンドラインから引数を与えることで、動作モードとパラメータを指定します。引数は順番に評価されます。プログラムには主に波動関数を計算する wf モードと、エネルギー依存の透過確率を計算する tr モードがあります。

python D:/git/sphinx/tkProg/source/Quantum/transfer_matrix.py [mode] [pottype] [nz] [Ez0] [Emin] [Emax] [nE]

各引数の意味は以下の通りです。

mode: 動作モード。wf（波動関数計算）または tr（透過確率計算）を指定します。
pottype: ポテンシャルの種類。mqw または入力用のExcelファイル名（例：potential.xlsx）を指定します。
nz: 空間座標（z軸）の分割数。
Ez0: 計算対象の基準エネルギー（eV）。
Emin: （tr モード時のみ）透過確率計算の最小エネルギー（eV）。
Emax: （tr モード時のみ）透過確率計算の最大エネルギー（eV）。
nE: （tr モード時のみ）透過確率を計算するエネルギーのステップ数。

コマンドラインでの具体的な使用例

単一エネルギーでの波動関数を計算し、ポテンシャルが多重量子井戸（mqw）である場合の使用例です。空間分割数を 201、対象エネルギーを 0.1 eV としています。

python D:/git/sphinx/tkProg/source/Quantum/transfer_matrix.py wf mqw 201 0.1

このコマンドを実行すると、エネルギー0.1 eVにおける波動関数が計算され、ポテンシャル形状、有効質量、波数の実部と虚部、振幅の絶対値、および波動関数の形状（実部、虚部、二乗絶対値）を示すグラフが表示されます。

次に、外部のExcelファイル（例：potential_defect.xlsx）を読み込んで透過確率を計算する場合の使用例です。分割数 201、基準エネルギー 0.1 eV、計算範囲を 0.01 eV から 9.5 eV、エネルギーステップ数を 1001 としています。

python D:/git/sphinx/tkProg/source/Quantum/transfer_matrix.py tr potential_defect.xlsx 201 0.1 0.01 9.5 1001

このコマンドを実行すると、指定したエネルギー範囲での透過確率のスペクトルが計算されます。計算中に進行状況が標準出力に表示され、完了するとポテンシャル情報とエネルギーに対する透過確率がそれぞれのExcelファイルとして保存され、結果がグラフに表示されます。