求解一类边界问题的最小曲面的数值并行算法
◎ 刘 婕
摘 要:本文给出了具有最小面积约束的一类边界问题的数值求法,同时实现了该算法的并行化。在算例中,介绍了利用Jacobi迭代求解曲顶柱体顶面面积最小值的一种并行算法,并阐述了解决这一问题的实际意义。算例结果表明,该并行算法的并行效率令人满意。
关键词:边界问题 Jacobi迭代法 并行算法 曲顶柱体 MPI
1. 引 言
在科学研究和生产实践中,经常会遇到这样一类问题,即在三维空间R3中给定封闭曲线 ,如何确定光滑曲面S,使该曲面以 为边界并具有最小曲面。
为方便叙述,设G是OXY平面上的一个矩形开区域,f是其边界 上的已知连续函数,则具有最小面积约束的边界问题可以叙述为:
求有界闭区域 上的连续函数U,使满足
(i) 函数U在边界 上满足: (1)
(ii) 函数U确定的闭曲面S具有最小面积,其中
本文将设计一种算法,求问题(1)的数值解.
随着并行计算机的问世,并行算法和并行系统不断的发展和完善,出现了大量的软件开发平台环境,其中著名的如,可移植的异构编程环境PVM(Parallel
Virtual machine)、消息传递标准平台MPI(Message Passing Interface)、支持数据并行程序设计的HPF(High
PerformanceFortran)等。其中MPI统一了互不兼容的用户界面,制定了消息传递界面的新标准,支持最先进的可移植平台,具备实用性、可移植性、效率高和可扩展性等特点,因此获得了来自研制比并行计算机的主要厂商、高校、实验室与工业界的广泛支持和热烈欢迎。
本文提出了一种应用MPI程序实现求解一类边界问题的最小曲面的并行Jacobi算法,分析了该并行Jacobi迭代算法的并行过程,给出了该算法数值测试结果,有效并合理的解决了曲顶柱体顶面面积的最小值求解这一现实问题。
2. 算法的设计
2.1 Jacobi迭代法的基本思想
Jacobi迭代是一种比较常见的迭代方法,其串行算法的实质是由原先旧值点的相邻四点值或八点值的某种运算而得到新值,通过多次迭代,内部点的值逐步稳定下来。
相邻四点Jacobi迭代可表示为
DO K=1,STEP
DO J=1,N
DO I=1,N
B(I,J)=0.25*(A(I-1,J)+A(I+1,J)+A(I,J+1)+A(I,J-1))
END DO
END DO
DO J=1,N
DO I=1,N
END DO
A(I,J)=B(I,J)
END DO
2.2推广的 Jacobi迭代法
2.1中所阐述的迭代法是运用的相邻四点的算术平均的运算,我们可以将Jacobi迭代进行推广,最重要的迭代过程可表示为:
B(I,J)=f(A(I-1,J-1),A(I,J-1),A(I+1,J-1),A(I-1,J),A(I+1,J),A(I-1,J+1),
A(I,J+1),A(I+1,J+1))
其中,迭代函数f应根据实际应用背景来确定.
2.3问题的描述及其Jacobi迭代法的具体设计
如图1所示,将矩形区域G分割成细网格:
取定沿x 方向和y方向的步长分别为a和b,从坐标原点(0,0)开始
作两族与坐标轴x和y平行的直线:
x=ia,i=0,1,2...m;
y=jb,j=0,1,2...n;
分割时,使区域 ={ }.两族直线的交点(ia,jb)称为网格结点,记为( , )。属于G的结点称为内点,用 ={( , ) }表示内点全体组成的集合。网线与边界
的交点称为界点,用 表示界点全体组成的集合。令 。
从数学理论上来说,问题(1)的解析解是存在的,设解为 。
将问题(1) 离散化,也即:
求 上的结点值 ,i=0,1,2...m, j=0,1,2...n,使满足
(i)在 上: ( , ), ( , ) (2)
(ii)在 上: ( , ), ( , )
下面导出 , )的近似表达式:
如图2所示,空间曲面ABCDEFGHO的面积可用8个小三角形ΔABO、ΔBCO、…、ΔHAO的面积S1、S2、…、S8的和 来近似表示。其中
、 、…、 。
事实上,当A的结点值 、B的结点值 、H的结点值 已经确定时,为使 达到最小,O的结点值 将被唯一确定,此时 可近似表示为
下面具体设计Jacobi迭代算法求解问题(2)
Step 1: , ,j=0,1,2…,n
, ,i=0,1,2…,m
i=1,2…,m-1; j=1,2…,n-1
Step 2:
i=1,2…,m-1; j=1,2…,n-1
Step 3: 若 , (其中 是允许误差限),则终止运算,
否则 i=1,2…,m-1; j=1,2…,n-1。并返回Step 2。
由于Jacobi迭代的局部性很好,可以取得很高的并行性,下节将进一步阐述Jacobi迭代算法的并行化过程。
3.算法的并行实现
3.1 Jacobi迭代的并行化
在给定边界问题的最小曲面面积的问题中,Jacobi迭代的边界值太多,如果采取串行算法肯定会花费很多的时间,这将给工作带来不便。由于Jacobi迭代的局部性很好,可以获得较高的并行化,那么如何将Jacobi迭代并行化呢?
为了叙述方便,假设处理机台数p=4,我们将从以下两个方面来实现Jacobi迭代的并行化:
1) 数据按行进行分割(结果如图3所示)
假设需要迭代的数据是4M×N的二维数组A(4M,N),按图示将数据按行划分,则分布在四个不同进程上的数据分别是:
进程0,A(1:M,N) 进程1,A(M+1:2M,N)
进程2,A(2M+1:3M, N) 进程3,A(3M+1:4M,N)
由于在迭代过程中,边界点新值的计算需要相邻边界其它块的数据,因此在每一个数据块的上下两侧又增加一行的数据空间,用于存放从相邻数据块通信得到的数据。这样原来每个数据块的大小从M×N扩大到(M+2)×N,进程0和3的数据块只需扩大一块即可满足通信的要求,进程1和2的数据块需要扩大两边数据块即可满足通信的要求,但这里为了编程的方便和形式的一致,在两边都增加了一行数据。进程0和进程3最上侧,最下侧的数据从何处而来,这就需要运用到MPI中的虚拟进程的概念.虚拟进程(MPI_PROC_NULL)是不存在的假想进程,当一个真实进程向虚拟进程发送消息时会立即成功返回;一个真实进程从虚拟进程接收消息时也会立即成功返回,并且对接收缓冲区没有任何改变。
2) 数据的通信
在Jacobi迭代中,每一个进程需要用到上一进程和下一进程边界数据,这就要求每一进程都要向相邻进程发送数据,同时从相邻的进程接收数据.以进程1为例说明边界数据的通信过程.事实上,进程1第一行数据是接收了进程0中的倒数第二行的数据;进程1还将其倒数第二行的数据向进程2发送.同时进程1还将其第二行数据向进程0发送,同时也接收进程2的第二行的数据.数据通信过程用图4表示.
在MPI中提供了捆绑发送和接收操作函数(MPI_SENDRECV),该操作把发送一个目的地和从另一个进程接收一个消息合并到一个调用中,全部通信都用MPI_SENDRECV语句实现,使得程序的设计变得较方便。
4.算 例
曲顶柱体这个名词我们并不陌生,其在数学
上的解释如下:设Z=f(x,y)是定义在平面有界闭
域D上的非负连续函数,其图形为曲面S。曲面
S位于坐标平面OXY上方,即f(x,y)>0,由曲面
S、坐标平面OXY上的区域D和以D的边界C
为准线且母线平行于Z轴的柱面所围成的立体称
为在区域D上以曲面S为顶面的曲顶柱体。
在高等数学求重积分的时候,我们已经研究过曲顶柱体体积的计算,是否曾想过曲顶柱体顶面面积的求解呢?如果在建筑设计上,我们为了追求艺术美观,同时又为了工程上节约顶面材料的用料,客观自然就需要用一种科学的方法来求解曲顶柱体顶面面积的最小值。即我们必须给出具有最小面积的顶面的解析式,并使之与既定的侧面边界吻合.
例:在曲顶柱体中,设
, 顶面的边界 , 其中
利用本文提出的并行算法求解上述算例,将求解得到的最小曲面S绘成图形.(见图5)
为了便于比较,我们又编写了相应的串行程序用于求解最小曲面面积,所用时间对比如表(1)所示.由表(1)可知,并行算法的效率令人满意.
处理机台数 串行(单PC) 并行(4PC)
运行时间 3.48 0.89
加速比 1.00 3.92
表1 并行算法的效率(时间单位:秒)
注:上表中得到的串行时间是在m=n=100时的运算时间。
5.结 束 语
以上给出了一类边界问题的最小曲面的求法,同时给出了相应的并行算法.运用Jacobi迭代基本原理介绍了该迭代法的MPI并行计算方法,解决了曲顶柱体顶面表面积最小值问题的求解问题。我们有理由相信,随着高性能计算机的飞速发展及其应用的日益广泛,并行计算必将在科学研究和工程应用中发挥巨大作用。事实上,就算法而言,仍有两处值得思考和改进,一是迭代函数的选择,二是网格的取法.其实,我们可以选取多重网格,令其逐渐由粗变细,这样可能在整体上加快算法的收敛速度.
(本文得到了上海市重点学科建设项目的资助)
参考文献:
[1]都志辉等编著,高性能计算并行编程技术-MPI并行程序设计,清华大学出版社,P52-63
[2]李晓梅等编著,可扩展并行算法的与分析,国防工业出版社,P157-159
[3]沈志宇等编著,并行编译算法,国防工业出版社,P157-159,P187-199
[4]王晓东等编著,计算机算法设计与分析,电子工业出版社,2001
作者简介:刘婕,1975年生,女,数学系2000级硕士研究生;刘灿文,1978年生,男,数学系2000级硕士研究生,师从赵书钦教授,主要研究计算数学、并行计算研究。
(责任编辑:陈俏巧)