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`.

kwant/kpm.py

@@ -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]))

kwant/tests/test_kpm.py

@@ -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.