Dear Kwant Authors and Users,
Greetings and thank you for a wonderful tool.
While trying to force Kwant to accept the unit cell of the lead in a
particular, connected form, I've encountered somewhat puzzling
situation illustrated in the attached piece of code. The described procedure
does not make sense physically, but its purpose is
to illustrate the behaviour which may possibly (?) lead to problems, at least
for newbies like me. ;)
Let's say I define two honeycomb lattices, differing in the choice of basis,
but entirely equivalent (i.e. they overlap)
prime=np.array([[sqrt(3)/2,sqrt(3)],[1.5,0.0]])*a
basis_at=np.array([[0.0,0.0],[-0.5,0.5]])*a
lat =
kwant.lattice.general([prime[:,0],prime[:,1]],[basis_at[:,0],basis_at[:,1]])
and
basis2=np.array([[0.0,sqrt(3)/2],[-0.5,-1.0]])*a
lat2 =
kwant.lattice.general([prime[:,0],prime[:,1]],[basis2[:,0],basis2[:,1]])
I'm using numpy so that vector algebra works. The rectangular scattering region
and the left lead (same width)
are populated with "lat". Now I'm adding a few atoms from lat2 to the right
edge of the scattering region,
e.g. a single hexagon
for i in range (0,1): # changing the range we can generate also the full
layer of lat2 sites
syst[lat2.sublattices[0](0+2*i,9-i)]=0.0
syst[lat2.sublattices[1](0+2*i,8-i)]=0.0
syst[lat2.sublattices[0](1+2*i,8-i)]=0.0
syst[lat2.sublattices[1](1+2*i,8-i)]=0.0
and then attach the lead, also populated with "lat2". I know, that's not how
it's described in tutorial. ;)
The filling algorithm generates the system which to naked eye looks like the
perfect nanowire.
The transmission however is not perfect i.e. T < N (number of modes). If the
full layer of lat2 sites is added to the right edge,
than everything is as it should be, i.e. perfect transmission.
I'm not sure if it's a bug - I'm doing things not in the manual - but I find it
strange. A perfect system, however generated, should
lead to perfect transmission. What gives?
Also I noticed that neighbours() method called for one of the lattices tends to
pick up the sites from the other lattice,
if they belong to he sublattice generated from (0,-0.5a), i.e. from the common
basis. Try commenting out one of the lines 53-54.
Is it intentional or not?
It's quite likely that the described situation never arises in real-world
situations especially if one does the right thing
and pads the interface with the full principal layer before attaching the
leads. I've stumbled upon it only thanks to my laziness.
But it might also illustrate some fragility of the algorithm, possibly worth
eliminating.
In any case I'd be grateful for some illumination. Why isn't the "perfect" wire
behaving as it should?
Where is the imperfection if it's not visible in the structure?
With Kind Regards
Maciej Zwierzycki
import kwant
import numpy as np
from numpy import sin, cos,tan,sqrt
from matplotlib import pyplot
def make_system(a, t, W, L1):
prime=np.array([[sqrt(3)/2,sqrt(3)],[1.5,0.0]])*a
basis_at=np.array([[0.0,0.0],[-0.5,0.5]])*a
lat =
kwant.lattice.general([prime[:,0],prime[:,1]],[basis_at[:,0],basis_at[:,1]])
syst = kwant.Builder()
#### Define the scattering region. ####
def scatt_shape(pos):
(x, y) = pos
return (-L1 < x < L1) and ( -W <y < W)
syst[lat.shape(scatt_shape,basis_at[:,0])]=0
syst[lat.neighbors()] = t
#### Left lead. ####
def lead_shape(pos):
(x, y) = pos
return (-W < y < W )
lsym=kwant.TranslationalSymmetry(lat.vec((0,-1)))
lsym.add_site_family(lat.sublattices[0], other_vectors=[(-2,1)])
lsym.add_site_family(lat.sublattices[1], other_vectors=[(-2,1)])
llead = kwant.Builder(lsym)
llead[lat.shape(lead_shape,basis_at[:,0])] = 0
llead[lat.neighbors()] = t
#### Right lead. ####
### Define a hexagonal lattice with different - but equivalent! - basis
basis2=np.array([[0.0,sqrt(3)/2],[-0.5,-1.0]])*a
lat2 =
kwant.lattice.general([prime[:,0],prime[:,1]],[basis2[:,0],basis2[:,1]])
#### Extra hexagon of lat2 on the right edge of the scattering region:
#### -- range(0,1) - single hexagon
##### -- range(-2,3) - full width
for i in range (0,1):
syst[lat2.sublattices[0](0+2*i,9-i)]=0.0
syst[lat2.sublattices[1](0+2*i,8-i)]=0.0
syst[lat2.sublattices[0](1+2*i,8-i)]=0.0
syst[lat2.sublattices[1](1+2*i,8-i)]=0.0
syst[lat2.neighbors()]=t
syst[lat.neighbors()]=t
kwant.plot(syst)
lsym=kwant.TranslationalSymmetry(lat2.vec((0,1)))
lsym.add_site_family(lat2.sublattices[0], other_vectors=[(2,-1)])
lsym.add_site_family(lat2.sublattices[1], other_vectors=[(2,-1)])
rlead= kwant.Builder(lsym)
rlead[lat2.shape(lead_shape,basis_at[:,0])] = 0
rlead[lat2.neighbors()] =t
syst.attach_lead(llead)
syst.attach_lead(rlead)
syst=syst.finalized()
return syst
syst = make_system(a=1.42,t=3.0,W=10.5,L1=20)
kwant.plot(syst)
data1 = []
data2 = []
params = dict(a=1.42, t=-3.0)
energies=[0.2*i+0.01 for i in range(20)]
for energy in energies:
print('En=',energy)
smatrix = kwant.smatrix(syst, energy,params=params)
data1.append(smatrix.transmission(1, 0))
data2.append(smatrix.num_propagating(0))
pyplot.figure()
pyplot.plot(energies, data1, 'r--', label='T')
pyplot.plot(energies, data2,'r>', label='N')
legend = pyplot.legend(loc='upper left')
pyplot.xlabel("energy")
pyplot.ylabel("conductance [e^2/h]")
pyplot.show()
