wavefunction2D.py 技術ドキュメント

プログラムの動作

wavefunction2D.py は、様々な物理モデルにおける2次元波動関数を計算し、その実部、虚部、および確率密度を可視化するPythonスクリプトです。このプログラムは、以下の3つの主要なモードで動作します。

自由電子 (Plain Wave, pw モード): 周期境界条件下での自由電子（または結晶中のブロッホ状態）の波動関数を計算します。
量子箱 (Quantum Box, qbox モード): 無限に高いポテンシャル障壁を持つ量子箱内の粒子の波動関数を計算します。これは、固体中の電子が有限の領域に閉じ込められている状況をモデル化するのに使えます。
水素様原子 (Hydrogen-like model, H モード): 水素原子（または水素様イオン）の電子軌道に対応する波動関数を計算します。ここでは、実数関数で表現される\(s, p_x, p_y, p_z, d_{xy}\) などの軌道を扱います。

プログラムは、指定されたx-y平面上の格子点における波動関数の値を計算し、その結果をグラフとして描画します。具体的には、以下のプロットを生成します。

x軸に沿った波動関数の実部・虚部とブロッホ因子の実部・虚部
y軸に沿った波動関数の実部・虚部とブロッホ因子の実部・虚部
x-y平面上での波動関数の実部・虚部の2次元カラーマップ（コンタープロット）
x-y平面上での波動関数の実部・虚部の3次元表面プロット

これらの可視化により、ユーザーは量子力学的な粒子の振る舞いを直感的に理解することができます。

原理

wavefunction2D.py は、選択された物理モデルに基づいて、以下の波動関数の計算原理を使用しています。

1. 共通のブロッホ因子

全てのモードにおいて、波動関数には共通のブロッホ因子 \(\phi(\mathbf{r})\) が乗じられます。これは、結晶の周期的なポテンシャル下での電子の振る舞いを模擬するため、または単に空間的な位相の傾きを導入するために使用されます。

\[\phi(\mathbf{r}) = \exp\left(i 2\pi \left( k_x \frac{x}{a_x} + k_y \frac{y}{a_y} + k_z \frac{z}{a_z} \right)\right)\]

ここで、\(i\) は虚数単位、\(k_x, k_y, k_z\) はブロッホ波数ベクトル成分（\(\pi/a\) 単位で指定されるが、コード内では \(2\pi\) が乗じられるため、実質的には \(1/a\) 単位として扱われる）、\(a_x, a_y, a_z\) は格子定数です。

2. 自由電子 (Plain Wave, `pw` モード)

自由電子の波動関数 \(\psi_{pw}(\mathbf{r})\) は、格子定数 \(a_x, a_y, a_z\) を持つ周期的な系における平面波として定義されます。量子数 \(n_1, n_2, n_3\) に対応する波数成分を持ち、ブロッホ因子 \(\phi(\mathbf{r})\) が乗じられます。

\[ \psi_{pw}(\mathbf{r}) = \exp\left(i 2\pi \left( n_1 \frac{x}{a_x} + n_2 \frac{y}{a_y} + n_3 \frac{z}{a_z} \right)\right) \cdot \phi(\mathbf{r}) \]

3. 量子箱 (Quantum Box, `qbox` モード)

無限ポテンシャル障壁を持つ1次元量子箱（長さ \(L\)）内の波動関数は、量子数 \(n\) に対して以下のように表されます。

\(n\) が奇数の場合: \(A \cos\left(\frac{n\pi x'}{L}\right)\)
\(n\) が偶数の場合: \(A \sin\left(\frac{n\pi x'}{L}\right)\)

ここで、\(x'\) は箱の中心を原点とする座標です。wavefunction2D.py では、reduce2unitcell 関数で座標が単位胞内にマッピングされ、kqbox 配列で指定された拡張比率が適用されます。psi_qbox 関数は、この量子箱の波動関数にブロッホ因子 \(\phi(\mathbf{r})\) を乗じ、さらに周期的な環境を考慮して複数のブロッホ波数ベクトルの成分を持つ重ね合わせ（\(k_x \pm 1, k_y \pm 1, k_z \pm 1\) など）を計算しています。

4. 水素様原子 (Hydrogen-like model, `H` モード)

水素様原子の波動関数 \(\Psi_{nlm}(r, \Theta, \Phi)\) は、動径部分 \(R_{nl}(r)\) と角度部分（球面調和関数）\(Y_{lm}(\Theta, \Phi)\) の積で表されます。

\[ \Psi_{nlm}(\mathbf{r}) = R_{nl}(r) Y_{lm}(\Theta, \Phi) \cdot \phi(\mathbf{r}) \]

動径部分 \(R_{nl}(r)\): 主量子数 \(n\) と方位量子数 \(l\) に依存し、ボーア半径 \(a_B\) と核電荷 \(Z\) を用いて計算されます。コードでは、kR 配列が半径を調整するために使用されます。
角度部分 \(Y_{lm}(\Theta, \Phi)\): 方位量子数 \(l\) と磁気量子数 \(m\) に依存し、極座標 \((\Theta, \Phi)\) で表現されます。通常、球面調和関数は複素数値ですが、このプログラムでは実数関数としての軌道（例: \(p_x, p_y, d_{x^2-y^2}, d_{xy}\) など）を表現するために、線形結合が計算されています。
- \(m=0\): \(Y_{l0} \propto \cos\Theta\) や \(\cos^2\Theta - 1\) など。
- \(m>0\): 実部として \(\frac{1}{\sqrt{2}} (Y_{lm} + Y_{l,-m}) \propto \sin^m\Theta \cos(m\Phi)\) のような形式。
- \(m<0\): 虚部として \(\frac{1}{i\sqrt{2}} (Y_{lm} - Y_{l,-m}) \propto \sin^m\Theta \sin(|m|\Phi)\) のような形式。これらの実数軌道は、wfname_H 関数で \(s, p_x, p_y, p_z, d_{z^2}, d_{xz}, d_{yz}, d_{x^2-y^2}, d_{xy}\) といった慣用的な名称に対応付けられています。ここでも、psi_H 関数は周期的な環境を考慮して複数のブロッホ波数成分を持つ重ね合わせを計算しています。

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

このプログラムの実行には、以下の非標準Pythonライブラリが必要です。

numpy: 数値計算を効率的に行うためのライブラリです。
matplotlib: グラフ描画のためのライブラリです。特に、3Dプロット機能 (mpl_toolkits.mplot3d) を使用します。

これらのライブラリは、Pythonのパッケージ管理システムである pip を用いてインストールできます。コマンドプロンプトやターミナルで以下のコマンドを実行してください。

pip install numpy matplotlib

必要な入力ファイル

このプログラムは、外部の入力ファイルを必要としません。全ての必要な設定値（モード、量子数、格子パラメータ、メッシュ数など）は、プログラムのソースコード内にデフォルト値として定義されているか、またはコマンドライン引数として与えられます。

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

wavefunction2D.py は、計算結果をファイルとして保存しません。プログラムの実行により、以下の情報が標準出力（コンソール）とグラフィカルなウィンドウに出力されます。

標準出力（コンソール）:
- プログラムの実行モード、格子パラメータ、プロット範囲、メッシュ数などの設定情報。
- x軸に沿った波動関数の実部、虚部、絶対値の2乗、およびブロッホ因子の実部、虚部の数値データ。
- y軸に沿った波動関数の実部、虚部、絶対値の2乗、およびブロッホ因子の実部、虚部の数値データ。
- x-y平面上での波動関数の実部、虚部、絶対値の2乗の数値データ（単位胞内の点のみ）。
- 使用された量子数または波数ベクトル。
グラフィカルなウィンドウ: matplotlib を使用して、以下のプロットを含む単一のウィンドウが表示されます。
- x軸に沿った波動関数の実部 (\(\psi_x_r\)) と虚部 (\(\psi_x_i\)) の1次元プロット。
- x軸に沿ったブロッホ因子の実部 (\(\phi_r\)) と虚部 (\(\phi_i\)) の1次元プロット。
- y軸に沿った波動関数の実部 (\(\psi_y_r\)) と虚部 (\(\psi_y_i\)) の1次元プロット。
- y軸に沿ったブロッホ因子の実部 (\(\phi_r\)) と虚部 (\(\phi_i\)) の1次元プロット。
- x-y平面上の波動関数の実部 (\(\Psi_r\)) の2次元カラーマップ（コンタープロット）。
- x-y平面上の波動関数の虚部 (\(\Psi_i\)) の2次元カラーマップ（コンタープロット）。
- x-y平面上の波動関数の実部 (\(\Psi_r\)) の3次元表面プロット。
- x-y平面上の波動関数の虚部 (\(\Psi_i\)) の3次元表面プロット。

グラフウィンドウを閉じるには、コンソールに表示される "Press ENTER to exit>>" のプロンプトでEnterキーを押す必要があります。

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

プログラムの基本的な使用方法は、引数としてモードと各種パラメータを指定することです。 wavefunction2D.py の usage() 関数に定義されている形式は以下の通りです。

usage: python wavefunction2D.py pw kx0 ky0 kz0 kx ky kz nmeshx nmeshy nmeshz
   For free electron (plain wave)
usage: python wavefunction2D.py qbox nx ny nz kx ky kz nmeshx nmeshy nmeshz
   For quantum dot separated by infite potential barrier
usage: python wavefunction2D.py H n l m kx ky kz nmeshx nmeshy nmeshz
   For hydrogen like model. Radius of wavefunctions are adjusted for clear visualization

引数の説明:

mode: 以下のいずれかを指定します。
- pw: 自由電子 (Plain Wave) モード
- qbox: 量子箱 (Quantum Box) モード
- H: 水素様原子 (Hydrogen-like) モード
kx0, ky0, kz0 (pwモードのみ): 自由電子の波数ベクトル成分（\(\pi/a\) 単位）。波動関数の形状を決定します。
nx, ny, nz (qboxモードのみ): 量子箱内の波動関数の量子数。
n, l, m (Hモードのみ): 水素様原子の主量子数 \(n\)、方位量子数 \(l\)、磁気量子数 \(m\)。
kx, ky, kz: ブロッホ波数ベクトル成分（\(\pi/a\) 単位）。全体の位相因子を決定します。
nmeshx, nmeshy, nmeshz: x, y, z方向のメッシュの数。プロットの解像度に影響します。

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

1. 自由電子 (Plain Wave) モードの実行例

目的: \(k_x = \pi/a\) (kx0=1.0) の自由電子波動関数を、ブロッホ波数 \(\mathbf{k}=(0,0,0)\) の下で可視化します。

python wavefunction2D.py pw 1.0 0.0 0.0 0.0 0.0 0.0 120 120 1

実行結果の説明: このコマンドは、pw モードで波動関数を計算します。 n1=1.0, n2=0.0, n3=0.0 は、波動関数の主要な波数成分がx方向に \(1.0 \times \pi/a\) であることを意味します。kx=0.0, ky=0.0, kz=0.0 は、ブロッホ因子が定数1であることを示します。 nmeshx=120, nmeshy=120, nmeshz=1 は、x-y平面で計算を行い、z方向のメッシュは1であることを示します。実行すると、x方向に周期的な正弦波/余弦波の形状を持つ波動関数の実部・虚部と、x-y平面でのその分布を示すグラフウィンドウが開きます。ブロッホ因子は定数なので、x方向の1次元プロットでは直線が表示されます。コンソールには計算された波動関数の数値データが出力されます。

2. 量子箱 (Quantum Box) モードの実行例

目的: 1,1,1状態（基底状態）の量子箱波動関数を、ブロッホ波数 \(\mathbf{k}=(0,0,0)\) の下で可視化します。

python wavefunction2D.py qbox 1 1 1 0 0 0 120 120 1

実行結果の説明: このコマンドは、qbox モードで波動関数を計算します。 n1=1, n2=1, n3=1 は、x, y, z方向の量子数がそれぞれ1であることを示し、これは量子箱の基底状態に対応します。kx=0, ky=0, kz=0 は、ブロッホ因子が定数1であることを示します。実行すると、中心付近でピークを持ち、端で0に近づく（cos 関数的な）形状を持つ波動関数の実部・虚部が描画されます。x-y平面のコンタープロットは、単位胞の中央で最大値を持つ分布を示します。コンソールには計算された波動関数の数値データが出力されます。

3. 水素様原子 (Hydrogen-like) モードの実行例

目的: 1s軌道 (\(n=1, l=0, m=0\)) の水素様原子波動関数を、ブロッホ波数 \(\mathbf{k}=(0,0,0)\) の下で可視化します。

python wavefunction2D.py H 1 0 0 0 0 0 120 120 1

実行結果の説明: このコマンドは、H モードで波動関数を計算します。 n=1, l=0, m=0 は、水素様原子の1s軌道に対応します。kx=0, ky=0, kz=0 は、ブロッホ因子が定数1であることを示します。実行すると、原点（0,0）で最大値を取り、中心から離れるほど指数関数的に減少する波動関数の実部が描画されます。1s軌道は球対称であるため、虚部は0になります。x-y平面のコンタープロットは、原点を中心とした円形の分布を示します。コンソールには計算された波動関数の数値データが出力されます。