From 37266338616bcc325eadae666fcc9e0189cd5a78 Mon Sep 17 00:00:00 2001 From: Christoph Groth <christoph.groth@cea.fr> Date: Thu, 5 Apr 2012 18:21:10 +0200 Subject: [PATCH] doc: replace the figure creation scripts ...by ones closely based on the actual tutorial examples. --- doc/source/images/tutorial2a.py | 183 ++++++++-------- doc/source/images/tutorial2b.py | 200 +++++++++-------- doc/source/images/tutorial2c.py | 365 +++++++++++++++----------------- doc/source/images/tutorial3a.py | 86 +++++--- doc/source/images/tutorial3b.py | 124 ++++++----- 5 files changed, 467 insertions(+), 491 deletions(-) diff --git a/doc/source/images/tutorial2a.py b/doc/source/images/tutorial2a.py index c078ae6e..942b6bdf 100644 --- a/doc/source/images/tutorial2a.py +++ b/doc/source/images/tutorial2a.py @@ -11,113 +11,96 @@ # - Numpy matrices as values in Builder import kwant + +# For plotting +import pylab + +# For matrix support import numpy import latex, html -# define sigma-matrices for convenience +# define Pauli-matrices for convenience sigma_0 = numpy.eye(2) sigma_x = numpy.array([[0, 1], [1, 0]]) sigma_y = numpy.array([[0, -1j], [1j, 0]]) sigma_z = numpy.array([[1, 0], [0, -1]]) -# First, define the tight-binding system - -sys = kwant.Builder() - -# Here, we are only working with square lattices - -# for simplicity, take lattice constant = 1 -a = 1 -lat = kwant.lattice.Square(a) - -t = 1.0 -alpha = 0.5 -e_z = 0.08 -W = 10 -L = 30 - -# Define the scattering region - -def rectangle(pos): - (x, y) = pos - return ( -0.5 < x < L - 0.5 ) and ( -0.5 < y < W - 0.5 ) - -sys[lat.shape(rectangle, (0, 0))] = 4 * t * sigma_0 + e_z * sigma_z -# hoppings in x-direction -sys[sys.possible_hoppings((1, 0), lat, lat)] = - t * sigma_0 - \ - 1j * alpha * sigma_y -# hoppings in y-directions -sys[sys.possible_hoppings((0, 1), lat, lat)] = - t * sigma_0 + \ - 1j * alpha * sigma_x - -# Then, define the leads: - -# First the lead to the left - -# (Note: in the current version, TranslationalSymmetry takes a -# realspace vector) -sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) -lead0 = kwant.Builder(sym_lead0) -lead0.default_site_group = lat - -def lead_shape(pos): - (x, y) = pos - return (-1 < x < 1) and ( -0.5 < y < W - 0.5 ) - -lead0[lat.shape(lead_shape, (0, 0))] = 4 * t * sigma_0 + e_z * sigma_z -# hoppings in x-direction -lead0[lead0.possible_hoppings((1, 0), lat, lat)] = - t * sigma_0 - \ - 1j * alpha * sigma_y -# hoppings in y-directions -lead0[lead0.possible_hoppings((0, 1), lat, lat)] = - t * sigma_0 + \ - 1j * alpha * sigma_x - -# Then the lead to the right -# there we can use a special function that simply reverses the direction - -lead1 = lead0.reversed() - -# Then attach the leads to the system - -sys.attach_lead(lead0) -sys.attach_lead(lead1) - -# finalize the system - -fsys = sys.finalized() - -# Now that we have the system, we can compute conductance - -energies = [] -data = [] -for ie in xrange(100): - energy = ie * 0.01 - 0.3 - - # compute the scattering matrix at energy energy - smatrix = kwant.solvers.sparse.solve(fsys, energy) - - # compute the transmission probability from lead 0 to - # lead 1 - energies.append(energy) - data.append(smatrix.transmission(1, 0)) - -# Use matplotlib to write output -# We should see conductance steps -import pylab -pylab.plot(energies, data) -pylab.xlabel("energy [in units of t]", - fontsize=latex.mpl_label_size) -pylab.ylabel("conductance [in units of e^2/h]", - fontsize=latex.mpl_label_size) -fig = pylab.gcf() -pylab.setp(fig.get_axes()[0].get_xticklabels(), - fontsize=latex.mpl_tick_size) -pylab.setp(fig.get_axes()[0].get_yticklabels(), - fontsize=latex.mpl_tick_size) -fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) -fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) -fig.savefig("tutorial2a_result.pdf") -fig.savefig("tutorial2a_result.png", - dpi=(html.figwidth_px/latex.mpl_width_in)) +def make_system(a=1, t=1.0, alpha=0.5, e_z=0.08, W=10, L=30): + # Start with an empty tight-binding system and a single square lattice. + # `a` is the lattice constant (by default set to 1 for simplicity). + lat = kwant.lattice.Square(a) + + sys = kwant.Builder() + sys.default_site_group = lat + + #### Define the scattering region. #### + sys[((x, y) for x in range(L) for y in range(W))] = 4 * t * sigma_0 + \ + e_z * sigma_z + # hoppings in x-direction + sys[sys.possible_hoppings((1, 0), lat, lat)] = - t * sigma_0 - \ + 1j * alpha * sigma_y + # hoppings in y-directions + sys[sys.possible_hoppings((0, 1), lat, lat)] = - t * sigma_0 + \ + 1j * alpha * sigma_x + + #### Define the leads. #### + # left lead + sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) + lead0 = kwant.Builder(sym_lead0) + lead0.default_site_group = lat + + lead0[((0, j) for j in xrange(W))] = 4 * t * sigma_0 + e_z * sigma_z + # hoppings in x-direction + lead0[lead0.possible_hoppings((1, 0), lat, lat)] = - t * sigma_0 - \ + 1j * alpha * sigma_y + # hoppings in y-directions + lead0[lead0.possible_hoppings((0, 1), lat, lat)] = - t * sigma_0 + \ + 1j * alpha * sigma_x + + # Then the lead to the right + # (again, obtained using reverse() + lead1 = lead0.reversed() + + #### Attach the leads and return the finalized system. #### + sys.attach_lead(lead0) + sys.attach_lead(lead1) + + return sys.finalized() + +def plot_conductance(fsys, energies): + # Compute conductance + data = [] + for energy in energies: + smatrix = kwant.solve(fsys, energy) + data.append(smatrix.transmission(1, 0)) + + pylab.plot(energies, data) + pylab.xlabel("energy [in units of t]", + fontsize=latex.mpl_label_size) + pylab.ylabel("conductance [in units of e^2/h]", + fontsize=latex.mpl_label_size) + fig = pylab.gcf() + pylab.setp(fig.get_axes()[0].get_xticklabels(), + fontsize=latex.mpl_tick_size) + pylab.setp(fig.get_axes()[0].get_yticklabels(), + fontsize=latex.mpl_tick_size) + fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) + fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) + fig.savefig("tutorial2a_result.pdf") + fig.savefig("tutorial2a_result.png", + dpi=(html.figwidth_px/latex.mpl_width_in)) + + +def main(): + fsys = make_system() + + # We should see non-monotonic conductance steps. + plot_conductance(fsys, energies=[0.01 * i - 0.3 for i in xrange(100)]) + + +# 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() diff --git a/doc/source/images/tutorial2b.py b/doc/source/images/tutorial2b.py index 6e00c7b7..de4a8c68 100644 --- a/doc/source/images/tutorial2b.py +++ b/doc/source/images/tutorial2b.py @@ -8,108 +8,102 @@ import kwant -import latex, html - -# First, define the tight-binding system - -sys = kwant.Builder() - -# Here, we are only working with square lattices - -# for simplicity, take lattice constant = 1 -a = 1 -lat = kwant.lattice.Square(a) - -t = 1.0 -alpha = 0.5 -e_z = 0.08 -W = 10 -L = 30 - -# Define the scattering region - -def rectangle(pos): - (x, y) = pos - return ( -0.5 < x < L - 0.5 ) and ( -0.5 < y < W - 0.5 ) - -def potential(site): - (x, y) = site.pos - if 10 < x < 20: - return pot - else: - return 0 - -def onsite(site): - return 4 * t + potential(site) - -sys[lat.shape(rectangle, (0, 0))] = onsite -for hopping in lat.nearest: - sys[sys.possible_hoppings(*hopping)] = - t - -# Then, define the leads: - -# First the lead to the left - -# (Note: in the current version, TranslationalSymmetry takes a -# realspace vector) -sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) -lead0 = kwant.Builder(sym_lead0) -lead0.default_site_group = lat - -def lead_shape(pos): - (x, y) = pos - return (-1 < x < 1) and ( -0.5 < y < W - 0.5 ) - -lead0[lat.shape(lead_shape, (0, 0))] = 4 * t -for hopping in lat.nearest: - lead0[lead0.possible_hoppings(*hopping)] = - t - -# Then the lead to the right -# there we can use a special function that simply reverses the direction - -lead1 = lead0.reversed() - -# Then attach the leads to the system - -sys.attach_lead(lead0) -sys.attach_lead(lead1) - -# finalize the system - -fsys = sys.finalized() - -# Now that we have the system, we can compute conductance - -energy = 0.2 -wellpot = [] -data = [] -for ipot in xrange(100): - pot = - ipot * 0.01 - - # compute the scattering matrix at energy energy - smatrix = kwant.solvers.sparse.solve(fsys, energy) - - # compute the transmission probability from lead 0 to - # lead 1 - wellpot.append(-pot) - data.append(smatrix.transmission(1, 0)) - -# Use matplotlib to write output -# We should see conductance steps +# For plotting import pylab -pylab.plot(wellpot, data) -pylab.xlabel("well depth [in units of t]", - fontsize=latex.mpl_label_size) -pylab.ylabel("conductance [in units of e^2/h]", - fontsize=latex.mpl_label_size) -fig = pylab.gcf() -pylab.setp(fig.get_axes()[0].get_xticklabels(), - fontsize=latex.mpl_tick_size) -pylab.setp(fig.get_axes()[0].get_yticklabels(), - fontsize=latex.mpl_tick_size) -fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) -fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) -fig.savefig("tutorial2b_result.pdf") -fig.savefig("tutorial2b_result.png", - dpi=(html.figwidth_px/latex.mpl_width_in)) +import latex, html + +# global variable governing the behavior of potential() in +# make_system() +pot = 0 + +def make_system(a=1, t=1.0, W=10, L=30, L_well=10): + # Start with an empty tight-binding system and a single square lattice. + # `a` is the lattice constant (by default set to 1 for simplicity. + lat = kwant.lattice.Square(a) + + sys = kwant.Builder() + sys.default_site_group = lat + + #### Define the scattering region. #### + # Potential profile + def potential(site): + (x, y) = site.pos + if (L - L_well) / 2 < x < (L + L_well) / 2: + # The potential value is provided using a global variable + return pot + else: + return 0 + + def onsite(site): + return 4 * t + potential(site) + + sys[((x, y) for x in range(L) for y in range(W))] = onsite + for hopping in lat.nearest: + sys[sys.possible_hoppings(*hopping)] = -t + + #### Define the leads. #### + # First the lead to the left, ... + # (Note: in the current version, TranslationalSymmetry takes a + # realspace vector) + sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) + lead0 = kwant.Builder(sym_lead0) + lead0.default_site_group = lat + + lead0[((0, j) for j in xrange(W))] = 4 * t + for hopping in lat.nearest: + lead0[lead0.possible_hoppings(*hopping)] = - t + + # ... then the lead to the right. We use a method that returns a copy of + # `lead0` with its direction reversed. + lead1 = lead0.reversed() + + #### Attach the leads and return the finalized system. #### + sys.attach_lead(lead0) + sys.attach_lead(lead1) + + return sys.finalized() + +def plot_conductance(fsys, energy, welldepths): + # We specify that we want to not only read, but also write to a + # global variable. + global pot + + # Compute conductance + data = [] + for welldepth in welldepths: + # Set the global variable that defines the potential well depth + pot = -welldepth + + smatrix = kwant.solve(fsys, energy) + data.append(smatrix.transmission(1, 0)) + + pylab.plot(welldepths, data) + pylab.xlabel("well depth [in units of t]", + fontsize=latex.mpl_label_size) + pylab.ylabel("conductance [in units of e^2/h]", + fontsize=latex.mpl_label_size) + fig = pylab.gcf() + pylab.setp(fig.get_axes()[0].get_xticklabels(), + fontsize=latex.mpl_tick_size) + pylab.setp(fig.get_axes()[0].get_yticklabels(), + fontsize=latex.mpl_tick_size) + fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) + fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) + fig.savefig("tutorial2b_result.pdf") + fig.savefig("tutorial2b_result.png", + dpi=(html.figwidth_px/latex.mpl_width_in)) + + +def main(): + fsys = make_system() + + # We should see conductance steps. + plot_conductance(fsys, energy=0.2, + welldepths=[0.01 * i for i in xrange(100)]) + + +# 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() diff --git a/doc/source/images/tutorial2c.py b/doc/source/images/tutorial2c.py index 2c1ea86a..2188d3d4 100644 --- a/doc/source/images/tutorial2c.py +++ b/doc/source/images/tutorial2c.py @@ -11,204 +11,177 @@ from cmath import exp from math import pi -import kwant - -import latex, html - -# First, define the tight-binding system - -sys = kwant.Builder() - -# Here, we are only working with square lattices - -# for simplicity, take lattice constant = 1 -a = 1 -lat = kwant.lattice.Square(a) - -t = 1.0 -W = 10 -r1 = 10 -r2 = 20 - -# Define the scattering region -# Now, we aim for a more compelx shape, namely a ring (or annulus) - -def ring(pos): - (x, y) = pos - rsq = x**2 + y**2 - return ( r1**2 < rsq < r2**2) - -sys[lat.shape(ring, (0, 11))] = 4 * t -for hopping in lat.nearest: - sys[sys.possible_hoppings(*hopping)] = - t - -# In order to introduce a flux through the ring, we introduce a phase -# on the hoppings on the line cut through one of the arms - -# since we want to change the flux without modifying Builder repeatedly, -# we define the modified hoppings as a function that takes the flux -# as a global variable. - -def fluxphase(site1, site2): - return exp(1j * phi) - -# Now go through all the hoppings and modify those in the lower -# arm of the ring that go from x=0 to x=1 - -for (site1, site2) in sys.hoppings(): - ix1, iy1 = site1.tag - ix2, iy2 = site2.tag - - hopx = tuple(sorted((ix1, ix2))) - - if hopx == (0, 1) and iy1 == iy2 and iy1 < 0: - sys[lat(hopx[1], iy1), lat(hopx[0], iy1)] = fluxphase - -# Then, define the leads: - -# First the lead to the left - -# (Note: in the current version, TranslationalSymmetry takes a -# realspace vector) -sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) -lead0 = kwant.Builder(sym_lead0) -lead0.default_site_group = lat - -def lead_shape(pos): - (x, y) = pos - return (-1 < x < 1) and ( -W/2 < y < W/2 ) - -lead0[lat.shape(lead_shape, (0, 0))] = 4 * t -for hopping in lat.nearest: - lead0[lead0.possible_hoppings(*hopping)] = - t - -# Then the lead to the right -# there we can use a special function that simply reverses the direction - -lead1 = lead0.reversed() - -# Then attach the leads to the system - -sys.attach_lead(lead0) -sys.attach_lead(lead1) - -# finalize the system - -fsys = sys.finalized() - -# and plot it, to make sure it's proper - -kwant.plot(fsys, "tutorial2c_sys.pdf", width=latex.figwidth_pt) -kwant.plot(fsys, "tutorial2c_sys.png", width=html.figwidth_px) -# Now that we have the system, we can compute conductance - -energy = 0.15 -phases = [] -data = [] -for iphi in xrange(100): - phi = iphi * 0.01 * 3 * 2 * pi - - # compute the scattering matrix at energy energy - smatrix = kwant.solve(fsys, energy) - - # compute the transmission probability from lead 0 to - # lead 1 - phases.append(phi / (2 * pi)) - data.append(smatrix.transmission(1, 0)) +import kwant -# Use matplotlib to write output -# We should see conductance steps +# For plotting import pylab -pylab.plot(phases, data) -pylab.xlabel("flux [in units of the flux quantum]", - fontsize=latex.mpl_label_size) -pylab.ylabel("conductance [in units of e^2/h]", - fontsize=latex.mpl_label_size) -fig = pylab.gcf() -pylab.setp(fig.get_axes()[0].get_xticklabels(), - fontsize=latex.mpl_tick_size) -pylab.setp(fig.get_axes()[0].get_yticklabels(), - fontsize=latex.mpl_tick_size) -fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) -fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) -fig.savefig("tutorial2c_result.pdf") -fig.savefig("tutorial2c_result.png", - dpi=(html.figwidth_px/latex.mpl_width_in)) - -# Finally, some plots needed for the notes - -sys = kwant.Builder() - -sys[lat.shape(ring, (0, 11))] = 4 * t -for hopping in lat.nearest: - sys[sys.possible_hoppings(*hopping)] = - t - -sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) -lead0 = kwant.Builder(sym_lead0) -lead0.default_site_group = lat - -def lead_shape(pos): - (x, y) = pos - return (-1 < x < 1) and ( 0.5 * W < y < 1.5 * W ) - -lead0[lat.shape(lead_shape, (0, W))] = 4 * t -for hopping in lat.nearest: - lead0[lead0.possible_hoppings(*hopping)] = - t - -# Then the lead to the right -# there we can use a special function that simply reverses the direction - -lead1 = lead0.reversed() - -# Then attach the leads to the system - -sys.attach_lead(lead0) -sys.attach_lead(lead1) - -# finalize the system - -fsys = sys.finalized() - -# and plot it, to make sure it's proper - -kwant.plot(fsys, "tutorial2c_note1.pdf", width=latex.figwidth_small_pt) -kwant.plot(fsys, "tutorial2c_note1.png", width=html.figwidth_small_px) - -sys = kwant.Builder() - -sys[lat.shape(ring, (0, 11))] = 4 * t -for hopping in lat.nearest: - sys[sys.possible_hoppings(*hopping)] = - t - -sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) -lead0 = kwant.Builder(sym_lead0) -lead0.default_site_group = lat - -def lead_shape(pos): - (x, y) = pos - return (-1 < x < 1) and ( -W/2 < y < W/2 ) - -lead0[lat.shape(lead_shape, (0, 0))] = 4 * t -for hopping in lat.nearest: - lead0[lead0.possible_hoppings(*hopping)] = - t - -# Then the lead to the right -# there we can use a special function that simply reverses the direction - -lead1 = lead0.reversed() - -# Then attach the leads to the system - -sys.attach_lead(lead0) -sys.attach_lead(lead1, lat(0, 0)) - -# finalize the system - -fsys = sys.finalized() - -# and plot it, to make sure it's proper +import latex, html -kwant.plot(fsys, "tutorial2c_note2.pdf", width=latex.figwidth_small_pt) -kwant.plot(fsys, "tutorial2c_note2.png", width=html.figwidth_small_px) +def make_system(a=1, t=1.0, W=10, r1=10, r2=20): + # Start with an empty tight-binding system and a single square lattice. + # `a` is the lattice constant (by default set to 1 for simplicity). + + lat = kwant.lattice.Square(a) + + sys = kwant.Builder() + + #### Define the scattering region. #### + # Now, we aim for a more complex shape, namely a ring (or annulus) + def ring(pos): + (x, y) = pos + rsq = x**2 + y**2 + return ( r1**2 < rsq < r2**2) + + # and add the corresponding lattice points using the `shape`-function + sys[lat.shape(ring, (0, 11))] = 4 * t + for hopping in lat.nearest: + sys[sys.possible_hoppings(*hopping)] = - t + + # In order to introduce a flux through the ring, we introduce a phase + # on the hoppings on the line cut through one of the arms + + # since we want to change the flux without modifying Builder repeatedly, + # we define the modified hoppings as a function that takes the flux + # through the global variable phi. + def fluxphase(site1, site2): + return exp(1j * phi) + + def crosses_branchcut(hop): + ix0, iy0 = hop[0].tag + + # possible_hoppings with the argument (1, 0) below + # returns hoppings ordered as ((i+1, j), (i, j)) + return iy0 < 0 and ix0 == 1 # ix1 == 0 then implied + + # Modify only those hopings in x-direction that cross the branch cut + sys[(hop for hop in sys.possible_hoppings((1,0), lat, lat) + if crosses_branchcut(hop))] = fluxphase + + #### Define the leads. #### + # left lead + sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) + lead0 = kwant.Builder(sym_lead0) + + def lead_shape(pos): + (x, y) = pos + return (-1 < x < 1) and ( -W/2 < y < W/2 ) + + lead0[lat.shape(lead_shape, (0, 0))] = 4 * t + for hopping in lat.nearest: + lead0[lead0.possible_hoppings(*hopping)] = - t + + # Then the lead to the right + # (again, obtained using reverse() + lead1 = lead0.reversed() + + #### Attach the leads and return the finalized system. #### + sys.attach_lead(lead0) + sys.attach_lead(lead1) + + return sys.finalized() + + +def make_system_note1(a=1, t=1.0, W=10, r1=10, r2=20): + lat = kwant.lattice.Square(a) + sys = kwant.Builder() + def ring(pos): + (x, y) = pos + rsq = x**2 + y**2 + return ( r1**2 < rsq < r2**2) + sys[lat.shape(ring, (0, 11))] = 4 * t + for hopping in lat.nearest: + sys[sys.possible_hoppings(*hopping)] = - t + sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) + lead0 = kwant.Builder(sym_lead0) + def lead_shape(pos): + (x, y) = pos + return (-1 < x < 1) and ( 0.5 * W < y < 1.5 * W ) + lead0[lat.shape(lead_shape, (0, W))] = 4 * t + for hopping in lat.nearest: + lead0[lead0.possible_hoppings(*hopping)] = - t + lead1 = lead0.reversed() + sys.attach_lead(lead0) + sys.attach_lead(lead1) + return sys.finalized() + + +def make_system_note2(a=1, t=1.0, W=10, r1=10, r2=20): + lat = kwant.lattice.Square(a) + sys = kwant.Builder() + def ring(pos): + (x, y) = pos + rsq = x**2 + y**2 + return ( r1**2 < rsq < r2**2) + sys[lat.shape(ring, (0, 11))] = 4 * t + for hopping in lat.nearest: + sys[sys.possible_hoppings(*hopping)] = - t + sym_lead0 = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) + lead0 = kwant.Builder(sym_lead0) + def lead_shape(pos): + (x, y) = pos + return (-1 < x < 1) and ( -W/2 < y < W/2 ) + lead0[lat.shape(lead_shape, (0, 0))] = 4 * t + for hopping in lat.nearest: + lead0[lead0.possible_hoppings(*hopping)] = - t + lead1 = lead0.reversed() + sys.attach_lead(lead0) + sys.attach_lead(lead1, lat(0, 0)) + return sys.finalized() + + +def plot_conductance(fsys, energy, fluxes): + # compute conductance + # global variable phi controls the flux + global phi + + normalized_fluxes = [flux/(2 * pi) for flux in fluxes] + data = [] + for flux in fluxes: + phi = flux + + smatrix = kwant.solve(fsys, energy) + data.append(smatrix.transmission(1, 0)) + + pylab.plot(normalized_fluxes, data) + pylab.xlabel("flux [in units of the flux quantum]", + fontsize=latex.mpl_label_size) + pylab.ylabel("conductance [in units of e^2/h]", + fontsize=latex.mpl_label_size) + fig = pylab.gcf() + pylab.setp(fig.get_axes()[0].get_xticklabels(), + fontsize=latex.mpl_tick_size) + pylab.setp(fig.get_axes()[0].get_yticklabels(), + fontsize=latex.mpl_tick_size) + fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) + fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) + fig.savefig("tutorial2c_result.pdf") + fig.savefig("tutorial2c_result.png", + dpi=(html.figwidth_px/latex.mpl_width_in)) + + +def main(): + fsys = make_system() + + # Check that the system looks as intended. + kwant.plot(fsys, "tutorial2c_sys.pdf", width=latex.figwidth_pt) + kwant.plot(fsys, "tutorial2c_sys.png", width=html.figwidth_px) + + # We should see a conductance that is periodic with the flux quantum + plot_conductance(fsys, energy=0.15, fluxes=[0.01 * i * 3 * 2 * pi + for i in xrange(100)]) + + # Finally, some plots needed for the notes + fsys = make_system_note1() + kwant.plot(fsys, "tutorial2c_note1.pdf", width=latex.figwidth_small_pt) + kwant.plot(fsys, "tutorial2c_note1.png", width=html.figwidth_small_px) + fsys = make_system_note2() + kwant.plot(fsys, "tutorial2c_note2.pdf", width=latex.figwidth_small_pt) + kwant.plot(fsys, "tutorial2c_note2.png", width=html.figwidth_small_px) + + +# 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() diff --git a/doc/source/images/tutorial3a.py b/doc/source/images/tutorial3a.py index 440a000a..2fa54f7a 100644 --- a/doc/source/images/tutorial3a.py +++ b/doc/source/images/tutorial3a.py @@ -7,52 +7,70 @@ # - Computing the band structure of a finalized lead. import kwant + import numpy as np from math import pi +# For plotting +import pylab + import latex, html -a = 1 -lat = kwant.lattice.Square(a) -t = 1.0 -W = 10 +def make_lead(a=1, t=1.0, W=10): + # Start with an empty lead with a single square lattice + lat = kwant.lattice.Square(a) -# Define a lead: + sym_lead = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) + lead = kwant.Builder(sym_lead) + lead.default_site_group = lat -sym_lead = kwant.TranslationalSymmetry([lat.vec((-1, 0))]) -lead = kwant.Builder(sym_lead) -lead.default_site_group = lat + # build up one unit cell of the lead, and add the hoppings + # to the next unit cell + for j in xrange(W): + lead[(0, j)] = 4 * t -for j in xrange(W): - lead[(0, j)] = 4 * t + if j > 0: + lead[(0, j), (0, j-1)] = - t - if j > 0: - lead[(0, j), (0, j-1)] = - t + lead[(1, j), (0, j)] = - t - lead[(1, j), (0, j)] = - t + # return a finalized lead + return lead.finalized() -# Now compute the band structure -# Only a finalized lead has the information about bandstructure -flead = lead.finalized() +def plot_bandstructure(flead, momenta): + # Use the method ``energies`` of the finalized lead to compute + # the bandstructure + energy_list = [flead.energies(k) for k in momenta] -momenta = np.arange(-pi, pi + .01, 0.02 * pi) -energy_list = [flead.energies(k) for k in momenta] + pylab.plot(momenta, energy_list) + pylab.xlabel("momentum [in untis of (lattice constant)^-1]", + fontsize=latex.mpl_label_size) + pylab.ylabel("energy [in units of t]", + fontsize=latex.mpl_label_size) + fig = pylab.gcf() + pylab.setp(fig.get_axes()[0].get_xticklabels(), + fontsize=latex.mpl_tick_size) + pylab.setp(fig.get_axes()[0].get_yticklabels(), + fontsize=latex.mpl_tick_size) + fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) + fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) + fig.savefig("tutorial3a_result.pdf") + fig.savefig("tutorial3a_result.png", + dpi=(html.figwidth_px/latex.mpl_width_in)) -import pylab -pylab.plot(momenta, energy_list) -pylab.xlabel("momentum [in untis of (lattice constant)^-1]", - fontsize=latex.mpl_label_size) -pylab.ylabel("energy [in units of t]", - fontsize=latex.mpl_label_size) -fig = pylab.gcf() -pylab.setp(fig.get_axes()[0].get_xticklabels(), - fontsize=latex.mpl_tick_size) -pylab.setp(fig.get_axes()[0].get_yticklabels(), - fontsize=latex.mpl_tick_size) -fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) -fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) -fig.savefig("tutorial3a_result.pdf") -fig.savefig("tutorial3a_result.png", - dpi=(html.figwidth_px/latex.mpl_width_in)) + +def main(): + flead = make_lead() + + # list of momenta at which the bands should be computed + momenta = np.arange(-pi, pi + .01, 0.02 * pi) + + plot_bandstructure(flead, momenta) + + +# 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() diff --git a/doc/source/images/tutorial3b.py b/doc/source/images/tutorial3b.py index c219a0f3..c0ac9db3 100644 --- a/doc/source/images/tutorial3b.py +++ b/doc/source/images/tutorial3b.py @@ -12,83 +12,91 @@ from cmath import exp import kwant -import latex, html +# For eigenvalue computation +import scipy.linalg as la -# First, define the tight-binding system +# For plotting +import pylab -sys = kwant.Builder() +import latex, html -# Here, we are only working with square lattices -# for simplicity, take lattice constant = 1 -a = 1 -lat = kwant.lattice.Square(a) +def make_system(a=1, t=1.0, r=10): + # Start with an empty tight-binding system and a single square lattice. + # `a` is the lattice constant (by default set to 1 for simplicity). -t = 1.0 -r = 10 + lat = kwant.lattice.Square(a) -# Define the quantum dot + sys = kwant.Builder() -def circle(pos): - (x, y) = pos - rsq = x**2 + y**2 - return rsq < r**2 + # Define the quantum dot + def circle(pos): + (x, y) = pos + rsq = x**2 + y**2 + return rsq < r**2 -def hopx(site1, site2): - y = site1.pos[1] - return - t * exp(-1j * B * y) + def hopx(site1, site2): + # The magnetic field is controlled by the global variable B + y = site1.pos[1] + return - t * exp(-1j * B * y) -sys[lat.shape(circle, (0, 0))] = 4 * t -# hoppings in x-direction -sys[sys.possible_hoppings((1, 0), lat, lat)] = hopx -# hoppings in y-directions -sys[sys.possible_hoppings((0, 1), lat, lat)] = - t + sys[lat.shape(circle, (0, 0))] = 4 * t + # hoppings in x-direction + sys[sys.possible_hoppings((1, 0), lat, lat)] = hopx + # hoppings in y-directions + sys[sys.possible_hoppings((0, 1), lat, lat)] = - t -# It's a closed system for a change, so no leads + # It's a closed system for a change, so no leads + return sys.finalized() -# finalize the system -fsys = sys.finalized() +def plot_spectrum(fsys, Bfields): + # global variable B controls the magnetic field + global B -# and plot it, to make sure it's proper + # In the following, we compute the spectrum of the quantum dot + # using dense matrix methods. This works in this toy example, as + # the system is tiny. In a real example, one would want to use + # sparse matrix methods -kwant.plot(fsys, "tutorial3b_sys.pdf") -kwant.plot(fsys, "tutorial3b_sys.png") + energies = [] + for Bfield in Bfields: + B = Bfield -# In the following, we compute the spectrum of the quantum dot -# using dense matrix methods. This works in this toy example, -# as the system is tiny. In a real example, one would want to use -# sparse matrix methods + # Obtain the Hamiltonian as a dense matrix + ham_mat = fsys.hamiltonian_submatrix()[0] -import scipy.linalg as la + ev = la.eigh(ham_mat, eigvals_only=True) -Bs = [] -energies = [] -for iB in xrange(100): - B = iB * 0.002 + # we only plot the 15 lowest eigenvalues + energies.append(ev[:15]) -# Obtain the Hamiltonian as a dense matrix - ham_mat = fsys.hamiltonian_submatrix()[0] + pylab.plot(Bfields, energies) + pylab.xlabel("magnetic field [some arbitrary units]", + fontsize=latex.mpl_label_size) + pylab.ylabel("energy [in units of t]", + fontsize=latex.mpl_label_size) + fig = pylab.gcf() + pylab.setp(fig.get_axes()[0].get_xticklabels(), + fontsize=latex.mpl_tick_size) + pylab.setp(fig.get_axes()[0].get_yticklabels(), + fontsize=latex.mpl_tick_size) + fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) + fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) + fig.savefig("tutorial3b_result.pdf") + fig.savefig("tutorial3b_result.png", + dpi=(html.figwidth_px/latex.mpl_width_in)) - ev = la.eigh(ham_mat, eigvals_only=True) - Bs.append(B) - energies.append(ev[:15]) +def main(): + fsys = make_system() + + # We should observe energy levels that flow towards Landau + # level energies with increasing magnetic field + plot_spectrum(fsys, [iB * 0.002 for iB in xrange(100)]) -import pylab -pylab.plot(Bs, energies) -pylab.xlabel("magnetic field [some arbitrary units]", - fontsize=latex.mpl_label_size) -pylab.ylabel("energy [in units of t]", - fontsize=latex.mpl_label_size) -fig = pylab.gcf() -pylab.setp(fig.get_axes()[0].get_xticklabels(), - fontsize=latex.mpl_tick_size) -pylab.setp(fig.get_axes()[0].get_yticklabels(), - fontsize=latex.mpl_tick_size) -fig.set_size_inches(latex.mpl_width_in, latex.mpl_width_in*3./4.) -fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15) -fig.savefig("tutorial3b_result.pdf") -fig.savefig("tutorial3b_result.png", - dpi=(html.figwidth_px/latex.mpl_width_in)) +# 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() -- GitLab