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

fix bug in passing `energy_resolution` instead of `num_moments`

  When passing `energy_resolution`, the number of moments is set using
  the yet undefined parameter `self._a`. Now the number of moments is
  set after defining `self._a`.
parent 939b57a1
......@@ -141,19 +141,6 @@ class SpectralDensity:
raise TypeError("either 'num_moments' or 'energy_resolution' "
"must be provided.")
if energy_resolution:
num_moments = math.ceil((1.6 * self._a) / energy_resolution)
elif num_moments is None:
num_moments = 100
must_be_positive_int = ['num_vectors', 'num_moments']
for var in must_be_positive_int:
val = locals()[var]
if val <= 0 or val != int(val):
raise ValueError('{} must be a positive integer'.format(var))
if eps <= 0:
raise ValueError('eps must be positive')
rng = ensure_rng(rng)
# self.eps ensures that the rescaled Hamiltonian has a
# spectrum strictly in the interval (-1,1).
......@@ -195,6 +182,19 @@ class SpectralDensity:
bounds=bounds)
self.bounds = (self._b - self._a, self._b + self._a)
if energy_resolution:
num_moments = math.ceil((1.6 * self._a) / energy_resolution)
elif num_moments is None:
num_moments = 100
must_be_positive_int = ['num_vectors', 'num_moments']
for var in must_be_positive_int:
val = locals()[var]
if val <= 0 or val != int(val):
raise ValueError('{} must be a positive integer'.format(var))
if eps <= 0:
raise ValueError('eps must be positive')
for r in range(num_vectors):
self._rand_vect_list.append(
self._vector_factory(self.hamiltonian.shape[0]))
......
......@@ -175,6 +175,24 @@ def test_api_single_eigenvalue_error():
SpectralDensity(np.identity(dim, dtype=complex))
def test_energy_resolution():
"""Check that energy resolution works and gives the same output as
setting the equivalent number of moments.
"""
ham = kwant.rmt.gaussian(dim)
rng = 1
energy_resolution = 0.05
sp1 = SpectralDensity(ham, rng=rng, energy_resolution=energy_resolution)
# number of moments for this energy resolution
num_moments = sp1.num_moments
# use the same number of moments but not passing energy resolution
sp2 = SpectralDensity(ham, rng=rng, num_moments=num_moments)
# Check bit for bit equality
assert np.all(sp1.densities == sp2.densities)
def test_operator_none():
"""Check operator=None gives the same results as operator=np.identity(),
with the same random vectors.
......
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