Hi Michael, Thanks for your help. I modified the hopping of the vertical leads from -t to hopping_Ax. Some LL-like flat bands are observed. However, they are asymmetry about k = 0 (As shown in the attached figures) and it looks a bit strange to me. Is it making sense? I suspect something is still wrong with my code. Thanks for your time. Best, Johnny On 24 August 2016 at 18:43, Michael Wimmer <m.t.wimmer@tudelft.nl> wrote:
Dear Johnny,
you simply did not include any magnetic field in the vertical leads, see lead_hop_vert(). Passing the parameters as you do is fine.
Best,
Michael
Dear Kwant users,
I would like to plot the band structure of my system under a magnetic field. Based on the code demostrated in http://topocondmat.org/w3_pump_QHE/Laughlinargument.html, I plot the band structure but I cannot observe any Landau level like spectrum. I suspect it has something to do with the line
p = SimpleNamespace(t=1.0, mu=0.3, mu_lead=0.3, B=0.5) bands = kwant.physics.Bands(flead,args = [p])
Can the B-field be passed to the Hamiltonian in this way? The complete code is attached below. Thanks for your help.
Best, Johnny
#------My code-----------# from __future__ import division import numpy as np import kwant from matplotlib import pyplot as plt from math import pi
class SimpleNamespace(object): def __init__(self, **kwargs): self.__dict__.update(kwargs)
def qhe_corbino(r_out=100, r_in=65, w_lead=10): """Create corbino disk.
Square lattice, one orbital per site. Returns kwant system.
Arguments required in onsite/hoppings: t, mu, mu_lead, B, phi """ # ring shape def ring(pos): (x, y) = pos rsq = x ** 2 + y ** 2 return (r_in ** 2 < rsq < r_out ** 2)
# Onsite and hoppings def onsite(site, p): (x, y) = site.pos return 4 * p.t - p.mu <http://p.mu>
def crosses_branchcut(hop): x1, y1 = hop[0].pos x2, y2 = hop[1].pos return y1 < 0 and x1 > 0.5 and x2 < 0.5
def hopping(site1, site2, p): x1, y1 = site1.pos x2, y2 = site2.pos # Check for correctness! return -p.t * np.exp(-0.5j * p.B * (x1 - x2) * (y1 + y2))
def branchcut_hopping(site1, site2, p): return hopping(site1, site2, p) * np.exp(1j * p.phi)
# Building system lat = kwant.lattice.square() sys = kwant.Builder()
sys[lat.shape(ring, (0, r_in + 1))] = onsite sys[lat.neighbors()] = hopping
# adding special hoppings def hops_across_cut(sys): for hop in kwant.builder.HoppingKind((1, 0), lat, lat)(sys): if crosses_branchcut(hop): yield hop sys[hops_across_cut] = branchcut_hopping
# Attaching leads sym_lead = kwant.TranslationalSymmetry((-1, 0)) lead = kwant.Builder(sym_lead)
def lead_shape(pos): (x, y) = pos return (-w_lead / 2 < y < w_lead / 2)
lead[lat.shape(lead_shape, (0, 0))] = lambda site, p: 4 * p.t - p.mu_lead lead[lat.neighbors()] = lambda site1, site2, p: -p.t
#### Attach the leads and return the system. #### sys.attach_lead(lead) sys.attach_lead(lead, origin=lat(0, 0))
return sys
def total_charge(value_array): determinants = [np.linalg.det(s) for s in value_array] charge = np.cumsum(np.angle(determinants / np.roll(determinants, 1))) return charge / (2 * np.pi)
def qhe_ribbon(W, periodic=False): """ Creates ribbon system
If we have periodic boundary conditions, the flux through a single unit cell is quantized. """ W = 2 * (W // 2)
def ribbon_shape(pos): (x, y) = pos return (-W / 2 <= y <= W / 2)
def onsite(site, p): (x, y) = site.pos return 4 * p.t - p.mu <http://p.mu>
def hopping(site1, site2, p): x1, y1 = site1.pos x2, y2 = site2.pos return -p.t * np.exp(-0.5j * p.B * (x1 - x2) * (y1 + y2))
def hopping_periodic(site1, site2, p): x1, y1 = site1.pos x2, y2 = site2.pos return -p.t * np.exp(-1j * np.pi / (W + 1) * np.round((W + 1) * p.B / (2 * np.pi)) * (x1 - x2) * (y1 + y2))
lat = kwant.lattice.square() sym_sys = kwant.TranslationalSymmetry((-1, 0)) sys = kwant.Builder(sym_sys)
sys[lat.shape(ribbon_shape, (0, 0))] = onsite
if periodic: sys[lat.neighbors()] = hopping_periodic sys[lat(0, - W / 2), lat(0, + W / 2)] = hopping_periodic else: sys[lat.neighbors()] = hopping
return sys
# Quantum hall bar codes
def qhe_hall_bar(L=50, W=10, w_lead=10, w_vert_lead=None): """Create a hall bar system.
Square lattice, one orbital per site. Returns kwant system.
Arguments required in onsite/hoppings: t, mu, mu_lead """
L = 2 * (L // 2) W = 2 * (W // 2) w_lead = 2 * (w_lead // 2) if w_vert_lead is None: w_vert_lead = w_lead else: w_vert_lead = 2 * (w_vert_lead // 2)
# bar shape def bar(pos): (x, y) = pos return (x >= -L / 2 and x <= L / 2) and (y >= -W / 2 and y <= W / 2)
# Onsite and hoppings def onsite(site, p): (x, y) = site.pos return 4 * p.t - p.mu <http://p.mu>
def hopping_Ax(site1, site2, p): x1, y1 = site1.pos x2, y2 = site2.pos return -p.t * np.exp(-0.5j * p.B * (x1 + x2) * (y1 - y2))
def make_lead_hop_y(x0): def hopping_Ay(site1, site2, p): x1, y1 = site1.pos x2, y2 = site2.pos return -p.t * np.exp(-1j * p.B * x0 * (y1 - y2)) return hopping_Ay
def lead_hop_vert(site1, site2, p): return -p.t
# Building system lat = kwant.lattice.square() sys = kwant.Builder()
sys[lat.shape(bar, (0, 0))] = onsite sys[lat.neighbors()] = hopping_Ax
# Attaching leads sym_lead = kwant.TranslationalSymmetry((-1, 0)) lead = kwant.Builder(sym_lead)
def lead_shape(pos): (x, y) = pos return (-w_lead / 2 <= y <= w_lead / 2)
lead_onsite = lambda site, p: 4 * p.t - p.mu_lead
sym_lead_vertical = kwant.TranslationalSymmetry((0, 1)) lead_vertical1 = kwant.Builder(sym_lead_vertical) lead_vertical2 = kwant.Builder(sym_lead_vertical)
def lead_shape_vertical1(pos): (x, y) = pos return (-L / 4 - w_vert_lead / 2 <= x <= -L / 4 + w_vert_lead /
On 24-08-16 07:02, T.C. Wu wrote: 2)
def lead_shape_vertical2(pos): (x, y) = pos return (+L / 4 - w_vert_lead / 2 <= x <= +L / 4 + w_vert_lead /
2)
lead_vertical1[lat.shape(lead_shape_vertical1, (-L / 4, 0))] = lead_onsite lead_vertical1[lat.neighbors()] = lead_hop_vert lead_vertical2[lat.shape(lead_shape_vertical2, (L / 4, 0))] = lead_onsite lead_vertical2[lat.neighbors()] = lead_hop_vert
sys.attach_lead(lead_vertical1) sys.attach_lead(lead_vertical2)
sys.attach_lead(lead_vertical1.reversed()) sys.attach_lead(lead_vertical2.reversed())
lead[lat.shape(lead_shape, (-1, 0))] = lead_onsite lead[lat.neighbors()] = make_lead_hop_y(-L / 2)
sys.attach_lead(lead)
lead = kwant.Builder(sym_lead) lead[lat.shape(lead_shape, (-1, 0))] = lead_onsite lead[lat.neighbors()] = make_lead_hop_y(L / 2)
sys.attach_lead(lead.reversed())
return sys
def calculate_sigmas(G): # reduce by one dimension G -> G[temp, temp] temp = [0, 1, 3, 4, 5] G = G[temp, :] G = G[:, temp] # invert R = G^-1 # find out whether it is a numpy object r = np.linalg.inv(G) # Voltages follow: V = R I[temp] V = r.dot(np.array([0, 0, 0, 1, -1])) # Completely solved the six terminal system. # Consider the 2x2 conductance now: Use I = sigma U E_x = V[1] - V[0] E_y = V[1] - V[3] # formula above sigma_xx = E_x / (E_x**2 + E_y**2) sigma_xy = E_y / (E_x**2 + E_y**2) return sigma_xx, sigma_xy
def plot_bandstructure(flead, momenta):
bands = kwant.physics.Bands(flead,args = [p]) energies = [bands(k) for k in momenta]
plt.figure() plt.plot(momenta, energies) plt.xlabel(r"$k[a^{-1}]$") plt.ylabel(r"$E [t]$") plt.show()
sys = qhe_hall_bar(L=60, W=100, w_lead=90, w_vert_lead=28).finalized() kwant.plot(sys) p = SimpleNamespace(t=1.0, mu=0.3, mu_lead=0.3, B=0.5) Bs = np.linspace(0.02, 0.15, 200) num_leads = len(sys.leads)
momenta = [-pi + 0.02 * pi * i for i in range(101)] plot_bandstructure(sys.leads[2], momenta)
