（计算数学专业论文）近似逆预条件子的研究.pdf

摘要 摘要 求解稀疏线性方程组是科学计算里的一个重要的课题。随着并行和分布式处 理器的出现与流行，使得寻求适合高性能计算机的可并行化预条件子变得越来越 重要。稀疏近似逆方法( s a i ) 因其良好的并行性受到特别的重视。 本文首先分类介绍了s a i 方法的主要理论和成功的算法。并结合数值实验， 对几种较成功的稀疏近似逆方法的优缺点( 强壮性，有效性，并行性等) 及其适 用的范围进行了分析，并与传统串行预条件子i l u t 进行了比较。得到结论如下： 多数情况下，s a i 预条件子没有i l u t 预条件子强壮、有效，但具有更好的并行 性。每一种具体的s a i 方法都有其优点和缺陷，换言之，在所有可比较的规则中， 没有哪种是绝对最优的。 本文将稀疏近似逆( s p a i ) 技术与多层块i l u 预条件子( b i l u m ) 结合起来， 得到一种新的预条件子，b i l u m s p a i ( 多层块不完全l u 稀疏近似逆) 。这种预 条件子保持了b i l u m 的强壮性，但比b i l u m 有更好的稀疏性和并行性。同时， 它比s a i 方法更强壮，更稳定。数值实验显示，b i l u m s p a i 是一种强壮有效的， 适用于求解一般稀疏线性方程组的可并行化预条件子。 本文还应用p h g u i l l a u m ee ta 1 给出的块常数矩阵模式，提出了一个并行性的 求稀疏矩阵近似逆的方法一块常数近似逆( b c a i ) 。该法与s a i 方法的区别在于用 低阶矩阵块( 块常数矩阵) 而不是稀疏矩阵，来近似矩阵的逆。块常数矩阵不要 求较多的存储，在有些情况下，甚至更少；而且无论是预条件的构造过程或是矩 阵与向量相乘都可以并行实现。文中给出了它的并行执行方法。通过数值实验可 以看出，b c a i 方法对于对称正定阵有着较好的效果，但是对于非对称j 下定阵可 能会失败，因为无法得到合适的块常数模式来近似矩阵的逆。 关键词：稀疏近似逆，预条件子，多层块不完全l u 分解，块常数矩阵，并行性 a b s t r a ( 了r a b s t r a c t s o l v i n gs p a r s el i n e a rs y s t e m si sa l li m p o r t a n ti s s u ei ns c i e n t i f i cc o m p u t i n g w i t h t h ea d v e n ta n dg r a d u a lp o p u l a r i t yo fp a r a l l e la n dd i s t r i b u t e da r c h i t e c t u r e s ，t h es e a r c h f o rp a r a l l e l i z a b l ep r e c o n d i t i o n e r sb e c o m e sa ni m p o r t a n ti s s u e o w i n gt ot h eg o o d p a r a l l e l i z a b i l i t y , s p a r s ea p p r o x i m a t ei n v e r s et e c h n i q u e sa r eb e c o m i n gp o p u l a r a tf i r s t ，t h ep a p e ri n t r o d u c e st h em a i nt h e o r yo fs a ia n ds o m e s p e c i f i ca l g o r i t h m s t h a tp r o v e ds u c c e s s f u l a c c o r d i n gt ot h en u m e r i c a lt e s t s ，w ec o m p a r ea n da n a l y z e t h e s ea l g o r i t h m sf r o ma s p e c t s ( r o b u s t n e s s ，e f f i c i e n c y , p a r a l l e l i s m ) a n dw ed r a wt h e c o n c l u s i o nt h 她i ng e n e r a l ，t h ee f f i c i e n c ya n dr o b u s t n e s so fs a ii sn o ta sg o o da si l u c l a s so fp r e c o n d i t i o n e r s ，b u ti th a sb e t t e rp a r a l l e l i z a b i l i t y t h en u m e r i c a lt e s t ss h o w e d t h a te a c hs a ia l g o r i t h mh a si t so w na d v a n t a g e sa n dd i s a d v a n t a g e s i no t h e rw o r d s ， n o n ec a l lb et h eb e s ta l g o r i t h mf o ra l lt h ea s p e c t ss h o u l db ec o n s i d e r e d b e s i d e s ，t h ep a p e ri n v e s t i g a t e st h eu s eo fs p a r s ea p p r o x i m a t ei n v e r s et e c h n i q u e s i nam u l t i l e v e lb l o c ki l up r e c o n d i t i o n e rt od e s i g nar o b u s ta n de f f i c i e n tp a r a l l e l i z a b l e p r e c o n d i t i o n e r f o rs o l v i n gg e n e r a l s p a r s em a t r i c e s t h ep r e c o n d i t i o n e r r e t a i n s r o b u s t n e s so ft h eb i l u ma n dh a sb e t t e ra b i l i t yo fc o n t r o ls p a r s i t ya n di n c r e a s e d p a r a l l e l i s m n u m e r i c a le x p e r i m e n t ss h o wt h ee f f e c t i v e n e s s a n de f f i c i e n c yo ft h e p r o p o s e dv a r i a n to fb i l u m t h ep a p e ra l s oe s t a b l i s h e sap a r a l l e l i z a b l ea p p r o x i m a t ei n v e r s ep r e c o n d i t i o n e r u s i n gt h eb l o c kc o n s t a n tp a t t e r ng i v e nb yp h g u i l l a u m ee ta 1 1 1 1 em a i ni d e ao ft h i s p r e c o n d i t i o nt e c h n i q u ec o n s i s t e di na p p r o x i m a t i n gt h ei n v e r s eo fam a t r i xb yab l o c k c o n s t a n tm a t r i xi n s t e a do fas p a r s em a t r i xl i k ei ns a im e t h o d s t h e yd on o tr e q u i r e m o r es t o r a g e ，o re v e nl e s s ，a n da l ew e l la d a p t e dt op a r a l l e l e dc o m p u t i n g , b o 付lf o rt h e c o n s t r u c t i o no ft h ep r e c o n d i t i o n e ra n df o rm a t r i x v e c t o rp r o d u c t s t h ep a r a l l e l i z a b l e i m p l e m e n t so fb c a ia r eg i v e ni nt h ep a p e r a c c o r d i n gt o t h en u m e r i c a lt e s t s ，a c o n c l u s i o nc o u l db er e a c h e dt h a tb c a ii ss u i t a b l ef o rs p d s ，a n db e h a v e sb a d l yf o r n o n - s p dm a t r i c e s k e y w o r d s ：s a i ，p r e c o n d i t i o n e r , b i l u m ，b l o c kc o n s t a n tm a t r i x ，p a r a l l e l i s m i i 独创性声明 本人声明所呈交的学位论文是本人在导师指导下进行的研究工 作及取得的研究成果。据我所知，除了文中特别加以标注和致谢的地 方外，论文中不包含其他人已经发表或撰写过的研究成果，也不包含 为获得电子科技大学或其它教育机构的学位或证书而使用过的材料。 与我一同工作的同志对本研究所做的任何贡献均已在论文中作了明 确的说明并表示谢意。 签名：垒i 邀 日期：7 年岁月2 7 日 关于论文使用授权的说明 本学位论文作者完全了解电子科技大学有关保留、使用学位论文 的规定，有权保留并向国家有关部门或机构送交论文的复印件和磁 盘，允许论文被查阅和借阅。本人授权电子科技大学可以将学位论文 的全部或部分内容编入有关数据库进行检索，可以采用影印、缩印或 扫描等复制手段保存、汇编学位论文。 ( 保密的学位论文在解密后应遵守此规定) 签名：垃导师签名： 日期：加7 年g - 月2 7 日 第一章引言 1 1 课题研究背景 第一章引言 求解稀疏线性方程组是科学计算里的一个重要的课题。大型的，稀疏的矩阵 在很多实际科学的领域里也有广泛的应用，例如模拟电路，生物科学，数据库检 索等。目前，一般采用预条件k r y l o v 子空间方法来求解大型稀疏线性代数方程 组，如【1 1 。但是需要用预条件子提高求解精度和速度。 对于线性系统a x = - b ，其中，么是 阶稀疏矩阵，可以应用左右预条件子求解 p a x = p b 或a p y = p b ，x = p y 。以下都仅考虑左预条件子的情况，右预条件子可以同 上参考。对于k r y l o v 子空间方法，迭代的速度与系数矩阵么的特征值的分布紧 密相关，特征值分布越集中在l 附近，收敛速度越快，否则，收敛速度慢，甚至 不收敛。所以，预条件子p 的选择标准之一是能使预条件后系统p a 有更紧凑的 特征值分布。 预条件子的作用在于提高求解线性系统的速度。因此，预条件子的构造和相 应的预条件过程希望时间耗费量越小越好。预条件过程的主要计算一般是矩阵与 向量的乘积，因此，预条件子p 与向量相乘应快速。 综上，预条件子的选取主要需考虑以下方面： p a 应该具有更聚集的特征值； p 应该能快速构造出来； 预条件子p 与向量相乘应快速。 目前，已经有多种方法可以得到预条件子尸，例如块j a c o b i ，多项式预条件 子，和最常用的i l u 技术。其中，i l u 技术以其强壮性和有效性得到了广泛的应 用，但是由于它内在的并行，使得它不可能并行处理。随着并行和分布式处理器 的应用和发展，构造可并行化预条件子越来越成为能否高效求解大型线性系统的 关键。 对可并行化的预条件子，近年来出现了稀疏近似逆技术，即选择稀疏矩阵m 作为预条件子，m a 一。稀疏近似逆预条件子m 可以很好的满足上面三个条件。 第一，膨近似于彳的逆，若是精确逆，则预条件系统m a 的特征值全部为1 ，当 然这是不肯能的，但是m a 的特征值会在1 附近聚集。第二，这个条件涉及到具 电子科技大学硕士学位论文 体的构造算法，但是稀疏近似逆预条件子m 具有良好的并行性，考虑到高性能并 行计算机的应用，肘在构造时间上会给出一个令人满意的结果。第三，由于m 也 是稀疏矩阵，m 与向量的乘积可以是比较快速的，而且m 是可并行的，即可以并 行处理m 与向量的乘积运算。综上，稀疏近似逆预条件子是非常有吸引力的和有 前途的预条件子。 近二十年，稀疏近似逆方法由于其良好的并行性得到了广泛的关注和快速的 发展。如今已经出现了多种构造方法。大量的数值实验也证明了这类预条件子可 提供较高的并行性且可以求解一些传统串行预条件子，如i l u t 1 0 难以求解的问 题。然而，s a i 方法发展尚未完善，对一般的稀疏矩阵，它们的健壮性并不强， 未达到可与传统串行预条件子( 如i l u t ) 相媲美的水平。 1 2 本文工作及创新 在这样的研究背景下，本文致力于近似逆方法的进一步研究，在保持稀疏近 似逆方法的并行性的前提下，希望能给出更强壮，更有效，更适用于一般稀疏矩 阵的近似逆预条件子。 本文的主要创新点包括如下方面： 结合数值实验，比较和分析现有的几种较成功的稀疏近似逆方法的优缺点( 强 壮性，稳定性，并行性等) 及其适用的范围； 将具有强壮性的块i l u t 方法和有良好并行性的s a i 预条件子( 其中的s p a i 算法) 结合起来，得到新的预条件子b i l u t - s p a i ，该预条件子在强壮性和并 行性中给出了一个良好的平衡； 应用p h g u i l l a u m ee ta 1 给出的块常数矩阵模式，建立一个新的求近似逆的方 法块常数稀疏近似逆( b c a d ，并给出它的并行执行方法。该预条件子虽然非 稀疏，但是保持了良好的并行性，并且在所需存储空间上更具有优势； 利用数值实验，定量的研究b c a i 预条件子。针对对称正定和非对称正定稀 疏矩阵两大类情况分别给出了与其它预条件子的比较和说明。 1 3 本文内容及组织 本文的主要内容是对稀疏近似逆方法进行了深入的研究，分析并比较了多种 较成功的稀疏近似逆预条件子。进一步的，将s a i 方法与前人给出的成功的算法 想结合，得到新的并行预条件子。并且结合数值实验来验证新方法的有效性，强 2 第一章引言 壮性和并行性。结合b i l u t 方法，给出了新的并行预条件子。 文章组织如下： 第一章引言 对课题有关背景知识做了简要的介绍 第二章稀疏近似逆方法介绍 介绍稀疏近似逆方法的主要理论和主要的成功的算法，并结合数值实 验，分析并比较了各种稀疏近似逆预条件子的效果。 第三章稀疏近似逆与多层块i l u 预条件方法 把s a i 方法( s p a i 算法) 和b i l u t 方法结合，给出了新的并行预条 件子，并给出数值实验。 第四章块常数近似逆方法 利用p h g u i l l a u m ee ta 1 给出的块常数矩阵模式，建立一个新的求近 似逆的方法b c a i ，并给出它的并行执行方法。并定量的研究b c a i 给出的近似逆的近似程度。 第五章结论和展望 对本文中所做的工作进行了总结，并提出了以后的发展方向。 电子科技大学硕士学位论文 第二章稀疏近似逆方法介绍 2 1 近似逆预条件子的发展背景 如何高效的求解大型线性系统a x = - b 是科学计算里一个重要的课题。目前， 一般采用迭代法，而其中k r y l o v 子空间方法是发展迅猛【l 】，应用广泛的一类。研 究证明【啦 3 】，对于k r y l o v 子空间方法，例如，c g 法，g m r e 法，有效的预条件 子的选取比迭代方法的选取更为重要。许多好的预条件予，如i l u 及s s o r ，以 其强壮性和好的收敛性被广泛应用。但是，随着向量计算机的发展，因算法的并 行性可以大大提高计算的效率，所以构造预条件子的算法的并行性变得很重要。 i l u 等算法固有的序列性很难通过改良实现并行执行，而近似逆预条件予以其天 然的并行性而受到重视，并在近二十年来得到迅速的发展。其基本思想如下： 应用左右预条件子求解 a x = - b p a x - = - p b 或么毋，叫砀，x - - - p y 选取p a ，直接计算出近似逆p 并应用于此线性系统，使得预条件的过程 只包含矩阵尸与向量的乘法运算。显然，近似逆预条件子具有良好的并行性。 除此之外，近似逆预条件子还有一个重要的优势，即也适用于非对称j 下定问 题，而在这种情况下，不完全分解方法效果常常很差【4 】，所以，即使不考虑并行 性，近似逆预条件子也有广泛的应用之处。 多项式求近似逆是一类重要的方法，它的主要思想是用低阶的系数矩阵彳的 多项式来近似彳。 近二十年，另一类近似逆方法被发展和重视起来，即稀疏近似逆方法。这种 方法最早于七十年代初被提出【5 一，这种方法主要思想在于用一个稀疏矩阵来近似 所求稀疏矩阵的逆，但是，这并不一定是显而易见的，因为稀疏矩阵的逆往往是 稠密的。事实上，不可约稀疏矩阵的逆已经被证明是稠密的【7 ，8 】。不过，一般情况 下，稠密的逆矩阵里会有很多的绝对值很小的元素，所以抓住重要的大的元素， 而令小的不重要的元素为0 ，即选择非零元模式成了一个关键的但很难的问题。 稀疏近似逆方法提出以来，正是因为缺乏自动选择非零元模式的算法，而没有引 起广泛的注意。随着许多自适应选择非零元模式的算法的出现，这类方法得到了 4 第二章稀疏近似逆方法的介绍 更多的重视。 现在存在着几种比较优秀的算法，在一些问题上可以与传统的不完全分解l u 方法相媲美，在一些情况下，可以优于i l u 方法。下面的内容致力于描述这些成 功的s a i 算法，并对一系列不同的稀疏矩阵给出数值实验，分析比较这些方法的 适用情况和效果。 下面将分类介绍稀疏近似逆方法的思想、发展过程和主要算法，和通过数值 实验的比较和分析。 2 2 稀疏近似逆( s a i ) 方法的介绍 稀疏近似逆方法目前主要分为二大类：基于f r o b e n i u s 范数最小化的方法，基 于矩阵分解的近似逆方法。每类都有一些经实验证明非常优秀的算法，例如基于 f r o b c n i u s 范数最小化的方法中，有s p a i 算法，m r 算法，s e l f - p r e c o n d i t i o n i n gm r 算法；第二类中有f s a i 算法，a i n v 算法，及i n v e r s ei l u 方法。下面的内容的 意图在于介绍这些成功的算法，并通过数值实验，深入研究和比较这些方法得到 的预条件子的效果，强壮性及并行性等。以便给出对具体问题应采用哪种算法最 优的指导。 以下都是针对线性大型系统a x = - b 进行讨论，其中，彳是以阶方阵，b 是，l 阶 列向量。 2 2 1 基于f r o b e n i u s 范数最小化的近似逆方法 这种方法是求近似逆中应用最广泛，并行性质最好的一类方法。它的主要思 想是： 求解最小值问题 m i n l l i - m a l l m，es i i ， 其中，s 是非零元模式。 这样，所得的m 就是彳的近似逆矩阵。显然， 孵卜彳m 炉嘉孵”圳， 其中，勺和唧分别表示，和m 的第歹列。这样，问题就变成了求解刀个独立 5 电子科技大学硕士学位论文 的阼阶线性方程的问题，具有良好的并行性。 如果稀疏模式s 在计算前确定了，而肘是一个满足模式s 的稀疏矩阵，贝, j j r n j 里有很多零元，且这些零元的位置是已知的，可以不计算。所以，将中的零元 素去掉得到砖，相应的，巳记为乞，a 中对应于的零元素的列也可以舍去， 记为五。于是，问题变成了求解”个独立的l e a s ts q u a r e 问题： a 叛j = 龟i 这显然是可以并行计算的。 这样，如果s 确定了，就可以直接计算膨。显然，选择s 成了至关重要的问 题。s 有多种静态的选择方法，例如直接采用彳或a r 的稀疏模式，采用彳”的稀 疏模式，m 是一个较小的整数。但是，这样近似的结果通常都是不理想的。这是 因为，稀疏近似逆方法建立在假设存在一个稀疏矩阵m ，可以较好的近似任意给 出的矩阵么的逆。但是，即使彳是一个稀疏矩阵，“- 1 通常都是稠密的，而且不 可约矩阵的逆已经被证明了是稠密的f 7 ，3 】。幸运的是，矩阵的逆中常常有大量元的 绝对值很小，如果能找到关键的元的位置，通过省略小的值就可以得到有效的稀 疏模式s 。t h o m a sh u c k l e 就提出了预测么- 1 的稀疏模式的方法【9 】。而且，有些预 测算法成功的应用于由积分方程数值解产生的稠密线性系统中，建立了有效的稀 疏近似逆预条件子【1 0 , 11 j 。 但是，对于一般的矩阵预测其逆的稀疏模式是很困难的。一个普遍的选择是 选择a 的稀疏模式作为么- 1 的稀疏模式，因为通过实验观察可以看出，在彳有非 零的位置，彳一的相同位置的元素值也比较大。但是，这个简单易行的方法并不 够强壮，对于一般的稀疏问题，彳一在么零元处也有很多相对大的值。还有一种 比较普遍的方法是选择的稀疏模式，k 2 ，k 是一个比较小的正整数。这个方 法在理论上有n e u m a n n 理论支撑。尽管的稀疏模式比直接选择a 的稀疏模式 更准确，但是仍然不够强壮。而且，相应的计算量，存储空间和预条件过程都随 着尼的增大而急速增加。 因为预测一般的矩阵逆的稀疏模式是很困难的，很多科学工作者致力于得到 自适应的动态的稀疏模式s 。这些方法通常猜测一个初始的s ，然后进行迭代，直 至满足停止条件，例如，i e ，一a m 川 0 ，或者是达到非零元个数的 l i ， oi 1 2 最大值。c o s g r o v ee la t 首先提出了这种方法【1 2 】。然后，g r o t e 和h u c k l e 1 3 】及 g o u l d 和s c o t t 1 4 】也提出了一些改进的方法。其中，g o u l d 和h u c k l e 1 3 】提出的 6 第二章稀疏近似逆方法的介绍 方法更为成功，后来被称为s p a i 预条件子。算法如下： 算法2 - 1 s p a i 算法 f 0 r 膨的每一列m ，： ( 1 ) 选择一个初始的稀疏模式zj = ii ( f ，) g ，g 是m 的稀疏模式。 ( 2 ) 决定行标集合，对应于非零元。求j = a ( z ，j ) 的q r 分解。然后计算l e a s t s q u a r e 问题二向= 乞的解砖，并计算其剩余向量，= 乞一a r h j 。 w h i l e 例2 g ： ( 3 ) 令l 等于z 的集合，r ( d 0 。 ( 4 ) 令歹等于所有包含于三但不包含于，的指标的集合。 ( 5 ) 对于每个k 歹，计算新剩余向量的范数露 霹= ：一( ，r 彳吼) 2 i i a p , i i ： 然后删除指标k ，除了使虞最小的m a x 个的k 。 ( 6 ) 决定新的指标集j ，并对a ( i u i ，j u 了) 计算q r 分解。然后解新的l e a s t s q u a r e s 问题，计算新的剩余向量，= 巳一彳，i = i u i ，j = j u - ，。 这个算法需要用户设置一个初始的稀疏模式，还有几个参数，包括最大残量 的阀值占及每次迭代的新的非零元的最大个数m a x ，同时m a x 也规定- j ( 3 ) 至1 j ( 6 ) 的最大迭代次数。s p a i 的具体实现中的问题和参见【删。 由s p a i 算法计算m 的代价相当大，c h o wa n ds a a d 1 5 1 提出了一个新的算法， 主要思想是利用几步最小残量法来近似求解刀个线性子方程。下面介绍【1 5 】中的基 本算法，被称为m r 算法： 算法2 - 2 m i n i m a lr e s i d u a l 算法 选择一个初始的m = m o = 【，z i ，m 2 ，】。 f o r 每列m ，j 2 1 2 “d o f o ri = 1 , 2 ，强，d o t j = e i a m j q = r f a r j ( ( a 5 ) r ( 以) ) m i2 m j 七q 产i 对m j 进行数值选元，即令较小的元变为零元。 7 电子科技大学硕士学位论文 e n dd o e n dd o 这里表示迭代次数。每次迭代，计算当前的剩余向量厂，并进行一维的残量 范数最小化。这里几个参数的选取详见【”】。 除这两个成功的算法之外，还有多种算法【”】。例如s e l f - p r e c o n d i t i o n i n g 。这些 方法的好处在于它们能给出比较有效的近似逆预条件子。s p a i 和m r 算法的缺点 是它们都有较多的参数且运算量大，而且，它们都不能保证m 是非奇异的，同样 地，即使彳是一个对称正定阵，m 也不一定能是对称正定的。这样，预处理过的 线性方程就不能用c g 法等迭代求解，只能用可以求解任意类型矩阵的迭代法来 解。而且，由于这些方法的自适应和不规则性，就很难并行执行。 虽然如此，这些方法也有很多值得注意的优点。例如，s p a i 和m r 预条件予 对于重排( r e c o r d e r i n g ) 不敏感。研究表明，不完全分解法产生的预条件子对重排很 敏感【1 6 ，1 7 1 。这样有一个好处就是，系数矩阵a 可以被随意分块，而不用为了得到 更好的收敛速度而考虑重排。还有一个优点是：数值实验证明了，s p a i 和m r 预 条件子可以解决一些其它预条件子，象i l u 1 8 】，解决不了的“难 问题。而且， 这类预条件子被证实更适用于有离散的内存电脑解决大型问题。由于这些好处， 这类方法仍需要得到更多的重视和发展。 2 2 2 基于矩阵分解的近似逆方法 这节将讨论基于矩阵分解的近似逆方法。这种方法的主要思想如下： 如果彳可以被分解为l d u , l 是一个单位下三角矩阵，d 是a 的对角阵，【厂 是单位上三角阵，则a 一= u - 1 d - 1 l - 1 = z d 叫形r ，其中，z = u 一，矽= l - 7 是单位上 三角阵。上面讨论过的，z 和形在上三角部分将是稠密的。所以，我们可以采用 稀疏近似逆方法来近似乞z ，矿形，则彳的分解近似逆是 膨= 历一矿7 a ，此处历d 计算一个非奇异阵的近似逆彳有几类方法。第一类方法不需要计算厶u 阵， 只需要有彳，预条件子就可直接得到。这类方法包括由k o l o t i l i n aa n dy e r e m i n 3 3 l 提出的f s a i 算法，k a p o r i n 【2 3 】提出的相关的方法，不完全双共轭法【2 5 , 2 4 1 等。第二 类方法与之对应的，要计算出a 的分解式，然后直接采用s a i 技术求其分解阵的 近似逆，得到预条件子。这类两阶段法也常常被单独分为一类，被称为不完全分 解因子近似逆方法( i n v e r s ei l ut e c h n i q u e s ) 。 第二章稀疏近似逆方法的介绍 下面大概介绍f s a i 法【2 6 1 。当a 是对称正定阵，一个下三角矩阵龟可有求解 下面的矩阵方程得到 印q b = 岛，( f ，j ) 其中，墨表示下三角阵的稀疏模式，包括对角元素。 0 ，可以一列- n 的求解，每列都要求解一个小的子正定对称线性系统，大小 是墨规定的每列非零元的个数。这里规定0 ：的对角元都大于零。定义 b ：( d i a g ( a 。) ) 一1 和q = 西7 2 龟，则预条件子为瓯彳q ，且为一个对角元都是1 的对称正定阵。 s ，的选取通常是直接取彳的下三角阵的稀疏模式，或者是取的下三角阵 的稀疏模式 2 7 1 ，k = 2 或3 。 f s a i 方法有良好的并行性，已经被成功应用于各种困难的有限元问题的求解 3 3 ，- 3 4 1 。它的最主要的缺点是需要预先定义稀疏模式s p 简单的采用上面提到的模 式往往对一些问题效果是不好的。对从有限元分析里出现的各种对称正定阵的近 似逆稀疏模式的选取， 2 8 ，2 9 给出了建议。 另一种方法被广泛的注意的方法是基于不完全共轭法的，首先由【1 3 】提出。 这类方法中最成功的是被称为a i n v 的算法【2 4 , 3 0 。a i n v 不要求提前规定稀疏模 式舅，以任何稀疏模式都可以应用。具体算法如下： 算法2 - 3 a i n v 算法 ( 1 ) 令z ：= q ；硝= q l ( 2 ) f o ri = 2 ，n 妒= q f o r 产1 扣l = 口；- i ) ， 乏力= z j - 0 - - ( 缔) e n d ， 砖卜d = 口_ 矿q e n d 基于矩阵分解的近似逆方法有很多优点。第一，它构造的预条件算子可以保 持矩阵彳的对称正定性，显然，这是一种好的值得保留的性质。而且，如果a 是 对称正定阵，则预处理过的线性系统可以用c g 法求解。第二，与其它近似逆方 法相比，它的计算量较小，且有较少的用户定义的参数，减小了不确定性。第三， 9 电子科技大学硕士学位论文 它对重排很敏感，也就意味着，可以通过适当的重排技术使得非零元填充的现象 减少，提高收敛效率【3 0 】。 另一方面，基于矩阵分解的近似逆方法也有它的缺点。首先，它的并行性不 如不考虑矩阵分解的方法。a i n v 法被证实难以并行实现，并难以在离散内存的 计算机上实现。而且，尽管对角变换【3 1 1 和对角补偿【3 2 】技术被用于保证近似逆的计 算，但是，却不能保证得到的预条件子是有效的，特别是需要进行大量的变换的 情况。 2 3 数值实验 这一章里，将给出相关的数值实验，来比较上面介绍的几种算法的强壮性， 有效性和适用范围。这里主要对第一类方法中的s p a i 算法，m g 算法，自适应的 m r 算法( m r p ) ，及第二类方法中的a i n v 算法，和不完全分解近似逆方法i i l u 方法进行了数值实验，实验矩阵来自各个方面得到的矩阵，分为对称正定和非对 称正定两大类。 虽然稀疏近似逆类的预条件子的最大优势是并行性，但是这罩的数值实验都 是用序列型计算机来运行的，虽然如此，仍然能够比较各种算法的优缺点和适用 矩阵类型。首先，进行对称正定矩阵的数值实验。 2 3 1 稀疏对称正定矩阵的数值实验 这里给出来自于各个领域的六个稀疏对称正定矩阵p o 。这些矩阵来源于不同 的领域和应用：偏微分方程有限元法，电力系统，电路分析，等等。这里用预条 件c g 法( p c g ) 来求解这类问题。表2 1 中，靠表示实验矩阵的阶数，行舷表示实 验矩阵中非零元数。 实验矩阵n o s * ，* b u s 从h a r w e l l b o e i n gc o l l e c t i o n 中选出【3 5 1 ，n a s a * 矩阵 从t i md 撕sc o l l e c t i o n 中选取【3 6 1 。这些问题中，n a s a * 的问题是最难的，虽然实 验矩阵并不多，但是都比较具有代表性。 对下列预条件子进行了数值实验：无填充的不完全c h o l e s k y 方法( i c ( o ) ) ， k o l o t i l i n a 和y e r c m i n 提出的矩阵分解稀疏近似逆算法f s a i ，及算法a i n v ，并给 出了分析和比较的结论。 i c ( 0 ) 和f s a i 预条件子没有参数，它们比较容易实现和使用，而a i n v 方法 需要设置参数，但是就更灵活，可以适应更难的情况。运用迭代法时，需要设置 1 0 第二章稀疏近似逆方法的介绍 最大迭代次数，这里设置的是m a x f f = i i l i n n ，1 0 0 0 ) 。下面的叙述里，如果预条件 子在小硎f 之内还不能收敛，就认为是迭代失败。迭代次数与最大残差值相关。 表2 - 1s p d 实验矩阵信息 实验矩阵 n，z ，王z 来源 n o s 3 9 6 0 8 4 0 2b i h a r m o n i ce q u a t i o n n o s 77 2 92 6 7 3d i f f u s i o ne q u a t i o n 6 8 5 b u s6 8 51 9 6 7p o w e rs y s t e mn e t w o r k 1 1 3 8 b u s1 1 3 82 5 9 6p o w e rs y s t e mn e t w o r k n a s a 2 1 4 62 1 4 63 7 1 9 8s t r u c t u r a la n a l y s i s n a s i a 2 9 1 02 9 1 08 8 6 0 3s t r u c t u r a la n a l y s i s i c ( 0 ) 是最强壮的方法，上面的六个实验矩阵都能够求解，并且收敛速度最快。 a i n v 算法在b c s s t k 2 3 失败，而且在其它问题上迭代收敛速度也稍慢于i c ( o ) ， 但快于f s a i 算法。 有效性定义为在最短的时间内得到在精度范围内的有效解，时间包括构造预 条件子的时间和迭代法计算的时间。i c ( o ) 仍是最有效的方法，因为预条件子的花 费远小于其它两种方法，总时间耗费比较小。而以f s a i 和a i n v 算法为代表的 s a i 方法由于预条件子构造花费巨大，所以在总时间上处于劣势。但是考虑到这 些数值实验是在序列性计算机上进行的，而s a i 方法具有良好的并行性，所以如 果在并行计算机上执行s a i 算法，很可能使其更加有效，并优于序列性的预条件 子。 如果只考察迭代部分的时间，即使在序列性计算机上执行，s a i 方法是很有 优势的。这是因为在迭代过程中，计算量主要由矩阵与向量相乘的运算组成，而 s a i 方法在矩阵向量乘积中特别高效，所以导致了较快的收敛速度。而在s a i 方 法中，a i n v 的速度优于f s a i 方法。但是如果在并行计算机上，由于f s a i 算法 是完全并行的，而且比较容易在分布式内存机器上执行，可能效果会有大的提高。 对于s p d 矩阵，给出了如下结论【3 0 】： i c ( o ) 是实验算法中最强壮的； 在非并行计算机上，i c ( o ) 也是最快的方法； f s a i 在s a i 方法中最强壮； 在一些问题上，a i n v 即使在非并行计算机上也是非常有效的。 因此，在并行计算机上，s a i 方法应该受到更多的重视，但是可能会在一些 电子科技大学硕士学位论文 比较难的问题上失败。s a i 方法在有限元问题和电力系统问题上表现良好，但在 体系分析问题上不理想。另外，f s a i 方法强壮性较好，而且并行性非常好，如果 能够采用更智能的稀疏模式选择方式会提高这种方法的可用度和效率。 2 3 2 稀疏非对称正定阵的数值实验 这里给出来自于各个领域的十个稀疏非对称正定矩阵【3 0 】。这些矩阵来源于不 同的领域和应用：石油开采模型，电路模型，半导体设备模拟，等等。这里用预 条件b i c g s t a b 来求解这类问题。表2 - 2 给出实验矩阵的名称，规模及来源。其 中，, 表示实验矩阵的阶数，, n z 表示实验矩阵中非零元数。 表2 - 2 非s p d 实验矩阵信息 实验矩阵 nn n z来源 o r s i r r l1 0 3 06 8 5 8o i lr e s e r v o i rs i m u l a t i o n o r s i r i t 28 8 65 9 7 0 o i lr e s e r v o i rs i m u l a t i o n s h e r m a n11 0 0 03 7 5 0o i lr e s e r v o i rs i m u l a t i o n s a y l r 31 0 0 03 7 5 0o i lr e s e r v o i rs i m u l a t i o n s h e rm a n 21 0 8 02 3 0 9 4o i lr e s e r v o i rs i m u l a t i o n p o r e s 35 3 23 4 7 4o i lr e s e r v o i rs i m u l a t i o n a d d 2 02 3 9 51 7 3 1 9c i r c u i tm o d e l i n g a d d 2 02 3 9 51 7 3 1 9c i r c u i tm o d e l i n g s 後n g l3 1 6 92 0 8 4 1s e m i c o n d u c t o rd e v i c es i m u l a t i o n 岍l1 8 5 61 1 3 6 0o i lr e s e r v o i rs i m u l a t i o n 所有的来源于石油开采模拟模型的实验矩阵都从h a r w e l l b o e i n gc o l l e c t i o n 得 到。来源与电路模拟和半导体设备模拟的实验矩阵从t i md a v i sc o l l e c t i o n 得到。 虽然矩阵大多来源于石油工业，但是这些问题包涵了不同的类型，无论从结构上 还是从数值上。 对以下算法进行了数值实验【蚓： 无填充的不完全l u 分解，i l u ( o ) ； g r o t e 和h u c k l e 的算法，s p a i ； c h o w 和s a a d 的算法，m r ； k o l o t i l i n a 和y e r e m i n 的算法，f s 触； 不完全双共轭法，a i n v 。 1 2 第二章稀疏近似逆方法的介绍 预条件后采用b i c g s t a b 求解，将所有预条件子的最大迭代次数设定为 m a x t = 5 0 0 ，当不采用预条件时，令m a x i t = n 。从计算结果上看，i l u ( 0 ) 仍是最健 壮的预条件子，可以求解所有问题。稀疏近似逆方法不够强壮：s p a i 在1 个实验 矩阵上失败了，m r 在2 个上失败，f s a i 在3 个上失败，a i n v 在5 个上失败。 所有s a i 方法的预条件子都不能求解s h e r m a n 2 。原因是不能寻找到有效的稀 疏模式来近似逆矩阵。而且在这些预条件子可以在最大迭代次数之内得到满足精 度需要的解的时候，迭代法的收敛速度也是非常慢的。没有预条件而直接采用 b i - c g s t a b 方法求解时，有3 个问题失败。 在减少迭代次数上，i l u 仍然比s a i 方法更优，至少在平均意义上。如果只 比较迭代时间，a i n v 对3 个实验矩阵最快，s p a i ，m r p 和无预条件子的情况各 对一种矩阵最快。总体平均来讲，a i n v 迭代时收敛最快，其次是s p a i 。而且， s p a i 方法在几乎所有问题上，无论从强壮性还是收敛速度上，都优于m r 方法。 而且，平均上来看，构造m r 的花费不比s p a i 小。但是，m r 算法有一个优势， 它占用的内存空间更小。s a i 技术在并行向量机上将会有更好的表现，而且有趣 的是，这类方法在电路模型问题上表现的比较好。 在预条件子构造花费上，a i n v 是最好的一个，而且它也是三个较快求解所 有实验问题的中的一个。s p a i 、m r p 和m r 都比较昂贵。当然，情况可能发生变 化，如果这些算法在并行性计算机上执行。s p a i 和m r 类的预条件子有很好的并 行性，而a i n v 相比较就差些。 通过这些实验，得到下面结论： i l u 类的预条件子强壮性最好； a i n v 是s a i 方法里最强壮的，最有效的预条件子； 基于矩阵分解的近似逆预条件子比没有分解的更有效； s p a i 比m r 方法更强壮，更有效。 通过上面两节的分析比较，可以看出，预条件子的强壮性有效性，在s p d 和 非s p d 矩阵上的差异并不大。总体来说，基于矩阵分解的，即第二类稀疏近似逆 预条件的效果比第一类的好，尤其在非s p d 的情况下。 2 4 总结 稀疏近似逆是一个稀疏矩阵膨，它是一般稀疏矩阵彳之逆彳一的一个好的近 似。寻求有效稀疏近似逆的主要动力是它们在并行计算中的潜在优势。而且，有证 电子科技大学硕士学位论文 据显示这种类型的预条件子可解一些i l u 预条件子难以求解的问题【7 j ，并且提供 了传统i l u 预条件子的一种替代。 当前的稀疏近似逆预条件技术大致有二种：基于f r o b e n i u s 范数极小化的稀 疏近似逆和基于矩阵分解的稀疏近似逆。每一种都有一些不同的构造方式且每一 种都有其优点和缺陷。换言之，在所有可比较的规则中( 构造花费、应用花费、健 壮性、有效性等) ，没有哪种是绝对最优的。 总体来说，在序列计算机上，i l u 类的预条件子普遍得到了较好的效果。在 s a i 方法中，第二类s a i