Hello Kwant community,
I want to calculate the local Chern marker in kwant, roughly a real-space
approach to topological numbers such as the Chern number. Basically, one needs
to calculate the expectation value of the commutator [x, y] of position
operators projected to the occupied states. The main reference is Eqs.(6-9) in
[ http://dx.doi.org/10.1103/PhysRevB.84.241106 ], which has been applied in
many other cases.
My understanding is that kwant can define x,y as ???density operators??? and
then this task should fit into the kwant.kpm.correlator framework, because it
is a correlation function at zero temperature between the two operators. So I
tried in the following code a few ways (some commented) for the simplest
square-lattice quantum Hall effect (QHE). But it does not seem to give a sign
of quantization, although I'm sure that the same system gives reasonable QHE in
conductivity and conductance calculation. E.g., for 1/B in [3, 12], it should
show the C=1, C=2 plateaus. I feel that I probably missed or misunderstood
something basic here. Is my use of operator and kpm.correlator correct here?
Thank you for your help!
from cmath import exp
import numpy as np
import matplotlib.pyplot as plt
import kwant
import scipy.sparse.linalg as sla
import scipy.sparse as sparse
import tinyarray
Nx, Ny = 101, 101
Lx, Ly = Nx-1, Ny-1
edge=18
eps = 1e-7
L_sys = [Lx, Ly]
Ef=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.
fsyst = syst.finalized()
return lat, fsyst
lat, fsyst = make_system()
sites = fsyst.sites
X = sparse.diags([site.pos[0] for site in sites])
Y = sparse.diags([site.pos[1] for site in sites])
def where(site):
pos=site.pos
x0_min, y0_min = edge-eps, edge-eps
x0_max, y0_max= Lx-edge+eps, Ly-edge+eps
return x0_min < pos[0] < x0_max and y0_min < pos[1] < y0_max
def op0(site):
return site.pos[0]
def op1(site):
return site.pos[1]
def chern(lat, fsyst):
# correlator=kwant.kpm.Correlator(fsyst, params=params,
operator1=kwant.operator.Density(fsyst, onsite=op0, where=where),\
# operator2=kwant.operator.Density(fsyst, onsite=op1,
where=where),\
# num_vectors=10, num_moments=400,
vector_factory=kwant.kpm.RandomVectors(fsyst, where=where),\
# bounds=(0, 8.), eps=0.01, rng=0, kernel=None,
mean=True, accumulate_vectors=False)
correlator=kwant.kpm.Correlator(fsyst, params=params,
operator1=kwant.operator.Density(fsyst, onsite=op0),\
operator2=kwant.operator.Density(fsyst, onsite=op1),\
num_vectors=10, num_moments=400,
vector_factory=kwant.kpm.RandomVectors(fsyst, where=where),\
bounds=(0, 8.), eps=0.01, rng=0, mean=True,
accumulate_vectors=False)
# correlator=kwant.kpm.Correlator(fsyst, params=params, operator1=X,\
# operator2=Y,\
# num_vectors=10, num_moments=400,
vector_factory=kwant.kpm.RandomVectors(fsyst, where=where),\
# bounds=(0, 8.), eps=0.01, rng=0, kernel=None,
mean=True, accumulate_vectors=False)
c = correlator(mu=Ef, temperature=0.0) / ((Lx-2*edge)*(Ly-2*edge))
print(2*np.pi*c)
print("Chern number: ", np.real(np.pi*1j*(c-c.conj()))) # Chern number
Bfields = []
for BInv in np.arange(3., 12., 2):
Bfields.append(1/BInv)
NB = len(Bfields)
params = dict(Bvec=0)
for B in Bfields:
params['Bvec'] = B
print("B: ", B)
chern(lat, fsyst)
