高斯消元 解方程组的预备工作 将矩阵化为上三角
不分主从线程 每个线程负责一个方程
#include "mpi.h"
#include <stdio.h>
#include <stdlib.h>
typedef struct{
float value;
int rank;
} MD;
int main(int argc,char *argv[]) {
int self,size;
MPI_Status s;
//equivalents
float equs[][3]={{1,2,3},
{2,3,4}};
MPI_Init(&argc,&argv);
MPI_Comm_rank(MPI_COMM_WORLD,&self);
MPI_Comm_size(MPI_COMM_WORLD,&size);
float *equivalent=(float*)malloc((size+1)*sizeof(float));
float *recv=(float*)malloc((size+1)*sizeof(float));
int marked=0;
MD me,max;
//distribute
if(0==self) {
for(int i=0;i<size+1;++i) {
equivalent[i]=equs[0][i];
}
for(int i=1;i<size;++i) {
MPI_Ssend(equs[i],size+1,MPI_FLOAT,i,0,MPI_COMM_WORLD);
}
} else {
MPI_Recv(equivalent,size+1,MPI_FLOAT,0,0,MPI_COMM_WORLD,&s);
}
//pivot
for(int i=0;i<size;++i) {
if(0==marked) {
me.value=equivalent[i];
} else {
me.value=0;
}
me.rank=self;
MPI_Allreduce(&me,&max,1,MPI_FLOAT_INT,MPI_MAXLOC,MPI_COMM_WORLD);
if(self==max.rank) {
for(int i=0;i<size+1;++i) {
recv[i]=equivalent[i];
}
marked=1;
}
MPI_Bcast(recv,size+1,MPI_FLOAT,max.rank,MPI_COMM_WORLD);
if(0==marked) {
float ratio=equivalent[i]/recv[i];
for(int j=i;j<size+1;++j) {
equivalent[j]-=recv[j]*ratio;
}
}
}
printf("%d %f %f %f \n",self,equivalent[0],equivalent[1],equivalent[2]);
free(recv);
free(equivalent);
MPI_Finalize();
return 0;
}