（计算数学专业论文）求解广义对称特征值问题的块jacobidavidson方法.pdf

，“ 1，u月*10 0j_o珀 一 l i 艟 ， b l o c k g e n e r a l ， o下j&吩飞， 承诺书 本人郑重声明：所呈交的学位论文，是本人在导师指导下，独立进 行研究工作所取得的成果。尽我所知，除文中已经注明引用的内容外， 本学位论文的研究成果不包含任何他人享有著作权的内容。对本论文所 涉及的研究工作做出贡献的其他个人和集体，均己在文中以明确方式标 明。 本人授权南京航空航天大学可以有权保留送交论文的复印件，允许 论文被查阅和借阅，可以将学位论文的全部或部分内容编入有关数据库 进行检索，可以采用影印、缩印或其他复制手段保存论文。 ( 保密的学位论文在解密后适用本承诺书) 作者签名： 日 期： 南京航空航天大学硕士学位论文 摘要 j a c o b i - d a v i d s o n 方法是求解广义对称特征值问题极端特征对的一种有效方法。本文将其进 行块推广并对其使用调和策略，提出了调和块j a c o b i - d a v i d s o n 方法，新方法可以有效计算广义 对称特征值问题的重或密集内部特征值。为了有效求解校正向量，对块j a c o b i - d a v i d s o n 方法的 校正方程进行了预处理。另外，为了减少块j a c o b i - d a v i d s o n 方法的计算量与存储量，本文还给 出了块j a c o b i - d a v i d s o n 方法的重新开始及收缩技术。数值实验表明：块j a e o b i d a v i d s o n 方法对于 计算重或密集极端特征值是有效的，调和块j a e o b i d a v i d s o n 方法对于计算谱内部特征值是有效 的。 关键词：广义对称特征值问题，j a c o b i - d a v i d s o n 方法，预处理，调和方法，重新开始及收缩 求解广义对称特征值问题的块j a c o b i - d a v i d s o n 方法 a b s t r a c t j a c o b i - d a v i d s o nm e t h o di s v e r ye f f i c i e n t f o rc o m p u t i n gt h ee x t r e m ee i g e n p a i r so ft h e g e n e r a l i z e ds y m m e t r i ce i g e n p r o b l e m s t h i sp a p e rp r e s e n t st h eb l o c kv e r s i o no ft h em e t h o da n d a p p l i e st h eh a r m o n i cs t r a t e g yt ot h eb l o c kj a c o b i - d a v i d s o nm e t h o d , t h e nd e r i v e st h eh a r m o n i cb l o c k j a c o b i d a v i d s o nm e t h o db yw h i c ht h em u l t i p l e0 1 c l u s t e r e di n t e r i o re i g e n v a l u e so ft h eg e n e r a l i z e d s y m m e t r i ce i g e n p r o b l e m sc a nb ec o m p u t e de f f i c i e n t l y t oa p p r o x i m a t et h ec o r r e c t i o nv e c t o r s e f f i c i e n t l y , t h ep r e c o n d i t i o n i n gt e c h n i q u eo ft h ec o r r e c t i o ne q u a t i o ni sc o n s i d e r e d m o r e o v e r , i no r d e r t or e d u c et h ec o m p u t a t i o n a lc o s ta n ds t o r a g eo ft h em e t h o d ，t h er e s t a r t i n ga n dd e f l a t i o nt e c h n i q u e sa r e d i s c u s s e d n u m e r i c a le x p e r i m e n t ss h o wt h a tt h eb l o c kj a c o b i - d a v i d s o nm e t h o di se f f i c i e n tf o r c o m p u t i n gt h em u l 石p l eo rc l u s t e r e de x t r e m ee i g e n p a i r sa n dt h eh a r m o n i cb l o c kj a c o b i - d a v i d s o n m e t h o di se f f i c i e n tf o rc o m p u t i n gt h ei n t e r i o re i g e r r p a i r so ft h eg e n e r a l i z e ds y m m e t r i ce i g e n p r o b l e r u s k e yw o r d s ：g e n e r a l i z e ds y m m e t r i ce i g e n p r o b l e m s ；b l o c kj a c o b i d a v i d s o nm e t h o d ；h a r m o n i cm e t h o d ； p r e c o n d i t i o n i n g ；r e s t a r t i n ga n d d e f l a t i o n i i 譬 南京航空航天大学硕士学位论文 目录 第一章绪论一1 1 1 背景介绍。1 1 2 记号和约定3 1 3 正交投影方法的理论知识3 第二章j a c o b i d a v i d s o n 方法及其块推广6 2 1 标准特征值问题的j a c o b i - d a v i d s o n 方法6 2 2 广义对称特征值问题的j a c o b i - d a v i d s o n 方法7 2 3 广义对称特征值问题的块j a c o b i - d a v i d s o n 方法9 第三章块j a e o b i - d a v i d s o n 方法的调和策略与预处理技术1 1 3 1 调和块j a e o b i d a v i d s o n 方法1 1 3 2 调和块j a e o b i d a v i d s o n 方法的预处理技术1 2 第四章块j a e o b i - d a v i d s o n 方法的重新开始与收缩技术1 7 4 1 块j a e o b i d a v i d s o n 方法的收缩技术1 7 4 2 动态稠密重新开始块j a e o b i - d a v i d s o n 方法1 8 第五章数值实验2 l 5 1 校正方程中a 1 9 l 曰的合理选取2 l 5 2j a c o l ：i i - d a v i d s o n 方法和块j a e o b i - d a v i d s o n 方法的比较2 2 5 3 块j a c o b i d a v i d s o n 方法使用不同预处理算子及调和策略的比较2 3 5 4 调和块j a e o b i - d a v i d s o n 方法的动态重新开始及收缩技术的应用2 7 第六章结论。2 9 参考文献3 0 j 醪c 谢3 3 在学期间的研究成果及发表的学术论文3 4 m 求解广义对称特征值问题的块j a c o b i d a v i d s o n 方法 i v 图表清单 图5 1j d 方法2 3 图5 2b j d 方法2 3 图5 3 三对角预处理b j d 方法2 6 图5 4 五对角预处理b j d 方法2 6 图5 5 三对角预处理耶j d 方法2 6 图5 6 五对角预处理 i b j d 方法2 6 表5 1 不同a 1 9 l b 对应的b j d 方法的数值结果2 1 表5 2j d 方法与b j d 方法的数值结果2 2 表5 3 使用不同预处理矩阵策略的b j d 方法的数值结果2 4 表5 4 使用三对角预处理的计算结果2 5 表5 5 使用五对角预处理的计算结果2 5 表5 6b i d 方法与h b j d 方法的比较2 7 表5 7d b j d 方法与d h b j d 方法的比较2 8 南京航空航天大学硕士学位论文 第一章绪论 1 1 背景介绍 矩阵特征值问题是数值代数领域的一个基本课题之一，其研究具有重要的理论意义和应用 价值。从数学角度来看，矩阵特征值问题大多来自数学物理方程、差分方程、m a r k o v 过程等； 从工程应用角度来看，机械、航空、航天、建筑、电子等大型结构的动力分析或稳定性分析都 涉及求解大规模矩阵特征值问题。 设彳为 阶实( 或复) 矩阵，所谓特征值问题就是求复数a 和非零向量x ，使得a x = a x 。 通常称a 为彳的特征值，工为么对应于特征值a 的特征向量。对予计算兄和x ，一个最直接 的方法是先求特征方程d e t ( a i 一彳) = 0 的根即特征值旯，再将其回代到方程a x = a x 中求相应 的特征向量。但在实际应用中，除一些特殊矩阵外，该方法对一般矩阵不可行。首先，若矩阵彳 的阶数较大，则计算特征多项式d e t ( a i 一彳) 的计算量非常大；其次，特征方程d e t ( 村- 毋= o 的根通常是不稳定的，即特征方程系数的微小扰动将导致方程根的较大误差。因此，人们一直 研究特征值问题的数值解法，经过一百多年的探索，人们提出了许多方法，促进了该领域的不 断发展。 目前，求解大型矩阵特征值问题的方法主要分为两类：变换方法和投影方法。变换方法利 用了矩阵特征值的相似不变性原理，直接对原矩阵进行处理，通过一系列变换( 如h o u s e h o l d e r 变换，g i v e n s 变换等) ，将其转化为一个容易求解的特殊矩阵( 如h e s s e n b e r g 矩阵，三角形矩阵， 拟三角形矩阵等) ，然后通过j a c o b i 方法，q r 方法等求解。变换方法收敛速度快，计算结果精 确可靠，但也存在一些不足，如变换方法存储量较大，且一般的变换方法会破坏矩阵的稀疏性。 因此，交换方法只适宜求解阶数较低的稠密矩阵特征值问题。投影方法是通过投影算子把阶数 较高的矩阵特征值问题投影到低维子空间中，将大规模矩阵特征值问题转化为中、小型矩阵特 征值问题，然后采用变换方法计算其特征对作为原矩阵特征对的近似。由于投影方法可采用压 缩存储技术，因而它是目前求解大型稀疏矩阵特征值问题的一类常用方法。对于大型稀疏对称 矩阵，目前常用的投影方法有子空间迭代法、l a n c z o s 方法、a m o l d i 方法等。子空间迭代法可 同时求解几个极端特征值和相应的特征向量，但它收敛速度较慢，运算量较大。l a n c z o s 方法 利用三项递推产生l a n c z o s 向量，并将矩阵彳投影为低阶三对角矩阵，然后求解低阶三对角矩 阵特征值问题，用r i t z 对逼近矩阵么的特征对。该方法存储量较小，计算速度较快，但是l a n c z o $ 向量之间容易失去正交性，通常需要使用一些重正交化技术。 1 9 7 5 年，d a v i d s o n 1 在计算量子化学中的特征值问题时提出了计算大型对称矩阵极端特征 值问题的一种新方法- - d a v i & o n 方法。最初，d a v i d s o n 方法是通过摄动方法得到的，其主要思 l 求解广义对称特征值问题的块j a c o b i d a v i d s o n 方法 想是用一个包含特征方向更多信息的向量扩充搜索子空间从而提高逼近效果。大量的数值试 验表明，d a v i d s o n 方法对于求解对角占优矩阵的特征值问题是非常有效的。后来许多学者对 d a v i d s o n 方法进行了研究和改进【2 6 】，如c r o u z e i x 等人【5 】将d a v i d s o n 方法推广到了块形式， 并给出了d a v i d s o n 方法的收敛性分析。d a v i d s o n 方法已发展成为求解对角占优h e r m i t i a n 矩阵 若干个最小或最大特征值的有效方法。然而，对非对角占优h e r m i t i a n 矩阵，d a v i d s o n 方法的 计算效果不尽如人意，而且精度很高的预条件子可能会导致收敛速度减慢甚至停滞。 1 9 9 6 年，s l e i j p e n 和v a nd e rv o r s t 7 将j a c o b i 方法的校正思想和d a v i d s o n 方法的内外迭代 思想相结合，提出了j a e o b i d a v i d s o n 方法：首先在近似特征向量的正交补空间上求解校正方程， 然后用求得的校正向量扩充搜索子空间，最后在搜索子空间中计算r i t z 对作为所求特征对的近 似。j a c o b i d a v i d s o n 方法的关键是校正方程的有效求解，其良好的收敛性质使得校正方程只需 少数迭代与适当的计算精度即可。另外，通过投影得到的校正向量可以避免d a v i d s o n 方法中可 能的停滞，且j a c o b i - d a v i d s o n 方法允许应用高精度的预条件子。目前，j a e o b i - d a v i d s o n 方法被 认为是计算特征值问题的最有效方法之一，特别是这一方法对于计算中间特征值有优势。另外， j a c o b i d a v i d s o n 方法还具有较好的稳定性，对于非对角占优、非正规矩阵也可能达到较快的收 敛速度。许多学者对j a c o b i - d a v i d s o n 方法进行了研究和改进，如s l e i j p e n 等人【8 ，9 】将其应用到 广义特征值问题以及多项式特征值问题。v a nd e ne s h o f f l o 讨论在对称情况下j a c o b i - d a v i d s o n 方法的收敛性问题。n o t a y 1 1 】将共轭梯度方法与j a c o b i - d a v i d s o n 方法有效结合，利用共轭梯 度法求解校正方程，并且根据共轭梯度法的性质建立了内部线性系统残量与外部过程收敛性的 关系。此外，戴小英等【1 2 】讨论了在j a c o b i - d a v i d s o n 方法中采用区域分解和多水平技术产生校 正方程的预条件子。当所求的特征值分布密集或矩阵阶数很高时，投影子空间的维数可能较大， 从而j a e o b i - d a v i d s o n 方法需要的计算量及存储量很大。为解决这个问题，一些学者研究了 d a v i d s o n 方法和j a c o b i d a v i d s o n 方法的重新开始技术，如隐式重新开始技术【1 3 】、动态稠密重 新开始技术【1 4 】、组合c g 递归与动态稠密重新开始技术【1 5 】等。 文献【8 ，9 】将j a c o b i - d a v i d s o n 方法应用于广义特征值问题a x = 2 b x ，其中矩阵彳，召是给 定的刀l 矩阵，且假定b 正定，并分析了如何使用合理而又有效的投影方法将原来阶数较高的 特征值问题转化为阶数较低的特征值问题，以及如何正确选取测试子空间和搜索子空间，此外 还分析了如何对校正方程进行合理的预处理以及如何构造增广校正方程等问题。数值结果表明， j a c o b i - d a v i d s o n 方法是求解广义特征值问题的有效方法之一。然而，当要求的特征值是重或密 集特征值时，j a c o b i - d a v i d s o n 方法的有效性和可靠性会大大降低，并且当所求特征值是谱内部 特征值时，j a e o b i - d a v i d s o n 方法的收敛性常常表现得没有规律，有时还会出现“假”特征值。为 解决j a c o b i - d a v i d s o n 方法计算重或密集特征值有效性和可靠性低的问题，本文研究广义对称特 征值问题j a c o b i - d a v i d s o n 方法的块格式，提出求解广义对称特征值问题的块j a c o b i - d a v i d s o n 2 方法。为了应用 略 1 6 1 9 1 ，提出 全文安排如 准特征值问题的 行块推广；第三章讨论了块j a c o b i - d a v i d s o n 方法的调和策略与预处理技术；第四章研究了块 j a e o b i d a v i d s o n 方法的重新开始与收缩技术；第五章给出数值实验；第六章对本文进行总结和 展望。 1 2 记号和约定 本文用r ”表示n 维实向量空间，歙“”表示刀咒实对称矩阵的全体，用大写字母表示矩 阵或者是线性算子，小写字母表示向量。 厶表示刀阶单位矩阵，若不产生混淆，用，表示。 设a s r “4 ，b s r “。，称a 和曰为广义对称特征值问题a x = a 擞的矩阵束，记为 ( 彳，b ) 。 ，五，以表示矩阵束( 彳，b ) 的特征值，x z ，屯，表示相应正交规范特征向量( 通 常意义下的正交规范) ，且满足 如丸。a x = l b x 的特征对表示为( a ，砷，其中名为 特征值，工为相应的特征向量。 本文涉及的特征向量均为正交规范特征向量，所给的范数若没有特别声明均为2 范数 ，和x r 分别表示向量x 和矩阵x 的转置。 列正交规范矩阵v r “是指满足矿。v = l 的矩阵。 m g s 表示m o d i f i e dg r a m s c h m i d t 正交化过程。 工上y 表示向量x 和y 正交，x 上u e r 蹦”表示向量x 和u 的列所张成的子空间正交。 如果m ，吃，是m 个刀维向量，那么v = 【m ，屹，】表示以m 矩阵，其中m 为y 的 第i 列。s p a n ( v ) 表示d 3 v 的列所张成的子空间。s p a n ( v ) r 表示子空间s p a n ( v ) 的正交补。 1 3 正交投影方法的理论知识 考虑如下标准特征值问题 a x = 3 , x ， ( 1 1 ) 其中a 是r “4 中的一个对称矩阵，令圪= h ，v 2 ，】是咒维空间中的m 个线性无关的向量 组成的刀m 矩阵，既= s p a n ( v ) 是由m ，屹，所张成的子空问。对子空间瓦的正交投 影方法就是寻求问题( 1 1 ) 的近似特征值a ( ”和近似特征向量伊( 驯，满足 笳- 瓦a , 恤。 m 2 ， 【彳伊辨”妒用上e 。 、7 3 求解广义对称特征值问题的块j a c o b i d a v i d s o n 方法 记缈所= v y 刖，其中y ”是一个朋维列向量，则( 1 2 ) 可改写为 ( 吃一兄“c 卅) y “) = 0 ， ( 1 3 ) 其中吃= 圪r a 圪，g = 圪r 圪。可见兄“和y ”是广义特征值问题( 1 3 ) 的特征值和特征向 量。 若匕是列正交规范矩阵，则问题( 1 3 ) 可转化为 ( 圪r 彳圪一旯4 d y 帕= 0 。 ( 1 4 ) 若表示子空间瓦上的正交投影算子，则( 1 2 ) 可写成 x , 。( a - 1 ) i ) t p ”) = 0 。 ( 1 5 ) 我们称问题( 1 2 ) ，( 1 3 ) ，( 1 5 ) 为( 1 1 ) 的近似问题，近似问题的所有解力伽) 称为矩阵彳对子 空间乜的r i t z 值，相应于旯”的向量妒”称为r i t z 向量。 为了便于理论分析，我们引入线性算子以= 彳。若圪是邑的一组规范正交基，则 矩阵彳相对于基圪的矩阵表示是圪r a 圪。因为缈”= 伊，于是有如下定理。 定理1 1 2 0 1 子空间已中的r i t z 值允( ”) 和r i t z 向量9 伽) 分别是4 ，l 的特征值和特征向量。 应用正交投影反方法求得的r i t z 值2 t )r i t z 向量缈4 是否收敛到特征值问题( 1 1 ) 的特征 值a 和特征向量x ? 事实上，当历= n 时，r i t z 值兄( “和r i t z 向量矽( 所就是特征值问题( 1 1 ) 的 特征值和特征向量。当m n 时，如果4 。的特征值允( ”) 和相应的特征向量缈( 帕接近彳的特征 值允和特征向量z ，则残量( 彳一九“，) 妒”是小量。 关于0 ( a 一, 无o x l l ，有如下估计式。 定理1 2 【2 0 】令= 0 彳u 一) l l ，则 i i ( 4 ，l 一旯，) 石忙h 2 + 2 一) x i l 。 ( 1 6 ) 由定理1 2 可见，特征向量工与子空间l 之间的距离i i ( ，一) 刮是一个非常重要的量， 通常9 ( ，一) 刈不是小量，当且仅当聊= 靠时，0 ( ，一) z i i = o 。而当所 以时，一是一 个正交投影算子，i i ( ，一刀) x 0 并不随朋的增大而必然减小。但是如果子空间瓦比较特殊，可 以保证8 ( ，一) 刘随着m 的增大而迅速减小。 定义1 1 【2 1 】对n 阶实对称矩阵a ，以m 列正交规范矩阵q ，令h = q 2a q ，假设日的 某一特征对为( 旯，y t m ) ) ，则称允”为矩阵a 在子空间s = s p a n ( q ) 上的r i t z 值，简称r i t z 值，称缈4 = ”为相应的r i t z 向量。用a ”和缈( “来逼近矩阵么的特征值和特征向量，称 此逼近过程为r a y l e i g h - r i t z 逼近过程。 事实上，正交投影方法就是r a y l e i g h r i t z 逼近过程的应用。 另外，假设子空间s 是么的一个m 维不变子空间，q 是由s 的正交规范基为列构成的 n m 矩阵，则必有m 扰矩阵日，使得 4 个不同的自然数，f 2 ，使得 i 九一段i 0 r ( 艿) 0 ，k = l ，2 ，m 。 ( 1 8 ) 定理1 3 说明了矩阵刀和彳的特征值之间的关系，只要0 r ( 曰) 0 充分小，则矩阵b 的特征 值就可以作为a 的特征值的较好近似。显然，i i r ( b ) i i 的大：j 、与矩阵曰的选取紧密相关。下面 的定理给出了当矩阵b 取何值时，l i r ( b ) i i 獭：j 、。 定理1 4 【2 2 】假设矩阵q 和b 满足定理1 3 的条件，记 r ( 曰) = 彳q 一妒，h = q r 彳q ， 则 i i r ( h ) i i - i i r ( b ) i i 。 对形如a q 一妒的矩阵，只有当b = q r a q 时，4 q q 召的范数才达到最小。因此，正 交投影方法是一种最佳逼近，在求解特征值问题时，若采用子空间方法进行逼近，一般都会采 用此方法。 5 求解广义对称特征值问题的块j a c o b i d a v i d s o n 方法 第二章j a c o b i d a v i d s o n 方法及其块推广 本章讨论标准和广义特征值问题的j a c o b i d a v i d s o n 方法。为了有效计算广义对称特征值问 题的重或密集特征值，文章提出求解广义对称特征值问题的块j a c o b i - d a v i d s o n 方法。 2 1 标准特征值问题的j a c o b i d a v i d s o n 方法 1 9 7 5 年，d a v i d s o n 【1 】在计算量子化学中的特征值问题时提出了计算大型稀疏对称矩阵极端 特征值问题的一种新方法d a v i d s o n 方法，该方法是r a y l e i g h - r i t z 过程和预处理过程的有效结 合，其中对于预处理矩阵和重新开始向量的合理选取非常重要。最初，d a v i d s o n 方法是通过摄 动方法得到的，其主要思想是用一个包含特征方向更多信息的向量扩充搜索子空间从而提高逼 近效果。若d a v i d s o n 方法没有预处理，则它与l a n c z o s 方法等价。大量数值试验表明，对于对 角占优矩阵，d a v i d s o n 方法的收敛性优于l a n c z o s 方法，因为d a v i d s o n 方法能够使得不需要的 特征向量的分量从所求近似特征予空间的基中逐渐减小，从而提高方法的收敛性；而对于非对 角占优矩阵，d a v i d s o n 方法的计算效果并不明显，且精度很高的预条件子常会导致收敛速度减 慢甚至停滞。s l e i j p e n 和v a nd e rv o r s t 【7 】于1 9 9 6 年将j a e o h i 方法的校正思想和d a v i d s o n 方法 的内外迭代思想相结合，提出了j a c o b i - d a v i d s o n 方法，该方法求解特征值问题的关键是校正方 程的有效求解，其良好的收敛性质使得校正方程只需少数迭代与适当精度即可。另外，由投影 得到的校正向量可以避免d a v i d s o n 方法中可能出现的停滞，j a e o b i - d a v i d s o n 方法还允许应用高 精度的预条件子。 对于标准特征值问题血= , t x ，将其投影到以标准正交基m ，v 2 ，吃张成的子空间 s p a n ( m ，1 2 ，咋) 上，用矩阵彳在该子空间上的r i t z 对作为所求特征对的近似。i a ( s i ，& ) 是 投影矩阵的特征对，圪= 【m ，k 】，则矩阵彳在子空间s p a n ( v 1 ，v 2 ，唯) 上的r i t z 对 ( 最，u i ) 满足r i t z - g a l e r k i n 条件 i ( a u i 一最) 上s p a n ( v l ，v 2 ，吃) l u t = 圪 、 即圪r a v k s k 一最& = 0 。 假设矩阵么在子空间删栉( ) 上上的投影矩阵为j = u 一蚝“k r ) a ( i 一咋r ) ，a p _- a = a + “i u l l a + a u i u l 2 一t g k t $ 七l l l 2 。( 2 1 ) 为了对r i t z 向量进行合理校正，假设在子空间印口，z ( ) 上中求得的校正向量气，满足 二4 ? t + 气) = 旯( 甜t + 气) 。 占2 ) 【气上 、7 6 南京航空航天大学硕士学位论文 将( 2 2 ) 式带入( 2 1 ) 式可得 ( a - a i ) t k = - r , + ( 允一b u 1 1 a u i ) “i ， 其中气= ( 么一, g k i ) u t ，则 “l r ( 7 4 一甜) = 心r 气+ ( 允一最一u k r a u t ) u k r “t 。 由向量和气与的正交性可得五一最- u , r a u k = 0 ，所以 ( 彳一九d 气= 一气。 因为所求特征值允未知，计算过程中用r i t z 值最代替允。由7 4 = ( i 一“k r ) a ( i 一心r ) 及 t k = ( i - u t u k t ) 可得j a c o b i - d a v i d s o n 方法的校正方程 f ( ，一“i u k r ) ( 彳一号t d ( j 一“t u k t ) 气= - r k k 上 9 由上述过程可得求解标准特征值问题的j a e o b i d a v i d s o n 算法。 算法2 1 标准特征值问题的j a c o b i - d a v i d s o n 方法 1 ) 选取最大迭代次数朋，允许误差拓，正交规范向量v ，令k - v 】，u = 1 ，； 2 ) 对k = 1 ，2 ，m - 1 ，执行 2 1 ) 计算q = k r 彳巧，1 9 = “r a u ，= 彳“一抛； 2 2 ) 求解校正方程 f u u u ，) ( 彳一聊一u u ，弦= 一， 1 z 上扰 。 2 3 )由校i f 向量z ，通过嬲正交化过程将搜索子空间圪扩充至圪+ l ； 2 4 ) 计算吼+ = 圪+ i r 4 圪的特征对( 口，y ) ，其中y 满足眇0 ：= 1 ； 2 5 ) 计算r i t z 向量“= 圪+ l y ，相应残量厂= a “一o a u ； 2 6 ) 检测收敛性，若收敛则算法停止，否则七= k + l ； 3 ) 重新开始：记巧= 【u 】，返回2 ) ； 若对最大迭代次数m 不加限制，则算法至多迭代n 次便可得到矩阵a 的特征值，因为当 m = n 时，矩阵峨= k r 彳圪与矩阵4 相似，从而矩阵以与矩阵彳特征值相同。但实际计算 时，由于计算量与存储量的限制，算法2 1 需采用重新开始策略，通常m n ( 一般取m 为所 求特征值个数的4 9 倍) ，即当迭代步数达到最大值m 时，用近似r i t z 向量作为新的初始迭代 向量，返回2 ) 。 2 2 广义对称特征值问题的j a c o b i d a v i d s o n 方法 考虑广义对称特征值问题 7 求解广义对称特征值问题的块j a c o b i d a v i d s o n 方法 a x = 2 b x ， ( 2 3 ) 其中矩阵彳，曰为给定的，z n 对称矩阵，且假定矩阵b 正定。假设( 口，“) 是( 2 3 ) 式某个特征 对( a ，x ) 的近似，则由r i t z - c r a l e r k i n 条件可得 厂暑au 1 9bu 上u 。 与标准特征值问题的j a e o b i - d a v i d s o n 方法推导过程类似，在子空间s p a n ( u ) 上中求u 的校 正向量z ，满足条件 a ( u + z ) = 2 b ( u + z ) 。( 2 4 ) 对( 2 。4 ) 式两边同时作用向量u 以及投影矩阵1 一u “7 ，可得 i 允：晕坠堕 u 1 b ( u + z ) ， ( 2 5 ) 【( ，一u u r ) ( 彳一a b ) z = 一( 口一肋) ，昆上“ 其中口=华，口柏一，=华，6=bu一肛，口=舌，一obuu uu 。 对于( 2 5 ) 式中第二式，由于所求特征值五未知，计算过程中用近似特征值汐代替旯。由 z 上u 可得z = ( ，一u u7 ) z ，所以( 2 5 ) 式中第二式可化为 ( i u ur ) ( 彳一ob ) ( ，一u ur ) z = 一，。( 2 6 ) 校正向量zj m j 立m g s 正交化过程可将投影子空间圪扩充至圪+ l ，然后在投影子空间圪+ l 中求解所求特征对( 旯，x ) 的近似。 综上所述，可得求解广义对称特征值问题的j a e o b i d a v i d s o n 方法。 算法2 2 广义对称特征值问题的j a c o b i - d a v i d s o n 方法 1 ) 选取正整数朋，允许误差幻z ，正交规范向量v l ，记k = 【v 1 】，u = h ； 2 ) 对k = 1 ，2 ，m - 1 ，执行 2 1 ) 计算口：r u r a u ，r ：a u - o b “； u 1b u 2 2 ) 求解校正方程 ( ，一u u r ) ( 彳一o b ) ( i u u r ) z = 一，： 2 3 ) 由校正向量z ，通过嬲正交化过程将投影子空间圪扩充至圪+ l ； 2 4 ) 计算圪+ 。r 彳，y 一彤+ 。r b 圪+ 。y = 0 的特征对( 口，y ) ，其中y 满足陟i i ：- - 1 ； 2 5 ) 计算r i t z 向量“= 圪+ 1 y 及相应残量，= a “一o b u ； 2 6 ) 检测收敛性，若收敛则算法停止，否则k = k + 1 ； 3 ) 重新开始：记k = 【u 】，返回2 ) ； 8 南京航空航天大学硕士学位论文 2 3 广义对称特征值问题的块j a c o b i d a v i d s o n 方法 在对j a e o b i - d a v i d s o n 方法的研究中发现当所求特征值为重或密集特征值时， j a e o b i - d a v i d s o n 方法的有效性会降低，可能出现下列情形， 1 )当近似特征对( 占，“) 与所求特征对 ，石) 接近时，校正方程的系数矩阵可能变得病态。 2 ) 对于重特征值问题，j a e o b i d a v i d s o n 方法只能有效求解部分特征对。 3 )校正方程的系数矩阵可能奇异，使得算法的收敛速度变慢，甚至出现较大偏差。 4 )对于矩阵良态，且其密集特征值也良态，采用块方法可达到较高的收敛性【2 1 】。 基于上述原因，下面导出广义对称特征值问题的块j a e o b i d a v i d s o n 方法。 对于广义对称特征值问题( 2 3 ) ，假设特征值为a ，乃，以，且a 如丸，相应的 正交规范特征向量( 通常意义下的正交规范) 为x a ，恐，假设( 逸，咋) o = l ，2 ，) 是 a x = 2 b x 的z 个所求特征对的近似，相应残量分别为= a u ，一 9 ，b u ，( f = 1 ，2 ，) 。记 u = 【均，u ，】，a = d i a g c & ，丑) ，r = k ，r t 】，a i = d i a g ( o i ，占，) ，则 r=a u b u a l ， 由r i t z - g a l e r k i n 条件可得 r 害彳u 一艿u 人1 上u 。 因为( 母，u ，) ( f = 1 ，2 ，) 是所求特征对的近似，所以希望求得矩阵u 的校正z 上u 满足 a ( u + z ) = 曰( u + z ) a 。( 2 7 ) 用子空间s p a n ( u ) 上上的正交投影矩阵i u ur 对( 2 7 ) 式进行投影处理，可得 ( ，一【，u2 ) ( 彳z 一曰za ) = - ( i u u1 ) ( 彳u bu 人) 。 ( 2 8 ) 由于所求特征值人未知，在计算过程中用近似特征值人，作为人的近似，从而( 2 8 ) 式可化为 ( ，一u u ) ( 彳一占，丑) z ，= 一 1 3 f ：1 ，2 ，) 。 ( 2 9 ) 【z l 上【， 又由矩阵z 上u ，所以乃= ( ，一u u l ) 毛，则( 2 9 ) 式可化为 ( ，一u u r ) ( 彳一, g i b ) ( ，一u u ，) z ，= 一 ( f ：1 ，2 ，) 。( 2 1 0 ) 【z f 上u 使用校正向量z ，( f = l ，2 ，z ) 可将投影子空间圪扩充至珞+ l ，然后在子空间圪+ l 中计算 新的近似特征对瓴，丑) ( f = 1 ，2 ，) 。 综上所述，可得求解广义对称特征值问题的块j a c o b i d a v i d s o n 方法。 算法2 3 广义对称特征值问题的块j a c o b i - d a v i d s o n 方法 1 ) 选取正整数朋，允许误差f d ，目标值f ，迭代块大小p ，列正交规范矩阵k ； 2 ) 对k = 1 ，2 ，聊- 1 ，执行 9 求解广义对称特征值问题的块j a c o b i - d a v i d s o n 方法 2 1 ) 计算圪r a v k y = 1 9 圪r b v k y 的，个最靠近f 的特征对( 母，咒) ( f = 1 ，2 ，) ， 并计算r i t z 向量“f = 圪y l 以及相应残量= a u f 一母舰f1 3 f = 1 ，2 ，1 ) ； 2 2 ) 检测收敛性，若收敛，则算法停止，否则进行下一步； 2 3 ) 求解1 个校正方程 ( ，_ u u r ) ( 彳一1 9 f 召) ( ，一u 【，r ) 毛2 一 o ：1 ，2 ，) ； 【z ，上u 、 2 4 ) 由校正向量刁( f = 1 ，2 ，z ) ，通过删正交化过程将k 扩充至吆+ l ； 3 ) 重新开始：记k = z l ，z 2 ，乃】，返回2 ) ： 迭代块大小p 可有多种选择，如果所求特征值为，重，则可取p 大于或等于，。然而，由 于一般z 事先未知，因此在实际计算中，p 的大小可以进行动态改变。另外，校正方程可以使 用m i n r e s 方法【2 3 】求解。 1 0 k 第三章块j a c o b i - 。a v i d s 。n 方法的调和策略与预处理技术 在块j a e o b i d a v i d s o n 方法中，近似特征对的选取对投影子空间的扩张具有重要意义。对于 p e t r o v - g a l e r k i n 投影方法，若计算极端特征值，则r i t z 值作为所求特征值的近似非常适合，而 计算内部特征值( 或者由于舍入误差等原因) ，算法的收敛性常常没有规律，还可能会产生“多 余”( s p u r i o u s ) r i t z 值。由块j a e o b i - d a v i d s o n 方法的校正方程可以看出，若近似特征对选取不合 适，则搜索子空间可能得不到有效扩展。调和策略【7 ，9 ，1 3 , 2 4 】对于内部特征值的求解具有较好性 质，可以有效克服上述问题，且可避免矩阵求逆及矩阵分解等问题。另外，校正方程( 2 1 0 ) 是块 j a e o b i - d a v i d s o n 方法的重要组成部分，合理的选取预处理算子可以减少块j a e o b i - d a v i d s o n 方法 的计算量与存储量。 本章首先对块j a e o b i d a v i d s o n 方法使用调和p e l r o v - g a l e r k i n 投影方法，其次，讨论如何使 用预处理算子求解校正方程。 3 1 调和块j a c o b i d a v i d s o n 方法 对于广义对称特征值问题( 2 3 ) ，假设求解，个靠近给定值f ( f 不是矩阵束( 么，占) 的特征值) 的内部特征值( 以，而) ( f = 1 ，2 ，z ) 由( 2 3 ) 式可得 ( a - r b ) x = ( 兄一r ) b x ， 则 ( 彳一f b ) b x ：_ l 石， 即靠近给定值f 的特征值对应于矩阵( 彳一r b ) - 1 b 的极端特征值。因为对于标准特征值问题， r a y l e i g h - r i t z 投影方法对于求解极端特征值效果明显，所以对矩阵( a - r b ) - 1 b 应用 r a y l e i g h - r i t z 投影方法，为避免使用矩阵求逆，测试子空间取为( a - r s ) v , ，则 f ( a - r b ) - l 觑一击z 上( 彳川彤， 【x 圪 即 肛么- r b ) ( a - 2