Accelerate eigs with GPU

Dear all,

I have implemented a numerical solver (of the Fokker-Planck equation) in MATLAB.

At some point, the algorithm needs to calculate

eigs(L, 1, 0)

of a very large sparse matrix L.

I would like to perform this calculation on a GPU to lower the computational costs. So, I created the array

Lgpu = gpuArray(L);

on the GPU and tried to calculate

eigs(Lgpu, 1, 0)

again.

Unfortunately, I receive the error message: "First argument must be a double matrix or a function."

I am wondering what the cause of this error might be and appreciate any help from you.

Thank you very much.

best,

Sven

Answer 1

[Walter Roberson]

There is no GPU support for sparse arrays.

[Edric Ellis]

This is no longer correct, support for sparse gpuArray was added in R2015a.

[Walter Roberson]

Ah!

Answer 2

[Christine Tobler]

The eigs function is not supported on the GPU. There is support for sparse matrices on the GPU, since R2015a: Release notes parallel computing toolbox.

You could pass a function handle to EIGS that would use computation on the GPU, but would need to accept and pass back out vectors on the CPU. I'm not sure how efficient that would be, but it could be worth a try.

[Sven Auschra]

Thanks very much for clarification.

TBH, I am not quite sure how to implement you suggestion (function handle to EIGS). Could you provide a minimal example, or provide a reference?

[Bruno Luong]

I posted it yesterday but I remove it due to Walter post.

But it's no complicated than this.

eigs( @(x) gather(Lgpu\gpuArray(x)), size(L,1), 1, 0);

or

eigs(@(x) gather(bicg(Lgpu,gpuArry(x))), size(L,1), 1, 0);

[Christine Tobler]

To compute the largest eigenvalue by absolute value, you would use

n = size(L, 1);
eigs(@(x) gather(Lgpu*x), n, 1)

which would apply the matrix-vector product on the GPU, and then move back to the CPU to pass that result back to EIGS. This could give you a first impression if there's anything to gain from using this as opposed to just calling

eigs(L, 1);

Since you're actually computing the eigenvalues of L that are closest to zero, things are a bit more complicated, and I'm still investigating if this can be done efficiently on the GPU at the moment.

[Christine Tobler]

Bruno, it looks like we were posting at the same time there. Using Lgpu\gpuArray(x) will definitely work, but since backslash is called many times here, to be efficient Lgpu should be factorized and this factorization used in the function handle. Otherwise I'd expect that the CPU version will still be faster.

Unfortunately, LU doesn't support sparse gpuArrays at this point, which means we have to call backslash which will compute this factorization underneath on every call.

Overall, it's worth trying out these ideas, but I kind of expect that while they will move some computation to the GPU, they won't actually make your EIGS call faster than the CPU version, because gpuArray doesn't yet support all tools to get the best possible performance out of EIGS.

[Sven Auschra]

Thank you so much for your support.

Unfortunately, calculations on the GPU (Kepler K20) are by a factor of 17 slower than on the CPU. So, I guess I will stick to calculating EIGS on the latter one.

Anyway, thank you so much. I really appreciate it.

[Bruno Luong]

Yeah I would also expect any speedup using GPU. Mostly EIGS on sparse is mainly an iterative process in double layers. Nothing really leans for GPU computation.