24
$\begingroup$

I came across Strassen's algorithm for matrix multiplication, which has time complexity $O(n^{2.81})$, significantly better than the naive $O(n^3)$. Of course, there have been several other improvements in matrix multiplication since Strassen, but my question is specific to this algorithm.

If you see the algorithm, you'll notice that 7 matrices $M_1$ to $M_7$ have been defined as intermediate computation steps, and the final matrix product can be expressed in terms of these. I understand how to verify this claim, and arrive at the expression for the desired time complexity, but I'm unable to grasp the intuition behind this algorithm, i.e. why are the matrices $M_1$ through $M_7$ defined the way they are?

Thank you!

$\endgroup$
4
  • $\begingroup$ It's a good question. It's not always possible to grasp the intuition, sometimes stuff is just a bit "accidental"(?) However, if you go through the steps of the algorithm you will see how the algorithms takes an 8 operation algorithm down to a 7 operation algorithm by reusing previous calculations. $\endgroup$
    – Pål GD
    Sep 11, 2020 at 7:15
  • 4
    $\begingroup$ Great question! I'd love to find ana answer too. I find the algorithm utterly mysterious - it works, but how on earth could anyone have come up with that? I'm not sure whether there is any systematic way to derive the algorithm, short of inspiration. $\endgroup$
    – D.W.
    Sep 11, 2020 at 7:27
  • 2
    $\begingroup$ cstheory.stackexchange.com/questions/17040/… $\endgroup$ Sep 11, 2020 at 8:33
  • 1
    $\begingroup$ Personally, I've found Strassen's particular algorithm to be unintuitive, especially due to the large amount of asymmetry. I doubt it was "accidentally discovered", but I've yet to understand it. In a similar line of thought, which may interest you since they don't require lots of theory, there are Winograd algorithms based on 4x4 or 6x6 matrix multiplication which are asymptotically faster and I find to be more intuitive, since they are symmetric, and can even be implemented with a few nice loops. $\endgroup$ Sep 11, 2020 at 21:50

2 Answers 2

15
$\begingroup$

The real answer to this question is that if you play with it long enough, you'll hit an algorithm requiring 7 multiplications – not necessarily the same as Strassen's, but an equivalent one, in a certain sense: it is known that all such algorithms are equivalent, as shown by de Groote in his 1978 paper, On varieties of optimal algorithms for the computation of bilinear mappings. II. Optimal algorithms for 2 × 2-matrix multiplication.

There are many attempts in the literature to explain how one could come up with such an algorithm, for example:

Gideon Yuval shows how you could come up with Strassen's algorithm. The starting point is to convert matrix multiplication to the problem of computing a matrix-vector product: computing $$ \begin{pmatrix} a & b \\ c & d \end{pmatrix} \times \begin{pmatrix} e & f \\ g & h \end{pmatrix} $$ is the same as computing $$ \begin{pmatrix} a & 0 & b & 0 \\ 0 & a & 0 & b \\ c & 0 & c & 0 \\ 0 & d & 0 & d \end{pmatrix} \times \begin{pmatrix} e \\ f \\ g \\ h \end{pmatrix} $$ Suppose that we could write the matrix on the left as a sum $\ell_1 M_1 + \cdots + \ell_7 M_7$, where $\ell_i$ is a linear combination of $a,b,c,d$ and $M_i$ is a rank one matrix, say $M_i = x_i y_i^T$. The product we are after is thus $$ \sum_{i=1}^7 \ell_i M_i \begin{pmatrix} c\\d\\e\\f \end{pmatrix} = \sum_{i=1}^7 \ell_i x_i y_i^T \begin{pmatrix} c\\d\\e\\f \end{pmatrix} = \sum_{i=1}^7 \ell_i r_i x_i, $$ where $r_i$ is a linear combination of $e,f,g,h$. This shows that each entry of the product matrix is some linear combination of the products $\ell_i,r_i$.

Let us now show how one could find the decomposition. We start by cancelling the top left and bottom right entries, in a way which avoids hitting zero entries: $$ \begin{pmatrix} a & 0 & b & 0 \\ 0 & a & 0 & b \\ c & 0 & d & 0 \\ 0 & c & 0 & d \end{pmatrix} = \begin{pmatrix} a & 0 & a & 0 \\ 0 & 0 & 0 & 0 \\ a & 0 & a & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & d & 0 & d \\ 0 & 0 & 0 & 0 \\ 0 & d & 0 & d \end{pmatrix} + \begin{pmatrix} 0 & 0 & b-a & 0 \\ 0 & a-d & 0 & b-d \\ c-a & 0 & d-a & 0 \\ 0 & c-d & 0 & 0 \end{pmatrix} $$ This results in a mess, which we try to fix by "flipping" the inner square: $$ \begin{pmatrix} 0 & 0 & b-a & 0 \\ 0 & a-d & 0 & b-d \\ c-a & 0 & d-a & 0 \\ 0 & c-d & 0 & 0 \end{pmatrix} = \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & a-d & a-d & 0 \\ 0 & d-a & d-a & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & b-a & 0 \\ 0 & 0 & d-a & b-d \\ c-a & a-d & 0 & 0 \\ 0 & c-d & 0 & 0 \end{pmatrix} $$ Since $d-a = (b-a)-(b-d)$ and $a-d = (c-d)-(c-a)$, it is easy to represent the last matrix as a sum of four rank one matrices: $$ \begin{pmatrix} 0 & 0 & b-a & 0 \\ 0 & 0 & d-a & b-d \\ c-a & a-d & 0 & 0 \\ 0 & c-d & 0 & 0 \end{pmatrix} = \begin{pmatrix} 0 & 0 & b-a & 0 \\ 0 & 0 & b-a & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & d-b & b-d \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & c-d & 0 & 0 \\ 0 & c-d & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ c-a & a-c & 0 & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} $$ In total, we obtain the following representation: $$ \begin{pmatrix} a & 0 & b & 0 \\ 0 & a & 0 & b \\ c & 0 & d & 0 \\ 0 & c & 0 & d \end{pmatrix} = \begin{pmatrix} a & 0 & a & 0 \\ 0 & 0 & 0 & 0 \\ a & 0 & a & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & d & 0 & d \\ 0 & 0 & 0 & 0 \\ 0 & d & 0 & d \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & a-d & a-d & 0 \\ 0 & d-a & d-a & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & b-a & 0 \\ 0 & 0 & b-a & 0 \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & d-b & b-d \\ 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ 0 & c-d & 0 & 0 \\ 0 & c-d & 0 & 0 \end{pmatrix} + \begin{pmatrix} 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 \\ c-a & a-c & 0 & 0 \\ 0 & 0 & 0 & 0 \end{pmatrix} $$

$\endgroup$
0
5
$\begingroup$

It’s reasonably obvious that if you can calculate a 2x2 matrix product with 7 multiplications and quite a few additions, you get an asymptotically faster algorithm. You need 8 products. But for example (a+b)*(c+d) gives you the sum of four products with one multiplication.

So it might be possible to calculate many products with seven multiplications in such away that all the unwanted products cancel each other out. I don’t know how many products he tried that in the end didn’t work out.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.