Issue 1:
Encountering a MemoryError when creating a large sparse random matrix.

The matrix is created with scipy.sparse.rand, and even if the density is low enough that only a few elements in the whole matrix should be nonzero (or even none at all), it still tries to allocate a large integer array.

I'm not sure if this really is a bug, or I'm just using the function wrong. I know it's possible to create an empty sparse matrix and fill it with random values manually. But I thought it a bit strange that it would try to allocate that much.

Reproducing code example:

from scipy.sparse import rand
rand(10**5, 10**5, density=0)

Error message:

---------------------------------------------------------------------------
MemoryError Traceback (most recent call last)
<ipython-input-68-4d3ba46fe0a7> in <module>
----> 1 rand(10**5, 10**5, density=0)

~/anaconda3/lib/python3.7/site-packages/scipy/sparse/construct.py in rand(m, n, density, format, dtype, random_state)
840
841 """
--> 842 return random(m, n, density, format, dtype, random_state)

~/anaconda3/lib/python3.7/site-packages/scipy/sparse/construct.py in random(m, n, density, format, dtype, random_state, data_rvs)
786 data_rvs = random_state.rand
787
--> 788 ind = random_state.choice(mn, size=k, replace=False)
789
790 j = np.floor(ind * 1. / m).astype(tp, copy=False)

mtrand.pyx in numpy.random.mtrand.RandomState.choice()

mtrand.pyx in numpy.random.mtrand.RandomState.permutation()

MemoryError: Unable to allocate 74.5 GiB for an array with shape (10000000000,) and data type int64

Scipy/Numpy/Python version information:

1.4.1 1.18.1 sys.version_info(major=3, minor=7, micro=6, releaselevel='final', serial=0)

Issue 2:
gh-9036 removed a special case in the scipy.sparse.random code, as the numpy implementation of random.choice(..., replace=False) had improved significantly.

I think we may need to reintroduce the special case, but only for more extreme cases. A quick test shows that np.random.choice(N, k, replace=False) scales linearly in N:

In [2]: np.__version__
Out[2]: '1.15.4'

In [3]: %timeit np.random.choice(100, 10, replace=False)
17 µs ± 62.3 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each)

In [4]: %timeit np.random.choice(1000, 10, replace=False)
28.9 µs ± 204 ns per loop (mean ± std. dev. of 7 runs, 10000 loops each)

In [5]: %timeit np.random.choice(10000, 10, replace=False)
152 µs ± 1.17 µs per loop (mean ± std. dev. of 7 runs, 10000 loops each)

In [6]: %timeit np.random.choice(100000, 10, replace=False)
1.45 ms ± 6.24 µs per loop (mean ± std. dev. of 7 runs, 1000 loops each)

In [7]: %timeit np.random.choice(1000000, 10, replace=False)
14.6 ms ± 101 µs per loop (mean ± std. dev. of 7 runs, 100 loops each)

The result I've observed is that scipy.sparse.random gets unbearably slow. For example, when running the example code in gh-9697 on scipy master, constructing the matrix takes over 2 minutes and completely locks my machine for the duration of the np.random.choice call.

Perhaps the new logic should be:

if mn < 1e8 or k > 1e4:
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