import os
import re
import csv
import numpy as np
from math import log, exp
import scipy.special
from scipy import optimize
from scipy.optimize import minimize
from matplotlib import pyplot as plt
from math import sin, cos, tan, exp, log, sinh, cosh, tanh, sqrt
"""
Plot band diagram of a Schottky junction
"""
# Physical constants
e = 1.60218e-19 # C";
me = 9.1093897e-31 # kg";
e0 = 8.854418782e-12; # C2N-1m-2";
e2_4pie0 = 2.30711e-28 # Nm2";
a0 = 5.29177e-11 # m";
kB = 1.380658e-23 # JK-1";
# Material and junction parameters
er = 11.8
Eg = 4.12 # eV
ND = 7.5e19 # cm-3
ED = 0.68 # eV
Vbi = 0.42 #eV
phi_SB = 0.97 #eV
# x range to plot band diagram
xmin = -1.0
xmax = 5.0
xstep = 0.02
nx = int((xmax - xmin) / xstep)
# V values to calculate band diagrams under V
xV = [0.0, -0.2, -0.4, -0.6, -0.8, -1.0, -2.0, -3.0, -5.0, -10.0, -20.0]
# V range to plot x(EC=EF) - V characteristics
xVx = [0.0 - i * 0.1 for i in range(200)]
# figure configration
fontsize = 16
legend_fontsize = 16
def main():
k = e * (ND * 1.0e6) / 2.0 / (er * e0)
print("k=", k)
# Width of depletion layer at V = 0
W0 = sqrt(Vbi * 2.0 * er * e0 / e / (ND * 1.0e6)) * 1.0e9
print("W(0)={} nm".format(W0))
ECs0 = phi_SB - Vbi
# Band diagram at V = 0
xx = [xmin + i * xstep for i in range(nx)]
xm = []
yEFm = []
yEFms = []
xs = []
yECs = []
yEVs = []
yEDs = []
for x in xx:
xm.append(x)
yEFm.append(0.0)
if x <= 0.0:
pass
elif x <= W0:
xs.append(x)
dphi = k * (x*x - 2.0 * W0 * x) * 1.0e-18
yECs.append(phi_SB + dphi)
yEVs.append(phi_SB - Eg + dphi)
yEDs.append(phi_SB - ED + dphi)
else:
xs.append(x)
yECs.append(phi_SB - Vbi)
yEVs.append(phi_SB - Vbi - Eg)
yEDs.append(phi_SB - Vbi - ED)
# Band diagram at finite V
xxs = []
yyECs = []
for V in xV:
WV = sqrt((Vbi - V) * 2.0 * er * e0 / e / (ND * 1.0e6)) * 1.0e9
EC = []
xxs = []
for x in xx:
if x < 0.0:
pass
elif x <= WV:
xxs.append(x)
dphi = k * (x*x - 2.0 * WV * x) * 1.0e-18
EC.append(phi_SB + dphi)
else:
xxs.append(x)
dphi = k * (- WV * WV) * 1.0e-18
EC.append(phi_SB + dphi)
pass
yyECs.append(EC)
# x at EC(V) = EF
print("")
print("x at EC(x) = EF")
print("V(V)\tx(nm)")
xEC0 = []
yEC0 = []
yEC0ex = []
for V in xVx:
WV = sqrt((Vbi - V) * 2.0 * er * e0 / e / (ND * 1.0e6)) # m
a = 1.0
b = -2.0 * WV
c = phi_SB / k
sq2 = b * b - 4.0 * a * c
if -V > ECs0 and sq2 >= 0.0:
xcross = phi_SB / (2.0 * k * WV) * 1.0e9
xcrossex = (-b - sqrt(sq2)) / (2.0 * a) * 1.0e9
xEC0.append(V)
yEC0.append(xcross)
yEC0ex.append(xcrossex)
print("{:6.3f}\t{:10.6f}".format(V, xcrossex))
# Plot graphs
fig = plt.figure(figsize = (12, 6))
plt.rcParams["font.size"] = fontsize
ax11 = fig.add_subplot(1, 2, 1)
ax12 = fig.add_subplot(1, 2, 2)
ax11.plot(xm, yEFm, label = 'EF(Pt)', linestyle = 'dashed', linewidth = 1.0, color = 'red')
ax11.plot(xs, yECs, label = 'EC(a-GO)', linestyle = '-', linewidth = 1.0, color = 'black')
# ax11.plot(xs, yEVs, label = 'EV(a-GO)', linestyle = '-', linewidth = 1.0, color = 'black')
ax11.plot(xs, yEDs, label = 'ED(a-GO)', linestyle = 'dashed', linewidth = 0.5, color = 'blue')
ax11.plot([0.0, 0.0], [0.0, phi_SB], linestyle = '-', linewidth = 1.0, color = 'black')
for i in range(len(xV)):
ax11.plot(xxs, yyECs[i], linestyle = 'dotted', linewidth = 0.5, color = 'black')
# ax11.tick_params(labelsize = fontsize)
ax11.set_xlabel("x (nm)")
ax11.set_ylabel("Energy (eV)")
ax11.legend(fontsize = legend_fontsize)
ax11.set_ylim([-1.0, 1.2])
# ax12.plot(xEC0, yEC0, label = 'approx', color = 'red', linewidth = 0.5)
ax12.plot(xEC0, yEC0ex, label = 'exact', color = 'black', linewidth = 1.0)
ax12.set_xlabel("V (V)")
ax12.set_ylabel("x(EC=EF) (nm)")
ax12.set_xlim([min(xEC0), 0.0])
# ax12.tick_params(labelsize = fontsize)
# ax12.legend(fontsize = legend_fontsize)
plt.tight_layout()
plt.pause(0.1)
print("Press ENTER to exit>>", end = '')
input()
if __name__ == '__main__':
main()
