請看一下有關我的程式的錯誤(ABOUT平行處理”MPI”)請

Question

以下是A*B的矩陣(皆4*4)分給四台電腦做的程式(可1~4)
#include
#include
#include
#include\"mpi.h\"
main(int argc, char *argv[])
{
int rank;
int size;
int tag=100;
int N1=4;
int N2;
int N3;
int N4;
int N5;
int i,j,k;
int m,n;
int a[N1][N1];
int b[N1][N1];
int c[N1][N1];
int d[N1][N1];
double time1,time2;
MPI_Status status;
MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
MPI_Comm_size(MPI_COMM_WORLD, &size);

srand(time(NULL));

N2=N1/size;

if(rank==0)
{
for(i=0;i {
for(j=0;j {
a[i][j]=rand()%10; //設定A矩陣(亂數)//
b[i][j]=rand()%10; //設定B矩陣(亂數)//
c[i][j]=0;
d[i][j]=0;
}
}
printf(\"This is a[i][j]\\n\");
for(i=0;i {
for(j=0;j {
printf(\"%d \",a[i][j]);
}
printf(\"\\n\");
}

printf(\"This is b[i][j]\\n\");
for(i=0;i {
for(j=0;j {
printf(\"%d \",b[i][j]);
}
printf(\"\\n\");
}
}
MPI_Barrier(MPI_COMM_WORLD);
time1=MPI_Wtime(); //計算時間用//

MPI_Bcast(a,N1*N1,MPI_INT,0,MPI_COMM_WORLD); //把A矩陣傳給每一台工作的電腦//
MPI_Bcast(b,N1*N1,MPI_INT,0,MPI_COMM_WORLD);//把B矩陣傳給每一台工作的電腦//


printf(\"This is %d\\n\",rank);//每台電腦各自計算分配到的區域//
for(k=rank*N2;k {
for(j=0;j {
for(i=0;i {

// printf(\"m=%d=%d n=%d=%d \",m,n,a[k][i],b[i][j]);
c[k][j]+=a[k][i]*b[i][j];
// printf(\"c[%d][%d]=%d*%d \",k,j,a[k][i],b[i][j]);
// printf(\"a[%d][%d]*b[%d][%d]=%d*%d \",k,i,i,j,a[k][i],b[i][j]);
// printf(\"%d \",c[k][j]);

// if(rank==2&&k==2&&j==3)
// printf(\"c[%d][%d]=a[%d][%d]*b[%d][%d]=%d=%d*%d=%d \",k,j,k,i,i,j,c[k][j],a[k][i],b[i][j],a[k][i]*b[i][j]);
// printf(\"\\n\");

// if(rank==3)
// printf(\"a[%d][%d]*b[%d][%d]=%d*%d \",k,i,i,j,a[k][i],b[i][j]);
}
printf(\"\\n\");
}
printf(\"\\n\");
}
MPI_Gather (c,N2*N1,MPI_INT,d,N2*N1,MPI_INT,0,MPI_COMM_WORLD);
//把每台電腦所得到的答案存到d陣列中，並傳給編號為0的電腦//
time2=MPI_Wtime()-time1;//計算時間//

if(rank==0)
{
// printf(\"Matrix %dx%d with %d processor(s)\\n\",N1,N1,size);

for(k=0;k {
for(j=0;j {
printf(\"%d \",d[k][j]);//輸出陣列//
}
printf(\"\\n\");
}

printf(\"Total time:%f\\n\",time2);

MPI_Finalize();
}
return(0);
}

基本上這個程式我是設計成"無論幾X幾的矩陣，幾台電腦皆可計算"
(前提是矩陣的邊長數可以被電腦的總數給整除)
但是這個程式是測試用的，
故矩陣邊長及工作電腦數量皆設為4；
但是這個程式的問題，
已以下的跑出來的例子來做說明

This is a[i][j]
9 9 4 2
6 9 0 7
5 6 2 2
5 8 2 0
This is b[i][j]
5 0 8 1
4 4 7 7
0 9 2 4
6 2 2 3
這是亂數跑出的兩個矩陣

93 76 147 94
120 0 0 0
120 0 0 0
120 0 0 0
這是答案

在測試之中，給4台電腦跑出來的結果，
是只有process0所負責的部份可以跑出正確的答案，
無論是用多少台電腦都是一樣，
之後再加上一些測試用的程式碼，
發現無論是哪一台電腦，
皆可接收到自己所負責的部份(a,b矩陣有傳過去)，
但是在做矩陣相乘的時候所得到的答案卻是錯誤的，
希望對平行處理有研究的大大可以幫忙解答一下，謝謝。
(不知道這樣說明有沒有人聽的懂........)

看到你的回答了，
不知道是我說的不清不楚，還是這隻程式真的寫爛了。

因為之前有參考到類似的程式而做出的簡易程式
a(4*4)乘b(4*1)的，同樣也是分給4台電腦做。

這隻程式就可以跑出正確的答案，
現在有問題的程式是改過擴大版的......
(大大要看的話我在放上來)

對於你的回答我想解釋一下，
我的構想是這樣的，
a*b=c的計算，分給4台電腦去做(process0~3)
process0負責(c00,c01,c02,c03)的答案計算
process1負責(c10,c11,c12,c13)的答案計算
......以此類推

以矩陣相乘來說的話
c00=a00*b00+a01*b10+a02*b20+a03*b30
c01=a00*b01+a01*b11+a02*b21+a03*b31
c02=a00*b02+a01*b12+a02*b22+a03*b32
c01=a00*b03+a01*b13+a02*b23+a03*b33
這是process0所負責計算的部份

c10=a10*b00+a11*b10+a12*b20+a13*b30
c11=a10*b01+a11*b11+a12*b21+a13*b31
c12=a10*b02+a11*b12+a12*b22+a13*b32
c11=a10*b03+a11*b13+a12*b23+a13*b33
這是process1所負責計算的部份
......以此類推

以process0為例，
因為只要做c00,c01,c02,c03的計算，
由上面可得知，
a矩陣的資料只需要a00,a01,a02,a03而已，
b矩陣的資料則是需要全部。

那一開始我就把這兩個矩陣的資料都傳給各個process，
(其實a的部份只要傳1/4而已，我想之後再改)
然後各個process接收到之後，計算完之後再把答案由c傳到d，
再由0輸出，
(因為由process0來說，他所獲得的資料只有c00,c01,c02,c03，
process1則是c10,c11,c12,c13，所以想用Gather把答案按照順序放在d中再輸出)

至於大大在"齊行"的地方我看不懂，因為我並沒有"全部做一樣的東東"，

程式不好請見諒，謝謝。

Cannon's Algorithm

Answer 1

我 100% 能幫你！
但，你的 code 我不能跑啊！
出 error message!!
我目前沒空幫你從頭找！
除非你能等上 4 天！！
我 4 天後要 Thesis defense
內容就是你要的東東！
更複雜的距陣平行運算加速。

你用啥 compiler? 哪個 MPI ?

2006-07-24 04:30:06 補充：
你是要啥？
1. 每台電腦跑出一樣的內容？（大家都算一樣的東東，並沒平行。）
2. 每台電腦分別跑部份內容，交回給 0 號統一印出？（大家分工合作。）

2006-07-24 04:52:04 補充：
我大概看懂你的東東了。
你想要分工，但，程式這樣寫是錯的！！
你要去看 Cannon's Algorithm
不然，就算你算出來，答案也是錯的！！

在 n^2 個 CPU 下，做 N * N 的矩陣乘法，
1. 你整個 矩陣 BroadCast 走，做啥？
　只要 BCast 1 / n^2 就可以啦！
　難道你真的要大家同樣把 N * N 都做出來？
　那是哪門子的平行？

2. 不管你要平行，要齊行，你的 for (k .. ) 都是錯的！
　平行：請先去看 Cannon's Algorithm，您目前的 code 差太遠了。
　　　　除非你真的放棄，要從頭我幫你寫；
　　　　不然，真的沒辦法〝改〞給你。
　　　　你若要懂，看程式大概是沒辦法懂的。
　　　　請去看 Cannon's Algorithm 。
　　　　看完，把程式骨架改成 Cannon's Algorithm 的樣子，還跑不出來的話，我再幫你改。
　齊行：既然是全部做一樣的東東，每個 CPU 要做的東東和它的 rank 就無關了！
　　　　那，為何要 for (k=rank* N2; ...) ，這非炸不可！！
　　　　以你的例子： N2 = 2，在 rank == 3 那個 process跑時：
　　　　　for (k=3*2; k<3*2+2; ...) => for (k=6; k < 8; ...)
　　　　　c[k][j] += a[k][i] * b[i][j] => c[6][j] += a[6][i] ...
　　　　炸了沒？

3. 去看 MPI_Gather 的功用，它不是你要的樣子！差很多！

加油！ ^_^
我要去趕我的 present 了！

有新的進度，〝不幸〞不能完成，要幫忙，再 Post！
^_^

2006-07-24 04:57:20 補充：
版大，有看到我的回答的話，話通報一聲。不然，我可能會比你還緊張（不好意思，我個性太雞婆！）謝謝。

2006-07-24 08:51:40 補充：
齊行是開玩笑的，沒有人真的這麼做啦！！Cannon Alg. 是標準的平行陣列演算法。課本沒有嗎？我這裡２、３本都有耶。你想要用你的演算法，也沒關係。速度會慢一點就是了。但，無論如何，我都得等到台灣時間：週五早上8:00以後，才能開始大手筆地幫你了。這之前，若沒人能幫你，你要自立自強一下。加油。 ^_^

2006-07-25 00:12:50 補充：
你的 alg. 是 columnwise decompositionCan. alg. 是 checker-board decomposition，且是 highly scalable！速度差很多的！我還在拚命，你也加油！^_^

2006-07-28 04:33:56 補充：
我忙完了。
你 parallel matrix multiplication 還需要幫忙嗎?
不然, 我要去忙別的了!

2006-08-05 22:40:57 補充：
#include#include#include#include"mpi.h"#define N 4#define FREE(A) { free(*A); free(A); }typedef int TYP;void **dim2(void ***A,int h,int w){ int i,w0; char *buf,**idx; *A=0; if (idx=(char**)malloc(sizeof(char*)*h))

2006-08-05 22:44:45 補充：
if (buf=(char*)malloc(sizeof(TYP)*h*w)){w*= sizeof(TYP);for (w0=i=0; i
2006-08-05 22:50:19 補充：
for(j=0;j
2006-08-05 22:56:23 補充：
{TYP **A;字數滿了，見〝意見〞

2006-08-05 22:58:06 補充：
TYP **randMat(int h, int w)
{ TYP **A;
int i, j;

if (dim2((void***)&A, h, w))
for (i=0; i for (j=0; j A[i][j] = rand()%10;
return A;
}

2006-08-05 23:03:34 補充：
main(int argc, char *argv[])
{double ti;
TYP **A, **B, **C, **lA, **Re;
int nPrc, rank;
int i, j, k, l, m, n, ll, lm, ln;

MPI_Init(&argc, &argv);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);

2006-08-05 23:04:29 補充：
MPI_Comm_size(MPI_COMM_WORLD, &nPrc);
srand(time(NULL));
l = m = n = N;
ll = N / nPrc; lm = ln = N;

if (rank) dim2( (void ***)&B, m, n);
else
{ A = randMat(l, m), B = randMat(m, n);
ti = -MPI_Wtime();
}
dim2((void***)&lA,ll,lm);
dim2((void***)&Re, l, n);

2006-08-05 23:05:24 補充：
MPI_Scatter(*A, ll*lm, MPI_INT, *lA, ll*lm, MPI_INT, 0, MPI_COMM_WORLD);
MPI_Bcast (*B, m* n, MPI_INT, 0, MPI_COMM_WORLD);
C = inerp(lA, B, ll, lm, ln);

MPI_Gather(*C, ll*ln, MPI_INT, *Re, ll*ln, MPI_INT, 0, MPI_COMM_WORLD);

2006-08-05 23:05:39 補充：
ti += MPI_Wtime();
if (!rank)
{ showMat(Re, l, n, "This is C");
printf("Total time:%f\n", ti);
}

FREE(C); FREE(Re); FREE(lA); FREE(B);
if (!rank) FREE(A);

MPI_Finalize();

return 0;
}

2006-08-05 23:09:50 補充：
程式確定能 run，結果正確。
測過 1, 4, 6, 8 process
N = 4, 6, 8, 14, 18, 20, 22
int 和 double 都試過了
限於字數，我把 double 的 #define #ifdef 都殺了
測試答案正礭與否的部份也．．．
註解也殺了(本來要留給你的)
你要試試看，看程式有沒有被我殺〝死〞了！

有問題再問。
^_^

2006-08-05 23:21:07 補充：
這是稍微根據你的 rowwise 的想法設計的(B 是全部傳給每個process)，它並不是標準的 rowwise Algorithm.

另外，它的乘法副程式只是標準作法（O(n^3)），不是什麼加速版（可以到約（O(n^2.2)，我只做過（O(n^2.7)））。
唯〝二〞的小小加速是 用了 tmp ！
把它拿掉，會慢約 20%
把 i, j, k 的順序換了，也會變慢！

Answer 2

這是大學程度的阿........

2006-07-24 04:45:15 補充：
我是要2（大家分工合作。） ，以這個例子的話，電腦編號分別是0,1,2,3，各自把自己負責的部份做完之後，把答案放在矩陣d之中，在交由編號0的電腦輸出答案，但是輸出的結果之中，只有編號0的電腦所負責的部份是正確的。

我是用putty連上學校的主機跑的，之前都是用VC在自己的電腦去模擬多台；至於程式是關於MPICH2。

2006-08-02 05:25:40 補充：
Jacob Lee大大還請幫忙一下，最近在忙其他的事情，所以沒有時間上來看，謝謝~~

2006-08-08 15:13:53 補充：
程式經測試可以執行，不過還是希望可以把註解寄給我，謝謝
[e-mail]

Answer 3

目前在知識＋的解答者大多是高中、大學程度，如果有一、二人能解答你的問題，那算很厲害了。