1.3.rst 9.37 KB
Newer Older
1
2
3
What's new in Kwant 1.3
=======================

Christoph Groth's avatar
Christoph Groth committed
4
5
This article explains the user-visible changes in Kwant 1.3.0,
released on 19 May 2017.
Christoph Groth's avatar
Christoph Groth committed
6
7
8
See also the `full list of changes up to the most recent bugfix
release of the 1.3 series
<https://gitlab.kwant-project.org/kwant/kwant/compare/v1.3.0...latest-1.3>`_.
9

Christoph Groth's avatar
Christoph Groth committed
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
Using high-symmetry builders as models
--------------------------------------
Builders now have a `~kwant.builder.Builder.fill` method that fills a builder
instance with copies of a template builder. This can be used to "cut out"
shapes from high-symmetry models, or to increase the symmetry period of a lead.

Thus Kwant gains the new concept of a "model".  Models may be created manually,
or with the new function `kwant.continuum.discretize` (see next paragraph).
There is also support for finalizing models and e.g. calculating their band
structure (see `Finalizing builders with multiple translational symmetries`_).

Tools for continuum Hamiltonians
--------------------------------
The new sub-package `~kwant.continuum` is a collection of tools for working
with continuum models and for discretizing them into tight-binding models. It
aims at providing a handy interface to convert symbolic Hamiltonians into a
builder with N-D translational symmetry that can be use to calculate
tight-binding band structures or construct systems with different/lower
28
29
symmetry. For example in just a few lines we can construct a two-band model that exhibits
a quantum anomalous spin Hall phase:
Christoph Groth's avatar
Christoph Groth committed
30

31
32
33
34
35
36
37
38
39
.. literalinclude:: ../../tutorial/plot_qahe.py
    :start-after: HIDDEN_BEGIN_model
    :end-before: HIDDEN_END_model

From: :download:`QAHE example script <../../tutorial/plot_qahe.py>`

See the tutorial: :doc:`../../tutorial/discretize`

See the reference documentation: :doc:`../../reference/kwant.continuum`
Christoph Groth's avatar
Christoph Groth committed
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61

Calculating charges and currents using the operator module
----------------------------------------------------------
Often one may wish to calculate quantities that are defined over sites of
the system (such as charge density, spin density along some axis etc),
or over hoppings of the system (such as current or spin current). With
the introduction of the ``operator`` module it has now become much easier
to calculate such quantities. To obtain the regular density and current
everywhere in a system due to a wavefunction ``psi``, one only needs to do
the following::

    syst = make_system().finalized()
    psi = kwant.wave_function(syst)(0)[0]

    # create the operators
    Q = kwant.physics.LocalOperator(syst)
    J = kwant.physics.Current(syst)

    # evaluate the expectation value with the wavefunction
    q = Q(psi)
    j = J(psi)

Joseph Weston's avatar
Joseph Weston committed
62
See the tutorial: :doc:`../../tutorial/operators`
Christoph Groth's avatar
Christoph Groth committed
63
64
65
66

Plotting of currents
--------------------
Quantities defined on system hoppings (e.g. currents calculated using
Joseph Weston's avatar
Joseph Weston committed
67
`~kwant.operator.Current`) can be directly plotted as a streamplot over the
Christoph Groth's avatar
Christoph Groth committed
68
69
system using `kwant.plotter.current`. This is similar to how
`kwant.plotter.map` can be used to plot quantities defined on sites.
70
71
72
73
74
75
76
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:: ../../tutorial/plot_qahe.py
    :start-after: HIDDEN_BEGIN_current
    :end-before: HIDDEN_END_current
Christoph Groth's avatar
Christoph Groth committed
77

78
.. image:: ../../images/plot_qahe_current.*
Christoph Groth's avatar
Christoph Groth committed
79

80
From: :download:`QAHE example script <../../tutorial/plot_qahe.py>`
Christoph Groth's avatar
Christoph Groth committed
81
82
83
84
85
86
87
88
89
90
91
92
93
94

Scattering states with discrete symmetries and conservation laws
----------------------------------------------------------------
Given a lead Hamiltonian that has a conservation law, it is now possible to
construct lead modes that have definite values of the conservation law. This
is done by declaring projectors that block diagonalize the Hamiltonian before
the modes are computed. For a Hamiltonian that has one or more of the three
fundamental discrete symmetries (time-reversal symmetry, particle-hole symmetry
and chiral symmetry), it is now possible to declare the symmetries in Kwant.
The symmetries are then used to construct scattering states that are properly
related by symmetry. The discrete symmetries may be combined with conservation
laws, such that if different blocks of the Hamiltonian are related by a discrete
symmetry, the lead modes are computed to reflect this.

Joseph Weston's avatar
Joseph Weston committed
95
See the updated tutorial: :doc:`../../tutorial/superconductors`
Christoph Groth's avatar
Christoph Groth committed
96
97
98

Named parameters for value functions
------------------------------------
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
In Kwant < 1.3 whenever Hamiltonian values were provided as functions,
they all had to take the same extra parameters (after the site(s))
regardless of whether or not they actually used them at all. For example,
if we had some onsite potential and a magnetic field that we
model using the Peierls substitution, we would have to define our value
functions like so::

    # formally depends on 'B', but 'B' is never used
    def onsite(site, V, B):
        return V

    # formally depends on 'V', but 'V' is never used
    def hopping(site_a, site_b, V, B):
        return (site_b.pos[1] - site_a.pos[1]) * B

This was because previously extra arguments were provided to the system
by passing them as a sequence via the ``args`` parameter to various Kwant
functions (e.g. ``kwant.smatrix`` or ``hamiltonian_submatrix``).

In Kwant 1.3 it is now possible for value functions to depend on different
parameters, e.g.::

    def onsite(site, V):
        return V

    def hopping(site_a, site_b, B):
        return (site_b.pos[1] - site_a.pos[1]) * B

If you make use of this feature then you must in addition pass your arguments
via the ``params`` parameter. The value provided to ``params`` must
be a ``dict`` that maps parameter names to values, e.g.::

    kwant.smatrix(syst, params=dict(B=0.1, V=2))

as opposed to the old way::

    kwant.smatrix(syst, args=(2, 0.1))

Passing a dictionary of parameters via ``params`` is now the recommended way
to provide parameters to the system.

Christoph Groth's avatar
Christoph Groth committed
140
141
142
143
144
Reference implementation of the kernel polynomial method
--------------------------------------------------------
The kernel polynomial method is now implemented within Kwant to obtain the
density of states or, more generally, the spectral density of a given operator
acting on a system or Hamiltonian.
Joseph Weston's avatar
Joseph Weston committed
145

Joseph Weston's avatar
Joseph Weston committed
146
147
148
See the tutorial: :doc:`../../tutorial/kpm`

See the reference documentation: :doc:`../../reference/kwant.kpm`
Joseph Weston's avatar
Joseph Weston committed
149

Christoph Groth's avatar
Christoph Groth committed
150
151
152
153
154
Finalizing builders with multiple translational symmetries
----------------------------------------------------------
While it remains impossible to finalize a builder with more than a single
direction of translational symmetry, the ``wraparound`` module has been added
as a temporary work-around until the above limitation gets lifted.
Joseph Weston's avatar
Joseph Weston committed
155

Christoph Groth's avatar
Christoph Groth committed
156
157
158
159
160
161
162
163
164
The function `~kwant.wraparound.wraparound` transforms all (or all but one)
translational symmetries of a given builder into named momentum parameters
`k_x`, `k_y`, etc.  This makes it easy to compute transport through systems
with periodic boundary conditions or across infinite planes.

Plotting the 2-d band structure of graphene is now as straightforward as::

    from matplotlib import pyplot
    import kwant
Joseph Weston's avatar
Joseph Weston committed
165

Christoph Groth's avatar
Christoph Groth committed
166
167
    lat = kwant.lattice.honeycomb()
    sym = kwant.TranslationalSymmetry(lat.vec((1, 0)), lat.vec((0, 1)))
Joseph Weston's avatar
Joseph Weston committed
168

Christoph Groth's avatar
Christoph Groth committed
169
170
171
172
173
174
175
176
    bulk = kwant.Builder(sym)
    bulk[ [lat.a(0, 0), lat.b(0, 0)] ] = 0
    bulk[lat.neighbors()] = 1
    wrapped = kwant.wraparound.wraparound(bulk).finalized()
    kwant.wraparound.plot_2d_bands(wrapped)

Consistent ordering of sites in finalized builders
--------------------------------------------------
177
178
In Python 3 the internal ordering of dictionaries is not deterministic. This
meant that running a Kwant script twice would produce systems with different
Bas Nijholt's avatar
Bas Nijholt committed
179
ordering of sites, which leads to non-reproducible calculations. Now, sites
180
181
182
183
184
185
186
187
188
189
in finalized builders are always ordered first by their site family, then by
their tag.

Coincidentally, this means that you can plot a wavefunction in a simple 1D
system by just saying::

    lattice_1D = chain()
    syst = make_system(lattice_1D)
    h = syst.hamiltonian_submatrix()
    pyplot.plot(np.eigs(h)[1][0])
190

Christoph Groth's avatar
Christoph Groth committed
191
192
193
194
195
196
197
198
attach_lead() can now handle leads with greater than nearest-neighbor hoppings
------------------------------------------------------------------------------
When attaching a lead with greater than nearest-neighbor hoppings, the symmetry
period of the finalized lead is suitably extended and the unit cell size is
increased.

Pickling support
----------------
Joseph Weston's avatar
Joseph Weston committed
199
200
It is now possible to pickle and unpickle `~kwant.builder.Builder` and
`~kwant.system.System` instances.
Christoph Groth's avatar
Christoph Groth committed
201

202
Improved build configuration
Joseph Weston's avatar
Joseph Weston committed
203
----------------------------
204
205
206
207
208
209
210
The name of the build configuration file, ``build.conf`` by default, is now
configurable with the ``--configfile=PATH`` option to ``setup.py``.  (This
makes build configuration usable with the ``pip`` tool.)  The build
configuration as specified in this file is now more general, allowing to
modify any build parameter for any of the compiled extensions contained in
Kwant.  See the :ref:`Installation instructions <build-configuration>` for
details.
211

212
213
Builder.neighbors() respects symmetries
---------------------------------------
Christoph Groth's avatar
Christoph Groth committed
214
Given a site, the method `~kwant.builder.Builder.neighbors` of
Christoph Groth's avatar
typo    
Christoph Groth committed
215
`~kwant.builder.Builder` returns an iterator over sites that are connected by a
Christoph Groth's avatar
Christoph Groth committed
216
217
218
hopping to the provided site.  This is in contrast to previous versions of
Kwant, where the neighbors were yielded not of the provided site, but of it's
image in the fundamental domain.
219
220
221

This change is documented here for completeness.  We expect that the vast
majority of users of Kwant will not be affected by it.