Dear Enrique,
What I see in your program is that you put a potential
''edge_potential=10'' to help you to color the sites at the edge where you
want to put disorder. In your case, the width is uniform. So you can
achieve your goal just by putting a constraint on the coordinates of your
sites: Ex: W-1>site.pos[1]>1.
The most interesting case, is when the shape of your system is arbitrary.
To choose the edge , you do:
for site in sys.expand(graphene.shape(rec, (0, 1))):
if sys.degree(site)<9 : sys[site]=edge_potential
Please, notice that I used <9 and not <=9 as you did in your program (this
is why everything is green in your plot).
To get rid of the sites at the edge, in front of the leads, just do, after
attaching the leads :
for site in sys.leads[0].interface+sys.leads[1].interface:
sys[site]=0
This way, your program will work for arbitrary shape, and you can think how
to put disorder after this.
I hope that this helps.
Adel
On Wed, Feb 8, 2017 at 5:31 AM, enrique colomes <[email protected]> wrote:
> Hello,
>
> My name is Enrique and first of all, thanks for the developping of such a
> useful tool.
>
> My question is related to that one of the edge potentials. I am trying to
> reproduce the Haldane model with disorder in a rectangular scattering
> region. My goal is to introduce edge disorder as well as compute the local
> current density.
>
> I have two problems:
>
> a) Regarding the edge disorder, I cannot remove the disorder in the edges
> with the leads (I saw your recent previous discussion)
>
> b) Is it possible already to compute directly from Kwant the local current
> density? If yes, how?
>
> Thanks again,
>
> Enrique
>
>
> from __future__ import division # so that 1/2 == 0.5, and not 0
>
> from scipy.sparse import spdiags
> from math import pi, sqrt, tanh, e
> import numpy as np
> import kwant
> import random
> import math
>
> # For computing eigenvalues
> import scipy.sparse.linalg as sla
>
> # For plotting
> from matplotlib import pyplot
>
> # Define the graphene lattice
> graphene = kwant.lattice.honeycomb()
> a, b = graphene.sublattices
>
>
> edge_potential=10
>
> def make_system(L=10, W=10, pot=0.001):
>
> t=1
> tt= 0.04 * t
> phi=1*pi/2
> edge=1 # 1 if there is edge disorder
> bulk=0 # 1 if there is bulk disorder
> nnn_hoppinga=tt*e**(1j*phi)
> nnn_hoppingb=tt*e**(1j*phi)
>
> #### Define the scattering region. ####
> def rec(pos):
> x, y = pos
> return 0 < x < L and 0 < y < W
>
> sys = kwant.Builder()
>
> # w: width and pot: potential maximum of the p-n junction
> def potential(site):
> (x, y) = site.pos
> d = y * 2 + x * 3
> return 0
>
> sys[graphene.shape(rec, (1, 1))] = potential
>
> # specify the hoppings of the graphene lattice
> nn_hoppings = (((0, 0), a, b), ((0, 1), a, b), ((-1, 1), a, b))
> nnn_hoppings_a = (((-1, 0), a, a), ((0, 1), a, a), ((1, -1), a, a))
> nnn_hoppings_b = (((1, 0), b, b), ((0, -1), b, b), ((-1, 1), b, b))
> sys[[kwant.builder.HoppingKind(*hopping) for hopping in nn_hoppings]] = -t
> sys[[kwant.builder.HoppingKind(*hopping) for hopping in nnn_hoppings_a]]
> = -nnn_hoppinga
> sys[[kwant.builder.HoppingKind(*hopping) for hopping in nnn_hoppings_b]]
> = -nnn_hoppingb
>
>
> # Modify the scattering region
> if edge==1:
> # change the values of the potential at the edge
> for site in sys.expand(graphene.shape(rec, (0, 1))):
> if sys.degree(site)<=9 and abs(site.pos[0])!=L-4:
> #abs(site.pos[0])!=19 execludes the interfaces with the leads
> sys[site]=edge_potential
>
> elif bulk==1:
> # change the values of the potential at the edge
> for site in sys.expand(graphene.shape(rec, (0, 1))):
> if sys.degree(site)==9 and abs(site.pos[0])!=19:
> #abs(site.pos[0])!=19 execludes the interfaces with the leads
> sys[site]=edge_potential
>
>
> #### Define the leads. ####
> # left lead
> sym0 = kwant.TranslationalSymmetry(graphene.vec((-1, 0)))
>
> sym0.add_site_family(graphene.sublattices[0], other_vectors=[(-1, 2)])
> sym0.add_site_family(graphene.sublattices[1], other_vectors=[(-1, 2)])
>
> def lead0_shape(pos):
> x, y = pos
> return (0 < y < W)
>
> lead0 = kwant.Builder(sym0)
> lead0[graphene.shape(lead0_shape, (0, 0))] = -pot
> lead0[[kwant.builder.HoppingKind(*hopping) for hopping in nn_hoppings]] =
> -t
> lead0[[kwant.builder.HoppingKind(*hopping) for hopping in
> nnn_hoppings_a]] = -nnn_hoppinga
> lead0[[kwant.builder.HoppingKind(*hopping) for hopping in
> nnn_hoppings_b]] = -nnn_hoppingb
>
>
> # right lead
> sym1 = kwant.TranslationalSymmetry(graphene.vec((1, 0)))
>
> sym1.add_site_family(graphene.sublattices[0], other_vectors=[(-1, 2)])
> sym1.add_site_family(graphene.sublattices[1], other_vectors=[(-1, 2)])
> def lead1_shape(pos):
> x, y = pos
> return (0 < y < W)
>
> lead1 = kwant.Builder(sym1)
> lead1[graphene.shape(lead1_shape, (0, 0))] = pot
> lead1[[kwant.builder.HoppingKind(*hopping) for hopping in nn_hoppings]] =
> -t
> lead1[[kwant.builder.HoppingKind(*hopping) for hopping in
> nnn_hoppings_a]] = -nnn_hoppinga
> lead1[[kwant.builder.HoppingKind(*hopping) for hopping in
> nnn_hoppings_b]] = -nnn_hoppingb
>
> sys.attach_lead(lead1)
> sys.attach_lead(lead0)
>
> return sys
>
>
> # Call the main function if the script gets executed (as opposed to imported).
>
> if __name__ == '__main__':
> main()t.show()
>
>
> def main():
>
> sys = make_system()
>
> # Check that the system looks as intended.
>
> def family_color(site):
> #print(sys[site])
> if sys[site]==edge_potential:
> return 'green'
> else: return 'black'
>
> def site_size(site):
> if sys[site]==edge_potential:
> return 0.35
> else:return 0.2
> kwant.plot(sys,site_color=family_color,site_size=site_size)
> # Finalize the system.
> sys = sys.finalized()
>
>
> # Compute number of neighbors
> i=sys.sites.index(graphene(5,6))
> all_the_neighbors=sys.graph.out_neighbors(i)
>
>
>
> # Call the main function if the script gets executed (as opposed to imported).
> # See <http://docs.python.org/library/__main__.html>
> <http://docs.python.org/library/__main__.html>.
> if __name__ == '__main__':
> main()
>
>
>
--
Abbout Adel
