进程ai内存不足

2003年高性能计算培训班材料迟学斌中国科学院计算机网络信息中心张林波中国科学院数学与系统科学研究院莫则尧北京应用物理与计算数学研究所2003年8月28日本文档及程序示例可从http://lssc.
cc.
ac.
cn/training2003下载目录第一部分并行计算基础1第一章并行计算基础知识3第二章并行计算机体系结构23第三章Linux机群37§3.
1引言37§3.
2构建Linux机群的要素37§3.
3几种典型的Linux机群结构38§3.
3.
1单台微机38§3.
3.
2由几台日常使用的微机构成的机群38§3.
3.
3专用并行机群38§3.
4在单机上安装、配置MPI并行环境39§3.
4.
1Linux的安装39§3.
4.
2MPICH的安装39§3.
5在联网的多台机器上安装、配置MPI并行环境40§3.
5.
1设置NFS41§3.
5.
2设置NIS42§3.
5.
3设置rsh43§3.
5.
4MPICH的安装43§3.
5.
5MPICH程序的编译、运行43§3.
6专用并行机群系统44§3.
7组建大型机群系统需要考虑的一些问题45第四章矩阵并行计算47§4.
1矩阵相乘的若干并行算法47§4.
1.
1行列划分算法47§4.
1.
2行行划分算法48§4.
1.
3列列划分算法48§4.
1.
4列行划分算法49§4.
1.
5Cannon算法49§4.
2线性方程组的解法50§4.
2.
1分布式系统的并行LU分解算法51§4.
2.
2具有共享存储系统的并行LU分解算法52§4.
2.
3三角方程组的并行解法53§4.
3对称正定线性方程组的并行解法55§4.
3.
1Cholesky分解列格式的并行计算55§4.
3.
2双曲变换Cholesky分解56iii目录§4.
3.
3修正的双曲变换Cholesky分解58§4.
4三对角方程组的并行解法59§4.
4.
1递推法59§4.
4.
2分裂法61§4.
5异步并行迭代法61§4.
5.
1异步并行迭代法基础62§4.
5.
2线性迭代的一般收敛性结果62第二部分MPI并行程序设计65第五章消息传递并行程序设计平台MPI67§5.
1MPI并行环境管理函数67§5.
2进程控制函数68§5.
3MPI进程组操作函数68§5.
4MPI通信子操作71§5.
5点到点通信函数73§5.
6阻塞式通信函数73§5.
7非阻塞式通信函数78§5.
8特殊的点到点通信函数82§5.
9MPI的通信模式84§5.
10用户定义的数据类型与打包85§5.
11用户定义的数据类型85§5.
12MPI的数据打包与拆包90§5.
13聚合通信93§5.
14障碍同步93§5.
15单点与多点通信函数93§5.
16多点与多点通信函数97§5.
17全局归约操作100§5.
18进程拓扑结构106§5.
18.
1迪卡尔拓扑结构106§5.
18.
2一般拓扑结构109§5.
18.
3底层支持函数110第六章文件输入输出(MPI–IO)113§6.
1基本术语113§6.
2基本文件操作114§6.
2.
1打开MPI文件114§6.
2.
2关闭MPI文件115§6.
2.
3删除文件115§6.
2.
4设定文件长度115§6.
2.
5为文件预留空间115目录iii§6.
2.
6查询文件长度116§6.
3查询文件参数116§6.
3.
1查询打开文件的进程组116§6.
3.
2查询文件访问模式116§6.
4设定文件视窗116§6.
4.
1文件中的数据表示格式117§6.
4.
2可移植数据类型117§6.
4.
3查询数据类型相应于文件数据表示格式的域118§6.
5文件读写操作118§6.
5.
1使用显式位移的阻塞型文件读写119§6.
5.
2使用独立文件指针的阻塞型文件读写119§6.
5.
3使用共享文件指针的阻塞型文件读写120§6.
5.
4非阻塞型文件读写函数120§6.
5.
5分裂型文件读写函数121§6.
6文件指针操作121§6.
6.
1独立文件指针操作121§6.
6.
2共享文件指针操作122§6.
6.
3文件位移在文件中的绝对地址122§6.
7不同进程对同一文件读写操作的相容性123§6.
7.
1设定文件访问的原子性123§6.
7.
2查询atomicity的当前值124§6.
7.
3文件读写与存储设备间的同步124§6.
8子数组数据类型创建函数124第七章MPI程序示例127§7.
1矩阵乘积127§7.
1.
1算法描述127§7.
1.
2MPI并行程序127§7.
1.
3MPI并行程序的改进131§7.
2Poisson方程求解133§7.
2.
1并行算法设计133§7.
2.
2MPI并行程序设计134§7.
2.
3MPI并行程序的改进138参考文献143附录145§A.
1Linpack性能测试145§A.
1.
1相关链接145§A.
1.
2HPL简介145§A.
1.
3适用于Linux的BLAS库145§A.
1.
4HPL的编译、运行与性能优化146iv目录第一部分并行计算基础第一章并行计算基础知识参见1.
3第二章并行计算机体系结构参见2.
2336第二章并行计算机体系结构第三章Linux机群本章介绍在微机环境中配置、安装MPI并行环境.
为方便叙述,这里仅给出最简单的配置方法而暂不考虑其它问题,如网络性能、网络安全等等.
本章叙述的方法以RedHat–9Linux为例,如果使用其它Linux系统在一些细节上会略有差异.
§3.
1引言Linux机群系统己成为最流行的高性能计算平台,在高性能计算机中占有越来越大的比重系统规模可从单机、少数几台联网的微机直到包括上千个结点的大规模并行系统既可作为廉价的并行程序调试环境,也可设计成真正的高性能并行机普及并行计算必不可少的工具用于高性能计算的机群系统在结构上、使用的软件工具上通常有别于用于提供网络、数据库服务的机群(后者亦称为服务器集群)TOP500SupercomputerSitesCluster@TOP500参考资料:用关键字"clusterhowto"在网上搜索相关材料.
.
.
§3.
2构建Linux机群的要素单台或联网的多台微机或服务器Linux系统:RedHat,Debian,SuSE,Mandrake,.
.
.
(可选)高速内联网络:千兆以太网,Myrinet,QsNet,DolphinSCI,Inniband,.
.
.
编译系统:gcc/g++/g77,PGI,Intel,.
.
.
MPI系统:MPICH,LAM-MPI,.
.
.
网络文件系统:NFS,PVFS,Lustre.
.
.
资源管理与作业调度:PBS,Condor,LSF,.
.
.
数学库:BLASMKL,ATLAS,KazushigeGoto'sBLAS(推荐)FFTWhttp://www.
fftw.
org/LAPACKhttp://www.
netlib.
org/lapack/ScaLAPACKhttp://www.
netlib.
org/scalapack/.
.
.
其它工具:PETSchttp://www-unix.
mcs.
anl.
gov/petsc/petsc-2/UGhttp://cox.
iwr.
uni-heidelberg.
de/~ug/.
.
.
3738第三章Linux机群§3.
3几种典型的Linux机群结构§3.
3.
1单台微机可以有一个或多个处理器安装Linux系统,C、Fortran编译器,以及MPICH或LAM-MPI模拟并行:一个处理器上运行多个进程或线程真实并行:多个处理器上运行多个进程或线程通过共享内存或TCP/IP进行通信§3.
3.
2由几台日常使用的微机构成的机群通常连接在同一个局域网,通信通过网络进行安装Linux,编译环境,MPICH或LAM-MPI为方便使用,最好设立一个共享的目录(NFS),以及一个NIS或LDAP服务器§3.
3.
3专用并行机群指专门建造的用于并行计算的机群.
通常单独形成一个局域网,再通过一个网关连接到Internet通常使用内部IP地址(如192.
168.
0.
x),对外部而言只有网关是可见的可根据需要设立一至数台服务器,分别承担网关、时钟同步(NTP)、NIS/LDAP、网络文件系统(NFS)、资源管理、用户登录、作业调度等服务通过IP伪装(IPMasquerading)或网络地址转换(NAT,NetworkAddressTranslation)使得内部结点能够直接访问Internet(ipchains或iptables)利用LVS(LinuxVirtualServer)将外部用户分配到登录结点大型机群通常采用安装在机柜中的机架式或刀片式服务器,并有专门配备的UPS及空调几个机群实例:CRack,LSSC–I,LSSC–II图3.
1:典型Linux机群结构示意图§3.
4在单机上安装、配置MPI并行环境39§3.
4在单机上安装、配置MPI并行环境§3.
4.
1Linux的安装可以安装任何Linux发布版,推荐RedHat-9一些必须安装的包(关于如何用RedHat的rpm命令安装软件包可参看MaximumRPM)gcc包GNUC,用于C程序的编译gcc-g77包GNUFortran77,用于Fortran程序的编译rsh包和rsh-server包用于启动MPI进程配置rsh-server,下面的操作必须以root身份执行.
–编辑文件/etc/hosts.
equiv,在其中加上本机主机名(单独占一行).
如果该文件不存在则创建一个.
–开启rsh服务:/sbin/chkconfigrshon–测试rshd的配置.
以普通用户(非root)登录并运行命令:rsh主机名/bin/hostname如果配置正确该命令应该显示本机主机名.
如果出错可查看/var/log/messages文件中的错误信息.
–注意保证shell初始化文件(.
cshrc,.
profile,.
bashrc等)不要往stdout和stderr输出任何信息,即上述命令的输出除了主机名外不应该有任何其它内容.
否则可能影响MPI进程的启动.
§3.
4.
2MPICH的安装MPICH的手册在源程序的doc目录中.
从网址http://www-unix.
mcs.
anl.
gov/mpi/mpich/处下载MPICH最新版本的源程序,文件名为mpich-1.
x.
x.
tar.
gz,其中1.
x.
x为MPICH的版本号(目前的最新版本为1.
2.
5).
在下面描述的过程中需用MPICH的实际版本号替换1.
x.
x展开MPICH源程序:tarxzpvfmpich-1.
x.
x.
tar.
gz配置、编译MPICH:cdmpich-1.
x.
x.
/configure-prefix=/usr/local/mpi--disable-weak-symbolsmake上述命令生成的MPICH库应该使用ch_p4进行底层通信.
选项--disable-weak-symbols的使用是为避免MPICH1.
2.
5的一个bug,否则产生的库中将没有MPI_File_xxxx等函数(只有PMPI_File_xxx等函数)安装MPICH:makeinstall这一步必须以root身份执行,它将MPICH的文件安装到目录/usr/local/mpi中.
安装完毕后可以删除MPICH源文件目录.
分别将"/usr/local/mpi/bin"和"/usr/local/mpi/man"加入到环境变量PATH和MANPATH中.
只需在目录/etc/profile.
d中创建两个文件mpich.
sh和mpich.
csh,它们分别对Bourneshell40第三章Linux机群和Cshell起作用,这两个文件的内容如下:–/etc/profile.
d/mpich.
sh:#!
/bin/bashexportMANPATH=${MANPATH}:/usr/local/mpi/manexportPATH=${PATH}:/usr/local/mpi/bin–/etc/profile.
d/mpich.
csh:#!
/bin/cshif($MANPATH==0)thensetenvMANPATH:/usr/local/mpi/manelsesetenvMANPATH${MANPATH}:/usr/local/mpi/manendifsetenvPATH${PATH}:/usr/local/mpi/bin测试MPICH:C程序cp/usr/local/mpi/examples/cpi.
c.
mpicc-ocpicpi.
c.
mpirun-np1cpimpirun-np2cpiFortran程序cp/usr/local/mpi/examples/pi3.
f.
mpif77-opi3pi3.
f.
mpirun-np1pi3mpirun-np2pi3(程序运行过程中依次输入10000和0).
§3.
5在联网的多台机器上安装、配置MPI并行环境Linux系统的安装要求与前一节单机环境一样.
此外,在开始下面的步骤之前还应该先配置好TCP/IP网络连接.
为避免额外麻烦,在配置网络时请不要开启任何防火墙设置.
此外,为方便相互访问,将所有机器的主机名放在/etc/hosts文件中.
可以在所有机器上使用同样的/etc/hosts文件,它包含如下形式的内容:127.
0.
0.
1localhost.
localdomainlocalhost10.
10.
10.
1node1.
mydomainnode110.
10.
10.
2node2.
mydomainnode210.
10.
10.
nnoden.
mydomainnoden(按真实情况替换其中的主机名和IP地址).
下面介绍的方法使用NIS(NetworkInformationService,也叫做SUNYellowPages)管理用户帐号,使用NFS(NetworkFileSystem)共享用户目录.
首先根据情况选择一台机器作为NIS和NFS服务器,我们将该机器称为服务结点或主结点机,而其它机器称为从结点机.
主结点机和从结点机的配置方式不同,下面将分别介绍.

§3.
5在联网的多台机器上安装、配置MPI并行环境41§3.
5.
1设置NFS主结点机创建目录:mkdir-p/home/local将/usr/local链接到/home/local:/bin/rm-rf/usr/localln-s/home/local/usr/.
(注:如果以前的/usr/local中安装有有用的文件,则在执行上述命令前应该将其拷贝或移走).
确认在主结点机上安装了nfs-utils包.
开启NFS服务:/sbin/chkconfignfson/sbin/chkconfignfslockon/etc/init.
d/nfslockrestart/etc/init.
d/nfsrestart在文件/etc/exports中加入下面一行:/home*(rw,no_root_squash)它将/home目录输出给所有机器.
(注:出于安全考虑可以限制仅将目录输出给指定的结点机,及将no_root_squash改成root_squash.
请用命令man5exports查看相关参数的写法).
输出指定目录(/home):exportfs-a(也可重新启动系统).
从结点机创建目录:mkdir-p/home在文件/etc/fstab中加入下面一行::/home/homenfsdefaults00(将替换成主结点机的主机名或IP地址).
运行命令:/sbin/chkconfignetfson它使得系统启动时自动挂接主结点机上的/home目录.
运行命令:mount/home(也可重新启动系统).
将/usr/local链接到/home/local:/bin/rm-rf/usr/localln-s/home/local/usr/.
(注:如果以前的/usr/local中有有用的文件,则在执行上述命令前应该将其拷贝或移走).
上述所有操作必须以root身份执行.
本步骤完成后,所有结点机上的/home和/usr/local目录的内容应该是一样的.
在从结点机上可用命令df检查挂接情况.
有关NFS的进一步材料可在网上搜索"NFS"等关键字.
42第三章Linux机群§3.
5.
2设置NIS以下描述中假设以'CLUSTER'作为NIS域名.
主结点机确认安装了下述包:ypserv,ypbind,yp-tools在文件/etc/sysconfig/network中加入下面一行:NISDOMAIN=CLUSTER开启NIS服务:/sbin/chkconfigypservon/etc/init.
d/ypservstart初始化NIS数据库:/usr/lib/yp/ypinit-m程序运行时按Ctrl-D,然后按"y"和回车.
该命令将生成NIS数据库.
可以忽略Noruletomaketarget.
.
.
之类的错误信息.
开启NIS客户程序:/sbin/chkconfigypbindon/etc/init.
d/ypbindstart验证NIS设置–命令"ypwhich"应该显示出主结点机主机名.
–命令"ypcatpasswd"应该显示出(主结点机上的)用户帐号.
从结点机确认安装了下述包:ypbind,yp-tools在文件/etc/sysconfig/network中加入下面一行:NISDOMAIN=CLUSTER开启NIS客户程序:/sbin/chkconfigypbindon/etc/init.
d/ypbindstart验证NIS设置–命令"ypwhich"应该显示出主结点机主机名.
–命令"ypcatpasswd"应该显示出主结点机上的用户帐号.
为了能够使用NIS用户登录,还需要修改/etc/nsswitch.
conf文件,使其包含下述设置:passwd:filesnisshadow:filesnisgroup:filesnishosts:filesnisdns注意在从结点机上不要重复定义用户信息,即在从结点机上应该将所有NIS用户(可用命令"ypcatpasswd"显示出来)从文件/etc/passwd和/etc/shadow中删除,将所有NIS组(可用命令"ypcatgroup"显示出来)从文件/etc/group中删除.
完成NIS配置后,创建新的用户帐号只需在主结点机上进行(注意将用户的主目录放到/home§3.
5在联网的多台机器上安装、配置MPI并行环境43下),然后运行命令"cd/var/yp;make"即可.
如果在主结点机上修改了一个用户帐号信息,也应该运行一次上述命令以刷新NIS数据库.
NIS用户在从结点机上不能用"passwd"命令修改用户口令,而必须用"yppasswd"命令来修改.
上述所有操作必须以root身份执行.
有关NIS的进一步材料可在网上搜索"NIS"等关键字.
§3.
5.
3设置rsh确认安装了3.
4.
1中列出的包.
将所主机名加到文件/etc/hosts.
equiv中.
开启rsh服务:/sbin/chkconfigrshon注:为使root用户能够用rsh执行远程命令,需将/etc/hosts.
equiv文件拷贝到/root/.
rhosts,并在文件/etc/securetty中加入"rsh".
上述操作必须在所有结点机上以root身份进行.
完成上述设置后,在任何一台结点机上应该可以在所有结点机(包括自己)上执行远程命令.
可按如下方法测试:以一个NIS用户登录到一个结点机上并运行命令:rsh另一台结点机名/bin/hostname如果配置正确该命令应该显示对方主机名.
如果出错可查看对方/var/log/messages文件中的错误信息.
注意保证shell初始化文件(.
cshrc,.
profile,.
bashrc等)不要往stdout和stderr输出任何信息,即上述命令的输出除了主机名外不应该有任何其它内容.
否则可能影响到MPI进程的启动.
§3.
5.
4MPICH的安装MPICH的安装与单机下的安装完全一样,只需在主结点机上进行,因为所有结点机/usr/local目录是共享的.
此外,需要将文件/etc/profile.
d/mpich.
sh/etc/profile.
d/mpich.
csh拷贝到所有结点机上.
§3.
5.
5MPICH程序的编译、运行MPICH程序的编译可以在任何一台结点机上用mpicc(C),mpif77(Fortran),mpiCC(C++)等命令进行.
它们是MPICH提供的shell脚本,用法与普通的C/C++/Fortran编译器一样.
MPICH程序的运行方式取决于编译MPICH系统时选择的底层driver.
这里介绍的编译方式使用ch_p4作为底层driver,这种情况下有两种方法选择运行一个MPI程序使用的结点机和进程数,即:mpirun-machinefile文件名-np4MPI程序名[MPI程序参数]文件文件名中列出希望使用的结点机名,一行一个.
mpirun将在给定的结点机上启动指定数目的进程(这里是4).
当进程数目大于结点机数目时mpirun会在一些结点机上启动两个或更多进程.
命令"mpirun-help"可以显示一个mpirun的简要使用说明.
.
/MPI程序名-p4pg文件名[MPI程序参数]这种方式可以精确控制在每台结点机上启动的MPI进程数与进程序号,并且允许在不同结点机44第三章Linux机群上启动不同的可执行文件(适用于Master/Slave模式的并行程序).
文件文件名中按下列格式列出各结点机上启动的程序名:结点机名10可执行文件名1结点机名21可执行文件名2结点机名31可执行文件名3.
.
.
.
.
.
结点机名n1可执行文件名n其中结点机名1必须是运行该命令时所在的结点机,可执行文件名1必须与命令行上的MPI程序名为同一文件.
所有可执行文件名必须使用绝对路径(如/home/zlb/test/cpi).
通常情况下,所有可执行文件名是一样的.
而当同一结点机名出现多次时表示在该结点机上启动多个进程.
例如,假设用户在结点机node1的/home/zlb/test目录下,该目录中有已经编译好的MPI程序cpi.
在该目录下创建一个名为p4file的文件,它包含如下内容:node10/home/zlb/test/cpinode21/home/zlb/test/cpinode11/home/zlb/test/cpinode21/home/zlb/test/cpi则命令".
/cpi-p4pgp4file"将在node1,node2上运行四个进程,其中进程0和进程2在node1上,进程1和进程3在node2上.
§3.
6专用并行机群系统指专门用于并行计算的机群系统.
这类系统的特点是结点机针对机群的用途经过专门选型,结点机配置比较统一.
结点机按照所负责的功能分成几类,如I/O结点、计算结点、登录结点、服务结点等,同一类结点机通常具有相同的硬件配置和文件系统结构一般单独构成一个局域网,使用内部网络地址配备资源调度和作业管理系统系统规模可以从数个结点到上千个结点图3.
2中给出一个包含5个结点的机群示例,它有1个主控结点,4个计算结点.
图3.
2:小型专用机群示例以图3.
2的机群为例,假设4台计算结点的硬件配置完全一样,则它们安装的系统文件也基本一样,唯一区别是主机名和IP地址.
在RedHatLinux中,主机名由文件§3.
7组建大型机群系统需要考虑的一些问题45/etc/sysconfig/network设定,而假设以太网卡设备名为eth0的话则IP地址由文件/etc/sysconfig/network-scripts/ifcfg-eth0设定,这是计算结点文件系统间的唯一区别.
为方便操作系统的安装与管理,我们希望所有计算结点上这两个文件也完全一样.
为达到这一目的,一个办法是利用DHCP协议从主控结点获取IP地址,具体作法可参看,例如:http://www.
linux.
org/docs/ldp/howto/mini/DHCP/我们在这里介绍另一个简单方法,它将所有MAC地址放在文件/etc/hosts中,计算结点启动时通过MAC地址得出自己的主机名和IP地址.
下面是具体做法.
/etc/hosts文件:127.
0.
0.
1localhost.
localdomainlocalhost192.
168.
0.
254node0node0.
cluster00:90:27:57:42:64192.
168.
0.
1node1node1.
cluster00:90:27:57:45:F2192.
168.
0.
2node2node2.
cluster00:90:27:57:42:B4192.
168.
0.
3node3node3.
cluster00:90:27:57:45:E1192.
168.
0.
4node4node4.
cluster00:90:27:57:3F:7E(最后一列应换成真实MAC地址,MAC地址可用命令/sbin/ifconfigeth0获得)/etc/sysconfig/network文件:NETWORKING=yeshwaddr='/sbin/ifconfigeth0|/bin/grep-iHWaddr|awk'{print$5}''HOSTNAME='grep'^[0-9]'/etc/hosts|grep-w$hwaddr|awk'{print$2}''unsethwaddrGATEWAY=192.
168.
0.
254NISDOMAIN=CLUSTER/etc/sysconfig/network-scripts/ifcfg-eth0文件:DEVICE=eth0BOOTPROTO=staticIPADDR='grep'^[0-9]'/etc/hosts|grep-w"$HOSTNAME"|awk'{print$1}''BROADCAST=192.
168.
0.
255NETMASK=255.
255.
255.
0NETWORK=192.
168.
0.
0ONBOOT=yes按上述方法配置,只需要安装好一台计算结点,其它计算结点可以简单地通过硬盘拷贝进行复制,也容易利用远程引导(PXEboot)实现自动安装.
§3.
7组建大型机群系统需要考虑的一些问题环境考虑:布局、电源、空调等结点机操作系统自动安装、恢复结点机状态监控:CPU负载、用户作业、内存使用等结点机硬件监控:电压、风扇、温度等(http://secure.
netroedge.
com/~lm78/)MPI系统进程数的限制rsh端口数限制并行程序的快速启动46第三章Linux机群第四章矩阵并行计算在科学与工程计算的许多问题中,矩阵运算,特别是矩阵相乘、求解线性方程组和矩阵特征值问题是最基本的内核.
随着MMP并行计算机和网络并行环境的不断发展,为了充分发挥分布式并行计算机的功效,有必要对计算方法进行深入的研究.
这里着重考虑矩阵乘积和求解线性方程组的多种并行算法,其向量化算法[2,3]本书将不作介绍.
代数特征值问题及其相关问题的计算方法在文献[4]中有较详细的讨论,这里将不作介绍.
为了讨论方便起见,先作一些约定,假设有p个处理机,Pj表示第j个处理机,Pmyid表示当前运行程序的处理机,send(x,j)和recv(x,j)分别表示在Pmyid中把x传送到Pj和从Pj中接收x,本章给出的算法都是在Pmyid处理机上的.
此外,用imodp表示i对p取模运算.
程序设计与机器实现是密不可分的,计算结果的好坏与编程技术有很大的关系,尤其是在并行计算机环境下,研制高质量的程序对发挥计算机的性能起着至关重要的作用.
我们结合算法研究给出在不同并行机上的一些重要例子来表明程序设计思想,并采用程序方式描述算法,以利于机器实现.
在并行算法的研究中,经常要用到加速比的概念,它是评价并行算法好坏的一个重要标准.
关于加速比的多种定义在第三章中已经给出,这里把要用到的再强调一下,用Sp表示p个处理机的加速比,其定义如下:Sp=单机求解问题的执行时间/p个处理机求解问题的执行时间这个定义在应用中经常受到单机解题规模的限制,也就是说,在多处理机上能求解的问题在单机上是做不了的,通常是问题规模增大时,单机内存不够的缘故.
另一个容易改进的加速比的定义为Sp=p个处理机求解问题的运算速度/单机求解问题的运算速度如此定义的加速比对工作量确定的问题是容易得到的,对于这里所考虑的问题如此定义的加速比是合适的,但对工作量不确定,单机又无法求解的问题这两种定义都不适合.
关于通信模式,数据传输方式有许多不同的假设,这样的假设在理论上讨论并行算法的加速比及效率是非常重要的.
对于消息传递型并行系统,通常考虑其为:Tc=α+β*N,其中α是启动时间,β是传输单位数据所需的时间,N是数据的传输量.
§4.
1矩阵相乘的若干并行算法矩阵乘积在实际应用中是经常要用到的,许多先进的计算机上都配有高效的串行程序库.
为了在并行计算环境上实现矩阵乘积,研究并行算法是非常必要的.
本节要考虑的计算问题是C=A*B(4.
1)其中A和B分别是m*k和k*n矩阵,C是m*n矩阵.
不失一般性,假设m=m*p,k=k*p和n=n*p,下面考虑基于对矩阵A和B的不同划分的并行计算方法.
§4.
1.
1行列划分算法这里将矩阵A和B分别划分为如下的行块子矩阵和列块子矩阵:A=AT0AT1.
.
.
ATp1T,B=B0B1.
.
.
Bp1(4.
2)4748第四章矩阵并行计算这时C=(Ci,j)=(Ai*Bj),其中Ci,j是m*n矩阵.
Ai,Bi和Ci,j,j=0,p1存放在Pi中,这种存放方式使数据在处理机中不重复.
在本节所考虑的算法中,都使用这种数据的存放方式.
由于使用p个处理机,每次每台处理机计算出一个Ci,j,计算C需要p次来完成.
Ci,j的计算是按对角线进行的,计算方法如下:算法4.
1.
1.
fori=0top1dol≡i+myidmodpCl=A*B,mp1≡myid+1modp,mm1≡myid-1modpifi=p1,send(B,mm1),recv(B,mp1)end{for}在这个算法中,Cl=Cmyid,l,A=Amyid,B在处理机中每次循环向前移动一个处理机,每次交换数据为k*n矩阵,交换次数为p1次.
如果用DTA1表示在算法4.
1.
1中数据的交换量,CA1表示在算法4.
1.
1中的计算量,则有DTA1=2*k*(nn),CA1=m*k*n/p.
§4.
1.
2行行划分算法这里将矩阵A和B均划分为行块子矩阵,矩阵A的划分同式(4.
2),B的划分如下:B=BT0BT1.
.
.
BTp1T(4.
3)Ci为和Ai相对应的C的第i个块,进一步把Ai按列分块与B的行分块相对应,记Ai=Ai0Ai1.
.
.
Ai,p1从而有Ci=Ai*B=p1j=0Ai,j*Bj初始数据A,B和C的存放方式与4.
1.
1相同,在结点Pmyid上的计算过程可归纳为算法4.
1.
2.
算法4.
1.
2.
fori=0top1dol≡i+myidmodpC=C+Al*B,mp1≡myid+1modp,mm1≡myid-1modpifi=p1,send(B,mm1),recv(B,mp1)end{for}这个算法中的数据交换量和计算量与算法4.
1.
1相同,所不同的只是计算C的方式,其中Al=Amyid,l.
§4.
1.
3列列划分算法这里将矩阵A和B均划分为列块子矩阵,B的划分与式(4.
2)相同,A划分为如下形式:A=A0A1.
.
.
Ap1(4.
4)§4.
1矩阵相乘的若干并行算法49这时C的划分与B的划分相对应,也是按列分块的,进一步把Bi按行分成与A的列分块相对应的行分块,记Bi=BTi0BTi1.
.
.
BTi,p1T,从而有下面计算C的方法.
Cj=A*Bj=p1i=0Ai*Bi,j这时的计算过程是传送矩阵A而不是B,具体的算法描述如下:算法4.
1.
3.
fori=0top1dol≡i+myidmodpC=C+A*Bl,mp1≡myid+1modp,mm1≡myid-1modpifi=p1,send(A,mm1),recv(A,mp1)end{for}算法4.
1.
3的计算量与算法4.
1.
1和算法4.
1.
2是相同的,算法4.
1.
3的数据交换量是DTA3=2*m*(kk).
当m=n时,DTA1=DTA3.
两种算法数据交换量的大小是由m和n决定的,即当m≤n时,DTA3≤DTA1.
由于它们的计算量是相同的,因此只要按通信量大小选择算法就可以得到好的并行效率.
§4.
1.
4列行划分算法这里将矩阵A和B分别划分为列和行块子矩阵,A的划分与式(4.
4)相同,B的划分与式(4.
3)相同.
由此得到C=A*B=p1i=0Ai*BiC的计算是通过p个规模和C相同的矩阵之和得到的,从对问题的划分可以看出,并行算法的关键是计算矩阵的和,设计有效地计算矩阵和的算法,对发挥分布式并行系统的效率起着重要作用.
为了使负载达到平衡,我们给出如下的一种简单计算方法,其中这里的C与其它算法中的C不同,它不是分块的,而是整个存放在所有的处理机中,但在计算C的和时采用与4.
1.
1相同的按行分块.
算法4.
1.
4.
C=A*Bfori=1top1dol≡i+myidmodpk≡myid-imodpsend(Cl,l),recv(T,k)Cmyid=Cmyid+Tend{for}当然在求和时有许多实现方法,这里给出的算法简洁易懂,其通信量DTA4=2*(mm)*n.
如果采用按列分块方式计算C,算法4.
1.
4也同样适合,这时的通信量为DTA4=2*m*(nn),由于m*n和m*n是相等的,因此选择何种方式计算C可根据需要而定.
§4.
1.
5Cannon算法假设矩阵A,B和C可以分成m*m块矩阵,也即,A=(Aij)m*m,B=(Bij)m*m,和C=50第四章矩阵并行计算(Cij)m*m,其中Aij,Bij和Cij是n*n矩阵,进一步假定有p=m*m个处理机.
为了讨论Cannon算法,引入块置换矩阵Q=(Qij)使得Qij=In,j≡i+1modm0n,其它情况其中In和0n分别是n阶单位矩阵和零矩阵.
定义对角块矩阵D(l)A=diag(D(l)i)=diag(Ai,i+lmodm),容易推导出A=m1l=0D(l)A*Ql.
因此C=A*B=m1l=0D(l)A*Ql*B=m1l=0D(l)A*B(l)其中B(l)=Ql*B=Q*B(l1).
利用这个递推关系式,并把处理机结点编号从一维映射到二维,即有Pmyid=Pmyrow,mycol,数据Aij,Bij,和Cij存放在Pij中,容易得到下面的在处理机Pmyid结点上的算法.
算法4.
1.
5.
C=0fori=0tom1dok≡myrow+imodm;mp1≡mycol+1modm;mm1≡mycol-1modm;ifmycol=kthensend(A,(myrow,mp1));copy(A,tmpA);elserecv(tmpA,(myrow,mm1));ifmycol=mp1,send(tmpA,(myrow,mp1));end{if}C=C+tmpA*B;mp1≡myrow+1modm;mm1≡myrow-1modm;ifi=m1thensend(B,(mm1,mycol));recv(B,(mp1,mycol));end{if}end{for}该算法具有很好的负载平衡,其特点是在同一行中广播A,计算出C的部分值之后,在同列中滚动B.
数据交换量DTA5=m*2*n2+(m1)*2*n2=2(2m1)n2=4m2n2/√p2m2n2/p.
由于计算量对每个处理机来说是相同的,因此在选择算法时只需考虑通信量.
从所给出的这五个并行计算矩阵乘积的算法可以看到,对于方阵的乘积,当p≥4时,Cannon算法具有优越性.
在理论上这些算法的加速比当问题规模增大时都趋向于处理机个数,它们可以由本章开始时给出的加速比的定义以及数据传送时间公式导出,感兴趣的读者可把它们当成练习来做.
§4.
2线性方程组的解法线性方程组是许多重要问题的核心,因此有效地求解线性方程组在科学与工程计算中是非常重要的.
并行计算机的问世,使求解问题的速度和解题规模大幅度地提高,同时也使计算方法产生了变化.
在传统的串行机上,Linpack是求解线性方程组的有效软件包.
然而在并行机上求解此问题,就需要设§4.
2线性方程组的解法51计出适合于该机的并行算法,算法的优劣会对并行机的效率产生很大的影响.
这里考虑的问题是Ax=b这里的任务可以分为两方面,一方面是并行计算矩阵A的LU分解,其中L,U分别是下三角和上三角矩阵,也即存在一排列矩阵Q,使AQ=LU.
另一方面是并行求解三角形方程组,即,求解方程组Ly=b和Ux=y.
下面我们描述有关的算法.
§4.
2.
1分布式系统的并行LU分解算法首先考虑n*n矩阵A=(aij)的串行LU分解法.
根据求解线性方程组的需要,采用部分选主元的Gauss消去法,进行列消去,使得L是单位下三角矩阵.
在算法中ak表示A的第k行.
算法4.
2.
1.
forj=0ton2dondl:|alj|=max{|aij|,i=j,n1}ifl=j,swapAjandAlifajj=0,Aissingularandreturnaij=aij/ajj,i=j+1,n1fork=j+1ton1doaik=aikaij*ajk,i=j+1,n1end{for}end{for}在算法4.
2.
1中,主要计算工作量是修正矩阵A,即做aikaij*ajk,因此并行计算的主要任务就是在多处理机上同时对矩阵A的不同部分做修正.
在多处理机上LU分解的重要工作是使载荷尽可能的平衡,我们采用卷帘(wrap)存储方式在各处理机上分配矩阵A,把矩阵A的第i列存放在Pimodp中.
假设n=p*m,在下面算法中A的第i列为原来A的第i*p+myid列,下面给出在Pmyid上的结点描述.
算法4.
2.
2.
icol=0forj=0ton2doifmyid=jmodpthenndl:|al,icol|=max{|ai,icol|,i=j,n1}ifl=j,swapaj,icolandal,icolifaj,icol=0,Aissingularandkillallprocessesai,icol=ai,icol/aj,icol,i=j+1,n1fij1=ai,icol,i=j+1,n1send(l,myid+1)andsend(f,myid+1)icol+1→icolelserecv(l,myid-1)andrecv(f,myid+1)ifmyid+1=jmodp,send(l,myid+1)andsend(f,myid+1)52第四章矩阵并行计算end{if}ifl=j,swapAjandAlfork=icoltom1doaik=aikfi*ajk,i=j+1,n1end{for}end{for}算法4.
2.
2是在分布式并行计算机上做LU分解的有效方法之一,我们在国产并行机上做了许多实验,效果很好.
如果采用分块卷帘方式存储,增大了算法的粒度,效果更好.
为减少在一个处理机上计算分解因子时其它处理机出现的等待,可以采用提前计算下一次要用到的因子.
也即,算法4.
2.
2中的条件判断由ifmyid=jmodp变成ifmyid=j+1modp,首先只修正计算下一个因子要用到的列,然后就计算分解因子,这对于具有发送消息不需要等待的并行系统来说是非常适合的.
如果并行系统不具备这样的性质,那就不要用这种改进的算法,否则会事倍功半.
§4.
2.
2具有共享存储系统的并行LU分解算法共享存储系统是并行计算机的另一代表,其上的算法设计与分布式系统的算法设计有很大的不同.
主要体现在通信方面,共享存储系统数据的交换是通过共享变量实现的.
在做LU分解时,我们已经知道计算工作分为两个部分,一是计算分解因子,另一是修正矩阵A.
由于需要交换的数据是每步的分解因子,因此我们把分解因子定义成共享变量,并用gfactor表示,与此相对应的选主元的变量记为gl,它们的局部复制品用lfactor和l表示,同时用ag变量作为gfactor被用过与否的标志.
算法如下:算法4.
2.
3.
icol=0forj=0ton2doifmyid=jmodpthenndl:|al,icol|=max{|ai,icol|,i=j,n1}ifl=j,swapaj,icolandal,icolifaj,icol=0,Aissingularandkillallprocessesai,icol=ai,icol/aj,icol,i=j+1,n1lfactorij1=ai,icol,i=j+1,n1ifp1k=0agk=p,waitcopy(lfactor,gfactor),copy(l,gl),agk=0,k=0,p1icol+1→icolelseifagmyid=0,waitcopy(gfactor,lfactor),copy(gl,l),agmyid=1end{if}ifl=j,swapAjandAlfork=icoltom1doaik=aiklfactori*ajk,i=j+1,n1end{for}end{for}§4.
2线性方程组的解法53算法4.
2.
3与算法4.
2.
2在许多方面是相同的,这里的copy代替了算法4.
2.
2中的send和recv.
在这个算法中,同样需要等待.
在这类系统上,可以设计出几乎不需要等待的算法.
在算法4.
2.
3的基础上,对于第0个处理机要首先计算出分解因子,之后的计算过程由如下的修改算法给出:算法4.
2.
4.
icol=0ifmyid=0,icol=1forj=0ton2docid=jmod2p,next=j+1modpifagcid=0,waitcopy(gfactor,lfactor),copy(gl,l)ifmyid=nextthenifl=j,swapal,icolandaj,icolai,icol=ai,icollfactori*aj,icol,i=j+1,n1ndgl:|agl,icol|=max{|ai,icol|,i=j+1,n1}ifgl=j+1,swapaj+1,icolandagl,icolifaj+1,icol=0,Aissingularandkillallprocessesai,icol=ai,icol/aj,icol,i=j+2,n1gfactorij2=ai,icol,i=j+2,n1agcid+1mod2p=0,agcid+2mod2p=1icol+1→icolend{if}ifl=j,swapAjandAlfork=icoltom1doaik=aiklfactori*ajk,i=j+1,n1end{for}end{for}这个算法的特点是计算分解因子时不是在修正全部矩阵之后进行的,而是把当前要计算因子的列做修正就计算分解因子,这样在下一次其它处理机做修正时就不需要等待,实现了计算分解因子和修正矩阵在时间上的重叠,这是在这类机器上做LU分解的一种理想算法,详见文献[5].
§4.
2.
3三角方程组的并行解法本节考虑三角方程组的并行计算方法,我们不妨仅讨论解下三角方程组Lx=b.
三角方程组的并行求解对有效的并行求解线性方程组是不可缺少的,它的并行效果的好坏对求解整个问题有直接的影响,这里给出两种并行实现方法.
首先给出一个串行算法.
算法4.
2.
5.
fori=0ton1doxi=bi/liiforj=i+1ton1dobj=bjlji*xi54第四章矩阵并行计算end{for}end{for}在这个算法中每次对b进行修正时用到L的一列,如果按这种方式并行修正b,则称之为列扫描方法.
对于列扫描算法,原始数据L适合于按行存放,当修正b的值时,就可以并行计算.
同时为使每个处理机的工作量尽可能均衡,要采取卷帘方式存放数据.
正如我们所描述的,为了实现并行计算,需要将每步计算出来的解的一个分量传送到所有其它处理机中,其通信次数是很多的,这对于消息传递型并行系统是不太适合的.
但是对于有共享存储的系统是可以采用这种计算方案的.
像前一节中做LU分解时一样,需要引进共享变量ag和gb,这时在共享存储系统上的算法可描述成如下形式:算法4.
2.
6.
k=0fori=0ton1doifmyid≡imodpthengbi=bk/lki,agi=1,k+1→kend{if}ifagi=1,waitforj=ktom1dobj=bjlji*gbiend{for}end{for}这是一种非常自然的并行实现方法,但是为减少等待,可以采用修正b的一个分量就计算下一个解的方式,这样其它的处理机在做下一次修正时就已经有了要用的新值,而不需要等待.
这些都是算法与程序设计时应该注意的技巧,反复试验就会有所体会.
在大部分现有的共享存储并行系统上,都提供了同步机制,在算法中条件判断的等待也可以由系统提供的库函数来实现,即在此同步.
但对采用异步方式计算下一个新解的算法就不适合了,这时就应用我们在算法中给出的条件判断的等待来实现.
下面介绍一种在分布式并行机上的下三角方程组的求解方法,该方法采用按列卷帘方式存放数据,每次传递的是部分修正的右端项,而不是新求出的解,通过叠加的方式计算下次的新解.
这个算法在分布式系统上被广泛应用,是非常有效的并行算法,其形式如下:算法4.
2.
7.
k=0ifmyid=0,thenui=bi,i=0,n1,vi=0,i=0,p2elseui=0,i=0,n1fori=0ton1doifi>0,recv(v,i1modp)xk=(ui+v0)/likvj=vj+1+ui+1+jli+1+j,k*xk,j=0,p3§4.
3对称正定线性方程组的并行解法55vp2=ui+p1li+p1,k*xksend(v,i+1modp)uj=ujljk*xk,j=i+p,n1k+1→kend{for}关于这个算法的详细讨论可参见文献[6].
§4.
3对称正定线性方程组的并行解法对于对称正定矩阵A的LU分解,我们采用Cholesky分解,也即A=RTR,其中R是上三角矩阵.
由于三角方程组的求解在前一节中已经给出,因此这里着重考虑对称正定矩阵的Cholesky分解.
一个方法是在传统的Cholesky分解列格式算法的基础上,对于发送数据不需等待的并行系统,提出的并行算法.
另一个方法是用双曲旋转变换的方式来做Cholesky分解.
下面就分别介绍这两种算法.
§4.
3.
1Cholesky分解列格式的并行计算这里给出的并行Cholesky分解算法,是在传统的Cholesky分解列格式算法的基础上,结合分布式并行计算机系统的特点给出的[7].
它是在分布式多处理机系统上求Cholesky分解的有效算法.
首先从串行算法出发算法4.
3.
1.
forj=0ton1doajj=ajjj1k=0ajk*ajk,ajj=√ajjfori=j+1ton1aij=(aijj1k=0ajk*aik)/ajjend{for}end{for}在这个算法中,矩阵RT存放在与矩阵A相对应的下三角位置,它对于j循环来说,每次计算出RT的一列,故称之为列格式算法.
由于第j列的计算用到前面的j列的值,因此在并行计算R时就要把它之前的列的信息传送到该列所在的结点上.
在该串行算法的基础上,引入分解因子变量F,它记录当前处理机之前的p1个处理机上的分解因子,是(p1)*n矩阵,原始矩阵A按行卷帘方式存放在处理机中,则在结点Pmyid上的算法如下:算法2fori=0tom1dok=i*p+myid,l=kp+1ifk>0then,recieveGfromPmyid1forj=0top2doai,j+l=(ai,j+lj+l1t=0ait*gjt)/gj,j+lFj=Gj+156第四章矩阵并行计算end{for}aik=aikk1t=0ait*aitaik=√aik,Fp2=AiSendFtoPmyid+1fore=i+1tom1doforj=0top2doae,j+l=(ae,j+lj+l1t=0aet*gjt)/gj,j+laek=(aekk1t=0aet*ait)/aikend{for}end{for}这个并行算法的特点是每步计算出RT的p列,在同一循环中,各处理机计算出的RT的p列是不相同的,从而实现计算与通信的异步进行,减少处理机的等待.
理论分析与数值计算结果可参见文献[7].
§4.
3.
2双曲变换Cholesky分解双曲变换Cholesky分解是指做Cholesky分解时使用如下的双曲变换:H=coshφsinhφsinhφcoshφ在我们的算法中,使用下面的形式:H=(1ρ2)121ρρ1其中ρ=tanh(φ).
通过这个变换把矩阵A化成RTR,首先我们从它的基本理论开始.
假设A=D+UT+U,其中D是对角矩阵,U是严格上三角矩阵,记W=D1/2U,和V=D1/2+W.
通过简单的推导有A=VTVWTW从而有下面的关系式成立:RT0I00IR0=VTWTI00IVW其中I是n*n单位矩阵,为做双曲变换Cholesky分解,要用到下面的一些定义与引理.
定义4.
3.
1.
如果一个2m*2m矩阵Θ满足下面的关系式:ΘTI00IΘ=I00I其中I是m*m单位矩阵,则称之为是伪正交的(pseudo-orthogonal)矩阵.
§4.
3对称正定线性方程组的并行解法57从这个定义可以看到,如果存在一个伪正交矩阵Q使得QVW=R0那么从前面的关系式中就可得出A=RTR,因此主要任务就是寻找伪正交矩阵Q.
下面不加证明地列出文献[8]中的引理.
引理4.
3.
1.
如果R和S都是n*n上三角矩阵,使得RTRSTS对称正定,则R是可逆的,并满足:|skkr1kk|0和K2>0.
令G=K1HK,其中K是2*2对角矩阵.
如果B=GB,则M=KB.
这里通过适当选取K,以达到减少运算次数和开方运算的目的.
在这里我们并不直接计算K,而是用它的平方形式.
假设L=K2,L=K2,则L的计算由下面的引理给出.
引理4.
3.
3.
假设α=L2L1,β=b21b11.
如果选取L=(1αβ2)1L,则G=1αββ1Proof.
从H的定义我们知道ρ=m21m111=K2b21K1b11.
因此就有HK=(1ρ2)12K1K22b21K1b11K2b21b11K2=(1αβ2)12K1αββ1引理得证.
引理4.
3.
4.
如果R和S都是n*n上三角矩阵,并且E和F是对角矩阵,使得RTERSTFS是正定的,并假设αk=FkEk,βk=skkrkk,如果RS=QRS其中Q=Q(n)Q(n1)···Q(1),Q(k),E和F的元素是如下定义的:q(k)ij=1,i=jαkβk,i=k,j=n+kβk,i=n+k,j=k0,其它和Ek=Ek1αkβ2k,Fk=Fk1αkβ2k,1≤k≤n则RTERSTFS=RTERSTFS,并且R是上三角矩阵,S是严格上三角矩阵.
这个引理的证明是容易的,故此略去.
假设A=VTEVWTFW,则修正的双曲Cholesky分解算法可被描述成如下形式:算法4.
3.
3.
V0=V,W0=W,E0=E,F0=Ffori=0ton1do§4.
4三对角方程组的并行解法59Vi+1Wi+1=I00PQ(i)ViWiEi+1k=Eik1αikβik2,Fi+1k=Fik1αikβik2,1≤k≤nFi+11=Fi+1n,Fi+1k+1=Fi+1k,1≤k≤n1end{for}这里Q(i)的定义与引理4.
3.
4中Q的定义相同.
这个算法与算法4.
3.
2一样都是易于并行实现的.
§4.
4三对角方程组的并行解法解三对角线性方程组在偏微分方程数值解中起着非常重要的作用,因此已经有了很多关于它的并行算法,这方面的工作可参见文献[9]–[10].
这里着重介绍几个重要的方法,为不太熟悉并行计算的读者提供一些实例,以便对进一步的并行计算打下基础.
在这里要求解的问题是Ax=d,其中A=T(c,a,b)=a0b0c1a1b1.
.
.
.
.
.
.
.
.
cn2an2bn2cn1an1§4.
4.
1递推法在本章的第二节中,我们已经介绍了如何做矩阵A的LU分解,这里由于A是三对角矩阵,其LU分解可表达成如下形式:A=1e11.
.
.
.
.
.
en21en11*f0b0f1b1.
.
.
.
.
.
fn2bn2fn1由此可以得出:f0=a0,fi1ei=ci,fi=aieibi1,i=1,n1归纳整理有如下的非线性递推关系式:fi=aici*bi1/fi1,i=1,n1如果令fi=ui/ui1,u0=f0,u1=1,则有下面的二阶线性递推关系:uiui1=aiti10ui1ui2,ti=cibi1,i=1,n1我们从一阶线性递推关系式开始进行线性递推关系式的并行计算,它具有如下形式:x0=0,xi=aixi1+bi,1≤i≤n60第四章矩阵并行计算66666666rrrrrrrr66666666rrrrrrrr66666666rrrrrrrrrrrrrrrr77777773333331111a1Q(1,1)a2Q(2,2)a3Q(3,3)a4Q(4,4)a5Q(5,5)a6Q(6,6)a7Q(7,7)a8Q(8,8)a2a1Q(2,1)a3a2Q(3,2)a4a3Q(4,3)a5a4Q(5,4)a6a5Q(6,5)a7a6Q(7,6)a8a7Q(8,7)a3a2a1Q(3,1)a4a3a2a1Q(4,1)a5a4a3a2Q(5,2)a6a5a4a3Q(6,3)a7a6a5a4Q(7,4)a8a7a6a5Q(8,7)Q(1,1)x1=Q(2,1)x2=Q(3,1)x3=Q(4,1)x4=Q(5,1)x5=Q(6,1)x6=Q(7,1)x7=Q(8,1)x8=图4.
1:递推法流程图为了并行计算的需要[11],引入函数Q(t,s)=tj=s(tr=j+1ar)bj其中假定tr=t+1ar=1,1≤s,t≤n,s≤t.
如此定义的函数Q(t,s)具有以下的性质:(1)Q(i,i)=bi,1≤i≤n(2)l≥0,使得s+l0,则A是M矩阵.
下面我们不加证明地引用文献[16]中的两个结论.
引理4.
5.
1.
设A=BC是A的弱正则分裂且A是M矩阵,则迭代形式(4.
5)是异步迭代收敛的.
这个定理可应用到许多迭代法中,比如Jacobi和Gauss-Seidel型迭代法等.
对于这种类型的迭代法,其矩阵A分解成A=DLU,其中D是对角矩阵,L和U分别是严格下和上三角矩阵.
基于这种分裂形式的迭代过程的异步迭代收敛性可用下面的定理给出.
引理4.
5.
2.
设A是严格或不可约对角占优矩阵,A=BC,其中B=DαL,C=(1α)L+U,0≤α≤1,则迭代形式(4.
5)是异步迭代收敛的.
这些结果在求解偏微分方程的差分离散的方程中有很大的作用,因为在采用区域分解法的时候,信息的交换可以不需要等待,所以可以提高并行处理机的效率.
64第四章矩阵并行计算第二部分MPI并行程序设计第五章消息传递并行程序设计平台MPIMPI是英文MessagePassingInterface的缩写,是基于消息传递编写并行程序的一种用户界面.
消息传递是目前并行计算机上广泛使用的一种程序设计模式,特别是对分布式存储的可扩展的并行计算机SPCs(ScalableParallelComputers)和工作站机群NOWs(NetworksofWorkstations)或COWs(ClustersofWorkstations).
尽管还有很多其它的程序实现方式,但是过程之间的通信采用消息传递已经是一种共识.
在MPI和PVM问世以前,并行程序设计与并行计算机系统是密切相关的,对不同的并行计算机就要编写不同的并行程序,给并行程序设计和应用带来了许多麻烦,广大并行计算机的用户迫切需要一些通用的消息传递用户界面,使并行程序具有和串行程序一样的可移植性.
在过去的4年中,国际上确定了MPI为消息传递用户界面标准,自从1994年6月推出MPI以来,它已被广泛接受和使用,目前国际上推出的所有并行计算机都支持MPI和PVM.
对于使用SPCs的用户来说,编写SPMD并行程序使用MPI可能更为方便.
在使用MPI函数与常数时要注意,MPI的所有函数与常数均以MPI_开头.
在C程序中,所有常数的定义除下划线外一律由大写字母组成,在函数和数据类型定义中,接MPI_之后的第一个字母大写,其余全部为小写字母,在所有函数调用之后都将返回一个错误信息码,对于Fortran程序,MPI函数全部以过程方式调用,其错误码以哑元参数返回.
此外,MPI是按进程组(ProcessGroup)方式工作的,所有MPI程序在开始时均被认为是在通信子MPI_COMM_WORLD所拥有的进程组中工作,之后用户可以根据自己的需要,建立其它的进程组.
此外还需注意的是,所有MPI的通信一定要在通信子(communicator)中进行.
此处我们选择了使用MPI编程中常用的函数作介绍,并结合并行计算的实际需要,给出了大量的Fortran程序示例,以便增强对MPI函数的理解,关于MPI的详细介绍可参见文献[1].
§5.
1MPI并行环境管理函数在编写MPI并行程序中,必不可少的二个函数是MPI_Init()和MPI_Finalize().
它们的详细定义如下:MPIINITCintMPI_Init(int*argc,char**argv)FortranMPI_INIT(IERROR)INTEGERIERROR由于ANSIC的主程序main接受参数argc和argv,如此使用MPI_Init可将这些参数传送到每个进程中.
在Fortran程序中,MPI_INIT只返回一个错误码.
这个函数初始化MPI并行程序的执行环境,它必须在调用所有其它MPI函数(除MPI_INITIALIZED)之前被调用,并且在一个MPI程序中,只能被调用一次.
MPIFINALIZECintMPI_Finalize(void)FortranMPI_FINALIZE(IERROR)INTEGERIERROR6768第五章消息传递并行程序设计平台MPI这个函数清除MPI环境的所有状态.
即一但它被调用,所有MPI函数都不能再调用,其中包括MPI_INIT.
在一个程序中,如果不清楚是否已经调用了MPI_INIT,可以使用MPI_INITIALIZED来检查,它是唯一的可以在调用MPI_INIT之前使用的函数.
MPIINITIALIZEDCintMPI_Initialized(intflag)FortranMPI_INITIALIZED(FLAG,IERROR)LOGICALFLAGINTEGERIERROR参数说明OUTFLAG,如果MPI_INIT被调用,返回值为TRUE,否则为FALSE.
另外在MPI中还提供了一个用来检查出何种错误的函数MPI_ERROR_STRING,它将给出错误信息.
MPIERRORSTRINGCintMPI_Error_string(interrorcode,char*string,int*len)FortranMPI_ERROR_STRING(ERRORCODE,STRING,LEN,IERROR)INTEGERERRORCODE,LEN,IERRORCHARACTER*(*)STRING参数说明INERRORCODE,由MPI函数返回的错误码.
OUTSTRING,对应ERRORCODE的错误信息.
OUTLEN,错误信息STRING的长度.
在使用此函数时,应注意STRING的长度最小应为MPI_MAX_ERROR_STRING.
§5.
2进程控制函数在这一章中,介绍与进程组有关的一些基本函数,其中包括如何建立进程组和通信子,灵活使用这些函数在实际程序设计中会带来巨大方便.
§5.
3MPI进程组操作函数MPICOMMGROUPCintMPI_Comm_group(MPI_Commcomm,MPI_Group*group)FortranMPI_COMM_GROUP(COMM,GROUP,IERROR)INTEGERCOMM,GROUP,IERROR参数说明INCOMM,通信子.
OUTGROUP,对应COMM的进程组.
§5.
3MPI进程组操作函数69MPI_COMM_GROUP这个函数是用来建立一个通信子对应的新进程组,之后就可以对此进程组进行需要的操作.
MPIGROUPFREECintMPI_group_free(MPI_Group*group)FortranMPI_GROUP_FREE(GROUP,IERROR)INTEGERGROUP,IERROR参数说明INOUTGROUP,释放进程组并返回MPI_GROUP_NULL.
当MPI_GROUP_FREE被调用之后,任何关于此进程组的操作将视为无效.
MPIGROUPSIZECintMPI_Group_size(MPI_Groupgroup,int*size)FortranMPI_GROUP_SIZE(GROUP,SIZE,IERROR)INTEGERGROUP,SIZE,IERROR参数说明INGROUP,进程组.
OUTSIZE,进程组中的进程个数.
如果进程组是MPI_GROUP_EMPTY,则返回值SIZE为0.
MPIGROUPRANKCintMPI_Group_rank(MPI_Groupgroup,int*rank)FortranMPI_GROUP_RANK(GROUP,RANK,IERROR)INTEGERGROUP,RANK,IERROR参数说明INGROUP,进程组.
OUTRANK,进程在进程组中的编号.
如果进程不是进程组中的成员,则返回值RANK为MPI_UNDEFINED.
MPIGROUPTRANSLATERANKSCintMPI_Group_translate_ranks(MPI_Groupgroup1,intn,int*ranks1,MPI_Groupgroup2,int*ranks2)FortranMPI_GROUP_TRANSLATE_RANKS(GROUP1,N,RANKS1,GROUP2,RANKS2,IERROR)INTEGERGROUP1,N,RANKS1(*),GROUP2,RANKS2(*),IERROR70第五章消息传递并行程序设计平台MPI参数说明INGROUP1,进程组1.
INN,RANKS1和RANKS2中数组元素个数.
INRANKS1,进程组1中有效编号组成的数组.
INGROUP2,进程组2.
OUTRANKS2,RANKS1中的元素在进程组2中的对应编号.
如果属于进程组1的某个进程可以在RANKS1中找到,而这个进程不属于进程组2,则在RANKS2中对应RANKS1的位置返回值为MPI_UNDEFINED.
MPIGROUPINCLCintMPI_Group_incl(MPI_Groupgroup,intn,int*ranks,MPI_Groupnewgroup)FortranMPI_GROUP_INCL(GROUP,N,RANKS,NEWGROUP,IERROR)INTEGERGROUP,N,RANKS(*),NEWGROUP,IERROR参数说明INGROUP,进程组.
INN,RANKS数组中元素的个数和新进程组的大小.
INRANKS,将在新进程组中出现的旧进程组中的编号.
OUTNEWGROUP,由RANKS定义的顺序导出的新进程组.
MPIGROUPEXCLCintMPI_Group_excl(MPI_Groupgroup,intn,int*ranks,MPI_Groupnewgroup)FortranMPI_GROUP_EXCL(GROUP,N,RANKS,NEWGROUP,IERROR)INTEGERGROUP,N,RANKS(*),NEWGROUP,IERROR参数说明INGROUP,进程组.
INN,RANKS数组中元素的个数.
INRANKS,在新进程组中不出现的旧进程组中的编号.
OUTNEWGROUP,旧进程组中不在RANKS里的元素组成的新进程组.
MPIGROUPUNIONCintMPI_Group_union(MPI_Groupgroup1,MPI_Groupgroup2,MPI_Groupnewgroup)FortranMPI_GROUP_UNION(GROUP1,GROUP2,NEWGROUP,IERROR)INTEGERGROUP1,GROUP2,NEWGROUP,IERROR§5.
4MPI通信子操作71参数说明INGROUP1,进程组1.
INGROUP2,进程组2.
OUTNEWGROUP,进程组1和进程组2的并.
MPIGROUPINTERSECTIONCintMPI_Group_intersection(MPI_Groupgroup1,MPI_Groupgroup2,MPI_Groupnewgroup)FortranMPI_GROUP_INTERSECTION(GROUP1,GROUP2,NEWGROUP,IERROR)INTEGERGROUP1,GROUP2,NEWGROUP,IERROR参数说明INGROUP1,进程组1.
INGROUP2,进程组2.
OUTNEWGROUP,进程组1和进程组2的交.
MPIGROUPDIFFERENCECintMPI_Group_difference(MPI_Groupgroup1,MPI_Groupgroup2,MPI_Groupnewgroup)FortranMPI_GROUP_DIFFERENCE(GROUP1,GROUP2,NEWGROUP,IERROR)INTEGERGROUP1,GROUP2,NEWGROUP,IERROR参数说明INGROUP1,进程组1.
INGROUP2,进程组2.
OUTNEWGROUP,进程组1和进程组2的差.
以上关于进程组操作的函数功能,和数学中集合运算相同.
§5.
4MPI通信子操作下面介绍的关于通信子的操作和上述关于进程组操作十分类似.
MPICOMMSIZECintMPI_Comm_size(MPI_Commcomm,int*size)FortranMPI_COMM_SIZE(COMM,SIZE,IERROR)INTEGERCOMM,SIZE,IERROR参数说明INCOMM,通信子.
OUTSIZE,通信子中的进程个数.
72第五章消息传递并行程序设计平台MPIMPICOMMRANKCintMPI_Comm_rank(MPI_Commcomm,int*rank)FortranMPI_COMM_RANK(COMM,RANK,IERROR)INTEGERCOMM,RANK,IERROR参数说明INCOMM,通信子.
OUTRANK,通信子中的进程编号.
MPICOMMDUPCintMPI_Comm_dup(MPI_Commcomm,MPI_Comm*newcomm)FortranMPI_COMM_DUP(COMM,NEWCOMM,IERROR)INTEGERCOMM,NEWCOMM,IERROR参数说明INCOMM,通信子.
OUTNEWCOMM,COMM通信子的复制.
MPICOMMCREATECintMPI_Comm_create(MPI_Commcomm,MPI_Groupgroup,MPI_Comm*newcomm)FortranMPI_COMM_CREATE(COMM,GROUP,NEWCOMM,IERROR)INTEGERCOMM,GROUP,NEWCOMM,IERROR参数说明INCOMM,通信子.
INGROUP,通信子COMM的一个子集.
OUTNEWCOMM,对应GROUP的新通信子.
MPICOMMSPLITCintMPI_Comm_split(MPI_Commcomm,intcolor,intkey,MPI_Comm*newcomm)FortranMPI_COMM_SPLIT(COMM,COLOR,KEY,NEWCOMM,IERROR)INTEGERCOMM,COLOR,KEY,NEWCOMM,IERROR参数说明INCOMM,通信子.
INCOLOR,子集控制值.
INKEY,子集中进程编号的顺序.
OUTNEWCOMM,由此产生的新通信子.
§5.
5点到点通信函数73这个函数划分COMM所对应的进程组为不相交的子进程组,每个子进程组有一个共同的值COLOR,就是说,每个子进程组中包含COLOR相同的所有进程.
MPICOMMFREECintMPI_Comm_free(MPI_Comm*comm)FortranMPI_COMM_FREE(COMM,IERROR)INTEGERCOMM,IERROR参数说明INOUTCOMM,通信子.
§5.
5点到点通信函数点到点通信是指一对进程之间的数据转换,也就是说,一边发送数据另一边接收数据.
点到点通信是MPI通信机制的基础,它分为同步通信和异步通信二种机制.
§5.
6阻塞式通信函数MPISENDCintMPI_Send(void*buf,intcount,MPI_Datatypedatatype,intdest,inttag,MPI_Commcomm)FortranMPI_SEND(BUF,COUNT,DATATYPE,DEST,TAG,COMM,IERROR)BUF(*)INTEGERCOUNT,DATATYPE,DEST,TAG,COMM,IERROR参数说明INBUF,所要发送消息数据的首地址.
INCOUNT,发送消息数组元素的个数.
INDATATYPE,发送消息的数据类型.
INDEST,接收消息的进程编号.
INTAG,消息标签.
INCOMM,通信子.
这里COUNT是BUF的元素个数而不是字节数.
MPIRECVCintMPI_Recv(void*buf,intcount,MPI_Datatypedatatype,intsource,inttag,MPI_Commcomm,MPI_Status*status)FortranMPI_RECV(BUF,COUNT,DATATYPE,SOURCE,TAG,COMM,STATUS,IERROR)BUF(*)INTEGERCOUNT,DATATYPE,SOURCE,TAG,COMM,STATUS(MPI_STATUS_SIZE),IERROR74第五章消息传递并行程序设计平台MPI参数说明OUTBUF,接收消息数据的首地址.
INCOUNT,接收消息数组元素的最大个数.
INDATATYPE,接收消息的数据类型.
INSOURCE,发送消息的进程编号.
INTAG,消息标签.
INCOMM,通信子.
OUTSTATUS,接收消息时返回的状态.
在这个接收函数中,可以不指定SOURCE和TAG,而分别用MPI_ANY_SOURCE和MPI_ANY_TAG来代替,用于接收任何进程发送的消息或任何编号的消息.
接收消息时返回的状态STATUS在C语言中是用结构定义的,在Fortran中是用数组定义的,其中包括MPI_SOURCE,MPI_TAG和MPI_ERROR.
此外STATUS还包含接收消息元素的个数,但它不是显式给出的,需要用到后面给出的函数MPI_GET_COUNT.
MPI的基本数据类型定义与相应的Fortran和C的数据类型对照关系如下:Fortran程序中可使用的MPI数据类型MPIdatatypeFortrandatatypeMPI_INTEGERINTEGERMPI_REALREALMPI_DOUBLE_PRECISIONDOUBLEPRECISIONMPI_COMPLEXCOMPLEXMPI_LOGICALLOGICALMPI_CHARACTERCHARACTER(1)MPI_BYTEMPI_PACKEDC程序中可使用的MPI数据类型MPIdatatypeFortrandatatypeMPI_CHARsignedcharMPI_SHORTsignedshortintMPI_INTsignedintMPI_LONGsignedlongintMPI_UNSIGNED_CHARunsignedcharMPI_UNSIGNED_SHORTunsignedshortintMPI_UNSIGNED_INTunsignedintMPI_UNSIGNED_LONGunsignedlongintMPI_FLOATfloatMPI_DOUBLEdoubleMPI_LONG_DOUBLElongdoubleMPI_BYTEMPI_PACKED这里MPI_BYTE和MPI_PACKED不对应Fortran或C中的任何数据类型,MPI_BYTE是由一个字节组成的,而MPI_PACKED将在后面介绍.
到目前为止,我们已经可以编写一些简单的MPI程序,下面给出一个MPI程序的例子,我们结合这个例子说明一个MPI程序的组成部分.
重点需要掌握的MPI函数有§5.
6阻塞式通信函数75初始化并行环境的函数MPI_INIT,得到本进程的编号函数MPI_COMM_RANK,得到全部进程个数的函数MPI_COMM_SIZE,退出MPI并行环境的函数MPI_FINALIZE,以及进行消息传递的二个函数MPI_SEND和MPI_RECV.
例5.
6.
1.
假设一共有p个进程,在进程编号为myid(myid=0,p1)的进程中有一个整数m,我们要把m传送到进程(myid+1)modp中.
programringccTheheaderfilempif.
hmustbeincludedwhenyouuseMPIfuctions.
cinclude'mpif.
h'integermyid,p,mycomm,ierr,m,status(mpi_status_size),&next,front,mod,nccCreateMPIparallelenvironmentandgetthenecessarydata.
ccallmpi_init(ierr)callmpi_comm_dup(mpi_comm_world,mycomm,ierr)callmpi_comm_rank(mycomm,myid,ierr)callmpi_comm_size(mycomm,p,ierr)ccBeginningthemainparallelworkoneachprocess.
cm=myidfront=mod(p+myid-1,p)next=mod(myid+1,p)ccCommunicationwitheachother.
cif(myid.
eq.
0)thencallmpi_recv(n,1,mpi_integer,front,1,mycomm,status,ierr)callmpi_send(m,1,mpi_integer,next,1,mycomm,ierr)m=nelsecallmpi_send(m,1,mpi_integer,next,1,mycomm,ierr)callmpi_recv(m,1,mpi_integer,front,1,mycomm,status,ierr)endifccEndingofparallelwork.
cprint*,'Thevalueofmis',m,'onProcess',myidcallmpi_comm_free(mycomm,ierr)ccRemoveMPIparallelenvironment.
ccallmpi_finalize(ierr)endMPIGETCOUNTCintMPI_Get_count(MPI_Statusstatus,MPI_Datatypedatatype,int*count)FortranMPI_GET_COUNT(STATUS,DATATYPE,COUNT,IERROR)INTEGERSTATUS(MPI_STATUS_SIZE),DATATYPE,COUNT,IERROR76第五章消息传递并行程序设计平台MPI参数说明INSTATUS,接收消息时返回的状态.
INDATATYPE,接收消息的数据类型.
OUTCOUNT,接收消息数组元素的个数.
由于在接收消息时使用的是最大个数,为了准确知道接收消息的个数,就要使用此函数.
MPISENDRECVCintMPI_Sendrecv(void*sendbuf,intsendcount,MPI_Datatypesendtype,intdest,intsendtag,void*recvbuf,intrecvcount,MPI_Datatyperecvtype,intsource,intrecvtag,MPI_Commcomm,MPI_Status*status)FortranMPI_SENDRECV(SENDBUF,SENDCOUNT,SENDTYPE,DEST,SENDTAG,RECVBUF,RECVCOUNT,RECVTYPE,SOURCE,RECVTAG,COMM,STATUS,IERROR)SENDBUF(*),RECVBUF(*)INTEGERSENDCOUNT,SENDTYPE,DEST,SENDTAG,RECVCOUNT,RECVTYPESOURCE,RECVTAG,COMM,STATUS(MPI_STATUS_SIZE),IERROR参数说明INSENDBUF,所要发送消息数据的首地址.
INSENDCOUNT,发送消息数组元素的个数.
INSENDTYPE,发送消息的数据类型.
INDEST,接收消息的进程编号.
INSENDTAG,发送消息标签.
OUTRECVBUF,接收消息数据的首地址.
INRECVCOUNT,接收消息数组元素的最大个数.
INRECVTYPE,接收消息的数据类型.
INSOURCE,发送消息的进程编号.
INRECVTAG,接收消息标签.
INCOMM,通信子.
OUTSTATUS,接收消息时返回的状态.
这是发送消息和接收消息组合在一起的一个函数,好处是不用考虑先发送还是先接收消息.
在示例5.
6.
1中从通信开始到通信结束部分可用如下的一个函数来完成:callmpi_sendrecv(m,1,mpi_integer,next,1,n,1,mpi_integer,&front,1,mycomm,status,ierr)m=n由此可见使用MPI_SENDRECV可以使程序简化,更重要的是在通信过程中不需要考虑哪个进程先发送还是先接收,从而可以避免消息传递过程中的死锁.
§5.
6阻塞式通信函数77MPISENDRECVREPLACECintMPI_Sendrecv_replace(void*sendbuf,intcount,MPI_Datatypedatatype,intdest,intsendtag,intsource,intrecvtag,MPI_Commcomm,MPI_Status*status)FortranMPI_SENDRECV_REPLACE(BUF,COUNT,DATATYPE,DEST,SENDTAG,SOURCE,RECVTAG,COMM,STATUS,IERROR)BUF(*)INTEGERCOUNT,DATATYPE,DEST,SENDTAG,SOURCE,RECVTAG,COMM,STATUS(MPI_STATUS_SIZE),IERROR参数说明OUTBUF,发送和接收消息数据的首地址.
INCOUNT,发送和接收消息数组元素的个数.
INDATATYPE,发送和接收消息的数据类型.
INDEST,接收消息的进程编号.
INSENDTAG,发送消息标签.
INSOURCE,发送消息的进程编号.
INRECVTAG,接收消息标签.
INCOMM,通信子.
OUTSTATUS,接收消息时返回的状态.
在这个函数中,发送和接收使用同一个消息缓冲区BUF.
因此示例5.
6.
1的通信部分可进一步改写成如下的形式:callmpi_sendrecv_replace(m,1,mpi_integer,next,1,&front,1,mycomm,status,ierr)在进行环形通信时,使用MPI_SENDRECV_REPLACE是非常方便的.
为了更方便使用此函数,MPI定义了一个空进程MPI_PROC_NULL,如果通信采用这个空进程将不起任何作用.
比如在示例5.
6.
1中,如果我们不需要从第p1个进程向第0个进程传送m,则可用如下方式来实现:.
.
.
.
.
.
if(myid.
eq.
0)front=MPI_PROC_NULLif(myid.
eq.
p-1)next=MPI_PROC_NULLcallmpi_sendrecv_replace(m,1,mpi_integer,next,1,&front,1,mycomm,status,ierr)MPIPROBECintMPI_Probe(intsource,inttag,MPI_Commcomm,MPI_Status*status)FortranMPI_PROBE(SOURCE,TAG,COMM,STATUS,IERROR)INTEGERSOURCE,TAG,COMM,STATUS(MPI_STATUS_SIZE),IERROR参数说明INSOURCE,发送消息进程的编号.
INTAG,接收消息的标签.
INCOMM,通信子.
OUTSTATUS,返回到达消息的状态.
78第五章消息传递并行程序设计平台MPI此函数与下面给出的MPI_IPROBE的功能基本类似,但是MPI_PROBE一定要等到消息才返回.
MPIIPROBECintMPI_Iprobe(intsource,inttag,MPI_Commcomm,int*flag,MPI_Status*status)FortranMPI_IPROBE(SOURCE,TAG,COMM,FLAG,STATUS,IERROR)LOGICALFLAGINTEGERSOURCE,TAG,COMM,STATUS(MPI_STATUS_SIZE),IERROR参数说明INSOURCE,发送消息进程的编号.
INTAG,接收消息的标签.
INCOMM,通信子.
OUTFLAG,如果指定消息已经达到,FLAG返回值为TRUE.
OUTSTATUS,返回到达消息的状态.
这个函数对由SOURCE,TAG和COMM确定的消息是否到达,它都将立即返回.
这里SOURCE和TAG同MPI_RECV中一样,可以用通配符(wildcard)MPI_ANY_SOURCE和MPI_ANY_TAG来代替.
§5.
7非阻塞式通信函数非阻塞式通信函数是指在通信过程中,不需要等待通信结束就返回,通常这种通信过程交由计算机的后台来处理.
如果计算机系统提供硬件支持非阻塞式通信函数,就可以使计算与通信在时间上的重叠,从而提高并行计算的效率.
MPIISENDCintMPI_Isend(void*buf,intcount,MPI_Datatypedatatype,intdest,inttag,MPI_Commcomm,MPI_Request*request)FortranMPI_ISEND(BUF,COUNT,DATATYPE,DEST,TAG,COMM,REQUEST,IERROR)BUF(*)INTEGERCOUNT,DATATYPE,DEST,TAG,COMM,REQUEST,IERROR参数说明INBUF,所要发送消息数据的首地址.
INCOUNT,发送消息数组元素的个数.
INDATATYPE,发送消息的数据类型.
INDEST,接收消息的进程编号.
INTAG,消息标签.
INCOMM,通信子.
OUTREQUEST,请求句柄以备将来查询.
§5.
7非阻塞式通信函数79MPIIRECVCintMPI_Irecv(void*buf,intcount,MPI_Datatypedatatype,intsource,inttag,MPI_Commcomm,MPI_Request*request)FortranMPI_IRECV(BUF,COUNT,DATATYPE,SOURCE,TAG,COMM,REQUEST,IERROR)BUF(*)INTEGERCOUNT,DATATYPE,SOURCE,TAG,COMM,REQUEST,IERROR参数说明OUTBUF,接收消息数据的首地址.
INCOUNT,接收消息数组元素的个数.
INDATATYPE,接收消息的数据类型.
INSOURCE,发送消息的进程编号.
INTAG,消息标签.
INCOMM,通信子.
OUTREQUEST,请求句柄以备将来查询.
MPI_ISEND和MPI_IRECV不需要等待发送或接收消息完成就可执行其余的任务,看发送或接收过程是否结束,有下面的一些函数:MPIWAITCintMPI_Wait(MPI_Request*request,MPI_Status*status)FortranMPI_WAIT(REQUEST,STATUS,IERROR)INTEGERREQUEST,STATUS(MPI_STATUS_SIZE),IERROR参数说明INOUTREQUEST,请求句柄.
OUTSTATUS,发送或接收消息的状态.
如果REQUEST所指的操作已经完成,MPI_WAIT将结束等待状态.
MPITESTCintMPI_Test(MPI_Request*request,int*flag,MPI_Status*status)FortranMPI_TEST(REQUEST,FLAG,STATUS,IERROR)LOGICALFLAGINTEGERREQUEST,STATUS(MPI_STATUS_SIZE),IERROR参数说明INOUTREQUEST,请求句柄.
OUTFLAG,REQUEST所指的操作已经完成返回值为TRUE.
OUTSTATUS,发送或接收消息的状态.
以上这二个函数中STATUS对发送消息操作是没有定义的,唯一可在发送消息操作中使用STATUS的是查询函数MPI_TEST_CANCELLED.
80第五章消息传递并行程序设计平台MPIMPIREQUESTFREECintMPI_Request_free(MPI_Request*request)FortranMPI_REQUEST_FREE(REQUEST,IERROR)INTEGERREQUEST,IERROR参数说明INOUTREQUEST,请求句柄,返回值为MPI_REQUEST_NULL.
以上是对单个REQUEST进行查询的函数,如果要对多个请求句柄进行查询,有下面的一些函数可以使用.
MPIWAITANYCintMPI_Waitany(intcount,MPI_Request*array_of_requests,int*index,MPI_Status*status)FortranMPI_WAITANY(COUNT,ARRAY_OF_REQUESTS,INDEX,STATUS,IERROR)INTEGERCOUNT,ARRAY_OF_REQUESTS(*),INDEX,STATUS(MPI_STATUS_SIZE),IERROR参数说明INCOUNT,请求句柄的个数.
INOUTARRAY_OF_REQUESTS,请求句柄数组.
OUTINDEX,已经完成通信操作的句柄指标.
OUTSTATUS,消息的状态.
这个函数当所有请求句柄中至少有一个已经完成通信操作,就返回,如果有多于一个请求句柄已经完成,MPI_WAITANY将随机选择其中的一个并立即返回.
MPITESTANYCintMPI_Testany(intcount,MPI_Request*array_of_requestS,int*index,int*flag,MPI_Status*status)FortranMPI_TESTANY(COUNT,ARRAY_OF_REQUESTS,INDEX,FLAG,STATUS,IERROR)LOGICALFLAGINTEGERCOUNT,ARRAY_OF_REQUESTS(*),INDEX,STATUS(MPI_STATUS_SIZE),IERROR参数说明INCOUNT,请求句柄的个数.
INOUTARRAY_OF_REQUESTS,请求句柄数组.
OUTINDEX,已经完成通信操作的句柄指标.
OUTFLAG,如果有一个已经完成,则FLAG=TRUE.
OUTSTATUS,消息的状态.
这个函数无论有没有通信操作完成都将立即返回.
§5.
7非阻塞式通信函数81MPIWAITALLCintMPI_Waitall(intcount,MPI_Request*array_of_requests,MPI_Status*array_of_statuses)FortranMPI_WAITALL(COUNT,ARRAY_OF_REQUESTS,ARRAY_OF_STATUSES,IERROR)INTEGERCOUNT,ARRAY_OF_REQUESTS(*),ARRAY_OF_STATUSES(MPI_STATUS_SIZE,*),IERROR参数说明INCOUNT,请求句柄的个数.
INOUTARRAY_OF_REQUESTS,请求句柄数组.
INOUTARRAY_OF_STATUSES,所有消息的状态数组.
这个函数当所有通信操作完成之后才返回,否则将一直等待.
MPITESTALLCintMPI_Testall(intcount,MPI_Request*array_of_requests,int*flag,MPI_Status*array_of_statuses)FortranMPI_TESTALL(COUNT,ARRAY_OF_REQUESTS,FLAG,ARRAY_OF_STATUSES,IERROR)LOGICALFLAGINTEGERCOUNT,ARRAY_OF_REQUESTS(*),ARRAY_OF_STATUSES(MPI_STATUS_SIZE,*),IERROR参数说明INCOUNT,请求句柄的个数.
INOUTARRAY_OF_REQUESTS,请求句柄数组.
OUTFLAG,如果有一个没完成,则FLAG=FALSE.
INOUTARRAY_OF_STATUSES,所有消息的状态数组.
这个函数无论所有通信操作是否完成都将立即返回.
MPICANCELCintMPI_Cancel(MPI_Request*request)FortranMPI_CANCEL(REQUEST,IERROR)INTEGERREQUEST,IERROR参数说明INOUTREQUEST,请求句柄.
这个函数用来取消一个发送或接收操作,其后要使用MPI_WAIT或MPI_TEST等函数,因此使用此函数会损失许多效率.
但是要知道关于REQUEST请求句柄的操作是否已经被取消,还要使用下面的检测函数:82第五章消息传递并行程序设计平台MPIMPITESTCANCELLEDCintMPI_Test_cancelled(MPI_Status*status,int*flag)FortranMPI_TEST_CANCELLED(STATUS,FLAG,IERROR)LOGICALFLAGINTEGERREQUEST,IERROR参数说明INSTATUS,消息的状态.
OUTFLAG,如果已经取消,则FLAG=TRUE.
§5.
8特殊的点到点通信函数MPISENDINITCintMPI_Send_init(void*buf,intcount,MPI_Datatypedatatype,intdest,inttag,MPI_Commcomm,MPI_Request*request)FortranMPI_SEND_INIT(BUF,COUNT,DATATYPE,DEST,TAG,COMM,REQUEST,IERROR)BUF(*)INTEGERCOUNT,DATATYPE,DEST,TAG,COMM,REQUEST,IERROR参数说明INBUF,所要发送消息数据的首地址.
INCOUNT,发送消息数组元素的个数.
INDATATYPE,发送消息的数据类型.
INDEST,接收消息的进程编号.
INTAG,消息标签.
INCOMM,通信子.
OUTREQUEST,请求句柄以备将来查询.
MPIRECVINITCintMPI_Recv_init(void*buf,intcount,MPI_Datatypedatatype,intsource,inttag,MPI_Commcomm,MPI_Request*request)FortranMPI_RECV_INIT(BUF,COUNT,DATATYPE,SOURCE,TAG,COMM,REQUEST,IERROR)BUF(*)INTEGERCOUNT,DATATYPE,SOURCE,TAG,COMM,REQUEST,IERROR参数说明OUTBUF,接收消息数据的首地址.
INCOUNT,接收消息数组元素的个数.
INDATATYPE,接收消息的数据类型.
INSOURCE,发送消息的进程编号.
INTAG,消息标签.
INCOMM,通信子.
OUTREQUEST,请求句柄以备将来查询.
§5.
8特殊的点到点通信函数83MPI_SEND_INIT和MPI_RECV_INIT是二个持久性通信函数,它们的参数和前面讲到的异步通信中的发送与接收相同,由此可以把通信以通道方式建立起来对应关系,从而提高通信效率.
但是仅有这二个函数还不能够达到通信的目的,它们还需要用MPI_START或MPI_STARTALL来激活.
MPISTARTCintMPI_Start(MPI_Request*request)FortranMPI_START(REQUEST,IERROR)INTEGERREQUEST,IERROR参数说明INOUTREQUEST,请求句柄.
在使用MPI_START之后,MPI_SEND_INIT和MPI_RECV_INIT就同异步发送与接收.
MPISTARTALLCintMPI_Startall(intcount,MPI_Request*array_of_request)FortranMPI_STARTALL(COUNT,ARRAY_OF_REQUEST,IERROR)INTEGERCOUNT,ARRAY_OF_REQUEST,IERROR参数说明INCOUNT,需要激活的请求句柄个数.
INOUTARRAY_OF_REQUEST,请求句柄数组.
例5.
8.
1.
在示例5.
6.
1中,整数m只是在进程环中向前移动一次,如果我们要移动多次,可以使用下面的程序实现:programmringinclude'mpif.
h'cintegeriterparameter(iter=3)integermyid,p,mycomm,ierr,m,status(mpi_status_size,2),&next,front,mod,n,i,requests(2)ccallmpi_init(ierr)callmpi_comm_dup(mpi_comm_world,mycomm,ierr)callmpi_comm_rank(mycomm,myid,ierr)callmpi_comm_size(mycomm,p,ierr)print*,'Process',myid,'of',p,'isrunning'cm=myidfront=mod(p+myid-1,p)next=mod(myid+1,p)callmpi_send_init(m,1,mpi_integer,next,1,mycomm,&requests(1),ierr)callmpi_recv_init(n,1,mpi_integer,front,1,mycomm,&requests(2),ierr)do10i=1,iter84第五章消息传递并行程序设计平台MPIcallmpi_startall(2,requests,ierr)callmpi_waitall(2,requests,status,ierr)m=n10continueprint*,'Thevalueofmis',m,'onProcess',myidccallmpi_comm_free(mycomm,ierr)callmpi_finalize(ierr)stopend§5.
9MPI的通信模式前面所提到的通信模式是标准通信模式,它是由MPI实现中确定的,可能选择的是缓冲区模式(buered-mode)或同步模式(synchronous-mode),此外MPI还提供了就绪模式(ready-mode).
这些模式仅对发送消息起作用,无论用任何模式发送的消息都可用一个接收函数来完成接收消息.
在下面给出的与发送消息有关的函数,其参数与对应标准通信模式时是完全一致的,因此这里我们仅给出函数名字及简要说明.
MPI_BSEND,使用缓冲区发送消息,只要缓冲区足够大,而不管接收进程是否开始接收它都将立即返回,也就是说,在发送消息时是不需要等待的.
MPI_SSEND,它可以在没有接收信号时就开始发送,但要等到接收完成之后才能结束通信操作,好处是无需用户自己申请缓冲区.
MPI_RSEND,只当已经有接收信号时才开始发送消息,否则将出现错误.
MPI_IBSEND,异步通信的MPI_BSEND.
MPI_ISSEND,异步通信的MPI_SSEND.
MPI_IRSEND,异步通信的MPI_RSEND.
MPI_BSEND_INIT,建立自己缓冲区的MPI_SEND_INIT.
MPI_SSEND_INIT,同步模式的MPI_SEND_INIT.
MPI_RSEND_INIT,就绪模式的MPI_SEND_INIT.
有关消息缓冲区的管理,可以使用如下的二个函数:MPIBUFFERATTACHCintMPI_Buffer_attach(void*buffer,intsize)FortranMPI_BUFFER_ATTACH(BUFFER,SIZE,IERROR)BUFFER(*)INTEGERSIZE,IERROR参数说明INBUFFER,初始化缓冲区的首地址.
INSIZE,缓冲区按字节的大小.
这个函数指明使用消息缓冲区模式通信时用哪一个缓冲区.
§5.
10用户定义的数据类型与打包85MPIBUFFERDETACHCintMPI_Buffer_detach(void*buffer,int*size)FortranMPI_BUFFER_DETACH(BUFFER,SIZE,IERROR)BUFFER(*)INTEGERSIZE,IERROR参数说明OUTBUFFER,要拆卸的缓冲区.
OUTSIZE,返回缓冲区的字节数.
§5.
10用户定义的数据类型与打包§5.
11用户定义的数据类型在MPI中用户定义的数据类型和C语言中的结构是非常类似的,它包括数据类型、数据长度和数据起始位置,这个新定义的数据类型同其它MPI定义的标准数据类型一样,可以用来作为发送和接收消息的数据类型.
对于这个新定义的数据类型,具体地说它的组成可以概括为如下的形象描述:它包含了一系列的原始数据类型;和一系列的称之为位移(displacement)的整数.
数据类型与位移是成对出现的,这一系列的数据对就是称之为类型映射(typemap).
用类型映射来表示一个新的数据类型为:Typemap={(type0,disp0)typen1,dispn1)}在这个新的数据类型中,包括位移是dispi的原始数据类型typei,其中i=0,n1.
如何构造新的数据类型我们将在以下进行详细介绍.
MPITYPECONTIGUOUSCintMPI_Type_contiguous(intcount,MPI_Datatypeoldtype,MPI_Datatype*newtype)FortranMPI_TYPE_CONTIGUOUS(COUNT,OLDTYPE,NEWTYPE,IERROR)INTEGERCOUNT,OLDTYPE,NEWTYPE,IERROR参数说明INCOUNT,重复放在一起的OLDTYPE的元素个数.
INOLDTYPE,原始数据类型.
OUTNEWTYPE,新构造的数据类型.
这是构造数据类型中最简单的一种,就是把原始数据类型重复COUNT次,我们举个简单的例子来说明其含义.
假设COUNT=5,OLDTYPE=INTEGER,则调用此函数之后产生的新的数据类型NEWTYPE为5个整数.
86第五章消息传递并行程序设计平台MPIMPITYPEVECTORCintMPI_Type_vector(intcount,intblocklength,intstride,MPI_Datatypeoldtype,MPI_Datatype*newtype)FortranMPI_TYPE_VECTOR(COUNT,BLOCKLENGTH,STRIDE,OLDTYPE,NEWTYPE,IERROR)INTEGERCOUNT,BLOCKLENGTH,STRIDE,OLDTYPE,NEWTYPE,IERROR参数说明INCOUNT,OLDTYPE的块的个数.
INBLOCKLENGTH,每块的长度.
INSTRIDE,用OLDTYPE的元素度量的每块起始的间距.
INOLDTYPE,原始数据类型.
OUTNEWTYPE,新构造的数据类型.
假设COUNT=3,BLOCKLENGTH=2,STRIDE=3,OLDTYPE=INTEGER,则调用此函数得到的新数据类型NEWTYPE为6个整数组成.
这6个整数是在原始数组中从开始取2个,隔一个再取2个,再隔一个再取2个构成的.
MPITYPEHVECTORCintMPI_Type_hvector(intcount,intblocklength,intstride,MPI_Datatypeoldtype,MPI_Datatype*newtype)FortranMPI_TYPE_HVECTOR(COUNT,BLOCKLENGTH,STRIDE,OLDTYPE,NEWTYPE,IERROR)INTEGERCOUNT,BLOCKLENGTH,STRIDE,OLDTYPE,NEWTYPE,IERROR参数说明INCOUNT,OLDTYPE的块的个数.
INBLOCKLENGTH,每块的长度.
INSTRIDE,用字节数度量的每块起始的间距.
INOLDTYPE,原始数据类型.
OUTNEWTYPE,新构造的数据类型.
这个函数要比MPI_TYPE_VECTOR更具普遍性,但是它的每一块是固定不变的,下面将介绍变块长度的数据类型如何构造.
MPITYPEINDEXEDCintMPI_Type_indexed(intcount,int*array_of_blocklengths,int*array_of_displacements,MPI_Datatypeoldtype,MPI_Datatype*newtype)FortranMPI_TYPE_INDEXED(COUNT,ARRAY_OF_BLOCKLENGTHS,ARRAY_OF_DISPLACEMENTS,OLDTYPE,NEWTYPE,IERROR)INTEGERCOUNT,ARRAY_OF_BLOCKLENGTHS(*),ARRAY_OF_DISPLACEMENTS(*),OLDTYPE,NEWTYPE,IERROR§5.
11用户定义的数据类型87参数说明INCOUNT,OLDTYPE的块的个数.
INARRAY_OF_BLOCKLENGTHS,每块长度数组.
INARRAY_OF_DISPLACEMENTS,用元素个数定义的每块的起始位移量数组.
INOLDTYPE,原始数据类型.
OUTNEWTYPE,新构造的数据类型.
这里ARRAY_OF_BLOCKLENGTHS和ARRAY_OF_DISPLACEMENTS是由OLDTYPE定义的元素个数来度量的,而且ARRAY_OF_DISPLACEMENTS是一个相对值.
前面给出的MPI_TYPE_VECTOR函数是此函数的特例,也即下面的二个调用产生相同的数据类型.
MPI_TYPE_VECTOR(COUNT,BLOCKLENGTH,STRIDE,OLDTYPE,NEWTYPE,IERROR)和MPI_TYPE_INDEXED(COUNT,LENGTHS,DISPLACEMENTS,OLDTYPE,NEWTYPE,IERROR)其中LENGTHS(i)=BLOCKLENGTH,DISPLACEMENTS(i)=i*STRIDE,i=0,.
.
.
,COUNT-1.
在这里我们给出的函数MPI_TYPE_INDEXED的位移量是由OLDTYPE定义的元素个数确定的,下面将要给出的函数的位移量是由字节数确定的.
MPITYPEHINDEXEDCintMPI_Type_hindexed(intcount,int*array_of_blocklengths,int*array_of_displacements,MPI_Datatypeoldtype,MPI_Datatype*newtype)FortranMPI_TYPE_HINDEXED(COUNT,ARRAY_OF_BLOCKLENGTHS,ARRAY_OF_DISPLACEMENTS,OLDTYPE,NEWTYPE,IERROR)INTEGERCOUNT,ARRAY_OF_BLOCKLENGTHS(*),ARRAY_OF_DISPLACEMENTS(*),OLDTYPE,NEWTYPE,IERROR参数说明INCOUNT,OLDTYPE的块的个数.
INARRAY_OF_BLOCKLENGTHS,每块长度数组.
INARRAY_OF_DISPLACEMENTS,用字节数定义的每块的起始位移量数组.
INOLDTYPE,原始数据类型.
OUTNEWTYPE,新构造的数据类型.
以上给出的关于数据类型构造函数都是单一的原始数据类型,下面将给出一个最一般的构造新数据类型的函数.
MPITYPESTRUCTCintMPI_Type_struct(intcount,int*array_of_blocklengths,int*array_of_displacements,MPI_Datatype*array_of_types,MPI_Datatype*newtype)FortranMPI_TYPE_STRUCT(COUNT,ARRAY_OF_BLOCKLENGTHS,ARRAY_OF_DISPLACEMENTS,ARRAY_OF_TYPES,NEWTYPE,IERROR)INTEGERCOUNT,ARRAY_OF_BLOCKLENGTHS(*),ARRAY_OF_DISPLACEMENTS(*),ARRAY_OF_TYPES(*),NEWTYPE,IERROR88第五章消息传递并行程序设计平台MPI参数说明INCOUNT,块的个数.
INARRAY_OF_BLOCKLENGTHS,每块长度数组.
INARRAY_OF_DISPLACEMENTS,用字节数定义的每块的起始位移量数组.
INARRAY_OF_TYPES,原始数据类型数组.
OUTNEWTYPE,新构造的数据类型.
在所有MPI的构造数据类型中,MPI_TYPE_STRUCT是使用最广泛的一个函数,正确使用此函数在实际应用中是非常重要的.
到目前为此,所有构造数据类型的函数都已经介绍了,但是要使新构造的数据类型真正能够使用,还需要有以下的一些函数配合才行,因此关于如何使用这些数据类型构造函数的例子将在其它一些有关的函数介绍完之后给出.
MPITYPECOMMITCintMPI_Type_commit(MPI_Datatype*datatype)FortranMPI_TYPE_COMMIT(DATATYPE,IERROR)INTEGERDATATYPE,IERROR参数说明INOUTDATATYPE,准备提交的数据类型.
在使用自定义的MPI数据类型之前,必须调用MPI_TYPE_COMMIT.
MPITYPEFREECintMPI_Type_free(MPI_Datatype*datatype)FortranMPI_TYPE_FREE(DATATYPE,IERROR)INTEGERDATATYPE,IERROR参数说明INOUTDATATYPE,释放不用的数据类型.
MPITYPEEXTENTCintMPI_Type_extent(MPI_Datatypedatatype,MPI_Aint*extent)FortranMPI_TYPE_EXTENT(DATATYPE,EXTENT,IERROR)INTEGERDATATYPE,EXTENT,IERROR参数说明INDATATYPE,数据类型.
OUTEXTENT,数据类型拓展后的字节数.
MPIADDRESSCintMPI_Address(void*location,MPI_Aint*address)FortranMPI_ADDRESS(LOCATION,ADDRESS,IERROR)LOCATION(*)INTEGERADDRESS,IERROR§5.
11用户定义的数据类型89参数说明INLOCATION,在存储器中的变量表示.
OUTADDRESS,变量在存储器中的地址.
为了说明如何构造数据类型和MPI_TYPE_EXTENT函数的含义以及MPI_ADDRESS的使用,我们给出一个具体的例子加以解释.
例5.
11.
1.
假设在C语言中定义了一个结构,其中包括一个双精度数和一个字符,我们要定义这样结构的MPI数据类型.
#include#include"mpi.
h"typedefstruct{doublevalue;charstr;}data;voidmain(argc,argv)intargc;char**argv;{intp,myid,lens[3]={1,1,1},i;MPI_Commmycomm;datatst[3];MPI_Datatypenew,type[3]={MPI_DOUBLE,MPI_CHAR,MPI_UB};MPI_Aintdisp[3],size;MPI_Statusstatus;MPI_Init(&argc,&argv);MPI_Comm_dup(MPI_COMM_WORLD,&mycomm);MPI_Comm_rank(mycomm,&myid);MPI_Comm_size(mycomm,&p);if(myid==0)for(i=0;i=0;i--)disp[i]-=disp[0];MPI_Type_struct(3,lens,disp,type,&new);MPI_Type_commit(&new);MPI_Type_extent(new,&size);if(myid==0)printf("\nTheExtentofStructDataTypeis%d.
\n",size);if(myid==0)MPI_Ssend(tst,3,new,1,1,mycomm);elseif(myid==1)MPI_Recv(tst,3,new,0,1,mycomm,&status);MPI_Type_free(&new);if(myid==1)printf("\nThevaluesare%fand%c\n",tst[1].
value,tst[1].
str);MPI_Comm_free(&mycomm);MPI_Finalize();}在这个C程序中,我们定义了一个新的MPI数据类型,而且还使用它进行消息传递.
这里我们使用了MPI中提供的伪数据类型(pseudodatatype)MPI_UB,在MPI中还有一个伪数据类型MPI_LB.
这样定义的MPI的新数据类型+new+在调用函数MPI_TYPE_EXTENT之后返回值与机器关于双精度数是按4位还是8位对齐有关,如果是按4位对齐,返回值为12,否则为16.
但是如果我们使用下面的函数MPI_TYPE_SIZE,它的返回值是9.
90第五章消息传递并行程序设计平台MPIMPITYPESIZECintMPI_Type_size(MPI_Datatypedatatype,int*size)FortranMPI_TYPE_SIZE(DATATYPE,SIZE,IERROR)INTEGERDATATYPE,SIZE,IERROR参数说明INDATATYPE,数据类型.
OUTSIZE,数据类型的字节数.
MPIGETELEMENTSCintMPI_Get_elements(MPI_Status*status,MPI_Datatypedatatype,int*count)FortranMPI_GET_ELEMENTS(STATUS,DATATYPE,COUNT,IERROR)INTEGERSTATUS(MPI_STATUS_SIZE),DATATYPE,COUNT,IERROR参数说明INSTATUS,接收消息的状态.
INDATATYPE,数据类型.
OUTCOUNT,MPI基本数据类型的个数.
该函数和MPI_GET_COUNT的含意是不同的,MPI_GET_COUNT是接收到的数据类型个数,而这里给出的MPI_GET_ELEMENTS得到的是MPI基本数据类型的个数.
如果这里的数据类型是MPI基本数据类型,如C语言中的MPI_INT(Fortran中MPI_INTEGER)等,则MPI_GET_COUNT和MPI_GET_ELEMENTS等价.
我们使用下面的一段程序来说明它们的差别.
.
.
.
.
.
.
callmpi_type_contiguous(2,mpi_integer,type,ierr)callmpi_type_commit(type,ierr)if(myid.
eq.
0.
and.
p.
gt.
1)thencallmpi_send(ia,3,mpi_integer,1,1,mycomm,ierr)elseif(myid.
eq.
1)thencallmpi_recv(ia,2,type,0,1,mycomm,status,ierr)callmpi_get_count(status,type,count1,ierr)callmpi_get_elements(status,type,count2,ierr)endif.
.
.
.
.
.
一般来说,在我们自己写的程序中很少会出现发送与接收时数据类型是不同的,因此在大部分实际应用中,只需会使用MPI_GET_COUNT就不会遇到任何困难.
在这一节中,我们已经介绍了如何构造自定义的MPI数据类型,通常在定义这些数据类型时使用的是相对位置移动量,如果使用绝对地址来构造新的MPI数据类型,在使用这个新的数据类型进行消息传递时,MPI定义了一个称之为MPI_BOTTOM的消息缓冲区常数,以供MPI的发送与接收函数使用,具体的使用方式将在以后的例子中给出.
§5.
12MPI的数据打包与拆包这里的数据打包与拆包和PVM(ParallelVirtualMachine)有相同的性质,在MPI的绝大多数应用中使用用户自定义的MPI数据类型就可以完成许多数据一起发送与接收的目的.
但是对于复杂的情况使用数据打包可能会收到好处.
§5.
12MPI的数据打包与拆包91MPIPACKCintMPI_Pack(void*inbuf,intincount,MPI_Datatypedatatype,void*outbuf,intoutsize,int*position,MPI_Commcomm)FortranMPI_PACK(INBUF,INCOUNT,DATATYPE,OUTBUF,OUTSIZE,POSITION,COMM,IERROR)INBUF(*),OUTBUF(*)INTEGERINCOUNT,DATATYPE,OUTSIZE,POSITION,COMM,IERROR参数说明ININBUF,准备打在包中的输入缓冲区.
ININCOUNT,输入元素的个数.
INDATATYPE,数据类型.
OUTOUTBUF,打包缓冲区.
INOUTSIZE,用字节数定义的打包缓冲区的大小.
INOUTPOSITION,打包缓冲区的当前位置.
INCOMM,打包消息的通信子.
MPIUNPACKCintMPI_Unpack(void*inbuf,intinsize,int*position,void*outbuf,intoutcount,MPI_Datatypedatatype,MPI_Commcomm)FortranMPI_UNPACK(INBUF,INSIZE,POSITION,OUTBUF,OUTCOUNT,DATATYPE,COMM,IERROR)INBUF(*),OUTBUF(*)INTEGERINSIZE,POSITION,OUTCOUNT,DATATYPE,COMM,IERROR参数说明ININBUF,拆包缓冲区.
ININSIZE,用字节数定义的拆包缓冲区的大小.
INOUTPOSITION,打包缓冲区的当前位置.
OUTOUTBUF,输出数据缓冲区.
INOUTCOUNT,输出数据缓冲区元素的个数.
INDATATYPE,数据类型.
INCOMM,通信子.
采用这种打包方式进行消息传递时,数据类型一定要使用MPI_PACKED.
为了能够知道一组数据在打包时需要多大的缓冲区,MPI提供了如下的函数:MPIPACKSIZECintMPI_Pack_size(intincount,MPI_Datatypedatatype,MPI_Commcomm,int*size)FortranMPI_PACK_SIZE(INCOUNT,DATATYPE,COMM,SIZE,IERROR)INTEGERINCOUNT,DATATYPE,COMM,SIZE,IERROR92第五章消息传递并行程序设计平台MPI参数说明ININCOUNT,准备打包的元素个数.
INDATATYPE,数据类型.
INCOMM,通信子.
OUTSIZE,按字节数定义的打包需要的缓冲区大小.
此函数返回的值SIZE要比原始数据本身所占有的字节数要大,这是因为在包中还有一些其它的信息.
以下我们给出使用MPI_PACK和MPI_UNPACK的一个例子,来结束关于用户自定义的数据类型部分.
例5.
12.
1.
假设我们要把一个整数数组和一个双精度数组从第0个进程传送到第1个进程,尽管我们可以使用用户自定义数据类型来达到一次传送的目的,但我们也可用打包和拆包的方式实现上述要求.
具体实现如下:programpackinclude'mpif.
h'cintegermaxbuf,lenparameter(maxbuf=200,len=10)integermyid,p,mycomm,ierr,status(mpi_status_size,2),&ia(len),count1,count2,i,posreal*8a(len)characterbuf(maxbuf)ccallmpi_init(ierr)callmpi_comm_dup(mpi_comm_world,mycomm,ierr)callmpi_comm_rank(mycomm,myid,ierr)callmpi_comm_size(mycomm,p,ierr)print*,'Process',myid,'of',p,'isrunning'cif(myid.
eq.
0)thendo10i=1,lenia(i)=ia(i)=dble(i)10continueendifcallmpi_pack_size(len,mpi_integer,mycomm,count1,ierr)callmpi_pack_size(len,mpi_double_precision,mycomm,count2,i)print*,'Thepacksizeof10integerandreal8are:',&count1,count2pos=1if(myid.
eq.
0.
and.
p.
gt.
1)thencallmpi_pack(ia,len,mpi_integer,buf,maxbuf,pos,&mycomm,ierr)callmpi_pack(a,len,mpi_double_precision,buf,maxbuf,pos,&mycomm,ierr)callmpi_send(buf,pos-1,mpi_packed,1,1,mycomm,ierr)elseif(myid.
eq.
1)thencallmpi_recv(buf,maxbuf,mpi_packed,0,1,mycomm,&status,ierr)callmpi_unpack(buf,maxbuf,pos,ia,len,mpi_integer,&mycomm,ierr)§5.
13聚合通信93callmpi_unpack(buf,maxbuf,pos,a,len,mpi_double_precision,&mycomm,ierr)print*,'Thereceivedvaluesare:',a(1),a(2),a(3)endifcallmpi_comm_free(mycomm,ierr)callmpi_finalize(ierr)stopend§5.
13聚合通信聚合通信(collectivecommunication)是指多个进程(通常大于2)之间的通信.
在这一节中,介绍三方面的内容,障碍同步(barriersynchronization)、全局通信函数(globalcommunicationfunctions)和全局归约操作(globalreductionoperations).
§5.
14障碍同步MPIBARRIERCintMPI_Barrier(MPI_Commcomm)FortranMPI_BARRIER(COMM,IERROR)INTEGERCOMM,IERROR参数说明INCOMM,通信子.
这是MPI提供的唯一的一个同步函数,当COMM中的所有进程都执行这个函数时才返回,如果有一个进程没有执行此函数,其余进程将处于等待状态.
在执行完这个函数之后,所有进程将同时执行其后的任务.
§5.
15单点与多点通信函数MPIBCASTCintMPI_Bcast(void*buffer,intcount,MPI_Datatypedatatype,introot,MPI_Commcomm)FortranMPI_BCAST(BUFFER,COUNT,DATATYPE,ROOT,COMM,IERROR)BUFFER(*)INTEGERCOUNT,DATATYPE,ROOT,COMM,IERROR参数说明INOUTBUFFER,缓冲区的首地址.
INCOUNT,缓冲区中元素的个数.
INDATATYPE,缓冲区的数据类型.
INROOT,以它为源进行广播的进程编号.
INCOMM,通信子.
94第五章消息传递并行程序设计平台MPI使用此函数时必须注意在COMM中的所有进程都执行此函数,如果进程编号ROOT,则表示把此进程中的BUFFER的内容广播到COMM所有其它的进程中.
此函数是并行程序中经常出现的,因此是个必须很好掌握的MPI通信函数.
MPIGATHERCintMPI_Gather(void*sendbuf,intsendcount,MPI_Datatypesendtype,void*recvbuf,intrecvcount,MPI_Datatyperecvtype,introot,MPI_Commcomm)FortranMPI_GATHER(SENDBUF,SENDCOUNT,SENDTYPE,RECVBUF,RECVCOUNT,RECVTYPE,ROOT,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERSENDCOUNT,SENDTYPE,RECVCOUNT,RECVTYPE,ROOT,COMM,IERROR参数说明INSENDBUF,发送缓冲区的首地址.
INSENDCOUNT,要发送的元素的个数.
INSENDTYPE,发送缓冲区的数据类型.
OUTRECVBUF,接收缓冲区的首地址.
INRECVCOUNT,接收每个进程中数据元素的个数.
INRECVTYPE,接收缓冲区的数据类型.
INROOT,接收进程的编号.
INCOMM,通信子.
这是MPI提供的聚集(gather)函数,其作用是把COMM中所有进程(包括ROOT)的数据聚集到ROOT进程中,并且自动按进程编号顺序存放在接收缓冲区RECVBUF中.
对于非ROOT进程,忽略接收缓冲区RECVBUF.
MPIGATHERVCintMPI_Gatherv(void*sendbuf,intsendcount,MPI_Datatypesendtype,void*recvbuf,int*recvcounts,int*displs,MPI_Datatyperecvtype,introot,MPI_Commcomm)FortranMPI_GATHERV(SENDBUF,SENDCOUNT,SENDTYPE,RECVBUF,RECVCOUNTS,DISPLS,RECVTYPE,ROOT,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERSENDCOUNT,SENDTYPE,RECVCOUNTS(*),DISPLS(*),RECVTYPE,ROOT,COMM,IERROR§5.
15单点与多点通信函数95参数说明INSENDBUF,发送缓冲区的首地址.
INSENDCOUNT,要发送的元素的个数.
INSENDTYPE,发送缓冲区的数据类型.
OUTRECVBUF,接收缓冲区的首地址.
INRECVCOUNTS,接收每个进程数据元素的个数整数数组.
INDISPLS,接收数据存放的位置数组.
INRECVTYPE,接收缓冲区的数据类型.
INROOT,接收进程的编号.
INCOMM,通信子.
这个函数是MPI_GATHER的扩充,它允许从不同的进程中接收不同长度的消息,而且接收到的消息可以存放在接收缓冲区的不同位置,因此使用要比MPI_GATHER更灵活.
下面我们给出一个例子来说明如何使用MPI的聚集函数:例5.
15.
1.
假设在每个进程上都有一组数据,现要将它们收集到进程编号为ROOT的进程中,并按进程编号的顺序存放,则有如下的程序实现此功能.
programgatherinclude'mpif.
h'cintegermaxbuf,len,mpparameter(maxbuf=200,len=10,mp=5)integermyid,p,mycomm,ierr,root,ia(len),iga(maxbuf),i,&displs(mp),counts(mp)ccallmpi_init(ierr)callmpi_comm_dup(mpi_comm_world,mycomm,ierr)callmpi_comm_rank(mycomm,myid,ierr)callmpi_comm_size(mycomm,p,ierr)print*,'Process',myid,'of',p,'isrunning'cdo10i=1,lenia(i)=i+myid*len10continuecdo20i=1,pdispls(i)=20*i-20counts(i)=len20continueroot=0ccallmpi_gather(ia,len,mpi_integer,iga,len,mpi_integer,&root,mycomm,ierr)cif(myid.
eq.
root)thenprint*,'Thegathervaluesare:',&iga(1),iga(len+1),iga(2*len+1)endifccallmpi_gatherv(ia,len,mpi_integer,iga,counts,displs,96第五章消息传递并行程序设计平台MPI&mpi_integer,root,mycomm,ierr)cif(myid.
eq.
root)thenprint*,'Thegathervvaluesare:',&iga(2),iga(2*len+2),iga(4*len+2)endifccallmpi_comm_free(mycomm,ierr)callmpi_finalize(ierr)stopend这个程序的输出值应分别为:1,11,21和2,12,22.
MPI_GATHER接收到的数据在接收缓冲区中是连续存放的,而MPI_GATHERV接收到的数据在接收缓冲区中是不连续存放的,在我们这里给出的例子中,由于COUNTS的所有分量是相同的,可以看成是每隔步长为20的位置存放一个新接收的数据.
MPISCATTERCintMPI_Scatter(void*sendbuf,intsendcount,MPI_Datatypesendtype,void*recvbuf,intrecvcount,MPI_Datatyperecvtype,introot,MPI_Commcomm)FortranMPI_SCATTER(SENDBUF,SENDCOUNT,SENDTYPE,RECVBUF,RECVCOUNT,RECVTYPE,ROOT,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERSENDCOUNT,SENDTYPE,RECVCOUNT,RECVTYPE,ROOT,COMM,IERROR参数说明INSENDBUF,发送缓冲区的首地址.
INSENDCOUNT,发送到每个进程中的元素的个数.
INSENDTYPE,发送缓冲区的数据类型.
OUTRECVBUF,接收缓冲区的首地址.
INRECVCOUNT,接收数据元素的个数.
INRECVTYPE,接收缓冲区的数据类型.
INROOT,发送进程的编号.
INCOMM,通信子.
这是MPI提供的散布(scatter)函数,此函数是MPI_GATHER的逆操作.
MPISCATTERVCintMPI_Scatterv(void*sendbuf,int*sendcounts,int*displs,MPI_Datatypesendtype,void*recvbuf,intrecvcount,MPI_Datatyperecvtype,introot,MPI_Commcomm)FortranMPI_SCATTERV(SENDBUF,SENDCOUNTS,DISPLS,SENDTYPE,RECVBUF,RECVCOUNT,RECVTYPE,ROOT,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERSENDCOUNTS(*),DISPLS(*),SENDTYPE,RECVCOUNT,RECVTYPE,ROOT,COMM,IERROR§5.
16多点与多点通信函数97参数说明INSENDBUF,发送缓冲区的首地址.
INSENDCOUNTS,发送到每个进程中的元素的个数数组.
INDISPLS,发送到每个进程中的数据起始位移.
INSENDTYPE,发送缓冲区的数据类型.
OUTRECVBUF,接收缓冲区的首地址.
INRECVCOUNT,接收数据元素的个数.
INRECVTYPE,接收缓冲区的数据类型.
INROOT,发送进程的编号.
INCOMM,通信子.
§5.
16多点与多点通信函数MPIALLGATHERCintMPI_Allgather(void*sendbuf,intsendcount,MPI_Datatypesendtype,void*recvbuf,intrecvcount,MPI_Datatyperecvtype,MPI_Commcomm)FortranMPI_ALLGATHER(SENDBUF,SENDCOUNT,SENDTYPE,RECVBUF,RECVCOUNT,RECVTYPE,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERSENDCOUNT,SENDTYPE,RECVCOUNT,RECVTYPE,COMM,IERROR参数说明INSENDBUF,发送缓冲区的首地址.
INSENDCOUNT,要发送的元素的个数.
INSENDTYPE,发送缓冲区的数据类型.
OUTRECVBUF,接收缓冲区的首地址.
INRECVCOUNT,接收每个进程中数据元素的个数.
INRECVTYPE,接收缓冲区的数据类型.
INCOMM,通信子.
MPIALLGATHERVCintMPI_Allgatherv(void*sendbuf,intsendcount,MPI_Datatypesendtype,void*recvbuf,intrecvcount,MPI_Datatyperecvtype,MPI_Commcomm)FortranMPI_ALLGATHERV(SENDBUF,SENDCOUNT,SENDTYPE,RECVBUF,RECVCOUNT,RECVTYPE,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERSENDCOUNT,SENDTYPE,RECVCOUNT,RECVTYPE,COMM,IERROR98第五章消息传递并行程序设计平台MPI参数说明INSENDBUF,发送缓冲区的首地址.
INSENDCOUNT,要发送的元素的个数.
INSENDTYPE,发送缓冲区的数据类型.
OUTRECVBUF,接收缓冲区的首地址.
INRECVCOUNT,接收每个进程中数据元素的个数.
INRECVTYPE,接收缓冲区的数据类型.
INCOMM,通信子.
此二函数是每个进程都进行数据聚集操作,而MPI_GATHER和MPI_GATHERV只是其中的一个进程进行数据的聚集.
这里每个进程都发送消息也都接收消息.
MPIALLTOALLCintMPI_Alltoall(void*sendbuf,intsendcount,MPI_Datatypesendtype,void*recvbuf,intrecvcount,MPI_Datatyperecvtype,MPI_Commcomm)FortranMPI_ALLTOALL(SENDBUF,SENDCOUNT,SENDTYPE,RECVBUF,RECVCOUNT,RECVTYPE,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERSENDCOUNT,SENDTYPE,RECVCOUNT,RECVTYPE,COMM,IERROR参数说明INSENDBUF,发送缓冲区的首地址.
INSENDCOUNT,要发送的元素的个数.
INSENDTYPE,发送缓冲区的数据类型.
OUTRECVBUF,接收缓冲区的首地址.
INRECVCOUNT,接收每个进程中数据元素的个数.
INRECVTYPE,接收缓冲区的数据类型.
INCOMM,通信子.
MPIALLTOALLVCintMPI_Alltoallv(void*sendbuf,int*sendcounts,int*sdispls,MPI_Datatypesendtype,void*recvbuf,int*recvcounts,int*rdispls,MPI_Datatyperecvtype,MPI_Commcomm)FortranMPI_ALLTOALLV(SENDBUF,SENDCOUNTS,SDISPLS,SENDTYPE,RECVBUF,RECVCOUNTS,RDISPLS,RECVTYPE,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERSENDCOUNTS(*),SDISPLS(*),SENDTYPE,RECVCOUNTS(*),RDISPLS(*),RECVTYPE,COMM,IERROR§5.
16多点与多点通信函数99参数说明INSENDBUF,发送缓冲区的首地址.
INSENDCOUNTS,发送到每个进程中的元素的个数数组.
INSDISPLS,发送到每个进程中的数据起始位移.
INSENDTYPE,发送缓冲区的数据类型.
OUTRECVBUF,接收缓冲区的首地址.
INRECVCOUNTS,接收数据元素的个数数组.
INRDISPLS,接收进程中存放数据的起始位移.
INRECVTYPE,接收缓冲区的数据类型.
INCOMM,通信子.
这二个函数等价于在所有进程中的数据进行散布.
为了便于理解和正确使用此函数,我们还是给出一个具体例子.
programallallinclude'mpif.
h'integermaxbuf,len,mpparameter(maxbuf=200,len=10,mp=5)integermyid,p,mycomm,ierr,igb(maxbuf),iga(maxbuf),i,&sdispls(mp),scounts(mp),rdispls(mp),rcounts(mp)ccallmpi_init(ierr)callmpi_comm_dup(mpi_comm_world,mycomm,ierr)callmpi_comm_rank(mycomm,myid,ierr)callmpi_comm_size(mycomm,p,ierr)print*,'Process',myid,'of',p,'isrunning'cdo10i=1,maxbufiga(i)=i+maxbuf*myid10continuedo20i=1,psdispls(i)=20*i-20rdispls(i)=15*i-15scounts(i)=lenrcounts(i)=len20continuecallmpi_alltoall(iga,len,mpi_integer,igb,len,mpi_integer,&mycomm,ierr)print*,'Thealltoallvaluesare:',&igb(1),igb(len+1),igb(2*len+1),'onProc.
',myidccallmpi_alltoallv(iga,scounts,sdispls,mpi_integer,igb,&rcounts,rdispls,mpi_integer,mycomm,ierr)print*,'Thealltoallvvaluesare:',&igb(1),igb(16),igb(31),'onProc.
',myidccallmpi_comm_free(mycomm,ierr)callmpi_finalize(ierr)stopend100第五章消息传递并行程序设计平台MPI§5.
17全局归约操作MPIREDUCECintMPI_Reduce(void*sendbuf,void*recvbuf,intcount,MPI_Datatypedatatype,MPI_Opop,introot,MPI_Commcomm)FortranMPI_REDUCE(SENDBUF,RECVBUF,COUNT,DATATYPE,OP,ROOT,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERCOUNT,DATATYPE,OP,ROOT,COMM,IERROR参数说明INSENDBUF,发送缓冲区的首地址.
OUTRECVBUF,接收缓冲区的首地址.
INCOUNT,缓冲区中元素的个数.
INDATATYPE,缓冲区的数据类型.
INOP,何种归约操作.
INROOT,接收进程编号.
INCOMM,通信子.
这个函数是按分量进行归约操作的,其运算由OP来确定,最后的结果放在ROOT进程中,其它进程中的RECVBUF不起作用.
在MPI中规定了一些允许的操作如下:操作名意义MPI_MAX求最大MPI_MIN求最小MPI_SUM求和MPI_PROD求积MPI_LAND逻辑与MPI_BAND按位与MPI_LOR逻辑或MPI_BOR按位或MPI_LXOR逻辑与或MPI_BXOR按位与或MPI_MAXLOC求最大和位置MPI_MINLOC求最小和位置这些运算是有数据类型要求的,首先对数据类型进行分类:Cinteger:MPI_INT,MPI_LONG,MPI_SHORT,MPI_UNSIGNED_SHORT,MPI_UNSIGNED,MPI_UNSIGNED_LONGFortraninteger:MPI_INTEGERFloatingpoint:MPI_DOUBLE,MPI_REAL,MPI_DOUBLE_PRECISION,MPI_LONG_DOUBLELogical:MPI_LOGICALComplex:MPI_COMPLEXByte:MPI_BYTE现在对每种操作允许的数据类型规定如下:§5.
17全局归约操作101OP允许的数据类型MPI_MAX,MPI_MINCinteger,Fortraninteger,FloatingpointMPI_SUM,MPI_PRODCinteger,Fortraninteger,Floatingpoint,ComplexMPI_LAND,MPI_LOR,MPI_LXORCinteger,LogicalMPI_BAND,MPI_BOR,MPI_BXORCinteger,Fortraninteger,Byte关于MPI_MAXLOC和MPI_MINLOC,MPI对Fortran程序和C程序使用的复合数据类型规定如下:Fortran程序复合数据类型类型描述MPI_2REALpairofREALsMPI_2DOUBLE_PRECISIONpairofDOUBLEPRECISIONsMPI_2INTEGERpairofINTEGERsC程序复合数据类型类型描述MPI_FLOAT_INTfloatandintMPI_DOUBLE_INTdoubleandintMPI_LONG_INTlongandintMPI_SHORT_INTshortandintMPI_LONG_DOUBLE_INTlongdoubleandintMPI_2INTpairofints例5.
17.
1.
假设在每个处理机中有一个数,我们要在这些数中找一个最大的,并确定这个最大数在哪个处理机中,则可用如下的程序:programreduceinclude'mpif.
h'cintegermyid,p,mycomm,ierr,m,n,root,pair(2),answer(2)ccallmpi_init(ierr)callmpi_comm_dup(mpi_comm_world,mycomm,ierr)callmpi_comm_rank(mycomm,myid,ierr)callmpi_comm_size(mycomm,p,ierr)print*,'Process',myid,'of',p,'isrunning'root=0m=myidcallmpi_reduce(m,n,1,mpi_integer,mpi_max,root,mycomm,ierr)cif(myid.
eq.
root)print*,'Themaxmumvalueis',npair(1)=mod(myid+1,p)pair(2)=myidcallmpi_reduce(pair,answer,1,mpi_2integer,mpi_maxloc,root,&mycomm,ierr)if(myid.
eq.
root)print*,'Themaxmumvalueis',answer(1),&'onprocess',answer(2)c102第五章消息传递并行程序设计平台MPIcallmpi_comm_free(mycomm,ierr)callmpi_finalize(ierr)stopendMPIALLREDUCECintMPI_Allreduce(void*sendbuf,void*recvbuf,intcount,MPI_Datatypedatatype,MPI_Opop,MPI_Commcomm)FortranMPI_ALLREDUCE(SENDBUF,RECVBUF,COUNT,DATATYPE,OP,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERCOUNT,DATATYPE,OP,COMM,IERROR参数说明INSENDBUF,发送缓冲区的首地址.
OUTRECVBUF,接收缓冲区的首地址.
INCOUNT,缓冲区中元素的个数.
INDATATYPE,缓冲区的数据类型.
INOP,何种归约操作.
INCOMM,通信子.
此函数和MPI_REDUCE的意思是相同的,只是最后结果在所有的进程中.
MPIREDUCESCATTERCintMPI_Reduce_scatter(void*sendbuf,void*recvbuf,int*recvcounts,MPI_Datatypedatatype,MPI_Opop,MPI_Commcomm)FortranMPI_REDUCE_SCATTER(SENDBUF,RECVBUF,RECVCOUNTS,DATATYPE,OP,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERRECVCOUNTS(*),DATATYPE,OP,COMM,IERROR参数说明INSENDBUF,发送缓冲区的首地址.
OUTRECVBUF,接收缓冲区的首地址.
INRECVCOUNTS,元素个数的数组.
INDATATYPE,缓冲区的数据类型.
INOP,何种归约操作.
INCOMM,通信子.
此函数相当于先做MPI_REDUCE,然后在做MPI_SCATTER.
例5.
17.
2.
设我们使用p个进程计算c=Ab,其中A是m阶矩阵,b是m-维向量.
假设m=np(如不能整除把剩余部分放在最后一个进程中),在每个进程中存放A的n列和对应的b的n个分量,现要在每个进程中得到与b对应的c或是c的全部,则我们可以用下面的程序实现:§5.
17全局归约操作103programredsctinclude'mpif.
h'cintegerlda,cols,maxnpparameter(lda=100,cols=100,maxnp=5)integermyid,p,mycomm,ierr,m,n,counts(maxnp)reala(lda,cols),b(cols),c(lda),sum(lda)logicalsctccallmpi_init(ierr)callmpi_comm_dup(mpi_comm_world,mycomm,ierr)callmpi_comm_rank(mycomm,myid,ierr)callmpi_comm_size(mycomm,p,ierr)print*,'Process',myid,'of',p,'isrunning'if(p.
gt.
maxnp)goto999m=100sct=.
false.
csct=.
true.
callinitab(m,myid,p,n,a,lda,b,sum)calllocmv(m,n,a,lda,b,sum)cifsct=true,callreduce_scatter,otherwizecallallreduceif(sct)thencallmpi_allgather(n,1,mpi_integer,counts,1,&mpi_integer,mycomm,ierr)callmpi_reduce_scatter(sum,c,counts,mpi_real,mpi_sum,&mycomm,ierr)elsecallmpi_allreduce(sum,c,m,mpi_real,mpi_sum,mycomm,ierr)endifcprint*,'Thevaluesofcare',c(1),c(2)c999callmpi_comm_free(mycomm,ierr)callmpi_finalize(ierr)stopendsubroutineinitab(m,myid,np,k,a,lda,b,sum)integerm,myid,np,k,lda,i,j,lreala(lda,*),b(*),sum(*)do10i=1,m10sum(i)=0.
0k=m/npl=myid*kif(myid.
eq.
np-1)k=m-ldo20i=1,k20b(i)=1.
0do40j=1,kdo40i=1,m40a(i,j)=real(i+l+j)returnend104第五章消息传递并行程序设计平台MPIsubroutinelocmv(m,n,a,lda,b,c)integerm,n,lda,i,jreala(lda,*),b(*),c(*)do20j=1,ndo20i=1,m20c(i)=c(i)+a(i,j)*b(j)returnendMPISCANCintMPI_Scan(void*sendbuf,void*recvbuf,intcount,MPI_Datatypedatatype,MPI_Opop,MPI_Commcomm)FortranMPI_SCAN(SENDBUF,RECVBUF,COUNT,DATATYPE,OP,COMM,IERROR)SENDBUF(*),RECVBUF(*)INTEGERCOUNT(*),DATATYPE,OP,COMM,IERROR参数说明INSENDBUF,发送缓冲区的首地址.
OUTRECVBUF,接收缓冲区的首地址.
INCOUNT,元素个数.
INDATATYPE,缓冲区的数据类型.
INOP,何种归约操作.
INCOMM,通信子.
此函数是MPI_REDUCE在不同进程中的重复使用,进程i中得到的是在由编号为{0,1,i}的进程构成的进程组中使用MPI_REDUCE结果,是一个并不常用的函数.
在MPI中提供了用户自定义的操作函数以满足不同需要,它由如下的几个函数来实现:MPIOPCREATECintMPI_Op_create(MPI_User_function*function,intcommute,MPI_Op*op)FortranMPI_OP_CREATE(FUNCTION,COMMUTE,OP,IERROR)EXTERNALFUNCTIONLOGICALCOMMUTEINTEGEROP,IERROR参数说明INFUNCTION,用户定义的函数.
INCOMMUTE,等于TRUE是可交换的,否则只是可结合的.
OUTOP,新的归约操作.
此函数定义了一个新的归约操作OP.
MPIOPFREECintMPI_Op_free(MPI_Op*op)FortranMPI_OP_FREE(OP,IERROR)INTEGEROP,IERROR§5.
17全局归约操作105参数说明INOP,归约操作.
此函数释放归约操作OP.
在MPI中,用户自定义操作的函数是有严格要求的,其形式如下:Cvoiduser_function(void*invec,void*inoutvec,int*len,MPI_Datatype*datatype)FortranFUNCTIONUSER_FUNCTION(INVEC,INOUTVEC,LEN,DATATYPE)INVEC(LEN),INOUTVEC(LEN)INTEGERLEN,DATATYPE下面的例子是用自定义的求和函数来实现全局操作的:programuserdefinclude'mpif.
h'integerlenparameter(len=100)integermyid,p,mycomm,ierr,m,root,myoprealx(len),y(len)externaluserfunccallmpi_init(ierr)callmpi_comm_dup(mpi_comm_world,mycomm,ierr)callmpi_comm_rank(mycomm,myid,ierr)callmpi_comm_size(mycomm,p,ierr)print*,'Process',myid,'of',p,'isrunning'm=100root=0callinitx(m,myid,x)callmpi_op_create(userfunc,.
true.
,myop,ierr)callmpi_reduce(x,y,m,mpi_real,myop,root,mycomm,ierr)callmpi_op_free(myop,ierr)if(myid.
eq.
root)&print*,'Thevaluesofanswerare',y(1),y(2),y(3)callmpi_comm_free(mycomm,ierr)callmpi_finalize(ierr)stopendsubroutineinitx(m,myid,x)integerm,myid,irealx(*)do10i=1,m10x(i)=real(i+myid)returnendsubroutineuserfunc(x,y,m,mpi_real)integermrealx(*),y(*)do10i=1,m10y(i)=y(i)+x(i)returnend106第五章消息传递并行程序设计平台MPI§5.
18进程拓扑结构进程拓扑结构是(域内)通信器的一个附加属性,它描述一个进程组各进程间的逻辑联接关系.
进程拓扑结构的使用一方面可以方便、简化一些并行程序的编制,另一方面可以帮助MPI系统更好地将进程映射到处理器以及组织通信的流向,从而获得更好的并行性能.
MPI的进程拓扑结构定义为一个无向图,图中结点(node)代表进程,而边(edge)代表进程间的联接.
MPI进程拓扑结构也被称为虚拟拓扑结构,因为它不一定对应处理器的物理联接.
MPI提供了一组函数用于创建各种进程拓扑结构.
应用问题中较为常见、也是较为简单的一类进程拓扑结构具有网格形式,这类结构中进程可以用迪卡尔坐标来标识,MPI中称这类拓扑结构为迪卡尔(Cartesian)拓扑结构,并且专门提供了一组函数对它们进行操作.
§5.
18.
1迪卡尔拓扑结构§5.
18.
1.
1创建迪卡尔拓扑结构CintMPI_Cart_create(MPI_Commcomm_old,intndims,int*dims,int*periods,intreorder,MPI_Comm*comm_cart)Fortran77MPI_CART_CREATE(COMM_OLD,NDIMS,DIMS,PERIODS,REORDER,+COMM_CART,IERR)INTEGERCOMM_OLD,NDIMS,DIMS(*),COMM_CART,IERRLOGICALPERIODS(*),REORDER该函数从一个通信器comm_old出发,创建一个具有迪卡尔拓扑结构的新通信器comm_cart.
ndims给出进程网格的维数.
数组dims给出每维中的进程数.
数组periods则说明进程在各个维上的联接是否具有周期性,即该维中第一个进程与最后一个进程是否相联,周期的迪卡尔拓扑结构也称为环面(torus)结构,periods[i]=true表明第i维是周期的,否则则是非周期的.
reorder指明是否允许在新通信器comm_cart中对进程进行重新排序.
在某些并行机上,根据处理器的物理联接方式及所要求的进程拓扑结构对进程重新排序有助于提高并行程序的性能.

comm_cart中各维的进程数之积必须不大于comm_old中的进程数,即:ndims1i=0dims[i]≤NPROCS其中NPROCS为comm_old的进程数.
当ndims1i=0dims[i]BUF(*)INTEGERFH,COUNT,DATATYPE,STATUS(MPI_STATUS_SIZE),+IERRINTEGER(KIND=MPI_OFFSET_KIND)OFFSETfh为文件句柄,offset为位移.
buf,count和datatype分别为数据的缓冲区地址、个数和类型.
status返回操作结果状态(与通信函数类似).
§6.
5.
2使用独立文件指针的阻塞型文件读写使用独立文件指针的阻塞型文件读写函数与使用显式位移的阻塞型文件读写函数的功能完全一样,只是文件位移由独立文件指针隐式设定.
这些函数的接口参数中比使用显式位移的函数少了一个offset参数,其它参数完全一样.
作为例子,我们给出MPI_FILE_READ的接口参数.
CintMPI_File_read(MPI_Filefh,void*buf,intcount,MPI_Datatypedatatype,MPI_Status*status)FortranMPI_FILE_READ(FH,BUF,COUNT,DATATYPE,STATUS,IERR)BUF(*)120第六章文件输入输出(MPI–IO)INTEGERFH,COUNT,DATATYPE,STATUS(MPI_STATUS_SIZE),+IERR§6.
5.
3使用共享文件指针的阻塞型文件读写使用共享文件指针的阻塞型文件读写函数的接口参数与使用独立文件指针的阻塞型文件读写函数的接口参数完全一样.
作为例子,这里给出MPI_FILE_READ_ORDERED的接口参数.
CintMPI_File_read_ordered(MPI_Filefh,void*buf,intcount,MPI_Datatypedatatype,MPI_Status*status)FortranMPI_FILE_READ_ORDERED(FH,BUF,COUNT,DATATYPE,STATUS,+IERR)BUF(*)INTEGERFH,COUNT,DATATYPE,STATUS(MPI_STATUS_SIZE),+IERR由于使用共享文件指针的文件操作函数中进程组的全部进程共同使用和修改同一个文件指针,因此这类操作非常类似于以文件为"根进程"的数据收集和散发,即它们相当于将进程组中各进程的数据块合并写入文件(收集)或读取文件中的数据并分发给各进程(散发).
当使用非聚合式函数MPI_FILE_READ_SHARED和MPI_FILE_WRITE_SHARED时,各进程从文件中读取或写入文件的数据块在文件中的相对位置是不确定的,而聚合式函数MPI_FILE_READ_ORDERED和MPI_FILE_WRITE_ORDERED则可确保这些数据块在文件中严格按进程序号排列.
§6.
5.
4非阻塞型文件读写函数每个阻塞型非聚合式文件读写函数都有一个对应的非阻塞型函数,由阻塞型函数的函数名中在READ或WRITE前面加I构成,如MPI_FILE_READ的非阻塞型函数为MPI_FILE_IREAD.
非阻塞型函数的接口参数中只需将对应的阻塞型函数的参数表中的status参数换成request,其它参数完全一样.
非阻塞型函数递交文件读或写的请求,在request中返回一个请求句柄,实际的读或写操作在后台进行.
非阻塞型文件读写函数返回的请求句柄与非阻塞型消息传递函数所返回的句柄的操作方式完全一样,即用户需在关闭文件前调用MPI_WAIT,MPI_TEST等函数来检查、等待操作的完成.
做为例子,下面列出MPI_IREAD_AT的参数.
CintMPI_File_iread_at(MPI_Filefh,MPI_Offsetoffset,void*buf,intcount,MPI_Datatypedatatype,MPI_Request*request)FortranMPI_FILE_IREAD_AT(FH,OFFSET,BUF,COUNT,DATATYPE,+REQUEST,IERR)BUF(*)INTEGERFH,COUNT,DATATYPE,REQUEST,IERRINTEGER(KIND=MPI_OFFSET_KIND)OFFSET§6.
6文件指针操作121§6.
5.
5分裂型文件读写函数MPI为每个阻塞型聚合式文件读写函数定义了一对分裂型函数,分别在阻塞型函数的函数名后面加_BEGIN和_END构成.
分裂型函数将文件读写操作分解成开始和结束两步,用户可以在开始和结束之间插入其它通信或计算,从而实现计算或通信与文件输入输出重叠进行.
这些函数的函数名及包含的接口参数如下:MPI_FILE_READ_AT_ALL_BEGIN(fh,offset,buf,count,datatype)MPI_FILE_READ_AT_ALL_END(fh,buf,status)MPI_FILE_WRITE_AT_ALL_BEGIN(fh,offset,buf,count,datatype)MPI_FILE_WRITE_AT_ALL_END(fh,buf,status)MPI_FILE_READ_ALL_BEGIN(fh,buf,count,datatype)MPI_FILE_READ_ALL_END(fh,buf,status)MPI_FILE_WRITE_ALL_BEGIN(fh,buf,count,datatype)MPI_FILE_WRITE_ALL_END(fh,buf,status)MPI_FILE_READ_ORDERED_BEGIN(fh,buf,count,datatype)MPI_FILE_READ_ORDERED_END(fh,buf,status)MPI_FILE_WRITE_ORDERED_BEGIN(fh,buf,count,datatype)MPI_FILE_WRITE_ORDERED_END(fh,buf,status)§6.
6文件指针操作§6.
6.
1独立文件指针操作§6.
6.
1.
1移动独立文件指针CintMPI_File_seek(MPI_Filefh,MPI_Offsetoffset,intwhence)FortranMPI_FILE_SEEK(FH,OFFSET,WHENCE,IERR)INTEGERFH,WHENCE,IERRINTEGER(KIND=MPI_OFFSET_KIND)OFFSET改变独立文件指针的位移.
参数whence可取为下列值:MPI_SEEK_SET—将指针的位移设为offsetMPI_SEEK_CUR—将指针的位移设为当前位移加上offsetMPI_SEEK_END—将指针的位移设为文件结尾加上offset§6.
6.
1.
2查询独立文件指针的当前位移CintMPI_File_get_position(MPI_Filefh,MPI_Offset*offset)122第六章文件输入输出(MPI–IO)FortranMPI_FILE_GET_POSITION(FH,OFFSET,IERR)INTEGERFH,IERRINTEGER(KIND=MPI_OFFSET_KIND)OFFSET在参数offset返回独立文件指针的位移.
§6.
6.
2共享文件指针操作§6.
6.
2.
1移动共享文件指针CintMPI_File_seek_shared(MPI_Filefh,MPI_Offsetoffset,intwhence)FortranMPI_FILE_SEEK_SHARED(FH,OFFSET,WHENCE,IERR)INTEGERFH,WHENCE,IERRINTEGER(KIND=MPI_OFFSET_KIND)OFFSET改变共享文件指针的位移.
参数whence可取为下列值:MPI_SEEK_SET—将指针的位移设为offsetMPI_SEEK_CUR—将指针的位移设为当前位移加上offsetMPI_SEEK_END—将指针的位移设为文件结尾加上offsetMPI_FILE_SEEK_SHARED是聚合型函数,进程组中所有进程必须同时调用并且提供相同的参数.
§6.
6.
2.
2查询共享文件指针的当前位移CintMPI_File_get_position_shared(MPI_Filefh,MPI_Offset*offset)FortranMPI_FILE_GET_POSITION_SHARED(FH,OFFSET,IERR)INTEGERFH,IERRINTEGER(KIND=MPI_OFFSET_KIND)OFFSET在参数offset返回共享文件指针的位移.
§6.
6.
3文件位移在文件中的绝对地址CintMPI_File_get_byte_offset(MPI_Filefh,MPI_Offsetoffset,MPI_Offset*disp)FortranMPI_FILE_GET_BYTE_OFFSET(FH,OFFSET,DISP,IERR)INTEGERFH,IERRINTEGER(KIND=MPI_OFFSET_KIND)OFFSET,DISP该函数将以etype为单位相对于当前文件视窗的位移(offset)换算成以字节为单位从文件开头计算的绝对地址(disp).
§6.
7不同进程对同一文件读写操作的相容性123§6.
7不同进程对同一文件读写操作的相容性当单个或多个进程同时对同一个文件进行访问时,MPI称这些访问是相容的,如果这些访问可以等效地被看成是以某种顺序依次进行的,即便它们的先后顺序是不确定的.
换言之,对同一个文件的多个访问是相容的,如果它们中的任一访问都不会在操作过程中间由于被另一个访问打断或干扰而影响到访问的结果.
MPI系统允许用户将对一个文件的访问设置成具有"原子性"(atomicity,意即"不可分的")来保证属于与该文件关联的进程组中的进程对该文件的访问的相容性.

但如果同一个文件分别被不同的进程组打开,并且两个进程组中对该文件的访问存在冲突,则用户必须通过在程序中调用MPI_FILE_SYNC函数以及同步函数(MPI_BARRIER)等来保证对文件访问的相容性与访问顺序.
§6.
7.
1设定文件访问的原子性CintMPI_File_set_atomicity(MPI_Filefh,intflag)FortranMPI_FILE_SET_ATOMICITY(FH,FLAG,IERR)INTEGERFH,IERRLOGICALFLAG该函数设定MPI是否需要保证打开文件的进程组中进程对该文件的访问的原子性.
当flag为true时,由MPI系统将保证文件访问的原子性从而保证属于(与该文件相关联的)同一进程组的进程对该文件的访问的相容性.
而当flag为false时,MPI不保证对文件访问的原子性,而需要用户通过其它途径来保证对文件的不同访问间的相容性.
MPI_FILE_SET_ATOMICITY是聚合型函数,进程组中所有进程必须同时调用并且提供相同的参数.
例6.
7.
1.
在文件的同一位置上一个进程写、另一个进程读.
INTEGERSTATUS(MPI_STATUS_SIZE),FH,A(10).
.
.
.
.
.
CALLMPI_FILE_OPEN(MPI_COMM_WORLD,"myfile",+MPI_MODE_RDWR+MPI_MODE_CREATE,+MPI_INFO_NULL,FH,IERR)CALLMPI_FILE_SET_VIEW(FH,0,MPI_INTEGER,MPI_INTEGER,+'native',MPI_INFO_NULL,IERR)CALLMPI_SET_ATOMICITY(FH,.
TRUE.
,IERR)IF(MYRANK.
EQ.
0)THENDOI=1,10A(I)=5ENDDOCALLMPI_FILE_WRITE_AT(FH,0,A,10,MPI_INTEGER,+STATUS,IERR)ELSEIF(MYRANK.
EQ.
1)THENCALLMPI_FILE_READ_AT(FH,0,A,10,MPI_INTEGER,+STATUS,IERR)ENDIF124第六章文件输入输出(MPI–IO)在该例中,因为atomicity被设为true,因此进程1将总是读到0个数或10个5.
如果改变上面的程序将atomicity设为false,则进程1读到的结果是不确定的,它与具体的MPI实现和程序运行过程有关.
§6.
7.
2查询atomicity的当前值CintMPI_File_get_atomicity(MPI_Filefh,int*flag)FortranMPI_FILE_GET_ATOMICITY(FH,FLAG,IERR)INTEGERFH,IERRLOGICALFLAG该函数在参数flag中返回atomicity的当前值.
§6.
7.
3文件读写与存储设备间的同步CintMPI_File_sync(MPI_Filefh)FortranMPI_FILE_SYNC(FH,IERR)INTEGERFH,IERR该函数确保将调用它的进程新近写入文件的数据写入存储设备.
如果文件在存储设备中的内容已被其它进程改变,则调用该函数可确保调用它的进程随后读该文件时得到的是改变后的数据.

调用该函数时不能有尚未完成的对该文件的非阻塞型或分裂型读写操作.
注意,如果打开文件的进程组的两个进程中一个进程往文件中写入一组数据,另一个进程希望从文件的同一位置读到这组数据,则各进程可能需要调用两次MPI_FILE_SYNC,并在两次调用间进行一次同步(MPI_BARRIER).
第一次对MPI_FILE_SYNC的调用确保第一个进程写的数据被写入存储设备,而第二次调用则可确保新写入存储设备的数据被另一个进程读到.
MPI_FILE_SYNC是聚合型函数,进程组中所有进程必须同时调用并且提供相同的参数.
§6.
8子数组数据类型创建函数CintMPI_Type_create_subarray(intndims,intarray_of_sizes[],intarray_of_subsizes[],intarray_of_starts[],intorder,MPI_Datatypeoldtype,MPI_Datatype*newtype)FortranMPI_TYPE_CREATE_SUBARRAY(NDIMS,ARRAY_OF_SIZES,+ARRAY_OF_SUBSIZES,ARRAY_OF_STARTS,ORDER,+OLDTYPE,NEWTYPE,IERR)INTEGERNDIMS,ARRAY_OF_SIZES(*),ARRAY_OF_SUBSIZES(*),+ARRAY_OF_STARTS(*),ORDER,OLDTYPE,NEWTYPE,IERR§6.
8子数组数据类型创建函数125这是一个辅助数据类型构造函数,可用于对分布式数组的读写操作.
该函数创建一个"子数组"数据类型,即描述一个n维(全局)数组中的一个n维子数组.
创建的新数据类型的域对应于全局数组,类型中数据的位移由子数组的元素在全局数组中的位移确定.
ndims给出数组的维数.
array_of_sizes[i]给出全局数组第i维的大小.
array_of_subsizes[i]给出子数组第i维的大小.
array_of_starts[i]给出子数组第i维在全局数组中的起始位置(不论C还是Fortran均用0表示从全局数组的第一个元素开始).
参数order给出数组元素的排列顺序,order=MPI_ORDER_C表示数组元素按C的数组顺序排列,order=MPI_ORDER_FORTRAN表示数组元素按Fortran的数组顺序排列.
oldtype给出数组元素的数据类型.
newtype返回所创建的子数组数据类型的句柄.
子数组各维的大小必须大于0并且小于或等于全局数组相应维的大小.
子数组的起始位置可以是全局数组中的任何位置,但必须确保子数组被包含在全局数组中,否则函数调用报错.
如果数据类型oldtype是可移植数据类型,则新数据类型newtype也是可移植数据类型.
例7.
2.
5中给出了一个使用MPI–IO函数的程序示例.
126第六章文件输入输出(MPI–IO)第七章MPI程序示例§7.
1矩阵乘积A为M*N阶矩阵,B为N*L阶矩阵,C为M*L阶矩阵.
计算矩阵乘积C=AB.
§7.
1.
1算法描述假设:使用NPROCS个MPI进程,为简单起见假定M和L均为NPROCS的倍数.
A和C按行等分成子块分别存储在不同的进程中,而B则按列等分成子块分别存储在不同的进程中.
A,B和C的子块大小分别为MLOC*N,N*LLOC和MLOC*L其中MLOC=M/NPROCS,LLOC=L/NPROCS.
具体存储方式为:{ai,j|k*MLOCb(可以通过在计算/通信中交替使用b/work来避免该操作)41doj=1,lloc42doi=1,n43b(i,j)=work(i,j)44enddo45enddo46*47jpos=jpos+lloc48if(jpos.
ge.
l)jpos=049*50enddo51*52doj=1,lloc53doi=1,mloc54sum=0.
d055dok=1,n56sum=sum+a(i,k)*b(k,j)57enddo58c(i,j+jpos)=sum59enddo60enddo61*62return63end例7.
1.
4.
调用BLAS库函数完成矩阵子块的乘积.
选用适当的BLAS库可以大幅度提高程序的实际运行性能.
源程序在文件matrix-blas.
f中.
注意,编译该程序时必须与BLAS库链接.
1*2*ParallelmultiplicationofmatricesusingMPI_Isend/MPI_IrecvandBLAS3*4subroutinematmul(m,n,l,myrank,nprocs,a,b,c,work)5implicitdoubleprecision(a-h,o-z)6include'mpif.
h'7dimensiona(m/nprocs,n),b(n,l/nprocs),c(m/nprocs,l),8&work(n,l/nprocs)9integersrc,dest,tag10integerstatus(MPI_STATUS_SIZE,2),request(2)11*12mloc=m/nprocs13lloc=l/nprocs14*§7.
2Poisson方程求解13315dest=mod(myrank-1+nprocs,nprocs)16src=mod(myrank+1,nprocs)17*18jpos=myrank*lloc19*20doip=1,nprocs-121tag=10000+ip22*23callMPI_Isend(b,n*lloc,MPI_DOUBLE_PRECISION,dest,tag,24&MPI_COMM_WORLD,request(1),ierr)25callMPI_Irecv(work,n*lloc,MPI_DOUBLE_PRECISION,src,tag,26&MPI_COMM_WORLD,request(2),ierr)27*28callDGEMM('n','n',mloc,lloc,n,1.
d0,a,mloc,b,n,0.
d0,29&c(1,1+jpos),mloc)30*31callMPI_Waitall(2,request,status,ierr)32*33*拷贝work->b(可以通过在计算/通信中交替使用b/work来避免该操作)34callDCOPY(n*lloc,work,1,b,1)35*36jpos=jpos+lloc37if(jpos.
ge.
l)jpos=038*39enddo40*41callDGEMM('n','n',mloc,lloc,n,1.
d0,a,mloc,b,n,0.
d0,42&c(1,1+jpos),mloc)43*44return45end作业7.
1.
1.
修改例7.
1.
4中的程序,以避免第34行的矩阵拷贝(DCOPY)操作.
作业7.
1.
2.
修改例7.
1.
4中的程序,要求在不使用中间数组work的情况下实现计算与通信的(部分)重叠.
修改后的程序对性能有什么影响(提示:将数组B分成两半,一半通信,一半计算).
§7.
2Poisson方程求解设计MPI程序,求解定义在二维规则区域上的Poisson方程:u(x,y)=f(x,y)(x,y)0,a)*(0,b)u(x,y)|=g(x,y)(7.
1)其中,f(x,y)和g(x,y)为已知函数,分别定义在区域的内部和边界.
§7.
2.
1并行算法设计沿坐标轴x和y方向,分别取步长hx=aIM,hy=bJM(7.
2)134第七章MPI程序示例将区域离散成规模为IM*JM的网格,其中IM和JM分别为沿坐标轴x和y方向的网格单元个数.
不妨设方程(7.
1)的近似解u(x,y)定义在所有网格结点上,且用如下未知量表示ui,j=u(i*hx,j*hy)1≤i≤IM1,1≤j≤JM1ui,j=gi,j=g(i*hx,j*hy)i=0或i=IM或j=0或j=JM(7.
3)用二阶中心差商近似导数:uxx(i*hx,j*hy)≈ui1,j2ui,j+ui+1,jh2x,uyy(i*hx,j*hy)≈ui,j12ui,j+ui,j+1h2y(7.
4)并记fi,j=f(i*hx,j*hy)(7.
5)将以上程式代入式(7.
1),则问题转化为求解稀疏线性代数方程组:2(h2x+h2y)ui,jh2y(ui1,j+ui+1,j)h2x(ui,j1+ui,j+1)=h2xh2yfi,j,1≤i≤IM1,1≤j≤JM1(7.
6)具体地,我们将采用众所周知的Jacobi点迭代算法求解方程(7.
6).
§7.
2.
2MPI并行程序设计设计求解方程(7.
1)的MPI并行程序必须考虑两个关键问题:第一,选择正确的区域分解策略.
将区域分解成多个子区域,分配给不同的进程,并保证进程间的负载平衡和最小的消息传递通信开销.
通常有两种方式:1)沿y方向的一维条分解策略,如图7.
2(a)所示;2)沿两个方向的二维块分解策略,如图7.
2(b)所示.
显然,如果取x方向的进程个数等于1,则二维块分解策略就退化为一维条分解策略.
无论哪种方式,都应该尽量保证每个子区域包含的网格结点个数相等,因为这样才能保证进程间的负载平衡.
(a)一维条分解(b)二维块分解图7.
2:两种区域分解策略第二,选择合适的通信数据结构.
由式(7.
6)可知,在任意网格结点上,执行Jacobi点迭代需要知道该结点上、下、左、右四个相邻结点上的近似解.
因此,在每次Jacobi迭代之前,每个进程必须与其相邻的进程交换边界结点上的近似解.
为了描述通信数据结构,不妨设方程近似解定义在网格单元的中心.
图7.
3给出了一个3*3的二维块区域分解,其中,各个子区域被分配给不同的进程,各个进程负责求解该子区域的近似解.
具体地,进程5将向其相邻的四个进程(进程2、进程4、进程6和进程8)输出""标示的网格单元的近§7.
2Poisson方程求解135似解;同时,从这四个进程接收用""标示的网格单元上的近似解.
因此,如何管理这些沿区域分解边界交换的网格单元上的量,将会直接影响到MPI程序的并行性能.
图7.
4给出了一个比较有效的办法,就是沿各个子区域的边界,向外增加一个宽度为1的辅助网格单元,用于存储相邻子区域在这些网格单元上的近似解.
图7.
3:3*3二维区域分解图7.
4:辅助网格单元示意图类似地,同样的数据结构也适应于近似解定义在网格结点上的情形,这里不再讨论.

下面,基于以上二维块区域分解策略和通信数据结构,给出求解差分方程(7.
6)的MPI并行程序.
其中,近似解定义在网格结点上.
为了简单,这里假设网格单元个数IM和JM能分别被沿x方向和y方向的进程个数NPX和NPY整除,进程的序号按自然序排列(先沿x方向,后沿y方向).
例7.
2.
1.
并行Jacobi点迭代MPI程序:二维块分解策略(源程序见文件poisson0.
f).
1!
Poisson方程求解:使用阻塞通信(可能死锁).
作者:莫则尧2INCLUDE"mpif.
h"3PARAMETER(DA=2.
0,DB=3.
0)!
问题求解区域沿X、Y方向的大小4PARAMETER(IM=30,JM=60)!
沿X、Y方向的全局网格规模5PARAMETER(NPX=1,NPY=1)!
沿X、Y方向的进程个数6PARAMETER(IML=IM/NPX,JML=JM/NPY)7!
各进程沿X、Y方向的局部网格规模,仅为全局网格规模的1/(NPX*NPY)8REALU(0:IML+1,0:JML+1)!
定义在网格结点的近似解9REALUS(0:IML+1,0:JML+1)!
定义在网格结点的精确解10REALU0(IML,JML)!
Jacobi迭代辅助变量11REALF(IML,JML)!
函数f(x,y)在网格结点上的值12INTEGERNPROC!
mpirun启动的进程个数,必须等于NPX*NPY13INTEGERMYRANK,MYLEFT,MYRIGHT,MYUPPER,MYLOWER14!
各进程自身的序号,4个相邻进程的序号15INTEGERMEPX,MEPY!
各进程自身的序号沿X,Y方向的坐标16REALXST,YST!
各进程拥有的子区域沿X,Y方向的起始坐标17REALHX,HY!
沿X,Y方向的网格离散步长18REALHX2,HY2,HXY2,RHXY19INTEGERIST,IEND,JST,JEND20!
各进程沿X,Y方向的内部网格结点的起始和终止坐标21INTEGERHTYPE,VTYPE22!
MPI用户自定义数据类型,表示各进程沿X,Y方向23!
与相邻进程交换的数据单元24INTEGERSTATUS(MPI_STATUS_SIZE)25!
In-linefunctions26solution(x,y)=x**2+y**2!
解析解136第七章MPI程序示例27rhs(x,y)=-4.
0!
Poisson方程源项(右端项)28!
29!
程序可执行语句开始30CALLMPI_INIT(IERR)31CALLMPI_COMM_SIZE(MPI_COMM_WORLD,NPROC,IERR)32IF(NPROC.
NE.
NPX*NPY.
OR.
MOD(IM,NPX).
NE.
0.
OR.
MOD(JM,NPY).
NE.
0)THEN33PRINT*,'+++mpirun-npxxxerrorORgridscaleerror,',34&'exitout+++'35CALLMPI_FINALIZE(IERR)36STOP37ENDIF38!
39!
按自然序(先沿X方向,后沿Y方向)确定各进程自身及其4个相邻进程的序号40CALLMPI_COMM_RANK(MPI_COMM_WORLD,MYRANK,IERR)41MYLEFT=MYRANK-142IF(MOD(MYRANK,NPX).
EQ.
0)MYLEFT=MPI_PROC_NULL43MYRIGHT=MYRANK+144IF(MOD(MYRIGHT,NPX).
EQ.
0)MYRIGHT=MPI_PROC_NULL45MYUPPER=MYRANK+NPX46IF(MYUPPER.
GE.
NPROC)MYUPPER=MPI_PROC_NULL47MYLOWER=MYRANK-NPX48IF(MYLOWER.
LT.
0)MYLOWER=MPI_PROC_NULL49MEPY=MYRANK/NPX50MEPX=MYRANK-MEPY*NPX51!
对应二维NPY(NPXCartesian行主序坐标为(MEPY,MEPX).
52!
53!
基本变量赋值,确定各进程负责的子区域54HX=DA/IM55HX2=HX*HX56HY=DB/JM57HY2=HY*HY58HXY2=HX2*HY259RHXY=0.
5/(HX2+HY2)60XST=MEPX*DA/NPX61YST=MEPY*DB/NPY62IST=163IEND=IML64IF(MEPX.
EQ.
NPX-1)IEND=IEND-1!
最右边的区域X方向少一个点65JST=166JEND=JML67IF(MEPY.
EQ.
NPY-1)JEND=JEND-1!
最上边的区域Y方向少一个点68!
69!
数据类型定义70CALLMPI_TYPE_CONTIGUOUS(IEND-IST+1,MPI_REAL,HTYPE,IERR)71CALLMPI_TYPE_COMMIT(HTYPE,IERR)72!
沿X方向的连续IEND-IST+1个MPI_REAL数据单元,73!
可用于表示该进程与其上、下进程交换的数据单元74CALLMPI_TYPE_VECTOR(JEND-JST+1,1,IML+2,MPI_REAL,VTYPE,IERR)75CALLMPI_TYPE_COMMIT(VTYPE,IERR)76!
沿Y方向的连续JEND-JST+1个MPI_REAL数据单元,77!
可用于表示该进程与其左、右进程交换的数据单元78!
79DOJ=JST-1,JEND+1§7.
2Poisson方程求解13780DOI=IST-1,IEND+181xx=(I+MEPX*IML)*HX!
xx=XST+I*HX82yy=(J+MEPY*JML)*HY!
yy=YST+J*HY83IF(I.
GE.
IST.
AND.
I.
LE.
IEND.
AND.
J.
GE.
JST.
AND.
J.
LE.
JEND)THEN84U(I,J)=0.
0!
近似解赋初值85US(I,J)=solution(xx,yy)!
解析解86F(I,J)=rhs(xx,yy)!
右端项87ELSEIF((I.
EQ.
IST-1.
AND.
MEPX.
EQ.
0).
OR.
88&(J.
EQ.
JST-1.
AND.
MEPY.
EQ.
0).
OR.
89&(I.
EQ.
IEND+1.
AND.
MEPX.
EQ.
NPX-1).
OR.
90&(J.
EQ.
JEND+1.
AND.
MEPY.
EQ.
NPY-1))THEN91U(I,J)=solution(xx,yy)!
边界值92ENDIF93ENDDO94ENDDO95!
96!
Jacobi迭代求解97NITER=098100CONTINUE99NITER=NITER+1100!
101!
交换定义在辅助网格结点上的近似解102CALLMPI_SEND(U(1,1),1,VTYPE,MYLEFT,NITER+100,103&MPI_COMM_WORLD,IERR)!
发送左边界104CALLMPI_SEND(U(IEND,1),1,VTYPE,MYRIGHT,NITER+100,105&MPI_COMM_WORLD,IERR)!
发送右边界106CALLMPI_SEND(U(1,1),1,HTYPE,MYLOWER,NITER+100,107&MPI_COMM_WORLD,IERR)!
发送下边界108CALLMPI_SEND(U(1,JEND),1,HTYPE,MYUPPER,NITER+100,109&MPI_COMM_WORLD,IERR)!
发送上边界110CALLMPI_RECV(U(IEND+1,1),1,VTYPE,MYRIGHT,NITER+100,111&MPI_COMM_WORLD,STATUS,IERR)!
接收右边界112CALLMPI_RECV(U(0,1),1,VTYPE,MYLEFT,NITER+100,113&MPI_COMM_WORLD,STATUS,IERR)!
接收左边界114CALLMPI_RECV(U(1,JEND+1),1,HTYPE,MYUPPER,NITER+100,115&MPI_COMM_WORLD,STATUS,IERR)!
接收上边界116CALLMPI_RECV(U(1,0),1,HTYPE,MYLOWER,NITER+100,117&MPI_COMM_WORLD,STATUS,IERR)!
接收下边界118DOJ=JST,JEND119DOI=IST,IEND120U0(I,J)=RHXY*(HXY2*F(I,J)+HX2*(U(I,J-1)+U(I,J+1))121&+HY2*(U(I-1,J)+U(I+1,J)))122ENDDO123ENDDO124!
125!
计算与精确解间的误差126ERR=0.
0127DOJ=JST,JEND128DOI=IST,IEND129U(I,J)=U0(I,J)130ERR=MAX(ERR,ABS(U(I,J)-US(I,J)))!
用L\infty模以使误差与NP无关131ENDDO132ENDDO138第七章MPI程序示例133ERR0=ERR134CALLMPI_ALLREDUCE(ERR0,ERR,1,MPI_REAL,MPI_MAX,135&MPI_COMM_WORLD,IERR)136IF(MYRANK.
EQ.
0.
AND.
MOD(NITER,100).
EQ.
0)THEN137PRINT*,'NITER=',NITER,',ERR=',ERR138ENDIF139!
140IF(ERR.
GT.
1.
E-3)THEN!
收敛性判断141GOTO100!
没有收敛,进入下次迭代142ENDIF143!
144IF(MYRANK.
EQ.
0)THEN145PRINTSuccessfullyconvergedafter',146&NITER,'iterations'147PRINTerror=',ERR148ENDIF149!
输出近似解(略)150!
151CALLMPI_FINALIZE(IERR)152END在上例中,我们固定了该MPI程序产生的进程个数NP=NPY*NPX,这样做的一个主要目的是为了使各进程的内存空间仅为相应串行程序的1/(NPY*NPX),从而使得原来串行程序在单机上由于内存不够而无法计算的问题通过多机并行计算而成为可能.
当然,这样做也带来一些不便,例如它要求MPI运行命令"mpirun-npxxx"中的参数xxx等于NPY*NPX,且如果参数NPY和NPX被改变后,必须重新编译该程序.
§7.
2.
3MPI并行程序的改进在例7.
2.
1中,交换定义在辅助网格结点近似解的8条消息传递语句(102–117行)是该MPI程序的关键.
但是,由MPI标准可知,它们是不安全的,因为在某些并行机上,当消息较长时(例如大于16KB),可能由于MPI系统缓存区大小的限制,而导致执行该MPI程序的进程的死锁.
因此,我们将它们替换成如下的非阻塞通信函数.
例7.
2.
2.
改进一:非阻塞通信(源程序见文件poisson1.
f)略)!
用下一行替换poisson0.
f第24行INTEGERREQ(8),STATUS(MPI_STATUS_SIZE,8)略)!
用下述内容替换poisson0.
f第102-117行CALLMPI_ISEND(U(1,1),1,VTYPE,MYLEFT,NITER+100,&MPI_COMM_WORLD,REQ(1),IERR)!
发送左边界CALLMPI_ISEND(U(IEND,1),1,VTYPE,MYRIGHT,NITER+100,&MPI_COMM_WORLD,REQ(2),IERR)!
发送右边界CALLMPI_ISEND(U(1,1),1,HTYPE,MYLOWER,NITER+100,&MPI_COMM_WORLD,REQ(3),IERR)!
发送下边界CALLMPI_ISEND(U(1,JEND),1,HTYPE,MYUPPER,NITER+100,&MPI_COMM_WORLD,REQ(4),IERR)!
发送上边界CALLMPI_IRECV(U(IEND+1,1),1,VTYPE,MYRIGHT,NITER+100,&MPI_COMM_WORLD,REQ(5),IERR)!
接收右边界§7.
2Poisson方程求解139CALLMPI_IRECV(U(0,1),1,VTYPE,MYLEFT,NITER+100,&MPI_COMM_WORLD,REQ(6),IERR)!
接收左边界CALLMPI_IRECV(U(1,JEND+1),1,HTYPE,MYUPPER,NITER+100,&MPI_COMM_WORLD,REQ(7),IERR)!
接收上边界CALLMPI_IRECV(U(1,0),1,HTYPE,MYLOWER,NITER+100,&MPI_COMM_WORLD,REQ(8),IERR)!
接收下边界CALLMPI_WAITALL(8,REQ,STATUS,IERR)!
阻塞式等待消息传递的结束略)在例7.
2.
1中,我们可以将118–123行的循环分裂成两个部分,其中一个部分需要辅助网格点上的近似解,而另一个部分不需要辅助网格点上的近似解.
这样,为了改进该MPI程序的并行性能,我们可以将后一个部分的计算与例7.
2.
2的非阻塞消息传递重叠起来,从而达到屏蔽网络延迟的目的.
具体改进如下.
例7.
2.
3.
改进二:重叠通信与计算(源程序见文件poisson2.
f)略)!
用下述内容替换poisson1.
f第119-124行DOJ=JST+1,JEND-1DOI=IST+1,IEND-1U0(I,J)=RHXY*(HXY2*F(I,J)+HX2*(U(I,J-1)+U(I,J+1))&+HY2*(U(I-1,J)+U(I+1,J)))ENDDOENDDOCALLMPI_WAITALL(8,REQ,STATUS,IERR)!
阻塞式等待消息传递的结束DOJ=JST,JEND,JEND-JSTDOI=IST,IENDU0(I,J)=RHXY*(HXY2*F(I,J)+HX2*(U(I,J-1)+U(I,J+1))&+HY2*(U(I-1,J)+U(I+1,J)))ENDDOENDDODOJ=JST,JENDDOI=IST,IEND,IEND-ISTU0(I,J)=RHXY*(HXY2*F(I,J)+HX2*(U(I,J-1)+U(I,J+1))&+HY2*(U(I-1,J)+U(I+1,J)))ENDDOENDDO略)在例7.
2.
1–7.
2.
3中,各进程按自然序(先沿x方向,后沿y方向)确定与它相邻的4个进程的序号(MYLEFT,MYRIGHT,MYLOWER,MYUPPER),以及它自己所处的行主序位置(MEPY,MEPX).
实际上,这些进程按区域分解策略可以很自然地映射到NPY*NPX的二维Cartesian拓扑结构(参看§5.
18.
1),而(MEPY,MEPX)就是各进程在该拓扑结构中的坐标.
因此,我们可以从通信器MPI_COMM_WORLD出发,建立二维Cartesian拓扑结构,从而方便地确定各进程的相邻关系,并使得之后的所有MPI消息传递均基于该拓扑结构进行.
例7.
2.
4.
改进三:二维Cartesian拓扑结构(源程序见文件poisson3.
f)略)!
在poisson[012].
f程序头(变量声明部分)加入下面二行INTEGERCOMM,DIMS(2),COORD(2)LOGICALPERIOD(2),REORDER140第七章MPI程序示例略)!
用下述内容替换poisson[012].
f第39-50行DIMS(1)=NPY!
拓扑结构中Y方向的进程个数DIMS(2)=NPX!
拓扑结构中X方向的进程个数PERIOD(1)=.
FALSE.
!
沿Y方向,拓扑结构非周期连接PERIOD(2)=.
FALSE.
!
沿X方向,拓扑结构非周期连接REORDER=.
TRUE.
!
在新通信器中,进程重新排列序号CALLMPI_CART_CREATE(MPI_COMM_WORLD,2,DIMS,PERIOD,REORDER,&COMM,IERR)CALLMPI_COMM_RANK(COMM,MYRANK,IERR)CALLMPI_CART_COORDS(COMM,MYRANK,2,COORD,IERR)MEPY=COORD(1)MEPX=COORD(2)CALLMPI_CART_SHIFT(COMM,0,1,MYLOWER,MYUPPER,IERR)!
Y方向CALLMPI_CART_SHIFT(COMM,1,1,MYLEFT,MYRIGHT,IERR)!
X方向略)在例7.
2.
1–7.
2.
4中,我们忽略了近似解的输出,这里,我们用第6章中介绍的MPI并行I/O函数实现近似解的并行输出.
特别地,我们要求输出的近似解按自然序排列,且包含物理边界结点.
例7.
2.
5.
改进四:并行I/O(源程序见文件poisson4.
f).
该程序使用了独立文件指针、聚合型函数MPI_FILE_WRITE_ALL(参看表6.
5).
略)!
在程序头(变量声明部分)加入下面内容INTEGERFH,FILETYPE,MEMTYPE,GSIZE(2),LSIZE(2),START(2)!
注意:从下面三种形式的变量声明中根据所使用的MPI系统选择一个正确的!
(可以参考文件mpiof.
h或mpif.
h中MPI_OFFSET_KIND的定义)!
INTEGER(kind=MPI_OFFSET_KIND)OFFSET!
适用于Fortran90!
INTEGER*8OFFSET!
适用于64位系统INTEGER*4OFFSET!
适用于某些32位系统略)!
在标有"输出近似解(略)"处(倒数第四行)加入下述内容GSIZE(1)=IM+1GSIZE(2)=JM+1LSIZE(1)=IEND-IST+1IF(MEPX.
EQ.
0)LSIZE(1)=LSIZE(1)+1IF(MEPX.
EQ.
NPX-1)LSIZE(1)=LSIZE(1)+1LSIZE(2)=JEND-JST+1IF(MEPY.
EQ.
0)LSIZE(2)=LSIZE(2)+1IF(MEPY.
EQ.
NPY-1)LSIZE(2)=LSIZE(2)+1START(1)=IML*MEPXIF(MEPX.
NE.
0)START(1)=START(1)+1START(2)=JML*MEPYIF(MEPY.
NE.
0)START(2)=START(2)+1!
!
定义局部子数组数据类型CALLMPI_TYPE_CREATE_SUBARRAY(2,GSIZE,LSIZE,START,&MPI_ORDER_FORTRAN,MPI_REAL,FILETYPE,IERR)CALLMPI_TYPE_COMMIT(FILETYPE,IERR)!
!
打开二进制文件CALLMPI_FILE_OPEN(COMM,'result.
dat',§7.
2Poisson方程求解141&MPI_MODE_CREATE+MPI_MODE_WRONLY,&MPI_INFO_NULL,FH,IERR)OFFSET=0!
注意使用正确的变量类型(INTEGER*4或INTEGER*8)!
(参考文件mpif.
h中MPI_OFFSET_KIND的定义)CALLMPI_FILE_SET_VIEW(FH,OFFSET,MPI_REAL,FILETYPE,&'native',MPI_INFO_NULL,IERR)!
!
定义数据类型,描述子数组在内存中的分布GSIZE(1)=IML+2GSIZE(2)=JML+2START(1)=1IF(MEPX.
EQ.
0)START(1)=0START(2)=1IF(MEPY.
EQ.
0)START(2)=0CALLMPI_TYPE_CREATE_SUBARRAY(2,GSIZE,LSIZE,START,&MPI_ORDER_FORTRAN,MPI_REAL,MEMTYPE,IERR)CALLMPI_TYPE_COMMIT(MEMTYPE,IERR)!
!
输出近似解(含物理边界结点)CALLMPI_FILE_WRITE_ALL(FH,U,1,MEMTYPE,STATUS,IERR)CALLMPI_FILE_CLOSE(FH,IERR)!
略)至止,例7.
2.
2例7.
2.
5分别从非阻塞通信、重叠通信与计算、拓扑结构和并行I/O四个方面,利用相应的MPI函数,依次改进了例7.
2.
1中MPI程序的功能和并行性能.
通过该应用示例,读者可以较好地将所介绍的MPI函数联系在一起,解决实际问题.
作业7.
2.
1.
在例7.
2.
1–7.
2.
5中用了一排辅助网格单元.
通过增加虚拟网格点的宽度可以提高通信粒度.
例如,如果使用两排辅助网格单元,则可以每两次迭代交换一次子区域边界附近的近似解,代价是子区域间增加了少量的重复计算.
试修改例7.
2.
4或例7.
2.
5中的程序,在程序中增加一个参数BW,它代表辅助网格单元的宽度(bw≥1),并比较不同BW的值对程序性能的影响.
作业7.
2.
2.
修改例7.
2.
5中的程序,将Jacobi迭代改为红黑顺序的Gauss–Seidel迭代.
142第七章MPI程序示例参考文献[1]莫则尧,袁国兴,《消息传递并行编程环境MPI》,科学出版社,2001.
[2]J.
M.
Ortega,IntroductiontoParallelandVectorSolutionofLinearSystems,PlenumPress,1988.
[3]J.
Dongarra,I.
Du,D.
SoresenandH.
vanderVorst,SolvingLinearSystemsonVectorandSharedMemoryComputers,SIAM,1991.
[4]G.
GolubandC.
vanLoan,MatrixComputation,TheJohnsHopkinsUniversityPress,1983.
(中译本:矩阵计算,廉庆荣、邓健新、刘秀兰译,大连理工大学出版社,1988年)[5]迟学斌,在具有局部内存与共享主存的并行机上并行求解线性方程组,计算数学,1995年第17卷第2期.
[6]G.
Y.
LiandT.
F.
Coleman,AParallelTriangularSolverforaHypercubeMultiprocessor,TR86–787,CornellUniversity,1986.
[7]迟学斌,Transputer上Cholesky分解的并行实现,计算数学,1993年第15卷第3期.
[8]J.
M.
DelosmeandI.
C.
F.
Ipsen,PositiveDeniteSystemswithHyperbolicRotations,LinearAlgebraAppl.
,77(1986),75–111.
[9]D.
H.
LawrieandA.
H.
Sameh,TheComputationandCommunicationComplexityofaParallelBandedSystemSolver,ACMTrans.
Math.
Soft.
,10(1984),185–195.
[10]迟学斌,Transputer上线性系统的并行求解,中国计算机用户,1991年第10期.
[11]陈景良,并行数值方法,清华大学出版社,1983年.
[12]张宝琳、袁国兴、刘兴平、陈劲,偏微分方程并行有限差分方法,科学出版社,1994年.
[13]D.
ChazanandW.
Miranker,ChaoticRelazation,J.
Lin.
Alg.
Appl.
,2(1969),199–222.
[14]G.
Baudet,AsynchronousIterativeMethodsforMultiprocessors,J.
ACM,25(1978),226–244.
[15]J.
M.
OrtegaandW.
C.
Rheinboldt,IterativeSolutionofNonlinearEquationsinSeveralVariables,AcademicPress,1970.
[16]迟学斌,线性方程组的异步迭代法,计算数学,1992年第14卷第3期.
[17]S.
Balay,W.
Gropp,L.
CurfmanMcInnes,andB.
Smith,PETSchomepage,http://www.
mcs.
anl.
gov/petsc.
[18]P.
Bridge,N.
Doss,W.
Gropp,etc.
,User'sGuidetoMPICH,aPortableImplementationofMPI,ANL/MCS-TM-000,1994[19]P.
Bridge,N.
Doss,W.
Gropp,etc.
,InstallationGuidetoMPICH,aPortableImplementationofMPI,ANL/MCS-TM-000,1994[20]D.
E.
Culler,J.
P.
Singh,A.
Gupta,ParallelComputerArchitecture:aHardware/SoftwareApproach,MorganKaufmannpublisher,1999.
[21]陈国良,并行计算:结构、算法与编程,北京:高等教育出版社,1999.
[22]J.
Dongarra,B.
Tourancheau,editors,ProdceedingsoftheWorkshoponEnvironmentsandToolsforParallelScienticComputing,SIAMpublications,1994.
[23]W.
Gropp,E.
Lusk,andA.
Skjellum,UsingMPI:PortableParallelProgrammingwiththeMessage–PassingInterface,MITPress,Cambridge,MA,1994.
[24]W.
Gropp,E.
Lusk,andA.
Skjellum,UsingMPI:PortableParallelProgrammingwiththeMessage–PassingInterface,2ndedition,MITPress,Cambridge,MA,1999.
[25]W.
Gropp,E.
Lusk,N.
Doss,andA.
Skjellum,AHigh–PerformancePortableImplementationoftheMPIMessage–PassingInterfaceStandard,ParallelComputing,Vol.
22,No.
6,1996,pp.
789–828.
[26]W.
Gropp,E.
Lusk,R.
Thakur,UsingMPI–2:AdvancedFeaturesoftheMessage–PassingInterface,MITPress,1999.
[27]A.
Geist,A.
Beguelin,J.
Dongarra,W.
Jiang,B.
Manchek,andV.
Sunderam,PVM:ParallelVirtualMachine143144参考文献—AUser'sGuideandTutorialforNetworkParallelComputing,MITPress,Cambridge,MA,1994.
[28]R.
Hempel,D.
W.
Walker,TheEmergenceofTheMPIMessagePassingStandardforParallelComputing,ComputerStandardandInterface,Vol.
21,1999,pp.
51–62.
[29]IMPISteeringCommittee.
IMPI—InteroperableMessage–PassingInterface,1998,http://impi.
nist.
gov/IMPI.
[30]李晓梅,莫则尧等,可扩展并行算法设计与分析,北京:国防工业出版社,2000.