《Introduction to Algorithm》notes:Strassen's Algorithm

发表于 | 分类于

You might at first think that any matrix multiplication algorithm must take o(n^3) time, since the natural definition of matrix multiplication requires that many multiplications. You would be incorrect, however: we have a way to multiply matrices in o(n^3) time.

在没学习这章之前，我确实是这么想的，我认为矩阵相乘这个问题本身的复杂度就是o(n^3)，不存在“捷径”，但是这里的Strassen’s algorithm居然能将矩阵相乘的复杂度降到o(n^2.81)，简直匪夷所思，是如何做到的？

normal solution

在提到矩阵相乘的时候我们一般都是用这种方法，很符合我们的思维习惯，也是我们最开始学习矩阵的时候会接触的方法。 0

伪代码是：

//Algorithm:SQUARE-MATRIX-MULTIPLY(A,B)
n=A.rows
for i=1 to n
	for j=1 to n
		c[i][j]=0 //save every element
		for cnt=1 to n
			c+=a[i][cnt]*b[cnt][j]
return C

不难分析，三层循环，时间复杂度是O(n^3)

a simple Divide-and-Conquer solution

有学习过线性代数的同学都知道矩阵有所谓的“分块乘法”，我们可以将大矩阵看成2x2矩阵，这就变成了原问题的一个子问题，递归解决即可。

//Algorithm:SQUARE-MATRIX-MULTIPLY-RECURSIVE(A,B)
n=A.rows
let C be a new n x n matrix
if n==1
	C[1][1]=A[1][1]*B[1][1]
else partition A,B to C11,C12,C21,C22
	C11=SQUARE-MATRIX-MULTIPLY-RECURSIVE(A11,B11)+SQUARE-MATRIX-MULTIPLY-RECURSIVE(A12,B21)
	C12=SQUARE-MATRIX-MULTIPLY-RECURSIVE(A11,A12)+SQUARE-MATRIX-MULTIPLY-RECURSIVE(A12,B22)
	C21=SQUARE-MATRIX-MULTIPLY-RECURSIVE(A21,B11)+SQUARE-MATRIX-MULTIPLY-RECURSIVE(A22,B21)
	C22=SQUARE-MATRIX-MULTIPLY-RECURSIVE(A21,B12)+SQUARE-MATRIX-MULTIPLY-RECURSIVE(A22,B22)
return C

它的recurrences是： 0

如果用后面会学到的master method，可以知道这个算法的复杂度是Θ(n^3)。

但还有一点需要注意的是，这里的矩阵分块具体如何实现呢？分别为A,B,C创建4个子矩阵也就是一共12个矩阵吗？这样会额外消耗Θ(n^2)时间复制矩阵，我们可以使用指定下标的方式给定子矩阵，具体请读者自行研究。 那么它recurrences的f(x)部分是Θ(n^2)，就是上面复制矩阵而消耗的额外时间吗？不是，无论是使用复制矩阵的Θ(n^2)方法还是后者的Θ(1)方法，recurrences的Θ(n^2)都是存在的，原因是这一部分其实是矩阵相加所带来的 消耗，每个子矩阵是n/2 x n/2 也就是n^2/4个元素，4个子矩阵一共有n^2次相加。

Strassen’s Algorithm

Strassen’s Algorithm的具体做法是：

  1. Divide the input matrices A and B and output matrix C into n/2 * n/2 submatrices。 This step takes Θ(1) time by index calculation, just as in SQUARE-MATRIX-MULTIPLY-RECURSIVE. 将矩阵分块，这一部分和前面的一样。
  2. Create 10 matrices S1,S2…S10, each of which is n/2 * n/2 and is the sum or difference of two matrices created in step 1. We can create all 10 matrices in Θ(n^2) time.创建10个新矩阵，这10个矩阵由4个子矩阵相乘或相减 所得，根据前面的经验，这一部分应为Θ(n^2)。 0
  3. Using the submatrices created in step 1 and the 10 matrices created in step 2,recursively compute seven matrix products P1,P2…P7 Each matrix Pi is n/2 * n/2. 用上面创建的10个矩阵和原来的4个子矩阵运算得到P1-P7。 0

  4. Compute the desired submatrices C11,C12,C21,C22 of the result matrix C by adding and subtracting various combinations of the Pi matrices. We can compute all four submatrices in Θ(n^2) time.用P1-P7运算得到C11,C12,C21,C22。

    C11=P5+P4-P2+P6
    C12=P1+P2
    C21=P3+P4
    C22=P5+P1-P3-P7

经过上述4步，Strassen’s Algorithm成功完成了矩阵相乘，它的recurrences是： 0

算法的复杂度用master method得Θ(n^log7)。

我原本以为，Strassen’s Algorithm是发现了传统的分块乘法有“信息的冗余”,也就是说在传统的分块乘法里有“重复计算”，才降低了复杂度，但是实际上不是这样，Strassen是把A*B+A*C变成了A*(B+C),前者先计算A*B再计算 A*C,再相加，经过了两次相乘一次相加，而后者先相加再相乘，经过了一次相乘和一次相加。

Exercise

4.2-3

How would you modify Strassen’s algorithm to multiply n*n matrices in which n is not an exact power of 2? Show that the resulting algorithm runs in time Θ(n^log7).

如果n不是2的n次方，那么就不能把原矩阵平均的分成4个子矩阵，我的解决办法是，损失一些运行效率，先将原矩阵分成4个大小不等的尽量是2的n次方的矩阵，再对这4个矩阵用Strassen’s Algorithm。

4.2-4

What is the largest k such that if you can multiply 33 matrices using k multiplications (not assuming commutativity of multiplication), then you can multiply nn matrices in time o(nlog 7)? What would the running time of this algorithm be?

之前书上的例子是：对一个2*2矩阵可以8次相乘得到相乘结果，也就是原问题分解成了8个(n/2)规模的子问题，而对于一个3*3矩阵进行k次相乘得到结果，也就是原问题分解成了k个(n/3)规模的子问题，这个 算法的复杂度，根据master method，取决于log3(k)和2的相对大小。 1.k=9,即log3(k)=2。算法复杂度为Θ(logn*n^2),T(n)=O(n^log7) 2.k>9,即log3(k)>2。算法复杂度为Θ(n^log3(k)),当log3(k)<log7的时候，T(n)=O(n^log7) 3.k<9,即log3(k)<2。算法复杂为Θ(n^2)，T(n)=O(n^log7)

4.2-5

V. Pan has discovered a way of multiplying 6868 matrices using 132,464 multiplications, a way of multiplying 7070 matrices using 143,640 multiplications, and a way of multiplying 72*72 matrices using 155,424 multiplications. Which method yields the best asymptotic running time when used in a divide-and-conquer matrix-multiplication algorithm? How does it compare to Strassen’s algorithm?

log68(132,464)=2.7951 log70(143,640)=2.7951 log72(155,424)=2.7951 log2(7)=2.8074 These three methods have the same asymptotic running time and they have a better one than Strassen’s algorithm.

4.2-7

Show how to multiply the complex numbers a + bi and c + di using only three multiplications of real numbers. The algorithm should take a, b, c, and d as input and produce the real component ac - bd and the imaginary component ad - bc separately.

这道题是道很好的题目，考了我们Strassen’s Algorithm的本质，如何只用3次相乘得到(ac-bd)和(ad-bc)？ 也就是说我们需要像Strassen’s Algorithm一样，去构造P1,P2…组合得到ac-bd和ad-bc。 我是这么做的,先构造:

S1=a+b
S2=c+d

对S1，S2进行运算得到P1,P2,P3

P1=cS1=ac+bc
P2=aS2=ac+ad
P3=S1*S2=ac+ad+bc+bd

令mP1+nP2+kP3=ac-bd/ad-bc，即解两个三元方程(如果熟悉线性代数，你会发现这其实是求一个4*3矩阵的列空间问题)，得到m,n,k的值。 最后得到

P1+P2-P3=ac+bc+ac+ad-ac-ad-bc-bd=ac-bd
-P1+P2=-ac-bc+ac+ad=ad-bc
即(a+bi)*(c+di)=(P1+P2+P3)+(-P1+P2)i

这道题其实还有挖掘的空间，有没有一个通用的方法解决此类问题？有没有这样一种算法去构造出了这种降低了复杂度的算法解决问题？其实前面S1,S2和P1,P2,P3都是我根据“经验(运气)”构造出来的，后面得到m,n,k 才算用了比较科学的方法。