kwant.kpm
– Kernel Polynomial Method¶kwant.kpm.
Correlator
(hamiltonian, operator1=None, operator2=None, **kwargs)[source]¶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
Internally, the correlation is approximated with a two dimensional KPM expansion,
with coefficients
and moments matrix \(µ_{n m}^{α β} = \mathrm{Tr}(O_α T_m(H) O_β T_n(H))\).
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\).
FiniteSystem
or matrix HamiltonianIf a system is passed, it should contain no leads.
Operators to be passed to two different instances of
SpectralDensity
.
Keyword arguments to pass to SpectralDensity
.
Notes
The operator1
must act to the right as \(O_α\rvert r\rangle\).
__call__
(mu=0, temperature=0)[source]¶Returns the linear response \(χ_{α β}(µ, T)\)
Chemical potential defined in the same units of energy as the Hamiltonian.
Temperature in units of energy, the same as defined in the Hamiltonian.
add_moments
(num_moments=None, *, energy_resolution=None)[source]¶Increase the number of Chebyshev moments
The number of Chebyshev moments to add. Mutually exclusive with ‘energy_resolution’.
Features wider than this resolution are visible
in the spectral density. Mutually exclusive with
num_moments
.
kwant.kpm.
LocalVectors
(syst, where=None, *args)[source]¶Generates a local vector iterator for the sites in ‘where’.
FiniteSystem
or matrix HamiltonianIf a system is passed, it should contain no leads.
Site
, or callable, optionalSpatial 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.
kwant.kpm.
RandomVectors
(syst, where=None, rng=None)[source]¶Returns a random phase vector iterator for the sites in ‘where’.
FiniteSystem
or matrix HamiltonianIf a system is passed, it should contain no leads.
Site
, or callable, optionalSpatial 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.
kwant.kpm.
SpectralDensity
(hamiltonian, params=None, operator=None, num_vectors=10, num_moments=None, energy_resolution=None, vector_factory=None, bounds=None, eps=0.05, rng=None, kernel=None, mean=True, accumulate_vectors=True)[source]¶Calculate the spectral density of an operator.
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
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\).
FiniteSystem
or matrix HamiltonianIf a system is passed, it should contain no leads.
Additional parameters to pass to the Hamiltonian and operator.
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.
Number of vectors used in the KPM expansion. If None
, the
number of vectors used equals the length of the ‘vector_factory’.
Number of moments, order of the KPM expansion. Mutually exclusive
with energy_resolution
.
The resolution in energy of the KPM approximation to the spectral
density. Mutually exclusive with num_moments
.
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].
Lower and upper bounds for the eigenvalue spectrum of the system. If not provided, they are computed.
Parameter to ensure that the rescaled spectrum lies in the
interval (-1, 1)
; required for stability.
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.
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
.
True
If True
, return the mean spectral density for the vectors
used, otherwise return an array of densities for each vector.
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
,
>>> import kwant
...
>>> def circle(pos):
... x, y = pos
... return x**2 + y**2 < 100
...
>>> lat = kwant.lattice.honeycomb()
>>> syst = kwant.Builder()
>>> syst[lat.shape(circle, (0, 0))] = 0
>>> syst[lat.neighbors()] = -1
and after finalizing the system, create an instance of
SpectralDensity
>>> fsyst = syst.finalized()
>>> rho = kwant.kpm.SpectralDensity(fsyst)
The energies
and densities
can be accessed with
>>> energies, densities = rho()
or
>>> energies, densities = rho.energies, rho.densities
Array of sampling points with length 2 * num_moments
in
the range of the spectrum.
Spectral density of the operator
evaluated at the energies.
__call__
(energy=None)[source]¶Return the spectral density evaluated at energy
.
Drawn from the nodes of the highest Chebyshev polynomial. Not returned if ‘energy’ was not provided
single float
if the energy
parameter is a single
float
, else an array of float
.
Notes
If energy
is not provided, then the densities are obtained
by Fast Fourier Transform of the Chebyshev moments.
add_moments
(num_moments=None, *, energy_resolution=None)[source]¶Increase the number of Chebyshev moments.
The number of Chebyshev moments to add. Mutually
exclusive with energy_resolution
.
Features wider than this resolution are visible
in the spectral density. Mutually exclusive with
num_moments
.
kwant.kpm.
conductivity
(hamiltonian, alpha='x', beta='x', positions=None, **kwargs)[source]¶Returns a callable object to obtain the elements of the conductivity tensor using the Kubo-Bastin approach.
A Correlator
instance is created to obtain the
correlation between two components of the current operator
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]
FiniteSystem
or matrix HamiltonianIf a system is passed, it should contain no leads.
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.
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’}.
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
,
>>> import kwant
...
>>> def circle(pos):
... x, y = pos
... return x**2 + y**2 < 100
...
>>> lat = kwant.lattice.honeycomb()
>>> syst = kwant.Builder()
>>> syst[lat.shape(circle, (0, 0))] = 0
>>> syst[lat.neighbors()] = -1
>>> syst[lat.a.neighbors()] = -0.5j
>>> syst[lat.b.neighbors()] = 0.5j
>>> fsyst = syst.finalized()
Now we can call conductivity
to calculate the transverse
conductivity at chemical potential 0 and temperature 0.01.
>>> cond = kwant.kpm.conductivity(fsyst, alpha='x', beta='y')
>>> cond(mu=0, temperature=0.01)
kwant.kpm.
fermi_distribution
(energy, mu, temperature)[source]¶Returns the Fermi distribution f(e, µ, T) evaluated at ‘e’.
Energy array where the Fermi distribution is evaluated.
Chemical potential defined in the same units of energy as the Hamiltonian.
Temperature in units of energy, the same as defined in the Hamiltonian.
kwant.kpm.
jackson_kernel
(moments)[source]¶Convolutes moments
with the Jackson kernel.
Taken from Eq. (71) of Rev. Mod. Phys., Vol. 78, No. 1 (2006).
kwant.kpm.
lorentz_kernel
(moments, l=4)[source]¶Convolutes moments
with the Lorentz kernel.
Taken from Eq. (71) of Rev. Mod. Phys., Vol. 78, No. 1 (2006).
The additional parameter l
controls the decay of the kernel.