Hello, I encountered an issue when dealing with a semi-infinite ribbon: when defining a 1D translationnal invariance, I end up with a unit cell with unwanted and unexpected additionnal degrees of freedom. The unit cell have additionnal segments of atoms (depending on the translationnal invariance vector). This additionnal segment is not visible with a "kwant.plot()" but appears in the Hilbert space dimensionnality. It then create unexpected states close to zero energy which are not visible on the band structure but are visible in the DOS. I copy the script at the end of the message, if you want to have a look. This is a Haldane model ribbon with armchair edges. I tried to change the translationnal invariance vector but I didn't manage to get around it. Do you have a solution ? Should I use wraparound ? My end goal is to have a large enough system with a kind of periodic boundary condition in one direction, but still preserving a plotable dispersion relation, so I though that using 1D translationnal invariance with a big unit cell was a good idea. Thank you in advance, Alexandre BERNARD ---------------------------------------------------- the script --------------------------------------------------------------- from types import SimpleNamespace import matplotlib import math from matplotlib import pyplot from mpl_toolkits import mplot3d import numpy as np import scipy.linalg as la import random import kwant def random_uni(param): return param*(random.random() - 0.5)*2 # returns a random number between +param and -param def onsite(site, p): # Here we affect a potential whose sign depends on which sublattice the atom lives in (breaking sublattice symmetry) if site.family == A: # sublattice A return p.m + random_uni(p.rand_os) else: # sublattice B return -p.m + random_uni(p.rand_os) def nn_hopping(site1, site2, params): return params.t_1 def nnn_hopping(site1, site2, params): return (math.cos(params.phi) + 1j*math.sin(params.phi)) * (params.t_2) def ribbon_shape_armchair(pos): return (-1 <= pos[0] < w) def hopping_color(site1, site2): if site1.family == A and site2.family == A: return 'b' elif site1.family == B and site2.family == B: return 'r' else: return '0.5' def site_color(site): if site.family == A: return 'b' else: return 'r' np.random.seed(271828) # Fixing the RNG seed for reproducibility # Graphene lattice "pointing up", with two sublattices graphene = kwant.lattice.general([[1, 0], [1/2, np.sqrt(3)/2]], # Here are the vectors of the Bravais lattice [[0, 0], [0, 1/np.sqrt(3)]]) # Here are the positions of the atoms in the unit cell A, B = graphene.sublattices w = 5 # width of the ribbon p = SimpleNamespace(t_1=1.0, t_2=0.1, m=0.0, phi=np.pi/2, rand_os=0.) # Creation of the list of (oriented) bonds that will have the complex hopping term nnn_hoppings_A = (((-1, 0), A, A), ((0, 1), A, A), ((1, -1), A, A)) # Bonds within the A sublattice nnn_hoppings_B = (((1, 0), B, B), ((0, -1), B, B), ((-1, 1), B, B)) # Bonds within the B sublattice nnn_hoppings = nnn_hoppings_A + nnn_hoppings_B haldane_armchair = kwant.Builder(kwant.TranslationalSymmetry((0, 1*np.sqrt(3)))) haldane_armchair[graphene.shape(ribbon_shape_armchair, (0, 0))] = onsite haldane_armchair[graphene.neighbors(1)] = nn_hopping haldane_armchair[[kwant.builder.HoppingKind(*hopping) for hopping in nnn_hoppings]] = nnn_hopping kwant.plot(haldane_armchair, site_color=site_color, hop_color=hopping_color, fig_size=(6, 9)); haldane_f2 = haldane_armchair.finalized() ham = haldane_f2.hamiltonian_submatrix(args = [p]) eigenVal, eigenVec = la.eigh(ham) nb_eVal = len(eigenVal[:]) print("Number of eigenvalues (i.e. unit cell DOF): ", nb_eVal) pyplot.hist(eigenVal[:], 300) pyplot.show() pyplot.hist(eigenVal[:], 100, cumulative=True) pyplot.show() kwant.plotter.bands(haldane_armchair.finalized(), args=[SimpleNamespace(t_1=p.t_1, t_2=p.t_2, m=p.m, phi=p.phi, rand_os=p.rand_os)], fig_size=(12, 8)); E_F = 0.0 for j in range(nb_eVal): if eigenVal[j] > E_F: break print("eigenVal[", j, "] = ", eigenVal[j]) pyplot.figure(figsize=(w,w/2)) for i, s in enumerate(haldane_f2.sites): x,y=s.pos prob = abs(eigenVec[i,j])**2 pyplot.scatter(x,y,c=prob,vmin=0,vmax=5e-3,s=60) pyplot.show()