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

プログラムの動作

本プログラム equation-bisection.py は、半導体における平衡状態のフェルミ準位 ($E_F$) を二分法 (Bisection Method) を用いて数値的に計算します。

目的: 半導体材料の重要な物理量であるフェルミ準位を、電荷中性条件に基づいて高精度に決定すること。
主な機能:
- 伝導帯の電子濃度、価電子帯の正孔濃度、イオン化されたアクセプター濃度、およびイオン化されたドナー濃度を、フェルミ・ディラック分布関数を用いて計算します。
- これらのキャリア濃度およびドーパント濃度から構成される電荷中性条件の関数 $f(E_F)$ を定義し、二分法によって $f(E_F) = 0$ となる $E_F$ を探索します。
- 探索の過程（各反復でのフェルミ準位候補、各キャリア濃度、電荷不均衡量）を標準出力に詳細に表示し、収束状況をユーザーに伝えます。
解決する課題: 半導体中のフェルミ準位は、材料の電気的特性を決定する上で不可欠なパラメータです。ドーピング条件などが複雑な場合、このフェルミ準位を解析的に求めることは困難であり、本プログラムはそのような場合に信頼性の高い数値解を提供します。

原理

本プログラムは、半導体中の電荷中性条件を満たすフェルミ準位 ($E_F$) を二分法によって探索します。

1. 電荷中性条件 半導体中では、全体として電荷的に中性であるという条件が成り立ちます。これは以下の式で表されます。

\[n + N_A^- = p + N_D^+\]

ここで、$n$ は電子濃度、$N_A^-$ はイオン化アクセプター濃度、$p$ は正孔濃度、$N_D^+$ はイオン化ドナー濃度です。この式は、フェルミ準位 $E_F$ の関数として次のように定義される関数 $f(E_F)$ がゼロとなる $E_F$ を探す問題に帰着されます。

\[f(E_F) = n(E_F) + N_A^-(E_F) - p(E_F) - N_D^+(E_F) = 0\]

2. 各電荷キャリア濃度の式 各濃度はフェルミ準位 $E_F$ の関数として、以下の式で計算されます。

電子濃度 ($n$): 伝導帯の有効状態密度 $N_c$ を用います。 $$n = N_c \exp\left(-\frac{E_c - E_F}{k_B T / e}\right)$$ プログラムでは Ne(EF, T) 関数に相当します。
正孔濃度 ($p$): 価電子帯の有効状態密度 $N_v$ を用います。 $$p = N_v \exp\left(-\frac{E_F - E_v}{k_B T / e}\right)$$ プログラムでは Nh(EF, T) 関数に相当します。
イオン化アクセプター濃度 ($N_A^-$): 全アクセプター濃度 $N_A$ と、アクセプター準位 $E_A$ に電子が占有されている確率を掛け合わせたものです。 $$N_A^- = N_A \cdot f_e(E_A, E_F, T)$$ ここで、$f_e(E, E_F, T)$ はフェルミ・ディラック分布関数であり、次式で与えられます。 $$f_e(E, E_F, T) = \frac{1}{\exp\left(\frac{E - E_F}{k_B T / e}\right) + 1}$$ プログラムでは NAm(EF, T) 関数に相当します。
イオン化ドナー濃度 ($N_D^+$): 全ドナー濃度 $N_D$ と、ドナー準位 $E_D$ に電子が占有されていない確率を掛け合わせたものです。 $$N_D^+ = N_D \cdot (1 - f_e(E_D, E_F, T))$$ プログラムでは NDp(EF, T) 関数に相当します。

上記の式において、$e$ は電気素量、$k_B$ はボルツマン定数、$T$ は絶対温度です。プログラム内の ekBT 変数は $e / (k_B T)$ と定義されており、これはエネルギーを電子ボルト (eV) 単位で扱う際に、指数関数の引数を無次元化するために使用されます。すなわち、$E$ がeV単位の場合、$\frac{E}{k_B T / e}$ は $\frac{E \cdot e}{k_B T}$ と等しくなります。$E_c$ は伝導帯の底のエネルギー、$E_v$ は価電子帯の頂上のエネルギーです。

3. 二分法 (Bisection Method) 二分法は、連続関数 $f(x)$ の根を数値的に求めるための反復法です。

根が区間 $[a, b]$ 内に存在することを確認します。具体的には、$f(a)$ と $f(b)$ の符号が異なること（つまり $f(a) \cdot f(b) < 0$）をチェックします。
区間の中点 $c = (a+b)/2$ を計算し、$f(c)$ の値を評価します。
もし $f(a) \cdot f(c) < 0$ ならば、根は区間 $[a, c]$ にあると判断し、$b = c$ と更新します。
そうでなければ、根は区間 $[c, b]$ にあると判断し、$a = c$ と更新します。
区間 $[a, b]$ の幅が許容誤差 eps 以下になるか、または最大反復回数 nmaxiter に達するまで、ステップ2から4を繰り返します。

本プログラムでは、探索変数 $x$ がフェルミ準位 $E_F$ に、初期区間の端点 $a, b$ がそれぞれ EFmin, EFmax に、そして中点 $c$ が EFhalf に対応します。

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

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

numpy: 数値計算、特に指数関数 exp のために使用されます。

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

pip install numpy

必要な入力ファイル

本プログラムは、実行時に外部からの入力ファイルを必要としません。すべての物理定数、半導体パラメータ（例: バンドギャップエネルギー、ドーパント濃度）、および二分法アルゴリズムのパラメータ（例: 初期探索範囲、許容誤差）は、プログラムのソースコード内に直接記述されています。

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

本プログラムは、いかなる出力ファイルも生成しません。計算の進捗状況、各反復におけるフェルミ準位の候補値、各キャリア濃度、および最終的な電荷不均衡量は、すべて標準出力（コンソール）にリアルタイムで表示されます。

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

本プログラムはコマンドライン引数を受け付けません。プログラムを実行するには、Pythonインタープリタを使用してスクリプトファイルを直接実行します。

python equation-bisection.py

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

以下のコマンドは、equation-bisection.py を実行する例です。

python equation-bisection.py

実行結果の例:

Solution of EF by bisection method

  EFmin =     0.00000000  dQmin =  1.0000e+15
  EFmax =     1.00000000  dQmax = -1.0000e+16
  Iter 0: EFhalf =     0.50000000  dQhalf = -1.0000e+16
     Ne=1.9216e+06  Nh=1.0022e+12  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 1: EFhalf =     0.25000000  dQhalf = -1.0000e+16
     Ne=6.5562e-01  Nh=2.2750e+15  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  Iter 2: EFhalf =     0.12500000  dQhalf = -1.0000e+16
     Ne=1.7618e-09  Nh=1.2014e+17  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-1.0000e+16
  ... (中略) ...
  Iter 22: EFhalf =     0.14925766  dQhalf = -2.7161e+12
     Ne=1.2384e+07  Nh=9.5532e+15  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-2.7161e+12
  Iter 23: EFhalf =     0.14926524  dQhalf = -2.7161e+12
     Ne=1.2464e+07  Nh=9.5606e+15  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-2.7161e+12
  Iter 24: EFhalf =     0.14927281  dQhalf = -2.7161e+12
     Ne=1.2543e+07  Nh=9.5679e+15  NA-=1.0000e+15  ND+=1.0000e+16  dQ=-2.7161e+12
  Success: Convergence reached at EF = 0.1492728066253662

実行結果の説明:

プログラムは「Solution of EF by bisection method」というタイトルメッセージで開始されます。
二分法の初期探索区間 [EFmin, EFmax] と、それぞれの端点における電荷不均衡量 dQmin および dQmax が表示されます。dQmin * dQmax < 0 という条件が満たされていることを確認し、二分法が適用可能であることが示されます。
各反復 (Iter) ごとに、現在の区間の中点であるフェルミ準位の候補 EFhalf と、その点における電荷不均衡量 dQhalf が表示されます。
さらに、その EFhalf における詳細なキャリア濃度（電子濃度 Ne、正孔濃度 Nh、イオン化アクセプター濃度 NA-、イオン化ドナー濃度 ND+）と、それらから算出される最終的な電荷不均衡量 dQ が表示されます。
区間の幅がプログラムで設定された許容誤差 eps 以下になった場合、プログラムは収束したと判断し、「Success: Convergence reached at EF = [収束したフェルミ準位]」というメッセージを出力して終了します。
もし、nmaxiter で指定された最大反復回数内に収束条件を満たせなかった場合は、「Failed: Convergence did not reach」というメッセージが表示されて終了します。上記の例では、24回の反復で収束に成功し、約 $0.1492728$ eV のフェルミ準位が計算されたことが示されています。