Oh, sorry I forgot the code in my post... Here is it.

from math import cos, sin, atan, atan2, pi
import numpy as np
import time
import kwant
from matplotlib import pyplot

mx=3
my=2
alpha=atan2(my,mx) # commensurate rotation angle
neg_direc = [(-mx,-my),(my,-mx)] # (negative) X-axis & Y-axis rotated from old 
x- & y-axis
newa=(mx**2+my**2)**0.5 # supercell length
a_sys = [newa, newa] # supercell is square-like
N_sys = [5, 3] # system size along the new XY-axes in units of newa
L_sys = [(N_sys[i]-1)*a_sys[i] for i in range(2)] # real-space system size 
along the new XY-axes 
Rinvmat = np.array([[cos(alpha),-sin(alpha)],[sin(alpha),cos(alpha)]]) # 
rotation matrix
Rmat = np.array([[cos(alpha),sin(alpha)],[-sin(alpha),cos(alpha)]]) # inverse 
rotation matrix

direcI = 0 # current/lead along which rotated axis: 0 for X, 1 for Y
eta = 1.0e-7 # numerical accuracy bound

def checkedge(a,b,c): # range check that deals with numeric error (when a site 
sits exactly at the edge)
    return ((a<=b<=c) or (abs(b-c)<eta) or (abs(b-a)<eta))

def central_reg(pos): # set the rectangular central region
    rot_pos = Rmat.dot(pos)
    for i in range(2):
        if not checkedge(0, rot_pos[i], L_sys[i]):
            return False
    return True

def lead_reg(pos): # confine the lead region
    rot_pos = Rmat.dot(pos)
    for i in range(2):
        if direcI != i:
            if not checkedge(0, rot_pos[i], L_sys[i]):
                return False
    return True

def make_system(a=1, t=1):
    lat = kwant.lattice.square(a)
    syst = kwant.Builder()

    syst[lat.shape(central_reg, (0,0))] = 4*t
    syst[kwant.builder.HoppingKind((1, 0), lat, lat)] = t
    syst[kwant.builder.HoppingKind((0, 1), lat, lat)] = t

    lead0 = kwant.Builder(kwant.TranslationalSymmetry(neg_direc[direcI]))
    lead0[lat.shape(lead_reg, (0,0))] = 4*t # (0,0) should work for I-leads as 
they fully cover the boundary surface
    lead0[kwant.builder.HoppingKind((1, 0), lat, lat)] = t
    lead0[kwant.builder.HoppingKind((0, 1), lat, lat)] = t

    #syst.attach_lead(lead0) # Kwant Lead 0
    #syst.attach_lead(lead0.reversed()) # Kwant Lead 1

    syst = syst.finalized()
    return syst

def main():
    syst = make_system()
    # Check that the system looks as intended.
    kwant.plot(syst, file="lattice_fig.pdf")

if __name__ == '__main__':
    main()

Reply via email to