sisl_toolbox.btd.DeviceGreen
- class sisl_toolbox.btd.DeviceGreen(H, elecs, pivot, eta=0.0)
Bases:
object
Block-tri-diagonal Green function calculator
This class enables the extraction and calculation of some important quantities not currently accessible in TBtrans.
For instance it may be used to calculate scattering states from the Green function. Once scattering states have been calculated one may also calculate the eigenchannels.
Both calculations are very efficient and uses very little memory compared to the full matrices normally used.
Consider a regular 2 electrode setup with transport direction along the 3rd lattice vector. Then the following example may be used to calculate the eigen-channels:
import sisl from sisl_toolbox.btd import * # First read in the required data H_elec = sisl.Hamiltonian.read("ELECTRODE.nc") H = sisl.Hamiltonian.read("DEVICE.nc") # remove couplings along the self-energy direction # to ensure no fake couplings. H.set_nsc(c=1) # Read in a single tbtrans output which contains the BTD matrices # and instructs this class how it should pivot the matrix to obtain # a BTD matrix. tbt = sisl.get_sile("siesta.TBT.nc") # Define the self-energy calculators which will downfold the # self-energies into the device region. # Since a downfolding will be done it requires the device Hamiltonian. H_elec.shift(tbt.mu("Left")) left = DownfoldSelfEnergy("Left", s.RecursiveSI(H_elec, "-C", eta=tbt.eta("Left"), tbt, H) H_elec.shift(tbt.mu("Right") - tbt.mu("Left")) left = DownfoldSelfEnergy("Right", s.RecursiveSI(H_elec, "+C", eta=tbt.eta("Right"), tbt, H) G = DeviceGreen(H, [left, right], tbt) # Calculate the scattering state from the left electrode # and then the eigen channels to the right electrode state = G.scattering_state("Left", E=0.1) eig_channel = G.eigenchannel(state, "Right")
To make this easier there exists a short-hand version that does the above:
G = DeviceGreen.from_fdf("RUN.fdf")
which reads all variables from the FDF file and parses them accordingly. This does not take all things into consideration, but should cover most problems.
Methods
Hk
(*args, **kwargs)Sk
(*args, **kwargs)eigenchannel
(state, elec_to[, ret_coeff])Calculate the eigenchannel from scattering states entering electrodes elec_to
from_fdf
(fdf[, prefix, use_tbt_se, eta])Return a new
DeviceGreen
using information gathered from the fdfgreen
(E[, k, format])Calculate the Green function for a given E and k point
reset
()Clean any memory used by this object
scattering_matrix
(elec_from, elec_to, E[, ...])Calculate the scattering matrix (S-matrix) between elec_from and elec_to
scattering_state
(elec, E[, k, cutoff, method])Calculate the scattering states for a given E and k point from a given electrode
spectral
(elec, E[, k, format, method, herm])Calculate the spectral function for a given E and k point from a given electrode
- eigenchannel(state, elec_to, ret_coeff=False)[source]
Calculate the eigenchannel from scattering states entering electrodes elec_to
The energy and k-point is inferred from the state object as returned from
scattering_state
.The eigenchannels are the eigen states of the transmission matrix in the DOS weighted scattering states:
\[\begin{split}\mathbf A_{\mathfrak e_{\mathrm{from}} }(E,\mathbf k) \mathbf u_i &= 2\pi a_i \mathbf u_i \\ \mathbf t_{\mathbf u} &= \sum \langle \mathbf u | \boldsymbol\Gamma_{ \mathfrak e_{\mathrm{to}} } | \mathbf u\rangle\end{split}\]where the eigenvectors of \(\mathbf t_{\mathbf u}\) are the coefficients of the DOS weighted scattering states (\(\sqrt{2\pi a_i} u_i\)) for the individual eigen channels. The eigenvalues are the transmission eigenvalues. Further details may be found in [8].
- Parameters:
state (
sisl.physics.StateCElectron
) – the scattering states as obtained fromscattering_state
elec_to (
str
orint (list
ornot)
) – which electrodes to consider for the transmission eigenchannel decomposition (the sum in the above equation)ret_coeff (
bool
, optional) – return also the scattering state coefficients
- Returns:
T_eig (
sisl.physics.StateCElectron
) – the transmission eigenchannels, theT_eig.c
contains the transmission eigenvalues.coeff (
sisl.physics.StateElectron
) – coefficients of state that creates the transmission eigenchannels Only returned if ret_coeff is True. There is a one-to-one correspondance betweencoeff
andT_eig
(with a prefactor of \(\sqrt{2\pi}\)). This is equivalent to theT_eig
states in the scattering state basis.
- classmethod from_fdf(fdf, prefix='TBT', use_tbt_se=False, eta=None)[source]
Return a new
DeviceGreen
using information gathered from the fdf- Parameters:
fdf (
str
orfdfSileSiesta
) – fdf file to read the parameters fromprefix (
{"TBT", "TS"}
) – which prefix to use, if TBT it will prefer TBT prefix, but fall back to TS prefixes. If TS, only these prefixes will be used.use_tbt_se (
bool
, optional) – whether to use the TBT.SE.nc files for self-energies or calculate them on the fly.
- green(E, k=(0, 0, 0), format='array')[source]
Calculate the Green function for a given E and k point
The Green function is calculated as:
\[\mathbf G(E,\mathbf k) = \big[\mathbf S(\mathbf k) E - \mathbf H(\mathbf k) - \sum \boldsymbol \Sigma(E,\mathbf k)\big]^{-1}\]- Parameters:
E (
float
) – the energy to calculate at, may be a complex value.k (
array_like
, optional) – k-point to calculate the Green function atformat (
{"array", "btd", "bm", "bd", "sparse"}
) – return the matrix in a specific formatarray: a regular numpy array (full matrix)
btd: a block-matrix object with only the diagonals and first off-diagonals
bm: a block-matrix object with diagonals and all off-diagonals
bd: a block-matrix object with only diagonals (no off-diagonals)
sparse: a sparse-csr matrix for the sparse elements as found in the Hamiltonian
- scattering_matrix(elec_from, elec_to, E, k=(0, 0, 0), cutoff=0.0)[source]
Calculate the scattering matrix (S-matrix) between elec_from and elec_to
The scattering matrix is calculated as
\[\mathbf S_{\mathfrak e_{\mathrm{to}}\mathfrak e_{\mathrm{from}} }(E, \mathbf) = -\delta_{\alpha\beta} + i \tilde{\boldsymbol\Gamma}_{\mathfrak e_{\mathrm{to}}} \mathbf G \tilde{\boldsymbol\Gamma}_{\mathfrak e_{\mathrm{from}}}\]Here \(\tilde{\boldsymbol\Gamma}\) is defined as:
\[\begin{split}\boldsymbol\Gamma(E,\mathbf k) \mathbf u_i &= \lambda_i \mathbf u_i \\ \tilde{\boldsymbol\Gamma}(E,\mathbf k) &= \operatorname{diag}\{ \sqrt{\boldsymbol\lambda} \} \mathbf u\end{split}\]Once the scattering matrices have been calculated one can calculate the full transmission function
\[T_{\mathfrak e_{\mathrm{from}}\mathfrak e_{\mathrm{to}} }(E, \mathbf k) = \operatorname{Tr}\big[ \mathbf S_{\mathfrak e_{\mathrm{to}}\mathfrak e_{\mathrm{from}} }^\dagger \mathbf S_{\mathfrak e_{\mathrm{to}}\mathfrak e_{\mathrm{from}} }\big]\]- Parameters:
elec_from (
str
orint
) – the electrode where the scattering matrix originates fromelec_to (
str
orint
orlist of
) – where the scattering matrix ends in.E (
float
) – the energy to calculate at, may be a complex value.k (
array_like
, optional) – k-point to calculate the scattering matrix atcutoff (
float
, optional) – cutoff the eigen states of the broadening matrix that are below this value. I.e. only \(\lambda\) values above this value will be used. A too high value will remove too many eigen states and results will be wrong. A small value improves precision at the cost of bigger matrices.
- Returns:
scat_matrix (
numpy.ndarray
ortuple
ofnumpy.ndarray
) – for each elec_to a scattering matrix will be returned. Its dimensions will be depending on the cutoff value at the cost of precision.
- scattering_state(elec, E, k=(0, 0, 0), cutoff=0.0, method: str = 'svd:gamma', *args, **kwargs)[source]
Calculate the scattering states for a given E and k point from a given electrode
The scattering states are the eigen states of the spectral function:
\[\mathbf A_{\mathfrak e}(E,\mathbf k) \mathbf u_i = 2\pi a_i \mathbf u_i\]where \(a_i\) is the DOS carried by the \(i\)’th scattering state.
- Parameters:
elec (
str
orint
) – the electrode to calculate the spectral function fromE (
float
) – the energy to calculate at, may be a complex value.k (
array_like
, optional) – k-point to calculate the spectral function atcutoff (
float
, optional) – cutoff the returned scattering states at some DOS value. Any scattering states with normalized eigenvalues lower than cutoff are discarded. The normalization is done by dividing the eigenvalue with the number of orbitals in the device region. This normalization ensures the same cutoff value has roughly the same meaning for different size devices. Values above or close to 1e-5 should be used with care.method (
{"svd:gamma", "svd:A", "full"}
) – which method to use for calculating the scattering states. Use only thefull
method for testing purposes as it is extremely slow and requires a substantial amount of memory. Thesvd:gamma
is the fastests while retaining complete precision. Thesvd:A
may be even faster for very large systems with very little loss of precision, since it diagonalizes \(\mathbf A\) in the subspace of the electrode elec and reduces the propagated part of the spectral matrix.cutoff_elec (
float
, optional) – Only used formethod=svd:A
. The initial block of the spectral function is diagonalized and only eigenvectors with eigenvalues>=cutoff_elec
are retained, thus reducing the initial propagated modes. The normalization explained for cutoff also applies here.
- Returns:
scat (
StateCElectron
) – the scattering states from the spectral function. Thescat.state
contains the scattering state vectors (eigenvectors of the spectral function).scat.c
contains the DOS of the scattering states scaled by \(1/(2\pi)\) so ensure correct density of states. One may recreate the spectral function withscat.outer(matrix=scat.c * 2 * pi)
.
- spectral(elec, E, k=(0, 0, 0), format: str = 'array', method: str = 'column', herm: bool = True)[source]
Calculate the spectral function for a given E and k point from a given electrode
The spectral function is calculated as:
\[\mathbf A_{\mathfrak{e}}(E,\mathbf k) = \mathbf G(E,\mathbf k)\boldsymbol\Gamma_{\mathfrak{e}}(E,\mathbf k) \mathbf G^\dagger(E,\mathbf k)\]- Parameters:
elec (
str
orint
) – the electrode to calculate the spectral function fromE (
float
) – the energy to calculate at, may be a complex value.k (
array_like
, optional) – k-point to calculate the spectral function atformat (
{"array", "btd", "bm", "bd"}
) – return the matrix in a specific formatarray: a regular numpy array (full matrix)
btd: a block-matrix object with only the diagonals and first off-diagonals
bm: a block-matrix object with diagonals and all off-diagonals
bd: same as btd, since they are already calculated
method (
{"column", "propagate"}
) – which method to use for calculating the spectral function. Depending on the size of the BTD blocks one may be faster than the other. For large systems you are recommended to time the different methods and stick with the fastest one, they are numerically identical.herm – The spectral function is a Hermitian matrix, by default (True), the methods that can utilize the Hermitian property only calculates the lower triangular part of \(\mathbf A\), and then copies the Hermitian to the upper part. By setting this to False the entire matrix is explicitly calculated.