Commit 11155d1d authored by Pablo Piskunow's avatar Pablo Piskunow Committed by Joseph Weston
Browse files

switch to argument passing via 'param' in kpm

Passing arguments via 'args' will be deprecated in a future
release, and will be removed in Kwant 2.0, so any new features
should only support 'params'.

We also now properly 'bind' parameters to the passed operator;
previously the parameters were not passed correctly.
parent 4c41b5df
Pipeline #3124 passed with stages
in 3 minutes and 40 seconds
......@@ -18,7 +18,7 @@ import scipy.fftpack as fft
from . import system
from ._common import ensure_isinstance, ensure_rng
from .operator import _LocalOperator
class SpectralDensity:
"""Calculate the spectral density of an operator.
......@@ -40,8 +40,8 @@ class SpectralDensity:
----------
syst_or_ham : `~kwant.system.FiniteSystem` or matrix Hamiltonian
If a system is passed, it should contain no leads.
args : tuple, optional
Positional arguments to pass to the system.
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
......@@ -132,8 +132,8 @@ class SpectralDensity:
Spectral density of the ``operator`` evaluated at the energies.
"""
def __init__(self, syst_or_ham, args=(), operator=None, num_rand_vecs=10,
num_moments=100, num_sampling_points=None,
def __init__(self, syst_or_ham, params=None, operator=None,
num_rand_vecs=10, num_moments=100, num_sampling_points=None,
vector_factory=None, bounds=None, epsilon=0.05, rng=None):
rng = ensure_rng(rng)
# self.epsilon ensures that the rescaled Hamiltonian has a
......@@ -146,24 +146,27 @@ class SpectralDensity:
try:
ensure_isinstance(syst_or_ham, system.System)
ham = scipy.sparse.csr_matrix(
syst_or_ham.hamiltonian_submatrix(args=args, sparse=True))
syst_or_ham.hamiltonian_submatrix(params=params,
sparse=True))
except TypeError:
raise ValueError('Parameter `syst_or_ham` is not a '
'Hamiltonian neither a (finalized) '
'`kwant.system`.')
# Check if operator is a sparse matrix or None.
# Normalize 'operator' to a common format.
if operator is None:
self.operator = None
elif isinstance(operator, _LocalOperator):
self.operator = operator.bind(params=params)
elif callable(operator):
self.operator = operator
elif hasattr(operator, 'dot'):
operator = scipy.sparse.csr_matrix(operator)
self.operator = lambda bra, ket: np.vdot(bra, operator.dot(ket))
else:
if callable(operator):
self.operator = operator
elif hasattr(operator, 'dot') and not hasattr(operator, 'act'):
operator = scipy.sparse.csr_matrix(operator)
self.operator = (lambda bra, ket:
np.vdot(bra, operator.dot(ket)))
else:
raise ValueError('Parameter `operator` has no `.dot` '
'attribute and is not callable.')
raise ValueError('Parameter `operator` has no `.dot` '
'attribute and is not callable.')
self.num_moments = num_moments
# Default number of sampling points
if num_sampling_points is None:
......
......@@ -43,7 +43,7 @@ def assert_allclose_sp(arr1, arr2):
np.testing.assert_allclose(arr1, arr2, rtol=0., atol=TOL_SP)
def make_spectrum(ham, p, operator=None, vector_factory=None, rng=None):
def make_spectrum(ham, p, operator=None, vector_factory=None, rng=None, params=None):
"""Create an instance of SpectralDensity class."""
return SpectralDensity(
ham,
......@@ -52,7 +52,8 @@ def make_spectrum(ham, p, operator=None, vector_factory=None, rng=None):
num_moments=p.num_moments,
num_rand_vecs=p.num_rand_vecs,
num_sampling_points=p.num_sampling_points,
rng=rng
rng=rng,
params=params
)
......@@ -65,6 +66,17 @@ def make_chain(r=dim, t=-1):
return syst.finalized()
def make_chain_with_params(r=10, pot=0, t=-1):
syst = kwant.Builder()
lat = kwant.lattice.chain(norbs=1)
pot = lambda site, pot: pot
hop = lambda site1, site2, t: t
for i in range(r):
syst[lat(i)] = pot
syst[lat.neighbors()] = hop
return syst.finalized()
def kpm_derivative(spectrum, e, order=1):
"""Calculate the coeficients for Chebyshev expansion and its derivates."""
e = np.asarray(e, dtype=complex)
......@@ -289,6 +301,32 @@ def test_kwant_op():
e = spectrum_syst.energies
assert_allclose_sp(spectrum_syst(e), spectrum_syst.densities)
def test_kwant_op_current():
"""Check that the kwant.operator.Density gives the same result as the
identity operator when the system has ``params``.
"""
params = {'r':dim, 'pot':1, 't':-2}
# build the system using parameters default values
# to be later rewritten
syst = make_chain_with_params()
current = kwant.operator.Current(syst)
# pass parameters to kpm to evaluate hamiltonian and operator
spectrum_syst = make_spectrum(syst, p, operator=current, rng=1,
params=params)
# compare with explicit hamiltonian and operator
ham = syst.hamiltonian_submatrix(params=params)
def my_current(bra, ket):
return current(bra, ket, params=params)
# do not pass parameters to kpm
spectrum = make_spectrum(ham, p, operator=my_current, rng=1)
# same algorithms are used so these arrays are equal up to TOL
assert_allclose(spectrum_syst.densities, spectrum.densities)
def test_kwant_op_average():
"""Check that the kwant.operator.Density gives the same result as the
identity operator.
......
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