Dear all, I make AGNR in Kwant and introduce the defect in both the scattering region and lead. When plotting the bandstructure of the lead there is a straight line that comes at E=0. I am confused that there is some error in the code or not. Please help
Thanking You code- ________________ import kwant import numpy as np import matplotlib.pyplot as plt pi = np.pi a = 1 d = 1.42; a1 = d*np.sqrt(3); lat = kwant.lattice.general( [[a1, 0], [a1/2, a1*np.sqrt(3)/2]],[[0, 0], [0, a1/np.sqrt(3)]] ) a,b = lat.sublattices def set_pub(): # for plotting the figures plt.rcParams.update({ "font.weight": "bold", # bold fonts "font.size":18, #"tick.labelsize": 15, # large tick labels "lines.linewidth": 2, # thick lines "lines.color": "k", # black lines "grid.color": "0.5", # gray gridlines "grid.linestyle": "-", # solid gridlines "grid.linewidth": 0.5, # thin gridlines "xtick.major.size":10, "xtick.minor.size":5, "ytick.major.size":10, "ytick.minor.size":5, "savefig.dpi": 600, # higher resolution output. }) def make_sys(w,l): def rect(pos): x,y = pos return -w*a1<x<w*a1 and -l*a1<y<l*a1 model1 = kwant.Builder() model1[lat.shape(rect,(1,1))] = 0 model1[lat.neighbors()] = 2.46 return model1 def lead_join(model1,w): def lead_shape(pos): x, y = pos return -w*a1<x<w*a1 Num = 8*(a1*np.sqrt(3)/2) sym = kwant.TranslationalSymmetry([0,Num]) lead = kwant.Builder(sym) lead[lat.shape(lead_shape, (1,1))] = 0 lead[lat.neighbors()] = 2.46 model1.attach_lead(lead) model1.attach_lead(lead.reversed()) return model1,lead def defect(model1,lead): def delet_lead(pos, shift = +a1*np.sqrt(3)): x, y = pos z=x**2+(y-shift)**2 return z<(a1)**2 def delet_str(pos, shift = a1*np.sqrt(3)): x, y = pos z = x**2 + (y-shift)**2 return z<(a1)**2#z<(2*a1)**2 del model1[lat.shape(delet_str, (1,1))] del lead[lat.shape(delet_lead, (2,2))] kwant.plot(model1) set_pub() plt.rc('axes', linewidth = 2.0) plt.show() return model1 def plot_bandstructure(flead, momenta): bands = kwant.physics.Bands(flead) energies = [bands(k) for k in momenta] band_gap = 100 for k in momenta: for band in bands(k): if(band <0): #Lower than fermi band min_band=band else: max_band=band if(max_band-min_band<band_gap): band_gap=max_band-min_band break print("Gap Energy= " + str(band_gap) +" ev ") kwant.plotter.bands(flead) # bandstructure of lead plt.ylim([-3,3]) set_pub() plt.xlabel("K $(\AA^{-1})$",fontweight='bold') plt.ylabel("Energy (eV)",fontweight='bold') plt.show() return band_gap def main(): w = 5 l = 8 sys = make_sys(w,l) syst,lead = lead_join(sys,w) syst = defect(syst,lead) fsys = syst.finalized() momenta = [-pi + 0.02 * pi * i for i in range(101)] Band_Gap = plot_bandstructure(lead.finalized(), momenta) if __name__ == '__main__': main()