《算法导论》学习（四）---- 矩阵乘法的Strassen（斯特拉森）算法

已于 2022-08-31 12:28:56 修改
阅读量4.2k 收藏 55
点赞数 6
分类专栏： 数据结构与算法 文章标签： 算法 矩阵 学习
于 2022-08-29 12:37:23 首次发布
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/weixin_52042488/article/details/126345027
版权

文章目录

前言

矩阵乘法可以采用分治的策略。
这里提供了两个分治策略的解决 n ∗ n n*n n∗n矩阵之间乘法的算法

1.矩阵乘法的普通递归方法
2.矩阵乘法的Strassen（斯特拉森）方法

但是着两个方法的缺点是只能是两个 n ∗ n n*n n∗n矩阵的乘法，同时n必须为2的幂

之后也对这两个算法进行了时间复杂度上的分析

一、矩阵乘法的普通递归方法

1.C语言代码实现

#include<stdio.h>
#include<stdlib.h>
#include<time.h>


//宏定义改变矩阵的大小
//此处size指的是n*n矩阵的宽度或者长度n
//由于Strassen算法本身的局限性
//两个相乘的矩阵只能是n*n,且n只能是2的幂
//即size为1,2,4,8,16...
#define size 4


//矩阵的合并
//就是将A11,A12,A21,A22合并为
//	|A11	A12| 
// 	|A21	A22|
//这样的形式 
void Merge_Matrix(int *a,int *b,int *c,int *d,int* c0,int rows)
{
	//rows是子矩阵的宽度，那么合并后矩阵的宽度就是2*rows 
	int i=0;//子矩阵遍历索引 
	int j=0;//合并矩阵遍历索引，此处先合并A11 
	for(i=0;i<(rows*rows);i++)
	{
		//如果执行了rows次就需要换行，要锁定到合并后矩阵的第二行，所以加上2*rows即可 
		if((i%rows==0)&&i!=0)
		{
			j=(rows*2)*(i/rows);//(i/rows)代表需要换多少行；(rows*2)就是行数 
			c0[j]=a[i];
			j++;
		}
		else
		{
			c0[j]=a[i];
			j++;
		}
	}
	//此处是A12的首元素，令索引等于子矩阵的那个函数
	j=rows;
	for(i=0;i<(rows*rows);i++)
	{
		//这里的条件和之前的一样的，其它四部分也是一样的，底层逻辑是一致的 
		if((i%rows==0)&&i!=0)
		{
			j=(rows*2)*(i/rows)+rows;//换行要换到第二部分，所以需要再加一个rows 
			c0[j]=b[i];
			j++;
		}
		else
		{
			c0[j]=b[i];
			j++;
		}
	}
	//此处是A21的首元素，令索引等于子矩阵的那个函数 
	j=rows*2*rows;
	for(i=0;i<(rows*rows);i++)
	{
		if((i%rows==0)&&i!=0)
		{
			j=rows*2*rows+(rows*2)*(i/rows);
			c0[j]=c[i];
			j++;
		}
		else
		{
			c0[j]=c[i];
			j++;
		}
	}
	//此处是A22的首元素，令索引等于子矩阵的那个函数
	j=rows*2*rows+rows;
	for(i=0;i<(rows*rows);i++)
	{
		if((i%rows==0)&&i!=0)
		{
			j=rows*2*rows+rows+(rows*2)*(i/rows);
			c0[j]=d[i];
			j++;
		}
		else
		{
			c0[j]=d[i];
			j++;
		}
	}
}


//两个大小相同的矩阵的加法 
void Matrix_Add(int* x,int* y,int* c0,int rows)
{
	int i=0;
	for(i=0;i<(rows*rows);i++)
	{
		c0[i]=x[i]+y[i];
	}
}


//矩阵的乘法递归函数 
void Matrix_Multiply(int* x,int* y,int* c0,int loca,int locb,int rows)
{
	//进行索引的申明
	//递归过程为了节约时间，没有进行原矩阵的切割赋值
	//直接采用索引的方法进行矩阵的伪切割 
	int loca11,loca12,loca21,loca22;
	int locb11,locb12,locb21,locb22;
	//结果矩阵的四分子矩阵申明 
	int c11[(rows/2)*(rows/2)];
	int c12[(rows/2)*(rows/2)];
	int c21[(rows/2)*(rows/2)];
	int c22[(rows/2)*(rows/2)];
	//中间暂存矩阵申明 
	int temp1[(rows/2)*(rows/2)];
	int temp2[(rows/2)*(rows/2)];
	//四分矩阵后的子矩阵宽度 
	int newrows=rows/2;
	//如果矩阵的宽度为1，那么就直接元素相乘，递归结束 
	if(rows==1)
	{
		c0[0]=x[(loca)]*y[(locb)];
	}
	else
	{
		//索引的计算，进行矩阵伪分割 
		//因为n为2的幂，所以直接二分即可 
		loca11=loca;
		locb11=locb;
		loca12=loca+(rows/2);
		locb12=locb+(rows/2);
		loca21=loca+(size*(rows/2));
		locb21=locb+(size*(rows/2));
		loca22=loca+(size*(rows/2))+(rows/2);
		locb22=locb+(size*(rows/2))+(rows/2);
		//根据算法公式进行运算
		//涉及矩阵相乘的进入递归 
		Matrix_Multiply(x,y,temp1,loca11,locb11,newrows);
		Matrix_Multiply(x,y,temp2,loca12,locb21,newrows);
		Matrix_Add(temp1,temp2,c11,newrows);
		Matrix_Multiply(x,y,temp1,loca11,locb12,newrows);
		Matrix_Multiply(x,y,temp2,loca12,locb22,newrows);
		Matrix_Add(temp1,temp2,c12,newrows);
		Matrix_Multiply(x,y,temp1,loca21,locb11,newrows);
		Matrix_Multiply(x,y,temp2,loca22,locb21,newrows);
		Matrix_Add(temp1,temp2,c21,newrows);		
		Matrix_Multiply(x,y,temp1,loca21,locb12,newrows);
		Matrix_Multiply(x,y,temp2,loca22,locb22,newrows);
		Matrix_Add(temp1,temp2,c22,newrows);
		//递归后，进行矩阵的合并 
		Merge_Matrix(c11,c12,c21,c22,c0,newrows);
	}
}


int main()
{
	int i=0;
	//动态生成指定大小的矩阵 
	int* a=NULL;
	int* b=NULL;
	int* c=NULL;
	a=(int *)malloc(size*size*sizeof(int));
	b=(int *)malloc(size*size*sizeof(int));
	c=(int *)malloc(size*size*sizeof(int));
	//时间随机数种子生成 
	srand((unsigned)time(NULL));
	//给A和B矩阵赋随机值 （0~10） 
	for(i=0;i<(size*size);i++)
	{
		a[i]=rand()%10;
		b[i]=rand()%10;
	}
	//打印矩阵A和B 
	printf("A matrix is:\n\n");
	for(i=0;i<(size*size);i++)
	{
		if(((i+1)%size)==0)
		{
			printf("%5d\n",a[i]);
		}
		else
		{
			printf("%5d",a[i]);
		}
	}
	printf("B matrix is:\n\n");
	for(i=0;i<(size*size);i++)
	{
		if(((i+1)%size)==0)
		{
			printf("%5d\n",b[i]);
		}
		else
		{
			printf("%5d",b[i]);
		}
	}
	//进行矩阵乘法 
	Matrix_Multiply(a,b,c,0,0,size);
	//打印结果矩阵C 
	printf("C matrix is:\n\n");
	for(i=0;i<(size*size);i++)
	{
		if(((i+1)%size)==0)
		{
			printf("%5d\n",c[i]);
		}
		else
		{
			printf("%5d",c[i]);
		}
	}
	//释放动态内存 
	free(a);
	free(b);
	free(c);
	return 0;
}

2.算法原理分析

对于 n ∗ n n*n n∗n矩阵A，B，C（n为2的幂）
我们可以进行如下矩阵分割：
A = [ A 11 A 12 A 21 A 22 ] , B = [ B 11 B 12 B 21 B 22 ] , C = [ C 11 C 12 C 21 C 22 ] A= \begin{bmatrix} A_{11}& A_{12}\\ A_{21}& A_{22}\\ \end{bmatrix} , B= \begin{bmatrix} B_{11}& B_{12}\\ B_{21}& B_{22}\\ \end{bmatrix} , C= \begin{bmatrix} C_{11}& C_{12}\\ C_{21}& C_{22}\\ \end{bmatrix} A=[A11​A21​​A12​A22​​],B=[B11​B21​​B12​B22​​],C=[C11​C21​​C12​C22​​]
对于矩阵运算式 C = A ⋅ B C=A\cdot B C=A⋅B，我们可以改写为：
[ C 11 C 12 C 21 C 22 ] = [ A 11 A 12 A 21 A 22 ] ⋅ [ B 11 B 12 B 21 B 22 ] \begin{bmatrix} C_{11}& C_{12}\\ C_{21}& C_{22}\\ \end{bmatrix}= \begin{bmatrix} A_{11}& A_{12}\\ A_{21}& A_{22}\\ \end{bmatrix} \cdot \begin{bmatrix} B_{11}& B_{12}\\ B_{21}& B_{22}\\ \end{bmatrix} [C11​C21​​C12​C22​​]=[A11​A21​​A12​A22​​]⋅[B11​B21​​B12​B22​​]
可以进一步写成：
C 11 = A 11 ⋅ B 11 + A 12 ⋅ B 21 C_{11}=A_{11}\cdot B_{11}+A_{12}\cdot B_{21} C11​=A11​⋅B11​+A12​⋅B21​
C 12 = A 11 ⋅ B 12 + A 12 ⋅ B 22 C_{12}=A_{11}\cdot B_{12}+A_{12}\cdot B_{22} C12​=A11​⋅B12​+A12​⋅B22​
C 21 = A 21 ⋅ B 11 + A 22 ⋅ B 21 C_{21}=A_{21}\cdot B_{11}+A_{22}\cdot B_{21} C21​=A21​⋅B11​+A22​⋅B21​
C 22 = A 21 ⋅ B 12 + A 22 ⋅ B 22 C_{22}=A_{21}\cdot B_{12}+A_{22}\cdot B_{22} C22​=A21​⋅B12​+A22​⋅B22​
这样我们就可以实现两个 n ∗ n n*n n∗n矩阵乘法的分割，分次运算。
再结合n为2的幂的条件，我们就可以进行算法的递归实现。

3.编程细节

（1）用索引的方式进行伪切割

根据算法，程序需要对原矩阵进行切割，将其分为四个大小一样的子矩阵，然后一直这样细分下去，直到递归结束。那么这个分割的操作就是很关键的。有两种可行的分割方式：

1.将矩阵实际分割，并且进行赋值，得到四个子矩阵
2.将矩阵进行伪分割，并没有去得到四个子矩阵，而是通过索引来让程序辨别子矩阵的位置，实际上一直对原矩阵操作

在进行索引的赋值时，关键就是要明确：该方法本质上是对原矩阵进行操作，那么换行时是原矩阵的宽度：即size

（2）编写递归结构

该方法的递归是针对矩阵乘法，那么我们只需要在乘法函数中明确递归结束条件，并且在函数实现时，将涉及矩阵乘法的部分，用调用函数本身的方法完成，即可完成递归结构了。

二、矩阵乘法的Strassen（斯特拉森）方法

1.C语言代码实现

#include<stdio.h>
#include<stdlib.h>
#include<time.h>


//宏定义改变矩阵的大小
//此处size指的是n*n矩阵的宽度或者长度n
//由于Strassen算法本身的局限性
//两个相乘的矩阵只能是n*n,且n只能是2的幂
//即size为1,2,4,8,16...
#define size 32


//矩阵的合并
//就是将A11,A12,A21,A22合并为
//	|A11	A12| 
// 	|A21	A22|
//这样的形式 
void Merge_Matrix(int *a,int *b,int *c,int *d,int* c0,int rows)
{
	//rows是子矩阵的宽度，那么合并后矩阵的宽度就是2*rows 
	int i=0;//子矩阵遍历索引 
	int j=0;//合并矩阵遍历索引，此处先合并A11 
	for(i=0;i<(rows*rows);i++)
	{
		//如果执行了rows次就需要换行，要锁定到合并后矩阵的第二行，所以加上2*rows即可 
		if((i%rows==0)&&i!=0)
		{
			j=(rows*2)*(i/rows);//(i/rows)代表需要换多少行；(rows*2)就是行数 
			c0[j]=a[i];
			j++;
		}
		else
		{
			c0[j]=a[i];
			j++;
		}
	}
	//此处是A12的首元素，令索引等于子矩阵的那个函数
	j=rows;
	for(i=0;i<(rows*rows);i++)
	{
		//这里的条件和之前的一样的，其它四部分也是一样的，底层逻辑是一致的 
		if((i%rows==0)&&i!=0)
		{
			j=(rows*2)*(i/rows)+rows;//换行要换到第二部分，所以需要再加一个rows 
			c0[j]=b[i];
			j++;
		}
		else
		{
			c0[j]=b[i];
			j++;
		}
	}
	//此处是A21的首元素，令索引等于子矩阵的那个函数 
	j=rows*2*rows;
	for(i=0;i<(rows*rows);i++)
	{
		if((i%rows==0)&&i!=0)
		{
			j=rows*2*rows+(rows*2)*(i/rows);
			c0[j]=c[i];
			j++;
		}
		else
		{
			c0[j]=c[i];
			j++;
		}
	}
	//此处是A22的首元素，令索引等于子矩阵的那个函数
	j=rows*2*rows+rows;
	for(i=0;i<(rows*rows);i++)
	{
		if((i%rows==0)&&i!=0)
		{
			j=rows*2*rows+rows+(rows*2)*(i/rows);
			c0[j]=d[i];
			j++;
		}
		else
		{
			c0[j]=d[i];
			j++;
		}
	}
}


//两个n*n矩阵的减法，是x-y 
void Matrix_SUB(int* x,int* y,int* c0,int rows)
{
	int i=0;
	for(i=0;i<(rows*rows);i++)
	{
		c0[i]=x[i]-y[i];
	}
}


//两个n*n矩阵的加法，是x+y 
void Matrix_ADD(int* x,int* y,int* c0,int rows)
{
	int i=0;
	for(i=0;i<(rows*rows);i++)
	{
		c0[i]=x[i]+y[i];
	}
}


//将一个矩阵 
//	|A11	A12|
//	|A21	A22|
//分解为四个矩阵
//A11,A12,A21,A22
//这个函数其实是Matrix_merge的一个逆运算
//逻辑是一致的，只需要交换赋值位置即可 
void Matrix_division(int* a,int* b,int* c,int* d,int* c0,int rows)
{
	int i=0;
	int j=0;
	for(i=0;i<(rows*rows);i++)
	{
		if((i%rows==0)&&i!=0)
		{
			j=(rows*2)*(i/rows);
			a[i]=c0[j];
			j++;
		}
		else
		{
			a[i]=c0[j];
			j++;
		}
	}
	j=rows;
	for(i=0;i<(rows*rows);i++)
	{
		if((i%rows==0)&&i!=0)
		{
			j=(rows*2)*(i/rows)+rows;
			b[i]=c0[j];
			j++;
		}
		else
		{
			b[i]=c0[j];
			j++;
		}
	}
	j=rows*2*rows;
	for(i=0;i<(rows*rows);i++)
	{
		if((i%rows==0)&&i!=0)
		{
			j=rows*2*rows+(rows*2)*(i/rows);
			c[i]=c0[j];
			j++;
		}
		else
		{
			c[i]=c0[j];
			j++;
		}
	}
	j=rows*2*rows+rows;
	for(i=0;i<(rows*rows);i++)
	{
		if((i%rows==0)&&i!=0)
		{
			j=rows*2*rows+rows+(rows*2)*(i/rows);
			d[i]=c0[j];
			j++;
		}
		else
		{
			d[i]=c0[j];
			j++;
		}
	}	
}


//矩阵乘法的函数 
void Matrix_Multiply(int* x,int* y,int* c0,int rows)
{
	//定义八个子矩阵，两个乘数矩阵各四个 
	int x11[(rows/2)*(rows/2)];
	int x12[(rows/2)*(rows/2)];
	int x21[(rows/2)*(rows/2)];
	int x22[(rows/2)*(rows/2)];
	int y11[(rows/2)*(rows/2)];
	int y12[(rows/2)*(rows/2)];
	int y21[(rows/2)*(rows/2)];
	int y22[(rows/2)*(rows/2)];
	
	//定义四个子矩阵，结果矩阵的子矩阵 
	int c11[(rows/2)*(rows/2)];
	int c12[(rows/2)*(rows/2)];
	int c21[(rows/2)*(rows/2)];
	int c22[(rows/2)*(rows/2)];
	
	//Strassen方法需要的中间矩阵 
	int s1[(rows/2)*(rows/2)];
	int s2[(rows/2)*(rows/2)];
	int s3[(rows/2)*(rows/2)];
	int s4[(rows/2)*(rows/2)];
	int s5[(rows/2)*(rows/2)];
	int s6[(rows/2)*(rows/2)];
	int s7[(rows/2)*(rows/2)];
	int s8[(rows/2)*(rows/2)];
	int s9[(rows/2)*(rows/2)];
	int s10[(rows/2)*(rows/2)];
	int p1[(rows/2)*(rows/2)];
	int p2[(rows/2)*(rows/2)];
	int p3[(rows/2)*(rows/2)];
	int p4[(rows/2)*(rows/2)];
	int p5[(rows/2)*(rows/2)];
	int p6[(rows/2)*(rows/2)];
	int p7[(rows/2)*(rows/2)];
	
	//代码实现需要的暂存矩阵 
	int temp1[(rows/2)*(rows/2)];
	int temp2[(rows/2)*(rows/2)];

	//rows==1说明矩阵只有一个元素，是递归的结束，递归树的终点 
	if(rows==1)
	{
		c0[0]=x[0]*y[0];
	}
	else
	{
		//得到新的大小，进入下次递归 
		int newrows=rows/2;
		
		//将矩阵分割 
		Matrix_division(y11,y12,y21,y22,y,newrows);
		Matrix_division(x11,x12,x21,x22,x,newrows);
		
		//按照Strassen方法，进行矩阵预处理		
		Matrix_SUB(y12,y22,s1,newrows);
		Matrix_ADD(x11,x12,s2,newrows);
		Matrix_ADD(x21,x22,s3,newrows);
		Matrix_SUB(y21,y11,s4,newrows);
		Matrix_ADD(x11,x22,s5,newrows);
		Matrix_ADD(y11,y22,s6,newrows);
		Matrix_SUB(x12,x22,s7,newrows);
		Matrix_ADD(y21,y22,s8,newrows);
		Matrix_SUB(x11,x21,s9,newrows);
		Matrix_ADD(y11,y12,s10,newrows);

		//按照Strassen方法，进行矩阵乘法上的递归 
		Matrix_Multiply(x11,s1,p1,newrows);
		Matrix_Multiply(s2,y22,p2,newrows);
		Matrix_Multiply(s3,y11,p3,newrows);
		Matrix_Multiply(x22,s4,p4,newrows);
		Matrix_Multiply(s5,s6,p5,newrows);
		Matrix_Multiply(s7,s8,p6,newrows);		
		Matrix_Multiply(s9,s10,p7,newrows);
		
		//按照Strassen方法，进行结果矩阵的计算
		//得到c11,c12,c21,c22 
		Matrix_ADD(p5,p6,temp1,newrows);
		Matrix_SUB(p4,p2,temp2,newrows);
		Matrix_ADD(temp1,temp2,c11,newrows);
		
		Matrix_ADD(p1,p2,c12,newrows);
		
		Matrix_ADD(p4,p3,c21,newrows);
				
		Matrix_SUB(p5,p3,temp1,newrows);
		Matrix_SUB(p1,p7,temp2,newrows);
		Matrix_ADD(temp1,temp2,c22,newrows);
		
		//最后将c11,c12,c21,c22合并为C 
		Merge_Matrix(c11,c12,c21,c22,c0,newrows);
	}
}


int main()
{
	int i=0;
	int* a=NULL;
	int* b=NULL;
	int* c=NULL;
	//动态生成三个n*n矩阵 
	a=(int *)malloc(size*size*sizeof(int));
	b=(int *)malloc(size*size*sizeof(int));
	c=(int *)malloc(size*size*sizeof(int));
	//随机数时间种子生成 
	srand((unsigned)time(NULL));
	//生成-10~10的随机数 
	for(i=0;i<(size*size);i++)
	{
		a[i]=(rand()%10)-(rand()%10);
		b[i]=rand()%10-(rand()%10);
	}
	//打印乘数矩阵A和B 
	printf("A matrix is:\n\n");
	for(i=0;i<(size*size);i++)
	{
		if(((i+1)%size)==0)
		{
			printf("%6d\n",a[i]);
		}
		else
		{
			printf("%6d",a[i]);
		}
	}
	printf("B matrix is:\n\n");
	for(i=0;i<(size*size);i++)
	{
		if(((i+1)%size)==0)
		{
			printf("%6d\n",b[i]);
		}
		else
		{
			printf("%6d",b[i]);
		}
	}
	//进行矩阵乘法 
	Matrix_Multiply(a,b,c,size);
	//打印结果矩阵 
	printf("C matrix is:\n\n");
	for(i=0;i<(size*size);i++)
	{
		if(((i+1)%size)==0)
		{
			printf("%6d\n",c[i]);
		}
		else
		{
			printf("%6d",c[i]);
		}
	}
	//释放动态生成的矩阵 
	free(a);
	free(b);
	free(c);
	return 0;
}

2.算法原理分析

对于 n ∗ n n*n n∗n矩阵A，B，C（n为2的幂）
我们还是可以进行如下矩阵分割：
A = [ A 11 A 12 A 21 A 22 ] , B = [ B 11 B 12 B 21 B 22 ] , C = [ C 11 C 12 C 21 C 22 ] A= \begin{bmatrix} A_{11}& A_{12}\\ A_{21}& A_{22}\\ \end{bmatrix} , B= \begin{bmatrix} B_{11}& B_{12}\\ B_{21}& B_{22}\\ \end{bmatrix} , C= \begin{bmatrix} C_{11}& C_{12}\\ C_{21}& C_{22}\\ \end{bmatrix} A=[A11​A21​​A12​A22​​],B=[B11​B21​​B12​B22​​],C=[C11​C21​​C12​C22​​]
分割之后我们先进行如下一些预处理运算
S 1 = B 12 − B 22 S_1=B_{12}-B_{22} S1​=B12​−B22​
S 2 = A 11 + A 12 S_2=A_{11}+A_{12} S2​=A11​+A12​
S 3 = A 21 + A 22 S_3=A_{21}+A_{22} S3​=A21​+A22​
S 4 = B 21 − B 11 S_4=B_{21}-B_{11} S4​=B21​−B11​
S 5 = A 11 + A 22 S_5=A_{11}+A_{22} S5​=A11​+A22​
S 6 = B 11 + B 22 S_6=B_{11}+B_{22} S6​=B11​+B22​
S 7 = A 12 − A 22 S_7=A_{12}-A_{22} S7​=A12​−A22​
S 8 = B 21 + B 22 S_8=B_{21}+B_{22} S8​=B21​+B22​
S 9 = A 11 − A 21 S_9=A_{11}-A_{21} S9​=A11​−A21​
S 10 = B 11 + B 12 S_{10}=B_{11}+B_{12} S10​=B11​+B12​
然后我们根据S1~S2的值，进行7次乘法运算，如下：
P 1 = A 11 ⋅ S 1 P_{1}=A_{11}\cdot S_1 P1​=A11​⋅S1​
P 2 = S 2 ⋅ B 22 P_{2}=S_{2}\cdot B_{22} P2​=S2​⋅B22​
P 3 = S 3 ⋅ B 11 P_{3}=S_3\cdot B_{11} P3​=S3​⋅B11​
P 4 = A 22 ⋅ S 4 P_{4}=A_{22}\cdot S_4 P4​=A22​⋅S4​
P 5 = S 5 ⋅ S 6 P_{5}=S_5\cdot S_6 P5​=S5​⋅S6​
P 6 = S 7 ⋅ S 8 P_{6}=S_7\cdot S_8 P6​=S7​⋅S8​
P 7 = S 9 ⋅ S 10 P_{7}=S_9\cdot S_{10} P7​=S9​⋅S10​
然后再通过简单的加法运算，就可以得到结果矩阵的四个子矩阵：
C 11 = P 5 + P 4 − P 2 + P 6 C_{11}=P_5+P_4-P_2+P_6 C11​=P5​+P4​−P2​+P6​
C 12 = P 1 + P 2 C_{12}=P_1+P_2 C12​=P1​+P2​
C 21 = P 3 + P 4 C_{21}=P_3+P_4 C21​=P3​+P4​
C 22 = P 5 + P 1 − P 3 − P 7 C_{22}=P_5+P_1-P_3-P_7 C22​=P5​+P1​−P3​−P7​
最后将子矩阵合并即可

3.编程细节

（1）分割矩阵

我们在普通的递归方法中，对于分割矩阵采用了利用索引的方法。但是在Strassen方法中，涉及了大量的中间矩阵，造成在递归调用的时候，原矩阵和中间矩阵会一起进入递归，如果用统一的索引方法，就会造成数据空间上的越界，故只能采用传统的对子矩阵赋值的方法。

当然我们可以对合并矩阵的代码进行逆向工程，稍微修改就可以得到正确的分割矩阵的代码。

三、算法的时间复杂度分析

1.两个方法的时间复杂度

1.普通递归方法的时间复杂度为： Ω ( n 3 ) \Omega(n^3) Ω(n3)

2.Strassen方法的时间复杂度为： O ( n 2.18 ) O(n^{2.18}) O(n2.18)

具体的计算过程将在下一篇文章中讲解。

链接: 《算法导论》学习（五）---- 分治策略的时间复杂度求解

2.两个方法时间上的比较

1.从渐进复杂性上，Strassen方法是要优于普通的递归方法的，归根结底是因为Strassen的递归数更小，因为只需要7次乘法，而普通递归方法需要8次乘法。

2.但是由于Strassen方法的中间过程更加繁琐，所以当输入规模小的时候，普通的递归方法反而要快不少。

SigmaBull
关注
专栏目录
[算法系列之十五]Strassen矩阵相乘算法
SmartSi
02-06 1120
引言 Strassen矩阵乘法是一种典型的分治算法。目前为止，我们已经见过一些分治策略的算法了，例如归并排序和Karatsuba大数快速乘法。现在，让我们看看分治策略的背后原理是什么。 同动态规划不同，在动态规划中，为了得到最终的答案，我们需要把一个大的问题“展开”为几个子问题（“expand” the solutions of sub-problems），但是在这里，我们会更多的谈到如何把...
strassen矩阵乘法
chenglinhust的专栏
07-02 5205
Strassen矩阵乘法简要解析 Strassen矩阵乘法具体描述如下： 两个n×n 阶的矩阵A与B的乘积是另一个n×n 阶矩阵C，C可表示为假如每一个C(i, j) 都用此公式计算，则计算C所需要的操作次数为n3 m+n2 (n- 1) a，其中m表示一次乘法，a 表示一次加法或减法。 为了使讨论简便，假设n 是2的幂（也就是说， n是1，2，4，8，1 6，...）。 首先，假
Strassen矩阵乘法
Ben的专栏
10-20 2061
Strassen矩阵乘法矩阵乘法是线性代数中最常见的运算之一，它在数值计算中有广泛的应用。若A和B是2个n×n的矩阵，则它们的乘积C=AB同样是一个n×n的矩阵。A和B的乘积矩阵C中的元素C[i,j]定义为:　若依此定义来计算A和B的乘积矩阵C，则每计算C的一个元素C[i,j]，需要做n个乘法和n-1次加法。因此，求出矩阵C的n2个元素所需的计算时间为0(n3)。60年代末
算法导论-矩阵乘法-strassen算法
weixin_30633405的博客
09-09 445
目录 1、矩阵相乘的朴素算法 2、矩阵相乘的strassen算法 3、完整测试代码c++ 4、性能分析 5、参考资料 内容 ...
【算法导论】矩阵乘法
青峰碧陋室
02-11 2142
离过年都不到十天了，还要等到这周五才能回家，想想也一年没回家了。从寒假开始到现在，已经有二十来天，这期间把2014年总结中的寒假计划也大多数完成了：The Element Of Style的阅读，三门数学课《随机过程》、《工程优化》、《数值分析》的算法实现。回家过年期间肯定不会写博客了，今天一看，这个月只写了三篇，于是乎今天必须再写一篇来完成这个月的基本工作量。言归正传，这篇文章写写选修课《算法设
算法导论（三）--分治法
so_vegetable的博客
07-17 718
算法导论（三）--分治法二分查找乘方问题斐波那契数列矩阵乘法VLSI(Very Large Scale Integration) Layout（超大规模集成电路布局问题） 分治法就是把一个大的问题分成若干个小的子问题，然后递归地解决小问题，最后再合并子问题的解。回顾一下归并排序，就是一个典型的分治法的例子：第一步把整个数组一分为二，第二步递归地对每个子数组排序，第三步合并两个有序数组。归并排序的复...
算法导论 英文版第三版 - Introduction to Algorithm
08-11
- 第二部分：分治法（Divide-and-Conquer），包括最大子数组问题（Maximum-Subarray Problem）、递归式解决方法（Solving recurrences）、斯特拉森算法（Strassen’s algorithm for matrix multiplication）等。...
算法导论第三版.英文版.高清非扫描版
01-11
书中详细讨论了递归的概念，包括递归式的设计与解决方法，特别是递归树方法和斯特拉森算法（Strassen’s algorithm）用于矩阵乘法的示例。另外，书中还讲述了主定理（Master Theorem）的证明和应用，这是求解递归式...
矩阵乘法的Strassen算法
03-16
很简洁的代码，对于初学矩阵算法的人很有帮助
矩阵乘法的strassen算法(C++实现)
04-23
一般情况下矩阵乘法需要三个for循环，时间复杂度为O(n^3)，现在我们将矩阵分块如图：( 来自MIT算法导论 ) 一般算法需要八次乘法 r = a * e + b * g ; s = a * f + b * h ; t = c * e + d * g; u = c * f + d * h; strassen将其变成7次乘法，因为大家都知道乘法比加减法消耗更多，所有时间复杂更高！ strassen的处理是： 令： p1 = a * ( f - h ) p2 = ( a + b ) * h p3 = ( c +d ) * e p4 = d * ( g - e ) p5 = ( a + d ) * ( e + h ) p6 = ( b - d ) * ( g + h ) p7 = ( a - c ) * ( e + f ) 那么我们可以知道： r = p5 + p4 + p6 - p2 s = p1 + p2 t = p3 + p4 u = p5 + p1 - p3 - p7
算法导论 英文第3版 清晰
04-21
其中，标准符号和常见函数（Standard Notations and Common Functions）以及斯特拉森算法（Strassen's Algorithm）被用于矩阵乘法的案例研究中，展示了分治法（Divide-and-Conquer）的算法设计策略。 在概率分析和...
斯特拉森矩阵乘法
aa663325的博客
10-04 2240
先上代码package Sort; public class Matrix { public static void main(String[] args) { int[][] cc={{1,1,1,1,1,1,1,1},{1,1,1,1,1,1,1,1},{1,1,1,1,1,1,1,1},{1,1,1,1,1,1,1,1},{1,1,1,1,1,1,1,1},{1,1,1,1,1,1,
矩阵乘法Strassen算法
死啦死啦滴的博客
06-06 889
这里主要是给出实现方法，至于算法的介绍，可以参考@金戈大王的介绍。下面算法是有bug的，虽然已经对一些非n*n的矩阵做出了处理，但是没有完善，当计算非n*n的矩阵是会出现数据越界的异常的。所以两个相乘的矩阵必须是n*n的。不啰嗦，下面直接贴出Java的实现方法：//调用入口：         public static int[][] StrassenMulti(int a [][], int b...
矩阵乘法-strassen算法
严霜九月
08-20 5290
 1.矩阵相乘的朴素算法 时间复杂度T(n)=Θ(n3)，朴素矩阵相乘算法，思想明了，编程实现简单。时间复杂度是Θ(n^3)。伪码如下 for i ← 1 to n     do for j ← 1 to n        do c[i][j] ← 0           for k ← 1 to n              do  c[i][j] ← c[i][j] + a[i][...
矩阵乘法——strassen算法
wbhua1984的博客
03-26 1206
strassen算法  是通过分治法进行快速矩阵乘法，运行时间达到O(n^lg7). 每一次进行18次加/减法，7次乘法。 分治过程如下： 而且这种运算方式是有局限的： 例如矩阵维数n必须是2的幂，用其他的数会导致分割矩阵时出现奇数。 每次分割计算生成的子矩阵消耗空间。 不像简单的方法那样有数值稳定性。 实现代码如下：
Strassen 算法
最新发布
07-28
Strassen算法是一种用于矩阵乘法的分治算法。它的基本思想是将两个大的矩阵分割成更小的子矩阵，然后通过一系列递归步骤来计算它们的乘积。相较于传统的矩阵乘法算法，Strassen算法在一定规模下具有更高的效率。 具体来说，Strassen算法将两个n × n的矩阵分割成四个n/2 × n/2的子矩阵，然后通过一系列计算来得到它们的乘积。这些计算包括7次乘法运算和18次加法或减法运算。通过递归地应用这些计算步骤，最终可以得到两个矩阵的乘积。 Strassen算法的时间复杂度为O(n^log2(7))，相较于传统矩阵乘法的时间复杂度O(n^3)有所改善。然而，由于Strassen算法在实际应用中存在常数因子较大和递归深度较大的问题，因此在一些情况下可能并不比传统算法更快。 总的来说，Strassen算法是一种用于矩阵乘法的高效分治算法，但在实际应用中需要综合考虑具体情况来选择是否使用。
评论 5
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

SigmaBull

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付：¥1
获取中
扫码支付

您的余额不足，请更换扫码支付或充值

打赏作者

实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值