Multiplying Sparse Matrices

Date: 12/01/2000 at 11:35:39
From: Inaki Garmendia
Subject: Matrix algebra and Computer programming

Hello,

I have to program a matrix product. The matrices are very big, but 
sparse (many zeros), and are stored in two arrays in what is called 
"skyline storage." Is there any reasonably fast way of doing this 
product?

Thank you very much in advance.

Inaki Garmendia

Date: 12/01/2000 at 13:38:31
From: Doctor Schwa
Subject: Re: Matrix algebra and Computer programming

Dear Inaki,

I'm not an expert on this sort of thing by any means, and I hope that 
another math doctor might find the time to give you a better answer to 
your question.

What springs to mind for me, though, is that you should start with the 
product (answer) matrix full of zeros. Search through the first 
matrix, and every time you find a nonzero entry, add the appropriate 
amounts to each location of the product matrix; that is, multiply that 
one entry by everything in the second matrix. Perhaps generating a 
lookup table first, of all the nonzero entries in the second matrix, 
would save time too.

For example (these matrices aren't very sparse; the savings from the 
algorithm I suggest get pretty big if the matrix is really sparse):

     [0 2 3]   [0 0 7]
     [4 0 0] * [0 8 0]
     [0 0 5]   [0 9 6]

So, first I generate a lookup table for the second matrix:

     (1, 3) 7
     (2, 2) 8
     (3, 2) 9
     (3, 3) 6

Then I look for nonzero entries in the first matrix. At (1,2) 2, I 
know I need to multiply that by each column, so looking for things in 
the second row of the second matrix, I find (2,2) 8, which means I add 
2*8 to the (1,2) spot of my final matrix.

At (1,3) 3, I look for things in row 3 of the second matrix, so 
(3,2) 9 means I add 3*9 to the (1,2) spot of the final matrix, and 
(3,3) 6 means I add 3*6 to the (1,3) spot of the final matrix.

At (2,1) 4, I look for things in row 1 of the second matrix, so 
(1, 3) 7 means I add 4*7 to the (2,3) spot of the final matrix.

At (3,3) 5, I look for things in row 3 of the second matrix, so 
(3,2) 9 means I add 5*9 to the (3,2) spot of the final matrix, and 
(3,3) 6 means I add 5*6 to the (3,3) spot of the final matrix.

So, in total, I have for the final matrix

     [0 2 3]   [0 0 7]   [0 43 18]
     [4 0 0] * [0 8 0] = [0  0 28]
     [0 0 5]   [0 9 6]   [0 45 30]

which, thankfully, is correct.

I think that method will be reasonably fast if you have big but very 
sparse matrices. It certainly saves you from ever wasting time 
multiplying by 0, and adding, anyway.

- Doctor Schwa, The Math Forum
  http://mathforum.org/dr.math/