Merge pull request #9051 from rwolst/sparse_rand_slow_fixe

MAINT: Fix slow sparse.rand for k < mn/3 (#9036).

benchmarks/benchmarks/sparse.py

@@ -423,3 +423,20 @@ def time_toarray(self, format, order):
 
     # Retain old benchmark results (remove this if changing the benchmark)
     time_toarray.version = "2fbf492ec800b982946a62785beda803460b913cc80080043a5d407025893b2b"
+
+
+class Random(Benchmark):
+    params = [
+        np.arange(0, 1.1, 0.1).tolist()
+    ]
+    param_names = ['density']
+
+    def setup(self, density):
+        warnings.simplefilter('ignore', SparseEfficiencyWarning)
+        self.nrows = 1000
+        self.ncols = 1000
+        self.format = 'csr'
+
+    def time_rand(self, density):
+        sparse.rand(self.nrows, self.ncols,
+                    format=self.format, density=density)

scipy/optimize/tests/test_lsq_linear.py

@@ -130,13 +130,13 @@ def test_sparse_bounds(self):
         lb = self.rnd.randn(n)
         ub = lb + 1
         res = lsq_linear(A, b, (lb, ub))
-        assert_allclose(res.optimality, 0.0, atol=1e-8)
+        assert_allclose(res.optimality, 0.0, atol=1e-6)
 
         res = lsq_linear(A, b, (lb, ub), lsmr_tol=1e-13)
-        assert_allclose(res.optimality, 0.0, atol=1e-8)
+        assert_allclose(res.optimality, 0.0, atol=1e-6)
 
         res = lsq_linear(A, b, (lb, ub), lsmr_tol='auto')
-        assert_allclose(res.optimality, 0.0, atol=1e-8)
+        assert_allclose(res.optimality, 0.0, atol=1e-6)
 
 
 class TestTRF(BaseMixin, SparseMixin):

scipy/sparse/construct.py

@@ -711,14 +711,20 @@ def random(m, n, density=0.01, format='coo', dtype=None,
     -------
     res : sparse matrix
 
+    Notes
+    -----
+    Only float types are supported for now.
+
     Examples
     --------
     >>> from scipy.sparse import random
     >>> from scipy import stats
-    >>> class CustomRandomState(object):
+
+    >>> class CustomRandomState(np.random.RandomState):
     ...     def randint(self, k):
     ...         i = np.random.randint(k)
     ...         return i - i % 2
+    >>> np.random.seed(12345)
     >>> rs = CustomRandomState()
     >>> rvs = stats.poisson(25, loc=10).rvs
     >>> S = random(3, 4, density=0.25, random_state=rs, data_rvs=rvs)
@@ -736,13 +742,10 @@ def random(m, n, density=0.01, format='coo', dtype=None,
     >>> Y = X()  # get a frozen version of the distribution
     >>> S = random(3, 4, density=0.25, random_state=2906, data_rvs=Y.rvs)
     >>> S.A
-    array([[ 0.        ,  1.9467163 ,  0.13569738, -0.81205367],
-           [ 0.        ,  0.        ,  0.        ,  0.        ],
+    array([[ 0.        ,  0.        ,  0.        ,  0.        ],
+           [ 0.13569738,  1.9467163 , -0.81205367,  0.        ],
            [ 0.        ,  0.        ,  0.        ,  0.        ]])
 
-    Notes
-    -----
-    Only float types are supported for now.
     """
     if density < 0 or density > 1:
         raise ValueError("density expected to be 0 <= density <= 1")
@@ -768,6 +771,7 @@ def random(m, n, density=0.01, format='coo', dtype=None,
         random_state = np.random
     elif isinstance(random_state, (int, np.integer)):
         random_state = np.random.RandomState(random_state)
+
     if data_rvs is None:
         if np.issubdtype(dtype, np.integer):
             randint = get_randint(random_state)
@@ -781,18 +785,7 @@ def data_rvs(n):
         else:
             data_rvs = random_state.rand
 
-    # Use the algorithm from python's random.sample for k < mn/3.
-    if mn < 3 * k:
-        ind = random_state.choice(mn, size=k, replace=False)
-    else:
-        ind = np.empty(k, dtype=tp)
-        selected = set()
-        for i in xrange(k):
-            j = random_state.randint(mn)
-            while j in selected:
-                j = random_state.randint(mn)
-            selected.add(j)
-            ind[i] = j
+    ind = random_state.choice(mn, size=k, replace=False)
 
     j = np.floor(ind * 1. / m).astype(tp, copy=False)
     i = (ind - j * m).astype(tp, copy=False)
@@ -841,8 +834,9 @@ def rand(m, n, density=0.01, format="coo", dtype=None, random_state=None):
     <3x4 sparse matrix of type '<class 'numpy.float64'>'
        with 3 stored elements in Compressed Sparse Row format>
     >>> matrix.todense()
-    matrix([[ 0.        ,  0.59685016,  0.779691  ,  0.        ],
-            [ 0.        ,  0.        ,  0.        ,  0.44583275],
-            [ 0.        ,  0.        ,  0.        ,  0.        ]])
+    matrix([[0.05641158, 0.        , 0.        , 0.65088847],
+            [0.        , 0.        , 0.        , 0.14286682],
+            [0.        , 0.        , 0.        , 0.        ]])
+
     """
     return random(m, n, density, format, dtype, random_state)