Computation Kernels
This essentially boils down to writing the functions for the RHS of the time-dependent Schrödinger equation. These will look mostly similar but will be slightly different to handle the following cases:
- We start in an eigenstate (we just need to evolve the difference)
- We start in an arbitrary state
- Finite system
- Infinite system (need a boundary-agnostic implementation)
These need to be written in Cython to be callable from C (they will be called by odeint).
We also need to write routines to handle the different kinds of boundary conditions that we want to support:
- fixed: a fixed number of extra layers in the lead, possibly with a complex absorbing potential applied
- adaptive: the number of extra layers in the lead grows as you request longer and longer times
- convolution: adds no extra lead layers, but solves the perfectly absorbing boundary condition, which has the form of a convolution integral involving the retarded self-energy
We also need to review different options for sparse matrix-vector product.
Designing against sparse BLAS seems a reasonable bet, there is already
librsb; we should benchmark this against
naive NIST sparse blas, against MKL (if we can
get it), and even against the scipy implementation. Need to test with realistic
matrices and also benchmark the time to go from a COO/CSR matrix given by
hamiltonian_submatrix
to the desired format for various different types of
sparse matrix (although this cost will be zero for the time-independent part of
the Hamiltonian, as we can do it once at the start).