专利名称：Sor的并行通信优化方法
技术领域：
本发明涉及SOR迭代方法并行化中的通信和同步幵销问题。该发明主要包括多层对称
S0R迭代方法的提出、通过交错条块法实现串行程序在高性能网络计算环境下的并行化等方 面的研究。
背景技术：
许多物理应用问题都归结为求微分方程数值解，其核心是高效地求解线性方程组，S0R 迭代方法以其高效性和实现简单的特性而成为大型线性方程组求解的重要方法，并且常用作 非定常迭代方法的高效预条件子，是许多应用软件的核心方法之一。
在并行处理系统中，通信和同步开销是许多并行方法程序设计必须优先考虑的问题。通 信和同步开销越大，程序的并行效率就越低。
S0R方法由于本质上的串行特性，其并行方法的实现一直是数值计算的重要研究对象。 关S0R迭代方法的并行化设计与实现，前人已经作了很多工作，Zhang通过使用基于区域分 解的多色排序方法实现了面向集群的并行SOR方法，但是由于多色排序的缺点，导致其通信 和同步开销很大。Xie通过改进区域分解方法，在保证相容序矩阵的情况下实现了 S0R的并 行化，但对优化同步开销缺乏考虑。Rohallah等人改进了区域分解方法，通过改变每一块数 据的执行顺序，实现SOR的并行化方法。但是，当数据量很大时，该方法同样没有考虑通信 优化的问题。Michelle提出一种空间条块划分方法，通过改变条块的执行顺序，得到了 S0R 的并行化方法，虽然提高了原有方法的数据局部性及并行度，但是该并行化方法只能应用于 共享内存机器。
发明内容
本发明的目的是提出一种解决S0R迭代方法并行化中的通信和同步开销问题的方法。该 发明提出了多层对称S0R迭代方法，通过交错条块方法实现串行方法在高性能网络计算环境 下的并行化，从而达到减少通信启动和同步次数的目的。
为实现上述目的，本发明包括提出多层对称SOR迭代方法，对其描述如下
设一维椭圆偏微分方程为
、=0, ;ceQ =
(1)
其中，Dirichet边界条件为
w(x) = x,xeSQ (2)对求解区域
使用网格划分，用空间间隔为/ = 1/(" + 1)的长度把求解区域网格化，其中，w
为正整数，网格点为jc, "=i7i， z' = 0,…，"十1}。我们用w,(!'-0,L ,")表示w(叫的有限差分近 似值。对(1)进行3点有限差分近似得<formula>formula see original document page 5</formula>
考虑边界值条件x。 = 0和、+1 = 1 ，将差分格式(3)使用矩阵形式表示为
<formula>formula see original document page 5</formula>其中，A
<formula>formula see original document page 5</formula>多层对称SOR方法是通过使用SOR迭代求解方法，对线性迭代方程求数值解。其中， 奇数K次([(p-lVC + l,;)*/:] ， p取奇数)迭代采用"向前"的SOR迭代方法，由(5)式计
算；偶数K次([(/M)《+ 1,"Z]， p取偶数)迭代采用"向后"的SOR迭代方法，由(6)式计算。
如图1所示，当迭代次数Ae[l，幻时，对所有网格点使用式(5)计算，网格点更新顺序如 箭头指向。类似的，当迭代次数"[《+ l，2i:]时，对所有网格点使用式(6)计算，网格点更新 顺序为由右向左。
传统的迭代方法执行特点是依据网格点的顺序依次对所有网格点进行迭代更新操作。
这种"串行"执行的特点阻碍了 SOR方法并行化执行。为了提高SOR方法的并行效率，本
发明包括改变以往传统的以迭代次序为执行序的特点，将迭代空间条块引入执行序中，通
过空间条块内部执行多次迭代，从而减少通信和同步开销。
该发明考虑将1维空间网格点和1维时间点组成1+1维的迭代空间网格系统
5^wc/ej卯ce(Z，r) ， Z为整数空间，r表示时间空间。其中网格点Pnode定义为 Pnode(x，f)^、(/OI"eZ,^:r"， x表示一维空间的网格点，Z表示时间维，、(0表示1 + 1维
迭代空间内的网格点。
以一维椭圆微分方程为例，第一阶段，迭代空间交错条块串行SOR方法过程如下-步骤l.网格块划分
用区域分解方法将网格区域《,。(^(2,|在7 = 0处进行等份划分。如图2所示，划分后
使得每块网格5^—j/wce(>U = 0)(其中，lSw^i:， w表示f = 0时的第附个网格块，
丄为/ = 0时的总块数)中的网格点数为R， R满足式(7):
i > 2《 (7) 其中，I为单向迭代次数，X的值为cache大小的一半。 步骤2.在网格块的基础上，沿时间轴对网格区域进行划分
划分方法采用时滞技术，对每层网格子区域修正其边界，图2中黑线表示修整后的边界。 修正方法描述如下
a) 对于奇数K层循环中的第一个网格块
Sub一Snode—space(m=l,t=k)=Sub—Snode—space(m=l，t=k-l)-Pynode( )， 其中 ,
Pynode(m,A:-l)为第m个网格块在第k-l次迭代时间步的右边界点；
b) 对于奇数X层循环中的最后 一 个网格块，即第L个网格块 Sub—Snode—space(m=L,t=k)=Sub—Snode_space(m=L,t=k-l)+Pynode( — l )， 其中 ，
Pynode(w-1,^-1 )为第m-l个网格块在第k-l次迭代时间步的右边界点；
c) 对于奇数《层循环中的除第一个和最后一个网格块外的其他网格块 Sub—Snode一space(m,t:k)-Sub—S加de一space(m,t-k-1)+ Pynode( m -1,_ 1 )-Pynode(m, — 1);
d) 对于偶数K层循环中的最后一个，即第L个网格块
Sub—Snode_space(m=L,t=k)=Sub—Snode—space(m=L,t=k-l)-Pznode( w，A: —1 )，其中 ，
Pznode(wJ-l)为第m个网格块在第k-l次迭代时间步的左边界点；
e) 对于偶数X层循环中的第 一 个网格块Sub—Snode—space(m=l,t=k)= Sub—Snode_space(m=l,t=k-l)+Pznode(w + l,yt —1)，其中，Pznode(m + l，& —1 )为第m个网格块 在第k-l次迭代时间步的左边界点；f)对于偶数AT层循环中的除第一个和最后一个网格块外的其他网格块 Sub—Snode_space(m，t=k)= Sub—Snode—space(m,t=k-1)+ Pznode(m + l,A:-1) - Pznode(/n,A: — 1);
以上修正中，1《Asi:。
边界修正后，在奇数k次迭代中，网格块和相应的边界点组成了下层空间网格条块，如
b—tilel b—tileL.在偶数《次迭代中，网格块和相应的边界点组成了上层网格条块，如 t—tilel t—tileL.
步骤3.按空间网格条块顺序执行SOR方法
网格条块生成后，按条块顺序执行。以条块内部网格点层作为内部执行序，条块序作为 外部执行序，其执行顺序描述如下
对于奇数k层循环中，依次执行块b—tilel b—tileL，每个块内，按时间步顺序执行，每个 时间步是按照从左到右的顺序进行计算。
对于偶数X层循环中，依次执行块t—tilel~t_tileL，每个块内，按时间步顺序执行，每个
时间步是按照从右到左的顺序进行计算。
为了实现网格条块的并行化，需要建立网格条块之间的依赖关系。条块SOR的并行度可 通过依赖关系表示。图3、图4给出了二维(空间维+时间维)迭代空间的两种块划分，对应 的网格条块依赖关系如图5和图6所示。图5显示网格条块按条块顺序 (b—tilel b—tile4 t一tile4 tjilel)串行执行更新，因此没有并行度。图6显示网格条块可以通过 两个进程并行执行迭代更新，因此并行度为2，进程1依次执行b一tilel、 b_tile3、 t_tile3、 t一tilel 条块，进程2依次执行b—tile2、 b_tile4、 t_tile4、 t—tile2条块.第一种划分方式的执行速度是第 二种划分方式的两倍。显然，这两种方式最大的不同在于初始化时的排序方式不同。排序方 式的不同代表了相邻条块间的执行顺序的差异，因此条块排序影响了条块间的数据依赖关系. 通过改变网格条块间数据依赖关系，提高条块执行的并行度。交错条块方法利用图着色原理 对网格块进行排序，相邻的网格块使用不同颜色，并对相同颜色使用连续的数字标示。这样， 同一颜色的条块不相邻，因此相邻标号的网格条块间没有数据依赖关系(如图6所示)。
为了方便说明，将迭代空间S^c^"划分为两个子空间sub—domainl和SUb_domain2。
为了保证CPU之间的负载均衡，sub—domainl的网格点数nl和sub—domain2的网格点数n2 必须满足式(8)和式(9):
<formula>formula see original document page 7</formula> (9) 网格条块的计算量可认为是条块内部网格点的函数，因此，b—tilel的计算量使用函数、 Y )表示，b_tile2的计算量使用函数/(A-(l +〖"iQ表示。s和t分别代表
b—tilel和b—tile2的网格点数。从负载均衡的角度考虑第一块与第二块计算量必须相等，所以 s和t必须满足一下关系
s"-(X十l)/2 (10)
空间条块的边界网格点必须发送给其它相邻子空间以维护数据依赖关系.子空间中的条 块分为三类b—tile2和t—tile3为"发送条块"。发送条块在执行更新后将自身边界网格 数据发送给"接收条块"。条块b—tile3与t—tile2为"接收条块"。接收条块必须在接受其 它子空间条块发送的边界网格数据后才执行更新。b—tilel、 b_Ule4、 t—tilel和t—tile4 为"非通信条块"。其自身计算不需要其它条块的边界值.在交错条块算法中，条块的执行序 为先执行"发送条块"，再执行"接收条块"。"发送条块"完成更新后才发送边界网格数据， 这样可以减少处理器之间的通信开销。
第二阶段，并行执行实现步骤描述如下
歩骤4.依据处理器数目g，将空间计算区域划分为《个子空间sub一domain。第一块子空 间的大小和最后一块子区域的大小通过添加和减少K/2个数据进行修正，其它子区域的大小 为"/《。
步骤5,该方法分奇数和偶数《次迭代，并且奇数K次迭代与偶数《次迭代执行方向相 反。例如，使用LR顺序(左到右)执行奇数X次迭代，使用RL顺序(右到左)执行偶数K 次迭代。
步骤6.各子区域通过时间轴进行划分，同串行执行的步骤2 步骤7.对网格条块重新排序
根据多色排序方法将所有的网格块进行重新排序，相邻的网格块使用不同颜色，并对相 同颜色使用连续的数字标示。
步骤8.以网格条块为单位执行奇数《次迭代数据更新
首先，执行奇数《次迭代更新，更新所有sub—domain中的"发送条块"，以及第一个条 块然后，将边界数据发送给相邻sub—domahi的"接收条块"；最后，更新所有sub一domain 中的"接收条块"，以及最后一个条块。
步骤9.以网格条块为单位执行偶数尺次迭代数据更新
首先，执行偶数《次迭代更新，更新所有sub—domain中的"发送条块"，以及最后一个 条块；然后，将边界数据发送给相邻sub—domain的"接收条块"；最后，更新所有sub—domain中的"接收条块"，以及第一个条块。
基于上述技术方案，本发明有效地减少了通信启动和同步次数。该并行方法比传统的区 域分解方法和红黑排序并行方法具有更好的并行效率和可扩展性。
图1为本发明的多层对称S0R方法的示意图。
图2为本发明的迭代空间交错条块S0R串行方法的示意图。
图3为条块原始执行顺序。 '
图4为木发明的网格条块重排后的执行顺序。
图5为图3的条块关系依赖图。
图6为本发明的图4的条块关系依赖图。
具体实施例方式
下面通过附图和实例，对本发明的技术方案做进一步的详细描述。
本发明包括三个方面提出多层对称SOR法，引入时滞技术，通过对网格条块重排，实 现S0K的并行化。
图1出示了多层对称SOR的执行过程的示意图。如图1所示，其中，奇数《次 (
，n取奇数)迭代采用"向前"的S0R迭代方法；偶数K次
([("+ l口"*幻，n取偶数)迭代采用"向后"的S0R迭代方法。例如当迭代次数A e [1，幻
时，网格点更新顺序如箭头指向。类似的，当迭代次数"[K + l，2jq时，网格点更新顺序为
由右向左。该方法改善了迭代间数据局部性，从而提高cache命中率。
图2出示了迭代空间交错条块S0R串行执行过程示意图。如图2所示，迭代空间交错条 块SOR串行执行与传统的迭代方法迥然不同。
传统的迭代方法执行特点是依据网格点的顺序依次对所有网格点进行迭代更新操作。 这种"串行"执行的特点阻碍了SOR并行化执行。为了提高SOR的并行效率，本发明提出新 的串行SOR执行方式。迭代空间交错条块方法的核心思想是改变以往传统的以迭代次序为 执行序的特点，将迭代空间条块引入执行序中，通过空间条块内部执行多次迭代，从而减少 通信和同步开销。考虑将1维空间网格点和1维时间点组成1 + 1维的迭代空间网格系统
S"Me—^ace(ZJ) ， Z为整数空间，r表示时间空间。其中网格点Pnode定义为 Pnode(x，f)={、(r)|(xsZ,fer)}， x表示一维空间的网格点，/表示时间维，w力)表示l + l维 迭代空间内的网格点。传统的并行SOR方法在一次更新未知数过程中，遍历并更新全部网格点。处理器需要通 信操作以得到必要的边界数据，迭代之间需要同步操作来维护计算的正确性。我们采用对迭 代空间进行时间轴方向划分的方法，实现对同一网格块进行递归式多个时间步更新，从而减 少通信和同步开销。
第一阶段，迭代空间交错条块串行S0R方法过程如下
步骤1.网格块划分
用区域分解方法将网格区域&<w(u)在T^处进行等份划分。如图2所示，划分后使
得每块网格Sub—Snode—space (m， t=0)(其中，lSm^Z)中的网格点数为R， R满足式(1):
i >2X (1)
其中，X为单向迭代次数。
步骤2.在网格块的基础上，沿时间轴对网格区域进行划分。
划分方法采用时滞技术，对每层网格子区域修正其边界，图2中黑线表示修整后的边界。 具体的修正方法如下描述-
a) 对于奇数K层循环中的第一个网格块
Sub—Snode—space(m=l，t=k)=Sub—Snode—space(m=l,t=k-l)-Pynode( )， 其中 ，
Pynode(附,yt-l)为第m个网格块在第k-l次迭代时间步的右边界点；
b) 对于奇数《层循环中的最后一个网格块，即第L个网格块 Sub—Snode_space(m=L，t=k)=Sub—Snode—space(m=L,t=k-l)+Pynode( 附-1,&_1)， 其中，
Pynode(m-lA-l)为第m-l个网格块在第k-l次迭代时间步的右边界点；
c) 对于奇数《层循环中的除第一个和最后一个网格块外的其他网格块 Sub—Snode_space(m，t=k)=Sub_Snode—space(m,t=k-1)+ Pynode( 7 -1, A: 一 1 )-Pynode( A: — 1);
d) 对于偶数K层循环中的最后一个，即第L个网格块 Sub—Snode—space(m=L,t=k)=Sub—Snode—space(m=L,t=k-l)-Pznode( w，t —1 )，其 中 ，
Pznode(m，A-l)为第m个网格块在第k-l次迭代时间步的左边界点；
e) 对于偶数《层循环屮的第 一 个网格块Sub_Snode_space(m=l，t=k)= Sub—Snode—space(m=l,t=k-l)+Pznode(> + l,A: —1)，其中，Pznode(m + l，/t-l)为第m个网格块
在第k-l次迭代时间步的左边界点；
f) 对于偶数《层循环中的除第一个和最后一个网格块外的其他网格块Sub—Snode—space(m，t=k)= Sub—Snode—space(m,t=k-1)+ Pznode( m +1, A: -1) - Pznode( m, & — 1);
边界修正后，在奇数《次迭代中，网格块和相应的边界点组成了下层空间网格条块，如 b—tilel b—tileL。在偶数尤次迭代中，网格块和相应的边界点组成了上层网格条块，如 t—tilert—tileL。
步骤3.按空间网格条块顺序执行SOR
网格条块生成后，按条块顺序执行SOR，以条块内部网格点层作为内部执行序，条块序 作为外部执行序。
为了实现网格条块的并行化，需要建立网格条块之间的依赖关系。条块SOR的并行度可 通过依赖关系表示。
图3和图4给出了二维(空间维+吋间维)迭代空间的两种交错条块划分，对应的网格条 块依赖关系如图5和图6所示。第 一 个依赖图显示网格条块按条块顺序 (b—tile广b一tile4、一tile4、—tilel)串行执行更新，因此没有并行度。第二个依赖图显示网 格条块可以通过两个进程并行执行迭代更新，因此并行度为2，进程1依次执行b—tilel 、 b—tile3、 t—tile3、 t—tilel条块，进程2依次执行b—tile2、 b—tile4、 t—tilc4、 t—tile2 条块。第一种划分方式的执行速度是第二种划分方式的两倍。显然，这两种方式最大的不同 在于初始化时的排序方式不同。排序方式的不同代表了相邻条块间的执行顺序的差异，因此 条块排序影响了条块间的数据依赖关系。通过改变网格条块间数据依赖关系，提高条块执行 的并行度。交错条块方法利用图着色原理对网格块进行排序，相邻的网格块使用不同颜色， 并对相同颜色使用连续的数字标示。这样，同 -颜色的条块不相邻，因此相邻标号的网格条 块间没有数据依赖关系(如图4所示)。
第二阶段，并行执行的具体实现步骤描述如下
步骤4.依据处理器数目《，将空间计算区域划分为《个子空间sub—domain。第一块子空
间的大小和最后一块子区域的大小通过添加和减少K/2个数据进行修正，其它子区域的大小 为n/q。其中，n代表迭代空间网格点数。
步骤5.交错条块S0R方法分奇数和偶数《次迭代，并且奇数《次迭代与偶数r次迭代 执行方向相反。例如，使用LR顺序(左到右)执行奇数〖次迭代，使用RL顺序(右到左) 执行偶数《次迭代。
步骤6.各子区域通过时间轴进行划分，同串行执行的步骤2
步骤7.对网格条块重新排序
根据多色排序方法将所有的网格块进行重新排序，相邻的网格块使用不同颜色，并对相同颜色使用连续的数字标示。
步骤s.以网格条块为单位执行奇数/a欠迭代数据更新
首先，执行奇数《次迭代更新，更新所有sub一domain中的"发送条块"，以及第一个条 块；然后，将边界数据发送给相邻sub—domain的"接收条块"；最后，更新所有sub—domain 中的"接收条块"，以及最后一个条块。
步骤9.以网格条块为单位执行偶数K次迭代数据更新
首先，执行偶数《次迭代更新，更新所有sub一domain中的"发送条块"，以及最后条块； 然后，将边界数据发送给相邻sub—domain的"接收条块"；最后，更新所有sub—domain中的 "接收条块"，以及第一个条块。
基于上述技术方案，本发明有效地减少了通信启动和同步次数。该并行方法比传统的区 域分解方法和红黑排序并行方法具有更好的并行效率和可扩展性。
权利要求
1、SOR的并行通信优化方法，其特征在于优化方法分为两个阶段顺序进行，包括串行执行阶段和并行执行阶段；将1维空间网格点和1维时间点组成1+1维的迭代空间网格系统Snode_space(Z，T)，Z为一维的网格空间，T表示时间空间；其中网格点Pnode定义为Pnode(x，t)＝{ux(t)|(x∈Z，t∈T)}，x表示一维网格点，t表示时间维，ux(t)表示1+1维迭代空间内的网格点；第一阶段串行执行阶段步骤1. 网格块划分用区域分解方法将网格区域Snode_space(Z，T)在T＝0处进行等份划分，划分后使得每块网格Sub_Snode_space(m，t＝0)中的网格点数为R，1≤m≤L，m表示t＝0时的第m个网格块，L为t＝0时的总块数，R满足R>2K，K为单向迭代次数，K的值为cache大小的一半；步骤2. 在网格块的基础上，沿时间轴对网格区域进行划分划分方法采用时滞技术，对每层网格子区域修正其边界，修正方法描述如下a)对于奇数K层循环中的第一个网格块Sub_Snode_space(m＝1，t＝k)＝Sub_Snode_space(m＝1，t＝k-1)-Pynode(m，k-1)，其中，Pynode(m，k-1)为第m个网格块在第k-1次迭代时间步的右边界点；b)对于奇数K层循环中的最后一个网格块，即第L个网格块Sub_Snode_space(m＝L，t＝k)＝Sub_Snode_space(m＝L，t＝k-1)+Pynode(m-1，k-1)，其中，Pynode(m-1，k-1)为第m-1个网格块在第k-1次迭代时间步的右边界点；c)对于奇数K层循环中的除第一个和最后一个网格块外的其他网格块Sub_Snode_space(m，t＝k)＝Sub_Snode_space(m，t＝k-1)+Pynode(m-1，k-1)-Pynode(m，k-1)；d)对于偶数K层循环中的最后一个，即第L个网格块Sub_Snode_space(m＝L，t＝k)＝Sub_Snode_space(m＝L，t＝k-1)-Pznode(m，k-1)，其中，Pznode(m，k-1)为第m个网格块在第k-1次迭代时间步的左边界点；e)对于偶数K层循环中的第一个网格块Sub_Snode_space(m＝1，t＝k)＝Sub_Snode_space(m＝1，t＝k-1)+Pznode(m+1，k-1)，其中，Pznode(m+1，k-1)为第m个网格块在第k-1次迭代时间步的左边界点；f)对于偶数K层循环中的除第一个和最后一个网格块外的其他网格块Sub_Snode_space(m，t＝k)＝Sub_Snode_space(m，t＝k-1)+Pznode(m+1，k-1)-Pznode(m，k-1)；以上修正中，1≤k≤K；边界修正后，在奇数K次迭代中，网格块和相应的边界点组成了下层空间网格条块b_tile1～b_tileL，在偶数K次迭代中，网格块和相应的边界点组成了上层网格条块t_tile1～t_tileL；步骤3. 按空间网格条块顺序执行SOR方法网格条块生成后，按条块顺序执行，以条块内部网格点层作为内部执行序，条块序作为外部执行序，其执行顺序描述如下对于奇数K层循环中，依次执行块b_tile1～b_tileL，每个块内，按时间步顺序执行，每个时间步是按照从左到右的顺序进行计算；对于偶数K层循环中，依次执行块t_tile1～t_tileL，每个块内，按时间步顺序执行，每个时间步是按照从右到左的顺序进行计算；第二阶段，并行执行实现步骤描述如下步骤4. 依据处理器数目q，将空间计算区域划分为q个子空间sub_domain；第一块子空间的大小和最后一块子区域的大小通过添加和减少K/2个数据进行修正，其它子区域的大小为n/q；步骤5. 该方法分奇数和偶数K次迭代，并且奇数K次迭代与偶数K次迭代执行方向相反；步骤6. 各子区域通过时间轴进行划分，同步骤2执行；步骤7. 对网格条块重新排序根据多色排序方法将所有的网格块进行重新排序，相邻的网格块使用不同颜色，并对相同颜色使用连续的数字标示；步骤8. 以网格条块为单位执行奇数K次迭代数据更新首先，执行奇数K次迭代更新，更新所有sub_domain中的“发送条块”，以及第一个条块；然后，将边界数据发送给相邻sub_domain的“接收条块”；最后，更新所有sub_domain中的“接收条块”，以及最后一个条块；步骤9. 以网格条块为单位执行偶数K次迭代数据更新首先，执行偶数K次迭代更新，更新所有sub_domain中的“发送条块”，以及最后一个条块；然后，将边界数据发送给相邻sub_domain的“接收条块”；最后，更新所有sub_domain中的“接收条块”，以及第一个条块。
全文摘要
SOR的并行通信优化方法，本发明针对传统的SOR迭代方法并行化中的通信和同步开销问题，提出了多层对称SOR迭代方法，并通过交错条块方法实现串行执行在高性能网络计算环境下的并行化。交错条块方法包括引入了时滞技术改进区域分解方法，实现空间网格的划分；并通过网格条块重排序，改变网格条块间的数据依赖关系，提高条块执行的并行度，较大幅度地减少了通信启动和同步次数。因此，本发明比传统的区域分解方法和红黑排序方法具有更好的并行效率和可扩展性。
文档编号G06F17/11GK101430654SQ20081023958
公开日2009年5月13日 申请日期2008年12月12日 优先权日2008年12月12日
发明者张继林, 李建江, 胡长军, 赵冲冲 申请人:北京科技大学