Strassen矩阵乘法分治与递归

最新推荐文章于 2023-04-15 16:54:54 发布

YbFree

最新推荐文章于 2023-04-15 16:54:54 发布

阅读量4.7k

点赞数 2

分类专栏：算法

算法专栏收录该内容

5 篇文章 0 订阅

订阅专栏

转自：http://blog.sina.com.cn/s/blog_7e9a88f70100zj2h.html

Strassen矩阵乘法

矩阵乘法是线性代数中最常见的运算之一，它在数值计算中有广泛的应用。若A和B是2个n×n的矩阵，则它们的乘积C=AB同样是一个n×n的矩阵。A和B的乘积矩阵C中的元素C[i,j]定义为:

分治与递归--strassen矩阵乘法

若依此定义来计算A和B的乘积矩阵C，则每计算C的一个元素C[i,j]，需要做n个乘法和n-1次加法。因此，求出矩阵C的n²个元素所需的计算时间为0(n³)。

60年代末，Strassen采用了类似于在大整数乘法中用过的分治技术，将计算2个n阶矩阵乘积所需的计算时间改进到O(n^log7)=O(n^2.18)。

首先，我们还是需要假设n是2的幂。将矩阵A，B和C中每一矩阵都分块成为4个大小相等的子矩阵，每个子矩阵都是n/2×n/2的方阵。由此可将方程C=AB重写为:

分治与递归--strassen矩阵乘法(1)

由此可得:

C₁₁=A₁₁B₁₁+A₁₂B₂₁ (2)

C₁₂=A₁₁B₁₂+A₁₂B₂₂ (3)

C₂₁=A₂₁B₁₁+A₂₂B₂₁ (4)

C₂₂=A₂₁B₁₂+A₂₂B₂₂ (5)

如果n=2，则2个2阶方阵的乘积可以直接用(2)-(3)式计算出来，共需8次乘法和4次加法。当子矩阵的阶大于2时，为求2个子矩阵的积，可以继续将子矩阵分块，直到子矩阵的阶降为2。这样，就产生了一个分治降阶的递归算法。依此算法，计算2个n阶方阵的乘积转化为计算8个n/2阶方阵的乘积和4个n/2阶方阵的加法。2个n/2×n/2矩阵的加法显然可以在c*n²/4时间内完成，这里c是一个常数。因此，上述分治法的计算时间耗费T(n)应该满足：

分治与递归--strassen矩阵乘法

这个递归方程的解仍然是T(n)=O(n³)。因此，该方法并不比用原始定义直接计算更有效。究其原因，乃是由于式(2)-(5)并没有减少矩阵的乘法次数。而矩阵乘法耗费的时间要比矩阵加减法耗费的时间多得多。要想改进矩阵乘法的计算时间复杂性，必须减少子矩阵乘法运算的次数。按照上述分治法的思想可以看出，要想减少乘法运算次数，关键在于计算2个2阶方阵的乘积时，能否用少于8次的乘法运算。Strassen提出了一种新的算法来计算2个2阶方阵的乘积。他的算法只用了7次乘法运算，但增加了加、减法的运算次数。这7次乘法是:

M₁=A₁₁(B₁₂-B₂₂)

M₂=(A₁₁+A₁₂)B₂₂

M₃=(A₂₁+A₂₂)B₁₁

M₄=A₂₂(B₂₁-B₁₁)

M₅=(A₁₁+A₂₂)(B₁₁+B₂₂)

M₆=(A₁₂-A₂₂)(B₂₁+B₂₂)

M₇=(A₁₁-A₂₁)(B₁₁+B₁₂)

做了这7次乘法后，再做若干次加、减法就可以得到:

C₁₁=M₅+M₄-M₂+M₆

C₁₂=M₁+M₂

C₂₁=M₃+M₄

C₂₂=M₅+M₁-M₃-M₇

以上计算的正确性很容易验证。例如:

C₂₂=M₅+M₁-M₃-M₇

=(A₁₁+A₂₂)(B₁₁+B₂₂)+A₁₁(B₁₂-B₂₂)-(A₂₁+A₂₂)B₁₁-(A₁₁-A₂₁)(B₁₁+B₁₂)

=A₁₁B₁₁+A₁₁B₂₂+A₂₂B₁₁+A₂₂B₂₂+A₁₁B₁₂

-A₁₁B₂₂-A₂₁B₁₁-A₂₂B₁₁-A₁₁B₁₁-A₁₁B₁₂+A₂₁B₁₁+A₂₁B₁₂

=A₂₁B₁₂+A₂₂B₂₂

由(2)式便知其正确性。

至此，我们可以得到完整的Strassen算法如下：

procedure STRASSEN(n,A,B,C);
begin
  if n=2 then MATRIX-MULTIPLY(A，B，C)
         else begin
                将矩阵A和B依(1)式分块;
                STRASSEN(n/2,A11,B12-B22,M1);
                STRASSEN(n/2,A11+A12,B22,M2);
                STRASSEN(n/2,A21+A22,B11,M3);
                STRASSEN(n/2,A22,B21-B11,M4);
                STRASSEN(n/2,A11+A22,B11+B22,M5);
                STRASSEN(n/2,A12-A22,B21+B22,M6);
                STRASSEN(n/2,A11-A21,B11+B12,M7);

                    分治与递归--strassen矩阵乘法;

               end;

end;

其中MATRIX-MULTIPLY(A，B，C)是按通常的矩阵乘法计算C=AB的子算法。

Strassen矩阵乘积分治算法中，用了7次对于n/2阶矩阵乘积的递归调用和18次n/2阶矩阵的加减运算。由此可知，该算法的所需的计算时间T(n)满足如下的递归方程:

分治与递归--strassen矩阵乘法

按照解递归方程的套用公式法，其解为T(n)=O(n^log7)≈O(n^2.81)。由此可见，Strassen矩阵乘法的计算时间复杂性比普通矩阵乘法有阶的改进。

有人曾列举了计算2个2阶矩阵乘法的36种不同方法。但所有的方法都要做7次乘法。除非能找到一种计算2阶方阵乘积的算法，使乘法的计算次数少于7次，按上述思路才有可能进一步改进矩阵乘积的计算时间的上界。但是Hopcroft和Kerr(197l)已经证明，计算2个2×2矩阵的乘积，7次乘法是必要的。因此，要想进一步改进矩阵乘法的时间复杂性，就不能再寄希望于计算2×2矩阵的乘法次数的减少。或许应当研究3×3或5×5矩阵的更好算法。在Strassen之后又有许多算法改进了矩阵乘法的计算时间复杂性。目前最好的计算时间上界是O(n^2.367)。而目前所知道的矩阵乘法的最好下界仍是它的平凡下界Ω(n²)。因此到目前为止还无法确切知道矩阵乘法的时间复杂性。关于这一研究课题还有许多工作可做。

importjava.util.*;

publicclass Strassen{

public Strassen(){

A = new int[NUMBER][NUMBER];

B = new int[NUMBER][NUMBER];

C = new int[NUMBER][NUMBER];

}

public void input(int a[][]){

Scanner scanner = new Scanner(System.in);

for(int i = 0; i < a.length; i++) {

for(int j = 0; j < a[i].length; j++) {

a[i][j] = scanner.nextInt();

}

public void output(int[][] resault){

for(int b[] : resault) {

for(int temp : b) {

System.out.print(temp + " ");

}

System.out.println();

}

public void Mul(int[][] first, int[][] second, int[][]resault){

for(int i = 0; i < 2; ++i) {

for(int j = 0; j < 2; ++j) {

resault[i][j] = 0;

for(int k = 0; k < 2; ++k) {

resault[i][j] += first[i][k] * second[k][j];

}

public void Add(int[][] first, int[][] second, int[][]resault){

for(int i = 0; i < first.length; i++) {

for(int j = 0; j < first[i].length; j++){

resault[i][j] = first[i][j] + second[i][j];

}

public void sub(int[][] first, int[][] second, int[][]resault){

for(int i = 0; i < first.length; i++) {

for(int j = 0; j < first[i].length; j++){

resault[i][j] = first[i][j] - second[i][j];

}

public void strassen(int[][] A, int[][] B, int[][] C){

//定义一些中间变量

int [][] M1=new int [NUMBER][NUMBER];

int [][] M2=new int [NUMBER][NUMBER];

int [][] M3=new int [NUMBER][NUMBER];

int [][] M4=new int [NUMBER][NUMBER];

int [][] M5=new int [NUMBER][NUMBER];

int [][] M6=new int [NUMBER][NUMBER];

int [][] M7=new int [NUMBER][NUMBER];

int [][] C11=new int [NUMBER][NUMBER];

int [][] C12=new int [NUMBER][NUMBER];

int [][] C21=new int [NUMBER][NUMBER];

int [][] C22=new int [NUMBER][NUMBER];

int [][] A11=new int [NUMBER][NUMBER];

int [][] A12=new int [NUMBER][NUMBER];

int [][] A21=new int [NUMBER][NUMBER];

int [][] A22=new int [NUMBER][NUMBER];

int [][] B11=new int [NUMBER][NUMBER];

int [][] B12=new int [NUMBER][NUMBER];

int [][] B21=new int [NUMBER][NUMBER];

int [][] B22=new int [NUMBER][NUMBER];

int [][] temp=new int [NUMBER][NUMBER];

int [][] temp1=new int [NUMBER][NUMBER];

if(A.length==2){

Mul(A, B, C);

}else{

//首先将矩阵A，B 分为4块

for(int i = 0; i < A.length/2; i++) {

for(int j = 0; j < A.length/2; j++) {

A11[i][j]=A[i][j];

A12[i][j]=A[i][j+A.length/2];

A21[i][j]=A[i+A.length/2][j];

A22[i][j]=A[i+A.length/2][j+A.length/2];

B11[i][j]=B[i][j];

B12[i][j]=B[i][j+A.length/2];

B21[i][j]=B[i+A.length/2][j];

B22[i][j]=B[i+A.length/2][j+A.length/2];

}

//计算M1

sub(B12, B22, temp);

Mul(A11, temp, M1);

//计算M2

Add(A11, A12, temp);

Mul(temp, B22, M2);

//计算M3

Add(A21, A22, temp);

Mul(temp, B11, M3);

//M4

sub(B21, B11, temp);

Mul(A22, temp, M4);

//M5

Add(A11, A22, temp1);

Add(B11, B22, temp);

Mul(temp1, temp, M5);

//M6

sub(A12, A22, temp1);

Add(B21, B22, temp);

Mul(temp1, temp, M6);

//M7

sub(A11, A21, temp1);

Add(B11, B12, temp);

Mul(temp1, temp, M7);

//计算C11

Add(M5, M4, temp1);

sub(temp1, M2, temp);

Add(temp, M6, C11);

//计算C12

Add(M1, M2, C12);

//C21

Add(M3, M4, C21);

//C22

Add(M5, M1, temp1);

sub(temp1, M3, temp);

sub(temp, M7, C22);

//结果送回C中

for(int i = 0; i < C.length/2; i++) {

for(int j = 0; j < C.length/2; j++) {

C[i][j]=C11[i][j];

C[i][j+C.length/2]=C12[i][j];

C[i+C.length/2][j]=C21[i][j];

C[i+C.length/2][j+C.length/2]=C22[i][j];

}

public static void main(String[] args){

Strassen demo=new Strassen();

System.out.println("输入矩阵A");

demo.input(A);

System.out.println("输入矩阵B");

demo.input(B);

demo.strassen(A, B, C);

demo.output(C);

}

private static int A[][];

private static int B[][];

private static int C[][];

private final static int NUMBER = 4;

}

YbFree

关注

2
点赞
踩
4

收藏

觉得还不错? 一键收藏
0
评论
Strassen矩阵乘法分治与递归

转自：http://blog.sina.com.cn/s/blog_7e9a88f70100zj2h.htmlStrassen矩阵乘法矩阵乘法是线性代数中最常见的运算之一，它在数值计算中有广泛的应用。若A和B是2个n×n的矩阵，则它们的乘积C=AB同样是一个n×n的矩阵。A和B的乘积矩阵C中的元素C[i,j]定义为:分治与递归--strassen矩阵乘法若依此定义来
复制链接

扫一扫