Dear all,
Using the Kwant package to calculate the spin-resolved density of states in a graphene nanoribbon system, where the central device region contains Rashba spin-orbit coupling (RSOC) in its Hamiltonian while the leads are free of spin-orbit coupling. We aim to obtain the total spin-up density of states distribution in the device region as the electron energy varies from 0 to 1.The attached files contain our program code. We found that when the wire contains a suspended key, we can run the result, but we don't know if the result is correct. When the wire does not have a hovering key, a strange value appears when executing, and an error is reported. We seek advice on the necessary modifications to correct these errors.
For reference and convenience, I attach the code I'm using at the end of this message.
Any help and guidance would be appreciated!
Best regards,
Rease
############################################
import numpy as np
import tinyarray as ta
import kwant
import matplotlib.pyplot as plt
sigma_0 = ta.array([[1, 0], [0, 1]])
sigma_x = ta.array([[0, 1], [1, 0]])
sigma_y = ta.array([[0, -1j], [1j, 0]])
sigma_z = ta.array([[1, 0], [0, -1]])
lat = kwant.lattice.honeycomb(norbs=2)
def make_system(t=1.0, lambda_rsoc=0.1):
syst = kwant.Builder()
def onsite(site):
return sigma_0
def hopping(site1, site2):
dx, dy = site2.pos - site1.pos
rsoc = 1j * lambda_rsoc * (sigma_x * dy - sigma_y * dx)
return -t * sigma_0 + rsoc
def rectangle(pos):
x, y = pos
return 0 <= x < 10 and 0.5 <= y < 10
syst[lat.shape(rectangle, (0, 0))] = onsite
syst[lat.neighbors()] = hopping
def lead_onsite(site):
return sigma_0
def lead_hopping(site1, site2):
return -t * sigma_0
sym_left = kwant.TranslationalSymmetry((-1, 0))
lead_left = kwant.Builder(sym_left, conservation_law=-sigma_z)
lead_left[lat.shape(rectangle, (0, 0))] = lead_onsite
lead_left[lat.neighbors()] = lead_hopping
syst.attach_lead(lead_left)
sym_right = kwant.TranslationalSymmetry((1, 0))
lead_right = kwant.Builder(sym_right, conservation_law=-sigma_z)
lead_right[lat.shape(rectangle, (0, 0))] = lead_onsite
lead_right[lat.neighbors()] = lead_hopping
syst.attach_lead(lead_right)
return syst
syst = make_system().finalized()
def spin_up_dos(syst):
spin_up = 0.5 * (sigma_0 + sigma_z)
return kwant.operator.Density(syst, spin_up, sum=True)
energies = np.linspace(0, 1, 100)
dos = []
for energy in energies:
ldos = kwant.ldos(syst, energy)
dos.append(np.sum(ldos[::2]))
plt.figure(figsize=(8, 6))
plt.plot(energies, dos, label="Spin-Up DOS")
plt.xlabel("Energy [t]")
plt.ylabel("Density of States")
plt.title("Spin-Up Density of States in Graphene Nanoribbon")
plt.legend()
plt.grid(True)
plt.show()
