矩阵相乘(分治法)
2016-04-21 22:30
633 查看
一个简单的分治算法求矩阵相乘
C=A * B ,假设三个矩阵均为n×n,n为2的幂。可以对其分解为4个n/2×n/2的子矩阵分别递归求解:
递归分治算法:
算法中一个重要的细节就是在分块的时候,采用的是下标的方式。
=========== 王杰 原创作品转载请注明出处==============
C=A * B ,假设三个矩阵均为n×n,n为2的幂。可以对其分解为4个n/2×n/2的子矩阵分别递归求解:
递归分治算法:
算法中一个重要的细节就是在分块的时候,采用的是下标的方式。
#include <stdio.h> #include <stdlib.h> #define ROW 16 //指定 行数 #define COL 16 //指定 列数 int a[ROW][COL],b[ROW][COL]; //矩阵a 和 矩阵b int **c; // c = a * b //保存一个矩阵的第一个元素的位置,即左上角元素的下标 //如果加上一个长度就可以知道整个矩阵了 typedef struct { //这里没有指定一个矩阵的长度,在分块时应该加入长度,否则不知道子块矩阵的大小 int str,stc; //str行下标 ; strc列下标 }subarr; // 两矩阵arr、brr相加减 保存在temp中 void operate(int **arr,int **brr,subarr te,char op,int **temp,int len); //分治法 求矩阵相乘 ,sa,sb分别为矩阵a,b参加运算的首元素 int ** square_recursive(subarr sa,subarr sb,subarr sc,int len){ int n=len; int **temp; int i; // 申请一个临时矩阵,用于保存a*b temp=(int**)malloc(sizeof(int *)*n); for ( i=0;i<n;++i){ temp[i]=(int *)malloc(sizeof(int)*n); } // 长度为1 则直接相乘 if (n==1) { temp[0][0]=a[sa.str][sa.stc]*b[sb.str][sb.stc]; }else{ // 这里都是对下标进行初始化 // sa,sb,sc代表输入矩阵A,B,temp参加运算的首元素下标,因为进行分块后只进行特定子块的运算 //标号1,2,3,4 分别代表第一、二、三、四个子块 subarr sa1,sb1, sc1; subarr sa2,sb2, sc2; subarr sa3, sb3,sc3; subarr sa4, sb4, sc4; // 矩阵A 进行分块后的各个子块下标 sa1.str=sa.str; sa1.stc=sa.stc; sa2.str=sa.str; sa2.stc=sa.stc+n/2; sa3.stc=sa.stc; sa3.str=sa.str+n/2; sa4.str=sa.str+n/2; sa4.stc=sa.stc+n/2; // 矩阵B 进行分块后的各个子块下标 sb1.str=sb.str; sb1.stc=sb.stc; sb2.str=sb.str; sb2.stc=sb.stc+n/2; sb3.stc=sb.stc; sb3.str=sb.str+n/2; sb4.str=sb.str+n/2; sb4.stc=sb.stc+n/2; // 矩阵temp 进行分块后的各个子块下标 sc1.str=sc1.stc=0; sc2.str=0; sc2.stc=n/2; sc3.stc=0; sc3.str=n/2; sc4.str=n/2; sc4.stc=n/2; // 将矩阵分为四块 分别求解。采用下标的方式进行分块,可以省去复制矩阵所产生的时间 // 若要复制矩阵则会产生 O(n*n)的时间复杂度 operate(square_recursive(sa1,sb1,sc1,n/2),square_recursive(sa2,sb3,sc1,n/2),sc1,'+',temp,n/2); operate(square_recursive(sa1,sb2,sc2,n/2),square_recursive(sa2,sb4,sc2,n/2),sc2,'+',temp,n/2); operate(square_recursive(sa3,sb1,sc3,n/2),square_recursive(sa4,sb3,sc3,n/2),sc3,'+',temp,n/2); operate(square_recursive(sa3,sb2,sc4,n/2),square_recursive(sa4,sb4,sc4,n/2),sc4,'+',temp,n/2); } return temp; } // temp矩阵的te位置(四个子块中的一个)=arr+brr // len表示arr,brr参加运算的长度 // op是运算符 ‘+’ void operate(int **arr,int **brr,subarr te,char op,int **temp,int len){ int i,j; switch(op){ case '+': for (i=0;i<len;++i){ for (j = 0; j < len; ++j) { temp[te.str+i][te.stc+j]=arr[i][j]+brr[i][j]; } } break; case '-': for (i=0;i<len;++i){ for (j = 0; j < len; ++j) { temp[te.str+i][te.stc+j]=arr[i][j]-brr[i][j]; } } break; } } //为矩阵初始化 即赋值 void createarr(int temp[][COL]){ int i,j; for (i = 0; i < ROW; ++i) { for (j = 0; j < COL; ++j) { temp[i][j]=(int)rand()%5; } } } // 打印C矩阵 void print(){ int i,j; printf("\n====================================\n"); for (i = 0; i < ROW; ++i) { for (j = 0; j < COL; ++j) { printf("%d\t", c[i][j]); } printf("\n"); } printf("===================================\n"); } // 打印矩阵 void printarray(int a[ROW][COL]){ int i,j; printf("-----------------------\n"); for (i = 0; i < ROW; ++i) { for (j = 0; j < COL; ++j) { printf("%d \t", a[i][j]); } printf("\n"); } printf("-----------------------\n"); } int main(){ int i,j; subarr sa,sb,sc; int len; //初始化各个下标 sa.str=sa.stc=0; sb.str=sb.stc=0; sc.str=sc.stc=0; // 长度赋值,因为在subarr结构里没有长度的定义 len=ROW; //申请空间 c=(int**)malloc(sizeof(int *)*len); for (i=0;i<len;++i){ c[i]=(int *)malloc(sizeof(int)*len); } // 给矩阵A,B 复制初始化 createarr(a); createarr(b); // 进行运算 c=square_recursive(sa,sb,sc,len); // 打印矩阵A,B,C printarray(a); printarray(b); print(); return 0; }
=========== 王杰 原创作品转载请注明出处==============
相关文章推荐
- 如何组织构建多文件 C 语言程序(二)
- 如何写好 C main 函数
- Lua和C语言的交互详解
- 有关数据库SQL递归查询在不同数据库中的实现方法
- C#中的递归APS和CPS模式详解
- WinForm实现按名称递归查找控件的方法
- 使用SqlServer CTE递归查询处理树、图和层次结构
- C#中的尾递归与Continuation详解
- 关于C语言中参数的传值问题
- 简要对比C语言中三个用于退出进程的函数
- 深入C++中API的问题详解
- 基于C语言string函数的详解
- C语言中fchdir()函数和rewinddir()函数的使用详解
- C语言内存对齐实例详解
- C语言编程中统计输入的行数以及单词个数的方法
- C语言自动生成enum值和名字映射代码
- 使用C语言判断英文字符大小写的方法
- c语言实现的带通配符匹配算法
- C语言实现顺序表基本操作汇总
- C语言中计算正弦的相关函数总结