Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
  • kwant/kwant
  • jbweston/kwant
  • anton-akhmerov/kwant
  • cwg/kwant
  • Mathieu/kwant
  • slavoutich/kwant
  • pacome/kwant
  • behrmann/kwant
  • michaelwimmer/kwant
  • albeercik/kwant
  • eunjongkim/kwant
  • basnijholt/kwant
  • r-j-skolasinski/kwant
  • sahmed95/kwant
  • pablopiskunow/kwant
  • mare/kwant
  • dvarjas/kwant
  • Paul/kwant
  • bbuijtendorp/kwant
  • tkloss/kwant
  • torosdahl/kwant
  • kel85uk/kwant
  • kpoyhonen/kwant
  • Fromeworld/kwant
  • quaeritis/kwant
  • marwahaha/kwant
  • fernandodfufrpe/kwant
  • oly/kwant
  • jiamingh/kwant
  • mehdi2369/kwant
  • ValFadeev/kwant
  • Kostas/kwant
  • chelseabaptiste03/kwant
33 results
Show changes
Hide whitespace changes
Inline Side-by-side
Showing
with 1318 additions and 899 deletions
doc/source/code/figure/quantum_well.py.diff deleted 100644 → 0
View file @ 6be96689
@@ -1,88 +1,95 @@
# Tutorial 2.3.2. Spatially dependent values through functions
# ============================================================
#
# Physics background
# ------------------
# transmission through a quantum well
#
# Kwant features highlighted
# --------------------------
# - Functions as values in Builder
+import _defs
import kwant
# For plotting
from matplotlib import pyplot
#HIDDEN_BEGIN_ehso
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)
syst = kwant.Builder()
#### Define the scattering region. ####
# Potential profile
def potential(site, pot):
(x, y) = site.pos
if (L - L_well) / 2 < x < (L + L_well) / 2:
return pot
else:
return 0
#HIDDEN_END_ehso
#HIDDEN_BEGIN_coid
def onsite(site, pot):
return 4 * t + potential(site, pot)
syst[(lat(x, y) for x in range(L) for y in range(W))] = onsite
syst[lat.neighbors()] = -t
#HIDDEN_END_coid
#### Define and attach the leads. ####
lead = kwant.Builder(kwant.TranslationalSymmetry((-a, 0)))
lead[(lat(0, j) for j in range(W))] = 4 * t
lead[lat.neighbors()] = -t
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
return syst
def plot_conductance(syst, energy, welldepths):
#HIDDEN_BEGIN_sqvr
# Compute conductance
data = []
for welldepth in welldepths:
smatrix = kwant.smatrix(syst, energy, params=dict(pot=-welldepth))
data.append(smatrix.transmission(1, 0))
- pyplot.figure()
+ fig = pyplot.figure()
pyplot.plot(welldepths, data)
- pyplot.xlabel("well depth [t]")
- pyplot.ylabel("conductance [e^2/h]")
- pyplot.show()
+ pyplot.xlabel("well depth [t]",
+ fontsize=_defs.mpl_label_size)
+ pyplot.ylabel("conductance [e^2/h]",
+ fontsize=_defs.mpl_label_size)
+ pyplot.setp(fig.get_axes()[0].get_xticklabels(),
+ fontsize=_defs.mpl_tick_size)
+ pyplot.setp(fig.get_axes()[0].get_yticklabels(),
+ fontsize=_defs.mpl_tick_size)
+ fig.set_size_inches(_defs.mpl_width_in, _defs.mpl_width_in * 3. / 4.)
+ fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15)
+ for extension in ('pdf', 'png'):
+ fig.savefig("quantum_well_result." + extension, dpi=_defs.dpi)
#HIDDEN_END_sqvr
def main():
syst = make_system()
- # Check that the system looks as intended.
- kwant.plot(syst)
-
# Finalize the system.
syst = syst.finalized()
# We should see conductance steps.
plot_conductance(syst, energy=0.2,
welldepths=[0.01 * i for i in range(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()
doc/source/code/figure/quantum_wire.py.diff deleted 100644 → 0
View file @ 6be96689
@@ -1,121 +1,130 @@
# Tutorial 2.2.2. Transport through a quantum wire
# ================================================
#
# Physics background
# ------------------
# Conductance of a quantum wire; subbands
#
# Kwant features highlighted
# --------------------------
# - Builder for setting up transport systems easily
# - Making scattering region and leads
# - Using the simple sparse solver for computing Landauer conductance
+import _defs
from matplotlib import pyplot
#HIDDEN_BEGIN_dwhx
import kwant
#HIDDEN_END_dwhx
# First, define the tight-binding system
#HIDDEN_BEGIN_goiq
syst = kwant.Builder()
#HIDDEN_END_goiq
# Here, we are only working with square lattices
#HIDDEN_BEGIN_suwo
a = 1
lat = kwant.lattice.square(a)
#HIDDEN_END_suwo
#HIDDEN_BEGIN_zfvr
t = 1.0
W = 10
L = 30
# Define the scattering region
for i in range(L):
for j in range(W):
# On-site Hamiltonian
syst[lat(i, j)] = 4 * t
# Hopping in y-direction
if j > 0:
syst[lat(i, j), lat(i, j - 1)] = -t
# Hopping in x-direction
if i > 0:
syst[lat(i, j), lat(i - 1, j)] = -t
#HIDDEN_END_zfvr
# Then, define and attach the leads:
# First the lead to the left
# (Note: TranslationalSymmetry takes a real-space vector)
#HIDDEN_BEGIN_xcmc
sym_left_lead = kwant.TranslationalSymmetry((-a, 0))
left_lead = kwant.Builder(sym_left_lead)
#HIDDEN_END_xcmc
#HIDDEN_BEGIN_ndez
for j in range(W):
left_lead[lat(0, j)] = 4 * t
if j > 0:
left_lead[lat(0, j), lat(0, j - 1)] = -t
left_lead[lat(1, j), lat(0, j)] = -t
#HIDDEN_END_ndez
#HIDDEN_BEGIN_fskr
syst.attach_lead(left_lead)
#HIDDEN_END_fskr
# Then the lead to the right
#HIDDEN_BEGIN_xhqc
sym_right_lead = kwant.TranslationalSymmetry((a, 0))
right_lead = kwant.Builder(sym_right_lead)
for j in range(W):
right_lead[lat(0, j)] = 4 * t
if j > 0:
right_lead[lat(0, j), lat(0, j - 1)] = -t
right_lead[lat(1, j), lat(0, j)] = -t
syst.attach_lead(right_lead)
#HIDDEN_END_xhqc
# Plot it, to make sure it's OK
#HIDDEN_BEGIN_wsgh
-kwant.plot(syst)
+size = (_defs.figwidth_in, 0.3 * _defs.figwidth_in)
+for extension in ('pdf', 'png'):
+ kwant.plot(syst, file="quantum_wire_syst." + extension,
+ fig_size=size, dpi=_defs.dpi)
#HIDDEN_END_wsgh
# Finalize the system
#HIDDEN_BEGIN_dngj
syst = syst.finalized()
#HIDDEN_END_dngj
# Now that we have the system, we can compute conductance
#HIDDEN_BEGIN_buzn
energies = []
data = []
for ie in range(100):
energy = ie * 0.01
# compute the scattering matrix at a given energy
smatrix = kwant.smatrix(syst, energy)
# compute the transmission probability from lead 0 to
# lead 1
energies.append(energy)
data.append(smatrix.transmission(1, 0))
#HIDDEN_END_buzn
# Use matplotlib to write output
# We should see conductance steps
#HIDDEN_BEGIN_lliv
-pyplot.figure()
+fig = pyplot.figure()
pyplot.plot(energies, data)
-pyplot.xlabel("energy [t]")
-pyplot.ylabel("conductance [e^2/h]")
-pyplot.show()
+pyplot.xlabel("energy [t]", fontsize=_defs.mpl_label_size)
+pyplot.ylabel("conductance [e^2/h]", fontsize=_defs.mpl_label_size)
+pyplot.setp(fig.get_axes()[0].get_xticklabels(), fontsize=_defs.mpl_tick_size)
+pyplot.setp(fig.get_axes()[0].get_yticklabels(), fontsize=_defs.mpl_tick_size)
+fig.set_size_inches(_defs.mpl_width_in, _defs.mpl_width_in * 3. / 4.)
+fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15)
+for extension in ('pdf', 'png'):
+ fig.savefig("quantum_wire_result." + extension, dpi=_defs.dpi)
#HIDDEN_END_lliv
doc/source/code/figure/spin_orbit.py.diff deleted 100644 → 0
View file @ 6be96689
@@ -1,104 +1,110 @@
# Tutorial 2.3.1. Matrix structure of on-site and hopping elements
# ================================================================
#
# Physics background
# ------------------
# Gaps in quantum wires with spin-orbit coupling and Zeeman splititng,
# as theoretically predicted in
# http://prl.aps.org/abstract/PRL/v90/i25/e256601
# and (supposedly) experimentally oberved in
# http://www.nature.com/nphys/journal/v6/n5/abs/nphys1626.html
#
# Kwant features highlighted
# --------------------------
# - Numpy matrices as values in Builder
+import _defs
import kwant
# For plotting
from matplotlib import pyplot
# For matrix support
#HIDDEN_BEGIN_xumz
import tinyarray
#HIDDEN_END_xumz
# define Pauli-matrices for convenience
#HIDDEN_BEGIN_hwbt
sigma_0 = tinyarray.array([[1, 0], [0, 1]])
sigma_x = tinyarray.array([[0, 1], [1, 0]])
sigma_y = tinyarray.array([[0, -1j], [1j, 0]])
sigma_z = tinyarray.array([[1, 0], [0, -1]])
#HIDDEN_END_hwbt
def make_system(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()
syst = kwant.Builder()
#### Define the scattering region. ####
#HIDDEN_BEGIN_uxrm
syst[(lat(x, y) for x in range(L) for y in range(W))] = \
4 * t * sigma_0 + e_z * sigma_z
# hoppings in x-direction
syst[kwant.builder.HoppingKind((1, 0), lat, lat)] = \
-t * sigma_0 + 1j * alpha * sigma_y / 2
# hoppings in y-directions
syst[kwant.builder.HoppingKind((0, 1), lat, lat)] = \
-t * sigma_0 - 1j * alpha * sigma_x / 2
#HIDDEN_END_uxrm
#### Define the left lead. ####
lead = kwant.Builder(kwant.TranslationalSymmetry((-1, 0)))
#HIDDEN_BEGIN_yliu
lead[(lat(0, j) for j in range(W))] = 4 * t * sigma_0 + e_z * sigma_z
# hoppings in x-direction
lead[kwant.builder.HoppingKind((1, 0), lat, lat)] = \
-t * sigma_0 + 1j * alpha * sigma_y / 2
# hoppings in y-directions
lead[kwant.builder.HoppingKind((0, 1), lat, lat)] = \
-t * sigma_0 - 1j * alpha * sigma_x / 2
#HIDDEN_END_yliu
#### Attach the leads and return the finalized system. ####
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
return syst
def plot_conductance(syst, energies):
# Compute conductance
data = []
for energy in energies:
smatrix = kwant.smatrix(syst, energy)
data.append(smatrix.transmission(1, 0))
- pyplot.figure()
+ fig = pyplot.figure()
pyplot.plot(energies, data)
- pyplot.xlabel("energy [t]")
- pyplot.ylabel("conductance [e^2/h]")
- pyplot.show()
+ pyplot.xlabel("energy [t]", fontsize=_defs.mpl_label_size)
+ pyplot.ylabel("conductance [e^2/h]",
+ fontsize=_defs.mpl_label_size)
+ pyplot.setp(fig.get_axes()[0].get_xticklabels(),
+ fontsize=_defs.mpl_tick_size)
+ pyplot.setp(fig.get_axes()[0].get_yticklabels(),
+ fontsize=_defs.mpl_tick_size)
+ fig.set_size_inches(_defs.mpl_width_in, _defs.mpl_width_in * 3. / 4.)
+ fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15)
+ for extension in ('pdf', 'png'):
+ fig.savefig("spin_orbit_result." + extension, dpi=_defs.dpi)
def main():
syst = make_system()
- # Check that the system looks as intended.
- kwant.plot(syst)
-
# Finalize the system.
syst = syst.finalized()
# We should see non-monotonic conductance steps.
plot_conductance(syst, energies=[0.01 * i - 0.3 for i in range(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()
doc/source/code/figure/superconductor.py.diff deleted 100644 → 0
View file @ 6be96689
@@ -1,132 +1,141 @@
# Tutorial 2.6. "Superconductors": orbitals, conservation laws and symmetries
# ===========================================================================
#
# Physics background
# ------------------
# - conductance of a NS-junction (Andreev reflection, superconducting gap)
#
# Kwant features highlighted
# --------------------------
# - Implementing electron and hole ("orbital") degrees of freedom
# using conservation laws.
# - Use of discrete symmetries to relate scattering states.
+import _defs
import kwant
import tinyarray
import numpy as np
# For plotting
from matplotlib import pyplot
+from contextlib import redirect_stdout
tau_x = tinyarray.array([[0, 1], [1, 0]])
tau_y = tinyarray.array([[0, -1j], [1j, 0]])
tau_z = tinyarray.array([[1, 0], [0, -1]])
#HIDDEN_BEGIN_nbvn
def make_system(a=1, W=10, L=10, barrier=1.5, barrierpos=(3, 4),
mu=0.4, Delta=0.1, Deltapos=4, t=1.0):
# Start with an empty tight-binding system. On each site, there
# are now electron and hole orbitals, so we must specify the
# number of orbitals per site. The orbital structure is the same
# as in the Hamiltonian.
lat = kwant.lattice.square(norbs=2)
syst = kwant.Builder()
#### Define the scattering region. ####
# The superconducting order parameter couples electron and hole orbitals
# on each site, and hence enters as an onsite potential.
# The pairing is only included beyond the point 'Deltapos' in the scattering region.
syst[(lat(x, y) for x in range(Deltapos) for y in range(W))] = (4 * t - mu) * tau_z
syst[(lat(x, y) for x in range(Deltapos, L) for y in range(W))] = (4 * t - mu) * tau_z + Delta * tau_x
# The tunnel barrier
syst[(lat(x, y) for x in range(barrierpos[0], barrierpos[1])
for y in range(W))] = (4 * t + barrier - mu) * tau_z
# Hoppings
syst[lat.neighbors()] = -t * tau_z
#HIDDEN_END_nbvn
#HIDDEN_BEGIN_ttth
#### Define the leads. ####
# Left lead - normal, so the order parameter is zero.
sym_left = kwant.TranslationalSymmetry((-a, 0))
# Specify the conservation law used to treat electrons and holes separately.
# We only do this in the left lead, where the pairing is zero.
lead0 = kwant.Builder(sym_left, conservation_law=-tau_z, particle_hole=tau_y)
lead0[(lat(0, j) for j in range(W))] = (4 * t - mu) * tau_z
lead0[lat.neighbors()] = -t * tau_z
#HIDDEN_END_ttth
#HIDDEN_BEGIN_zuuw
# Right lead - superconducting, so the order parameter is included.
sym_right = kwant.TranslationalSymmetry((a, 0))
lead1 = kwant.Builder(sym_right)
lead1[(lat(0, j) for j in range(W))] = (4 * t - mu) * tau_z + Delta * tau_x
lead1[lat.neighbors()] = -t * tau_z
#### Attach the leads and return the system. ####
syst.attach_lead(lead0)
syst.attach_lead(lead1)
return syst
#HIDDEN_END_zuuw
#HIDDEN_BEGIN_jbjt
def plot_conductance(syst, energies):
# Compute conductance
data = []
for energy in energies:
smatrix = kwant.smatrix(syst, energy)
# Conductance is N - R_ee + R_he
data.append(smatrix.submatrix((0, 0), (0, 0)).shape[0] -
smatrix.transmission((0, 0), (0, 0)) +
smatrix.transmission((0, 1), (0, 0)))
#HIDDEN_END_jbjt
- pyplot.figure()
+ fig = pyplot.figure()
pyplot.plot(energies, data)
pyplot.xlabel("energy [t]")
pyplot.ylabel("conductance [e^2/h]")
- pyplot.show()
+ pyplot.setp(fig.get_axes()[0].get_xticklabels(),
+ fontsize=_defs.mpl_tick_size)
+ pyplot.setp(fig.get_axes()[0].get_yticklabels(),
+ fontsize=_defs.mpl_tick_size)
+ fig.set_size_inches(_defs.mpl_width_in, _defs.mpl_width_in * 3. / 4.)
+ fig.subplots_adjust(left=0.15, right=0.95, top=0.95, bottom=0.15)
+ for extension in ('pdf', 'png'):
+ fig.savefig("superconductor_transport_result." + extension,
+ dpi=_defs.dpi)
#HIDDEN_BEGIN_pqmp
def check_PHS(syst):
# Scattering matrix
s = kwant.smatrix(syst, energy=0)
# Electron to electron block
s_ee = s.submatrix((0,0), (0,0))
# Hole to hole block
s_hh = s.submatrix((0,1), (0,1))
print('s_ee: \n', np.round(s_ee, 3))
print('s_hh: \n', np.round(s_hh[::-1, ::-1], 3))
print('s_ee - s_hh^*: \n',
np.round(s_ee - s_hh[::-1, ::-1].conj(), 3), '\n')
# Electron to hole block
s_he = s.submatrix((0,1), (0,0))
# Hole to electron block
s_eh = s.submatrix((0,0), (0,1))
print('s_he: \n', np.round(s_he, 3))
print('s_eh: \n', np.round(s_eh[::-1, ::-1], 3))
print('s_he + s_eh^*: \n',
np.round(s_he + s_eh[::-1, ::-1].conj(), 3))
#HIDDEN_END_pqmp
def main():
syst = make_system(W=10)
- # Check that the system looks as intended.
- kwant.plot(syst)
-
# Finalize the system.
syst = syst.finalized()
# Check particle-hole symmetry of the scattering matrix
- check_PHS(syst)
+ with open('check_PHS_out.txt', 'w') as f:
+ with redirect_stdout(f):
+ check_PHS(syst)
# Compute and plot the conductance
plot_conductance(syst, energies=[0.002 * i for i in range(-10, 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()
doc/source/code/include/quantum_wire_revisited.py deleted 100644 → 0
View file @ 6be96689
# Tutorial 2.2.3. Building the same system with less code
# =======================================================
#
# Physics background
# ------------------
# Conductance of a quantum wire; subbands
#
# Kwant features highlighted
# --------------------------
# - Using iterables and builder.HoppingKind for making systems
# - introducing `reversed()` for the leads
#
# Note: Does the same as tutorial1a.py, but using other features of Kwant.
#HIDDEN_BEGIN_xkzy
import kwant
# For plotting
from matplotlib import pyplot
def make_system(a=1, t=1.0, 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)
syst = kwant.Builder()
#HIDDEN_END_xkzy
#### Define the scattering region. ####
#HIDDEN_BEGIN_vvjt
syst[(lat(x, y) for x in range(L) for y in range(W))] = 4 * t
#HIDDEN_END_vvjt
#HIDDEN_BEGIN_nooi
syst[lat.neighbors()] = -t
#HIDDEN_END_nooi
#### Define and attach the leads. ####
# Construct the left lead.
#HIDDEN_BEGIN_iepx
lead = kwant.Builder(kwant.TranslationalSymmetry((-a, 0)))
lead[(lat(0, j) for j in range(W))] = 4 * t
lead[lat.neighbors()] = -t
#HIDDEN_END_iepx
# Attach the left lead and its reversed copy.
#HIDDEN_BEGIN_yxot
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
return syst
#HIDDEN_END_yxot
#HIDDEN_BEGIN_ayuk
def plot_conductance(syst, energies):
# Compute conductance
data = []
for energy in energies:
smatrix = kwant.smatrix(syst, energy)
data.append(smatrix.transmission(1, 0))
pyplot.figure()
pyplot.plot(energies, data)
pyplot.xlabel("energy [t]")
pyplot.ylabel("conductance [e^2/h]")
pyplot.show()
#HIDDEN_END_ayuk
#HIDDEN_BEGIN_cjel
def main():
syst = make_system()
# Check that the system looks as intended.
kwant.plot(syst)
# Finalize the system.
syst = syst.finalized()
# We should see conductance steps.
plot_conductance(syst, energies=[0.01 * i for i in range(100)])
#HIDDEN_END_cjel
# Call the main function if the script gets executed (as opposed to imported).
# See <http://docs.python.org/library/__main__.html>.
#HIDDEN_BEGIN_ypbj
if __name__ == '__main__':
main()
#HIDDEN_END_ypbj
doc/source/conf.py
View file @ 39b65acd
......@@ -15,8 +15,15 @@ import sys
import os
import string
from distutils.util import get_platform
sys.path.insert(0, "../../build/lib.{0}-{1}.{2}".format(
get_platform(), *sys.version_info[:2]))
package_path = os.path.abspath(
"../../build/lib.{0}-{1}.{2}"
.format(get_platform(), *sys.version_info[:2]))
# Insert into sys.path so that we can import kwant here
sys.path.insert(0, package_path)
# Insert into PYTHONPATH so that jupyter-sphinx will pick it up
os.environ['PYTHONPATH'] = ':'.join((package_path, os.environ.get('PYTHONPATH','')))
import kwant
import kwant.qsymm
......@@ -31,7 +38,8 @@ sys.path.insert(0, os.path.abspath('../sphinxext'))
extensions = ['sphinx.ext.autodoc', 'sphinx.ext.autosummary',
'sphinx.ext.todo', 'sphinx.ext.mathjax', 'numpydoc',
'kwantdoc', 'sphinx.ext.linkcode']
'kwantdoc', 'sphinx.ext.linkcode', 'jupyter_sphinx.execute',
'sphinxcontrib.rsvgconverter']
# Add any paths that contain templates here, relative to this directory.
templates_path = ['../templates']
......@@ -99,6 +107,15 @@ pygments_style = 'sphinx'
# Do not show all class members automatically in the class documentation
numpydoc_show_class_members = False
# Jupyter Sphinx config
jupyter_sphinx_thebelab_config = {
"binderOptions": {
"repo": "kwant-project/binder",
"ref": "master",
}
}
# -- Options for HTML output ---------------------------------------------------
# The theme to use for HTML and HTML Help pages. Major themes that come with
......@@ -218,7 +235,7 @@ latex_documents = [
latex_engine = 'xelatex'
latex_use_xindy = False
latex_use_xindy = False # Xindy not installable in CI environment
# The name of an image file (relative to this directory) to place at the top of
# the title page.
......@@ -239,8 +256,9 @@ latex_domain_indices = False
autosummary_generate = True
autoclass_content = "both"
autodoc_default_flags = ['show-inheritance']
autodoc_default_options = {
'show-inheritance': True,
}
# -- Teach Sphinx to document bound methods like functions ---------------------
import types
......
doc/source/pre/whatsnew/1.3.rst
View file @ 39b65acd
......@@ -28,11 +28,43 @@ tight-binding band structures or construct systems with different/lower
symmetry. For example in just a few lines we can construct a two-band model that exhibits
a quantum anomalous spin Hall phase:
.. literalinclude:: ../../code/include/plot_qahe.py
:start-after: HIDDEN_BEGIN_model
:end-before: HIDDEN_END_model
.. jupyter-kernel::
:id: plot_qahe
.. jupyter-execute::
:hide-code:
# Comprehensive example: quantum anomalous Hall effect
# ====================================================
#
# Physics background
# ------------------
# + Quantum anomalous Hall effect
#
# Features highlighted
# --------------------
# + Use of `kwant.continuum` to discretize a continuum Hamiltonian
# + Use of `kwant.operator` to compute local current
# + Use of `kwant.plotter.current` to plot local current
import math
import matplotlib.pyplot
import kwant
import kwant.continuum
.. jupyter-execute:: ../../tutorial/boilerplate.py
:hide-code:
.. jupyter-execute::
From: :download:`QAHE example script <../../code/download/plot_qahe.py>`
def make_model(a):
ham = ("alpha * (k_x * sigma_x - k_y * sigma_y)"
"+ (m + beta * kk) * sigma_z"
"+ (gamma * kk + U) * sigma_0")
subs = {"kk": "k_x**2 + k_y**2"}
return kwant.continuum.discretize(ham, locals=subs, grid=a)
From: :jupyter-download:script:`plot_qahe`
See the tutorial: :doc:`../../tutorial/discretize`
......@@ -71,13 +103,47 @@ The example below shows edge states of a quantum anomalous Hall phase
of the two-band model shown in the `above section
<#tools-for-continuum-hamiltonians>`_:
.. literalinclude:: ../../code/include/plot_qahe.py
:start-after: HIDDEN_BEGIN_current
:end-before: HIDDEN_END_current
.. jupyter-execute::
:hide-code:
def make_system(model, L):
def lead_shape(site):
x, y = site.pos / L
return abs(y) < 0.5
# QPC shape: a rectangle with 2 gaussians
# etched out of the top and bottom edge.
def central_shape(site):
x, y = site.pos / L
return abs(x) < 3/5 and abs(y) < 0.5 - 0.4 * math.exp(-40 * x**2)
lead = kwant.Builder(kwant.TranslationalSymmetry(
model.lattice.vec((-1, 0))))
lead.fill(model, lead_shape, (0, 0))
syst = kwant.Builder()
syst.fill(model, central_shape, (0, 0))
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
return syst.finalized()
# Set up our model and system, and define the model parameters.
params = dict(alpha=0.365, beta=0.686, gamma=0.512, m=-0.01, U=0)
model = make_model(1)
syst = make_system(model, 70)
# Calculate the scattering states at energy 'm' coming from the left
# lead, and the associated particle current.
psi = kwant.wave_function(syst, energy=params['m'], params=params)(0)
.. jupyter-execute::
.. image:: ../../code/figure/plot_qahe_current.*
J = kwant.operator.Current(syst).bind(params=params)
current = sum(J(p) for p in psi)
kwant.plotter.current(syst, current);
From: :download:`QAHE example script <../../code/download/plot_qahe.py>`
From: :jupyter-download:script:`plot_qahe`
Scattering states with discrete symmetries and conservation laws
----------------------------------------------------------------
......
doc/source/pre/whatsnew/1.5.rst 0 → 100644
View file @ 39b65acd
What's new in Kwant 1.5
=======================
This article explains the user-visible changes in Kwant 1.5.0.
Deprecation of leaving 'norbs' unset for site families
------------------------------------------------------
When constructing site families (e.g. lattices) it is now deprecated to
leave the 'norbs' parameter unset. This will now raise a
KwantDeprecationWarning and will be disallowed in a future version of
Kwant. For example, when constructing a square lattice with 1 orbital
per site, use::
kwant.lattice.square(norbs=1)
rather than::
kwant.lattice.square()
Automatic addition of Peierls phase terms to Builders
-----------------------------------------------------
Kwant 1.4 introduced `kwant.physics.magnetic_gauge` for computing Peierls
phases for arbitrary geometries and for systems with leads. Using this
functionality requires that the system value functions are equipped to
take the required Peierls phase parameters, which is not possible when
you are not in direct control of the value functions (e.g. discretized
systems). In Kwant 1.5 we have added the missing piece of functionality,
`kwant.builder.add_peierls_phase`, which properly adds the Peierls phase
parameters to a Builder's value functions::
syst = make_system()
lead = make_lead()
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
fsyst, phase = kwant.builder.add_peierls_phase(syst)
def B_syst(pos):
return np.exp(-np.sum(pos * pos))
kwant.smatrix(fsyst, parameters=dict(t=-1, **phase(B_syst, 0, 0)))
Improved tutorial building
--------------------------
Improving or adding to Kwant's tutorial is now much simpler. Now
the text and code for each tutorial is kept in the same file, making
it easy to see where changes need to be made, and images generated by
the code are inserted directly into the document thanks to the magic of
`jupyter-sphinx <https://github.com/jupyter-widgets/jupyter-sphinx/>`_.
It has never been easier to get started contributing to Kwant by
helping us improve our documentation.
Discretization onto a Landau level basis
----------------------------------------
The Hamiltonian for a system infinite in at least two dimensions and with
a constant applied magnetic field may be expressed in a basis of Landau levels.
The momenta in the plane perpendicular to the magnetic field direction are
written in terms of the Landau level ladder operators:
.. math::
k_x = \sqrt{\frac{B}{2}} (a + a^\dagger) \quad\quad
k_y = i \sqrt{\frac{B}{2}} (a - a^\dagger)
The Hamiltonian is then expressed in terms of these ladder operators, which
allows for a straight-forward discretization in the basis of Landau levels,
provided that the basis is truncated after $N$ levels.
`kwant.continuum.discretize_landau` makes this procedure simple::
syst = kwant.continuum.discretize_landau("k_x**2 + k_y**2", N)
syst.finalized().hamiltonian_submatrix(params=dict(B=0.5))
`~kwant.continuum.discretize_landau` produces a regular Kwant Builder that
can be inspected or finalized as usual. 3D Hamiltonians for systems that
extend into the direction perpendicular to the magnetic field are also
possible::
template = kwant.continuum.discretize_landau("k_x**2 + k_y**2 + k_z**2 + V(z)", N)
This will produce a Builder with a single translational symmetry, which can be
directly finalized, or can be used as a template for e.g. a heterostructure stack
in the direction of the magnetic field::
def stack(site):
z, = site.pos
return 0 <= z < 10
template = kwant.continuum.discretize_landau("k_x**2 + k_y**2 + k_z**2 + V(z)", N)
syst = kwant.Builder()
syst.fill(template, stack, (0,))
doc/source/pre/whatsnew/index.rst
View file @ 39b65acd
......@@ -2,6 +2,7 @@ What's new in Kwant
===================
.. toctree::
1.5
1.4
1.3
1.2
......
doc/source/reference/kwant.builder.rst
View file @ 39b65acd
......@@ -26,3 +26,10 @@ Abstract base classes
SiteFamily
Symmetry
Lead
Functions
---------
.. autosummary::
:toctree: generated/
add_peierls_phase
doc/source/reference/kwant.continuum.rst
View file @ 39b65acd
......@@ -11,6 +11,7 @@ Discretizer
discretize
discretize_symbolic
build_discretized
discretize_landau
Symbolic helpers
----------------
......@@ -19,3 +20,11 @@ Symbolic helpers
sympify
lambdify
to_landau_basis
Other
-----
.. autosummary::
:toctree: generated/
LandauLattice
doc/source/tutorial/boilerplate.py 0 → 100644
View file @ 39b65acd
import matplotlib
import matplotlib.pyplot
from IPython.display import set_matplotlib_formats
matplotlib.rcParams['figure.figsize'] = matplotlib.pyplot.figaspect(1) * 2
set_matplotlib_formats('svg')
doc/source/tutorial/discretize.rst
View file @ 39b65acd
......@@ -16,10 +16,40 @@ The `~kwant.continuum` sub-package aims to be a solution to this problem.
It is a collection of tools for working with
continuum models and for discretizing them into tight-binding models.
.. seealso::
You can execute the code examples live in your browser by
activating thebelab:
.. thebe-button:: Activate Thebelab
.. seealso::
The complete source code of this tutorial can be found in
:download:`discretize.py </code/download/discretize.py>`
:jupyter-download:script:`discretize`
.. jupyter-kernel::
:id: discretize
.. jupyter-execute::
:hide-code:
# Tutorial 2.9. Processing continuum Hamiltonians with discretize
# ===============================================================
#
# Physics background
# ------------------
# - tight-binding approximation of continuous Hamiltonians
#
# Kwant features highlighted
# --------------------------
# - kwant.continuum.discretize
import matplotlib as mpl
from matplotlib import pyplot
import kwant
.. jupyter-execute:: boilerplate.py
:hide-code:
.. _tutorial_discretizer_introduction:
......@@ -55,9 +85,16 @@ Using `~kwant.continuum.discretize` to obtain a template
........................................................
First we must explicitly import the `kwant.continuum` package:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_import
:end-before: #HIDDEN_END_import
.. jupyter-execute::
import kwant.continuum
.. jupyter-execute::
:hide-code:
import scipy.sparse.linalg
import scipy.linalg
import numpy as np
The function `kwant.continuum.discretize` takes a symbolic Hamiltonian and
turns it into a `~kwant.builder.Builder` instance with appropriate spatial
......@@ -65,17 +102,16 @@ symmetry that serves as a template.
(We will see how to use the template to build systems with a particular
shape later).
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_symbolic_discretization
:end-before: #HIDDEN_END_symbolic_discretization
.. jupyter-execute::
template = kwant.continuum.discretize('k_x * A(x) * k_x')
print(template)
It is worth noting that ``discretize`` treats ``k_x`` and ``x`` as
non-commuting operators, and so their order is preserved during the
discretization process.
The builder produced by ``discretize`` may be printed to show the source code of its onsite and hopping functions (this is a special feature of builders returned by ``discretize``):
.. literalinclude:: /code/figure/discretizer_intro_verbose.txt
Printing the Builder produced by ``discretize`` shows the source code of its onsite and hopping functions (this is a special feature of builders returned by ``discretize``).
.. specialnote:: Technical details
......@@ -139,26 +175,45 @@ where :math:`V(x, y)` is some arbitrary potential.
First, use ``discretize`` to obtain a
builder that we will use as a template:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_template
:end-before: #HIDDEN_END_template
.. jupyter-execute::
hamiltonian = "k_x**2 + k_y**2 + V(x, y)"
template = kwant.continuum.discretize(hamiltonian)
print(template)
We now use this system with the `~kwant.builder.Builder.fill`
method of `~kwant.builder.Builder` to construct the system we
want to investigate:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_fill
:end-before: #HIDDEN_END_fill
.. jupyter-execute::
def stadium(site):
(x, y) = site.pos
x = max(abs(x) - 20, 0)
return x**2 + y**2 < 30**2
syst = kwant.Builder()
syst.fill(template, stadium, (0, 0));
syst = syst.finalized()
After finalizing this system, we can plot one of the system's
energy eigenstates:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_plot_eigenstate
:end-before: #HIDDEN_END_plot_eigenstate
.. jupyter-execute::
def plot_eigenstate(syst, n=2, Vx=.0003, Vy=.0005):
def potential(x, y):
return Vx * x + Vy * y
ham = syst.hamiltonian_submatrix(params=dict(V=potential), sparse=True)
evecs = scipy.sparse.linalg.eigsh(ham, k=10, which='SM')[1]
kwant.plotter.density(syst, abs(evecs[:, n])**2, show=False)
.. jupyter-execute::
:hide-code:
.. image:: /code/figure/discretizer_gs.*
plot_eigenstate(syst)
Note in the above that we pass the spatially varying potential *function*
to our system via a parameter called ``V``, because the symbol :math:`V`
......@@ -168,6 +223,8 @@ In addition, the function passed as ``V`` expects two input parameters ``x``
and ``y``, the same as in the initial continuum Hamiltonian.
.. _discretize-bhz-model:
Models with more structure: Bernevig-Hughes-Zhang
.................................................
......@@ -176,32 +233,95 @@ model [1]_ [2]_, one can provide matrix input to `~kwant.continuum.discretize`
using ``identity`` and ``kron``. For example, the definition of the BHZ model can be
written succinctly as:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_define_qsh
:end-before: #HIDDEN_END_define_qsh
.. jupyter-execute::
hamiltonian = """
+ C * identity(4) + M * kron(sigma_0, sigma_z)
- B * (k_x**2 + k_y**2) * kron(sigma_0, sigma_z)
- D * (k_x**2 + k_y**2) * kron(sigma_0, sigma_0)
+ A * k_x * kron(sigma_z, sigma_x)
- A * k_y * kron(sigma_0, sigma_y)
"""
a = 20
template = kwant.continuum.discretize(hamiltonian, grid=a)
We can then make a ribbon out of this template system:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_define_qsh_build
:end-before: #HIDDEN_END_define_qsh_build
.. jupyter-execute::
L, W = 2000, 1000
def shape(site):
(x, y) = site.pos
return (0 <= y < W and 0 <= x < L)
def lead_shape(site):
(x, y) = site.pos
return (0 <= y < W)
syst = kwant.Builder()
syst.fill(template, shape, (0, 0))
lead = kwant.Builder(kwant.TranslationalSymmetry([-a, 0]))
lead.fill(template, lead_shape, (0, 0))
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
syst = syst.finalized()
and plot its dispersion using `kwant.plotter.bands`:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_plot_qsh_band
:end-before: #HIDDEN_END_plot_qsh_band
.. jupyter-execute::
params = dict(A=3.65, B=-68.6, D=-51.1, M=-0.01, C=0)
kwant.plotter.bands(syst.leads[0], params=params,
momenta=np.linspace(-0.3, 0.3, 201), show=False)
.. image:: /code/figure/discretizer_qsh_band.*
pyplot.grid()
pyplot.xlim(-.3, 0.3)
pyplot.ylim(-0.05, 0.05)
pyplot.xlabel('momentum [1/A]')
pyplot.ylabel('energy [eV]')
pyplot.show()
In the above we see the edge states of the quantum spin Hall effect, which
we can visualize using `kwant.plotter.map`:
we can visualize using `kwant.plotter.density`:
.. jupyter-execute::
# get scattering wave functions at E=0
wf = kwant.wave_function(syst, energy=0, params=params)
# prepare density operators
sigma_z = np.array([[1, 0], [0, -1]])
prob_density = kwant.operator.Density(syst, np.eye(4))
spin_density = kwant.operator.Density(syst, np.kron(sigma_z, np.eye(2)))
# calculate expectation values and plot them
wf_sqr = sum(prob_density(psi) for psi in wf(0)) # states from left lead
rho_sz = sum(spin_density(psi) for psi in wf(0)) # states from left lead
fig, (ax1, ax2) = pyplot.subplots(1, 2, sharey=True, figsize=(16, 4))
kwant.plotter.density(syst, wf_sqr, ax=ax1)
kwant.plotter.density(syst, rho_sz, ax=ax2)
ax = ax1
im = [obj for obj in ax.get_children()
if isinstance(obj, mpl.image.AxesImage)][0]
fig.colorbar(im, ax=ax)
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_plot_qsh_wf
:end-before: #HIDDEN_END_plot_qsh_wf
ax = ax2
im = [obj for obj in ax.get_children()
if isinstance(obj, mpl.image.AxesImage)][0]
fig.colorbar(im, ax=ax)
.. image:: /code/figure/discretizer_qsh_wf.*
ax1.set_title('Probability density')
ax2.set_title('Spin density')
pyplot.show()
Limitations of discretization
......@@ -238,29 +358,62 @@ example. Let us start from the continuum Hamiltonian
We start by defining this model as a string and setting the value of the
:math:`α` parameter:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_ls_def
:end-before: #HIDDEN_END_ls_def
.. jupyter-execute::
hamiltonian = "k_x**2 * identity(2) + alpha * k_x * sigma_y"
params = dict(alpha=.5)
Now we can use `kwant.continuum.lambdify` to obtain a function that computes
:math:`H(k)`:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_ls_hk_cont
:end-before: #HIDDEN_END_ls_hk_cont
.. jupyter-execute::
h_k = kwant.continuum.lambdify(hamiltonian, locals=params)
k_cont = np.linspace(-4, 4, 201)
e_cont = [scipy.linalg.eigvalsh(h_k(k_x=ki)) for ki in k_cont]
We can also construct a discretized approximation using
`kwant.continuum.discretize`, in a similar manner to previous examples:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_ls_hk_tb
:end-before: #HIDDEN_END_ls_hk_tb
.. jupyter-execute::
def plot(ax, a=1):
template = kwant.continuum.discretize(hamiltonian, grid=a)
syst = kwant.wraparound.wraparound(template).finalized()
def h_k(k_x):
p = dict(k_x=k_x, **params)
return syst.hamiltonian_submatrix(params=p)
k_tb = np.linspace(-np.pi/a, np.pi/a, 201)
e_tb = [scipy.linalg.eigvalsh(h_k(k_x=a*ki)) for ki in k_tb]
ax.plot(k_cont, e_cont, 'r-')
ax.plot(k_tb, e_tb, 'k-')
ax.plot([], [], 'r-', label='continuum')
ax.plot([], [], 'k-', label='tight-binding')
ax.set_xlim(-4, 4)
ax.set_ylim(-1, 14)
ax.set_title('a={}'.format(a))
ax.set_xlabel('momentum [a.u.]')
ax.set_ylabel('energy [a.u.]')
ax.grid()
ax.legend()
Below we can see the continuum and tight-binding dispersions for two
different values of the discretization grid spacing :math:`a`:
.. image:: /code/figure/discretizer_lattice_spacing.*
.. jupyter-execute::
:hide-code:
_, (ax1, ax2) = pyplot.subplots(1, 2, sharey=True, figsize=(12, 4))
plot(ax1, a=1)
plot(ax2, a=.25)
pyplot.show()
We clearly see that the smaller grid spacing is, the better we approximate
the original continuous dispersion. It is also worth remembering that the
......@@ -284,20 +437,26 @@ It is possible to use ``identity`` (for identity matrix), ``kron`` (for Kronecke
expressions involving matrices. Matrices can also be provided explicitly using
square ``[]`` brackets. For example, all following expressions are equivalent:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_subs_1
:end-before: #HIDDEN_END_subs_1
.. jupyter-execute::
sympify = kwant.continuum.sympify
subs = {'sx': [[0, 1], [1, 0]], 'sz': [[1, 0], [0, -1]]}
e = (
sympify('[[k_x**2, alpha * k_x], [k_x * alpha, -k_x**2]]'),
sympify('k_x**2 * sigma_z + alpha * k_x * sigma_x'),
sympify('k_x**2 * sz + alpha * k_x * sx', locals=subs),
)
.. literalinclude:: /code/figure/discretizer_subs_1.txt
print(e[0] == e[1] == e[2])
We can use the ``locals`` keyword parameter to substitute expressions
and numerical values:
.. literalinclude:: /code/include/discretize.py
:start-after: #HIDDEN_BEGIN_subs_2
:end-before: #HIDDEN_END_subs_2
.. jupyter-execute::
.. literalinclude:: /code/figure/discretizer_subs_2.txt
subs = {'A': 'A(x) + B', 'V': 'V(x) + V_0', 'C': 5}
print(sympify('k_x * A * k_x + V + C', locals=subs))
Symbolic expressions obtained in this way can be directly passed to all
``discretizer`` functions.
......
doc/source/tutorial/faq.rst
View file @ 39b65acd
......@@ -5,6 +5,23 @@ questions that are discussed there. The `Kwant paper
<https://downloads.kwant-project.org/doc/kwant-paper.pdf>`_ also digs deeper
into Kwant's structure.
.. seealso::
You can execute the code examples live in your browser by
activating thebelab:
.. thebe-button:: Activate Thebelab
.. jupyter-execute::
:hide-code:
import kwant
import numpy as np
import tinyarray
import matplotlib
from matplotlib import pyplot as plt
.. jupyter-execute:: boilerplate.py
:hide-code:
What is a system, and what is a builder?
========================================
......@@ -48,11 +65,16 @@ combination of family and tag uniquely defines a site.
For example let us create an empty tight binding system and add two sites:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_site
:end-before: #HIDDEN_END_site
.. jupyter-execute::
a = 1
lat = kwant.lattice.square(a, norbs=1)
syst = kwant.Builder()
syst[lat(1, 0)] = 4
syst[lat(1, 1)] = 4
.. image:: /code/figure/faq_site.*
kwant.plot(syst);
In the above snippet we added 2 sites: ``lat(1, 0)`` and ``lat(0, 1)``. Both
of these sites belong to the same family, ``lat``, but have different tags:
......@@ -79,11 +101,11 @@ tuples, for example lists, are not treated as hoppings.
Starting from the example code from `What is a site?`_, we can add a hopping
to our system in the following way:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_hopping
:end-before: #HIDDEN_END_hopping
.. jupyter-execute::
.. image:: /code/figure/faq_hopping.*
syst[(lat(1, 0), lat(1, 1))] = 1j
kwant.plot(syst);
Visually, a hopping is represented as a line that joins two sites.
......@@ -139,18 +161,29 @@ Let's create two monatomic lattices (``lat_a`` and ``lat_b``). ``(1, 0)``
and ``(0, 1)`` will be the primitive vectors and ``(0, 0)`` and ``(0.5, 0.5)``
the origins of the two lattices:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_lattice_monatomic
:end-before: #HIDDEN_END_lattice_monatomic
.. jupyter-execute::
# Two monatomic lattices
primitive_vectors = [(1, 0), (0, 1)]
lat_a = kwant.lattice.Monatomic(primitive_vectors, offset=(0, 0), norbs=1)
lat_b = kwant.lattice.Monatomic(primitive_vectors, offset=(0.5, 0.5), norbs=1)
# lat1 is equivalent to kwant.lattice.square(norbs=1)
syst = kwant.Builder()
.. image:: /code/figure/faq_lattice.*
syst[lat_a(0, 0)] = 4
syst[lat_b(0, 0)] = 4
kwant.plot(syst);
We can also create a ``Polyatomic`` lattice with the same primitive vectors and
two sites in the basis:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_lattice_polyatomic
:end-before: #HIDDEN_END_lattice_polyatomic
.. jupyter-execute::
# One polyatomic lattice containing two sublattices
lat = kwant.lattice.Polyatomic([(1, 0), (0, 1)], [(0, 0), (0.5, 0.5)], norbs=1)
sub_a, sub_b = lat.sublattices
The two sublattices ``sub_a`` and ``sub_b`` are nothing else than ``Monatomic``
instances, and are equivalent to ``lat_a`` and ``lat_b`` that we created
......@@ -169,17 +202,40 @@ When plotting, how to color the different sublattices differently?
==================================================================
In the following example we shall use a kagome lattice, which has three sublattices.
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_colors1
:end-before: #HIDDEN_END_colors1
.. jupyter-execute::
lat = kwant.lattice.kagome(norbs=1)
syst = kwant.Builder()
a, b, c = lat.sublattices # The kagome lattice has 3 sublattices
As we can see below, we create a new plotting function that assigns a color for each family, and a different size for the hoppings depending on the family of the two sites. Finally we add sites and hoppings to our system and plot it with the new function.
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_colors2
:end-before: #HIDDEN_END_colors2
.. jupyter-execute::
# Plot sites from different families in different colors
def family_color(site):
if site.family == a:
return 'red'
if site.family == b:
return 'green'
else:
return 'blue'
def plot_system(syst):
kwant.plot(syst, site_lw=0.1, site_color=family_color)
## Add sites and hoppings.
for i in range(4):
for j in range (4):
syst[a(i, j)] = 4
syst[b(i, j)] = 4
syst[c(i, j)] = 4
syst[lat.neighbors()] = -1
.. image:: /code/figure/faq_colors.*
## Plot the system.
plot_system(syst)
How to create many similar hoppings in one go?
......@@ -197,31 +253,64 @@ integers).
The following example shows how this can be used:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_direction1
:end-before: #HIDDEN_END_direction1
.. jupyter-execute::
.. image:: /code/figure/faq_direction1.*
# Create hopping between neighbors with HoppingKind
a = 1
syst = kwant.Builder()
lat = kwant.lattice.square(a, norbs=1)
syst[ (lat(i, j) for i in range(5) for j in range(5)) ] = 4
syst[kwant.builder.HoppingKind((1, 0), lat)] = -1
kwant.plot(syst);
Note that ``HoppingKind`` only works with site families so you cannot use
them directly with ``Polyatomic`` lattices; you have to explicitly specify
the sublattices when creating a ``HoppingKind``:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_direction2
:end-before: #HIDDEN_END_direction2
.. jupyter-execute::
:hide-code:
lat = kwant.lattice.kagome(norbs=1)
syst = kwant.Builder()
a, b, c = lat.sublattices # The kagome lattice has 3 sublattices
for i in range(4):
for j in range (4):
syst[a(i, j)] = 4
syst[b(i, j)] = 4
syst[c(i, j)] = 4
.. jupyter-execute::
# equivalent to syst[kwant.builder.HoppingKind((0, 1), b)] = -1
syst[kwant.builder.HoppingKind((0, 1), b, b)] = -1
Here, we want the hoppings between the sites from sublattice b with a direction of (0,1) in the lattice coordinates.
.. image:: /code/figure/faq_direction2.*
.. jupyter-execute::
:hide-code:
plot_system(syst)
.. jupyter-execute::
:hide-code:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_direction3
:end-before: #HIDDEN_END_direction3
# Delete the hoppings previously created
del syst[kwant.builder.HoppingKind((0, 1), b, b)]
.. jupyter-execute::
syst[kwant.builder.HoppingKind((0, 0), a, b)] = -1
syst[kwant.builder.HoppingKind((0, 0), a, c)] = -1
syst[kwant.builder.HoppingKind((0, 0), c, b)] = -1
Here, we create hoppings between the sites of the same lattice coordinates but from different families.
.. image:: /code/figure/faq_direction3.*
.. jupyter-execute::
plot_system(syst)
How to set the hoppings between adjacent sites?
......@@ -230,32 +319,53 @@ How to set the hoppings between adjacent sites?
that returns a list of ``HoppingKind`` instances that connect sites with their
(n-nearest) neighors:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_adjacent1
:end-before: #HIDDEN_END_adjacent1
.. jupyter-execute::
.. image:: /code/figure/faq_adjacent1.*
.. image:: /code/figure/faq_adjacent2.*
# Create hoppings with lat.neighbors()
syst = kwant.Builder()
lat = kwant.lattice.square(norbs=1)
syst[(lat(i, j) for i in range(3) for j in range(3))] = 4
As we can see in the figure above, ``lat.neighbors()`` (on the left) returns the hoppings between the first nearest neighbors and ``lat.neighbors(2)`` (on the right) returns the hoppings between the second nearest neighbors.
syst[lat.neighbors()] = -1 # Equivalent to lat.neighbors(1)
kwant.plot(syst);
.. jupyter-execute::
del syst[lat.neighbors()] # Delete all nearest-neighbor hoppings
syst[lat.neighbors(2)] = -1
kwant.plot(syst);
As we can see in the figures above, ``lat.neighbors()`` returns the hoppings between the first nearest neighbors and ``lat.neighbors(2)`` returns the hoppings between the second nearest neighbors.
When using a ``Polyatomic`` lattice ``neighbors()`` knows about the different
sublattices:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_adjacent2
:end-before: #HIDDEN_END_adjacent2
.. jupyter-execute::
# Create the system
lat = kwant.lattice.kagome(norbs=1)
syst = kwant.Builder()
a, b, c = lat.sublattices # The kagome lattice has 3 sublattices
.. image:: /code/figure/faq_adjacent3.*
for i in range(4):
for j in range (4):
syst[a(i, j)] = 4 # red
syst[b(i, j)] = 4 # green
syst[c(i, j)] = 4 # blue
syst[lat.neighbors()] = -1
plot_system(syst)
However, if we use the ``neighbors()`` method of a single sublattice, we will
only get the neighbors *on that sublattice*:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_adjacent3
:end-before: #HIDDEN_END_adjacent3
.. jupyter-execute::
.. image:: /code/figure/faq_adjacent4.*
del syst[lat.neighbors()] # Delete the hoppings previously created
syst[a.neighbors()] = -1
plot_system(syst)
Note in the above that there are *only* hoppings between the red sites. This
is an artifact of the visualisation: the blue and green sites just happen to lie
......@@ -268,12 +378,35 @@ To make a hole in the system, use ``del syst[site]``, just like with any other
mapping. In the following example we remove all sites inside some "hole"
region:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_hole
:end-before: #HIDDEN_END_hole
.. jupyter-execute::
# Define the lattice and the (empty) system
a = 2
lat = kwant.lattice.cubic(a, norbs=1)
syst = kwant.Builder()
L = 10
W = 10
H = 2
# Add sites to the system in a cuboid
syst[(lat(i, j, k) for i in range(L) for j in range(W) for k in range(H))] = 4
kwant.plot(syst);
.. jupyter-execute::
# Delete sites to create a hole
def in_hole(site):
x, y, z = site.pos / a - (L/2, W/2, H/2) # position relative to centre
return abs(x) < L / 4 and abs(y) < W / 4
for site in filter(in_hole, list(syst.sites())):
del syst[site]
kwant.plot(syst);
.. image:: /code/figure/faq_hole1.*
.. image:: /code/figure/faq_hole2.*
``del syst[site]`` also works after hoppings have been added to the system.
If a site is deleted, then all the hoppings to/from that site are also deleted.
......@@ -287,9 +420,18 @@ what is a builder?`_ if you do not know the difference).
We can access the sites of a ``Builder`` by using its `~kwant.builder.Builder.sites` method:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_sites1
:end-before: #HIDDEN_END_sites1
.. jupyter-execute::
:hide-code:
builder = kwant.Builder()
lat = kwant.lattice.square(norbs=1)
builder[(lat(i, j) for i in range(3) for j in range(3))] = 4
.. jupyter-execute::
# Before finalizing the system
sites = list(builder.sites()) # sites() doe *not* return a list
The ``sites()`` method returns an *iterator* over the system sites, and in the
above example we create a list from the contents of this iterator, which
......@@ -300,9 +442,11 @@ well be returned in a different order.
After finalization, when we are dealing with a ``System``, the sites themselves
are stored in a list, which can be accessed via the ``sites`` attribute:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_sites2
:end-before: #HIDDEN_END_sites2
.. jupyter-execute::
# After finalizing the system
syst = builder.finalized()
sites = syst.sites # syst.sites is an actual list
The order of sites in a ``System`` is fixed, and also defines the ordering of
the system Hamiltonian, system wavefunctions etc. (see `How does Kwant order components of an individual wavefunction?`_ for details).
......@@ -310,9 +454,9 @@ the system Hamiltonian, system wavefunctions etc. (see `How does Kwant order com
``System`` also contains the inverse mapping, ``id_by_site`` which gives us
the index of a given site within the system:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_sites3
:end-before: #HIDDEN_END_sites3
.. jupyter-execute::
i = syst.id_by_site[lat(0, 2)] # we want the id of the site lat(0, 2)
How to use different lattices for the scattering region and a lead?
......@@ -322,19 +466,35 @@ which we want to connect to leads that contain sites from a square lattice.
First we construct the central system:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_different_lattice1
:end-before: #HIDDEN_END_different_lattice1
.. jupyter-execute::
# Define the scattering Region
L = 5
W = 5
lat = kwant.lattice.honeycomb(norbs=1)
subA, subB = lat.sublattices
.. image:: /code/figure/faq_different_lattice1.*
syst = kwant.Builder()
syst[(subA(i, j) for i in range(L) for j in range(W))] = 4
syst[(subB(i, j) for i in range(L) for j in range(W))] = 4
syst[lat.neighbors()] = -1
kwant.plot(syst);
and the lead:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_different_lattice2
:end-before: #HIDDEN_END_different_lattice2
.. jupyter-execute::
# Create a lead
lat_lead = kwant.lattice.square(norbs=1)
sym_lead1 = kwant.TranslationalSymmetry((0, 1))
.. image:: /code/figure/faq_different_lattice2.*
lead1 = kwant.Builder(sym_lead1)
lead1[(lat_lead(i, 0) for i in range(2, 7))] = 4
lead1[lat_lead.neighbors()] = -1
kwant.plot(syst);
We cannot simply use `~kwant.builder.Builder.attach_lead` to attach this lead to the
system with the honeycomb lattice because Kwant does not know how sites from
......@@ -343,19 +503,22 @@ these two lattices should be connected.
We must first add a layer of sites from the square lattice to the system and manually
add the hoppings from these sites to the sites from the honeycomb lattice:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_different_lattice3
:end-before: #HIDDEN_END_different_lattice3
.. jupyter-execute::
syst[(lat_lead(i, 5) for i in range(2, 7))] = 4
syst[lat_lead.neighbors()] = -1
# Manually attach sites from graphene to square lattice
syst[((lat_lead(i+2, 5), subB(i, 4)) for i in range(5))] = -1
.. image:: /code/figure/faq_different_lattice3.*
kwant.plot(syst);
``attach_lead()`` will now be able to attach the lead:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_different_lattice4
:end-before: #HIDDEN_END_different_lattice4
.. jupyter-execute::
.. image:: /code/figure/faq_different_lattice4.*
syst.attach_lead(lead1)
kwant.plot(syst);
How to cut a finite system out of a system with translational symmetries?
......@@ -369,30 +532,50 @@ will be repeated in the desired shape to create the final system.
For example, say we want to create a simple model on a cubic lattice:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_fill1
:end-before: #HIDDEN_END_fill1
.. jupyter-execute::
# Create 3d model.
cubic = kwant.lattice.cubic(norbs=1)
sym_3d = kwant.TranslationalSymmetry([1, 0, 0], [0, 1, 0], [0, 0, 1])
model = kwant.Builder(sym_3d)
model[cubic(0, 0, 0)] = 4
model[cubic.neighbors()] = -1
We have now created our "template" ``Builder`` which has 3 translational
symmetries. Next we will fill a system with no translational symmetries with
sites and hoppings from the template inside a cuboid:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_fill2
:end-before: #HIDDEN_END_fill2
.. jupyter-execute::
# Build scattering region (white).
def cuboid_shape(site):
x, y, z = abs(site.pos)
return x < 4 and y < 10 and z < 3
.. image:: /code/figure/faq_fill2.*
cuboid = kwant.Builder()
cuboid.fill(model, cuboid_shape, (0, 0, 0));
kwant.plot(cuboid);
We can then use the original template to create a lead, which has 1 translational
symmetry. We can then use this lead as a template to fill another section of
the system with a cylinder of sites and hoppings:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_fill3
:end-before: #HIDDEN_END_fill3
.. jupyter-execute::
# Build electrode (black).
def electrode_shape(site):
x, y, z = site.pos - (0, 5, 2)
return y**2 + z**2 < 2.3**2
.. image:: /code/figure/faq_fill3.*
electrode = kwant.Builder(kwant.TranslationalSymmetry([1, 0, 0]))
electrode.fill(model, electrode_shape, (0, 5, 2)) # lead
# Scattering region
cuboid.fill(electrode, lambda s: abs(s.pos[0]) < 7, (0, 5, 4))
cuboid.attach_lead(electrode)
kwant.plot(cuboid);
How does Kwant order the propagating modes of a lead?
=====================================================
......@@ -402,16 +585,44 @@ achieved with the `~kwant.system.InfiniteSystem.modes` method, which returns a
pair of objects, the first of which contains the propagating modes of the
system in a `~kwant.physics.PropagatingModes` object:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_pm
:end-before: #HIDDEN_END_pm
.. jupyter-execute::
lat = kwant.lattice.square(norbs=1)
lead = kwant.Builder(kwant.TranslationalSymmetry((-1, 0)))
lead[(lat(0, i) for i in range(3))] = 4
lead[lat.neighbors()] = -1
flead = lead.finalized()
E = 2.5
prop_modes, _ = flead.modes(energy=E)
``PropagatingModes`` contains the wavefunctions, velocities and momenta of the
modes at the requested energy (2.5 in this example). In order to understand
the order in which these quantities are returned it is often useful to look at
the a section of the band structure for the system in question:
.. image:: /code/figure/faq_pm1.*
.. jupyter-execute::
:hide-code:
def plot_and_label_modes(lead, E):
# Plot the different modes
pmodes, _ = lead.modes(energy=E)
kwant.plotter.bands(lead, show=False)
for i, k in enumerate(pmodes.momenta):
plt.plot(k, E, 'ko')
plt.annotate(str(i), xy=(k, E), xytext=(-5, 8),
textcoords='offset points',
bbox=dict(boxstyle='round,pad=0.1',fc='white', alpha=0.7))
plt.plot([-3, 3], [E, E], 'r--')
plt.ylim(E-1, E+1)
plt.xlim(-2, 2)
plt.xlabel("momentum")
plt.ylabel("energy")
plt.show()
plot_and_label_modes(flead, E)
On the above band structure we have labelled the 4 modes in the order
that they appear in the output of ``modes()`` at energy 2.5. Note that
......@@ -425,7 +636,22 @@ the modes are sorted in the following way:
For more complicated systems and band structures this can lead to some
unintuitive orderings:
.. image:: /code/figure/faq_pm2.*
.. jupyter-execute::
:hide-code:
lat = kwant.lattice.square(norbs=2)
s0 = np.eye(2)
sz = np.array([[1, 0], [0, -1]])
lead2 = kwant.Builder(kwant.TranslationalSymmetry((-1, 0)))
lead2[(lat(0, i) for i in range(2))] = np.diag([1.8, -1])
lead2[lat.neighbors()] = -1 * sz
flead2 = lead2.finalized()
plot_and_label_modes(flead2, 1)
How does Kwant order scattering states?
......@@ -452,10 +678,38 @@ When all the site families present in a system have only 1 degree of freedom
per site (i.e. the all the onsites are scalars) then the index into a
wavefunction defined over the system is exactly the site index:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_ord1
:end-before: #HIDDEN_END_ord1
.. literalinclude:: /code/figure/faq_ord1.txt
.. jupyter-execute::
:hide-code:
def circle(R):
return lambda r: np.linalg.norm(r) < R
def make_system(lat):
norbs = lat.norbs
syst = kwant.Builder()
syst[lat.shape(circle(3), (0, 0))] = 4 * np.eye(norbs)
syst[lat.neighbors()] = -1 * np.eye(norbs)
lead = kwant.Builder(kwant.TranslationalSymmetry((-1, 0)))
lead[(lat(0, i) for i in range(-1, 2))] = 4 * np.eye(norbs)
lead[lat.neighbors()] = -1 * np.eye(norbs)
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
return syst.finalized()
.. jupyter-execute::
lat = kwant.lattice.square(norbs=1)
syst = make_system(lat)
scattering_states = kwant.wave_function(syst, energy=1)
wf = scattering_states(0)[0] # scattering state from lead 0 incoming in mode 0
idx = syst.id_by_site[lat(0, 0)] # look up index of site
print('wavefunction on lat(0, 0): ', wf[idx])
We see that the wavefunction on a single site is a single complex number, as
expected.
......@@ -467,10 +721,20 @@ to one another. In the case where all site families in the system have the
wavefunction into a matrix, where the row number indexes the site, and the
column number indexes the degree of freedom on that site:
.. literalinclude:: /code/include/faq.py
:start-after: #HIDDEN_BEGIN_ord2
:end-before: #HIDDEN_END_ord2
.. literalinclude:: /code/figure/faq_ord2.txt
.. jupyter-execute::
lat = kwant.lattice.square(norbs=2)
syst = make_system(lat)
scattering_states = kwant.wave_function(syst, energy=1)
wf = scattering_states(0)[0] # scattering state from lead 0 incoming in mode 0
idx = syst.id_by_site[lat(0, 0)] # look up index of site
# Group consecutive degrees of freedom from 'wf' together; these correspond
# to degrees of freedom on the same site.
wf = wf.reshape(-1, 2)
print('wavefunction on lat(0, 0): ', wf[idx])
We see that the wavefunction on a single site is a *vector* of 2 complex numbers,
as we expect.
......
doc/source/tutorial/first_steps.rst
View file @ 39b65acd
......@@ -13,7 +13,7 @@ Schrödinger equation
.. math::
H = \frac{-\hbar^2}{2m}(\partial_x^2 + \partial_y^2) + V(y)
with a hard-wall confinement :math:`V(y)` in y-direction.
with a hard-wall confinement :math:`V(y)` in the y-direction.
To be able to implement the quantum wire with Kwant, the continuous Hamiltonian
:math:`H` has to be discretized thus turning it into a tight-binding
......@@ -60,25 +60,58 @@ simplicity, we choose to work in such units that :math:`t = a = 1`.
Transport through a quantum wire
................................
.. seealso::
You can execute the code examples live in your browser by
activating thebelab:
.. thebe-button:: Activate Thebelab
.. seealso::
The complete source code of this example can be found in
:download:`quantum_wire.py </code/download/quantum_wire.py>`
:jupyter-download:script:`quantum_wire`
In order to use Kwant, we need to import it:
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_dwhx
:end-before: #HIDDEN_END_dwhx
.. jupyter-kernel::
:id: quantum_wire
.. jupyter-execute::
:hide-code:
# Tutorial 2.2.2. Transport through a quantum wire
# ================================================
#
# Physics background
# ------------------
# Conductance of a quantum wire; subbands
#
# Kwant features highlighted
# --------------------------
# - Builder for setting up transport systems easily
# - Making scattering region and leads
# - Using the simple sparse solver for computing Landauer conductance
# For plotting
from matplotlib import pyplot
.. jupyter-execute:: boilerplate.py
:hide-code:
.. jupyter-execute::
import kwant
Enabling Kwant is as easy as this [#]_ !
The first step is now the definition of the system with scattering region and
The first step is now to define the system with scattering region and
leads. For this we make use of the `~kwant.builder.Builder` type that allows to
define a system in a convenient way. We need to create an instance of it:
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_goiq
:end-before: #HIDDEN_END_goiq
.. jupyter-execute::
# First define the tight-binding system
syst = kwant.Builder()
Observe that we just accessed `~kwant.builder.Builder` by the name
``kwant.Builder``. We could have just as well written
......@@ -92,16 +125,19 @@ Apart from `~kwant.builder.Builder` we also need to specify
what kind of sites we want to add to the system. Here we work with
a square lattice. For simplicity, we set the lattice constant to unity:
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_suwo
:end-before: #HIDDEN_END_suwo
.. jupyter-execute::
a = 1
lat = kwant.lattice.square(a, norbs=1)
Since we work with a square lattice, we label the points with two
integer coordinates `(i, j)`. `~kwant.builder.Builder` then
allows us to add matrix elements corresponding to lattice points:
``syst[lat(i, j)] = ...`` sets the on-site energy for the point `(i, j)`,
and ``syst[lat(i1, j1), lat(i2, j2)] = ...`` the hopping matrix element
**from** point `(i2, j2)` **to** point `(i1, j1)`.
**from** point `(i2, j2)` **to** point `(i1, j1)`. In specifying ``norbs=1``
in the definition of the lattice we tell Kwant that there is 1 degree
of freedom per lattice site.
Note that we need to specify sites for `~kwant.builder.Builder`
in the form ``lat(i, j)``. The lattice object `lat` does the
......@@ -111,9 +147,25 @@ needed in Builder (more about that in the technical details below).
We now build a rectangular scattering region that is `W`
lattice points wide and `L` lattice points long:
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_zfvr
:end-before: #HIDDEN_END_zfvr
.. jupyter-execute::
t = 1.0
W, L = 10, 30
# Define the scattering region
for i in range(L):
for j in range(W):
# On-site Hamiltonian
syst[lat(i, j)] = 4 * t
# Hopping in y-direction
if j > 0:
syst[lat(i, j), lat(i, j - 1)] = -t
# Hopping in x-direction
if i > 0:
syst[lat(i, j), lat(i - 1, j)] = -t
Observe how the above code corresponds directly to the terms of the discretized
Hamiltonian:
......@@ -140,9 +192,14 @@ Next, we define the leads. Leads are also constructed using
`~kwant.builder.Builder`, but in this case, the
system must have a translational symmetry:
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_xcmc
:end-before: #HIDDEN_END_xcmc
.. jupyter-execute::
# Then, define and attach the leads:
# First the lead to the left
# (Note: TranslationalSymmetry takes a real-space vector)
sym_left_lead = kwant.TranslationalSymmetry((-a, 0))
left_lead = kwant.Builder(sym_left_lead)
Here, the `~kwant.builder.Builder` takes a
`~kwant.lattice.TranslationalSymmetry` as the optional parameter. Note that the
......@@ -155,17 +212,22 @@ as the hoppings inside one unit cell and to the next unit cell of the lead.
For a square lattice, and a lead in y-direction the unit cell is
simply a vertical line of points:
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_ndez
:end-before: #HIDDEN_END_ndez
.. jupyter-execute::
for j in range(W):
left_lead[lat(0, j)] = 4 * t
if j > 0:
left_lead[lat(0, j), lat(0, j - 1)] = -t
left_lead[lat(1, j), lat(0, j)] = -t
Note that here it doesn't matter if you add the hoppings to the next or the
previous unit cell -- the translational symmetry takes care of that. The
isolated, infinite is attached at the correct position using
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_fskr
:end-before: #HIDDEN_END_fskr
.. jupyter-execute::
:hide-output:
syst.attach_lead(left_lead)
This call returns the lead number which will be used to refer to the lead when
computing transmissions (further down in this tutorial). More details about
......@@ -175,9 +237,21 @@ We also want to add a lead on the right side. The only difference to
the left lead is that the vector of the translational
symmetry must point to the right, the remaining code is the same:
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_xhqc
:end-before: #HIDDEN_END_xhqc
.. jupyter-execute::
:hide-output:
# Then the lead to the right
sym_right_lead = kwant.TranslationalSymmetry((a, 0))
right_lead = kwant.Builder(sym_right_lead)
for j in range(W):
right_lead[lat(0, j)] = 4 * t
if j > 0:
right_lead[lat(0, j), lat(0, j - 1)] = -t
right_lead[lat(1, j), lat(0, j)] = -t
syst.attach_lead(right_lead)
Note that here we added points with x-coordinate 0, just as for the left lead.
You might object that the right lead should be placed `L`
......@@ -187,13 +261,9 @@ you do not need to worry about that.
Now we have finished building our system! We plot it, to make sure we didn't
make any mistakes:
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_wsgh
:end-before: #HIDDEN_END_wsgh
This should bring up this picture:
.. jupyter-execute::
.. image:: /code/figure/quantum_wire_syst.*
kwant.plot(syst);
The system is represented in the usual way for tight-binding systems:
dots represent the lattice points `(i, j)`, and for every
......@@ -203,16 +273,29 @@ fading color.
In order to use our system for a transport calculation, we need to finalize it
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_dngj
:end-before: #HIDDEN_END_dngj
.. jupyter-execute::
# Finalize the system
syst = syst.finalized()
Having successfully created a system, we now can immediately start to compute
its conductance as a function of energy:
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_buzn
:end-before: #HIDDEN_END_buzn
.. jupyter-execute::
# Now that we have the system, we can compute conductance
energies = []
data = []
for ie in range(100):
energy = ie * 0.01
# compute the scattering matrix at a given energy
smatrix = kwant.smatrix(syst, energy)
# compute the transmission probability from lead 0 to
# lead 1
energies.append(energy)
data.append(smatrix.transmission(1, 0))
We use ``kwant.smatrix`` which is a short name for
`kwant.solvers.default.smatrix` of the default solver module
......@@ -227,13 +310,15 @@ Finally we can use ``matplotlib`` to make a plot of the computed data
(although writing to file and using an external viewer such as
gnuplot or xmgrace is just as viable)
.. literalinclude:: /code/include/quantum_wire.py
:start-after: #HIDDEN_BEGIN_lliv
:end-before: #HIDDEN_END_lliv
.. jupyter-execute::
This should yield the result
.. image:: /code/figure/quantum_wire_result.*
# Use matplotlib to write output
# We should see conductance steps
pyplot.figure()
pyplot.plot(energies, data)
pyplot.xlabel("energy [t]")
pyplot.ylabel("conductance [e^2/h]")
pyplot.show()
We see a conductance quantized in units of :math:`e^2/h`,
increasing in steps as the energy is increased. The
......@@ -333,7 +418,12 @@ Building the same system with less code
.. seealso::
The complete source code of this example can be found in
:download:`quantum_wire_revisited.py </code/download/quantum_wire_revisited.py>`
:jupyter-download:script:`quantum_wire_revisited`
.. jupyter-kernel::
:id: quantum_wire_revisited
Kwant allows for more than one way to build a system. The reason is that
`~kwant.builder.Builder` is essentially just a container that can be filled in
......@@ -347,19 +437,57 @@ the code into separate entities. In this example we therefore also aim at more
structure.
We begin the program collecting all imports in the beginning of the
file and put the build-up of the system into a separate function
`make_system`:
file and defining the a square lattice and empty scattering region.
.. jupyter-execute::
:hide-code:
# Tutorial 2.2.3. Building the same system with less code
# =======================================================
#
# Physics background
# ------------------
# Conductance of a quantum wire; subbands
#
# Kwant features highlighted
# --------------------------
# - Using iterables and builder.HoppingKind for making systems
# - introducing `reversed()` for the leads
#
# Note: Does the same as tutorial1a.py, but using other features of Kwant.
.. jupyter-execute::
import kwant
# For plotting
from matplotlib import pyplot
.. jupyter-execute:: boilerplate.py
:hide-code:
.. jupyter-execute::
a = 1
t = 1.0
W, L = 10, 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).
# Each lattice site has 1 degree of freedom, hence norbs=1.
lat = kwant.lattice.square(a, norbs=1)
syst = kwant.Builder()
.. literalinclude:: /code/include/quantum_wire_revisited.py
:start-after: #HIDDEN_BEGIN_xkzy
:end-before: #HIDDEN_END_xkzy
Previously, the scattering region was build using two ``for``-loops.
Instead, we now write:
.. literalinclude:: /code/include/quantum_wire_revisited.py
:start-after: #HIDDEN_BEGIN_vvjt
:end-before: #HIDDEN_END_vvjt
.. jupyter-execute::
syst[(lat(x, y) for x in range(L) for y in range(W))] = 4 * t
Here, all lattice points are added at once in the first line. The
construct ``((i, j) for i in range(L) for j in range(W))`` is a
......@@ -375,9 +503,9 @@ hoppings. In this case, an iterable like for the lattice
points becomes a bit cumbersome, and we use instead another
feature of Kwant:
.. literalinclude:: /code/include/quantum_wire_revisited.py
:start-after: #HIDDEN_BEGIN_nooi
:end-before: #HIDDEN_END_nooi
.. jupyter-execute::
syst[lat.neighbors()] = -t
In regular lattices, hoppings form large groups such that hoppings within a
group can be transformed into one another by lattice translations. In order to
......@@ -397,49 +525,48 @@ values for all the nth-nearest neighbors at once, one can similarly use
The left lead is constructed in an analogous way:
.. literalinclude:: /code/include/quantum_wire_revisited.py
:start-after: #HIDDEN_BEGIN_iepx
:end-before: #HIDDEN_END_iepx
.. jupyter-execute::
lead = kwant.Builder(kwant.TranslationalSymmetry((-a, 0)))
lead[(lat(0, j) for j in range(W))] = 4 * t
lead[lat.neighbors()] = -t
The previous example duplicated almost identical code for the left and
the right lead. The only difference was the direction of the translational
symmetry vector. Here, we only construct the left lead, and use the method
`~kwant.builder.Builder.reversed` of `~kwant.builder.Builder` to obtain a copy
of a lead pointing in the opposite direction. Both leads are attached as
before and the finished system returned:
before:
.. literalinclude:: /code/include/quantum_wire_revisited.py
:start-after: #HIDDEN_BEGIN_yxot
:end-before: #HIDDEN_END_yxot
.. jupyter-execute::
:hide-output:
The remainder of the script has been organized into two functions. One for the
plotting of the conductance.
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
.. literalinclude:: /code/include/quantum_wire_revisited.py
:start-after: #HIDDEN_BEGIN_ayuk
:end-before: #HIDDEN_END_ayuk
The remainder of the script proceeds identically. We first finalize the system:
And one ``main`` function.
.. jupyter-execute::
.. literalinclude:: /code/include/quantum_wire_revisited.py
:start-after: #HIDDEN_BEGIN_cjel
:end-before: #HIDDEN_END_cjel
syst = syst.finalized()
Finally, we use the following standard Python construct [#]_ to execute
``main`` if the program is used as a script (i.e. executed as
``python quantum_wire_revisited.py``):
and then calculate the transmission and plot:
.. literalinclude:: /code/include/quantum_wire_revisited.py
:start-after: #HIDDEN_BEGIN_ypbj
:end-before: #HIDDEN_END_ypbj
.. jupyter-execute::
If the example, however, is imported inside Python using ``import
quantum_wire_revisted as qw``, ``main`` is not executed automatically.
Instead, you can execute it manually using ``qw.main()``. On the other
hand, you also have access to the other functions, ``make_system`` and
``plot_conductance``, and can thus play with the parameters.
energies = []
data = []
for ie in range(100):
energy = ie * 0.01
smatrix = kwant.smatrix(syst, energy)
energies.append(energy)
data.append(smatrix.transmission(1, 0))
The result of the example should be identical to the previous one.
pyplot.figure()
pyplot.plot(energies, data)
pyplot.xlabel("energy [t]")
pyplot.ylabel("conductance [e^2/h]")
pyplot.show()
.. specialnote:: Technical details
......@@ -453,11 +580,9 @@ The result of the example should be identical to the previous one.
For technical reasons it is not possible to add several points
using a tuple of sites. Hence it is worth noting
a subtle detail in
a subtle detail in:
.. literalinclude:: /code/include/quantum_wire_revisited.py
:start-after: #HIDDEN_BEGIN_vvjt
:end-before: #HIDDEN_END_vvjt
>>> syst[(lat(x, y) for x in range(L) for y in range(W))] = 4 * t
Note that ``(lat(x, y) for x in range(L) for y in range(W))`` is not
a tuple, but a generator.
......@@ -512,6 +637,132 @@ The result of the example should be identical to the previous one.
it would be impossible to distinguish whether one would like to add two
separate sites, or one hopping.
Tips for organizing your simulation scripts
...........................................
.. seealso::
The complete source code of this example can be found in
:jupyter-download:script:`quantum_wire_organized`
.. jupyter-kernel::
:id: quantum_wire_organized
.. jupyter-execute::
:hide-code:
# Tutorial 2.2.4. Organizing a simulation script
# ==============================================
#
# Physics background
# ------------------
# Conductance of a quantum wire; subbands
#
# Note: Does the same as quantum_write_revisited.py, but features
# better code organization
The above two examples illustrate some of the core features of Kwant, however
the code was presented in a style which is good for exposition, but which is
bad for making your code understandable and reusable. In this example we will
lay out some best practices for writing your own simulation scripts.
In the above examples we constructed a single Kwant system, using global variables
for parameters such as the lattice constant and the length and width of the system.
Instead, it is preferable to create a *function* that you can call, and which will
return a Kwant ``Builder``:
.. jupyter-execute::
from matplotlib import pyplot
import kwant
.. jupyter-execute:: boilerplate.py
:hide-code:
.. jupyter-execute::
def make_system(L, W, a=1, t=1.0):
lat = kwant.lattice.square(a, norbs=1)
syst = kwant.Builder()
syst[(lat(i, j) for i in range(L) for j in range(W))] = 4 * t
syst[lat.neighbors()] = -t
lead = kwant.Builder(kwant.TranslationalSymmetry((-a, 0)))
lead[(lat(0, j) for j in range(W))] = 4 * t
lead[lat.neighbors()] = -t
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
return syst
By encapsulating system creation within ``make_system`` we *document* our code
by telling readers that *this* is how we create a system, and that creating a system
depends on *these* parameters (the length and width of the system, in this case, as well
as the lattice constant and the value for the hopping parameter). By defining a function
we also ensure that we can consistently create different systems (e.g. of different sizes)
of the same type (rectangular slab).
We similarly encapsulate the part of the script that does computation and plotting into
a function ``plot_conductance``:
.. jupyter-execute::
def plot_conductance(syst, energies):
# Compute conductance
data = []
for energy in energies:
smatrix = kwant.smatrix(syst, energy)
data.append(smatrix.transmission(1, 0))
pyplot.figure()
pyplot.plot(energies, data)
pyplot.xlabel("energy [t]")
pyplot.ylabel("conductance [e^2/h]")
pyplot.show()
And the ``main`` function that glues together the components that we previously defined:
.. jupyter-execute::
def main():
syst = make_system(W=10, L=30)
# Check that the system looks as intended.
kwant.plot(syst)
# Finalize the system.
fsyst = syst.finalized()
# We should see conductance steps.
plot_conductance(fsyst, energies=[0.01 * i for i in range(100)])
Finally, we use the following standard Python construct [#]_ to execute
``main`` if the program is used as a script (i.e. executed as
``python quantum_wire_organized.py``):
.. jupyter-execute::
# 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()
If the example, however, is imported inside Python using ``import
quantum_wire_organized as qw``, ``main`` is not executed automatically.
Instead, you can execute it manually using ``qw.main()``. On the other
hand, you also have access to the other functions, ``make_system`` and
``plot_conductance``, and can thus play with the parameters.
The result of this example should be identical to the previous one.
.. rubric:: Footnotes
.. [#] https://docs.python.org/3/library/__main__.html
doc/source/tutorial/graphene.rst
View file @ 39b65acd
......@@ -3,9 +3,15 @@
Beyond square lattices: graphene
--------------------------------
.. seealso::
You can execute the code examples live in your browser by
activating thebelab:
.. thebe-button:: Activate Thebelab
.. seealso::
The complete source code of this example can be found in
:download:`graphene.py </code/download/graphene.py>`
:jupyter-download:script:`graphene`
In the following example, we are going to calculate the
conductance through a graphene quantum dot with a p-n junction
......@@ -18,9 +24,44 @@ We begin by defining the honeycomb lattice of graphene. This is
in principle already done in `kwant.lattice.honeycomb`, but we do it
explicitly here to show how to define a new lattice:
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_hnla
:end-before: #HIDDEN_END_hnla
.. jupyter-kernel::
:id: graphene
.. jupyter-execute::
:hide-code:
# Tutorial 2.5. Beyond square lattices: graphene
# ==============================================
#
# Physics background
# ------------------
# Transport through a graphene quantum dot with a pn-junction
#
# Kwant features highlighted
# --------------------------
# - Application of all the aspects of tutorials 1-3 to a more complicated
# lattice, namely graphene
from math import pi, sqrt, tanh
from matplotlib import pyplot
import kwant
# For computing eigenvalues
import scipy.sparse.linalg as sla
sin_30, cos_30 = (1 / 2, sqrt(3) / 2)
.. jupyter-execute:: boilerplate.py
:hide-code:
.. jupyter-execute::
graphene = kwant.lattice.general([(1, 0), (sin_30, cos_30)],
[(0, 0), (0, 1 / sqrt(3))],
norbs=1)
a, b = graphene.sublattices
The first argument to the `~kwant.lattice.general` function is the list of
primitive vectors of the lattice; the second one is the coordinates of basis
......@@ -31,9 +72,30 @@ itself forms a regular lattice of the same type as well, and those
In the next step we define the shape of the scattering region (circle again)
and add all lattice points using the ``shape``-functionality:
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_shzy
:end-before: #HIDDEN_END_shzy
.. jupyter-execute::
:hide-code:
r = 10
w = 2.0
pot = 0.1
.. jupyter-execute::
#### Define the scattering region. ####
# circular scattering region
def circle(pos):
x, y = pos
return x ** 2 + y ** 2 < r ** 2
syst = kwant.Builder()
# w: width and pot: potential maximum of the p-n junction
def potential(site):
(x, y) = site.pos
d = y * cos_30 + x * sin_30
return pot * tanh(d / w)
syst[graphene.shape(circle, (0, 0))] = potential
As you can see, this works exactly the same for any kind of lattice.
We add the onsite energies using a function describing the p-n junction;
......@@ -45,9 +107,11 @@ As a next step we add the hoppings, making use of
`~kwant.builder.HoppingKind`. For illustration purposes we define
the hoppings ourselves instead of using ``graphene.neighbors()``:
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_hsmc
:end-before: #HIDDEN_END_hsmc
.. jupyter-execute::
# 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))
The nearest-neighbor model for graphene contains only
hoppings between different basis atoms. For this type of
......@@ -63,27 +127,50 @@ respect to the two primitive vectors ``[(1, 0), (sin_30, cos_30)]``.
Adding the hoppings however still works the same way:
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_bfwb
:end-before: #HIDDEN_END_bfwb
.. jupyter-execute::
syst[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = -1
Modifying the scattering region is also possible as before. Let's
do something crazy, and remove an atom in sublattice A
(which removes also the hoppings from/to this site) as well
as add an additional link:
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_efut
:end-before: #HIDDEN_END_efut
.. jupyter-execute::
# Modify the scattering region
del syst[a(0, 0)]
syst[a(-2, 1), b(2, 2)] = -1
Note again that the conversion from a tuple `(i,j)` to site
is done by the sublattices `a` and `b`.
The leads are defined almost as before:
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_aakh
:end-before: #HIDDEN_END_aakh
.. jupyter-execute::
#### 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))] = -pot
lead0[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = -1
# 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))] = pot
lead1[[kwant.builder.HoppingKind(*hopping) for hopping in hoppings]] = -1
Note the method `~kwant.lattice.Polyatomic.vec` used in calculating the
parameter for `~kwant.lattice.TranslationalSymmetry`. The latter expects a
......@@ -98,28 +185,63 @@ in a square lattice -- they follow the non-orthogonal primitive vectors defined
in the beginning.
Later, we will compute some eigenvalues of the closed scattering region without
leads. This is why we postpone attaching the leads to the system. Instead,
we return the scattering region and the leads separately.
leads. This is why we postpone attaching the leads to the system.
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_kmmw
:end-before: #HIDDEN_END_kmmw
The computation of some eigenvalues of the closed system is done
in the following piece of code:
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_zydk
:end-before: #HIDDEN_END_zydk
.. jupyter-execute::
def compute_evs(syst):
# Compute some eigenvalues of the closed system
sparse_mat = syst.hamiltonian_submatrix(sparse=True)
evs = sla.eigs(sparse_mat, 2)[0]
print(evs.real)
The code for computing the band structure and the conductance is identical
to the previous examples, and needs not be further explained here.
Finally, in the `main` function we make and plot the system:
Finally, we plot the system:
.. jupyter-execute::
:hide-code:
def plot_conductance(syst, energies):
# Compute transmission as a function of energy
data = []
for energy in energies:
smatrix = kwant.smatrix(syst, energy)
data.append(smatrix.transmission(0, 1))
pyplot.figure()
pyplot.plot(energies, data)
pyplot.xlabel("energy [t]")
pyplot.ylabel("conductance [e^2/h]")
pyplot.show()
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_itkk
:end-before: #HIDDEN_END_itkk
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()
.. jupyter-execute::
# 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(syst, site_color=family_colors, site_lw=0.1, colorbar=False);
We customize the plotting: we set the `site_colors` argument of
`~kwant.plotter.plot` to a function which returns 0 for
......@@ -132,27 +254,43 @@ The function `~kwant.plotter.plot` shows these values using a color scale
(grayscale by default). The symbol `size` is specified in points, and is
independent on the overall figure size.
Plotting the closed system gives this result:
.. image:: /code/figure/graphene_syst1.*
Computing the eigenvalues of largest magnitude,
.. literalinclude:: /code/include/graphene.py
:start-after: #HIDDEN_BEGIN_jmbi
:end-before: #HIDDEN_END_jmbi
.. jupyter-execute::
should yield two eigenvalues equal to ``[ 3.07869311,
compute_evs(syst.finalized())
yields two eigenvalues equal to ``[ 3.07869311,
-3.06233144]``.
The remaining code of `main` attaches the leads to the system and plots it
The remaining code attaches the leads to the system and plots it
again:
.. image:: /code/figure/graphene_syst2.*
.. jupyter-execute::
# Attach the leads to the system.
for lead in [lead0, lead1]:
syst.attach_lead(lead)
# Then, plot the system with leads.
kwant.plot(syst, site_color=family_colors, site_lw=0.1,
lead_site_lw=0, colorbar=False);
It computes the band structure of one of lead 0:
We then finalize the system:
.. image:: /code/figure/graphene_bs.*
.. jupyter-execute::
syst = syst.finalized()
and compute the band structure of one of lead 0:
.. jupyter-execute::
# Compute the band structure of lead 0.
momenta = [-pi + 0.02 * pi * i for i in range(101)]
plot_bandstructure(syst.leads[0], momenta)
showing all the features of a zigzag lead, including the flat
edge state bands (note that the band structure is not symmetric around
......@@ -160,7 +298,11 @@ zero energy, due to a potential in the leads).
Finally the transmission through the system is computed,
.. image:: /code/figure/graphene_result.*
.. jupyter-execute::
# Plot conductance.
energies = [-2 * pot + 4. / 50. * pot * i for i in range(51)]
plot_conductance(syst, energies)
showing the typical resonance-like transmission probability through
an open quantum dot
......
doc/source/tutorial/index.rst
View file @ 39b65acd
......@@ -12,4 +12,5 @@ Tutorial: learning Kwant through examples
plotting
kpm
discretize
magnetic_field
faq