算法分析与设计实验报告一12页

1、华中科技大学 算法分析与设计实验报告学生姓名：庞 亮 系别：软件学院 专业与班号：软件工程0805学号：U200818042实验时间：第 二 周，星期三 ，晚上 实验房间号：软件学院五楼机房实验名称Strassens 矩阵乘法和最近点对算法实验目的1、理解“分治法”算法设计思想及其实现步骤2、掌握分治算法效率递归分析方法3、掌握主方式求解递归式方法实验条件硬件：计算机软件：计算机程序语言开发平台，如C、C+、Java、Matlab。学生：至少掌握一门计算机程序设计语言，如C、C+、Java、Matlab。实验内容及要求1、利用计算机程序设计语言，实现教材第28.2 章介绍的“Strassens

2、 矩阵乘法算法”，自主生成两个88 的矩阵，检验算法的正确性并输出算法结果。2、比较 Strassens 矩阵乘法算法和数学定义的矩阵乘法算法效率之间的区别，并用直观的表达方式把两种不同矩阵乘法的效率随矩阵维数的变化趋势。3、利用计算机程序设计语言，实现教材第33.4 章介绍的“最近点对算法”，在拟定的二维空间点集上检验算法的正确性并输出算法结果。实验原理1. Strassens矩阵乘法简介Strassens算法是将矩阵分成了如图所示的均等的四块。分后的每一块儿任然还是方阵。所以可以由大问题分解成若干子问题进行解决。为了能使子问题能够返回到原始问题。 Strassens算法提出了如下的计算公式

3、，可以用矩阵的子矩阵计算出S1-S7，然后又由S1-S7合成原始矩阵。而S1-S7的计算又是方阵的乘法。由此使用分治算法便可以解决问题。2. 最近点对问题(Closest Pair Problems)算法简介首先这个问题也是采用了分治的思想，将空间内的距离分成三类，分界线左边的点之间的距离，分界线右边的点之间的距离，还有分界线两边距离为D的区域内的两点间距离。算法具体代码1.矩阵相乘问题/ Strassen.cpp : 定义控制台应用程序的入口点。/*#include stdafx.h#include stdlib.h#define AM Copy(A,0,0,A.v/2)#define BM

4、 Copy(A,A.v/2,0,A.v/2)#define CM Copy(A,0,A.v/2,A.v/2)#define DM Copy(A,A.v/2,A.v/2,A.v/2)#define EM Copy(B,0,0,B.v/2)#define FM Copy(B,B.v/2,0,B.v/2)#define GM Copy(B,0,B.v/2,B.v/2)#define HM Copy(B,B.v/2,B.v/2,B.v/2)#define V 2/*/矩阵结构typedef struct matrixint v;int x1616;Matrix; /*/输入输出文件FILE *fout

5、;FILE *fin;/*/矩阵打印（文件）void fPrint(Matrix A)for(int j=0;jA.v;j+)for(int i=0;iA.v;i+)fprintf(fout,%d ,A.xij);fprintf(fout,n);/*/矩阵打印（屏幕）void Print(Matrix A)for(int j=0;jA.v;j+)for(int i=0;iA.v;i+)printf(%d ,A.xij);printf(n);/*/矩阵截取Matrix Copy(Matrix X,int x,int y,int v)Matrix temp;temp.v=v;for(int i=x

6、;ix+v;i+)for(int j=y;jy+v;j+)temp.xi-xj-y=X.xij;return temp;/*/矩阵相减Matrix Minus(Matrix A,Matrix B)Matrix temp;temp.v=A.v;for(int j=0;jA.v;j+)for(int i=0;iA.v;i+)temp.xij=A.xij-B.xij;return temp;/*/矩阵相加Matrix Plus(Matrix A,Matrix B)Matrix temp;temp.v=A.v;for(int j=0;jA.v;j+)for(int i=0;iA.v;i+)temp.x

7、ij=A.xij+B.xij;return temp;/A: Copy(A,0,0,A.v/2)/B: Copy(A,A.v/2,0,A.v/2)/C: Copy(A,0,A.v/2,A.v/2)/D: Copy(A,A.v/2,A.v/2,A.v/2);/E: Copy(B,0,0,B.v/2)/F: Copy(B,B.v/2,0,B.v/2)/G: Copy(B,0,B.v/2,B.v/2)/H: Copy(B,B.v/2,B.v/2,B.v/2);/*/矩阵合并Matrix Merge4(Matrix A,Matrix B,Matrix C,Matrix D)Matrix temp;te

8、mp.v=A.v*2;for(int i=0;iA.v;i+)for(int j=0;jA.v;j+)temp.xij=A.xij;for(int i=0;iB.v;i+)for(int j=0;jB.v;j+)temp.xi+B.vj=B.xij;for(int i=0;iC.v;i+)for(int j=0;jC.v;j+)temp.xij+C.v=C.xij;for(int i=0;iD.v;i+)for(int j=0;jD.v;j+)temp.xi+D.vj+D.v=D.xij;return temp;/*/矩阵相乘Matrix Mutiply(Matrix A,Matrix B)i

9、f(A.v=1 & B.v=1)A.x00=A.x00*B.x00;return A;Matrix s1,s2,s3,s4,s5,s6,s7;s1=Mutiply(AM,Minus(FM,HM);s2=Mutiply(Plus(AM,BM),HM);s3=Mutiply(Plus(CM,DM),EM);s4=Mutiply(DM,Minus(GM,EM);s5=Mutiply(Plus(AM,DM),Plus(EM,HM);s6=Mutiply(Minus(BM,DM),Plus(GM,HM);s7=Mutiply(Minus(AM,CM),Plus(EM,FM);fPrint(Plus(Pl

10、us(s5,s6),Minus(s4,s2);fPrint(Plus(s1,s2);fPrint(Plus(s3,s4);fPrint(Minus(Minus(s1,s7),Minus(s3,s5);return Merge4(Plus(Plus(s5,s6),Minus(s4,s2),Plus(s1,s2),Plus(s3,s4),Minus(Minus(s1,s7),Minus(s3,s5);/*/数据输入Matrix RnMatrix(int v)Matrix temp;temp.v=v;for(int j=0;jV;j+)for(int i=0;ix (Point *)p2)-x ?

11、1 : -1;int Compare_y(const void * p1,const void * p2)return (Point *)p1)-y (Point *)p2)-y ? 1 : -1;void Sort_xy()qsort(p_set_x,n,sizeof(Point),Compare_x);qsort(p_set_y,n,sizeof(Point),Compare_y);void data_from_file()fin=fopen(in.txt,r);fscanf(fin,%d,&n);for(int i=0;in;i+)fscanf(fin,%d %d,&p_set_xi.x

12、,&p_set_xi.y);p_set_yi.x=p_set_xi.x;p_set_yi.y=p_set_xi.y;fclose(fin);double distance(int a,int b)return (p_set_xa.x-p_set_xb.x)*(p_set_xa.x-p_set_xb.x)+(p_set_xa.y-p_set_xb.y)*(p_set_xa.y-p_set_xb.y);double ClosestOf3(int s,int e)double min=65536e10;for(int i=s;i=e;i+)for(int j=i+1;jdistance(i,j)mi

13、n=distance(i,j);return min;double FindClosest(int s,int e)double dis1=65536e10,dis2=65536e10,dis;if(e-s+1=3)return ClosestOf3(s,e);elseif(s(s+e)/2)dis2=FindClosest(s+e)/2+1,e);dis=dis1dis2 ? dis2:dis1;if(dis0)printf(%d %d,s,e);int p,q;p=(s+e)/2;q=(s+e)/2;while(p_set_x(s+e)/2.x-p_set_xp.xs)p-;elsebre

14、ak;while(p_set_xq.x-p_set_x(s+e)/2.xdis) if(qe)q+;elsebreak;dis1=ClosestOf3(p,q);if(dis1dis ? dis:dis1;void ShowOut()for(int i=0;in;i+)printf(%d %d n,p_set_xi.x,p_set_xi.y);printf(%lf,FindClosest(0,n-1);int _tmain(int argc, _TCHAR* argv)data_from_file();Sort_xy();ShowOut();getchar();return 0;思考题1、 分

15、治法算法设计思想的三个基本步骤是什么？如何证明分治算法的正确性？三个基本步骤：分解、解决、合并。用数学归纳发可以证明分治算法的正确性。2、 利用主方式求解 Strassens 矩阵乘法和最近点对算法效率的递归分析结果。Strassens矩阵乘法的递归式为Tn=7Tn2+(n2)主定理的公式为：所以由主定理得Tn=nlg73、解释怎样修改 Strassens 矩阵乘法算法，使得它也可以用于大小不必为2 的幂的矩阵？假设矩阵的大小为N，找一个2M使它是大于N的最小M值，此时将矩阵其他的部分填充0，使其变为一个2M大小的矩阵，然后就可以利用Strassens矩阵算法了。计算的结果取左上角的N的大小的矩阵，那就是矩阵乘法的结果。实验体会首先，说说Strassens矩阵乘法，在数学理论上，这个算法是十分优秀的，计算出来的界限明显比一般的靠定义计算要简单的多。但是，有一个问题，我想只要是写了这个程序的同学都知道，也就是这个程序如果计算很大的矩阵栈会溢