Skip to content
GitLab
Menu
Projects
Groups
Snippets
Loading...
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in / Register
Toggle navigation
Menu
Open sidebar
SPEC-GMT
tkwant-energy-transport
Commits
bc63427c
Commit
bc63427c
authored
Feb 13, 2020
by
Kloss
Browse files
move examples in separate dir and symlink
parent
1d0a9f8d
Changes
30
Hide whitespace changes
Inline
Side-by-side
doc/
source/tutorial
/restarting.py
→
doc/
examples
/restarting.py
View file @
bc63427c
File moved
doc/
source/tutorial
/restarting.rst
→
doc/
examples
/restarting.rst
View file @
bc63427c
File moved
doc/
source/tutorial
/voltage_raise.py
→
doc/
examples
/voltage_raise.py
View file @
bc63427c
...
...
@@ -67,14 +67,13 @@ def tkwant_calculation(syst, operator, tmax, chemical_potential, V_static, V_dyn
params
=
{
'V_static'
:
V_static
,
'V_dynamic'
:
V_dynamic
}
# occupation -- for each lead
occupation_l
=
tkwant
.
manybody
.
make
_occupation
(
chemical_potential
+
V_static
)
occupation_r
=
tkwant
.
manybody
.
make
_occupation
(
chemical_potential
)
occupation
=
[
occupation_l
,
occupation_r
]
occupation_l
=
tkwant
.
manybody
.
lead
_occupation
(
chemical_potential
+
V_static
)
occupation_r
=
tkwant
.
manybody
.
lead
_occupation
(
chemical_potential
)
occupation
s
=
[
occupation_l
,
occupation_r
]
# Create the solver
solver
=
tkwant
.
manybody
.
State
(
syst
,
tmax
,
occupation
,
params
=
params
)
# TODO: add when doc has changed
# solver.refine_intervals(rtol=1E-2, atol=1E-2) # require only low accuracy
# TODO: remove refine=False when doc needs less time to build
solver
=
tkwant
.
manybody
.
State
(
syst
,
tmax
,
occupations
,
params
=
params
,
refine
=
False
)
# loop over time, calculating the current as we go
expectation_values
=
[]
...
...
doc/
source/tutorial
/voltage_raise.rst
→
doc/
examples
/voltage_raise.rst
View file @
bc63427c
File moved
doc/source/conf.py
View file @
bc63427c
...
...
@@ -351,5 +351,7 @@ intersphinx_mapping = {
'python'
:
(
'https://docs.python.org/3.4/'
,
None
),
'numpy'
:
(
'https://docs.scipy.org/doc/numpy/'
,
None
),
'scipy'
:
(
'https://docs.scipy.org/doc/scipy/reference/'
,
None
),
'kwant'
:
(
'http://kwant-project.org/doc/1/'
,
None
),
'mpi4py'
:
(
'https://mpi4py.readthedocs.io/en/stable/'
,
None
),
'kwant'
:
(
'https://kwant-project.org/doc/1/'
,
None
),
'kwant_spectrum'
:
(
'https://kwant-project.org/extensions/spectrum/'
,
None
)
}
doc/source/index.rst
View file @
bc63427c
...
...
@@ -8,6 +8,7 @@ tkwant |release| documentation
pre/index
tutorial/index
examples/index
reference/index
* :ref:`genindex`
doc/source/reference/tkwant.manybody.rst
View file @
bc63427c
...
...
@@ -19,11 +19,12 @@ Helper functions
.. autosummary::
:toctree: generated
make
_occupation
lead
_occupation
calc_intervals
split_intervals
combine_intervals
calc_tasks
calc_initial_
manybody_
state
calc_initial_state
calc_energy_cutoffs
make_boundstates_time_dependent
add_boundstates
...
...
doc/source/tutorial/index.rst
View file @
bc63427c
...
...
@@ -9,4 +9,3 @@ Tutorial
boundary_condition
logging
mpi
examples
doc/source/tutorial/manybody_advanced.rst
View file @
bc63427c
...
...
@@ -19,7 +19,7 @@ Advanced manybody settings
from tkwant import manybody, leads
def make_system(a=1,
gamma=1.0
,
W
=2
, L=3
):
def make_system(a=1,
W=1
,
L
=2):
lat = kwant.lattice.square(a=a, norbs=1)
syst = kwant.Builder()
...
...
@@ -32,11 +32,11 @@ Advanced manybody settings
# time dependent coupling with gaussian pulse
def coupling_nn(site1, site2, time):
return -
gamma *
cmath.exp(- 1j * gaussian(time))
return - cmath.exp(- 1j * gaussian(time))
#### Define the scattering region. ####
syst[(lat(x, y) for x in range(L) for y in range(W))] =
2 * gamma
syst[lat.neighbors()] = -
gamma
syst[(lat(x, y) for x in range(L) for y in range(W))] =
1
syst[lat.neighbors()] = -
1
# time dependent coupling between two sites 0 and 1
for y in range(W):
syst[lat(0, y), lat(1, y)] = coupling_nn
...
...
@@ -44,24 +44,24 @@ Advanced manybody settings
#### Define and attach the leads. ####
# Construct the left lead.
lead = kwant.Builder(kwant.TranslationalSymmetry((-a, 0)))
lead[(lat(0, j) for j in range(W))] =
2 * gamma
lead[lat.neighbors()] = -
gamma
lead[(lat(0, j) for j in range(W))] =
1
lead[lat.neighbors()] = -
1
# Attach the left lead and its reversed copy.
syst.attach_lead(lead)
syst.attach_lead(lead.reversed())
#
syst.attach_lead(lead.reversed())
return syst
syst = make_system().finalized()
density_operator = kwant.operator.Density(syst)
state = manybody.State(syst, tmax=
50
0)
state = manybody.State(syst, tmax=
1
0)
def calc_my_boundstates(): # a routine to mimic boundstates
tasks = {100: tkwant.onebody.Task(weight=0.1, energy=0.5, lead=None, mode=None)}
psi = {key: tkwant.onebody.ScatteringStates(syst, energy=task.energy,
lead=0, tmax=
50
0)[0]
lead=0, tmax=
1
0)[0]
for key, task in tasks.items()}
return psi, tasks
...
...
@@ -88,7 +88,7 @@ intervals into subintervals: ``manybody.split_intervals()``.
that form the manybody state. Default method: ``manybody.calc_tasks()``.
**initial state** Calculates the initial manybody state. Default method:
``manybody.calc_initial_
manybody_
state()``.
``manybody.calc_initial_state()``.
**manybody** MPI scheduler to evolve the manybody state and evaluate
observables. Default method: ``manybody.WaveFunction()``.
...
...
@@ -129,16 +129,16 @@ operator, each one-body state contributes according to its weighting
factor, that is tabulated in the ``tasks`` dictionary. The basic setup
to set up the state is as follows:
.. jupyter-execute::
#
.. jupyter-execute::
spectra = kwant_spectrum.spectra(syst.leads)
occupation = manybody.
make
_occupation()
occupation = manybody.
lead
_occupation()
intervals = manybody.calc_intervals(spectra, occupation)
emin, emax = manybody.calc_energy_cutoffs(occupation)
boundaries = leads.automatic_boundary(spectra, tmax=
50
0, emin=emin, emax=emax)
boundaries = leads.automatic_boundary(spectra, tmax=
1
0, emin=emin, emax=emax)
tasks = manybody.calc_tasks(intervals, spectra, occupation)
psi_init = manybody.calc_initial_
manybody_
state(syst, tasks, boundaries)
psi_init = manybody.calc_initial_state(syst, tasks, boundaries)
state = manybody.WaveFunction(psi_init, tasks)
This way offers maximal flexibility, but is in most cases not accurate
...
...
@@ -156,7 +156,7 @@ maximal simulation time:
.. jupyter-execute::
state = manybody.State(syst, tmax=
50
0)
state = manybody.State(syst, tmax=
1
0)
In addition to the two previous state shown before, the adaptive
state can refine the momentum intervals in a semi-automatic way. The
...
...
@@ -188,10 +188,10 @@ Statistical average
The occupation depends on the temperature :math:`T`, the chemical
potential :math:`\mu`, as well as the distribution function :math:`f`,
over which the thermal average takes place. The function
``manybody.
make
_occupation()`` offers an easy way to set the lead occupation.
``manybody.
lead
_occupation()`` offers an easy way to set the lead occupation.
By default, the average is taken at zero-temperature using
the noninteracting Fermi-Dirac distribution function. The manybody
state uses the output of ``manybody.
make
_occupation()`` to extract the
state uses the output of ``manybody.
lead
_occupation()`` to extract the
intervals to perform the momentum or respectively energy integration.
In the case that several leads are attached to the system, the occupation
can be a sequence of ``occupation`` objects, one for each lead.
...
...
@@ -203,20 +203,20 @@ system:
.. jupyter-execute::
occupation = manybody.
make
_occupation(chemical_potential=1)
occupation = manybody.
lead
_occupation(chemical_potential=1)
``manybody.
make
_occupation()`` returns a ``manybody.Occupation`` instance with the physical
``manybody.
lead
_occupation()`` returns a ``manybody.Occupation`` instance with the physical
information on how to perform the statistical average. The object looks
like this:
.. jupyter-execute::
occupation = manybody.
make
_occupation()
occupation = manybody.
lead
_occupation()
print(occupation)
By default, the thermal average is taken at zero-temperature using the
noninteracting Fermi-Dirac distribution function. In the following we
show how to change the default behavior of ``manybody.
make
_occupation()``.
show how to change the default behavior of ``manybody.
lead
_occupation()``.
Finite temperature
~~~~~~~~~~~~~~~~~~
...
...
@@ -227,7 +227,7 @@ example, we set the temperature of the lead to the finite value
.. jupyter-execute::
occupation = manybody.
make
_occupation(temperature=0.5)
occupation = manybody.
lead
_occupation(temperature=0.5)
Distribution function
~~~~~~~~~~~~~~~~~~~~~
...
...
@@ -241,7 +241,7 @@ at finite temperature:
def bose_function(mu, T, energy):
return 1 / (exp((energy - mu) / T) - 1)
occupation = manybody.
make
_occupation(temperature=0.5, distribution=bose_function)
occupation = manybody.
lead
_occupation(temperature=0.5, distribution=bose_function)
Note that the distribution function :math:`f` must have the calling
signature ``(chemical_potential, temperature, energy)``.
...
...
@@ -257,7 +257,7 @@ only the modes with energies :math:`0.4 \leq E_n \leq 0.6`:
.. jupyter-execute::
occupation = manybody.
make
_occupation(energy_range=[(0.4, 0.6)])
occupation = manybody.
lead
_occupation(energy_range=[(0.4, 0.6)])
Include / exclude bands
~~~~~~~~~~~~~~~~~~~~~~~
...
...
@@ -271,13 +271,14 @@ over the two lowest energy :math:`E_n` bands with band index
.. jupyter-execute::
occupation = manybody.
make
_occupation(bands=[0, 1])
occupation = manybody.
lead
_occupation(bands=[0, 1])
Note that the provided band indicees must be in the physically valid
range, that is, they must not exceed the maximum number of bands. A
realistic example with band selection is shown in example
`[ex2] <#examples>`__.
Include / exclude leads
~~~~~~~~~~~~~~~~~~~~~~~
...
...
@@ -285,17 +286,18 @@ For system with several leads, the occupation passed to
``manybody.State()`` can be a sequence with one element per leads. If a
lead should not contribute to the statistial average, the corresponding
element of the ``occupation`` sequence must be set to ``None``.
In the following example, our system has two leads,
but only contribution of the lead with index 0 is taken into account:
In the following example, our system
is expected to
has two leads,
but only
the
contribution of the lead with index 0 is taken into account:
.. jupyter-execute::
occup = manybody.
make
_occupation()
occupation = [occup, None]
occup = manybody.
lead
_occupation()
occupation
_0
= [occup, None]
A realistic example with lead selection is shown in example
`[ex2] <#examples>`__.
Numerical integration
---------------------
...
...
@@ -320,8 +322,10 @@ with ``kwant_spectrum.spectra()``.
.. jupyter-execute::
spectra = kwant_spectrum.spectra(syst.leads)
occupation = manybody.lead_occupation()
intervals = manybody.calc_intervals(spectra, occupation)
While the above intervals would lead to the default result when passed to the manybody state,
one can now easily modify individual elements of the ``intervals`` list.
An alternative approach to precalculate an interval sequence with changed
...
...
@@ -338,7 +342,7 @@ Finally, the intervals are passed to the state with the keyword ``intervals``:
.. jupyter-execute::
state = manybody.State(syst, tmax=
50
0, intervals=intervals)
state = manybody.State(syst, tmax=
1
0, intervals=intervals)
For the second way, the modified ``manybody.Interval`` data class can passed
also directly to the manybody state.
...
...
@@ -347,11 +351,12 @@ The analog second approach to the above example with changed quadrature order is
.. jupyter-execute::
Interval = functools.partial(manybody.Interval, order=10)
state = manybody.State(syst, tmax=
50
0, intervals=Interval)
state = manybody.State(syst, tmax=
1
0, intervals=Interval)
In the following, both ways are presented for different default parameters.
Momentum intervals
~~~~~~~~~~~~~~~~~~
...
...
@@ -374,6 +379,7 @@ concrete quadrature rule that will be applied. ``integration_variable`` is an
additional information to switch between energy and momentum
integration.
Integration order
~~~~~~~~~~~~~~~~~
...
...
@@ -394,7 +400,7 @@ Alternatively, the interval order used by the manybody state is changed via
.. jupyter-execute::
state = manybody.State(syst, tmax=
50
0, intervals=Interval)
state = manybody.State(syst, tmax=
1
0, intervals=Interval)
The order is usually taken between 10 and 20. If the integration is not
accurate enough one should rather divide each interval into subintervals
...
...
@@ -413,6 +419,7 @@ subintervals to increase the numerical accuracy. Here we divide into
# (kmin, kmax) -> [(kmin, k_0), (k_0, k_1).. (k_{n-1}, kmax)]
intervals = manybody.split_intervals(intervals, number_subintervals=10)
**Energy split**
Splitting intervals in energy space is possible by using `energy_range`.
The snippet below is to reproduce older calculations performed with
...
...
@@ -432,7 +439,7 @@ energy integrations on an energy interval splitted equidistantly into subinterva
# energy interval split (0, 1) -> [(0, 0.1), (0.1, 0.2).. (0.9, 1)]
energy_range = [i for i in pairwise(np.linspace(0, chemical_potential, 11))]
occupation = manybody.
make
_occupation(chemical_potential,
occupation = manybody.
lead
_occupation(chemical_potential,
energy_range=energy_range)
Interval = functools.partial(manybody.Interval, integration_variable='energy')
...
...
@@ -443,6 +450,7 @@ with ``manybody.split_intervals`` is needed. If
``integration_variable="momentum"``, the split is still performed in the energy
range, but the quadrature will be integrate over momentum.
Quadrature error
~~~~~~~~~~~~~~~~
...
...
@@ -466,7 +474,9 @@ calculated with the higher (:math:`2 n +1`) quadrature rule.
intervals = manybody.calc_intervals(spectra, occupation)
tasks = manybody.calc_tasks(intervals, spectra, occupation)
psi_init = manybody.calc_initial_manybody_state(syst, tasks, boundaries)
emin, emax = manybody.calc_energy_cutoffs(occupation)
boundaries = leads.automatic_boundary(spectra, tmax=10, emin=emin, emax=emax)
psi_init = manybody.calc_initial_state(syst, tasks, boundaries)
state = manybody.WaveFunction(psi_init, tasks)
density = state.evaluate(density_operator)
print('integration error delta= {}'.format(maximal_absolute_error(density)))
...
...
@@ -475,6 +485,7 @@ The quadrature rule applied to each interval can be accessed using
the keyword ``quadrature``. Note that the adaptive state ``manybody.State``
needs a Gauss-Kronrod like rule for the error estimate.
Integral type
~~~~~~~~~~~~~
...
...
@@ -525,14 +536,14 @@ Alternatively, the integral type used by the manybody state is changed via
.. jupyter-execute::
state = manybody.State(syst, tmax=
50
0, intervals=Interval)
state = manybody.State(syst, tmax=
1
0, intervals=Interval)
Modes and weights
~~~~~~~~~~~~~~~~~
The low-level state ``manybody.WaveFunction()`` and the initial state
(via ``calc_initial_
manybody_
state()``) work directly with the ``kwant``
(via ``calc_initial_state()``) work directly with the ``kwant``
modes defined at a certain energy with a corresponding weight. The
routine ``calc_tasks()`` calculates the modes and weights for all
intervals:
...
...
@@ -553,7 +564,7 @@ and looks like this:
print(task)
For each set of ``(lead, mode, energy)`` attributes per element, the function
``manybody.calc_initial_
manybody_
state()`` calculates the corresponding
``manybody.calc_initial_state()`` calculates the corresponding
one-body state :math:`| \psi_{\alpha E_l}(t=0) \rangle` that form the
initial manybody state. The ``weights`` attribute correspond to the weighting
factors :math:`w` in the manybody average. Evaluating an observable over
...
...
@@ -565,7 +576,7 @@ task list as arguments:
.. jupyter-execute::
psi_init = manybody.calc_initial_
manybody_
state(syst, tasks, boundaries)
psi_init = manybody.calc_initial_state(syst, tasks, boundaries)
state = manybody.WaveFunction(psi_init, tasks)
Adaptive refinement
...
...
@@ -578,7 +589,7 @@ can change to a different observable as follows:
.. jupyter-execute::
current_operator = kwant.operator.Current(syst)
state = manybody.State(syst, tmax=
50
0, error_estimate_operator=current_operator)
state = manybody.State(syst, tmax=
1
0, error_estimate_operator=current_operator)
.. jupyter-execute::
...
...
@@ -602,7 +613,7 @@ example for ``manybody.State()``:
solver_type = functools.partial(tkwant.onebody.solvers.default, rtol=1E-5)
onebody_wavefunction = functools.partial(tkwant.onebody.WaveFunction.from_kwant, solver_type=solver_type)
state = manybody.State(syst, tmax=
50
0, onebody_wavefunction_type=onebody_wavefunction)
state = manybody.State(syst, tmax=
1
0, onebody_wavefunction_type=onebody_wavefunction)
A similar strategy is possible to change the onebody kernels
``onebody.kernels`` that evaluate the right-hand-side of the one-body
...
...
@@ -705,7 +716,7 @@ the keyword ``spectra``:
.. jupyter-execute::
spectra = kwant_spectrum.spectra(syst.leads)
state = manybody.State(syst, tmax=
50
0, spectra=spectra)
state = manybody.State(syst, tmax=
1
0, spectra=spectra)
Boundary conditions
~~~~~~~~~~~~~~~~~~~
...
...
@@ -716,7 +727,7 @@ precalculate and transferred to the state with keyword ``boundaries``:
.. jupyter-execute::
boundaries = [tkwant.leads.SimpleBoundary(tmax=
50
0)] * len(syst.leads)
boundaries = [tkwant.leads.SimpleBoundary(tmax=
1
0)] * len(syst.leads)
state = manybody.State(syst, boundaries=boundaries)
Initial state
...
...
@@ -754,7 +765,7 @@ here called ``calc_my_boundstates()``:
``boundstate_psi`` is a dictionary of all boundstate wavefunctions.
It has basically the same form as an initial manybody state calculated from
``manybody.calc_initial_
manybody_
state()``.
``manybody.calc_initial_state()``.
``boundstate_tasks`` contain the weighting factor and the energy of each
boundstate and is similar to the output of ``manybody.calc_tasks()``.
...
...
@@ -792,3 +803,5 @@ devices <https://link.springer.com/article/10.1007%2Fs10825-016-0855-9>`__, J. C
Electron. **15**, 1148 (2016).
`[arXiv] <https://arxiv.org/abs/1604.01198>`__
doc/source/tutorial/quantum_dynamic_simulations.rst
View file @
bc63427c
...
...
@@ -199,11 +199,10 @@ The tkwant code to perform above task is
lead_left[lat(0, 0)] = 1
lead_left[lat.neighbors()] = -1
syst.attach_lead(lead_left)
syst.attach_lead(lead_left.reversed())
return syst
syst = create_system(
5
).finalized()
syst = create_system(
2
).finalized()
.. jupyter-execute::
...
...
@@ -308,7 +307,7 @@ The manybody state is initialized as:
.. jupyter-execute::
state = manybody.State(syst, tmax=
50
0)
state = manybody.State(syst, tmax=
1
0)
The ``state`` instance has an ``evolve()`` and an ``evaluate()``
method, to evolve it forward in time and to evaluate an operator. We will
...
...
@@ -323,7 +322,7 @@ method ``evolve()`` with the time argument:
.. jupyter-execute::
state.evolve(time=
10
)
state.evolve(time=
5
)
Note that the evolution time :math:`t` must be less than the maximal
time ``tmax`` passed to ``manybody.State()`` during instantiation.
...
...
@@ -394,7 +393,7 @@ The two pieces of pseudo code compare the standard setup for both states:
.. jupyter-execute::
psi = onebody.ScatteringStates(syst, energy=1, lead=0, tmax=
50
0)[0]
psi = onebody.ScatteringStates(syst, energy=1, lead=0, tmax=
1
0)[0]
psi.evolve(time=5)
density = psi.evaluate(density_operator)
...
...
@@ -402,7 +401,7 @@ The two pieces of pseudo code compare the standard setup for both states:
.. jupyter-execute::
state = manybody.State(syst, tmax=
50
0)
state = manybody.State(syst, tmax=
1
0)
state.evolve(time=5)
density = state.evaluate(density_operator)
...
...
Prev
1
2
Next
Write
Preview
Markdown
is supported
0%
Try again
or
attach a new file
.
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment