并 行 计 算中国科学技术大学计算机科学与技术系国家高性能计算中心 (合肥 )
2003年 9月第三篇 并行数值算法第八章 基本通讯操作第九章 稠密矩阵运算第十章 线性方程组的求解第十一章 快速傅里叶变换第十章 线性方程组的求解
10.1 三角形方程组的求解
10.2 三对角方程组的求解
10.3 稠密线性方程组的求解
10.4 稀疏线性方程组的求解
10.1 三角形方程组的求解
10.1.1 基本术语
10.1.2 上三角方程组的求解国家高性能计算中心（合肥） 52009-7-24
基本术语
线性方程组的定义和符号
a1,1x1 + a1,2x2 + … + a1,nxn = b1
a2,1x1 + a2,1x2 + … + a2,nxn = b2
an,1x1 + an,1x2 + … + an,nxn = bn
记为 AX=b
nnnnnn
n
n
b
b
b
b
x
x
x
x
aaa
aaa
aaa
A
2
1
2
1
21
22221
11211
,,
10.1 三角形方程组的求解
10.1.1 基本术语
10.1.2 上三角方程组的求解国家高性能计算中心（合肥） 72009-7-24
上三角方程组的求解
上三角方程组的回代解法并行化
(1)SISD上的回代算法
Begin
(1)for i=n downto 1 do
(1.1)xi=bi/aii
(1.2)for j=1 to i-1 do
bj=bj-ajixi
aji=0
endfor
endfor
End
可并行化国家高性能计算中心（合肥） 82009-7-24
上三角方程组的求解
上三角方程组的回代解法并行化
(2)SIMD-CREW上的并行回代算法
-划分,p个处理器行循环带状划分
-算法
Begin
for i=n downto 1 do
xi=bi/aii
for all Pj,where 1≤j≤p do
for k=j to i-1 step p do
bk=bk-akixi
aki=0
endfor
endfor
endfor
End // p(n)=n,t(n)=n
P
1
P
j
P
1
P
j
行
1
j
1 + p
j + p
第十章 线性方程组的求解
10.1 三角形方程组的求解
10.2 三对角方程组的求解
10.3 稠密线性方程组的求解
10.4 稀疏线性方程组的求解
10.2 三对角方程组的求解
10.2.1 直接求解法
10.2.2 奇偶规约法国家高性能计算中心（合肥） 112009-7-24
i
i
待 改 变 的 元 素不 变 化 的 元 素待 消 为 0 的 元 素已 消 为 0 的 元 素主 元 a
ii
三对角方程组的直接求解法
Gauss消去法 (难以并行化 )
① 消元
② 回代注：由于三对角的方程组的特殊性，
一次消元或一次回代，只涉及邻近一个方程，故难以并行化。
10.2 三对角方程组的求解
10.2.1 直接求解法
10.2.2 奇偶规约法国家高性能计算中心（合肥） 132009-7-24
三对角方程组的直接求解法
奇偶规约求解法 (可并行化 )
三对角方程可以写成如下形式
fixi-1+gixi+hixi+1=bi i=1~n
f1=hn=0
串行算法描述
① 利用上下相邻方程消去偶序号方程中的奇下标变量,
f2i-1x2i-2+g2i-1x2i-1+h2i-1x2i =b2i-1
f2ix2i-1 + g2ix2i + h2ix2i+1 =b2i
f2i+1x2i +g2i+1x2i+1+h2i+1x2i+2 = b2i+1
2i-1方程乘上某个数消去 2i方程中的 f2ix2i-1项,2i+1方程乘上某个数消去 2i方程中的 h2ix2i+1项,使 2i方程变为
α ix2i-2+β ix2i+γ ix2i+2=η i i=1,2,…,n/2
国家高性能计算中心（合肥） 142009-7-24
三对角方程组的直接求解法
② 重复 ① 最终可得,
case 1,case 2,
g1x1+h1x2 =b1,
f2x1+g2x2+h2x3 =b2,
f3x2 +g3x3+h3x4=b3,
f4x3 +g4x4 =b4,
可以分别得到
g1x1+h1x2 =b1 或 g1x1+h1x2 =b1
f2x1+g2x2 =b2 f2x1+g2x2+h2x3 =b2
f3x2+g3x3 =b3
解得 x1,x2 或 x1,x2,x3
③ 回代求解 x
并行化分析,①,② 消去奇下标可以并行化；
③ 回代求解可以并行化第十章 线性方程组的求解
10.1 三角形方程组的求解
10.2 三对角方程组的求解
10.3 稠密线性方程组的求解
10.4 稀疏线性方程组的求解
10.3 稠密线性方程组的求解
10.3.1 有回代的高斯消去法
10.3.2 无回代的高斯 -约旦法
10.3.3 迭代求解的高斯 -赛德尔法国家高性能计算中心（合肥） 172009-7-24
有回代的高斯消去法
算法基本原理
求解过程分为消元和回代两个阶段，消元是将系数矩阵
A化为上三角阵 T，然后对 TX=c进行回代求解。
消元过程中可以应用选主元方法，增加算法的数值稳定性。
下面是消元过程图：
i
i
主 元 行图 1 0,1
待 改 变 的 元 素不 变 化 的 元 素待 消 为 0 的 元 素已 消 为 0 的 元 素国家高性能计算中心（合肥） 182009-7-24
有回代的高斯消去法
并行化分析
消元和回代均可以并行化；
选主元也可以并行化；
消元过程的并行化图示：处理器按行划分
i
i
主 元 行图 1 0,1
待 改 变 的 元 素不 变 化 的 元 素待 消 为 0 的 元 素已 消 为 0 的 元 素可 并 行 化 部 分
10.3 稠密线性方程组的求解
10.3.1 有回代的高斯消去法
10.3.2 无回代的高斯 -约旦法
10.3.3 迭代求解的高斯 -赛德尔法国家高性能计算中心（合肥） 202009-7-24
无回代的高斯 -约旦法
串行算法原理
① 消元,通过初等行变换,将 (A,b)化为主对角线矩阵,
(方便起见,记 b为 A的第 n+1列 )
② 求解,xj=a’j,n+1/a’jj
1,
1,222
1,111
1,
1,2
1,1
21
22221
11211
nnnn
n
n
nn
n
n
nnnn
n
n
aa
aa
aa
a
a
a
aaa
aaa
aaa
国家高性能计算中心（合肥） 212009-7-24
无回代的高斯 -约旦法
SIMD-CREW上的并行算法
(1)处理器,n2+n个处理器,这些处理器排成 n× (n+1)的矩阵,
处理器编号为 Pik,i=1~n,k=1~n+1
(2)并行化分析
① 消元的并行化,// O(n)
for j=1 to n-1,each Pik Par-do //第 j次消元
Pij(i<>j),aij <— 0
Pik(i<>j,k=j+1~n+1),aik <— aik-ajk(aij/ajj)
end for
② 求解,for each Pjj(j=1~n) Par-do,xj <— aj,n+1/ajj //O(1)
(3)时间分析,t(n)=O(n),p(n)=O(n2),c(n)=O(n3) 成本最优？
国家高性能计算中心（合肥） 222009-7-24
无回代的高斯 -约旦法
成本最优？
串行算法的最优时间,由于 x=A-1b
① A-1b(假设已有 A-1),O(n2)
② 求 A-1,
∴ 求 A-1需要,2次 n/2× n/2矩阵的逆 i(n/2)
6次 n/2× n/2矩阵的乘 m(n/2)
2次 n/2× n/2矩阵的加 a(n/2)
i(n)=i(n/2)+6m(n/2)+2a(n/2)
a(n/2)=n2/2,m(n/2)=O((n/2)x) 2<x<2.5
=> i(n)=O(nx) 综上，串行算法的最优时间为 O(nx) 2<x<2.5
12
1
1121221
1121
1
1
12
1
111
12
1
1111
1
11212221
1211
0
0
0
0
00
00
11 AAAAB
IAA
I
B
A
I
AAI
A
I
AAI
B
A
IAA
I
AA
AA
A
其中，
10.3 稠密线性方程组的求解
10.3.1 有回代的高斯消去法
10.3.2 无回代的高斯 -约旦法
10.3.3 迭代求解的高斯 -赛德尔法国家高性能计算中心（合肥） 242009-7-24
迭代求解的高斯 -赛德尔法
串行算法原理如果对某个 k,给定的误差允许值 c有则认为迭代是收敛的。
并行化分析由于每次迭代需要使用本次迭代的前面部分值，因而难以得到同步的并行算法，下面给出一个异步的并行算法。
i
ij
k
jij
ij
k
jij
ii
k
i bxaxaax
11 1
cxx
n
i
k
i
k
i
1
1
国家高性能计算中心（合肥） 252009-7-24
迭代求解的高斯 -赛德尔法
MIMD异步并行算法
N个处理器 (N≤n) 生成 n个进程,每个进程计算 x的一个分量
算法
Begin
(1)oldi? xi0,newi?xi0
(2)生成进程 i
(3)进程 i
repeat
(i) oldi?newi
(ii)newi? (bi-∑ k<iaik× oldk-∑ k>iaik× oldk)/aii
until ∑ i=1~n| oldi - newi |<c
xi? newi
End
第十章 线性方程组的求解
10.1 三角形方程组的求解
10.2 三对角方程组的求解
10.3 稠密线性方程组的求解
10.4 稀疏线性方程组的求解
10.4 稀疏线性方程组的求解
10.4.1 线性方程组的并行化方法
10.4.2 稀疏线性方程组的迭代解法
10.4.3 高斯 -赛德尔迭代法的并行化国家高性能计算中心（合肥） 282009-7-24
线性方程方程的并行化方法
线性方程组选择算法的考虑因素
系数矩阵 A的结构
dense Gaussian elimination,etc
Sparse iterative method
triangular substitution,odd-even reduction
certain PDEs multigrid,etc
计算精度要求
Gaussian elimination,more accurate,more expensive
Conjugate gradients,less accurate,less expensive
计算环境要求
architecture,available languages,compiler quality
libraries?
国家高性能计算中心（合肥） 292009-7-24
线性方程方程的并行化方法
求解方法的并行化
(1)直接解法的并行化 (用于稠密线性方程组 )
- Gauss消去法 (包括选主元的 Gauss消去法 )
- Gauss-Jordan消去法
- LU分解法
(2)迭代法的并行化 (用于稠密和稀疏线性方程组 )
- Jacobi
- Gauss-Seidel(可异步并行化 )
- Jacobi OverRelaxation(JOR)
- Gauss-Seidel OverRelaxation(SOR)
- Conjugate Gradient
10.4 稀疏线性方程组的求解
10.4.1 线性方程组的并行化方法
10.4.2 稀疏线性方程组的迭代解法
10.4.3 高斯 -赛德尔迭代法的并行化国家高性能计算中心（合肥） 312009-7-24
稀疏线性方程方程的迭代解法
迭代解法
G r a d i e n tC o n j u g a t e
axaxabxxSOR
axabxxJ O R
axaxabxS e i d e lG a u s s
axabxJ a c o b i
ii
ij
k
jij
ij
k
jiji
k
i
k
i
ii
ij
k
jiji
k
i
k
i
ii
ij
k
jij
ij
k
jiji
k
i
ii
ij
k
jiji
k
i
)5(
/)()1(:)4(
/)()1(:)3(
/)(:)2(
/)(:)1(
11
1
11
1
10.4 稀疏线性方程组的求解
10.4.1 线性方程组的并行化方法
10.4.2 稀疏线性方程组的迭代解法
10.4.3 高斯 -赛德尔迭代法的并行化国家高性能计算中心（合肥） 332009-7-24
高斯 -赛德尔迭代法的并行化
由 PDE离散产生的稀疏线性方程组
(1)Laplace方程
),(),(
]1,0[]1,0[),(),(),(),( 2
2
2
2
yxyxu
Dyxyxg
y
yxu
x
yxu
且满足边值满足上面的方程，求
u =?
0 1 2 3
4 5 6 7
8 9 1 0 1 1
1 2 1 3 1 4 1 5
u = φ
h 步 长 ＝ 1 / N
国家高性能计算中心（合肥） 342009-7-24
高斯 -赛德尔迭代法的并行化
(2)由五点格式的离散化 (假设 g(x,y)=0)
为单位阵方程代入用二阶中心差分逼近：
，记
ER
RE
ERE
ERE
ER
A
NjiuuuuuL a p l a c e
hyxuyxuhyxu
hy
u
yhxuyxuyhxu
hx
u
Njiuu
jijijijiij
N
i
N
i
ij
,
41
141
141
14
,
,0][
4
1
)],(),(2),([
1
)],(),(2),([
1
,0),(
1,1,,1,1
2
2
2
2
2
2
国家高性能计算中心（合肥） 352009-7-24
高斯 -赛德尔迭代法的并行化
A为稀疏的块三对角矩阵国家高性能计算中心（合肥） 362009-7-24
高斯 -赛德尔迭代法的并行化
Gauss-Seidel迭代解法的并行化
(1)两种串行算法的计算顺序及其并行化顺序 1 顺序 2
注,顺序 1难以并行化；顺序 2可以小规模并行化
0 1 2 3
4 5 6 7
8 9 1 0 1 1
1 2 1 3 1 4 1 5
0 1 2 3
4 5 6 7
8 9 1 0 1 1
1 2 1 3 1 4 1 5
国家高性能计算中心（合肥） 372009-7-24
高斯 -赛德尔迭代法的并行化
(2) 红黑着色并行算法
① 并行计算所有的红点
② 并行计算所有的黑点
③ 重复 ①,② 直至满足精度要求
0 1 2 3
4 5 6 7
8 9 1 0 1 1
1 2 1 3 1 4 1 5
2003年 9月第三篇 并行数值算法第八章 基本通讯操作第九章 稠密矩阵运算第十章 线性方程组的求解第十一章 快速傅里叶变换第十章 线性方程组的求解
10.1 三角形方程组的求解
10.2 三对角方程组的求解
10.3 稠密线性方程组的求解
10.4 稀疏线性方程组的求解
10.1 三角形方程组的求解
10.1.1 基本术语
10.1.2 上三角方程组的求解国家高性能计算中心（合肥） 52009-7-24
基本术语
线性方程组的定义和符号
a1,1x1 + a1,2x2 + … + a1,nxn = b1
a2,1x1 + a2,1x2 + … + a2,nxn = b2
an,1x1 + an,1x2 + … + an,nxn = bn
记为 AX=b
nnnnnn
n
n
b
b
b
b
x
x
x
x
aaa
aaa
aaa
A
2
1
2
1
21
22221
11211
,,
10.1 三角形方程组的求解
10.1.1 基本术语
10.1.2 上三角方程组的求解国家高性能计算中心（合肥） 72009-7-24
上三角方程组的求解
上三角方程组的回代解法并行化
(1)SISD上的回代算法
Begin
(1)for i=n downto 1 do
(1.1)xi=bi/aii
(1.2)for j=1 to i-1 do
bj=bj-ajixi
aji=0
endfor
endfor
End
可并行化国家高性能计算中心（合肥） 82009-7-24
上三角方程组的求解
上三角方程组的回代解法并行化
(2)SIMD-CREW上的并行回代算法
-划分,p个处理器行循环带状划分
-算法
Begin
for i=n downto 1 do
xi=bi/aii
for all Pj,where 1≤j≤p do
for k=j to i-1 step p do
bk=bk-akixi
aki=0
endfor
endfor
endfor
End // p(n)=n,t(n)=n
P
1
P
j
P
1
P
j
行
1
j
1 + p
j + p
第十章 线性方程组的求解
10.1 三角形方程组的求解
10.2 三对角方程组的求解
10.3 稠密线性方程组的求解
10.4 稀疏线性方程组的求解
10.2 三对角方程组的求解
10.2.1 直接求解法
10.2.2 奇偶规约法国家高性能计算中心（合肥） 112009-7-24
i
i
待 改 变 的 元 素不 变 化 的 元 素待 消 为 0 的 元 素已 消 为 0 的 元 素主 元 a
ii
三对角方程组的直接求解法
Gauss消去法 (难以并行化 )
① 消元
② 回代注：由于三对角的方程组的特殊性，
一次消元或一次回代，只涉及邻近一个方程，故难以并行化。
10.2 三对角方程组的求解
10.2.1 直接求解法
10.2.2 奇偶规约法国家高性能计算中心（合肥） 132009-7-24
三对角方程组的直接求解法
奇偶规约求解法 (可并行化 )
三对角方程可以写成如下形式
fixi-1+gixi+hixi+1=bi i=1~n
f1=hn=0
串行算法描述
① 利用上下相邻方程消去偶序号方程中的奇下标变量,
f2i-1x2i-2+g2i-1x2i-1+h2i-1x2i =b2i-1
f2ix2i-1 + g2ix2i + h2ix2i+1 =b2i
f2i+1x2i +g2i+1x2i+1+h2i+1x2i+2 = b2i+1
2i-1方程乘上某个数消去 2i方程中的 f2ix2i-1项,2i+1方程乘上某个数消去 2i方程中的 h2ix2i+1项,使 2i方程变为
α ix2i-2+β ix2i+γ ix2i+2=η i i=1,2,…,n/2
国家高性能计算中心（合肥） 142009-7-24
三对角方程组的直接求解法
② 重复 ① 最终可得,
case 1,case 2,
g1x1+h1x2 =b1,
f2x1+g2x2+h2x3 =b2,
f3x2 +g3x3+h3x4=b3,
f4x3 +g4x4 =b4,
可以分别得到
g1x1+h1x2 =b1 或 g1x1+h1x2 =b1
f2x1+g2x2 =b2 f2x1+g2x2+h2x3 =b2
f3x2+g3x3 =b3
解得 x1,x2 或 x1,x2,x3
③ 回代求解 x
并行化分析,①,② 消去奇下标可以并行化；
③ 回代求解可以并行化第十章 线性方程组的求解
10.1 三角形方程组的求解
10.2 三对角方程组的求解
10.3 稠密线性方程组的求解
10.4 稀疏线性方程组的求解
10.3 稠密线性方程组的求解
10.3.1 有回代的高斯消去法
10.3.2 无回代的高斯 -约旦法
10.3.3 迭代求解的高斯 -赛德尔法国家高性能计算中心（合肥） 172009-7-24
有回代的高斯消去法
算法基本原理
求解过程分为消元和回代两个阶段，消元是将系数矩阵
A化为上三角阵 T，然后对 TX=c进行回代求解。
消元过程中可以应用选主元方法，增加算法的数值稳定性。
下面是消元过程图：
i
i
主 元 行图 1 0,1
待 改 变 的 元 素不 变 化 的 元 素待 消 为 0 的 元 素已 消 为 0 的 元 素国家高性能计算中心（合肥） 182009-7-24
有回代的高斯消去法
并行化分析
消元和回代均可以并行化；
选主元也可以并行化；
消元过程的并行化图示：处理器按行划分
i
i
主 元 行图 1 0,1
待 改 变 的 元 素不 变 化 的 元 素待 消 为 0 的 元 素已 消 为 0 的 元 素可 并 行 化 部 分
10.3 稠密线性方程组的求解
10.3.1 有回代的高斯消去法
10.3.2 无回代的高斯 -约旦法
10.3.3 迭代求解的高斯 -赛德尔法国家高性能计算中心（合肥） 202009-7-24
无回代的高斯 -约旦法
串行算法原理
① 消元,通过初等行变换,将 (A,b)化为主对角线矩阵,
(方便起见,记 b为 A的第 n+1列 )
② 求解,xj=a’j,n+1/a’jj
1,
1,222
1,111
1,
1,2
1,1
21
22221
11211
nnnn
n
n
nn
n
n
nnnn
n
n
aa
aa
aa
a
a
a
aaa
aaa
aaa
国家高性能计算中心（合肥） 212009-7-24
无回代的高斯 -约旦法
SIMD-CREW上的并行算法
(1)处理器,n2+n个处理器,这些处理器排成 n× (n+1)的矩阵,
处理器编号为 Pik,i=1~n,k=1~n+1
(2)并行化分析
① 消元的并行化,// O(n)
for j=1 to n-1,each Pik Par-do //第 j次消元
Pij(i<>j),aij <— 0
Pik(i<>j,k=j+1~n+1),aik <— aik-ajk(aij/ajj)
end for
② 求解,for each Pjj(j=1~n) Par-do,xj <— aj,n+1/ajj //O(1)
(3)时间分析,t(n)=O(n),p(n)=O(n2),c(n)=O(n3) 成本最优？
国家高性能计算中心（合肥） 222009-7-24
无回代的高斯 -约旦法
成本最优？
串行算法的最优时间,由于 x=A-1b
① A-1b(假设已有 A-1),O(n2)
② 求 A-1,
∴ 求 A-1需要,2次 n/2× n/2矩阵的逆 i(n/2)
6次 n/2× n/2矩阵的乘 m(n/2)
2次 n/2× n/2矩阵的加 a(n/2)
i(n)=i(n/2)+6m(n/2)+2a(n/2)
a(n/2)=n2/2,m(n/2)=O((n/2)x) 2<x<2.5
=> i(n)=O(nx) 综上，串行算法的最优时间为 O(nx) 2<x<2.5
12
1
1121221
1121
1
1
12
1
111
12
1
1111
1
11212221
1211
0
0
0
0
00
00
11 AAAAB
IAA
I
B
A
I
AAI
A
I
AAI
B
A
IAA
I
AA
AA
A
其中，
10.3 稠密线性方程组的求解
10.3.1 有回代的高斯消去法
10.3.2 无回代的高斯 -约旦法
10.3.3 迭代求解的高斯 -赛德尔法国家高性能计算中心（合肥） 242009-7-24
迭代求解的高斯 -赛德尔法
串行算法原理如果对某个 k,给定的误差允许值 c有则认为迭代是收敛的。
并行化分析由于每次迭代需要使用本次迭代的前面部分值，因而难以得到同步的并行算法，下面给出一个异步的并行算法。
i
ij
k
jij
ij
k
jij
ii
k
i bxaxaax
11 1
cxx
n
i
k
i
k
i
1
1
国家高性能计算中心（合肥） 252009-7-24
迭代求解的高斯 -赛德尔法
MIMD异步并行算法
N个处理器 (N≤n) 生成 n个进程,每个进程计算 x的一个分量
算法
Begin
(1)oldi? xi0,newi?xi0
(2)生成进程 i
(3)进程 i
repeat
(i) oldi?newi
(ii)newi? (bi-∑ k<iaik× oldk-∑ k>iaik× oldk)/aii
until ∑ i=1~n| oldi - newi |<c
xi? newi
End
第十章 线性方程组的求解
10.1 三角形方程组的求解
10.2 三对角方程组的求解
10.3 稠密线性方程组的求解
10.4 稀疏线性方程组的求解
10.4 稀疏线性方程组的求解
10.4.1 线性方程组的并行化方法
10.4.2 稀疏线性方程组的迭代解法
10.4.3 高斯 -赛德尔迭代法的并行化国家高性能计算中心（合肥） 282009-7-24
线性方程方程的并行化方法
线性方程组选择算法的考虑因素
系数矩阵 A的结构
dense Gaussian elimination,etc
Sparse iterative method
triangular substitution,odd-even reduction
certain PDEs multigrid,etc
计算精度要求
Gaussian elimination,more accurate,more expensive
Conjugate gradients,less accurate,less expensive
计算环境要求
architecture,available languages,compiler quality
libraries?
国家高性能计算中心（合肥） 292009-7-24
线性方程方程的并行化方法
求解方法的并行化
(1)直接解法的并行化 (用于稠密线性方程组 )
- Gauss消去法 (包括选主元的 Gauss消去法 )
- Gauss-Jordan消去法
- LU分解法
(2)迭代法的并行化 (用于稠密和稀疏线性方程组 )
- Jacobi
- Gauss-Seidel(可异步并行化 )
- Jacobi OverRelaxation(JOR)
- Gauss-Seidel OverRelaxation(SOR)
- Conjugate Gradient
10.4 稀疏线性方程组的求解
10.4.1 线性方程组的并行化方法
10.4.2 稀疏线性方程组的迭代解法
10.4.3 高斯 -赛德尔迭代法的并行化国家高性能计算中心（合肥） 312009-7-24
稀疏线性方程方程的迭代解法
迭代解法
G r a d i e n tC o n j u g a t e
axaxabxxSOR
axabxxJ O R
axaxabxS e i d e lG a u s s
axabxJ a c o b i
ii
ij
k
jij
ij
k
jiji
k
i
k
i
ii
ij
k
jiji
k
i
k
i
ii
ij
k
jij
ij
k
jiji
k
i
ii
ij
k
jiji
k
i
)5(
/)()1(:)4(
/)()1(:)3(
/)(:)2(
/)(:)1(
11
1
11
1
10.4 稀疏线性方程组的求解
10.4.1 线性方程组的并行化方法
10.4.2 稀疏线性方程组的迭代解法
10.4.3 高斯 -赛德尔迭代法的并行化国家高性能计算中心（合肥） 332009-7-24
高斯 -赛德尔迭代法的并行化
由 PDE离散产生的稀疏线性方程组
(1)Laplace方程
),(),(
]1,0[]1,0[),(),(),(),( 2
2
2
2
yxyxu
Dyxyxg
y
yxu
x
yxu
且满足边值满足上面的方程，求
u =?
0 1 2 3
4 5 6 7
8 9 1 0 1 1
1 2 1 3 1 4 1 5
u = φ
h 步 长 ＝ 1 / N
国家高性能计算中心（合肥） 342009-7-24
高斯 -赛德尔迭代法的并行化
(2)由五点格式的离散化 (假设 g(x,y)=0)
为单位阵方程代入用二阶中心差分逼近：
，记
ER
RE
ERE
ERE
ER
A
NjiuuuuuL a p l a c e
hyxuyxuhyxu
hy
u
yhxuyxuyhxu
hx
u
Njiuu
jijijijiij
N
i
N
i
ij
,
41
141
141
14
,
,0][
4
1
)],(),(2),([
1
)],(),(2),([
1
,0),(
1,1,,1,1
2
2
2
2
2
2
国家高性能计算中心（合肥） 352009-7-24
高斯 -赛德尔迭代法的并行化
A为稀疏的块三对角矩阵国家高性能计算中心（合肥） 362009-7-24
高斯 -赛德尔迭代法的并行化
Gauss-Seidel迭代解法的并行化
(1)两种串行算法的计算顺序及其并行化顺序 1 顺序 2
注,顺序 1难以并行化；顺序 2可以小规模并行化
0 1 2 3
4 5 6 7
8 9 1 0 1 1
1 2 1 3 1 4 1 5
0 1 2 3
4 5 6 7
8 9 1 0 1 1
1 2 1 3 1 4 1 5
国家高性能计算中心（合肥） 372009-7-24
高斯 -赛德尔迭代法的并行化
(2) 红黑着色并行算法
① 并行计算所有的红点
② 并行计算所有的黑点
③ 重复 ①,② 直至满足精度要求
0 1 2 3
4 5 6 7
8 9 1 0 1 1
1 2 1 3 1 4 1 5
