Hello Kwant community,
I am using the KPM method to calculate conductivities in the simplest QHE
system -- square lattice nearest hopping gas. The goal is to see well-quantized
plateaus up to n=5.
I use local vectors at the center of the system. But it seems not available at
all as far as I've tried increasing system size (like 80*80) and number of
moments (like 500). I doubt if it really needs so demanding computational
resources in order to get n=5. Hopefully I missed something basic here.
Thank you!
from cmath import exp
import numpy as np
import time
import sys
import kwant
#from matplotlib import pyplot
from matplotlib import pyplot as plt
L_sys = [50, 50]
E_F=0.4
def make_system(a=1, t=1):
lat = kwant.lattice.square(a, norbs=1)
syst = kwant.Builder()
def fluxx(site1, site2, Bvec):
pos = [(site1.pos[i]+site2.pos[i]-L_sys[i])/2.0 for i in range(2)]
return exp(-1j * Bvec * pos[1] )
def fluxy(site1, site2, Bvec):
pos = [(site1.pos[i]+site2.pos[i]-L_sys[i])/2.0 for i in range(2)]
return 1 #exp(-1j * (-Bvec * pos[0]))
def hopx(site1, site2, Bvec):
# The magnetic field is controlled by the parameter Bvec
return fluxx(site1, site2, Bvec) * t
def hopy(site1, site2, Bvec):
# The magnetic field is controlled by the parameter Bvec
return fluxy(site1, site2, Bvec) * t
def onsite(site):
return 4*t
#### Define the scattering region. ####
syst[(lat(x, y) for x in range(L_sys[0]) for y in range(L_sys[1]))] = onsite
# hoppings in x-direction
syst[kwant.builder.HoppingKind((1, 0), lat, lat)] = hopx
# hoppings in y-directions
syst[kwant.builder.HoppingKind((0, 1), lat, lat)] = hopy
# Finalize the system.
syst = syst.finalized()
return lat, syst
def cond(lat, syst, random_vecs_flag=False):
center_pos = [x/2 for x in L_sys]
center0 = lambda s: np.linalg.norm(s.pos-center_pos) < 2
num_mmts = 800;
if random_vecs_flag:
center_vecs = kwant.kpm.RandomVectors(syst, where=center0)
one_vec = next(center_vecs)
num_vecs = 10 # len(center_vecs)
else:
center_vecs = np.array(list(kwant.kpm.LocalVectors(syst,
where=center0)))
one_vec = center_vecs[0]
num_vecs = len(center_vecs)
area_per_orb = np.abs(np.cross(*lat.prim_vecs))
norm = area_per_orb * np.linalg.norm(one_vec) ** 2
Binvfields = np.arange(3.01, 26.2, 1)
params = dict(Bvec=0)
dataxx = []
dataxy = []
for Binv in Binvfields:
params['Bvec'] = 1/Binv
cond_xx = kwant.kpm.conductivity(syst, params=params, alpha='x',
beta='x', mean=True,
num_moments=num_mmts, num_vectors=num_vecs,
vector_factory=center_vecs)
cond_xy = kwant.kpm.conductivity(syst, params=params, alpha='x',
beta='y', mean=True,
num_moments=num_mmts, num_vectors=num_vecs,
vector_factory=center_vecs)
dataxx.append(5*cond_xx(mu=E_F, temperature=0.0)/norm)
dataxy.append(-cond_xy(mu=E_F, temperature=0.0)/norm)
print(dataxx)
print(dataxy)
plot_curves(
[
(r'$\sigma_{xx}\times5$',(Binvfields,dataxx)),
(r'$\sigma_{xy}$',(Binvfields,dataxy)),
]
)
def plot_curves(labels_to_data):
plt.figure(figsize=(12,10))
for label, (x, y) in labels_to_data:
plt.plot(x, y, label=label, linewidth=2)
plt.legend(loc=2, framealpha=0.5)
plt.xlabel("1/B")
plt.ylabel(r'$\sigma [e^2/h]$')
plt.show()
plt.clf()
def main():
lat, syst = make_system()
cond(lat, syst)
if __name__ == '__main__':
main()
