Fork me on GitHub

Some applications have extremely large matrices with very few nonzero entries.
Such matrices occur naturally when solving partial differential equations
or as the adjacency matrix of large graphs. The \(n \times n\) matrix
$$A = \begin{bmatrix}
2 & -1 \\
-1 & 2 & -1 \\
& -1 & 2 & -1 \\
& & \ddots & \ddots & \ddots \\
& & & -1 & 2 \end{bmatrix}$$
is an example of a sparse matrix. Click here to see this example in the Workshop.

When stored as a dense matrix, it occupies \(n^2\) memory locations:

IN> A = numeric.diag([2,2,2,2,2,2]);
    for(i=0;i<5;++i) A[i+1][i] = A[i][i+1] = -1;
    A
OUT> [[     2,    -1,     0,     0,     0,     0],
      [    -1,     2,    -1,     0,     0,     0],
      [     0,    -1,     2,    -1,     0,     0],
      [     0,     0,    -1,     2,    -1,     0],
      [     0,     0,     0,    -1,     2,    -1],
      [     0,     0,     0,     0,    -1,     2]]

Instead of explicitly storing all the entries of the matrix (including the nonzero entries), we can use the Column Compressed Storage format:

IN> ccsA = numeric.ccsSparse(A);
OUT> [[  0,  2,  5,  8, 11, 14, 16],
      [  0,  1,  0,  1,  2,  1,  2,  3,  2,  3,  4,  3,  4,  5,  4,  5],
      [  2, -1, -1,  2, -1, -1,  2, -1, -1,  2, -1, -1,  2, -1, -1,  2]]

This is a bit hard to read, but ccsA[2] contains the nonzero entries of A, starting with column 0 (ie. entries 2,-1), followed by column 1, etc… Meanwhile, the array ccsA[1] contains the corresponding row numbers. The array ccsA[0] to the beginning of each column in ccsA[1] and ccsA[2]. The function numeric.ccsFull(ccsA) can be used put the entries of ccsA into a dense matrix form.

The main advantage of column compressed storage is that the matrix can be stored using much less than \(n^2\) storage. As a result, sparse matrices are not usually converted into dense matrices.

Given a vector \(b\) and a sparse matrix \(A\), one can solve the linear problem \(Ax=b\) for the unknown \(x\) by using the functions numeric.ccsLUPSolve() and numeric.ccsLUP(), as follows:

IN> numeric.ccsLUPSolve(numeric.ccsLUP(ccsA),[1,1,1,1,1,1])
OUT> [    3,    5,    6,    6,     5,    3]

The algorithm used is a sparse LUP decomposition \(PA=LU\), which may be better than computing the matrix inverse \(A^{-1}\) since the inverse of a sparse matrix is usually dense.

If you want to solve large sparse linear problems then you should use a sparse matrix and never convert it to a dense or full matrix.

September 19, 2012 at 1:37 am |

2 Responses to “Sparse Matrices in Column Compressed Storage format”

  1. Darcy Parker says:

    I think you made a small typo in this sentence:

    “Given a vector b and a sparse matrix A, one can solve the linear problem Ax=b for the unknown x by using the functions numeric.ccsLUPSolve() and numeric.LUP(), as follows”

    You probably meant numeric.cssLUP() and not numeric.LUP().

Leave a Reply