MPI 高斯消元

高斯消元 解方程组的预备工作 将矩阵化为上三角

不分主从线程 每个线程负责一个方程

#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;
}