基于MPI的PSRS并行排序算法的实现
摘 要 本文介绍MPI并行编程环境和讨论和分析PSRS的并行排序算法,并在MP并行编程环境下实现PSRS的并行排序算法。
关键词 PSRS;MPI; 并行计算;消息传递接口;并行算法
一、 引言
并行计算是与多学科相联系的交叉学科, 研究范围十分广泛。经过二十多年的研究与发展, 硬件系统已经发展为第四代分布主存大规模并行机。并行计算现已牢固确立了它在高科技关键技术中的地位, 其主要目标是解决具有重大挑战意义的科学工程计算问题, 并行计算研究与应用在国内外均是最引人注目的前沿课题之一。随着并行计算机种类的迅速发展, 并行机的应用正从传统的科学工程计算逐步渗透到金融管理等多个部门, 并且由于网络技术的迅猛发展, 免费的公开并行软件环境(如PVM、MPI等) 以及多任务操作系统的应用与推广, 使国内许多科研、教育部门只需要拥有少量工作站甚至微机, 就能建立自己的网络并行计算环境。
据统计,在大型计算中心,仅排序工作往往就要占去大约1/4的计算时间。对排序算法的研究无论是在理论上还是在实践上都具有重大的意义。随着并行处理技术的发展,并行排序已经成为并行算法中的一个重要研究领域,目前已有大量关于对并行排序的讨论和研究的文章问世。本文讨论PSRS的并行排序算法,通过基于消息传递接口MPI (Message Passing Interface) 编程实现该并行算法及其运行结果的分析。
二、 MPI 并行编程环境
1、MPI 简介
MPI 是为开发基于消息传递模型的并行程序而制定的工业标准,其目的是为了提高并行程序的可移植性和易用性,用户必须显式地通过发送和接受消息来实现处理器之间的数据交换。为了开发一个高标准并具有可移植性的消息传递库, 由厂商(如IBM , Intel 等) 、软件开发商(Parasoft , KAI 等) 、研究中心(ANL , GMP 等) 和大学(Endinburgh ,Maryland 等) 共同于1992 年成立了MPI 论坛。经过共同研究和使用,MPI 论坛先后于1994 年和1997 年后提出了MPI-1 和MPI- 2 。MPI 具有移植性好、功能强大、效率高等多种优点,而且有多种不同的免费、高效、实用的实现版本,几乎所有的并行计算机厂商都提供对它的支持,这是其它并行环境所无法比拟的。MPI只是一个并行编程语言标准,要编写基于MPI 的并行程序,还必须借助MPI的某种具体实现。MPICH 是最重要的一种MPI 实现,是一个与MPI 规范同步发展的版本。每当MPI 推出新的版本时,就会有相应的MPICH 的实现版本,目前常用的版本是mpich- 1.2.6 ,在Windows下的一个常见版本是mpich.nt.1.2.5 ,文中的所有实验均采用此版本作为实验环境。CHIMP 是Edinburgh 开发的另一个免费MPI 实现,是在EPCC( Edinburgh Parallel Computing Centre) 的支持下进行的。LAM(Local Area Multicomputer) 是Linux平台下另一个免费的MPI实现。它由Ohio州立大学开发,主要用于异构的网格计算并行系统。
2、基于MPI 的程序设计方法
并行程序设计可以通过多种方法实现,它们各有特点,并适合于不同的场合。最常用的是共享变量方法和消息传递方法。共享变量方法的主要特点是利用共享存储器中共享变量来实现处理机间的通信,它主要用在共享存储结构中。消息传递方法的主要特点是不同处理机间必须通过网络传送消息来相互通信,并达到任务间需要的同步操作,它主要应用于分布式存储结构中。MPI 是基于消息传递模型的,在这种模型中,把任务分成若干部分(进程) ,让每个处理器(节点) 独自运行不同或者相同的部分(进程) 。驻留在不同处理器(节点) 上的进程可以通过网络传递消息来互相通信,从而达到并行计算的效果,减少计算的时间,其执行过程如图1 所示。其中进程的分配、结果的收集均通过消息传递来完成。在MPI 中,通过指定通信因子和组来对进程进行一种逻辑上的划分,通讯因子定义了进程组内或组间通讯的上下文。
图1 消息传递模型的执行过程
消息传递的多处理机编程可以通过下面的方法实现:
(1) 设计全新的并行编程模型;
(2) 扩展现有串行高级编程语言的语法和保留字来实现消息传递;
(3) 使用现有的高级串行编程语言并提供消息传递的扩展链接库。
MPI 是运用方法(3)来实现的。采用MPI 构建一个并行计算环境主要用到6个函数。这些函数的函数名MPI_Init,MPI_Comm_Size,MPI_Comm_Rank,MPI_Send,
MPI_Recv,MPI_Finalize。MPI_ Init 和MPI_Finalize实现MPI 环境的初始化和结束。MPI_COMMWORLD通信因子是在MPI环境初始化过程中创建的。MPI_Comm_Size获取缺省组( Group )的大小; MPI_Comm_Rank获取本进程〔调用进程〕在缺省组中的逻辑号(从0 开始) ;MPI_Send和MPI_Recv实现消息的传递。一个基本的MPI程序的流程图如图2 所示。
图2 MPI 程序设计流程图
三、 PSRS算法原理
并行正则采样排序,简称PSRS(Parallel Sorting by Regular Sampling)它是一种基于均匀划分(Uniform Partition)原理的负载均衡的并行排序算法。假定待排序的元素有n个,系统中有p个处理器,那么系统首先将n个元素均匀地分割成p段,每段含有n/p个元素,每段指派一个处理器,然后各个处理器同时施行局部排序。为了使各段中诸局部有序的元素在整个序列中也能占据正确的位置,那么就首先从各段中抽取几个代表元素,再从他们产生出p-1个主元,然后按这些主元与原各局部有序中的元素之间的偏序关系,将各个局部有序段划分成p段,接着通过全局交换将各个段中的对应部分集合在一起,最后将这些集合在一起的各部分采用多路归并的方法进行排序,这些有序段汇合起来就自然成为全局有序序列了。
设有n个数据,P个处理器,以及均匀分布在P个处理器上的n个数据。 则PSRS算法可描述如下:
算法1 PSRS排序算法
输入:n个待排序的数据,均匀地分布在P个处理器上
输出:分布在各个处理器上,得到全局有序的数据序列
Begin
(1) 每个处理器将自己的n/P个数据用串行快速排序(Quicksort),得到一个排好序的序列;
(2) 每个处理器从排好序的序列中选取第w,2w,3w,…,(P-1)w个共P-1个数据作为代表元素,其中w=n/P2;
(3) 每个处理器将选好的代表元素送到处理器P0中,并将送来的P段有序的数据序列做P路归并,再选择排序后的第P-1,2(P-1),…,(P-1)(P-1)个共P-1个主元;
(4) 处理器P0将这P-1个主元播送到所有处理器中;
(5) 每个处理器根据上步送来的P-1个主元把自己的n/P个数据分成P段,记 为处理器Pi的第j+1段,其中i=0,…,P-1,j=0,…,P-1;
(6) 每个处理器送它的第i+1段给处理器Pi,从而使得第i个处理器含有所有处理器的第i段数据(i=0,…,P-1);
(7) 每个处理器再通过P路归并排序将上一步的到的数据排序;从而这n个数据便是有序的。
End
在该算法中,针对其中的计算部分(1)中的快速排序的代价为O(klogk),其中k=n/p;第(2)步中各个处理器选择P个主元的代价为O(P);(3)中对P2个主元进行归并并选取新的主元所需代价为O(P2+logP);(5)中对数据划分的代价为O(P+n/p);最后第(7)步中各个处理器进行并行归并的代价为O(n/p+logP)。针对通信部分⑶中处理器P0收集P2个主元的代价为O(P2);(4)中播送新选择的P个主元的代价为O(P);最后第(6)步的多对多播送的通信复杂度与具体的算法实现相关,其最大不会超过O(n)。(MPI源程序请附录)
四、 结论
MPI 提供了良好的并行程序接口,通过调用MPI的库程序,可以达到并行化程序的目的,但并行程序的设计相对于串行程序的设计来说要复杂得多。本文主要介绍和分析了PSRS的并行排序算法和用MPI实现了PSRS的并行排序算法。在此,只给出简单用MPI实现了PSRS的并行排序算法思想和实例,希望能对从事并行算法的工作者有所启发。
参考文献
[1]. 陈国良 编著. 并行算法的设计与分析(修订版). 高等教育出版社, 2002.11
[2]. 黄凯,徐志伟. 可扩展并行计算——技术、结构与编程[M] 机械工业出版社,2000
[3]. Ananth Grama,Anshul Gupta,George Karypis,Vipin Kumar著 并行计算导论[M] 机械工业出版社,2005.01
[4]. Khalid Al - Tawil , Csaba Andras Moritz. Performance Modeling and Evaluation of MPI[J ] . Journal of Parallel and Distributed Computing , 2001 ,61 (2) :202 - 223.
[5]. Dale Shires ,Ram Mohan. An Evaluation of HPF and MPI Approaches and Performance in Unstructured Finite Element Simulations[J] . Journal of Mathematical Modelling and Algorithms ,2001 ,1(3) : 153 - 167.
[6]. Shi H, Schaeffer J. Parallel Sorting by Regular Sampling. Journal of Parallel and Distributed Computing, 14(4), 1992
[7]. Li X, Lu P, Schaeffer J, Shillington J, Wong P S, Shi H. On the Versatility of Parallel Sorting by Regular Sampling. Parallel Computing, 19, 1993
附录
#define MPICH_SKIP_MPICXX 1
#include <stdio.h>
#include <stdlib.h>
#include "mpi.h"
#define INIT_TYPE 10
#define ALLTOONE_TYPE 100
#define ONETOALL_TYPE 200
#define MULTI_TYPE 300
#define RESULT_TYPE 400
#define RESULT_LEN 500
#define MULTI_LEN 600
int Spt;
long DataSize;
int *arr,*arr1;
int mylength;
int *index;
int *temp1;
Psrs_Main( );
/*串行多路归并算法*/
multimerge(int *data1,int * ind ,int *data,int *iter,int SumID);
/*串行快速排序算法*/
quicksort(int *datas,int bb,int ee);
merge(int *data1,int s1,int s2,int *data2);
/*输出错误信息*/
void merror(char* ch);
main(int argc,char* argv[])
{
long BaseNum = 1;
int PlusNum;
int MyID;
MPI_Init(&argc,&argv);
MPI_Comm_rank(MPI_COMM_WORLD,&MyID);
PlusNum=60;
DataSize = BaseNum*PlusNum;
if (MyID==0)
printf("The DataSize is : %lu/n",PlusNum);
Psrs_Main();
if (MyID==0)
printf("/n");
MPI_Finalize();
return 0;
}
Psrs_Main( )
{
int i,j;
int MyID,SumID;
int n,c1,c2,c3,c4,k,l;
//FILE * fp;
//int ready;
MPI_Status status[32*32*2];
MPI_Request request[32*32*2];
MPI_Comm_rank(MPI_COMM_WORLD,
&MyID);
MPI_Comm_size(MPI_COMM_WORLD,
&SumID);
Spt = SumID-1;
/*初始化参数*/
arr = (int *)malloc(2*DataSize*sizeof(int));
if (arr==0)
merror("malloc memory for arr error!");
arr1 = &arr[DataSize];
if (SumID>1)
{
temp1 = (int *)malloc(sizeof(int)*
SumID*Spt);
if (temp1==0) merror("malloc memory for temp1 error!");
index = (int *)malloc(sizeof(int)*
2*SumID);
if (index==0) merror("malloc memory for index error!");
}
MPI_Barrier( MPI_COMM_WORLD);
mylength = DataSize / SumID;
srand(MyID);
printf("This is node %d /n",MyID);
printf("On node %d the input data is:/n",MyID);
for (i=0;i<mylength;i++)
{
arr[i] = (int)rand();
printf("%d : ",arr[i]);
}
printf("/n");
/*每个处理器将自己的n/P个数据用串行快速
排序(Quicksort),得到一个排好序的序
列,对应于算法1步骤(1)*/
MPI_Barrier( MPI_COMM_WORLD);
quicksort(arr,0,mylength - 1);
MPI_Barrier( MPI_COMM_WORLD);
/*每个处理器从排好序的序列中选取第w,2w,
3w,…,(P-1)w个共P-1个数据作为代
表元素,其中w=n/P*P,对应于算法1
步骤(2)*/
if (SumID>1)
{
MPI_Barrier(MPI_COMM_WORLD);
n = (int)(mylength/(Spt+1));
for (i=0;i<Spt;i++)
temp1[i] = arr[(i+1)*n-1];
MPI_Barrier(MPI_COMM_WORLD);
if (MyID==0)
{
/*每个处理器将选好的代表元素送
到处理器P0中,对应于算法
1步骤(3) */
j = 0;
for (i=1;i<SumID;i++)
MPI_Irecv(&temp1[i*Spt],
sizeof(int)*Spt,
MPI_CHAR, i,ALLTOONE_TYPE+i, MPI_COMM_WORLD, &request[j++]);
MPI_Waitall(SumID-1, request,
status);
/* 处理器P0将上一步送来的P段
有序的数据序列做P路归并,
再选择排序后的第P-1,2(P-1),…,(P-1)(P-1)个共P-1个主元,,对应于算法1步骤(3)*/
MPI_Barrier(
MPI_COMM_WORLD);
quicksort(temp1,0,SumID*Spt-1);
MPI_Barrier(
MPI_COMM_WORLD);
for (i=1;i<Spt+1;i++)
temp1[i] = temp1[i*Spt-1];
/*处理器P0将这P-1个主元播送到
所有处理器中,对应于算法1步骤(4)*/
MPI_Bcast(temp1,
sizeof(int)*(1+Spt), MPI_CHAR, 0,
MPI_COMM_WORLD);
MPI_Barrier(
MPI_COMM_WORLD);
}
else
{ MPI_Send(temp1,sizeof(int)*Spt,
MPI_CHAR,0,
ALLTOONE_TYPE+MyID, MPI_COMM_WORLD);
MPI_Barrier(
MPI_COMM_WORLD);
MPI_Barrier(
MPI_COMM_WORLD);
MPI_Bcast(temp1,
sizeof(int)*(1+Spt),
MPI_CHAR, 0,
MPI_COMM_WORLD);
MPI_Barrier(
MPI_COMM_WORLD);
}
/*每个处理器根据上步送来的P-1个主
元把自己的n/P个数据分成P段,
记为处理器Pi的第j+1段,其中i=0,…,P-1,j=0,…,P-1,对应于算法1步骤(5)*/
n = mylength;
index[0] = 0;
i = 1;
while ((arr[0]>=temp1[i])&&(i<SumID))
{
index[2*i-1] = 0;
index[2*i] = 0;
i++;
}
if (i==SumID)
index[2*i-1] = n;
c1 = 0;
while (i<SumID)
{
c4 = temp1[i];
c3 = n;
c2 = (int)((c1+c3)/2);
while ((arr[c2]!=c4)&&(c1<c3))
{
if (arr[c2]>c4)
{
c3 = c2-1;
c2 = (int)((c1+c3)/2);
}
else
{
c1 = c2+1;
c2 = (int)((c1+c3)/2);
}
}
while ((arr[c2]<=c4)&&(c2<n))
c2++;
if (c2==n)
{
index[2*i-1] = n;
for (k=i;k<SumID;k++)
{
index[2*k] = 0;
index[2*k+1] = 0;
}
i = SumID;
}
else
{
index[2*i] = c2;
index[2*i-1] = c2;
}
c1 = c2;
c2 = (int)((c1+c3)/2);
i++;
}
if (i==SumID)
index[2*i-1] = n;
MPI_Barrier( MPI_COMM_WORLD);
/*每个处理器送它的第i+1段给处理器
Pi,从而使得第i个处理器含有所
有处理器的第i段数据
(i=0,…,P-1),,对应于算法1步
骤(6)*/
j = 0;
for (i=0;i<SumID;i++)
{
if (i==MyID)
{
temp1[i] = index[2*i+1]-
index[2*i];
for (n=0;n<SumID;n++)
if (n!=MyID)
{
k = index[2*n+1]-
index[2*n];
MPI_Send(&k,
sizeof(int), MPI_CHAR, n, MULTI_LEN+MyID,
MPI_COMM_WORLD);
}
}
else
{
MPI_Recv(&temp1[i],
sizeof(int), MPI_CHAR, i,MULTI_LEN+i,
MPI_COMM_WORLD,
&status[j++]);
}
}
MPI_Barrier(MPI_COMM_WORLD);
j = 0;
k = 0;
l = 0;
for (i=0;i<SumID;i++)
{
MPI_Barrier(
MPI_COMM_WORLD);
if (i==MyID)
{
for (n=index[2*i];
n<index[2*i+1];
n++)
arr1[k++] = arr[n];
}
MPI_Barrier(
MPI_COMM_WORLD);
if (i==MyID)
{
for (n=0;n<SumID;n++)
if (n!=MyID)
{ MPI_Send(&arr[
index[2*n]], sizeof(int)*
(index[2*n+1]-index[2*n]),MPI_CHAR, n, MULTI_TYPE+MyID, MPI_COMM_WORLD);
}
}
else
{
l = temp1[i]; MPI_Recv(&arr1[k],
l*sizeof(int),
MPI_CHAR,
i ,MULTI_TYPE+i,
MPI_COMM_WORLD,
&status[j++]);
k=k+l;
}
MPI_Barrier(
MPI_COMM_WORLD);
}
mylength = k;
MPI_Barrier(MPI_COMM_WORLD);
/*每个处理器再通过P路归并排序将上
一步的到的数据排序;从而这n个
数据便是有序的,,对应于算法1
步骤(7) */
k = 0;
multimerge(arr1,temp1,arr,&k,SumID);
MPI_Barrier(MPI_COMM_WORLD);
}
printf("On node %d the sorted data is :/n",MyID);
for (i=0;i<mylength;i++)
printf("%d : ",arr[i]);
printf("/n");
return 0;
}
/*输出错误信息*/
void merror(char* ch)
{
printf("%s/n",ch);
exit(1);
}
/*串行快速排序算法*/
quicksort(int *datas,int bb,int ee)
{
int tt,i,j;
tt = datas[bb];
i= bb;
j = ee;
if (i<j)
{
while(i<j)
{
while ((i<j)&&(tt<=datas[j]))
j--;
if (i<j)
{
datas[i] = datas[j];
i++;
while ((i<j)&&
(tt>datas[i]))
i++;
if (i<j)
{
datas[j] = datas[i];
j--;
if (i==j)
datas[i] = tt;
}
else
datas[j] = tt;
}
else
datas[i] = tt;
}
quicksort(datas,bb,i-1);
quicksort(datas,i+1,ee);
}
return 0;
}
/*串行多路归并算法*/
multimerge(int *data1,int * ind ,int *data,int *iter,int
SumID)
{
int i,j,n;
j = 0;
for (i=0;i<SumID;i++)
if ( ind [i]>0)
{
ind [j++] = ind [i];
if (j<i+1) ind [i] = 0;
}
if ( j>1 )
{
n = 0;
for (i =0;i<j,i+1<j;i=i+2)
{ merge(&(data1[n]), ind [i],
ind [i+1],&(data[n]));
ind [i] += ind [i+1];
ind [i+1] = 0;
n += ind [i];
}
if (j%2==1 )
for (i=0;i< ind [j-1];i++)
data[n++]=data1[n];
(*iter)++;
multimerge(data, ind ,data1,iter,
SumID);
}
return 0;
}
merge(int *data1,int s1,int s2,int *data2)
{
int i,l,m;
l = 0;
m = s1;
for (i=0;i<s1+s2;i++)
{
if (l==s1)
data2[i]=data1[m++];
else
if (m==s2+s1)
data2[i]=data1[l++];
else
if (data1[l]>data1[m])
data2[i]=data1[m++];
else
data2[i]=data1[l++];
}
return 0;
}