Dear Henry, The problem you experience comes from the fact, that if you specify `norbs` name of the lattice changes from '<Monatomic lattice e0>' to '<Monatomic lattice e0 with 1 orbitals>' while you conditionally insert hoppings, comparing site family string representation with the first variant. The best approach to fix this is to compare site families directly and not with their string representation, this will be faster and more reliable. With best regards, Viacheslav Ostroukh P.S. @Abbout Adel: special thanks for cleaning up the code, helped a lot! On 16/04/2021 02:52, Abbout Adel wrote:
Dear Christoph,
I simplified the posted program in order to track the error. It seems that it comes from this part:
###################################### for site in lead2.sites(): (x,y) = site.tag if str(site.family) == "<Monatomic lattice e0>": lead2[c(x,y),a(x,y)] =1 if str(site.family) == "<Monatomic lattice e1>": lead2[d(x,y),b(x,y)] = 1 ######################################
I could not go deeper to find the source of the error. The simplified program which shows the difference between"norbs=None" and "norbs=1" is pasted below:
import kwant from matplotlib import pyplot from math import pi, sqrt, tanh, cos, ceil from cmath import exp import numpy as np import time import textwrap as tw sin_30, cos_30, tan_30 = (1 / 2, sqrt(3) / 2, 1 / sqrt(3))
def create_system(N=31, Norbs=None ):
def rectangle(pos): x,y = pos if (0 < x <= 0.4*6.5) and (-3<= y <= 3): return True return False
def lead_shape(pos): x, y = pos return - 3 <= y <= 3 graphene_e = kwant.lattice.honeycomb(a=1,name='e', norbs=Norbs) graphene_h = kwant.lattice.honeycomb(a=1,name='h', norbs=Norbs) a, b = graphene_e.sublattices c, d = graphene_h.sublattices
sys = kwant.Builder()
sys[graphene_e.shape(rectangle, (0.5, 0))] = 0.3 sys[graphene_h.shape(rectangle, (0.5, 0))] = -0.3 sys[graphene_e.neighbors()] = -1 sys[graphene_h.neighbors()] = 1
# Symmetry of the left and right leads sym_l = kwant.TranslationalSymmetry((-1, 0)) sym_r = kwant.TranslationalSymmetry((1, 0))
lead0 = kwant.Builder(sym_l) # Left electron lead0[graphene_e.shape(lead_shape, (1,1))] = -1 lead0[graphene_e.neighbors()] = -1
lead1 = kwant.Builder(sym_l) # Left hole lead1[graphene_h.shape(lead_shape, (1,1))] = 1 lead1[graphene_h.neighbors()] = 1
lead3 = kwant.Builder(sym_r) # Right electron lead3[graphene_e.shape(lead_shape, (1,1))] = -0.3 lead3[graphene_e.neighbors()] = -1
lead4 = kwant.Builder(sym_r) # Right hole lead4[graphene_h.shape(lead_shape, (1,1))] = 0.3 lead4[graphene_h.neighbors()] = 1
lead2 = kwant.Builder(sym_r) # Right superconducting lead2.update(lead3) lead2.update(lead4)
# Hopping between electron and hole lead lattices in superconducting lead for site in lead2.sites(): (x,y) = site.tag
if str(site.family) == "<Monatomic lattice e0>": lead2[c(x,y),a(x,y)] =1 if str(site.family) == "<Monatomic lattice e1>": lead2[d(x,y),b(x,y)] = 1 # Attaching the left electron, left hole and right superconducting lead sys.attach_lead(lead0) sys.attach_lead(lead1) sys.attach_lead(lead2)
return sys, [lead0, lead1, lead2]
def plot_conductance(syst, energies, lead_in, lead_out):
data = [] for energy in energies: smatrix = kwant.smatrix(syst, energy) data.append(smatrix.transmission(lead_out, lead_in)) pyplot.figure() pyplot.plot(energies, data)
return data
def main():
#case 1 sys, leads = create_system(Norbs=None )
sys = sys.finalized() kwant.plot(sys,fig_size=(20,20))
Es = np.linspace(0.0001, 2.00, 10) data = plot_conductance(sys, Es, 0,1) #case 2 sys, leads = create_system(Norbs=1 )
sys = sys.finalized() kwant.plot(sys,fig_size=(20,20))
Es = np.linspace(0.0001, 2.00, 10) data = plot_conductance(sys, Es, 0,1)
if __name__ == "__main__": main()
On Mon, Apr 12, 2021 at 11:38 PM Christoph Groth <christoph.groth@cea.fr <mailto:christoph.groth@cea.fr>> wrote:
Henry, Adel, thanks a lot for reporting this problem and finding a workaround.
Setting norbs may be required in some cases, but computation results of Kwant should not depend on whether norbs is set or not! This could be a rather nasty and serious bug, I would like to look closely into it.
Could you show me the exact workaround that you have found? That will hopefully help me to track down this problem.
Cheers Christoph
-- Abbout Adel
-- With best regards, Slava Ostroukh