Hello All,

I am trying to simulate black phosphorous quantum dots of circular shape, but with their ends (that will be attached to the leads) will be broad. The shape would be something like a circle with two thick lines on either side of the circle. I have the following code for it, but the problem I face is that the two thick lines are disconnected from the circle. I am not able to understand the reason behind it and rectify it. If any of you could help me, it would be of great value to me.

------ CODE -------

import kwant

from matplotlib import pyplot

import numpy as np

h = 2.707 * np.cos(13.69*np.pi/180)

dx = 4.43

dy = 3.27

shift_x = 2.164 * np.cos(0.5*98.15*np.pi/180)

shift_y = 2.164 * np.sin(0.5*98.15*np.pi/180)

r_const = 20

lat0 = kwant.lattice.general([(dx, 0, 0), (0, dy, 0)], [(0, 0, 0), (shift_x, shift_y, 0)])

lat1 = kwant.lattice.general([(dx, 0, 0), (0, dy, 0)], [(0.5*dx, 0.5*dy, h), (0.5*dx + shift_x, 0.5*dy + shift_y, h)])

def cuboid_shape(pos):

x, y, z = pos

if -1.2*r_const <= x < -r_const:

return abs(y) < r_const

if -r_const <= x <= r_const:

return x**2 + y**2 <= r_const**2

if r_const < x <= 1.2*r_const:

return abs(y) < r_const

def potential(site):

x, y, z = site.pos

return (-tanh(5*(x-r_const)) - tanh(5*(x+r_const))) * pot

sys = kwant.Builder()

sys[lat0.shape(cuboid_shape, (0, 0, 0))] = potential

sys[lat0.neighbors(1)] = 1

sys[lat0.neighbors(2)] = 1

sys[lat1.shape(cuboid_shape, (0.5*dx, 0.5*dy, h))] = potential

sys[lat1.neighbors(1)] = 1

sys[lat1.neighbors(2)] = 1

lat0_0 = [ ]

lat0_1 = [ ]

lat1_0 = [ ]

lat1_1 = [ ]

for site in kwant.Builder.sites(sys):

if site.family == lat0.sublattices[0]:

lat0_0.append(site[1])

elif site.family == lat0.sublattices[1]:

lat0_1.append(site[1])

elif site.family == lat1.sublattices[0]:

lat1_0.append(site[1])

elif site.family == lat1.sublattices[1]:

lat1_1.append(site[1])

for w in lat0_0:

for v in lat0_1:

for k in lat1_0:

for l in lat1_1:

if v == k:

sys[lat0.sublattices[1](v[0],v[1]), lat1.sublattices[0](k[0],k[1])] = 2

if w == l:

if (w[0]+1,w[1]+1) in lat0_0:

sys[lat0.sublattices[0](w[0]+1,w[1]+1), lat1.sublattices[1](l[0],l[1])] = 2

if w == l:

if (l[0]-1,l[1]-1) in lat1_1:

sys[lat0.sublattices[0](w[0],w[1]), lat1.sublattices[1](l[0]-1,l[1]-1)] = 2

kwant.plot(sys, site_size=0.2, hop_lw=0.1)

Any help would be greatly appreciated.

Regards,

Shivang Agarwal

Shivang Agarwal

Senior Undergraduate

Coordinator - Academic Discussion Hours

Discipline of Electrical Engineering

IIT Gandhinagar

Contact: +91-9869321451