Strassen矩阵乘法
矩阵乘法数量积_概述说明以及解释
矩阵乘法数量积概述说明以及解释1. 引言1.1 概述矩阵乘法数量积是线性代数中的一个重要概念,它用于计算两个矩阵之间的相乘结果。
通过对每个元素按一定规则进行乘法和求和运算,数量积可以得到一个新的矩阵。
这种操作在各个学科领域有广泛的应用,包括数学、物理和工程等。
1.2 文章结构本文将从以下几个方面对矩阵乘法数量积进行详细说明。
首先,我们将介绍矩阵乘法的基本概念,包括定义和性质。
然后,我们将解释矩阵乘法数量积的原理,并说明其实现过程。
接下来,我们将探讨矩阵乘法数量积在不同领域中的应用情况,包括数学、物理和工程等方面。
此外,本文还将介绍一些常见的算法和计算优化技巧,以提高矩阵乘法数量积的效率。
最后,在结论部分,我们会总结以上内容,并展望未来矩阵乘法数量积的发展趋势并给出相关建议。
1.3 目的本文旨在深入探讨矩阵乘法数量积的概念和原理,以及其在不同领域中的应用。
通过介绍常用的算法和计算优化技巧,我们希望读者能够了解到如何提高矩阵乘法数量积的计算效率。
同时,本文还旨在为未来研究者提供一些思考点,并展望矩阵乘法数量积在未来可能的发展方向。
2. 矩阵乘法数量积的定义与原理2.1 矩阵乘法的基本概念矩阵乘法是指将两个矩阵相乘得到一个新的矩阵的操作。
如果矩阵A是一个m 行n列的矩阵,而矩阵B是一个n行p列的矩阵,那么它们的乘积C将是一个m行p列的矩阵。
在此过程中,对应位置上两个矩阵元素的相乘并求和得到结果矩阵C中对应位置上的元素。
2.2 数量积的定义与性质数量积也被称为内积、点积或标量积。
对于两个向量a和b,它们之间的数量积表示为a∙b。
数量积满足以下性质:- 若a和b平行(夹角为0度),则a∙b = |a|*|b|- 若a和b垂直(夹角为90度),则a∙b = 0- 对任意向量c和标量k,有(kc)∙(kc) = k^2 * (c∙c)2.3 矩阵乘法数量积的原理解释矩阵乘法数量积可视作将两个向量进行投影、放缩和重新组合的过程。
矩阵相乘法则
矩阵相乘法则矩阵相乘法则是线性代数中的重要内容。
它描述了如何将两个矩阵相乘,并且提供了一些非常有用的解决问题的方法。
在本文中,我们将介绍矩阵相乘法则的各个方面。
1. 矩阵的乘法矩阵的乘法是线性代数中一个基本概念。
如果有两个矩阵$A$和$B$,它们可以相乘当且仅当第一个矩阵的列数等于第二个矩阵的行数。
如果$A$是$m×n$的矩阵,$B$是$n×p$的矩阵,那么它们的乘积为 $C=AB$,结果矩阵$C$是$m×p$的矩阵。
在矩阵$C$中,元素$c_{ij}$的值是矩阵$A$的第$i$行和矩阵$B$的第$j$列的乘积之和,即:$${\displaystyle c_{ij}=\sum_{k=1}^{n}a_{ik}b_{kj}}$$以下是矩阵乘法的一个例子:$${\displaystyle \begin{pmatrix}1 & 2 & 3\\4 & 5 & 6\end{pmatrix}\begin{pmatrix}7 & 8\\9 & 10\\11 & 12\end{pmatrix}=\begin{pmatrix}58 & 64\\139 & 154\end{pmatrix}}$$2. 矩阵相乘的性质矩阵相乘具有以下性质:(1)结合律:$(AB)C=A(BC)$(2)分配律:$A(B+C)=AB+AC$;$(A+B)C=AC+BC$(3)不满足交换律:$AB\neq BA$。
可以看到,矩阵相乘的结合律和分配律与实数的运算性质相似。
但是,矩阵相乘不满足交换律,即矩阵的乘积与乘法的顺序有关。
这是因为在矩阵相乘时,乘法的顺序会影响结果矩阵中元素的计算方式。
3. 矩阵乘法的应用矩阵相乘法则不仅仅是线性代数的基本内容,还被广泛应用于其他领域,如计算机科学、物理学、经济学、统计学等。
以下是一些矩阵相乘的应用:(1)图像处理图像可以表示为像素矩阵,矩阵相乘可以实现图像的旋转、缩放等变换。
strassen算法-现实2次幂,偶数奇数
Strassen算法
姓名: 王军袖 学号: 31609059
1
目录
1
算法思想
2
3
算法代码
运行结果
2
算法思想
分治算法
其核心思想是将大的问题分解成若干小的子问题进行 处理,从而使复杂的问题变得简单。表现在算法中就 是递归算法。
m1 m2 m3 m4 m5 m6 m7 = = = = = = = (a11 + a22) * (b11 + b22); (a21 + a22) * b11; a11 * (b12 - b22); a22 * (b21 - b11); (a11 + a12) * b22; (a21 - a11) * (b11 + b12); (a12 - a22) * (b21 + b22); 再做若干次加、减法: c11 = m1 + m4 - m5 + m7 c12 = m3 + m5 c21 = m2 + m4 c22 = m1 +m3 -m2 + m6
16
Matrix **cr=new Matrix*[n];//初始化cr,储存ar*br结果 //阶为偶数矩阵想乘 for ( i=0;i<n;i++){ Matrix MutiEven(Matrix A,Matrix B){ cr[i]=new Matrix[n]; int n=A.n,m=1; for ( j=0;j<n;j++) int i,j; cr[i][j]=Matrix(m); //将矩阵阶拆分成m*n,n为奇数 } m为2的幂次 for( i=0;i<n;i++)//用传统方法将ar和br相乘,储存在cr中 while(!(n&1)){ for( j=0;j<n;j++) n>>=1; for(int k=0;k<n;k++) m<<=1; cr[i][j]=Plus(cr[i][j],Muti2n(ar[i][k],br[k][j])); Matrix C(m*n);//声明C } for( i=0;i<n;i++)//将cr合并成一个矩阵C //将A拆分为n^2个m阶矩阵,在ar中储存。 for( j=0;j<n;j++) Matrix **ar=new Matrix*[n]; for(int x=0;x<m;x++) for ( i=0;i<n;i++){ for(int y=0;y<m;y++) ar[i]=new Matrix[n]; for ( j=0;j<n;j++) C.m[i*m+x][j*m+y]=cr[i][j].m[x][y]; ar[i][j]=Split(A,i*m,j*m,m); return C;} } Matrix **br=new Matrix*[n];//同A for (i=0;i<n;i++){ br[i]=new Matrix[n]; for ( j=0;j<n;j++) 17 br[i][j]=Split(B,i*m,j*m,m);}
矩阵相乘的快速算法
矩阵相乘的快速算法算法介绍矩阵相乘在进行3D变换的时候是经常用到的。
在应用中常用矩阵相乘的定义算法对其进行计算。
这个算法用到了大量的循环和相乘运算,这使得算法效率不高。
而矩阵相乘的计算效率很大程度上的影响了整个程序的运行速度,所以对矩阵相乘算法进行一些改进是必要的。
这里要介绍的矩阵算法称为斯特拉森方法,它是由v.斯特拉森在1969年提出的一个方法。
我们先讨论二阶矩阵的计算方法。
对于二阶矩阵a11 a12 b11 b12A = a21 a22B = b21 b22先计算下面7个量(1)x1 = (a11 + a22) * (b11 + b22);x2 = (a21 + a22) * b11;x3 = a11 * (b12 - b22);x4 = a22 * (b21 - b11);x5 = (a11 + a12) * b22;x6 = (a21 - a11) * (b11 + b12);x7 = (a12 - a22) * (b21 + b22);再设C = AB。
根据矩阵相乘的规则,C的各元素为(2)c11 = a11 * b11 + a12 * b21c12 = a11 * b12 + a12 * b22c21 = a21 * b11 + a22 * b21c22 = a21 * b12 + a22 * b22比较(1)(2),C的各元素可以表示为(3)c11 = x1 + x4 - x5 + x7c12 = x3 + x5c21 = x2 + x4c22 = x1 + x3 - x2 + x6根据以上的方法,我们就可以计算4阶矩阵了,先将4阶矩阵A和B划分成四块2阶矩阵,分别利用公式计算它们的乘积,再使用(1)(3)来计算出最后结果。
ma11 ma12 mb11 mb12A4 = ma21 ma22 B4 = mb21 mb22其中a11 a12 a13 a14 b11 b12 b13 b14ma11 = a21 a22 ma12 = a23 a24 mb11 = b21 b22 mb12 = b23 b24a31 a32 a33 a34 b31 b32 b33 b34ma21 = a41 a42 ma22 = a43 a44 mb21 = b41 b42 mb22 = b43 b44实现// 计算2X2矩阵void Multiply2X2(float& fOut_11, float& fOut_12, float& fOut_21, float& fOut_22,float f1_11, float f1_12, float f1_21, float f1_22,float f2_11, float f2_12, float f2_21, float f2_22){const float x1((f1_11 + f1_22) * (f2_11 + f2_22));const float x2((f1_21 + f1_22) * f2_11);const float x3(f1_11 * (f2_12 - f2_22));const float x4(f1_22 * (f2_21 - f2_11));const float x5((f1_11 + f1_12) * f2_22);const float x6((f1_21 - f1_11) * (f2_11 + f2_12));const float x7((f1_12 - f1_22) * (f2_21 + f2_22));fOut_11 = x1 + x4 - x5 + x7;fOut_12 = x3 + x5;fOut_21 = x2 + x4;fOut_22 = x1 - x2 + x3 + x6;}// 计算4X4矩阵void Multiply(CLAYMATRIX& mOut, const CLAYMATRIX& m1, const CLAYMATRIX& m2){float fTmp[7][4];// (ma11 + ma22) * (mb11 + mb22)Multiply2X2(fTmp[0][0], fTmp[0][1], fTmp[0][2], fTmp[0][3],m1._11 + m1._33, m1._12 + m1._34, m1._21 + m1._43, m1._22 + m1._44,m2._11 + m2._33, m2._12 + m2._34, m2._21 + m2._43, m2._22 + m2._44);// (ma21 + ma22) * mb11Multiply2X2(fTmp[1][0], fTmp[1][1], fTmp[1][2], fTmp[1][3],m1._31 + m1._33, m1._32 + m1._34, m1._41 + m1._43, m1._42 + m1._44,m2._11, m2._12, m2._21, m2._22);// ma11 * (mb12 - mb22)Multiply2X2(fTmp[2][0], fTmp[2][1], fTmp[2][2], fTmp[2][3],m1._11, m1._12, m1._21, m1._22,m2._13 - m2._33, m2._14 - m2._34, m2._23 - m2._43, m2._24 - m2._44);// ma22 * (mb21 - mb11)Multiply2X2(fTmp[3][0], fTmp[3][1], fTmp[3][2], fTmp[3][3],m1._33, m1._34, m1._43, m1._44,m2._31 - m2._11, m2._32 - m2._12, m2._41 - m2._21, m2._42 - m2._22);// (ma11 + ma12) * mb22Multiply2X2(fTmp[4][0], fTmp[4][1], fTmp[4][2], fTmp[4][3],m1._11 + m1._13, m1._12 + m1._14, m1._21 + m1._23, m1._22 + m1._24,m2._33, m2._34, m2._43, m2._44);// (ma21 - ma11) * (mb11 + mb12)Multiply2X2(fTmp[5][0], fTmp[5][1], fTmp[5][2], fTmp[5][3],m1._31 - m1._11, m1._32 - m1._12, m1._41 - m1._21, m1._42 - m1._22,m2._11 + m2._13, m2._12 + m2._14, m2._21 + m2._23, m2._22 + m2._24);// (ma12 - ma22) * (mb21 + mb22)Multiply2X2(fTmp[6][0], fTmp[6][1], fTmp[6][2], fTmp[6][3],m1._13 - m1._33, m1._14 - m1._34, m1._23 - m1._43, m1._24 - m1._44,m2._31 + m2._33, m2._32 + m2._34, m2._41 + m2._43, m2._42 + m2._44);// 第一块mOut._11 = fTmp[0][0] + fTmp[3][0] - fTmp[4][0] + fTmp[6][0];mOut._12 = fTmp[0][1] + fTmp[3][1] - fTmp[4][1] + fTmp[6][1];mOut._21 = fTmp[0][2] + fTmp[3][2] - fTmp[4][2] + fTmp[6][2];mOut._22 = fTmp[0][3] + fTmp[3][3] - fTmp[4][3] + fTmp[6][3];// 第二块mOut._13 = fTmp[2][0] + fTmp[4][0];mOut._14 = fTmp[2][1] + fTmp[4][1];mOut._23 = fTmp[2][2] + fTmp[4][2];mOut._24 = fTmp[2][3] + fTmp[4][3];// 第三块mOut._31 = fTmp[1][0] + fTmp[3][0];mOut._32 = fTmp[1][1] + fTmp[3][1];mOut._41 = fTmp[1][2] + fTmp[3][2];mOut._42 = fTmp[1][3] + fTmp[3][3];// 第四块mOut._33 = fTmp[0][0] - fTmp[1][0] + fTmp[2][0] + fTmp[5][0];mOut._34 = fTmp[0][1] - fTmp[1][1] + fTmp[2][1] + fTmp[5][1];mOut._43 = fTmp[0][2] - fTmp[1][2] + fTmp[2][2] + fTmp[5][2];mOut._44 = fTmp[0][3] - fTmp[1][3] + fTmp[2][3] + fTmp[5][3];}比较在标准的定义算法中我们需要进行n * n * n次乘法运算,新算法中我们需要进行7log2n次乘法,对于最常用的4阶矩阵:原算法新算法加法次数 48 72(48次加法,24次减法)乘法次数 64 49需要额外空间 16 * sizeof(float) 28 * sizeof(float)新算法要比原算法多了24次减法运算,少了15次乘法。
4-2.矩阵乘法的Strassen算法详解
4-2.矩阵乘法的Strassen 算法详解题⽬描述请编程实现矩阵乘法,并考虑当矩阵规模较⼤时的优化⽅法。
思路分析根据wikipedia 上的介绍:两个矩阵的乘法仅当第⼀个矩阵B 的列数和另⼀个矩阵A 的⾏数相等时才能定义。
如A 是m×n 矩阵和B 是n×p 矩阵,它们的乘积AB 是⼀个m×p 矩阵,它的⼀个元素其中 1 ≤ i ≤ m,1 ≤ j ≤ p 。
值得⼀提的是,矩阵乘法满⾜结合律和分配率,但并不满⾜交换律,如下图所⽰的这个例⼦,两个矩阵交换相乘后,结果变了:下⾯咱们来具体解决这个矩阵相乘的问题。
解法⼀、暴⼒解法其实,通过前⾯的分析,我们已经很明显的看出,两个具有相同维数的矩阵相乘,其复杂度为O (n^3),参考代码如下:解法⼆、Strassen 算法在解法⼀中,我们⽤了3个for 循环搞定矩阵乘法,但当两个矩阵的维度变得很⼤时,O (n^3)的时间复杂度将会变得很⼤,于是,我们需要找到⼀种更优的解法。
⼀般说来,当数据量⼀⼤时,我们往往会把⼤的数据分割成⼩的数据,各个分别处理。
遵此思路,如果丢给我们⼀个很⼤的两个矩阵呢,是否可以考虑分治的⽅法循序渐进处理各个⼩矩阵的相乘,因为我们知道⼀个矩阵是可以分成更多⼩的矩阵的。
如下图,当给定⼀个两个⼆维矩阵A B 时:这两个矩阵A B 相乘时,我们发现在相乘的过程中,有8次乘法运算,4次加法运算:矩阵乘法的复杂度主要就是体现在相乘上,⽽多⼀两次的加法并不会让复杂度上升太多。
故此,我们思考,是否可以让矩阵乘法的运算过程中乘法的运算次数减少,从⽽达到降低矩阵乘法的复杂度呢?答案是肯定的。
1969年,德国的⼀位数学家Strassen 证明O (N^3)的解法并不是矩阵乘法的最优算法,他做了⼀系列⼯作使得最终的时间复杂度降低到了O(n^2.80)。
他是怎么做到的呢?还是⽤上⽂A B 两个矩阵相乘的例⼦,他定义了7个变量:如此,Strassen算法的流程如下:;表⾯上看,Strassen 算法仅仅⽐通⽤矩阵相乘算法好⼀点,因为通⽤矩阵相乘算法时间复杂度是,⽽Strassen 算法复杂度只是。
大整数乘法和strassen矩阵乘法
大整数乘法和strassen矩阵乘法大整数乘法和Strassen矩阵乘法是计算机科学中两个非常重要的算法。
在我们的日常生活中,我们经常需要比较大的数字,例如,考虑到我们的身份证号码或者信用卡号码中的位数就很大。
对于这些大数字的计算,实现乘法运算的标准方法导致了效率低下。
这就要依靠大整数乘法的算法来完成。
同样的,矩阵乘法是人们常用的数据分析和机器学习等领域的基础算法之一,Strassen矩阵乘法则是一种可以在更短时间内完成的矩阵乘法算法。
在接下来的文档中,我将详细讲解大整数乘法和Strassen矩阵乘法。
一、大整数乘法大整数乘法是指对于两个比较大的数字,我们如何快速且准确的计算它们的乘积。
在介绍大整数乘法的算法之前,先考虑乘法的基本方法。
在日常乘法中,乘法运算是通过对乘数和被乘数的每一位进行相乘并将结果相加而得到的。
例如,计算96 ×57,我们将乘数96 和被乘数57 的每一位相乘,去得到:``` 96 × 57 ------- 672 (6 x 7) 480 (9 x 5) <<1> +57600 (9 x 5 << 2> ) ------- 5472 ```我们在这个过程中需要进行至Less 2次的延时计算,6次的加法,这在数字比较小时的时候是可行的。
但是,如果数字的位数变得更大,传统的乘法算法就会非常不切实际,执行效率非常慢。
在大整数乘法中,我们需要考虑实现优化的算法以处理大量位数的数字,其中最流行和普遍使用的算法有以下两种:1.传统的分治算法2.卡拉茨巴乘法1.传统的分治算法传统的分治算法涉及将大的数字分解为更小的数字部分进行乘法运算,并且在计算此过程之间能够快速地进行组合。
该算法需要依靠递归算法的思想,在整个运算过程中采用分治策略。
如果给定两个长度为n的数字,我们可以将这些数字分成两个较小,长度为 n/2 的数字,然后将它们相乘。
在递归调用中,相同的方法会被重复递归调用。
用Strassen算法实现矩阵乘法
用Strassen算法实现矩阵乘法#coding=utf-8'''第10题:用Strassen的思想实现两个n×n矩阵的乘法'''import mathfrom numpy import *#check函数用来检查num是否是2的整数次幂def check(num):i = 1;while (True):if (i > num):return Falseif (i == num):return Truei = i * 2#Multiply函数即为实现方法def Matrix_Multiply(A,B):#由于本题只考虑两个n×n矩阵的乘法,且n为2的整数次幂的情况,故做一些判断去掉不合法的输入if A.shape!=B.shape:print '本题只考虑两个n×n矩阵的乘法,且n为2的整数次幂的情况'return Noneif A.shape[0]!=A.shape[1]:print '本题只考虑两个n×n矩阵的乘法,且n为2的整数次幂的情况'return Noneif check(A.shape[0])==False:print '本题只考虑两个n×n矩阵的乘法,且n为2的整数次幂的情况'return Nonen = A.shape[0]result = zeros([n,n])if(n == 2):#最基本的情况A11=A[0,0]A12=A[0,1]A21=A[1,0]A22=A[1,1]B11=B[0,0]B12=B[0,1]B21=B[1,0]B22=B[1,1]P1 = A11*(B12-B22)P2 = (A11+A12)*B22P3 = (A21+A22)*B11P4 = A22*(B21-B11)P5 = (A11+A22)*(B11+B22)P6 = (A12-A22)*(B21+B22)P7 = (A11-A21)*(B11+B12)else:#输入复杂的情况就采取divide-and-conquer策略A11=A[:n/2,:n/2]A12=A[:n/2,n/2:]A21=A[n/2:,:n/2]A22=A[n/2:,n/2:]B11=B[:n/2,:n/2]B12=B[:n/2,n/2:]B21=B[n/2:,:n/2]B22=B[n/2:,n/2:]#combineP1 = Matrix_Multiply(A11,B12-B22)P2 = Matrix_Multiply(A11+A12,B22)P3 = Matrix_Multiply(A21+A22,B11)P4 = Matrix_Multiply(A22,B21-B11)P5 = Matrix_Multiply(A11+A22,B11+B22)P6 = Matrix_Multiply(A12-A22,B21+B22)P7 = Matrix_Multiply(A11-A21,B11+B12)result[:n/2,:n/2] = P4 + P5 + P6 - P2 #C11result[:n/2,n/2:] = P1 + P2 #C12result[n/2:,:n/2] = P3 + P4 #C21result[n/2:,n/2:] = P1 + P5 - P3 - P7 #C22return result'''demo'''print '矩阵A='a = array([[1.1,2,3,4],[5,6,7,8],[9,10,11,12],[13,14,15,16]]) print aprint '计算A×A'print ''result = Matrix_Multiply(a,a)print '本题所实现的方法结果是:' print resultprint '正确结果应为:'print dot(a,a)。
斯特拉森矩阵相乘算法(c语言实现)
斯特拉森矩阵相乘算法(c语⾔实现)我们所要介绍的斯特拉森矩阵相乘算法是德国数学家沃尔克·施特拉森 (Volker Strassen) 于1969年提出的,该算法的主要思想是⼀种分治思想,即将⼀个2n的⽅阵分解成4个2n-1的⼩⽅阵。
借助这种办法,任何有穷⽅阵都可以简化为有限个2×2⽅阵,所以今天我们主要介绍斯特拉森算法在2×2⽅阵上的应⽤。
⾸先我们假设有两个2×2矩阵,A,B.A = a11 a12B = b11 b12a21 a22 b21 b22我们设矩阵相乘的结果矩阵为C。
则C = c11 c12c21 c22由斯特拉森算法可得7个⼩矩阵1. m1 = (a12 - a22)(b21 + b22)2. m2 = (a11 + a22)(b11 + b22)3. m3 = (a21-a11)(b11+b12)4. m4 = (a11+a12) * b225. m5 = a11 * (b12 - b22)6. m6 = a22 * (b21 - b11)7. m7 = (a21 + a22) * b11易得1. c11 = m6 + m1 + m2 - m42. c12 = m5 + m43. c21 = m6 + m74. c22 = m5 + m3 + m2 - m7应⽤该种⽅法可将花费的时间由最开始的Ω(n3 )变成了O(n lg7)⼜因为lg7在2.80和2.81之间,所以其运⾏时间为O(n2.81),相⽐⼀开始有较为明显的优化(尤其是在处理较⼤的矩阵时)。
⽰例代码如下:1int a[2][2] = { 2,3,4,5 },2 b[2][2] = { 1,7,2,6 };//输⼊矩阵3int c[2][2];//建⽴结果矩阵4int m1 = (a[0][1] - a[1][1]) * (b[1][0] + b[1][1]),5 m2 = (a[0][0] + a[1][1]) * (b[0][0] + b[1][1]),6 m3 = (a[1][0] - a[0][0]) * (b[0][0] + b[0][1]),7 m4 = (a[0][0] + a[0][1]) * b[1][1],8 m5 = a[0][0] * (b[0][1] - b[1][1]),9 m6 = a[1][1] * (b[1][0] - b[0][0]),10 m7 = (a[1][0] + a[1][1]) * b[0][0];11 c[0][0] = m6 + m2 + m1 - m4;12 c[0][1] = m5 + m4;13 c[1][0] = m6 + m7;14 c[1][1] = m5 + m3 + m2 - m7;。
strassen矩阵乘法
1.将大小为n×n的矩阵A和B划分为nቤተ መጻሕፍቲ ባይዱ2×n/2的小矩阵。
2.计算7个n/2×n/2的小矩阵:
M1=(A11+A22)(B11+B22) M2=(A21+A22)B11 M3=A11(B12-B22) M4=A22(B21-B11) M5=(A11+A12)B22 M6=(A21-A11)(B11+B12) M7=(A12-A22)(B21+B22)
Strassen矩阵乘法是一种快速矩阵乘法算法,由Volker Strassen在1969年提出。它的时间复杂度为O(n^2.81),比传统的矩阵乘法算法(时间复杂度为O(n^3))要快得多。
Strassen矩阵乘法的基本原理是将大小为n×n的矩阵A和B划分为n/2×n/2的小矩阵,并通过递归调用Strassen矩阵乘法来求解。
3.计算结果
继续
矩阵C:
C11=M1+M4-M5+M7 C12=M3+M5 C21=M2+M4 C22=M1-M2+M3+M6
4.返回结果矩阵C。
由于Strassen矩阵乘法的时间复杂度较低,因此在计算机科学中常常使用它来优化矩阵乘法的运算效率。但是,由于Strassen矩阵乘法的常数较大,因此在实际应用中并不总是比传统的矩阵乘法算法要快。
实验二Strassen矩阵乘法
实验2 Strassen矩阵乘法一、 实验目的1.理解Strassen矩阵乘法的分治思想Strassen矩阵乘法的分治法设计模式是:半分+混合2.改进Strassen矩阵乘法对内存的需求若按Strassen矩阵乘法的直接表述实现,则空间复杂度将是O(3n2),本实验将试图改进这个方面。
3.Strassen矩阵乘法的性能问题改进Strassen矩阵乘法的内存需求,并不一定能改进Strassen矩阵乘法的效率,本实验将试图测试一些较大规模(n>=1024)的n阶方阵的Strassen矩阵乘,探讨其效率问题。
二、 实验环境C/C++编程环境或任何编程语言环境三、 实验内容1. Strassen矩阵乘法描述尽管Strassen矩阵乘法的实用价值在当今的多核计算环境下可能不是那么显著,但其理论价值仍值得我们研究。
Strassen矩阵乘法体现了一类重要的分治算法设计模式,即半分+混合,同样具有这种算法设计模式的是FFT(Fast Fourier Transform)-“由于FFT这个卓越算法在实践上的重要意义,有些人把它看作是有史以来人们发明的最重要的算法之一。
”[1]Strassen 矩阵乘法的基本思想,可由下述矩阵乘法概括:输入:两个n=2k 维方阵A 和B (若A 和B 的维度不是2k ,则通过增加元素为0的行和列,以使A 和B 均为2k 维的方阵)输出:n 维方阵C(1)1+41+443-11+24134.12-43+4113.123-11+224.3424.1314.11.2412.4-4344+=A B =A B =A B =A B =A B M =A B M =A M M B M M M(2)为方便表示,这里采用与书本不同的下标表示法,如对于1个n/2维矩阵14.14M ,下标14.14表示其由两个矩阵A 1+4和B 1+4乘积而成,A 1+4表示A1+A4,同理B 1+4表示B1+B4,同理12.4M 表示A1+2与B4的乘积。
矩阵乘法(分治法)
算法设计与分析实验报告二、模型拟制、算法设计和正确性证明:设A和B是两个n*n阶矩阵,求他们两的成绩矩阵C。
这里假设n是2的幂次方;A和B是两个n*n的矩阵,他们的乘积C=AB也是一个n*n的矩阵,矩阵C中的元素C[i][j]定义为C[i][j]= ,则每计算一个C[i][j],需要做n次乘法和n-1次加法。
因此计算C的n2个元素需要n3次乘法和n3- n2次加法。
因此,所需的时间复杂度是O(n3)。
但是使用分治法可以改进算法的时间复杂度。
这里,假设n是2的幂。
将矩阵A,B,C中每一矩阵都分成4个大小相等的子矩阵,每个子矩阵是(n/2)*(n/2)的方阵。
由此,可将方阵C=AB重写为因此可得:C11=A11B11+A12B21C12=A11B12+A12B22C21=A21B11+A22B22C22=A21B12+A22B22这样就将2个n阶矩阵的乘积变成计算8个n/2阶矩阵的乘积和4个n/2阶矩阵的加法。
当n=1时,2个1阶方阵的乘积可直接算出,只需要做一次乘法。
当子矩阵阶n>1时,为求两个子矩阵的乘积,可继续对两个子矩阵分块,直到子矩阵的阶为1。
由此,便产生了分治降阶的递归算法。
但是这个算法并没有降低算法的时间复杂度。
由strassen矩阵乘法,M1=A11(B12-B22)M2=(A11+A12)B22M3=(A21+A22)B11M4=A22(B21-B11)M5=(A11+A22)(B11+B22)M6=(A12-A22)(B21+B22)M7=(A11-A21)(B11+B12)C11=M5+M4-M2+M6C12=M1+M2C21=M3+M4C22=M5+M1-M3-M7四、程序实现和测试过程:程序测试过程(1)测试过程(2)源程序:#include<iostream>#include<math.h>#include<fstream>using namespace std;ifstream infile("123.txt",ios::in);void Input(int n,int **A){//infile>>n;for(int i=0;i<n;i++)for(int j=0;j<n;j++)infile>>A[i][j];}void Output(int n,int **A){for(int i=0;i<n;i++){for(int j=0;j<n;j++)cout<<A[i][j]<<'\t';cout<<endl;}cout<<endl;}void Divide(int n,int **A,int **A11,int **A12,int **A21,int **A22) {int i,j;for(i=0;i<n;i++)for(j=0;j<n;j++){A11[i][j]=A[i][j];A12[i][j]=A[i][j+n];A21[i][j]=A[i+n][j];A22[i][j]=A[i+n][j+n];}}void Unit(int n,int **A,int **A11,int **A12,int **A21,int **A22) {int i,j;for(i=0;i<n;i++)for(j=0;j<n;j++){A[i][j]=A11[i][j];A[i][j+n]=A12[i][j];A[i+n][j]=A21[i][j];A[i+n][j+n]=A22[i][j];}}void Sub(int n,int **A,int **B,int **C){int i,j;for(i=0;i<n;i++)for(j=0;j<n;j++)C[i][j]=A[i][j]-B[i][j];}void Add(int n,int **A,int **B,int **C){int i,j;for(i=0;i<n;i++)for(j=0;j<n;j++)C[i][j]=A[i][j]+B[i][j];}void Mul(int n,int **A,int **B,int **M){if(n==1)M[0][0]=A[0][0]*B[0][0];else{n=n/2;int **A11,**A12,**A21,**A22;int **B11,**B12,**B21,**B22;int **M11,**M12,**M21,**M22;int **M1,**M2,**M3,**M4,**M5,**M6,**M7;int **T1,**T2;A11=new int*[n];A12=new int*[n];A21=new int*[n];A22=new int*[n];B11=new int*[n];B12=new int*[n];B21=new int*[n];B22=new int*[n];M11=new int*[n];M12=new int*[n];M21=new int*[n];M22=new int*[n];M1=new int*[n];M2=new int*[n];M3=new int*[n];M4=new int*[n];M5=new int*[n];M6=new int*[n];M7=new int*[n];T1=new int*[n];T2=new int*[n];int i;for(i=0;i<n;i++){A11[i]=new int[n];A12[i]=new int[n];A21[i]=new int[n];A22[i]=new int[n];B11[i]=new int[n];B12[i]=new int[n];B21[i]=new int[n];B22[i]=new int[n];M11[i]=new int[n];M12[i]=new int[n];M21[i]=new int[n];M22[i]=new int[n];M1[i]=new int[n];M2[i]=new int[n];M3[i]=new int[n];M4[i]=new int[n];M5[i]=new int[n];M6[i]=new int[n];M7[i]=new int[n];T1[i]=new int[n];T2[i]=new int[n];}Divide(n,A,A11,A12,A21,A22);Divide(n,B,B11,B12,B21,B22);// cout<<"A11,A12,A21,A22"<<endl;// Output(n,A11);Output(n,A12);Output(n,A21);Output(n,A22); Sub(n,B12,B22,T1);// cout<<"B12-B22"<<endl;// Output(n,T1);Mul(n,A11,T1,M1);Add(n,A11,A12,T2);Mul(n,T2,B22,M2);Add(n,A21,A22,T1);Mul(n,T1,B11,M3);Sub(n,B21,B11,T1);Mul(n,A22,T1,M4);Add(n,A11,A22,T1);Add(n,B11,B22,T2);Mul(n,T1,T2,M5);Sub(n,A12,A22,T1);Add(n,B21,B22,T2);Mul(n,T1,T2,M6);Sub(n,A11,A21,T1);Add(n,B11,B12,T2);Mul(n,T1,T2,M7);Add(n,M5,M4,T1);Sub(n,T1,M2,T2);Add(n,T2,M6,M11);Add(n,M1,M2,M12);Add(n,M3,M4,M21);Add(n,M5,M1,T1);Sub(n,T1,M3,T2);Sub(n,T2,M7,M22);Unit(n,M,M11,M12,M21,M22);}}int main(){int n;cout<<"please input number n"<<endl;cin>>n;int **A,**B,**C;A=new int*[n];B=new int*[n];C=new int*[n];for(int i=0;i<n;i++){A[i]=new int[n];B[i]=new int[n];C[i]=new int[n];}Input(n,A);cout<<"A Matrix is"<<endl;Output(n,A);Input(n,B);cout<<"B Matrix is"<<endl;Output(n,B);Input(n,C);//Output(n,C);Mul(n,A,B,C);cout<<"The Product of A and B is"<<endl;Output(n,C);// cout<<n<<endl;in();return 0;}1 / 1。
Strassen矩阵乘法
// test_project.cpp : 定义控制台应用程序的入口点。
// Strassen矩阵乘法是通过递归实现的,它将一般情况下二阶矩阵乘法//(可扩展到n阶,但Strassen矩阵乘法要求n是的幂)所需的次乘法降低为次#include"stdafx.h"#include<iostream>using namespace std;const int N=4;//设置矩阵的大小,仅为测试用//const int N=4;//Matrix_Siz/2template<typename T>//矩阵加法void Matrix_Add(int n,T X[][N],T Y[][N],T Z[][N]);template<typename T>//矩阵减法void Matrix_Sub(int n,T X[][N],T Y[][N],T Z[][N]);template<typename T>//矩阵数据输入void input(int n,T p[][N]);template<typename T>void output(int n,T C[][N]);template<typename T>//矩阵乘法void Strassen_Matrix(int n,T A[][N],T B[][N],T C[][N]);void main(){/*double**X;double**Y;double**Z;X=new double*[N];for(int i=0;i<N;i++)X[i]=new double[N];Y=new double*[N];for(int i=0;i<N;i++)Y[i]=new doubleN];Y=new double*[N];for(int i=0;i<N;i++)Y[i]=new double[N];*/int X[N][N]={0},Y[N][N]={0},Z[N][N]={0};cout<<"请输入第一个矩阵的值:"<<endl;input(N,X);cout<<"请输入第二个矩阵的值:"<<endl;input(N,Y);Strassen_Matrix(N,X,Y,Z);output(N,Z);}template<typename T>//矩阵加法void Matrix_Add(int n,T X[][N],T Y[][N],T Z[][N]) {for(int i=0;i<n;i++)for(int j=0;j<n;j++)Z[i][j]=X[i][j]+Y[i][j];}template<typename T>//矩阵减法void Matrix_Sub(int n,T X[][N],T Y[][N],T Z[][N]){for(int i=0;i<n;i++)for(int j=0;j<n;j++)Z[i][j]=X[i][j]-Y[i][j];}template<typename T>//矩阵数据输入void input(int n,T p[][N]){for(int i=0;i<n;i++){cout<<"请输入矩阵"<<i+1<<"行的"<<n<<"个数"<<endl;for(int j=0;j<n;j++){cin>>p[i][j];}}}template<typename T>void output(int n,T C[][N]){cout<<"输出矩阵是:"<<endl;for(int i=0;i<n;i++){for(int j=0;j<n;j++){cout<<C[i][j]<<" ";}cout<<endl;}}template<typename T>void Strassen_Matrix(int n,T A[][N],T B[][N],T C[][N]) {if(n==2)//当为阶方阵时for(int i=0;i<2;i++){for(int j=0;j<2;j++){C[i][j]=0;for(int t=0;t<2;t++){C[i][j]=C[i][j]+A[i][t]*B[t][j];}}}else{/*intA11[n/2][n/2],A12[n/2][n/2],A21[n/2][n/2],A22[n/2][n/2] ,B11[n/2][n/2],B12[n/2][n/2],B21[n/2][n/2],B22[n/2][n /2],M1[n/2][n/2],M2[n/2][n/2],M3[n/2][n/2],M4[n/2][n/2], M5[n/2][n/2],M6[n/2][n/2],M7[n/2][n/2],TMP[n/2][n/2],TMP1[n/2][n/ 2];//注:TMP,TMP1只用于存放中间变量*/int A11[N][N],A12[N][N],A21[N][N],A22[N][N], B11[N][N],B12[N][N],B21[N][N],B22[N][N],M1[N][N],M2[N][N],M3[N][N],M4[N][N],M5[N][N],M6[N][N],M7[N][N],TMP[N][N],TMP1[N][N],C11[N][N],C12[N][N],C21[N][N],C22[N][N];for(int i=0;i<n/2;i++)for(int j=0;j<n/2;j++){A11[i][j]=A[i][j];A12[i][j]=A[i][n/2+j];A21[i][j]=A[n/2+i][j];A22[i][j]=A[n/2+i][n/2+j];B11[i][j]=B[i][j];B12[i][j]=B[i][n/2+j];B21[i][j]=B[n/2+i][j];B22[i][j]=B[n/2+i][n/2+j];}//计算M1Matrix_Sub(n/2,B12,B22,TMP);Strassen_Matrix(n/2,A11,TMP,M1);Matrix_Add(n/2,A11,A12,TMP);Strassen_Matrix(n/2,TMP,B22,M2);//计算M3Matrix_Add(n/2,A21,A22,TMP);Strassen_Matrix(n/2,TMP,B11,M3);//计算M4Matrix_Sub(n/2,B21,B11,TMP);Strassen_Matrix(n/2,A22,TMP,M4);//计算M5;Matrix_Add(n/2,A11,A22,TMP);Matrix_Add(n/2,B11,B22,TMP1);Strassen_Matrix(n/2,TMP,TMP1,M5);//计算M6Matrix_Sub(n/2,A12,A22,TMP);Matrix_Add(n/2,B21,B22,TMP1);Strassen_Matrix(n/2,TMP,TMP1,M6);Matrix_Sub(n/2,A11,A21,TMP);Matrix_Add(n/2,B11,B12,TMP1);Strassen_Matrix(n/2,TMP,TMP1,M7);//计算C11,Matrix_Add(n/2,M5,M4,TMP);Matrix_Sub(n/2,TMP,M2,TMP1);Matrix_Add(n/2,TMP1,M6,C11);//就算C22,Matrix_Add(n/2,M1,M2,C12);//计算C21Matrix_Add(n/2,M3,M4,C21);//计算C22Matrix_Add(n/2,M5,M1,TMP);Matrix_Sub(n/2,TMP,M3,TMP1);Matrix_Sub(n/2,TMP1,M7,C22);for(int i=0;i<n/2;i++){for(int j=0;j<n/2;j++){C[i][j]=C11[i][j];C[i][j+n/2]=C12[i][j]; C[i+n/2][j]=C21[i][j]; C[i+n/2][j+n/2]=C22[i][j];}}}}。
使用递归与分治法求解3.strassen矩阵乘法
一、概述1.介绍矩阵乘法的概念和意义2.引出递归与分治法在矩阵乘法中的应用二、传统矩阵乘法算法1.介绍传统的矩阵乘法算法原理2.分析传统算法的时间复杂度和空间复杂度3.讨论传统算法在大规模矩阵计算中的局限性三、Strassen矩阵乘法算法原理1.介绍Strassen算法的基本思想和原理2.引出递归与分治法在Strassen算法中的运用3.分析Strassen算法的时间复杂度和空间复杂度四、递归与分治法在Strassen算法中的运用1.详细解释递归与分治法在Strassen算法中的具体应用过程2.分析递归与分治法对算法性能的影响3.讨论递归与分治法在其他算法中的推广应用五、实例分析1.通过具体实例演示Strassen算法和传统算法的计算过程2.对比分析两种算法的计算效率和精度3.总结实例分析结果,展示递归与分治法在Strassen算法中的优势六、改进和优化1.讨论现有Strassen算法的局限性和不足2.提出改进和优化的方案,探讨递归与分治法在算法优化中的作用3.展望递归与分治法在矩阵计算领域的未来发展方向七、结论1.总结文中讨论的内容,强调递归与分治法在Strassen算法中的重要性和价值2.展望递归与分治法在矩阵计算领域的广阔应用前景3.对读者提出建议,鼓励更多的研究者投身于这一领域的研究和探索。
六、改进和优化1. Strassen算法的局限性和不足尽管Strassen算法在理论上具有较低的时间复杂度,但实际应用中也存在一些局限性和不足。
Strassen算法中涉及到的矩阵分块操作会引入额外的运算开销和存储开销,使得在小规模矩阵计算中,并不能体现出明显的优势。
Strassen算法要求矩阵的维度必须为2的幂次方,而实际场景中的矩阵往往难以满足这一条件,限制了算法的适用范围。
另外,由于Strassen算法引入了额外的递归调用,对于小规模矩阵,递归调用会使得算法的性能反而不如传统的矩阵乘法算法。
2. 改进和优化的方案针对Strassen算法的局限性和不足,可以考虑一些改进和优化的方案。
Strassen算法
问题分析
一般的矩阵乘法:
c[i, j ] k 1 a[i, k ] * b[k , j ]
k n
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,k]*b[k,j] 一般算法: T(n)=n3 stressen算法:T(n)=n3*log7
函数实现
主函数调用下面六个函数可以实现2的N次幂阶数的矩阵相乘 : int judgement(int n) /*判断用户输入的矩阵的阶数是 否为2的N次方*/ void divide(int *d,int *d11,int *d12,int *d21,int *d22,int n) /*分割矩阵*/ int *merge(int *d11,int *d12,int *d21,int *d22,int n) /*合并矩阵*/ int *minus(int *a,int *b,int n) /*两个矩阵相减*/ int *plus(int *a,int *b,int n) /*两个矩阵相加*/ int *strassen(int *a,int *b,int n) /*strassen算法 */
Strassen算法
学生:解思南 专业:计算机技术 LOGO 学号:31209027
算法回顾
1、用擅长的语言编程实现Strassen算法,并 进行演示。 2、Strassen要求阶n是2的幂,但这样的情况 很少。对任一偶数n,总有 n=m*2k, 将矩阵 分为m*m个的2k阶矩阵。小矩阵用Strassen相 乘,大矩阵用传统算法。设计出偶数阶矩阵相 乘的算法。 3、进一步:对任意的正整数n,如何求n阶矩 阵相乘?
Stranssen矩阵乘法
Strassen矩阵乘法A和B都是n*n的方阵,计算C=A* B。
这个通常做法的乘法计算量是0(exp(n,3)),为了把这个计算数量级从3降到 log2(7),log2(7)略小于2.81。
我们使用Strassen乘法,方法如下:输入方阵的阶,和方阵A,B,输出C。
测试数据:Sample In (in.txt)91 32 5 8 93 8 35 3 2 7 7 3 2 2 65 6 7 8 4 2 1 7 32 5 4 7 4 8 1 5 78 3 2 6 3 8 6 5 82 6 4 83 2 1 8 93 54 75 46 3 02 5 8 5 93 1 6 83 5 2 1 6 7 5 8 311 13 12 25 38 49 33 28 335 3 2 7 7 3 2 2 2352 6 7 8 43 2 1 47 4424 5 4 7 4 8 1 5 281 3 24 6 36 8 36 5 02 62 34 83 2 1 8 23 5 34 75 456 3 212 35 8 5 943 15 6 10 34 33 23 22 1 2 5 3Sample Out (out.txt)1031 1038 815 347 1313 346 487 337 268 1015 612 660 439 916 478 489 382 254 1166 653 506 430 1444 446 455 617 518 926 1026 763 447 1108 292 334 441 374 728 1137 992 612 1243 786 521 521 364 898 832 631 468 1391 294 344 423 396 1021 511 586 310 853 553 385 410 451 1483 814 774 483 1553 318 528 615 540 867 930 779 369 1319 488 491 357 350 Code:(strassen.cpp)#include<iostream>#include<fstream>#define N 65using namespace std;void show(int n,int array[][N]){for(int i=0;i<n;i++){for(int j=0;j<n;j++){cout<<array[i][j]<<" ";}cout<<endl;}cout<<endl;}/*strassen矩?阵¨®计?算?函¡¥数ºy,ê?申¦¨º明¡Â不?定¡§义°?*/ void strassen(int n,int A[][N],int B[][N],int C[][N]);/*矩?阵¨®加¨®法¤¡§*/void matMin(int a[][N],int b[][N],int c[][N],int n){for(int i=0;i<n;i++){for(int j=0;j<n;j++){c[i][j]=a[i][j]-b[i][j];}}}/*方¤?阵¨®减?法¤¡§*/void matPlu(int a[][N],int b[][N],int c[][N],int n){for(int i=0;i<n;i++){for(int j=0;j<n;j++){c[i][j]=a[i][j]+b[i][j];}}}/*m函¡¥数ºy的Ì?计?算?*/void funcM1(int m1[][N],int n,int a11[][N],int a12[][N],int a21[][N],int a22[][N],int b11[][N],int b12[][N],int b21[][N],int b22[][N]){int temp1[N][N];matPlu(a21,a22,temp1,n);matMin(temp1,a11,temp1,n);int temp2[N][N];matMin(b22,b12,temp2,n);matPlu(temp2,b11,temp2,n);strassen(n,temp1,temp2,m1);}void funcM2(int m2[][N],int n,int a11[][N],int a12[][N],int a21[][N],int a22[][N],int b11[][N],int b12[][N],int b21[][N],int b22[][N]){strassen(n,a11,b11,m2);}void funcM3(int m3[][N],int n,int a11[][N],int a12[][N],int a21[][N],int a22[][N],int b11[][N],int b12[][N],int b21[][N],int b22[][N]){strassen(n,a12,b21,m3);}void funcM4(int m4[][N],int n,int a11[][N],int a12[][N],int a21[][N],int a22[][N],int b11[][N],int b12[][N],int b21[][N],int b22[][N]){int temp1[N][N];matMin(a11,a21,temp1,n);int temp2[N][N];matMin(b22,b12,temp2,n);strassen(n,temp1,temp2,m4);}void funcM5(int m5[][N],int n,int a11[][N],int a12[][N],int a21[][N],int a22[][N],int b11[][N],int b12[][N],int b21[][N],int b22[][N]){int temp1[N][N];matPlu(a21,a22,temp1,n);int temp2[N][N];matMin(b12,b11,temp2,n);strassen(n,temp1,temp2,m5);}void funcM6(int m6[][N],int n,int a11[][N],int a12[][N],int a21[][N],int a22[][N],int b11[][N],intb12[][N],int b21[][N],int b22[][N]){int temp1[N][N];matMin(a12,a21,temp1,n);matPlu(temp1,a11,temp1,n);matMin(temp1,a22,temp1,n);strassen(n,temp1,b22,m6);}void funcM7(int m7[][N],int n,int a11[][N],int a12[][N],int a21[][N],int a22[][N],int b11[][N],intb12[][N],int b21[][N],int b22[][N]){int temp1[N][N];matPlu(b11,b22,temp1,n);matMin(temp1,b12,temp1,n);matMin(temp1,b21,temp1,n);strassen(n,a22,temp1,m7);}/*把ã?A【?】?B【?】?方¤?阵¨®分¤?别Àe拆e开a为aa11,a12,a21,a22和¨ªb11,b12,b21,b22*/void part(int n,int A[][N],int B[][N],int a11[][N],int a12[][N],int a21[][N],int a22[][N],int b11[][N],int b12[][N],int b21[][N],int b22[][N]){int i,j;for(i=0;i<n/2;i++)for(j=0;j<n/2;j++){a11[i][j]=A[i][j];a12[i][j]=A[i][j+n/2];a21[i][j]=A[i+n/2][j];a22[i][j]=A[i+n/2][j+n/2];b11[i][j]=B[i][j];b12[i][j]=B[i][j+n/2];b21[i][j]=B[i+n/2][j];b22[i][j]=B[i+n/2][j+n/2];}}/*用®?mi组Á¨¦建¡§方¤?阵¨®乘?积yC*/void comb(int n,int C[][N],int m[][N][N]){int c[2][2][N][N];//c00matPlu(m[2],m[3],c[0][0],n/2);//c01matPlu(m[1],m[2],c[0][1],n/2);matPlu(c[0][1],m[5],c[0][1],n/2);matPlu(c[0][1],m[6],c[0][1],n/2);//c10matPlu(m[1],m[2],c[1][0],n);matPlu(c[1][0],m[4],c[1][0],n/2);matMin(c[1][0],m[7],c[1][0],n/2);//c11matPlu(m[1],m[2],c[1][1],n/2);matPlu(c[1][1],m[4],c[1][1],n/2);matPlu(c[1][1],m[5],c[1][1],n/2);//combfor(int i=0;i<n/2;i++)for(int j=0;j<n/2;j++){C[i][j]=c[0][0][i][j];C[i][j+n/2]=c[0][1][i][j];C[i+n/2][j]=c[1][0][i][j];C[i+n/2][j+n/2]=c[1][1][i][j];}}void strassen(int n,int A[][N],int B[][N],int C[][N]){if(n==1){C[0][0]=A[0][0]*B[0][0];}else{int a11[N][N],a12[N][N],a21[N][N],a22[N][N];int b11[N][N],b12[N][N],b21[N][N],b22[N][N];int m[8][N][N];part(n,A,B,a11,a12,a21,a22,b11,b12,b21,b22);funcM1(m[1],n/2,a11,a12,a21,a22,b11,b12,b21,b22);funcM2(m[2],n/2,a11,a12,a21,a22,b11,b12,b21,b22);funcM3(m[3],n/2,a11,a12,a21,a22,b11,b12,b21,b22);funcM4(m[4],n/2,a11,a12,a21,a22,b11,b12,b21,b22);funcM5(m[5],n/2,a11,a12,a21,a22,b11,b12,b21,b22);funcM6(m[6],n/2,a11,a12,a21,a22,b11,b12,b21,b22);funcM7(m[7],n/2,a11,a12,a21,a22,b11,b12,b21,b22);comb(n,C,m);}}int main(){int i,j,n;int A[N][N];int B[N][N];int C[N][N];ifstream fin("in2.txt");ofstream fout("out22.txt");//输º?入¨?数ºy据Yfin>>n;for(i=0;i<n;i++)for(j=0;j<n;j++)fin>>A[i][j];for(i=0;i<n;i++)for(j=0;j<n;j++)fin>>B[i][j];//填¬?充?int m=1;while(m<n){m=m<<1;}fout<<m<<endl;for(i=n;i<m;i++)for(j=n;j<m;j++)A[i][j]=B[i][j]=0;//strassen矩?阵¨®乘?法¤¡§计?算?strassen(m,A,B,C);//结¨¢果?显?示º?for(i=0;i<n;i++){for(j=0;j<n;j++){fout<<C[i][j]<<" ";}fout<<endl;}}。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
实验报告课程名称:算法设计与分析班级:12软工一班实验成绩:实验名称:Strassen矩阵乘法学号:1242159101 批阅教师签字:实验编号:实验一姓名:陈双飞实验日期:2014年12月14日指导教师:陈平组号:实验时间:时分-时分一、实验目的通过算法的程序实现和执行时间测试、并与理论上的结论进行对比分析,深入理解算法时间复杂度分析中对于输入数据考虑其等价类的意义,理解算法时间复杂度渐进性态和和增长率的概念,为后续学习和实验奠定基础,同时也学习程序效率测试的基本思路。
二、实验内容与实验步骤(1)了解分治的基本思想并编写算法解决Strassen矩阵乘法问题。
(2)打开一台装有MyEclipse-10.7.1的PC。
(3)把已经写好的代码在Java环境下运行并调试。
(4)记录运行结果。
三、实验环境Windows 7家庭普通版,MyEclipse-10.7.1四、阐述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采用了类似于在大整数乘法中用过的分治技术,将计算2个n阶矩阵乘积所需的计算时间改进到O(nlog7)=O(n2.18)。
五、问题分析首先,我们还是需要假设n是2的幂。
将矩阵A,B和C中每一矩阵都分块成为4个大小相等的子矩阵,每个子矩阵都是n/2×n/2的方阵。
由此可将方程C=AB重写为:(1)由此可得:C11=A11B11+A12B21 (2)C12=A11B12+A12B22(3)C21=A21B11+A22B21(4)C22=A21B12+A22B22 (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*n2/4时间内完成,这里c是一个常数。
因此,上述分治法的计算时间耗费T(n)应该满足:这个递归方程的解仍然是T(n)=O(n3)。
因此,该方法并不比用原始定义直接计算更有效。
究其原因,乃是由于式(2)-(5)并没有减少矩阵的乘法次数。
而矩阵乘法耗费的时间要比矩阵加减法耗费的时间多得多。
要想改进矩阵乘法的计算时间复杂性,必须减少子矩阵乘法运算的次数。
按照上述分治法的思想可以看出,要想减少乘法运算次数,关键在于计算2个2阶方阵的乘积时,能否用少于8次的乘法运算。
Strassen提出了一种新的算法来计算2个2阶方阵的乘积。
他的算法只用了7次乘法运算,但增加了加、减法的运算次数。
这7次乘法是:M1=A11(B12-B22)M2=(A11+A12)B22M3=(A21+A22)B11M4=A22(B21-B11)M5=(A11+A22)(B11+B22)M6=(A12-A22)(B21+B22)M7=(A11-A21)(B11+B12)做了这7次乘法后,再做若干次加、减法就可以得到:C11=M5+M4-M2+M6C12=M1+M2C21=M3+M4C22=M5+M1-M3-M7以上计算的正确性很容易验证。
例如:C22=M5+M1-M3-M7=(A11+A22)(B11+B22)+A11(B12-B22)-(A21+A22)B11-(A11-A21)(B11+B12)=A11B11+A11B22+A22B11+A22B22+A11B12-A11B22-A21B11-A22B11-A11B11-A11B12+A21B11+A21B12=A21B12+A22B22由(2)式便知其正确性。
六、问题解决我们可以得到完整的Strassen算法如下:procedure STRASSEN(n,A,B,C);beginif 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);; end;end;其中MATRIX-MULTIPLY(A,B,C)是按通常的矩阵乘法计算C=AB的子算法。
Strassen矩阵乘积分治算法中,用了7次对于n/2阶矩阵乘积的递归调用和18次n/2阶矩阵的加减运算。
由此可知,该算法的所需的计算时间T(n)满足如下的递归方程:按照解递归方程的套用公式法,其解为T(n)=O(nlog7)≈O(n2.81)。
由此可见,Strassen矩阵乘法的计算时间复杂性比普通矩阵乘法有阶的改进。
有人曾列举了计算2个2阶矩阵乘法的36种不同方法。
但所有的方法都要做7次乘法。
除非能找到一种计算2阶方阵乘积的算法,使乘法的计算次数少于7次,按上述思路才有可能进一步改进矩阵乘积的计算时间的上界。
但是Hopcroft 和Kerr(197l)已经证明,计算2个2×2矩阵的乘积,7次乘法是必要的。
因此,要想进一步改进矩阵乘法的时间复杂性,就不能再寄希望于计算2×2矩阵的乘法次数的减少。
或许应当研究3×3或5×5矩阵的更好算法。
在Strassen之后又有许多算法改进了矩阵乘法的计算时间复杂性。
目前最好的计算时间上界是O(n2.367)。
而目前所知道的矩阵乘法的最好下界仍是它的平凡下界Ω(n2)。
因此到目前为止还无法确切知道矩阵乘法的时间复杂性。
关于这一研究课题还有许多工作可做。
七、实验结果总结实现算法代码:/*** Strassen矩阵乘法* */import java.util.*;public class 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();}}/*** 矩阵乘法,此处只是定义了2*2矩阵的乘法* */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];}}}/*** strassen矩阵算法* */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];}}//计算M1sub(B12, B22, temp);Mul(A11, temp, M1);//计算M2Add(A11, A12, temp);Mul(temp, B22, M2);//计算M3Add(A21, A22, temp);Mul(temp, B11, M3);//M4sub(B21, B11, temp);Mul(A22, temp, M4);//M5Add(A11, A22, temp1);Add(B11, B22, temp);Mul(temp1, temp, M5);//M6sub(A12, A22, temp1);Add(B21, B22, temp);Mul(temp1, temp, M6);//M7sub(A11, A21, temp1);Add(B11, B12, temp);Mul(temp1, temp, M7);//计算C11Add(M5, M4, temp1);sub(temp1, M2, temp);Add(temp, M6, C11);//计算C12Add(M1, M2, C12);//C21Add(M3, M4, C21);//C22Add(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;}运行截图:心得体会:这个关于用分治法解决的Strassen矩阵问题,是把2个n*n阶矩阵相乘问题的计算时间进行了较大的改进。