Commit ee5d281f authored by Kloss's avatar Kloss
Browse files

add boundstates routine as methods to the manybody solver

parent c2ae2014
......@@ -38,5 +38,6 @@ Data types
Occupation
Interval
Task
BoundState
......@@ -618,6 +618,51 @@ The time attribute ``new_solve.time`` is set to the time of the initial
state. The user must ensure that all one-body states forming the initial
state have the same time.
Boundstates
~~~~~~~~~~~
The following pseudocode illustrates how to include boundstates in the manybody
dynamics. First one has to precalculate the boundstates:
.. code:: python
boundstates = calc_my_boundstates(...)
In above pseudo-code, ``boundstates`` is a sequence of named tuples of type ``manybody.BoundState``:
.. jupyter-execute::
:hide-code:
print(manybody.BoundState.__doc__);
Apart from ``energy``, ``mode`` distinguishes degenerate states, and
``wavefunction`` is a single boundstate that can be evolved and evaluated like
a ``onebody.Wavefunction`` instance.
One also has to provide a callable distribution function that returns a
scalar numerical weight factor for the boundstates:
.. code:: python
def my_boundstate_distribution(energy):
...
return weight
If ``solver`` is a manybody solver instance created with ``manybody.Solver``
or ``manybody.fixed_quadrature_solver``, boundstates are easily added by
the method ``add_boundstates()``:
.. code:: python
solver.add_boundstates(boundstates, distribution=my_boundstate_distribution)
For a manybody solver instance created with ``manybody.ManybodyTaskScheduler``,
the boundstates are added via the function ``manybody.add_boundstates()``:
.. code:: python
manybody.add_boundstates(solver, boundstates, my_boundstate_distribution)
Examples
--------
......
......@@ -19,7 +19,7 @@ import numpy as np
from . import onebody, leads, integration
__all__ = ['Occupation', 'Interval', 'Task',
__all__ = ['Occupation', 'Interval', 'Task', 'BondState',
'occupation', 'calc_intervals', 'split_intervals',
'get_distribution', 'calc_tasks', 'calc_initial_manybody_state',
'ManybodyTaskScheduler', 'fixed_quadrature_solver', 'Solver',
......@@ -37,6 +37,9 @@ Interval = namedtuple('Interval', ['lead', 'band', 'kmin', 'kmax',
Task = namedtuple('Task', ['lead', 'mode', 'energy', 'weight'])
"""Data format for a single one-body state and its weight in the manybody sum."""
BoundState = namedtuple('BoundState', ['energy', 'mode', 'wavefunction'])
"""Data format for a single boundstate."""
# helper routines
......@@ -1072,10 +1075,33 @@ def fixed_quadrature_solver(syst, intervals, distribution, spectrum,
solver : `tkwant.manybody.ManybodyTaskScheduler`
Manybody solver instance
"""
def _add_boundstates(boundstates, distribution):
"""Add a sequence of boundstates to the manybody solver
Parameters
----------
boundstates : `tkwant.manybody.BoundState` or sequence thereof
Boundstates to be added to the solver
distribution : callable
Distribution function for the boundstates.
Calling signature: `(energy)`.
Returns
-------
keys : list
List with all dict keys of the added boundstates.
The list has length `len(boundstates)` and the same ordering as
the `boundstates` sequence.
"""
keys = add_boundstates(solver, boundstates, distribution)
return keys
tasks = calc_tasks(intervals, distribution, spectrum)
psi_init = calc_initial_manybody_state(syst, tasks, boundaries, args,
onebody_solver, mpi_distribute, comm)
return manybody_solver(psi_init, tasks, comm)
solver = manybody_solver(psi_init, tasks, comm)
solver.add_boundstates = _add_boundstates
return solver
def get_adaptive_solver(base=ManybodyTaskScheduler):
......@@ -1422,29 +1448,70 @@ def get_adaptive_solver(base=ManybodyTaskScheduler):
return super().evaluate(observable)[return_element]
return super().evaluate(observable)
def add_boundstates(self, boundstates, distribution):
"""Add a sequence of boundstates to the manybody solver
Parameters
----------
boundstates : `tkwant.manybody.BoundState` or sequence thereof
Boundstates to be added to the solver
distribution : callable
Distribution function for the boundstates.
Calling signature: `(energy)`.
Returns
-------
keys : list
List with all dict keys of the added boundstates.
The list has length `len(boundstates)` and the same ordering as
the `boundstates` sequence.
"""
keys = add_boundstates(self, boundstates, distribution)
return keys
return Solver
Solver = get_adaptive_solver()
def add_boundstates(solver, boundstates):
def add_boundstates(solver, boundstates, distribution):
"""Add a sequence of boundstates to the manybody solver
Parameters
----------
solver : `tkwant.manybody.ManybodyTaskScheduler`
Manybody solver
boundstates : sequence of boundstates
Manybody solver to which the boundstates should be added.
boundstates : `tkwant.manybody.BoundState` or sequence thereof
Boundstates to be added to the solver
distribution : callable
Distribution function for the boundstates. Calling signature: `(energy)`.
Must return a real scalar weight factor.
Returns
-------
solver : `tkwant.manybody.ManybodyTaskScheduler`
Manybody solver with boundstates
keys : list
List with all dict keys of the added boundstates.
The list has length `len(boundstates)` and the same ordering as
the `boundstates` sequence.
Notes
-----
`solver` is modified upon return.
If any exceptions are raised in this function there are no
guarantees that the `solver` was not modified.
"""
Task = namedtuple('task', 'weight')
# all weights are 1 for boundstates
task = Task(np.ones(shape=solver._get_weight_shape()))
if isinstance(boundstates, BoundState):
boundstates = [boundstates]
if not callable(distribution):
raise TypeError('boundstate distribution function is not callable')
weight_shape = np.ones(shape=solver._get_weight_shape())
keys = []
for state in boundstates:
solver.add_onebody_state(state, task)
return solver
# weights are 1 for boundstates but include only thermal distribution
weight = distribution(state.energy) * weight_shape
task = Task(lead=None, mode=state.mode, energy=state.energy, weight=weight)
key = solver.add_onebody_state(state=state.wavefunction, task=task)
keys.append(key)
return keys
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment