EF-bisection.py の技術ドキュメント

プログラムの動作

EF-bisection.py は、Python言語で記述されたプログラムであり、半導体のフェルミエネルギー (EF) を2分法 (bisection method) を用いて数値的に求めることを目的としています。

このプログラムの主な機能は以下の通りです。

物理定数と半導体パラメータの定義: 電気素量、ボルツマン定数、温度、バンド構造（伝導帯下端エネルギー、価電子帯上端エネルギー）、有効状態密度、ドナー準位・密度、アクセプター準位・密度など、計算に必要な物理量および半導体固有のパラメータをコード内で定義しています。
電荷キャリア密度の計算: フェルミ・ディラック分布関数を基に、伝導帯の電子密度、価電子帯の正孔密度、イオン化したドナー密度、およびイオン化したアクセプター密度を計算する関数を提供します。
フェルミエネルギーの探索: 半導体の電荷中性条件（電子とイオン化アクセプターの総電荷が、正孔とイオン化ドナーの総電荷と等しいという条件）を満たすフェルミエネルギーを、2分法アルゴリズムによって探索します。
計算進行状況と結果の表示: 探索の各ステップにおけるフェルミエネルギーの候補値と、その時点での電荷中性条件の評価値、および各キャリア密度の内訳を標準出力に表示します。計算が収束した場合、最終的なフェルミエネルギーを出力します。

本プログラムは、特定のドーピングレベルや温度における半導体の電子状態、特にフェルミエネルギーの位置を決定するという課題を解決します。これは半導体デバイスの設計や解析において基礎となる重要な情報です。

原理

EF-bisection.py は、半導体のフェルミエネルギー $E_F$ が、熱平衡状態における電荷中性条件によって決定されるという物理原理に基づいています。このプログラムでは、関数 $dQ(E_F)$ がゼロとなる $E_F$ を2分法で探索します。

電荷中性条件: 半導体中の電荷中性条件は、以下の式で表されます。 $$ n + N_A^- = p + N_D^+ $$ ここで、$n$ は電子密度、$N_A^-$ はイオン化アクセプター密度、$p$ は正孔密度、$N_D^+$ はイオン化ドナー密度です。この条件は、以下の関数 $dQ(E_F)$ がゼロとなる $E_F$ を求める問題として定式化されます。 $$ dQ(E_F) = n(E_F) + N_A^-(E_F) - p(E_F) - N_D^+(E_F) = 0 $$
各電荷キャリア密度の計算式: プログラムでは、以下の近似式を用いて各キャリア密度を計算します。
- 電子密度 ($n$): 伝導帯の電子密度は、ボルツマン近似により次式で与えられます。 $$ n = N_c \exp \left( -\frac{E_c - E_F}{k_B T / e} \right) $$ ここで、$N_c$ は伝導帯有効状態密度、$E_c$ は伝導帯下端エネルギー、$k_B$ はボルツマン定数、$T$ は絶対温度、$e$ は電気素量です。
- 正孔密度 ($p$): 価電子帯の正孔密度は、ボルツマン近似により次式で与えられます。 $$ p = N_v \exp \left( -\frac{E_F - E_v}{k_B T / e} \right) $$ ここで、$N_v$ は価電子帯有効状態密度、$E_v$ は価電子帯上端エネルギーです。
- イオン化アクセプター密度 ($N_A^-$): アクセプター準位のイオン化割合はフェルミ・ディラック分布関数 $f(E)$ によって決まります。 $$ N_A^- = N_A f(E_A) = N_A \frac{1}{1 + \exp\left(\frac{E_A - E_F}{k_B T / e}\right)} $$ ここで、$N_A$ はアクセプター密度、$E_A$ はアクセプター準位です。
- イオン化ドナー密度 ($N_D^+$): ドナー準位のイオン化割合は、フェルミ・ディラック分布関数の補数によって決まります。 $$ N_D^+ = N_D (1 - f(E_D)) = N_D \left(1 - \frac{1}{1 + \exp\left(\frac{E_D - E_F}{k_B T / e}\right)}\right) $$ ここで、$N_D$ はドナー密度、$E_D$ はドナー準位です。
2分法 (Bisection Method): 関数 $f(x)=0$ の根を求める数値解析手法の一つです。初期区間 $[a, b]$ において $f(a)$ と $f(b)$ が異なる符号を持つ場合、その区間内に必ず根が存在するという中間値の定理に基づいています。 EF-bisection.py では、$f(E_F)$ を $dQ(E_F)$ として以下の手順でフェルミエネルギーを探索します。
1. 初期区間 $[E_{Fmin}, E_{Fmax}]$ を設定し、$dQ(E_{Fmin})$ と $dQ(E_{Fmax})$ の符号が異なることを確認します。異なる符号でなければ、根が区間に存在しないか、複数存在するため、エラーとして終了します。
2. 区間の中点 $E_{Fhalf} = (E_{Fmin} + E_{Fmax}) / 2$ を計算します。
3. $dQ(E_{Fhalf})$ を評価します。
4. $dQ(E_{Fmin}) \cdot dQ(E_{Fhalf}) < 0$ であれば、根は $[E_{Fmin}, E_{Fhalf}]$ に存在するため、$E_{Fmax} = E_{Fhalf}$ とします。同時に $dQ_{max}$ も $dQ_{half}$ に更新します。
5. そうでなければ、根は $[E_{Fhalf}, E_{Fmax}]$ に存在するため、$E_{Fmin} = E_{Fhalf}$ とします。同時に $dQ_{min}$ も $dQ_{half}$ に更新します。
6. 区間幅 $|E_{Fmax} - E_{Fmin}|$ が許容誤差 eps 以下になるか、最大繰り返し回数 nmaxiter に達するまで上記2-5を繰り返します。収束条件を満たした時点で、現在の $E_{Fhalf}$ を根として出力します。

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

このプログラムは以下の非標準ライブラリを使用します。

numpy: 高度な数値計算、特に数学関数（exp など）の効率的な利用のために必要です。

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

pip install numpy

必要な入力ファイル

EF-bisection.py は、プログラムの実行に際して外部の入力ファイルを必要としません。すべての物理定数、半導体パラメータ、および2分法の計算パラメータは、プログラムのソースコード内に直接ハードコードされています。

パラメータを変更したい場合は、テキストエディタでソースコード EF-bisection.py を開き、以下のセクションを直接編集する必要があります。

T: 計算対象の温度（ケルビン）。
Ev, Ec, Nv, Nc: 半導体の価電子帯上端エネルギー、伝導帯下端エネルギー、価電子帯有効状態密度、伝導帯有効状態密度。
EA, NA: アクセプター準位のエネルギーとアクセプター密度。
ED, ND: ドナー準位のエネルギーとドナー密度。
EFmin, EFmax: 2分法の初期探索範囲。デフォルトでは価電子帯上端と伝導帯下端に設定されています。
eps: 収束判定に使用されるフェルミエネルギーの許容誤差。
nmaxiter: 2分法の最大繰り返し回数。

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

EF-bisection.py は、外部の出力ファイルを生成しません。すべての計算結果と途中経過は、プログラムの実行中に標準出力（コンソール）に直接表示されます。

出力される主な情報は以下の通りです。

初期条件の表示: 2分法の初期区間として設定された EFmin, EFmax、およびそれぞれの点における電荷中性条件の評価値 dQmin, dQmax。
イテレーションごとの詳細: 各繰り返し (イテレーション) において、以下が表示されます。
- イテレーション回数 (Iter)
- その時点でのフェルミエネルギーの候補値 (EFhalf)
- その EFhalf における電荷中性条件の評価値 (dQhalf)
- 電子密度 (Ne), 正孔密度 (Nh), イオン化アクセプター密度 (NA-), イオン化ドナー密度 (ND+) の各値。
収束結果: 計算が指定された許容誤差内で収束した場合、最終的に決定されたフェルミエネルギーの値が表示されます。
失敗メッセージ: 最大繰り返し回数に達しても収束しなかった場合、その旨が通知されます。

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

EF-bisection.py プログラムは、特にコマンドライン引数を必要としません。Pythonインタープリタを使用して直接実行します。

python EF-bisection.py

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

以下のコマンドを実行することで、プログラムのソースコード内にハードコードされた物理定数と半導体パラメータに基づいたフェルミエネルギーの計算が開始され、その結果が標準出力に表示されます。

python EF-bisection.py

実行結果の例:

（以下の出力は、プログラムにハードコードされたデフォルトパラメータで実行した場合の典型的な出力例です。実際の値は、パラメータ設定によって変動します。）

Solution of EF by bisection method

  EFmin =   0.00000000  dQmin =  1.200e+19
  EFmax =   1.00000000  dQmax = -1.000e+16
  Iter 0: EFhalf =   0.50000000  dQhalf = -1.000e+16
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 1: EFhalf =   0.25000000  dQhalf = -1.000e+16
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 2: EFhalf =   0.12500000  dQhalf = -1.000e+16
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 3: EFhalf =   0.06250000  dQhalf = -1.000e+16
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 4: EFhalf =   0.03125000  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 5: EFhalf =   0.04687500  dQhalf = -1.000e+16
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 6: EFhalf =   0.03906250  dQhalf = -1.000e+16
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 7: EFhalf =   0.03515625  dQhalf = -1.000e+16
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 8: EFhalf =   0.03320312  dQhalf = -1.000e+16
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 9: EFhalf =   0.03222656  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 10: EFhalf =   0.03271484  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 11: EFhalf =   0.03247070  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 12: EFhalf =   0.03234863  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 13: EFhalf =   0.03228760  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 14: EFhalf =   0.03225708  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 15: EFhalf =   0.03224182  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 16: EFhalf =   0.03223419  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 17: EFhalf =   0.03223038  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 18: EFhalf =   0.03222847  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 19: EFhalf =   0.03222751  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 20: EFhalf =   0.03222703  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 21: EFhalf =   0.03222679  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 22: EFhalf =   0.03222667  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Iter 23: EFhalf =   0.03222661  dQhalf =  1.200e+19
     Ne=1.0000e-09  Nh=1.2000e+19  NA-=1.0000e+15  ND+=1.0000e+16  dQ=1.2000e+19
  Success: Convergence reached at EF = 0.032226607000000005