I’ve been playing around with Kwant and come across some situations where the transmission between two leads is identically equal to zero when I wouldn’t expect that result.
I’ve come up with a simple working example: a 1D step function with a step height that’s large relative to the lattice spacing (code below), altho this problem seems to sometimes crop up in other, somewhat less extreme situations.
I understand that the numerical result should become less accurate as the step-height-to-lattice-spacing ratio increases, but why does the transmission become identically equal to zero at some point? Are there well-defined conditions for when this happens?
Is there some way to know that the transmission is zero because of numerical issues rather than the underlying physics?
(Ultimately, I’m interested in modeling some systems where the potential mostly varies gradually but has a few small regions with abrupt changes in potential. Moving to a finer mesh (smaller lattice constant) everywhere is cost-prohibitive. Having some
tool to easily refine the mesh in a region would be very useful.)
Thanks.
# In[1]:
get_ipython().magic('load_ext autoreload')
get_ipython().magic('autoreload 2')
from numpy import *
import matplotlib.pyplot as plt
get_ipython().magic('matplotlib inline')
import tqdm
import kwant
# In[2]:
m0 = 9.10938215e-31 # Electron mass, [kg]
hbar = 1.054571726e-34 # hbar in [J] [s]
q = 1.602176565e-19 # Elementary charge, [C]
mt = 0.19
ml = 0.92
m = mt*m0
# In[3]:
V0 = 100 # step height
x = linspace(0,100,30) # thirty grid points
U = zeros_like(x)
U[len(x)//2:] = V0
plt.plot(x,U)
# In[4]:
a = x[1]-x[0] # grid spacing [nm]
t = hbar**2/(2.*m*(a*1e-9)**2)/q*1e3 #hopping parameter [meV]
lat = kwant.lattice.chain(a) # Set up the transport simulation on a 1D latice
sys = kwant.Builder() # initialize the transport simulation
for i in range(len(U)): # populate based on the potential landscape
sys[lat(i)]=U[i]+2*t
sys[lat.neighbors()] = -t # set the finite-difference hopping parameters
leftLead = kwant.Builder(kwant.TranslationalSymmetry((-a,))) # the lead to the left
leftLead[lat(0)] = 2*t + U[0]
leftLead[lat.neighbors()] = -t
sys.attach_lead(leftLead) # attach it
rightLead = kwant.Builder(kwant.TranslationalSymmetry((a,))) # the lead to the right
rightLead[lat(0)] = 2*t + U[-1]
rightLead[lat.neighbors()] = -t
sys.attach_lead(rightLead) # attach it
sys = sys.finalized()
# In[5]:
def plot_conductance(sys, energies):
# Compute transmission numerically
data = []
for energy in tqdm.tqdm(energies,leave=True):
smatrix = kwant.smatrix(sys, energy)
data.append(smatrix.transmission(1, 0))
# Compute exact conductance
k1 = sqrt(2*m*energies/hbar**2)
k2 = sqrt(2*m*(energies - V0)/hbar**2)
T = 4*k1*k2/(k1+k2)**2
T[energies <= V0] = 0
plt.figure()
plt.plot(energies, data, energies, T)
plt.legend(('numerical','exact'), loc=4)
plt.xlabel("energy [V0]")
plt.ylabel("Transmission")
plt.show()
return data, T
# In[6]:
stepNumerical, stepExact = plot_conductance(sys,linspace(1e-9,2*V0,201))
# In[7]:
print(stepNumerical)