Calculates the response of the correlation between two operators.
The response tensor \(χ_{α β}\) of an operator \(O_α\)
to a perturbation in an operator \(O_β\), is defined here based
on [3]_, and [4]_, and takes the form
\[ \begin{align}\begin{aligned}Γ_{m n}(E) =
(E - i n \sqrt{1 - E^2}) e^{i n \arccos(E)} T_m(E)\\+ (E + i m \sqrt{1 - E^2}) e^{-i m \arccos(E)} T_n(E),\end{aligned}\end{align} \]
The trace is calculated creating two instances of
SpectralDensity, and saving the vectors
\(Ψ_{n r} = O_β T_n(H)\rvert r\rangle\),
and \(Ω_{m r} = T_m(H) O_α\rvert r\rangle\) , where
\(\rvert r\rangle\) is a vector provided by the
vector_factory.
The moments matrix is formed with the product
\(µ_{m n} = \langle Ω_{m r} \rvert Ψ_{n r}\rangle\) for
every \(\rvert r\rangle\).
Parameters:
hamiltonian (FiniteSystem or matrix Hamiltonian) – If a system is passed, it should contain no leads.
operator1 (operators, dense matrix, or sparse matrix, optional) – Operators to be passed to two different instances of
SpectralDensity.
operator2 (operators, dense matrix, or sparse matrix, optional) – Operators to be passed to two different instances of
SpectralDensity.
**kwargs (dict) – Keyword arguments to pass to SpectralDensity.
Notes
The operator1 must act to the right as \(O_α\rvert r\rangle\).
num_moments (positive int, optional) – The number of Chebyshev moments to add. Mutually
exclusive with ‘energy_resolution’.
energy_resolution (positive float, optional) – Features wider than this resolution are visible
in the spectral density. Mutually exclusive with
num_moments.
Generates a local vector iterator for the sites in ‘where’.
Parameters:
syst (FiniteSystem or matrix Hamiltonian) – If a system is passed, it should contain no leads.
where (sequence of int or Site, or callable, optional) – Spatial range of the vectors produced. If syst is a
FiniteSystem, where behaves as in
Density. If syst is a matrix, where
must be a list of integers with the indices where column vectors
are nonzero.
Returns a random phase vector iterator for the sites in ‘where’.
Parameters:
syst (FiniteSystem or matrix Hamiltonian) – If a system is passed, it should contain no leads.
where (sequence of int or Site, or callable, optional) – Spatial range of the vectors produced. If syst is a
FiniteSystem, where behaves as in
Density. If syst is a matrix, where
must be a list of integers with the indices where column vectors
are nonzero.
This class makes use of the kernel polynomial
method (KPM), presented in [1], to obtain the spectral density
\(ρ_A(e)\), as a function of the energy \(e\), of some
operator \(A\) that acts on a kwant system or a Hamiltonian.
In general
\[ρ_A(E) = ρ(E) A(E),\]
where \(ρ(E) = \sum_{k=0}^{D-1} δ(E-E_k)\) is the density of
states, and \(A(E)\) is the expectation value of \(A\) for
all the eigenstates with energy \(E\).
Parameters:
hamiltonian (FiniteSystem or matrix Hamiltonian) – If a system is passed, it should contain no leads.
params (dict, optional) – Additional parameters to pass to the Hamiltonian and operator.
operator (operator, dense matrix, or sparse matrix, optional) – Operator for which the spectral density will be evaluated. If
it is callable, the densities at each energy will have the
dimension of the result of operator(bra,ket). If it has a
dot method, such as numpy.ndarray and
scipy.sparse.matrices, the densities will be scalars.
num_vectors (positive int, or None, default: 10) – Number of vectors used in the KPM expansion. If None, the
number of vectors used equals the length of the ‘vector_factory’.
num_moments (positive int, default: 100) – Number of moments, order of the KPM expansion. Mutually exclusive
with energy_resolution.
energy_resolution (positive float, optional) – The resolution in energy of the KPM approximation to the spectral
density. Mutually exclusive with num_moments.
vector_factory (iterable, optional) – If provided, it should contain (or yield) vectors of the size of
the system. If not provided, random phase vectors are used.
The default random vectors are optimal for most cases, see the
discussions in [1] and [2].
bounds (pair of floats, optional) – Lower and upper bounds for the eigenvalue spectrum of the system.
If not provided, they are computed.
eps (positive float, default: 0.05) – Parameter to ensure that the rescaled spectrum lies in the
interval (-1,1); required for stability.
rng (seed, or random number generator, optional) – Random number generator used for the calculation of the spectral
bounds, and to generate random vectors (if vector_factory is
not provided). If not provided, numpy’s rng will be used; if it
is an Integer, it will be used to seed numpy’s rng, and if it is
a random number generator, this is the one used.
kernel (callable, optional) – Callable that takes moments and returns stabilized moments.
By default, the jackson_kernel is used.
The Lorentz kernel is also accesible by passing
lorentz_kernel.
mean (bool, default: True) – If True, return the mean spectral density for the vectors
used, otherwise return an array of densities for each vector.
accumulate_vectors (bool, default: True) – Whether to save or discard each vector produced by the vector
factory. If it is set to False, it is not possible to
increase the number of moments, but less memory is used.
Notes
When passing an operator defined in operator, the
result of operator(bra,ket) depends on the attribute sum
of such operator. If sum=True, densities will be scalars, that
is, total densities of the system. If sum=False the densities
will be arrays of the length of the system, that is, local
densities.
Examples
In the following example, we will obtain the density of states of a
graphene sheet, defined as a honeycomb lattice with first nearest
neighbors coupling.
We start by importing kwant and defining a
FiniteSystem,
num_moments (positive int) – The number of Chebyshev moments to add. Mutually
exclusive with energy_resolution.
energy_resolution (positive float, optional) – Features wider than this resolution are visible
in the spectral density. Mutually exclusive with
num_moments.
Returns the integral over the whole spectrum with an optional
distribution function. distribution_function should be able
to take arrays as input. Defined using Gauss-Chebyshev
integration.
where \(V\) is the volume where the conductivity is sampled.
In this implementation it is assumed that the vectors are normalized
and \(V=1\), otherwise the result of this calculation must be
normalized with the corresponding volume.
The equations used here are based on [3]_ and [4]_
Parameters:
hamiltonian (FiniteSystem or matrix Hamiltonian) – If a system is passed, it should contain no leads.
alpha (str, or operators) – If hamiltonian is a kwant system, or if the positions
are provided, alpha and beta can be the directions of the
velocities as strings {‘x’, ‘y’, ‘z’}.
Otherwise alpha and beta should be the proper velocity
operators, which can be members of operator or matrices.
beta (str, or operators) – If hamiltonian is a kwant system, or if the positions
are provided, alpha and beta can be the directions of the
velocities as strings {‘x’, ‘y’, ‘z’}.
Otherwise alpha and beta should be the proper velocity
operators, which can be members of operator or matrices.
positions (array of float, optioinal) – If hamiltonian is a matrix, the velocities can be calculated
internally by passing the positions of each orbital in the system
when alpha or beta are one of the directions {‘x’, ‘y’, ‘z’}.
**kwargs (dict) – Keyword arguments to pass to Correlator.
Examples
We will obtain the conductivity of the Haldane model, defined as a
honeycomb lattice with first nearest neighbors coupling, and
imaginary second nearest neighbors coupling.
We start by importing kwant and defining a
FiniteSystem,