# Algorithm for Eigenvalue Problem of a Real Symmetric nxn Matrix

I have a nxn covariance matrix (so, real, symmetric, dense, nxn... I mean positive semi-definite). 'n' may be very very very big! I'd like to solve partial (3 largest) eigenvalue (+eigenvectors) problem for this matrix. Could somebody tell me what the fastest algorithm to do it? Is it power iteration method (because it have O(n^2) complexity)?

P.S. I'd like to make GPGPU implementation using OpenCL.

• How big is "very very big"? Commented Feb 2, 2014 at 15:22
• 100000x100000 for example Commented Feb 2, 2014 at 16:04
• If even reading your input data costs $O(n^2)$, it's difficult to think one can do better... Commented Jul 2, 2016 at 14:11

A 100,0000 by 100,000 symmetric dense matrix in single precision requires 20 gigabytes of memory (storing only the upper triangle) or 40 gigabytes of memory for double precision. Thus it is too large to fit within the memory of available GPU's.

In order to solve this problem using GPU acceleration you'd have to develop an algorithm that sends smaller chunks of work to the GPU. Copying parts of the matrix into the GPU in each iteration is unlikely to work well. Rather, you could try to use multiple GPU's in parallel, with each GPU storing part of the matrix throughout the algorithm. For example, you could use four of the NVIDIA Tesla K20X's to store your single precision matrix.

The appropriate algorithms to consider for this kind of problem are iterative schemes that use matrix-vector multiplications. Start by looking at the algorithms implemented in ARPACK. Since ARPACK can use BLAS routines for the matrix-vector multiplications, you should be able to simply link ARPACK with some GPU implementation of the BLAS to get a working code.

• Thank you for your answer. But my question was "is power iteration method fastest algorithm to solve partial eigenvalue problem"? I need to find only first 2-3 eigenvalue with corresponding eigenvectors Commented Feb 2, 2014 at 17:33
• ARPACK uses the implicitly restarted Lanczos method- This is a bit more sophisticated than simple power iteration but also uses matrix-vector multiplications as its basic operations. Commented Feb 2, 2014 at 17:46

There are many methods to solve your problem, and choosing the appropriate one is sometimes more of an art than anything else. A few things are definitely true, though:

• $10^5 \times 10^5$ is not unthinkably large in the world of scientific computing, there are many tools to deal with problems of this size
• GPU-accelerated eigenvalue solvers are still relatively new and typically not the preferred tool
• It is too large to fit onto a single compute node or GPU, which brings about the problem of scalability
• You should distribute your problem across multiple nodes to deal with the memory storage issue, many libraries deal with this and the gold standard is SLEPc which has several routines that you can experiment with

The crux: it's not clear what method is quickest for finding those eigenvalues, often depending on the matrix structure and conditioning, but the appropriate tools are available to find out.