MPI编程实例（PI、AllGather、排序）

1.MPI 简介

MPI（Message Passing Interface）是目前最重要的一个基于消息传递的并行编程工具，它具有移植性好、功能强大、效率高等许多优点，而且有多种不同的免费、高效、实用的实现版本，几乎所有的并行计算机厂商都提供对它的支持，成为了事实上的并行编程标准。

MPI是一个库，而不是一门语言，因此对MPI的使用必须和特定的语言结合起来进行。MPI不是一个独立的自包含系统，而是建立在本地并行程序设计环境之上，其进程管理和I/O均由本地并行程序设计环境提供。例如，MPI可以建立在IBM SP2的POE/MPL之上，也可以建立在Intel Paragon的OSF/NX。除了这些商业版本的MPI实现，还有一些免费版的MPI实现，主要有MPICH，LAM和CHIMP。

2.实验环境

本文是在KD机群上跑的实验，在Windows下也可搭建MPI运行环境，本文在这里不做详细介绍。

编译mpi程序： mpicc demo.c –o demo.o；运行mpi程序： mpirun -np 4 ./demo.o

3.实例

(1)用MPI编程实现PI的计算。

算法描述：键盘输入步骤次数 n，并把 n 广播给本通信环境中的所有进程，通过MPI_Bcast(&n, 1, MPI_INT, 0, MPI_COMM_WORLD); 各个进程计算自己的mypi = width * sum; 通过 MPI_Reduce(&mypi,  &pi,  1,  MPI_DOUBLE,  MPI_SUM,  0, MPI_COMM_WORLD);由进程 0 进行归约，把每个进程计算出来的 mypi 进行相加（MPI_SUM），赋给pi 。

#include<stdio.h>
#include<mpi.h>
#include<math.h>

int main(int argc, char *argv[]){
     int my_rank, num_procs;
     int i, n = 0;
     double sum, width, local, mypi, pi;
     double start = 0.0, stop = 0.0;
     int proc_len;
     char processor_name[MPI_MAX_PROCESSOR_NAME];
     MPI_Init(&argc, &argv);
     MPI_Comm_size(MPI_COMM_WORLD,  &num_procs);
      MPI_Comm_rank(MPI_COMM_WORLD,  &my_rank);
     MPI_Get_processor_name(processor_name,  &proc_len);
     printf("Process %d of %d\n", my_rank, num_procs);
     if(my_rank == 0){
         printf("please give step number n:");
         scanf("%d", &n);
         printf("\n");
         start = MPI_Wtime();
     }
 //  printf("Process %d of %d\n", my_rank, num_procs);
     MPI_Bcast(&n, 1, MPI_INT, 0, MPI_COMM_WORLD);
     width = 1.0 / n;
      sum = 0.0;
     for(i = my_rank; i < n; i += num_procs){
         local = width * ((double)i + 0.5);
         sum += 4.0 / (1.0 + local * local);
     }
     mypi = width * sum;
     MPI_Reduce(&mypi,  &pi,  1,  MPI_DOUBLE,  MPI_SUM,  0,
             MPI_COMM_WORLD);
     if(my_rank == 0){
         printf("PI is %.20f\n", pi);
         stop = MPI_Wtime();
         printf("Time: %f on %s\n", stop-start, processor_name);
         fflush(stdout);
     }
     MPI_Finalize();
     return 0;
}

(2)采用MPI_SEND和MPI_RECV编写代码来实现MPI_Allgather的功能，并写一个完整的MPI程序测试并对比实现和MPI原有实现的性能。

算法描述：对每一个进程pi，使用MPI_Send向其他所有进程发送消息，同时使用MPI_Recv从其他所有进程接收消息。

#include<stdio.h>
#include<stdlib.h>
#include"mpi.h"
#include<memory.h>

void MPI_Allgather_my(int * senddata, int sendcount, MPI_Datatype senddatatype, int * recvdata, int recvcount, MPI_Datatype recvdatatype, MPI_Comm comm)
{
	int rank, size, i;
	MPI_Request request;
	MPI_Status status;
	MPI_Comm_rank(MPI_COMM_WORLD, &rank);
	MPI_Comm_size(MPI_COMM_WORLD, &size);
	for(i = 0; i < size; i++)
	{
		if(i != rank)
		{
			MPI_Send( senddata, sendcount, senddatatype, i, 1, MPI_COMM_WORLD );
			MPI_Recv( recvdata + i*recvcount, recvcount, recvdatatype, i, 1, MPI_COMM_WORLD, &status );
		}
	}
	memcpy(recvdata+rank*recvcount,senddata,sizeof(senddatatype)*sendcount);
}

int main(int argc, char* argv[])
{
	int i, rank, size, tag = 1;
	int senddata[500] = {1}, recvdata[500*32];
	double start_time, end_time, s_t, e_t;
	int count = 500;
	MPI_Init(&argc, &argv);
	MPI_Comm_rank(MPI_COMM_WORLD, &rank);
	MPI_Comm_size(MPI_COMM_WORLD, &size);
	start_time = MPI_Wtime();
	MPI_Allgather_my(senddata, count, MPI_INT, recvdata, count, MPI_INT, MPI_COMM_WORLD);
	end_time = MPI_Wtime();	
	MPI_Reduce(&start_time, &s_t, 1, MPI_DOUBLE, MPI_MIN, 0, MPI_COMM_WORLD);
	MPI_Reduce(&end_time, &e_t, 1, MPI_DOUBLE, MPI_MAX, 0, MPI_COMM_WORLD);
	if(rank == 0)
	{
		printf("myallgather : count = %d total time = %f\n", count, e_t - s_t);
	}
	MPI_Barrier(MPI_COMM_WORLD);
	
	start_time = MPI_Wtime();
	MPI_Allgather(senddata, count, MPI_INT, recvdata, count, MPI_INT, MPI_COMM_WORLD);
	end_time = MPI_Wtime();	
	MPI_Reduce(&start_time, &s_t, 1, MPI_DOUBLE, MPI_MIN, 0, MPI_COMM_WORLD);
	MPI_Reduce(&end_time, &e_t, 1, MPI_DOUBLE, MPI_MAX, 0, MPI_COMM_WORLD);
	if(rank == 0)
	{
		printf("allgather : count = %d total time = %f\n", count, e_t - s_t);
	}
	MPI_Finalize();
	return 0;
}

(3)利用MPI实现一个大数据量的排序算法。

算法描述：使用MPI_Send 和MPI_Receive 来进行数据的通信，使用MPI_Send 函数将数据分成n 个部分传送给n 个进程（n 为通信域Comm 包括的进程个数），然后各个
进程对各自所含数据进行快速排序，再使用MPI_Receiv 函数将各个进程都排好的数据传送给进程0，最后在进程0 中进行快排序。

 #include <stdio.h>
 #include <sys/time.h>
 #include <stdlib.h>
 #include <time.h>
 #include <mpi.h>
 #define length 1000000

 void swap(int *data, int i, int j) {
  int temp = data[i];
  data[i] = data[j];
  data[j] = temp;
 }

 int partition(int *data, int start, int end) {
     if (start >= end) return 0;
     int pivotValue = data[start];
     int low = start;
     int high = end - 1;
     while (low < high) {
         while (data[low] <= pivotValue && low < end) low++;
         while (data[high] > pivotValue && high > start) high--;
         if (low < high) swap(data, low, high);
     }
     swap(data, start, high);
     return high;
 }

 void quicksort(int *data, int start, int end) {
     if (end-start+1 < 2) return;
     int pivot = partition(data, start, end);
     quicksort(data, start, pivot);
     quicksort(data, pivot+1, end);
 }

 int main(int argc, char *argv[]) {
     MPI_Init(&argc, &argv);
     int rank, size;
     MPI_Comm_rank (MPI_COMM_WORLD, &rank);
     MPI_Comm_size (MPI_COMM_WORLD, &size);
     srand(time(0));
     int *data = (int*)malloc(sizeof(int)*length);
     int i;
     for (i=0; i<length/size; i++)
         data[i] = rand();
     MPI_Status status;
     if (rank == 0) {
         for (i=1; i<size; i++)
             MPI_Recv(data+i*length/size, length/size, MPI_INT, i, MPI_ANY_TAG, MPI_COMM_WORLD, &status);
     }
     else {
         MPI_Send(data,  length/size,  MPI_INT,  0,  0,  MPI_COMM_WORLD);
     }

     struct timeval start, end;
     gettimeofday(&start, 0);
     int s;
     int localDataSize =  length;
     int pivot;
     for (s=size; s > 1; s /= 2) {
         if (rank % s == 0) {
             pivot = partition(data, 0, localDataSize);
             MPI_Send(data+pivot,  localDataSize  -  pivot,MPI_INT, rank + s/2, 0, MPI_COMM_WORLD);
             localDataSize = pivot;
         }
         else if (rank % s == s/2) {
             MPI_Recv(data,  length,  MPI_INT,  rank  -  s/2,
                     MPI_ANY_TAG, MPI_COMM_WORLD, &status);
             MPI_Get_count(&status,  MPI_INT,
                     &localDataSize);
         }
     }
     quicksort(data, 0, localDataSize);
     gettimeofday(&end, 0);
     if (rank == 0)
     {
         float  time  =  (end.tv_sec  -  start.tv_sec)  +
             0.000001*(end.tv_usec - start.tv_usec);
         printf("Time: %f s\n", time);
     }
     MPI_Finalize();
     return 0;
}

4.小结

陈国良等．并行算法实践［M］．北京：高等教育出版社，2004.1．