# Multiply Multiple Sparse Matrices

When we calculate products of multiple matrices, e.g., $ABC$, do you think it can be done in a cheaper way than as two consecutive multiplications? Note that I'm not talking about applying matrices to vectors but exactly matrix multiplication.

• I think the answer to your question depends strongly on what you will use the result $M = ABC$ for. Can you give some more information on this? – Pedro May 18 '12 at 18:46
• This question is for general purpose. Say, one can do LU, QR factorization for the resulting matrix $M$. Or, one uses $M$ to multiply so many vectors that to have $M$ is better than do $A$, $B$, $C$ multiplying vectors consecutively. – Hui Zhang May 18 '12 at 19:07
• Are all three matrices the same size? – Dan May 18 '12 at 19:40
• This question seems interesting even in the case of 1x1 matrices. I.e. multiply 3 floating point numbers fast. Or prove 2 multiplications is fastest. And then how about the case where they aren't native machine floats? – Andrew Moylan May 18 '12 at 19:49
• Even if the matrices are sparse, it depends very much on their structure. For example, if any two of those are arrowhead matrices, then the product is dense and you very likely do not want to form it. – Jed Brown May 19 '12 at 3:05

The answer completely depends upon how you've chosen to represent your matrices (e.g., unstructured dense, sparse, analytic, low-rank, hierarchically low-rank, etc.), whether or not there is a priori knowledge of properties of the product, and whether or not you are happy with forming an approximation to $M=ABC$.

Since most people seem to be thinking of the sparse case, I'll provide a different viewpoint. Suppose that it is known that $M$ is low-rank, with rank $r$, but that $A$, $B$, and $C$ are dense matrices, each of size $n \times n$, where $n \gg r$.

One could form $M=A(BC)$ with $4n^3$ flops, but a smarter method would start by generating a random matrix $\Omega \in \mathbb{F}^{n \times f(r)}$, where $f(r)$ is slightly greater than $r$, and recognize that $M \Omega$ almost certainly spans the column space of $M$. The $QR$ factorization can be computed in $O(nr^2)$ work (since $f(r)$ is $O(r)$), yielding a unitary matrix $Q$, and thus we can approximate

$$Mx \approx QQ'Mx = QQ'ABCx$$

for any vector $x$. Then one can compute $Q'M$ such that $Q(Q'M)$ is a low-rank approximation to $M$. Thus, once $M \Omega$ has been formed, only $O(n r^2)$ work is required to form $M$ with high probability.

In the case where $A$, $B$, and $C$ are just dense matrices, forming $M \Omega$ requires $O(n^2 r)$ work, which dominates the $O(n r^2)$ work required afterward, but this is still much better than the $O(n^3)$ naive method. Since your title claims that $A$, $B$, and $C$ are sparse, $M \Omega$ can be formed in $O(nr)$ work, so the overall complexity would be $O(nr^2)$.

See the review paper by Halko et al. for a more complete description of the ideas.

If $A$,$B$, and $C$ are the same size and don't have interactions between their sparsity patterns, then it should take the same amount of time to multiply them in any order.

If, however, the matrices are different sizes, then you can gain something my multiplying in a specific order.

After multiplying one of the matrices by one of the the others, you will have to store an intermediate matrix. The multiplication will be faster if you pick the order of multiplications such that the size of intermediate matrices is minimized.

Different sparsity patterns can interact, though, so watch out. In general you want to minimize the number of entries in the intermediate matrices, not the overall size.

It is likely that some algorithms for matrix multiply could benefit by considering a single 3-multiply rather than two 2-multiplies. For example one might be able to better optimize memory access if one considers all of the work at once.

However, this will be a substantially more complex and less general algorithm. There is a trade-off between writing your own fine-tuned algorithm and restricting yourself to a set of well understood and established operations.

Your question was "can we improve upon the standard solution in this specific case?" The answer to this question is almost always yes. A good followup question is "is it worth it?"