Thanks for your report and interest in contributing to CuPy. I'd like to see the output you're reporint, so would you share some input to reproduce the result with us?

Right now, for example with this code:

import cupy as cp
import cupyx.scipy.sparse as cupy_sparse
import cupyx.scipy.sparse.linalg as cl
import numpy as np

a = np.random.rand(80, 100)

a_cupy = cp.array(a)

a_sparse = cupy_sparse.csr_matrix(a_cupy)

u, s, vt = cl.svds(a_sparse, k=77, return_singular_vectors=True)

print(u, s, vt)

The max value I can give k for getting result is k=77 which gives me at max 77 singular values. This is the error for k=78, 79:

---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
[<ipython-input-33-cda5c70b8c85>](https://localhost:8080/#) in <cell line: 11>()
      9 a_sparse = cupy_sparse.csr_matrix(a_cupy)
     10 
---> 11 u, s, vt = cl.svds(a_sparse, k=79, return_singular_vectors=True)
     12 
     13 print(u, s, vt)

1 frames
[/usr/local/lib/python3.10/dist-packages/cupyx/scipy/sparse/linalg/_eigen.py](https://localhost:8080/#) in eigsh(a, k, which, ncv, maxiter, tol, return_eigenvectors)
    112         cublas.gemv(_cublas.CUBLAS_OP_C, 1, V[:k].T, u, 0, uu)
    113         cublas.gemv(_cublas.CUBLAS_OP_N, -1, V[:k].T, uu, 1, u)
--> 114         V[k] = u / cublas.nrm2(u)
    115 
    116         u[...] = a @ V[k]

cupy/_core/core.pyx in cupy._core.core._ndarray_base.__setitem__()

cupy/_core/_routines_indexing.pyx in cupy._core._routines_indexing._ndarray_setitem()

cupy/_core/_routines_indexing.pyx in cupy._core._routines_indexing._scatter_op()

cupy/_core/_routines_indexing.pyx in cupy._core._routines_indexing._view_getitem()

IndexError: Index 79 is out of bounds for axis 0 with size 79

and for k=80 the code is not compatible.

I want to make it work for k=80 (which is min(80, 100)) for a particular use case in a project. Please guide me how should I proceed with this?

Thanks for sharing the reproducer. I can see the result you reported.

Seems that the default value of ncv, i.e. the number of generated Lanczos vectors, is slightly different from SciPy's and it makes the out of bounds error.

CuPy:

SciPy:
https://github.com/scipy/scipy/blob/297f6cd82d196cca73e11c8153c9fad34c0fcf1b/scipy/sparse/linalg/_eigen/arpack/arpack.py#L1448

If you're interested, would you make a PR that fixes this issue as well as tests the case? We have tests for svds in tests/cupyx_tests/scipy_tests/sparse_tests/test_linalg.py

There are some pointers to contribute to CuPy:

• https://docs.cupy.dev/en/stable/install.html#installing-cupy-from-source
• https://docs.cupy.dev/en/stable/contribution.html#unit-testing

Thanks for this information, I have rectified it in my fork and will create a PR for it by today.

Can you also look into this:

I want to make it work for k=80 (which is min(80, 100)) for a particular use case in a project. Please guide me how should I proceed with this?

scipy has 3 solvers implemented and for solver='propack' we are able to calculate all the singular values (k=min(n,m)) of the matrice, however in cupy k is always lesser then min(n,m). Can we add a functionality to maybe allow k upto min(n,m)?

propack seems to be a fortran library. We could add the solver to CuPy, but it would require more effort than naive SciPy to CuPy porting in Python layer. Another point I'm not sure is if the solver's algorithm well suits on GPU or not.

Hi @takagi

Could you please share some resources to help me get started on this? It's currently a major bottleneck in our pipeline.

You can find proback's code below:
https://github.com/scipy/PROPACK

Porting proback to cupy may not be an easy task, and I'm not sure if proback on GPU have enough performance improvement you need in your lab. Are the sparse matrices you're handling too large for computing on CPU with Scipy?

@takagi Actually, we were using it for some CFD simulations, and in our use case number of singular values should be equal to min(m,n) for further processing.

Thanks for your help.