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 > <[email protected] <mailto:[email protected]>> 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
