Use scikit-cuda to compute singular value decomposition with cuSOLVER

Question

I am trying to use scikit-cuda's wrappers for the cuSOLVER functions, in particular I want to execute cusolverDnSgesvd to compute full-matrix single precision SVD on a matrix of real numbers.

Using the code here and here as a reference, I managed to get this far:

import pycuda.autoinit
import pycuda.driver as drv
import pycuda.gpuarray as gpuarray
import numpy as np

from skcuda import cusolver


handle = cusolver.cusolverDnCreate()

m = 50
n = 25

a = np.asarray(np.random.random((m, n)))
a_gpu = gpuarray.to_gpu(a)

ldu = m
ldvt = n

s_gpu = gpuarray.empty(min(m, n), np.float32)
u_gpu = gpuarray.empty((ldu, m), np.float32)
vh_gpu = gpuarray.empty((n, n), np.float32)

work_size = cusolver.cusolverDnSgesvd_bufferSize(handle, m, n)

work = gpuarray.empty((m,n), np.float32)

u_gpu, s_gpu, vh_gpu = cusolver.cusolverDnSgesvd(
    handle=handle,
    jobu='A',
    jobvt='A',
    m=m,
    n=n,
    A=a,
    lda=m,
    S=s_gpu,
    U=u_gpu,
    ldu=ldu,
    VT=vh_gpu,
    ldvt=ldvt,
    Work=work,
    Lwork=work_size,
    rwork=None,
    devInfo=0
)

But the code isn't working, probably because I'm messing up with types.

Traceback (most recent call last):
  File "/home/vektor/PycharmProjects/yancut/test_svd.py", line 44, in <module>
    devInfo=0
  File "/home/vektor/anaconda3/lib/python3.4/site-packages/skcuda/cusolver.py", line 577, in cusolverDnSgesvd
    int(A), lda, int(S), int(U),
TypeError: only length-1 arrays can be converted to Python scalars

How should I provide all the arguments so that the SVD is executed in a proper way?

UPDATE1: After using this question as reference, I edited my code and I'm getting a new error.

import pycuda.autoinit
import pycuda.driver as drv
import pycuda.gpuarray as gpuarray
import numpy as np

import ctypes

from skcuda import cusolver

rows = 20
cols = 10

a = np.asarray(np.random.random((rows, cols)))
a_gpu = gpuarray.to_gpu(a.copy())
lda = rows

u_gpu = gpuarray.empty((rows, rows), np.float32)
v_gpu = gpuarray.empty((cols, cols), np.float32)
s_gpu = gpuarray.empty(cols, np.float32)
devInfo = gpuarray.zeros(1, np.int32)

handle = cusolver.cusolverDnCreate()

worksize = cusolver.cusolverDnSgesvd_bufferSize(handle, rows, cols)
print("SIZE", worksize)

Workspace = gpuarray.empty(worksize, np.float32)

svd_status = cusolver.cusolverDnSgesvd(
    handle=handle,
    jobu='A',
    jobvt='A',
    m=rows,
    n=cols,
    A=a_gpu.ptr,
    lda=rows,
    S=s_gpu.ptr,
    U=u_gpu.ptr,
    ldu=rows,
    VT=v_gpu.ptr,
    ldvt=cols,
    Work=Workspace.ptr,
    Lwork=worksize,
    rwork=Workspace.ptr,
    devInfo=devInfo.ptr
)
status = cusolver.cusolverDnDestroy(handle)

And I'm getting a new error

Traceback (most recent call last):
  File "/home/vektor/PycharmProjects/yancut/test_svd.py", line 53, in <module>
    devInfo=devInfo.ptr
  File "/home/vektor/anaconda3/lib/python3.4/site-packages/skcuda/cusolver.py", line 579, in cusolverDnSgesvd
    Lwork, int(rwork), int(devInfo))
ctypes.ArgumentError: argument 2: <class 'TypeError'>: wrong type

It now seems that I'm doing something wrong with devInfo

*"But the code isn't working"*: what is your desired output, which are you getting instead? please edit your question — m.s., Oct 30 '15 at 10:53
@m.s. sorry, I was distracted and forgot to add further details. Basically I'm providing the - probably - right arguments but in a wrong format, the desired output are the full matrices U,S,VT of SVD. — Vektor88, Oct 30 '15 at 11:01
You haven't copied all of the error message text, it is partially cut off. Could you please edit the complete error text into the question? — talonmies, Oct 30 '15 at 12:10

score 2 · Answer 1 · answered Oct 30 '15 at 12:22

2

From the documentation it looks like each of the matrices (so A, S, U, VT) need to be passed as device pointers. So for PyCUDA gpuarrays, pass A.ptr rather than A. etc and it should work.

answered Oct 30 '15 at 12:22

talonmies

70,661
34
192
269

Thanks for the directions, I was already trying something similar and my code is close to compile, I've edited my question with the new code. – Vektor88 Oct 30 '15 at 12:34

Use scikit-cuda to compute singular value decomposition with cuSOLVER

1 Answers1