广义Hermitian对称特征问题标准化转换的有效并行块算法

关键词：特征问题；块算法；Cholesky 分解；Hermitian 矩阵；并行计算；可扩展性
Efficient parallel blocked algorithms for generalized Hermitian/symmetric eigenproblem
Zhao Yonghua Chi Xuebin Cheng Qiang Zhao Tao
2 广义特征问题标准化的串行块算法
2.1 广义特征问题标准化的一般求解方法
广义特征问题标准化转换主要涉及正定矩阵 B 的 Cholesky 分解 B=LLH 和 C=L-1AL-H 的求解， 这里 L 为下三角矩阵。而求矩阵 C 的有效方法是对 C=L-1A (L-1) H 的等价形式 LCLH=A （1） 进行处理。为此我们将矩阵 A、B 和 L 划分为下面的 块形式
* A11 A21
* B11 B21
* * A22 * * B22
enddo end
3 广义特征问题标准化的传统并行块 算法
在传统的广义特征问题标准化并行算法中，通 常采用（1）首先并行计算 B 的 cholesky 分解， （2） 然后并行求解矩阵 C。 下面为 ScaLAPACK 软件包中 广义特征问题标准化例程 pdsygst 对应的部分并行化 代码，其中矩阵 A，B 是作为下三角矩阵处理的。 * Compute inv(L-1)*sub(A)*inv(L-T) K=1 NB = DESCA( MB_ ) KB = MIN( ICEIL( IA, NB )*NB, IA+N-1 ) - IA + 1 20 CONTINUE * Update the lower triangle of A(ia+k-1:ia+n-1, * ja+k-1:ja+n-1) CALL PDSYGS2( IBTYPE, UPLO, KB, A, IA+K-1, $ JA+K-1, DESCA, B, IB+K-1, JB+K-1, $ DESCB, INFO ) IF( K+KB.LE.N ) THEN CALL PDTRSM( 'Right', UPLO, 'Transpose', $ 'Non-unit', N-K-KB+1, KB, ONE, B, IB+K-1, $ B+K-1,DESCB,A,IA+K+KB-1, JA+K-1, DESCA ) CALL PDSYMM( 'Right', UPLO, N-K-KB+1, KB, $ -HALF,A,IA+K-1,JA+K-1,DESCA,B, $ IB+K+KB-1,JB+K-1, DESCB, ONE, $ A, IA+K+KB-1, JA+K-1, DESCA ) CALL PDSYR2K( UPLO, 'No transpose', $ N-K-KB+1, KB, -ONE, A, IA+K+KB-1, JA+K-1, $ DESCA, B, IB+K+KB-1,JB+K-1, DESCB, $ ONE, A, IA+K+KB-1, JA+K+KB-1, DESCA ) CALL PDSYMM( 'Right', UPLO, N-K-KB+1, KB, $ -HALF, A, IA+K-1, JA+K-1, DESCA, B, $ IB+K+KB-1,JB+K-1, DESCB, ONE, A, $ IA+K+KB-1, JA+K-1, DESCA ) CALL PDTRSM( 'Left', UPLO, 'No transpose', $ 'Non-unit',N-K-KB+1, KB, ONE, B, $ IB+K+KB-1, JB+K+KB-1,DESCB, A, $ IA+K+KB-1, JA+K-1, DESCA ) END IF K = K + KB
ATL ABL * BTL , B= BBL A BR * 。 B BR
其中 ATL 和 BTL 是 0×0 矩阵。块大小为 nb do until ABR 是 0×0 矩阵 确定块大小: b=min(nb, ABR 阶) 重新划分矩阵
m L= m L11 n - m L21 m m A11 A= n - m A21
Keywords: eigenproblem; blocking algorithm; Cholesky decomposition; Hermitian matrix; parallel computing; scalability [s1]1
引言
Ax=λBX
广义 Hermitian（对称）特征问题
影响， 并行算法在求解大规模矩阵的特征问题 时存在一定性能和扩展性问题[7]，从而影响到 并行求解特征问题时的整体性能。
在这些转换中，标准化转换是一个重要的转换。
基金项目：国家自然科学基金（60373060） ；国家 973 项目 （2005CB321702）
该转换主要涉及正定矩阵B的Cholesky分解B=LLH和 C=L-1AL-H 的求 解 ，这 里 L 为 下 三角矩 阵 。在 包括 ScaLAPACK[1]、PeIGS和HJS[3]等在内的并行软件包 中，所采用的并行计算方法是首先进行B的Cholesky 分解，然后求解矩阵C。 在并行计算过程中，由于矩阵被分布到不同的 处理器上，这必然引起大量的通信开销以及同一问 题在多个处理器上的并行处理。而采用上面算法思 想设计的并行算法增加了许多不必要的消息传递， 并使并行求解过程中失去许多并行的机会。 本文给出了一个新的广义特征问题标准化转换 的并行块算法。该算法将 Cholesky 分解和矩阵 C 的 求解结合在一起，避免了不必要的通信开销。另外， 并行算法通过将存在于各迭代步的具有多个右端项 的三角矩阵方程的并行求解延迟，最后通过一次并 行求解完成，增加了算法的并行性，并改善了并行 算法的负载平衡。并行算法采用了矩阵处理器的二 维块-循环方式和算法块并行策略。本文同时给出了 一个有效地并行求解具有多个右端项的三角矩阵方 程的并行块算法。
We present a new parallelization which combines the Cholesky into the transformation from generalized to standard form. By reducing the communication cost and extending the parallelism, the new algorithm can obviously improve the performance and scalability of the original algorithm. Moreover, an efficient parallel algorithm is proposed to compute a triangular AX=B with multiple right hand sides. From the tests using our parallel software PSEPS, the speed of the parallel algorithm is about two times of that of the classical parallel algorithms, and it has better performance and scalability than the classical parallel algorithms.
（Supercomputing Center, Computer Network Information Center, Chinese Academy of Sciences, Beijing 100080） Abstract: To solve generalized eigensystem problem, we need to reduce the generalized eigenproblem to an ordinary eigenproblem.
n-m
m
n-m
H 0 m C11 C21 C = L 22 n - m C21 C 22
ATL A= ABL BTL B= BBL
∗
n-m
H A21 A 22
(3)
A00 = A10 ABR A 20 B00 * = B10 B BR B20
广义 Hermitian/对称特征问题标准化转换的有效并行块算法
赵永华 迟学斌 程强 赵涛
（中科院网络信息中心超级计算中心，北京 100080）
摘要：广义 Hermitian（对称）特征问题并行求解的性能依赖特征求解器所选择的并行算法和矩阵的分布策略等诸多方面。在
本文，一个新的标准化转化的并行算法被给出，该并行算法将 Cholesky 分解结合到广义特征问题标准化转换中, 从而降低了 已有并行算法的通信开销，增加了算法的并行性。新算法可显著改善已有并行算法的性能和可扩展性。另外本文给出了一个有 效求解具有多个右端项的三角矩阵方程 AX=B 的并行块算法。通过自主开发的特征问题并行软件包 PSEPS 的测试结果表明，文 中的并行算法比传统的并行算法快大约 1 倍，并具有较好的可扩展性。
1 L21C11 2 (6) B22 = B22 − L21LT 21 (5) A21 := A21 T − A21LT (7) A22 := A22 - L21 A21 21
由关系式（3）和矩阵 C 的对称性，我们可以得到下 列关系式：
(8) A21 = A21 -
Z=
1 1 −H L21C11 + L22C21 = A21L11 − L21C11 2 2
先得到 B 的 Cholesky 分解，然后求解矩阵 C。为了 更好的描述本文给出的并行算法，算法 1 给出了将 正定矩阵 B 的 Cholesky 分解结合矩阵 C 的计算过程 中的算法。在算法 1 中，带♦的部分对应着矩阵 B 的 Cholesky 分解。算法中矩阵 B 的下三角元素在求解 过程中被矩阵 C 元素覆盖。 算法 1 包含 Cholesky 分解的串行块算法，块大小 为 nb begin 划分矩阵： A =
-T ) B21:=L21=A21L11 （求得下三角 L21）
* A11 = A Z 21
* (3) A22
−T (4) A21 : =A21L11
ห้องสมุดไป่ตู้
其中
H H H H Z= L21C11LH 21 + L22C21L21 + L21C21 L22 + L22C22 L22 。
的求解存在于许多科学领域中，并在电子结构、
有限元模型和量子化学 [s2] 等中具有重要作 用。 因为具有重要意义和广泛的应用, 矩阵特 征值问题成为高性能计算机的主要计算任务 之一。 广义特征问题的并行求解涉及到 5 个转 换[8]，许多人对它们的并行算法深入的研究［2,
3, 4, 5, 6］
。但受到消息传递和负载不平衡问题的
(5) 对下列的矩阵形式继续上述过程
ATL A= ABL
(6)
∗
ABR
2.2 结合 Cholesky 分解的串行块算法
广义 Hermitian 特征问题标准化的通常方法是首
BTL B= BBL
A00 = A10 A 20 B * 00 = B10 BBR B 20
1 L21C11 2
1 (9) A21 := C21 = L− 22 A21 （求得矩阵 C21）
(4)
H H L22C22 LH 22 = A22 − L21Z − ZL21
1 1 C21 = L− 22 ( Z − L21C11 ) 2
这样由(4)式的右端项可以求得矩阵 Z，进而通过(6) 式求得 C21。而(5)式具有（1）的形式，因而上面的 过程可以继续下去。
该并行算法已用于我们开发的 PSEPS（Parallel Symmetric Eigenproblem Solver）软件包中。基于该 并行软件包的测试结果表明，本文给出的并行算法 比传统的并行算法在性能上高出大约 2 倍左右，并 具有较好的可扩展性。 在本文中我们也给出了 PSEPS 软件包的整体性能测试结果。
* A11 A21 * B11 B21
* * A22 * * B22
这样由（2）式矩阵的分解关系，我们可容易的得到 等式
H L11C11L11 LCLH = H （L21C11 + L22C21 ）L11
♦(1) B11:=L11=Chol(B11)（由 B11 得到 L11） (2)A11:=C11=sygst(A11, L11)(由 A11 和 L11 形成矩阵 C11) ♦