Quantum Hall effect in graphene
Dear all, I am trying to use kwant to reproduce the QHE in graphene. I combine the example in qhe.py and graphene.py provided in the official website(see the code below).However, I got the following error when I plot out the wavefunction at energy = 0: RuntimeError: Numbers of left and rightpropagating modes differ, possibly due to a numerical instability. Moreover, if I don't choose energy = 0, say energy = 0.1, then the error disappear. However, I cannot recover the plateaus in the plot of conductance. Lastly, it is correct to change t > t*sigma_0 in order to extract the distribution of spin(i.e. choosing the odd /even elements in the wavefunction by [1::2]/[0::2]) in this case? Thank you. Best, Johnny Code: # Tutorial 2.5. Beyond square lattices: graphene # ============================================== # # Physics background #  # Transport through a graphene quantum dot with a pnjunction # # Kwant features highlighted #  #  Application of all the aspects of tutorials 13 to a more complicated # lattice, namely graphene #execfile('C:/Users/wtc/Google Drive/Linnian/plot_graphene.py') from __future__ import division # so that 1/2 == 0.5, and not 0 from math import pi, sqrt, tanh import random import kwant import numpy as np # For computing eigenvalues import scipy.sparse.linalg as sla # For plotting from matplotlib import pyplot # Define the graphene lattice sin_30, cos_30 = (1 / 2, sqrt(3) / 2) graphene = kwant.lattice.general([(1, 0), (sin_30, cos_30)], [(0, 0), (0, 1 / sqrt(3))]) a, b = graphene.sublattices Bs = np.linspace(4,50,200) t = 1 def make_system(r=10, w=2.0, pot=0.1): #### Define the scattering region. #### # circular scattering region def circle(pos): x, y = pos return x ** 2 + y ** 2 < r ** 2 #return y > (abs(x)*np.tan(np.pi/3) size) and y < size sys = kwant.Builder() # w: width and pot: potential maximum of the pn junction def potential(site): (x, y) = site.pos d = y * cos_30 + x * sin_30 return pot * tanh(d / w) def hopping2(site1, site2, phi = 0,salt=0): xt, yt = site1.pos xs, ys = site2.pos return t * np.exp(1j * phi * (xtxs) *(yt+ys)) sys[graphene.shape(circle, (0, 0))] = 0 # specify the hoppings of the graphene lattice in the # format expected by builder.HoppingKind hoppings = (((0, 0), a, b), ((0, 1), a, b), ((1, 1), a, b)) sys[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2 # Modify the scattering region # del sys[a(0, 0)] # del sys[a(5, 5)] # del sys[a(5,5)] # sys[a(2, 1), b(2, 2)] = 1 #### Define the leads. #### # left lead sym0 = kwant.TranslationalSymmetry(graphene.vec((1, 0))) def lead0_shape(pos): x, y = pos return (0.4 * r < y < 0.4 * r) lead0 = kwant.Builder(sym0) lead0[graphene.shape(lead0_shape, (0, 0))] = 0 lead0[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2 # The second lead, going to the top right sym1 = kwant.TranslationalSymmetry(graphene.vec((0, 1))) def lead1_shape(pos): v = pos[1] * sin_30  pos[0] * cos_30 return (0.4 * r < v < 0.4 * r) lead1 = kwant.Builder(sym1) lead1[graphene.shape(lead1_shape, (0, 0))] = 0 lead1[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2 return sys, [lead0, lead1] def compute_evs(sys): # Compute some eigenvalues of the closed system sparse_mat = sys.hamiltonian_submatrix(sparse=True) evs = sla.eigs(sparse_mat, 2)[0] print evs.real def plot_conductance(sys, energy): # Compute transmission as a function of energy data = [] for phi in Bs: smatrix = kwant.smatrix(sys, energy,args =[phi,""]) data.append(smatrix.transmission(0, 1)) pyplot.figure() pyplot.plot(Bs, data) pyplot.xlabel("energy [t]") pyplot.ylabel("conductance [e^2/h]") pyplot.show() def plot_bandstructure(flead, momenta): bands = kwant.physics.Bands(flead) energies = [bands(k) for k in momenta] pyplot.figure() pyplot.plot(momenta, energies) pyplot.xlabel("momentum [(lattice constant)^1]") pyplot.ylabel("energy [t]") pyplot.show() def density(sys, energy, args, lead_nr): wf = kwant.wave_function(sys, energy, args) return (abs(wf(lead_nr))**2).sum(axis=0) def main(): pot = 0.1 energy = 0.0 sys, leads = make_system(pot=pot) # To highlight the two sublattices of graphene, we plot one with # a filled, and the other one with an open circle: def family_colors(site): return 0 if site.family == a else 1 # Plot the closed system without leads. kwant.plot(sys, site_color=family_colors, site_lw=0.1, colorbar=False) # Compute some eigenvalues. compute_evs(sys.finalized()) # Attach the leads to the system. for lead in leads: sys.attach_lead(lead) # Then, plot the system with leads. kwant.plot(sys, site_color=family_colors, site_lw=0.1, lead_site_lw=0, colorbar=False) # Finalize the system. sys = sys.finalized() # Compute the band structure of lead 0. momenta = [4*pi + 0.02 * 4*pi * i for i in xrange(101)] plot_bandstructure(sys.leads[0], momenta) # plot wf d = density(sys, energy, [1/40.0, ""], 0) kwant.plotter.map(sys, d) reciprocal_phis = np.linspace(4, 50, 200) conductances = [] for phi in 1 / reciprocal_phis: smatrix = kwant.smatrix(sys, energy, args=[phi, ""]) conductances.append(smatrix.transmission(0, 1)) # from lead 1 to lead 0 pyplot.plot(reciprocal_phis, conductances) pyplot.show() # Plot conductance. # energies = [2 * pot + 4. / 50. * pot * i for i in xrange(51)] # plot_conductance(sys, 0) # Call the main function if the script gets executed (as opposed to imported). # See http://docs.python.org/library/__main__.html. if __name__ == '__main__': main()
Dear Johnny,
The error message you're seeing isn't lying. You probably have a lead
with zigzag boundary conditions that has a bad singularity at E=0.
Choosing even E=1e4 would resolve the problem. And yes, you're right
about the spin density (but of course if you don't add any
spindependent terms, there won't be any spin physics happening.
Best regards,
Anton Akhmerov
On Sun, Nov 29, 2015 at 11:28 AM, T.C. Wu
Dear all,
I am trying to use kwant to reproduce the QHE in graphene. I combine the example in qhe.py and graphene.py provided in the official website(see the code below).However, I got the following error when I plot out the wavefunction at energy = 0: RuntimeError: Numbers of left and rightpropagating modes differ, possibly due to a numerical instability.
Moreover, if I don't choose energy = 0, say energy = 0.1, then the error disappear. However, I cannot recover the plateaus in the plot of conductance.
Lastly, it is correct to change t > t*sigma_0 in order to extract the distribution of spin(i.e. choosing the odd /even elements in the wavefunction by [1::2]/[0::2]) in this case?
Thank you.
Best, Johnny
Code: # Tutorial 2.5. Beyond square lattices: graphene # ============================================== # # Physics background #  # Transport through a graphene quantum dot with a pnjunction # # Kwant features highlighted #  #  Application of all the aspects of tutorials 13 to a more complicated # lattice, namely graphene #execfile('C:/Users/wtc/Google Drive/Linnian/plot_graphene.py') from __future__ import division # so that 1/2 == 0.5, and not 0 from math import pi, sqrt, tanh import random import kwant import numpy as np # For computing eigenvalues import scipy.sparse.linalg as sla
# For plotting from matplotlib import pyplot
# Define the graphene lattice sin_30, cos_30 = (1 / 2, sqrt(3) / 2) graphene = kwant.lattice.general([(1, 0), (sin_30, cos_30)], [(0, 0), (0, 1 / sqrt(3))]) a, b = graphene.sublattices Bs = np.linspace(4,50,200) t = 1
def make_system(r=10, w=2.0, pot=0.1):
#### Define the scattering region. #### # circular scattering region def circle(pos): x, y = pos return x ** 2 + y ** 2 < r ** 2 #return y > (abs(x)*np.tan(np.pi/3) size) and y < size sys = kwant.Builder()
# w: width and pot: potential maximum of the pn junction def potential(site): (x, y) = site.pos d = y * cos_30 + x * sin_30 return pot * tanh(d / w)
def hopping2(site1, site2, phi = 0,salt=0): xt, yt = site1.pos xs, ys = site2.pos return t * np.exp(1j * phi * (xtxs) *(yt+ys))
sys[graphene.shape(circle, (0, 0))] = 0
# specify the hoppings of the graphene lattice in the # format expected by builder.HoppingKind hoppings = (((0, 0), a, b), ((0, 1), a, b), ((1, 1), a, b)) sys[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2
# Modify the scattering region # del sys[a(0, 0)] # del sys[a(5, 5)] # del sys[a(5,5)] # sys[a(2, 1), b(2, 2)] = 1
#### Define the leads. #### # left lead sym0 = kwant.TranslationalSymmetry(graphene.vec((1, 0)))
def lead0_shape(pos): x, y = pos return (0.4 * r < y < 0.4 * r)
lead0 = kwant.Builder(sym0) lead0[graphene.shape(lead0_shape, (0, 0))] = 0 lead0[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2
# The second lead, going to the top right sym1 = kwant.TranslationalSymmetry(graphene.vec((0, 1)))
def lead1_shape(pos): v = pos[1] * sin_30  pos[0] * cos_30 return (0.4 * r < v < 0.4 * r)
lead1 = kwant.Builder(sym1) lead1[graphene.shape(lead1_shape, (0, 0))] = 0 lead1[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2
return sys, [lead0, lead1]
def compute_evs(sys): # Compute some eigenvalues of the closed system sparse_mat = sys.hamiltonian_submatrix(sparse=True)
evs = sla.eigs(sparse_mat, 2)[0] print evs.real
def plot_conductance(sys, energy): # Compute transmission as a function of energy data = [] for phi in Bs: smatrix = kwant.smatrix(sys, energy,args =[phi,""]) data.append(smatrix.transmission(0, 1))
pyplot.figure() pyplot.plot(Bs, data) pyplot.xlabel("energy [t]") pyplot.ylabel("conductance [e^2/h]") pyplot.show()
def plot_bandstructure(flead, momenta): bands = kwant.physics.Bands(flead) energies = [bands(k) for k in momenta]
pyplot.figure() pyplot.plot(momenta, energies) pyplot.xlabel("momentum [(lattice constant)^1]") pyplot.ylabel("energy [t]") pyplot.show() def density(sys, energy, args, lead_nr): wf = kwant.wave_function(sys, energy, args) return (abs(wf(lead_nr))**2).sum(axis=0)
def main(): pot = 0.1 energy = 0.0 sys, leads = make_system(pot=pot)
# To highlight the two sublattices of graphene, we plot one with # a filled, and the other one with an open circle: def family_colors(site): return 0 if site.family == a else 1
# Plot the closed system without leads. kwant.plot(sys, site_color=family_colors, site_lw=0.1, colorbar=False)
# Compute some eigenvalues. compute_evs(sys.finalized())
# Attach the leads to the system. for lead in leads: sys.attach_lead(lead)
# Then, plot the system with leads. kwant.plot(sys, site_color=family_colors, site_lw=0.1, lead_site_lw=0, colorbar=False)
# Finalize the system. sys = sys.finalized()
# Compute the band structure of lead 0. momenta = [4*pi + 0.02 * 4*pi * i for i in xrange(101)] plot_bandstructure(sys.leads[0], momenta)
# plot wf d = density(sys, energy, [1/40.0, ""], 0) kwant.plotter.map(sys, d) reciprocal_phis = np.linspace(4, 50, 200) conductances = [] for phi in 1 / reciprocal_phis: smatrix = kwant.smatrix(sys, energy, args=[phi, ""]) conductances.append(smatrix.transmission(0, 1)) # from lead 1 to lead 0 pyplot.plot(reciprocal_phis, conductances) pyplot.show() # Plot conductance. # energies = [2 * pot + 4. / 50. * pot * i for i in xrange(51)] # plot_conductance(sys, 0)
# Call the main function if the script gets executed (as opposed to imported). # See http://docs.python.org/library/__main__.html. if __name__ == '__main__': main()
Dear Anton,
Thanks for your prompt reply. However, how should I amend the code (which
works for square lattice in the example) in order to recover the
conductance plateaus? Thank you again.
Best,
Johnny
On 29 November 2015 at 18:35, Anton Akhmerov
Dear Johnny,
The error message you're seeing isn't lying. You probably have a lead with zigzag boundary conditions that has a bad singularity at E=0. Choosing even E=1e4 would resolve the problem. And yes, you're right about the spin density (but of course if you don't add any spindependent terms, there won't be any spin physics happening.
Best regards, Anton Akhmerov
Dear all,
I am trying to use kwant to reproduce the QHE in graphene. I combine the example in qhe.py and graphene.py provided in the official website(see
code below).However, I got the following error when I plot out the wavefunction at energy = 0: RuntimeError: Numbers of left and rightpropagating modes differ,
On Sun, Nov 29, 2015 at 11:28 AM, T.C. Wu
wrote: the possibly due to a numerical instability.
Moreover, if I don't choose energy = 0, say energy = 0.1, then the error disappear. However, I cannot recover the plateaus in the plot of conductance.
Lastly, it is correct to change t > t*sigma_0 in order to extract the distribution of spin(i.e. choosing the odd /even elements in the wavefunction by [1::2]/[0::2]) in this case?
Thank you.
Best, Johnny
Code: # Tutorial 2.5. Beyond square lattices: graphene # ============================================== # # Physics background #  # Transport through a graphene quantum dot with a pnjunction # # Kwant features highlighted #  #  Application of all the aspects of tutorials 13 to a more complicated # lattice, namely graphene #execfile('C:/Users/wtc/Google Drive/Linnian/plot_graphene.py') from __future__ import division # so that 1/2 == 0.5, and not 0 from math import pi, sqrt, tanh import random import kwant import numpy as np # For computing eigenvalues import scipy.sparse.linalg as sla
# For plotting from matplotlib import pyplot
# Define the graphene lattice sin_30, cos_30 = (1 / 2, sqrt(3) / 2) graphene = kwant.lattice.general([(1, 0), (sin_30, cos_30)], [(0, 0), (0, 1 / sqrt(3))]) a, b = graphene.sublattices Bs = np.linspace(4,50,200) t = 1
def make_system(r=10, w=2.0, pot=0.1):
#### Define the scattering region. #### # circular scattering region def circle(pos): x, y = pos return x ** 2 + y ** 2 < r ** 2 #return y > (abs(x)*np.tan(np.pi/3) size) and y < size sys = kwant.Builder()
# w: width and pot: potential maximum of the pn junction def potential(site): (x, y) = site.pos d = y * cos_30 + x * sin_30 return pot * tanh(d / w)
def hopping2(site1, site2, phi = 0,salt=0): xt, yt = site1.pos xs, ys = site2.pos return t * np.exp(1j * phi * (xtxs) *(yt+ys))
sys[graphene.shape(circle, (0, 0))] = 0
# specify the hoppings of the graphene lattice in the # format expected by builder.HoppingKind hoppings = (((0, 0), a, b), ((0, 1), a, b), ((1, 1), a, b)) sys[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2
# Modify the scattering region # del sys[a(0, 0)] # del sys[a(5, 5)] # del sys[a(5,5)] # sys[a(2, 1), b(2, 2)] = 1
#### Define the leads. #### # left lead sym0 = kwant.TranslationalSymmetry(graphene.vec((1, 0)))
def lead0_shape(pos): x, y = pos return (0.4 * r < y < 0.4 * r)
lead0 = kwant.Builder(sym0) lead0[graphene.shape(lead0_shape, (0, 0))] = 0 lead0[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2
# The second lead, going to the top right sym1 = kwant.TranslationalSymmetry(graphene.vec((0, 1)))
def lead1_shape(pos): v = pos[1] * sin_30  pos[0] * cos_30 return (0.4 * r < v < 0.4 * r)
lead1 = kwant.Builder(sym1) lead1[graphene.shape(lead1_shape, (0, 0))] = 0 lead1[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2
return sys, [lead0, lead1]
def compute_evs(sys): # Compute some eigenvalues of the closed system sparse_mat = sys.hamiltonian_submatrix(sparse=True)
evs = sla.eigs(sparse_mat, 2)[0] print evs.real
def plot_conductance(sys, energy): # Compute transmission as a function of energy data = [] for phi in Bs: smatrix = kwant.smatrix(sys, energy,args =[phi,""]) data.append(smatrix.transmission(0, 1))
pyplot.figure() pyplot.plot(Bs, data) pyplot.xlabel("energy [t]") pyplot.ylabel("conductance [e^2/h]") pyplot.show()
def plot_bandstructure(flead, momenta): bands = kwant.physics.Bands(flead) energies = [bands(k) for k in momenta]
pyplot.figure() pyplot.plot(momenta, energies) pyplot.xlabel("momentum [(lattice constant)^1]") pyplot.ylabel("energy [t]") pyplot.show() def density(sys, energy, args, lead_nr): wf = kwant.wave_function(sys, energy, args) return (abs(wf(lead_nr))**2).sum(axis=0)
def main(): pot = 0.1 energy = 0.0 sys, leads = make_system(pot=pot)
# To highlight the two sublattices of graphene, we plot one with # a filled, and the other one with an open circle: def family_colors(site): return 0 if site.family == a else 1
# Plot the closed system without leads. kwant.plot(sys, site_color=family_colors, site_lw=0.1,
colorbar=False) > > # Compute some eigenvalues. > compute_evs(sys.finalized()) > > # Attach the leads to the system. > for lead in leads: > sys.attach_lead(lead) > > # Then, plot the system with leads. > kwant.plot(sys, site_color=family_colors, site_lw=0.1, > lead_site_lw=0, colorbar=False) > > # Finalize the system. > sys = sys.finalized() > > # Compute the band structure of lead 0. > momenta = [4*pi + 0.02 * 4*pi * i for i in xrange(101)] > plot_bandstructure(sys.leads[0], momenta) > > > # plot wf > d = density(sys, energy, [1/40.0, ""], 0) > kwant.plotter.map(sys, d) > reciprocal_phis = np.linspace(4, 50, 200) > conductances = [] > for phi in 1 / reciprocal_phis: > smatrix = kwant.smatrix(sys, energy, args=[phi, ""]) > conductances.append(smatrix.transmission(0, 1)) # from lead 1 to > lead 0 > pyplot.plot(reciprocal_phis, conductances) > pyplot.show() > # Plot conductance. > # energies = [2 * pot + 4. / 50. * pot * i for i in xrange(51)] > # plot_conductance(sys, 0) > > > # Call the main function if the script gets executed (as opposed to > imported). > # See <http://docs.python.org/library/__main__.html>. > if __name__ == '__main__': > main() >
Dear all,
For honeycomb lattice, is it still correct to use the hopping t *
np.exp(1j * phi * (xtxs) *(yt+ys)) in the presence of a magnetic field? I
am just wondering the vaildity of it because I cannot recover the
conductance steps in QHE.
Any suggestions will be appreciated. Thanks.
Best,
Johnny
On 29 November 2015 at 19:13, T.C. Wu
Dear Anton,
Thanks for your prompt reply. However, how should I amend the code (which works for square lattice in the example) in order to recover the conductance plateaus? Thank you again.
Best, Johnny
On 29 November 2015 at 18:35, Anton Akhmerov
wrote: Dear Johnny,
The error message you're seeing isn't lying. You probably have a lead with zigzag boundary conditions that has a bad singularity at E=0. Choosing even E=1e4 would resolve the problem. And yes, you're right about the spin density (but of course if you don't add any spindependent terms, there won't be any spin physics happening.
Best regards, Anton Akhmerov
Dear all,
I am trying to use kwant to reproduce the QHE in graphene. I combine the example in qhe.py and graphene.py provided in the official website(see
code below).However, I got the following error when I plot out the wavefunction at energy = 0: RuntimeError: Numbers of left and rightpropagating modes differ,
On Sun, Nov 29, 2015 at 11:28 AM, T.C. Wu
wrote: the possibly due to a numerical instability.
Moreover, if I don't choose energy = 0, say energy = 0.1, then the error disappear. However, I cannot recover the plateaus in the plot of conductance.
Lastly, it is correct to change t > t*sigma_0 in order to extract the distribution of spin(i.e. choosing the odd /even elements in the wavefunction by [1::2]/[0::2]) in this case?
Thank you.
Best, Johnny
Code: # Tutorial 2.5. Beyond square lattices: graphene # ============================================== # # Physics background #  # Transport through a graphene quantum dot with a pnjunction # # Kwant features highlighted #  #  Application of all the aspects of tutorials 13 to a more complicated # lattice, namely graphene #execfile('C:/Users/wtc/Google Drive/Linnian/plot_graphene.py') from __future__ import division # so that 1/2 == 0.5, and not 0 from math import pi, sqrt, tanh import random import kwant import numpy as np # For computing eigenvalues import scipy.sparse.linalg as sla
# For plotting from matplotlib import pyplot
# Define the graphene lattice sin_30, cos_30 = (1 / 2, sqrt(3) / 2) graphene = kwant.lattice.general([(1, 0), (sin_30, cos_30)], [(0, 0), (0, 1 / sqrt(3))]) a, b = graphene.sublattices Bs = np.linspace(4,50,200) t = 1
def make_system(r=10, w=2.0, pot=0.1):
#### Define the scattering region. #### # circular scattering region def circle(pos): x, y = pos return x ** 2 + y ** 2 < r ** 2 #return y > (abs(x)*np.tan(np.pi/3) size) and y < size sys = kwant.Builder()
# w: width and pot: potential maximum of the pn junction def potential(site): (x, y) = site.pos d = y * cos_30 + x * sin_30 return pot * tanh(d / w)
def hopping2(site1, site2, phi = 0,salt=0): xt, yt = site1.pos xs, ys = site2.pos return t * np.exp(1j * phi * (xtxs) *(yt+ys))
sys[graphene.shape(circle, (0, 0))] = 0
# specify the hoppings of the graphene lattice in the # format expected by builder.HoppingKind hoppings = (((0, 0), a, b), ((0, 1), a, b), ((1, 1), a, b)) sys[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2
# Modify the scattering region # del sys[a(0, 0)] # del sys[a(5, 5)] # del sys[a(5,5)] # sys[a(2, 1), b(2, 2)] = 1
#### Define the leads. #### # left lead sym0 = kwant.TranslationalSymmetry(graphene.vec((1, 0)))
def lead0_shape(pos): x, y = pos return (0.4 * r < y < 0.4 * r)
lead0 = kwant.Builder(sym0) lead0[graphene.shape(lead0_shape, (0, 0))] = 0 lead0[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2
# The second lead, going to the top right sym1 = kwant.TranslationalSymmetry(graphene.vec((0, 1)))
def lead1_shape(pos): v = pos[1] * sin_30  pos[0] * cos_30 return (0.4 * r < v < 0.4 * r)
lead1 = kwant.Builder(sym1) lead1[graphene.shape(lead1_shape, (0, 0))] = 0 lead1[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = hopping2
return sys, [lead0, lead1]
def compute_evs(sys): # Compute some eigenvalues of the closed system sparse_mat = sys.hamiltonian_submatrix(sparse=True)
evs = sla.eigs(sparse_mat, 2)[0] print evs.real
def plot_conductance(sys, energy): # Compute transmission as a function of energy data = [] for phi in Bs: smatrix = kwant.smatrix(sys, energy,args =[phi,""]) data.append(smatrix.transmission(0, 1))
pyplot.figure() pyplot.plot(Bs, data) pyplot.xlabel("energy [t]") pyplot.ylabel("conductance [e^2/h]") pyplot.show()
def plot_bandstructure(flead, momenta): bands = kwant.physics.Bands(flead) energies = [bands(k) for k in momenta]
pyplot.figure() pyplot.plot(momenta, energies) pyplot.xlabel("momentum [(lattice constant)^1]") pyplot.ylabel("energy [t]") pyplot.show() def density(sys, energy, args, lead_nr): wf = kwant.wave_function(sys, energy, args) return (abs(wf(lead_nr))**2).sum(axis=0)
def main(): pot = 0.1 energy = 0.0 sys, leads = make_system(pot=pot)
# To highlight the two sublattices of graphene, we plot one with # a filled, and the other one with an open circle: def family_colors(site): return 0 if site.family == a else 1
# Plot the closed system without leads. kwant.plot(sys, site_color=family_colors, site_lw=0.1,
colorbar=False) > > # Compute some eigenvalues. > compute_evs(sys.finalized()) > > # Attach the leads to the system. > for lead in leads: > sys.attach_lead(lead) > > # Then, plot the system with leads. > kwant.plot(sys, site_color=family_colors, site_lw=0.1, > lead_site_lw=0, colorbar=False) > > # Finalize the system. > sys = sys.finalized() > > # Compute the band structure of lead 0. > momenta = [4*pi + 0.02 * 4*pi * i for i in xrange(101)] > plot_bandstructure(sys.leads[0], momenta) > > > # plot wf > d = density(sys, energy, [1/40.0, ""], 0) > kwant.plotter.map(sys, d) > reciprocal_phis = np.linspace(4, 50, 200) > conductances = [] > for phi in 1 / reciprocal_phis: > smatrix = kwant.smatrix(sys, energy, args=[phi, ""]) > conductances.append(smatrix.transmission(0, 1)) # from lead 1 to > lead 0 > pyplot.plot(reciprocal_phis, conductances) > pyplot.show() > # Plot conductance. > # energies = [2 * pot + 4. / 50. * pot * i for i in xrange(51)] > # plot_conductance(sys, 0) > > > # Call the main function if the script gets executed (as opposed to > imported). > # See <http://docs.python.org/library/__main__.html>. > if __name__ == '__main__': > main() >
Christoph Groth wrote:
T.C. Wu wrote:
(...) because I cannot recover the conductance steps in QHE.
Please check the QHE example in the Kwant paper. That example calculates conductance steps.
I attach a honeycomblattice version of qhe.py that does show conductance steps.
participants (3)

Anton Akhmerov

Christoph Groth

T.C. Wu