Hi everybody,
In the following code I want to delete few random sites from the scattering
region.
Please help me in this case.
Best regards,
Hossein
import kwant
from math import sqrt
import matplotlib.pyplot as plt
import tinyarray
import numpy as np
import math
import matplotlib
#Constants.........................................................
#Parameter.........................................................
a0=2.617*sqrt(3);
d0=2.617;
dz0=1.59047;
D0=sqrt((d0)**2+(dz0)**2);
#st=0;
st=-0.0;
dz=1.59047+st;
d=-0.3*dz+3.09;
a1=-0.3*sqrt(3)*dz+5.35;
D=sqrt((d)**2+(dz)**2);
scale=(D/D0)**(-8);
qo=[];
N=dz/D;
delta1=0j;
delta2=0;
R=(2*a1)**2;
#on-site energy...................................................
Es=-10.906;
Ep=-0.486;
Ep1=-0.486;
Ep2=-0.486;
#nearest neighbor................................................
sssigma=-0.608*scale;
spsigma=1.32*scale;
ppsigma=1.854*scale;
pppay=-0.6*scale;
#next nearest neighbor.........................................
ppsigma2=0.156*scale;
pppay2=0;
sssigma2=0;
spsigma2=0;
#onsite potentials ...........................................
pot_a=10;
pot_b=0;
#nearest neighbors.............................................................
l1=(math.cos(math.pi/6))*(d/D);
m1=(math.sin(math.pi/6))*(d/D);
N1=dz/D;
#second neighbors..............................................................
l2=(math.cos(math.pi/3))
m2=(math.sin(math.pi/3))
latt = kwant.lattice.general([(a1,0),(a1*0.5,a1*math.sqrt(3)/2)],
[(a1/2,-d/2),(a1/2,d/2)])
a,b = latt.sublattices
syst= kwant.Builder()
# Onsite matrice energt
Onsite_a=tinyarray.array([[Es,0,0,0],[0,Ep,delta1,0],[0,-delta1,Ep,0],
[0,0,0,Ep1]])
Onsite_b=tinyarray.array([[Es,0,0,0],[0,Ep,delta1,0],[0,-delta1,Ep,0],
[0,0,0,Ep2]])
pot_edge_a=tinyarray.array([[Es+pot_a,0,0,0],[0,Ep+pot_a,delta1,0],[0,-delta1,Ep+pot_a,0],
[0,0,0,Ep1+pot_a]])
pot_edge_b=tinyarray.array([[Es+pot_b,0,0,0],[0,Ep+pot_b,delta1,0],[0,-delta1,Ep+pot_b,0],
[0,0,0,Ep2+pot_b]])
#nearest neighbors.............................................................
first_up=tinyarray.array([[sssigma, l1*spsigma,m1*spsigma,N1*spsigma],
[-l1*spsigma,
l1**2*ppsigma+(1-l1**2)*pppay,l1*m1*(ppsigma-pppay),l1*N1*(ppsigma-pppay)],
[-m1*spsigma ,l1*m1*(ppsigma-pppay),
m1**2*ppsigma+(1-m1**2)*pppay,N1*m1*(ppsigma-pppay)],
[-N1*spsigma,l1*N1*(ppsigma-pppay),N1*m1*(ppsigma-pppay),N1**2*ppsigma+(1-N1**2)*pppay]])
first_left_up=tinyarray.array([[sssigma,-l1*spsigma,m1*spsigma,N1*spsigma],
[l1*spsigma, l1**2*ppsigma+(1-l1**2)*pppay,-
l1*m1*(ppsigma-pppay),-l1*N1*(ppsigma-pppay)],
[-m1*spsigma ,-l1*m1*(ppsigma-pppay),
m1**2*ppsigma+(1-m1**2)*pppay,N1*m1*(ppsigma-pppay)],
[-N1*spsigma,-l1*N1*(ppsigma-pppay),N1*m1*(ppsigma-pppay),N1**2*ppsigma+(1-N1**2)*pppay]])
first=tinyarray.array([[sssigma,0,-(d/D)*spsigma,N1*spsigma],
[0,pppay,0,0],
[(d/D)*spsigma ,0,
(d/D)**2*ppsigma+(1-(d/D)**2)*pppay,-N1*(d/D)*(ppsigma-pppay)],
[-N1*spsigma,0,-N1*(d/D)*(ppsigma-pppay),N1**2*ppsigma+(1-N1**2)*pppay]])
#second neighbors..............................................................
second_up=tinyarray.array([[sssigma2,l2*spsigma2,m2*spsigma2,0],
[-l2*spsigma2,l2**2*ppsigma2+(1-l2**2)*pppay2,l2*m2*(ppsigma2-pppay2),0],
[-m2*spsigma2
,l2*m2*(ppsigma2-pppay2),m2**2*ppsigma2+(1-m2**2)*pppay2,0],
[0,0,0,pppay2]])
second_left_up=tinyarray.array([[sssigma2,-l2*spsigma2,m2*spsigma2,0],
[l2*spsigma2,l2**2*ppsigma2+(1-l2**2)*pppay2,-l2*m2*(ppsigma2-pppay2),0],
[-m2*spsigma2
,-l2*m2*(ppsigma2-pppay2),m2**2*ppsigma2+(1-m2**2)*pppay2,0],
[0,0,0,pppay2]])
second_right_down=tinyarray.array([[sssigma2,l2*spsigma2,-m2*spsigma2,0],
[-l2*spsigma2,l2**2*ppsigma2+(1-l2**2)*pppay2,-l2*m2*(ppsigma2-pppay2),0],
[m2*spsigma2,-l2*m2*(ppsigma2-pppay2),m2**2*ppsigma2+(1-m2**2)*pppay2,0],
[0,0,0,pppay2]])
second_down=tinyarray.array([[sssigma2,-l2*spsigma2,-m2*spsigma2,0],
[l2*spsigma2,l2**2*ppsigma2+(1-l2**2)*pppay2,l2*m2*(ppsigma2-pppay2),0],
[m2*spsigma2
,l2*m2*(ppsigma2-pppay2),m2**2*ppsigma2+(1-m2**2)*pppay2,0],
[0,0,0,pppay2]])
second_right=tinyarray.array([[sssigma2,spsigma2,0,0],
[-spsigma2,ppsigma2,0,0],[0 ,0,pppay2,0],
[0,0,0,pppay2]])
second_left=tinyarray.array([[sssigma2,-spsigma2,0,0],
[spsigma2,ppsigma2,0,0],[0 ,0,pppay2,0],
[0,0,0,pppay2]])
#...................................................................................
def rectangle(pos):
x, y = pos
# z=x**2+y**2
return -9.2*a1<x<9.2*a1 and -12*d<=y<12*d
syst[a.shape(rectangle, (1,1))] = Onsite_a
syst[b.shape(rectangle, (1,1))] = Onsite_b
#nearest neighbors.............................................................
syst[[kwant.builder.HoppingKind((0,0),a,b)]] = first
syst[[kwant.builder.HoppingKind((0,1),a,b)]] = first_up
syst[[kwant.builder.HoppingKind((-1,1),a,b)]] = first_left_up
#second neighbors..............................................................
syst[[kwant.builder.HoppingKind((0,1),a,a)]]=second_up
syst[[kwant.builder.HoppingKind((-1,1),a,a)]]=second_left_up
syst[[kwant.builder.HoppingKind((1,-1),a,a)]]=second_right_down
syst[[kwant.builder.HoppingKind((0,-1),a,a)]]=second_down
syst[[kwant.builder.HoppingKind((1,0),a,a)]]=second_right
syst[[kwant.builder.HoppingKind((-1,0),a,a)]]=second_left
syst[[kwant.builder.HoppingKind((0,1),b,b)]]=second_up
syst[[kwant.builder.HoppingKind((-1,1),b,b)]]=second_left_up
syst[[kwant.builder.HoppingKind((1,-1),b,b)]]=second_right_down
syst[[kwant.builder.HoppingKind((0,-1),b,b)]]=second_down
syst[[kwant.builder.HoppingKind((1,0),b,b)]]=second_right
syst[[kwant.builder.HoppingKind((-1,0),b,b)]]=second_left
sym = kwant.TranslationalSymmetry(latt.vec((1,0)))
sym.add_site_family(latt.sublattices[0], other_vectors=[(-1, 2)])
sym.add_site_family(latt.sublattices[1], other_vectors=[(-1, 2)])
lead = kwant.Builder(sym)
# ---- random impurity
# --- leads
def lead_shape(pos):
x, y = pos
return -12*d<=y<12*d and -5.2*a1<x<5.2*a1
lead[a.shape(lead_shape, (1,1))] = Onsite_a
lead[b.shape(lead_shape, (1,1))] = Onsite_b
lead[[kwant.builder.HoppingKind((0,0),a,b)]] =first
lead[[kwant.builder.HoppingKind((0,1),a,b)]] =first_up
lead[[kwant.builder.HoppingKind((-1,1),a,b)]] =first_left_up
lead[[kwant.builder.HoppingKind((0,1),a,a)]]=second_up
lead[[kwant.builder.HoppingKind((-1,1),a,a)]]=second_left_up
lead[[kwant.builder.HoppingKind((1,-1),a,a)]]=second_right_down
lead[[kwant.builder.HoppingKind((0,-1),a,a)]]=second_down
lead[[kwant.builder.HoppingKind((1,0),a,a)]]=second_right
lead[[kwant.builder.HoppingKind((-1,0),a,a)]]=second_left
lead[[kwant.builder.HoppingKind((0,1),b,b)]]=second_up
lead[[kwant.builder.HoppingKind((-1,1),b,b)]]=second_left_up
lead[[kwant.builder.HoppingKind((1,-1),b,b)]]=second_right_down
lead[[kwant.builder.HoppingKind((0,-1),b,b)]]=second_down
lead[[kwant.builder.HoppingKind((1,0),b,b)]]=second_right
lead[[kwant.builder.HoppingKind((-1,0),b,b)]]=second_left
syst.attach_lead(lead,add_cells=0)
syst.attach_lead(lead.reversed(),add_cells=0)
fsys = syst.finalized()
def family_color(site):
#print(sys[site])
if syst[site]==Onsite_a:
return 'black'
elif syst[site]==Onsite_b:
return 'green'
elif syst[site]==pot_edge_a:
return 'yellow'
else:
return 'blue'
ax=kwant.plot(syst,site_color=family_color);
ax.savefig('syst.pdf')
