2013-10-26 58 views
1

我有一個程序下面應該根據不同的索引來分析矩陣乘法的時間。但是一些函數給出了分段錯誤,其中兩個是matmul_kji(C ,A,B)和matmul_jki(C,A,B)。有人可以向我解釋我做錯了什麼和可能的改進。此代碼應檢查矩陣的時間從(1000 * 1000)到(10000 * 10000)。C++中的分段錯誤在做矩陣乘法時

#include <iostream> 
#include <stdio.h> 
#include <stdlib.h> 
#include <unistd.h> 
#include <math.h> 
#include <sys/time.h> 

using namespace std; 
# define MATSIZE 1000 


double get_wall_time(){ 
    struct timeval time; 
    if (gettimeofday(&time,NULL)){ 
     // Handle error 
     return 0; 
    } 
    return (double)time.tv_sec + (double)time.tv_usec * .000001; 
} 
double get_cpu_time(){ 
    return (double)clock()/CLOCKS_PER_SEC; 
} 

void init_mat (double M[MATSIZE][MATSIZE]) 
{ 
int i, j; 
for (i=0; i<MATSIZE; i++) { 
for (j=0; j<MATSIZE; j++) { 
M[i][j] = sin(double(i)); 
} 
} 
} 

void matmul_ijk (double C[MATSIZE][MATSIZE], 
double A[MATSIZE][MATSIZE], double B[MATSIZE][MATSIZE]) //1 
{ 
int i, j, k;   
double aux; 

double cpu0 = get_cpu_time(); 
for (i=0; i<MATSIZE; i++) { 
for (j=0; j<MATSIZE; j++) { 
aux = 0; 
for (k=0; k<MATSIZE; k++) { 
aux += A[i][k] *B[k][j]; 
} 
C[i][j] += aux; 
} 
} 
double cpu1 = get_cpu_time(); 
cout << "CPU Time Of IKJ Ordering = "<< cpu1-cpu0<<endl; 
}      

void matmul_ikj (double C[MATSIZE][MATSIZE], 
double A[MATSIZE][MATSIZE], double B[MATSIZE][MATSIZE]) //2 
{ 
int i, j, k;   
double aux; 

double cpu0 = get_cpu_time(); 
for (i=0; i<MATSIZE; i++) { 
for (k=0; k<MATSIZE; k++) { 
aux = 0; 
for (j=0; j<MATSIZE; j++) { 
aux += A[i][k] *B[k][j]; 
} 
C[i][j] += aux; 
} 
} 
double cpu1 = get_cpu_time(); 
cout << "CPU Time Of IKJ Ordering = "<< cpu1-cpu0<<endl; 
}   

void matmul_kij (double C[MATSIZE][MATSIZE], 
double A[MATSIZE][MATSIZE], double B[MATSIZE][MATSIZE]) //3 
{ 
int i, j, k;   
double aux; 

double cpu0 = get_cpu_time(); 
for (k=0; k<MATSIZE; k++) { 
for (i=0; i<MATSIZE; i++) { 
aux = 0; 
for (j=0; j<MATSIZE; j++) { 
aux += A[i][k] *B[k][j]; 
} 
C[i][j] += aux; 
} 
} 
double cpu1 = get_cpu_time(); 
cout << "CPU Time of KIJ ordering = "<< cpu1-cpu0<<endl; 
}      

void matmul_kji(double C[MATSIZE][MATSIZE], 
double A[MATSIZE][MATSIZE], double B[MATSIZE][MATSIZE]) //4 
{ 
int i, j, k;   
double aux; 

double cpu0 = get_cpu_time(); 
for (k=0; k<MATSIZE; k++) { 
for (j=0; j<MATSIZE; j++) { 
aux = 0; 
for (i=0; i<MATSIZE; i++) { 
aux += A[i][k] *B[k][j]; 
} 
C[i][j] += aux; 
} 
} 
double cpu1 = get_cpu_time(); 
cout << "CPU Time Of IKJ Ordering = "<< cpu1-cpu0<<endl; 
}      


void matmul_jik(double C[MATSIZE][MATSIZE], 
double A[MATSIZE][MATSIZE], double B[MATSIZE][MATSIZE]) //5 
{ 
int i, j, k;   
double aux; 

double cpu0 = get_cpu_time(); 
for (j=0; j<MATSIZE; j++) { 
for (i=0; i<MATSIZE; i++) { 
aux = 0; 
for (k=0; k<MATSIZE; k++) { 
aux += A[i][k] *B[k][j]; 
} 
C[i][j] += aux; 
} 
} 
double cpu1 = get_cpu_time(); 
cout << "CPU Time Of KJI Ordering = "<< cpu1-cpu0<<endl; 
}  

void matmul_jki (double C[MATSIZE][MATSIZE], 
double A[MATSIZE][MATSIZE], double B[MATSIZE][MATSIZE]) //6 
{ 
int i, j, k;   
double aux; 

double cpu0 = get_cpu_time(); 
for (j=0; j<MATSIZE; j++) { 
for (k=0; k<MATSIZE; k++) { 
aux = 0; 
for (i=0; i<MATSIZE; i++) { 
aux += A[i][k] *B[k][j]; 
} 
C[i][j] += aux; 
} 
} 
double cpu1 = get_cpu_time(); 
cout << "CPU Time Of JKI Ordering = "<< cpu1-cpu0<<endl; 
}     

/*void matmul_ikj(double C[MATSIZE][MATSIZE],    
double A[MATSIZE][MATSIZE], double B[MATSIZE][MATSIZE]) 
{ 
int i, j, k;   
double aux; 

double cpu0 = get_cpu_time(); 

for (i=0; i<MATSIZE; i++) { 
for (k=0; j<MATSIZE; j++) { 
aux = 0; 
for (j=0; k<MATSIZE; k++) { 
aux += A[i][k]*B[k][j]; 
} 
C[i][j] += aux; 
} 
} 
double cpu1 = get_cpu_time(); 
cout << "CPU Time = " << cpu1 - cpu0 << endl; 
} */   

int main() 
{ 
static double A[MATSIZE][MATSIZE]; 
static double B[MATSIZE][MATSIZE]; 
static double C[MATSIZE][MATSIZE]= {{0}}; 

init_mat (A); 
init_mat (B); 

//matmul_ijk(C,A,B); 
//matmul_ikj(C,A,B); 
//matmul_kij(C,A,B); 
matmul_kji(C,A,B);  // Gives a segmentation error 
//matmul_jik(C,A,B); 
matmul_jki(C,A,B); // Gives a segmentation error 

/*for (int i=0; i<10; i++){ 
    for(int j=0; j < 10; j++){ 
     cout << C[i][j]<<" "; 
    } 
    cout <<endl; } */  

return 0; 
}     

回答

2

這裏的問題是:

C[i][j] += aux; // Inside matmul_kji 

這是內部的兩個嵌套循環執行與變量jk,不ij,所以它應該是

C[j][k] = aux; 

注意您需要用=替換+=,因爲所有添加都是由t hird嵌套循環。

還要注意的是,編譯器可以幫助您避免這樣的問題,如果你使用的定義循環變量的現代*方式:而不是定義他們的前期,定義它們的循環水箱內,像這樣:

for (int k=0; k<MATSIZE; k++) { 
    for (int j=0; j<MATSIZE; j++) { 
     double aux = 0; 
     for (int i=0; i<MATSIZE; i++) { 
      aux += A[i][k] *B[k][j]; 
     } 
     C[i][j] += aux; // <<== Compile error! 
    } 
} 


*「現代」是非常相對的:聲明循環變量的這種方式已經出現了自八十年代中期,使得它很老。

+0

可能要還需要注意的是'AUX + = A [1] [ k] * B [k] [j];'無論如何可能是錯誤的 - 'B [k] [j]'應該依賴於'i',即最內層的循環變量。由於我不確定第一個或第二個索引是代表行/列還是代表「A」還是「B」是左邊的哪一個,所以不能100%清楚地知道*是正確的。乘法中的右矩陣。 –

0

你可能得到的是不要在外環設置i兩種情況分段故障,由於

C[i][j] += aux;