2017-02-21 79 views
0

我在写一个代码,其目的是使用Strassen算法乘以两个稠密矩阵。主代码是一个递归函数,其基本情况是当两个矩阵的大小是2×2时,因此矩阵的乘法只是实数的乘法。如果不满足基本情况,则通过再次调用递归函数计算用于计算最终乘法结果的七个矩阵,矩阵的大小除以一半。然而,当试图运行这段代码时,我遇到了段错误,但我无法弄清楚为什么。我使用gdb调试器,似乎在某些时候我得到NULL temp1和C1,但我不知道为什么会发生这种情况。另外,for循环中使用的计数器超出了限制(它们应该被限制在n/2以内)。最后,这是递归执行Strassen算法的正确方法吗?这是我的代码。关于在矩阵乘法中执行Strassen算法的C编程中的分段错误

#include "assignment2.h" 

void denseMatrixMult(int ** A, int ** B, int *** resultMatrix, int n) 
{ 
    if(n==2) 
    { 
     int M0,M1,M2,M3,M4,M5,M6; 
     M0=(A[0][0]+A[1][1])*(B[0][0]+B[1][1]); 
     M1=(A[1][0]+A[1][1])*B[0][0]; 
     M2=A[0][0]*(B[0][1]-B[1][1]); 
     M3=A[1][1]*(B[1][0]-B[0][0]); 
     M4=(A[0][0]+A[0][1])*B[1][1]; 
     M5=(A[1][0]-A[0][0])*(B[0][0]+B[0][1]); 
     M6=(A[0][1]-A[1][1])*(B[1][0]+B[1][1]); 
     int** temp; 
     initMatrix(&temp,2); 
     resultMatrix=&temp; 
     *(resultMatrix)[0][0]=M0+M3-M4+M6; 
     *(resultMatrix)[0][1]=M2+M4; 
     *(resultMatrix)[1][0]=M1+M3; 
     *(resultMatrix)[1][1]=M0-M1+M2+M5; 
     /*free(freeMatrix(temp);*/ 
     return; 
    } 
    else 
    { 
     int a,b,c,d,e,f,g,h; 
     int** N0; 
     int** N1; 
     int** N2; 
     int** N3; 
     int** N4; 
     int** N5; 
     int** N6; 
     int** zero; 
     int** C1; 
     int** C2; 
     int** C3; 
     int** C4; 
     initMatrix(&N0,n/2); 
     initMatrix(&N1,n/2); 
     initMatrix(&N2,n/2); 
     initMatrix(&N3,n/2); 
     initMatrix(&N4,n/2); 
     initMatrix(&N5,n/2); 
     initMatrix(&N6,n/2); 
     initMatrix(&zero,n/2); 
     denseMatrixMult(sum(A,A,0,0,n/2,n/2,n/2),sum(B,B,0,0,n/2,n/2,n/2),&N0,n/2); 
     denseMatrixMult(sum(A,A,n/2,0,n/2,n/2,n/2),sum(B,zero,0,0,0,0,n/2),&N1,n/2); 
     denseMatrixMult(sum(A,zero,0,0,0,0,n/2),sub(B,B,0,n/2,n/2,n/2,n/2),&N2,n/2); 
     denseMatrixMult(sum(A,zero,n/2,n/2,0,0,n/2),sub(B,B,n/2,0,0,0,n/2),&N3,n/2); 
     denseMatrixMult(sum(A,A,0,0,0,n/2,n/2),sum(B,zero,n/2,n/2,0,0,n/2),&N4,n/2); 
     denseMatrixMult(sub(A,A,n/2,0,0,0,n/2),sum(B,B,0,0,0,n/2,n/2),&N5,n/2); 
     denseMatrixMult(sub(A,A,0,n/2,n/2,n/2,n/2),sum(B,B,n/2,0,n/2,n/2,n/2),&N6,n/2); 
     C1=sum(sub(sum(N0,N3,0,0,0,0,n/2),N4,0,0,0,0,n/2),N6,0,0,0,0,n/2); 
     C2=sum(N2,N4,0,0,0,0,n/2); 
     C3=sum(N1,N3,0,0,0,0,n/2); 
     C4=sum(sum(sub(N0,N1,0,0,0,0,n/2),N2,0,0,0,0,n/2),N5,0,0,0,0,n/2); 
     int** temp1; 
     initMatrix(&temp1,n); 
     resultMatrix=&temp1; 
     for(a=0;a<n/2;a++) 
     { 
      for(b=0;b<n/2;b++) 
      { 
       (*resultMatrix)[a][b]=C1[a][b]; 
      } 
     } 
     for(c=n/2;c<n;c++) 
     { 
      for(d=0;d<n/2;d++) 
      { 
       (*resultMatrix)[c][d]=C3[c-n/2][d]; 
      } 
     } 
     for(e=0;e<n/2;e++) 
     { 
      for(f=n/2;f<n;f++) 
      { 
       (*resultMatrix)[e][f]=C2[e][f-n/2]; 
      } 
     } 
     for(g=n/2;g<n;g++) 
     { 
      for(h=n/2;h<n;h++) 
      { 
       (*resultMatrix)[g][h]=C4[g-n/2][h-n/2]; 
      } 
      } 
      /*freeMatrix(N0); 
      freeMatrix(N1); 
      freeMatrix(N2); 
      freeMatrix(N3); 
      freeMatrix(N4); 
      freeMatrix(N5); 
      freeMatrix(N6); 
      freeMatrix(C1); 
      freeMatrix(C2); 
      freeMatrix(C3); 
      freeMatrix(C4);*/ 


    } 
} 
int **sum(int ** A, int ** B, int x1, int y1, int x2, int y2, int n) 
{ 
    int i,j,k; 

    int ** res=(int**)malloc(n*sizeof(int*)); 
    if(res!=NULL) 
    { 
     for(i=0;i<n;i++) 
     { 
      res[i]=(int*)malloc(n*sizeof(int)); 
     } 

     for(j=0;j<n;j++) 
     { 
      for(k=0;k<n;k++) 
      { 
       res[j][k]=A[x1+j][y1+k]+B[x2+j][y2+k]; 
      } 
     } 
    } 
    return res; 

} 
int **sub(int **A, int **B, int x1, int y1, int x2, int y2, int n) 
{ 
    int i,j,k; 

    int ** res=(int**)malloc(n*sizeof(int*)); 
    if(res!=NULL) 
    { 
     for(i=0;i<n;i++) 
     { 
      res[i]=(int*)malloc(n*sizeof(int)); 
     } 

     for(j=0;j<n;j++) 
     { 
      for(k=0;k<n;k++) 
      { 
       res[j][k]=A[x1+j][y1+k]-B[x2+j][y2+k]; 
      } 
     } 
    } 
    return res; 
} 
void initMatrix(int ***mat,int n) 
{ 
    int i,j,k; 
    int ** Mat=(int**)malloc(n*sizeof(int*)); 
    for(i=0;i<n;i++) 
    { 
     Mat[i]=(int*)malloc(n*sizeof(int)); 
    } 
    for(j=0;i<n;i++) 
    { 
     for(k=0;k<n;i++) 
     { 
      Mat[j][k]=0; 
     } 
    } 
    *mat=Mat; 
} 

int **readMatrix(char * filename) 
{ 
    int i,j,k; 
    int **mat=(int**)malloc(MATSIZE*sizeof(int*)); 
    for(i=0;i<MATSIZE;i++) 
    { 
     mat[i]=(int*)malloc(MATSIZE*sizeof(int)); 
    } 
    FILE *fp=fopen(filename,"r"); 
    for(j=0;j<MATSIZE;j++) 
    { 
     for(k=0;k<MATSIZE;k++) 
     { 
      fscanf(fp,"%d",&mat[j][k]); 
     } 
    } 
    fclose(fp); 
    return mat; 

} 
void freeMatrix(int n, int ** matrix) 
{ 
    int i; 
    for(i=0;i<n;i++) 
    { 
     free(matrix[i]); 
    } 
    free(matrix); 
} 

void printMatrix(int n, int ** A) 
{ 
    int i,j; 
    for(i=0;i<n;i++) 
    { 
     for(j=0;j<n;j++) 
     { 
      printf("%d",A[i][j]); 
     } 
    } 
} 
#include "assignment2.h" 

void p1(void) 
{ 
    int **matrix; 
    initMatrix(&matrix,MATSIZE); 
    printMatrix(MATSIZE,matrix); 
    freeMatrix(MATSIZE, matrix); 
} 

void p2(void) 
{ 
    int ** matrix1=readMatrix("matrix1.txt"); 
    printMatrix(MATSIZE,matrix1); 
    freeMatrix(MATSIZE, matrix1); 
} 

void p3a(void) 
{ 
    int ** matrix1=readMatrix("matrix1.txt"); 
    int ** matrix2=readMatrix("matrix2.txt"); 
    int ** sumMatrix = sum(matrix1, matrix2, 1, 1, 0, 1, 3); 
    printMatrix(MATSIZE,matrix1); 
    printMatrix(MATSIZE,matrix2); 
    printMatrix(3,sumMatrix); 
    freeMatrix(MATSIZE, matrix1); 
    freeMatrix(MATSIZE, matrix2); 
    freeMatrix(3, sumMatrix); 
} 

void p3b(void) 
{ 
    int ** matrix1=readMatrix("matrix1.txt"); 
    int ** matrix2=readMatrix("matrix2.txt"); 
    int ** subMatrix = sub(matrix1, matrix2, 1, 1, 0, 1, 3); 
    printMatrix(MATSIZE,matrix1); 
    printMatrix(MATSIZE,matrix2); 
    printMatrix(3,subMatrix); 
    freeMatrix(MATSIZE, matrix1); 
    freeMatrix(MATSIZE, matrix2); 
    freeMatrix(3, subMatrix); 
} 

void p4(void) 
{ 
    char dataFileMat1[]="matrix1.txt"; 
    char dataFileMat2[]="matrix2.txt"; 
    int ** matrix1=readMatrix(dataFileMat1); 
    int ** matrix2=readMatrix(dataFileMat2); 
    int ** resultingMatrix; 
    denseMatrixMult(matrix1, matrix2, &resultingMatrix, MATSIZE); 
    printMatrix(MATSIZE,resultingMatrix); 
    freeMatrix(MATSIZE,resultingMatrix); 
    freeMatrix(MATSIZE,matrix1); 
    freeMatrix(MATSIZE,matrix2); 
} 

int main(int argc, char *argv[]) 
{ 
    if(argc < 2) 
    { 
     printf("Expecting at least one argument. Please try again\n"); 
    } 
    else if(argc==2) 
    { 
     if(atoi(argv[1])==3) 
     { 
      printf("Expecting two arguments for this part. Please try again.\n"); 
     } 
     else 
     { 
      if(atoi(argv[1])==1) 
      { 
       p1(); 
      } 
      else if(atoi(argv[1])==2) 
      { 
       p2(); 
      } 
      else if(atoi(argv[1])==4) 
      { 
       p4(); 
      } 
      else 
      { 
       printf("Incorrect argument supplied.\n"); 
      } 
     } 
    } 
    else if(argc==3) 
    { 
     if(atoi(argv[1])!=3) 
     { 
      printf("Expecting two arguments only for Part 3. Please try again.\n"); 
     } 
     else 
     { 
      if(atoi(argv[2])==1) 
      { 
       p3a(); 
      } 
      else if(atoi(argv[2])==2) 
      { 
       p3b(); 
      } 
     } 
    } 
    else 
    { 
     printf("The argument supplied is %s\n", argv[1]); 
    } 
} 
+2

Stackoverflow不是“转储代码并要求其他人为您调试它”。特别是不是300多行代码,其中有一两个字母变量。建议您使用[valgrind](http://valgrind.org)这样的工具来帮助您缩小根本原因。 – kaylum

+0

似乎问题出现在第20行和第47行。 –

回答

0

的失败free让这堆管理已损坏的嫌疑。可能的原因是超出范围写入访问分配的内存。 为什么不简化矩阵的分配?

只有一个malloc()需要分配。为此,您可以分配适当大小的一维数组。

  1. 您可以使用分配的数组连续存储矩阵元素(行的行)。读取和写入访问必须以i * N + ji ...行索引,j ...列索引,N ...列数)的形式完成。
  2. 你可以使用一些C-cast-magic,甚至可以用它作为二维数组。

予制备的小样本来演示此:

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

enum { N = 4 }; 

int main(int argc, char **argv) 
{ 
    int i, n; 
    /* allocation of an array with size NxN: */ 
    double *arr = (double*)malloc(N * N * sizeof (double)); 
    /* initialization of array (with some illustrative values) */ 
    for (i = 0, n = N * N; i < n; ++i) { 
    int row = i/N, col = i % N; 
    arr[i] = (double)(row + col * 0.1); 
    } 
    /* show contents */ 
    printf("arr[]:\n"); 
    for (i = 0; i < N; ++i) { 
    int j; 
    for (j = 0; j < N; ++j) { 
     printf("%s%.1f", j ? "\t" : "", arr[i * N + j]); 
    } 
    printf("\n"); 
    } 
    /* how this is used as two-dimensional array */ 
    { double (*mat)[N] = (double(*)[N])arr; 
    printf("mat[][]:\n"); 
    for (i = 0; i < N; ++i) { 
     int j; 
     for (j = 0; j < N; ++j) { 
     printf("%s%.1f", j ? "\t" : "", mat[i][j]); 
     } 
     printf("\n"); 
    } 
    } 
    /* clean-up */ 
    free(arr); 
    /* done */ 
    return 0; 
} 

编写和用gcc上cygwin的测试:

$ gcc -o test-mat test-mat.c 

$ ./test-mat 
arr[]: 
0.0  0.1  0.2  0.3 
1.0  1.1  1.2  1.3 
2.0  2.1  2.2  2.3 
3.0  3.1  3.2  3.3 
mat[][]: 
0.0  0.1  0.2  0.3 
1.0  1.1  1.2  1.3 
2.0  2.1  2.2  2.3 
3.0  3.1  3.2  3.3 

类型double (*)[N]必须被理解为“一个一个的地址大小为N的三维数组“。括号可能看起来令人惊讶,但是必要的。 (没有,编译器会将其读作“一维地址数组”)什么不是我的意图。)

+0

谢谢先生,我会检查一下。 –