import csv
import numpy as np
from numpy import sin, cos, tan, pi, exp
import sys
"""
2分法を用いて半導体のフェルミエネルギーを計算するスクリプト。
詳細説明:
キャリア密度、イオン化不純物密度を考慮した電荷中性条件を2分法で解き、
半導体におけるフェルミエネルギーを決定します。
関連リンク:
:doc:`EF-bisection_usage`
"""
# physical parameeters
e = 1.60218e-19 # C";
kB = 1.380658e-23 # JK<sup>-1</sup>";
T = 300.0
ekBT = e / kB / T
# semiconductor parameters
#価電子帯上端エネルギー eV
Ev = 0.0
#伝導帯下端エネルギー
Ec = 1.0
#価電子帯有効状態密度
Nv = 1.2e19
#伝導帯有効状態密度
Nc = 2.1e18
#アクセプター準位
EA = 0.05
#アクセプター密度 cm-3
NA = 1.0e15
#ドナー準位
ED = Ec - 0.05
#ドナー密度
ND = 1.0e16
# parameters
#初期範囲として、価電子帯上端と伝導帯下端エネルギーを設定する
EFmin = Ev
#EFmax = Ec
EFmax = Ec
#EFがepsより小さくなったら敬さん終了
eps = 1.0e-5
#2分法の最大繰り返し数
nmaxiter = 200
#繰り返し中に途中経過を出力するサイクル数
iprintiterval = 1
# Fermi-Dirac function
[ドキュメント]
def fe(E, EF, T):
"""
フェルミ・ディラック分布関数を計算します。
詳細説明:
電子または正孔が特定のエネルギー準位を占有する確率を返します。
グローバル変数 `ekBT` を使用します。
:param E: エネルギー準位 (eV)。
:type E: float
:param EF: フェルミエネルギー (eV)。
:type EF: float
:param T: 温度 (K)。
:type T: float
:returns: フェルミ・ディラック分布関数の値。
:rtype: float
"""
global Nc, Ec, kB # これらのグローバル変数は使用されていませんが、既存のコードに合わせます。
return 1.0 / (exp((E - EF) * ekBT) + 1.0)
# electron density
[ドキュメント]
def Ne(EF, T):
"""
伝導帯中の電子密度を計算します。
詳細説明:
伝導帯下端エネルギー `Ec` と伝導帯有効状態密度 `Nc` を用いて、
ボルツマン近似に基づく電子密度を計算します。
グローバル変数 `Nc`, `Ec`, `ekBT` を使用します。
:param EF: フェルミエネルギー (eV)。
:type EF: float
:param T: 温度 (K)。
:type T: float
:returns: 電子密度 (cm^-3)。
:rtype: float
"""
global Nc, Ec, kB # これらのグローバル変数は使用されていませんが、既存のコードに合わせます。
return Nc * exp(-(Ec - EF) * ekBT)
# hole density
[ドキュメント]
def Nh(EF, T):
"""
価電子帯中の正孔密度を計算します。
詳細説明:
価電子帯上端エネルギー `Ev` と価電子帯有効状態密度 `Nv` を用いて、
ボルツマン近似に基づく正孔密度を計算します。
グローバル変数 `Nv`, `Ev`, `ekBT` を使用します。
:param EF: フェルミエネルギー (eV)。
:type EF: float
:param T: 温度 (K)。
:type T: float
:returns: 正孔密度 (cm^-3)。
:rtype: float
"""
global Nv, Ev, kB # これらのグローバル変数は使用されていませんが、既存のコードに合わせます。
return Nv * exp(-(EF - Ev) * ekBT)
# ionized donor density
[ドキュメント]
def NDp(EF, T):
"""
ドナー準位がイオン化したドナー密度を計算します。
詳細説明:
ドナー密度 `ND` とドナー準位 `ED` を用いて、
フェルミ・ディラック分布関数に基づいてイオン化されたドナーの密度を計算します。
グローバル変数 `ND`, `ED` および `fe` 関数を使用します。
:param EF: フェルミエネルギー (eV)。
:type EF: float
:param T: 温度 (K)。
:type T: float
:returns: イオン化ドナー密度 (cm^-3)。
:rtype: float
"""
global ND, ED, kB # これらのグローバル変数は使用されていませんが、既存のコードに合わせます。
return ND * (1.0 - fe(ED, EF, T))
# ionized acceptor density
[ドキュメント]
def NAm(EF, T):
"""
アクセプター準位がイオン化したアクセプター密度を計算します。
詳細説明:
アクセプター密度 `NA` とアクセプター準位 `EA` を用いて、
フェルミ・ディラック分布関数に基づいてイオン化されたアクセプターの密度を計算します。
グローバル変数 `NA`, `EA` および `fe` 関数を使用します。
:param EF: フェルミエネルギー (eV)。
:type EF: float
:param T: 温度 (K)。
:type T: float
:returns: イオン化アクセプター密度 (cm^-3)。
:rtype: float
"""
global NA, EA, kB # これらのグローバル変数は使用されていませんが、既存のコードに合わせます。
return NA * fe(EA, EF, T)
[ドキュメント]
def main():
"""
半導体のフェルミエネルギーを2分法で計算します。
詳細説明:
電荷中性条件 (Ne + NAm - Nh - NDp = 0) を満たすフェルミエネルギー (EF) を探索します。
初期範囲 [EFmin, EFmax] が適切であることを確認し、
指定された精度 `eps` または最大繰り返し数 `nmaxiter` に達するまで計算を繰り返します。
グローバル変数 `EFmin`, `EFmax`, `eps`, `nmaxiter`, `iprintinterval`
および `T` を使用し、各密度計算関数を呼び出します。
:returns: 計算が成功した場合は1、失敗した場合は0を返します。
:rtype: int
"""
global EFmin, EFmax, eps, nmaxiter, iprintiterval # 既存のロジックに合わせます。
print("Solution of EF by bisection method")
print("")
# まず、EFmin,EFmaxにおけるΔQを計算し、それぞれが正・負あるいは負・生となっていることを確認する
dQmin = Ne(EFmin, T) + NAm(EFmin, T) - Nh(EFmin, T) - NDp(EFmin, T)
dQmax = Ne(EFmax, T) + NAm(EFmax, T) - Nh(EFmax, T) - NDp(EFmax, T)
print(" EFmin = {:12.8f} dQmin = {:12.4g}".format(EFmin, dQmin))
print(" EFmax = {:12.8f} dQmax = {:12.4g}".format(EFmax, dQmax))
if dQmin * dQmax > 0.0:
print("Error: Initial Emin and Emax should be chosen as dQmin * dQmax < 0")
return 0
# 2分法開始
for i in range(nmaxiter):
EFhalf = (EFmin + EFmax) / 2.0
Neh = Ne(EFhalf, T)
NAmh = NAm(EFhalf, T)
Nhh = Nh(EFhalf, T)
NDph = NDp(EFhalf, T)
dQhalf = Neh + NAmh - Nhh - NDph
print(" Iter {}: EFhalf = {:12.8f} dQhalf = {:12.4g}".format(i, EFhalf, dQhalf))
print(" Ne={:10.4e} Nh={:10.4e} NA-={:10.4e} ND+={:10.4e} dQ={:10.4e}".format(Neh, Nhh, NAmh, NDph, dQhalf))
# EFの精度がepsより小さくなったら敬さん終了
if abs(EFmin - EFhalf) < eps and abs(EFmax - EFhalf) < eps:
print(" Success: Convergence reached at EF = {}".format(EFhalf))
return 1
if dQmin * dQhalf < 0.0:
EFmax = EFhalf
dQmax = dQhalf
else:
EFmin = EFhalf
dQmin = dQhalf
else:
print(" Failed: Convergence did not reach")
return 0
if __name__ == "__main__":
main()