有两个矩阵A和B，计算两矩阵的乘积存于C中；
方法一：每个进程分配到矩阵A的部分行，并且有矩阵B的所有元素，每个进程计算矩阵A的部分行与矩阵B的乘积，得到矩阵C的部分行，所有进程计算完毕后，将结果传递给根进程进行汇总，得到最后的结果；
#include <stdio.h>
#include <mpi.h>
#include <malloc.h>
#include <time.h>
#include <stdlib.h>
int main(int argc, char**argv){
int size, rank, i, j, k;
int N, len, indexC = 0;
int *A, *B, *C, *localA, *localC;//原矩阵和部分矩阵
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &size);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
if(rank == 0){
printf("the Matrix of N:");
scanf("%d", &N);//输入矩阵的阶数
C = (int*)malloc(sizeof(int)*N*N);
A = (int*)malloc(sizeof(int)*N*N);
for(i = 0; i < N*N; i++){
srand((unsigned int)i);
A[i] = rand() % 10;
}//赋值
}
MPI_Bcast(&N, 1, MPI_INT, 0, MPI_COMM_WORLD);
B = (int*)malloc(sizeof(int)*N*N);
for(i = 0; i < N*N; i++){
srand((unsigned int)i);
B[i] = rand() % 10;
}
len = N / size;//每个进程所得到的矩阵A的行数
localA = (int*)malloc(sizeof(int)*len*N);
localC = (int*)malloc(sizeof(int)*len*N);
MPI_Scatter(A, len*N, MPI_INT, localA, len*N, MPI_INT, 0, MPI_COMM_WORLD);
for(i = 0; i < len; i++){
for(j = 0; j < N; j++){
localC[indexC] = 0;
for(k = 0; k < N; k++){
localC[indexC] += localA[N*i+k] * B[N*k+j];
}
indexC++;
}
}//计算矩阵C的部分行
MPI_Gather(localC, len*N, MPI_INT, C, len*N, MPI_INT, 0, MPI_COMM_WORLD);
if(rank == 0){
printf("the Matrix of A:\n");
for(i = 0; i < N*N; i++)
printf("%3d%s", A[i], (i+1) % N == 0 ? "\n" : " ");
printf("the Matrix of B;\n");
for(i = 0; i < N*N; i++)
printf("%3d%s", B[i], (i+1) % N == 0 ? "\n" : " ");
printf("the Matrix of C:\n");
for(i = 0; i < N*N; i++)
printf("%3d%s", C[i], (i+1) % N == 0 ? "\n" : " ");
}
MPI_Finalize();
return 0;
}
方法二：每个进程分配到矩阵A的部分行和矩阵B的部分列，每一次的计算都将得到矩阵C的部分元素；每次计算完之后，一次将每个进程存储的矩阵A的部分行传递给下一个进程，再进行同样的计算，循环的次数为进程的数量；所有进程计算完毕，将结果发送给跟进程，得到最后的结果；
#include <stdio.h>
#include <mpi.h>
#include <malloc.h>
#include <stdlib.h>
#include <string.h>
int main(int argc, char**argv){
int size, rank, i, j, k;
int *A, *B, *C, *localA, *localB, *localC, *temp;
int N, len, circle = 0, index, cnt;//circle为循环的次数
int source, dest;
MPI_Init(&argc, &argv);
MPI_Comm_size(MPI_COMM_WORLD, &size);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
if(rank == 0){
printf("the Matrix of N:");
scanf("%d", &N);
A = (int*)malloc(sizeof(int)*N*N);
B = (int*)malloc(sizeof(int)*N*N);
C = (int*)malloc(sizeof(int)*N*N);
for(i = 0; i < N*N; i++){
srand((unsigned int)i);
A[i] = rand() % 10;
srand((unsigned int)i+1);
B[i] = rand() % 10;
}
/*
for(i = 0; i < N; i++){
for(j = 0; j<=i; j++){
k = B[N*i+j];
B[N*i+j] = B[N*j+i];
B[N*j+i] = k;
}
}
*/
}
MPI_Bcast(&N, 1, MPI_INT, 0, MPI_COMM_WORLD);
len = N / size;//每个进程分配的矩阵A的行数
localA = (int*)malloc(sizeof(int)*len*N);
localB = (int*)malloc(sizeof(int)*len*N);
localC = (int*)malloc(sizeof(int)*len*N);
temp = (int*)malloc(sizeof(int)*len*N);
MPI_Scatter(A, len*N, MPI_INT, localA, len*N, MPI_INT, 0, MPI_COMM_WORLD);
MPI_Scatter(B, len*N, MPI_INT, localB, len*N, MPI_INT, 0, MPI_COMM_WORLD);
index = rank;//指明当前进程实际所在进程号
//cnt = 0;
dest = rank==0 ? size-1 : rank-1;//发送的目标进程号
source = (rank+1) % size;//接收的目标进程号
do{
for(i = 0; i < len; i++){
for(j = 0; j < len; j++){
cnt = j*N + (index*len+i);
localC[cnt] = 0;
for(k = 0; k <N; k++)
localC[cnt] += localA[N*i+k] * localB[N*j+k];
}
//cnt++;
}//计算矩阵C中的部分元素
//一次传递矩阵A中的数据
if(rank == 0){
MPI_Ssend(localA, len*N, MPI_INT, dest, 99, MPI_COMM_WORLD);
MPI_Recv(localA, len*N, MPI_INT, source, 99, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
/*printf("Process %d:\n", rank);
for(i = 0; i < len*N; i++)
printf("%d%s", localA[i], (i+1)%N == 0?"\n":" ");*/
} else{
memcpy(temp, localA, sizeof(int)*len*N);
MPI_Recv(localA, len*N, MPI_INT, source, 99, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
MPI_Ssend(temp, len*N, MPI_INT, dest, 99, MPI_COMM_WORLD);
/*printf("Process %d:\n", rank);
for(i = 0; i < len*N; i++)
printf("%d%s", localA[i], (i+1)%N == 0?"\n":" ");*/
}
circle++;
index = (index+1) % size;
} while(circle < size);
MPI_Gather(localC, len*N, MPI_INT, C, len*N, MPI_INT, 0, MPI_COMM_WORLD);
if(rank == 0){
printf("the Matrix of A:\n");
for(i = 0; i < N*N; i++)
printf("%3d%s", A[i], (i+1) % N == 0 ? "\n" : " ");
printf("the Matirx of B:\n");
for(i = 0; i < N; i++)
for(j = 0; j < N; j++)
printf("%3d%s", B[N*j+i], (j+1) == N ? "\n" : " ");//矩阵B按列存
printf("the Matrix of C:\n");
for(i = 0; i < N; i++)
for(j = 0; j < N; j++)
printf("%3d%s", C[N*j+i], (j+1) == N ? "\n" : " ");
}
MPI_Finalize();
return 0;
}
