华东师范大学：《高等数值分析（高性能计算/并行计算）》课程教学资源（讲义）09 线性方程组并行直接法（基于 MPI）

《并行计算》 线性方程组并行直接法 (基于MPI) 选主元LU分解 三角矩阵求解

http://math.ecnu.edu.cn/~jypan 1 线性方程组并行直接法 （基于 MPI） 《并行计算》 —— 选主元 LU 分解 —— 三角矩阵求解

线性方程组求解 Linear algebra-in particular,the solution of linear systems of equations-lies at the heart of most calculations in scientific computing. 一Dongarra&Eijkhout J.J.Dongarra and V.Eijkhout,Numerical linear algebra algorithms and software,/CAM,123(2000),489-514 口线性方程组是许多重要问题的核心，有效地求解线性方程组在科学与工程计算中非常重要 口并行计算机的问世，使求解问题的速度和规模大幅提高，同时也使计算方法产生了变化 口软件包：Lapack,ScaLapack http://math.ecnu.edu.cn/~jypan

http://math.ecnu.edu.cn/~jypan 线性方程组求解 Linear algebra — in particular, the solution of linear systems of equations — lies at the heart of most calculations in scientific computing. — Dongarra & Eijkhout J. J. Dongarra and V. Eijkhout, Numerical linear algebra algorithms and software, JCAM, 123 (2000), 489–514  线性方程组是许多重要问题的核心，有效地求解线性方程组在科学与工程计算中非常重要  并行计算机的问世，使求解问题的速度和规模大幅提高，同时也使计算方法产生了变化  软件包：Lapack，ScaLapack

华东师范大学数学科学学院 目录页 School of Mathematical Sciences,ECNU Contents 案 LU分解 PA=LU http://math.ecnu.edu.cn/~jypan

http://math.ecnu.edu.cn/~jypan 目录页 Contents 华东师范大学 数学科学学院 School of Mathematical Sciences, ECNU LU 分解 𝑷𝑷𝑷𝑷 = 𝑳𝑳𝑳𝑳

数据存储方式 ▣1 U分解的主要计算量是更新矩阵A。 ▣ 根据算法计算过程可知，如果是按列（或按行）连续分块存储在各个结点，则会出 现越往后计算越多结点空闲的情况，因此建议采用卷帘方式存储。 u12u13u14u15u16u17u18 u12u13u14u15u16u17u18 u1zu13u14u15u16u17u18 u12u13u14u15u16u17u18 21 22 u23u24u25426u27u28 l22u23u24u2su26u27u28 31l32 131132 u34u35u36u37u38 l332 u34u35u36u3u38 4 4l42 l41l42l43 L4l4243 151152 L51l52l53 151152 153 l61l62 l61l62l63 161162163 l71l2 1l22l23 1l23 L81182 181 182 183 la1 l82 l83 http://math.ecnu.edu.cn/~jypan

http://math.ecnu.edu.cn/~jypan 数据存储方式  LU 分解的主要计算量是更新矩阵 A。  根据算法计算过程可知，如果是按列（或按行）连续分块存储在各个结点，则会出 现越往后计算越多结点空闲的情况，因此建议采用卷帘方式存储。 ...

数据存储方式 口可以采用一维划分，也可以采用二维划分。在实际应用中，通常采用二维划分，即 在两个方向上都进行循环划分，然后存储到按二维排列的结点上。 A00 A01 A02 A A A2 A10 Au A12 A20 A21 A22 为简单起见，这里介绍一维列循环划分的并行算法 http://math.ecnu.edu.cn/~jypan

http://math.ecnu.edu.cn/~jypan 数据存储方式  可以采用一维划分，也可以采用二维划分。在实际应用中，通常采用二维划分，即 在两个方向上都进行循环划分，然后存储到按二维排列的结点上。 ► 为简单起见，这里介绍一维列循环划分的并行算法 A00 A01 A02 A10 A11 A12 A20 A21 A22 A0 A1 A2

串行算法 算法1.8部分选主元LU分解 1:p=1:n %用于记录置换矩阵 2:for k 1 to n -1 do 3: [amax,)=max<i≤nak%选列主元，其中I表示主元所在的行 4 ifl≠k then 5 forj=1 to n do 6: tmp=akj,ak=a,aj=tmp%交换第k行与第1行 7: end for 8: tmp=p(k),p()=p(),p(l)=tmp%更新置换矩阵 9: end if 10: for i=k 1 to n do 11: aikaik/akk %计算L的第k列 12: end for 13: for i=k+1 to n do 14: for i=k+1 to n do 15 a)=a-ak*aj%更新A(k+1:n,k+1:n) 16: end for 17: end for 18:end for http://math.ecnu.edu.cn/~jypan

http://math.ecnu.edu.cn/~jypan 串行算法

并行算法 计算过程 每次选主元时，只涉及到同一列，因此没有数据通信 ·确定主元后，计算L的第k列 ·将1（主元所在的行）和L的第k列传给其他结点 ·各处理器（进程）了 更新自己的数据 http://math.ecnu.edu.cn/~jypan

http://math.ecnu.edu.cn/~jypan 并行算法 计算过程 ►每次选主元时，只涉及到同一列，因此没有数据通信 ►确定主元后，计算 L 的第 k 列 ►将 l（主元所在的行）和 L 的第 k 列传给其他结点 ►各处理器（进程）更新自己的数据

并行算法 1:icol =0 2: for j=0 to n-2 do 3: if myid=j mod p then 4: find l laLicoll maxflai,icoll,i=j,j+1,...,n-1} if al.icol=0 then kill all processes A is singular and return 6 if Ij then swap aj,icol and aiicol 作 aiicol diicol/ajricol,i=j+1,...,n-1 fi-j-1=ai,icol;i=j+1,...:n-1 9 send(l,myid+1 mod p)and send(f,myid+1 mod p) 10: icol icol+1 11: else 12 recv(l,myid-1modp)and recv(f,myid-1modp)%广播1和f 13: if myid+1 mod pthen send(l,myid+1 mod p)and send(f,myid +1 mod p) 14 end if 15: if l j then swap Aj.:and AL.: 16: for k icol to m-1 do 17: aik aik fi-j-lajk,i=j+1,...n-1 18: end for 19:end for http://math.ecnu.edu.cn/~jypan

http://math.ecnu.edu.cn/~jypan 并行算法

华东师范大学数学科学学院 目录页 School of Mathematical Sciences,ECNU Contents 三角方程并行求解 一以Lx=b为例 http://math.ecnu.edu.cn/~jypan

http://math.ecnu.edu.cn/~jypan 目录页 Contents 华东师范大学 数学科学学院 School of Mathematical Sciences, ECNU 三角方程并行求解 —— 以 𝑳𝑳𝑳𝑳 = 𝒃𝒃 为例

三角方程组 Lo … n-1,0 列扫描算法 x=b/0 x=(6-lox)/. x2=(色，-lx-lx)/2 七3=(b-l0七。-lg1-l252)/3 x=(b-LXo-IxIn http://math.ecnu.edu.cn/~jypan

http://math.ecnu.edu.cn/~jypan 三角方程组 00 0 0 `0 11 1 1 n n 1,0 1,1 n n 1, 1 n n 1 1 l x b l l x b ll l − − − − x b − −              =                 列扫描算法 ( ) ( ) ( ) 0 0 00 1 1 10 0 11 2 2 20 0 21 1 22 3 3 30 0 31 1 32 2 33 x bl x b lx l x b lx lx l x b lx lx lx l = = − =− − =− − − x b lx lx l x l i i i i i i i ii = − − −− ( 00 11 ,1 1  − − )  