（计算数学专业论文）计算最小奇异三元组的隐式重启精化lanczos双对角化方法.pdf

计算最小奇异三元组的隐式重新启动精化l a n c z o s 双对角方法 摘要 随着科学和技术的发展，越来越多的应用和计算问题需要求解大规模矩阵的少 数，1 1 个最小奇异值及其相应的左( 右漪异向量众所周知，求矩阵的奇异值分解湎常 用l o n c z o sb i d i a g o n a l i z a t i o n ( l b d ) 先把矩阵转化为双对角阵，然后用g o l u b - k a h a n s v d 求出双对角阵的奇异值分解【9 1 ，从而求出矩阵的奇异值分解最早l b d 算法 是g o l u b 等人在文献【1 0 中提出的。后来s o r e n s e n 在文献【2 5 1 提出了种有效的 隐式重新启动策略并成功应用到a m o l d i 方法和l a n c z o s 三对角化方法中，这种方 法可看作是种截断q r 算法( c u r t a i l e dq r ) 再后来助6 哦等人把这种策略应用到 l b d 并得到相应的递推公式 2 1 ，接着贾仲孝通过添加嘴化向量”提出了精化算法 【l 砑最近k o k i o p o u l o u 等人应胃l b d 的递推公式并添加。精化向量”来求解_ 些 奇异值及相应的奇异向量f 1 5 】本文就是在k o k i o p o u l o u 等人的工作匕进行改进的 本文的主要工作是给出种计算少数1 个最小奇异三元组的隐式重新启动精化 l a l l c o z 8 双对角方法，与e 述方法区别在于，我们求的是小的奇异值，因此我们用调 和r i t z 值作为位移，能有效地逼近大规模矩阵的最小奇异值的奇异三元组在算法 中我们还用精化残量，精化奇异向量和精化r a y l e i g h 商，同时采取压缩技术压缩掉 已经求出的小的奇异三元组最后，我们用数值实验表明，我们的算法可以成功地求 解大规模矩阵的小的奇异三元组，收敛速度也较快 本文的结构和内容如下t 在第章，我们在第节对奇异值的背景及数值解法 的研究历史给出个比较详细的概述，第二节介绍了怎样用低维的奇异三元组来逼 近大规模矩阵的奇异三元组在第二章，我仃了首先介绍l b d ，然看给出隐式重启l b d 及讨论位移的选择问题在第三章中我衍证明正交压缩变换( o d t ) 的，1 个结论，给 出精化向量、精化残量及精化p 隰y l e i g h 商的求解第四章，我们给出完整的算法， 最后章是数值实验 关键词，隐式重新启动l a n c z o s 双对角化；调和r i t z 值；正交压缩变换；精化奇异 向量；精化残量 计算最小奇异三元组的隐式重新启动精化l a n c 2 0 s 双对角方法 i v a b s t r a c t m a n ys c i e n t i f i ca p p l i c a t i o n sl e a dt ol a r g e - s c a l es i n g u l a rt r i p l e t sp r o b l e m s ， w h e r et y p i c a l l yo n l yaf e ws m a l l e s ts i n g u l a rv a l u e sa r eo fi n t e r e s t f o rs u c h p r o m b l e m s ，w ea l w a y sl l s el b d t oc h a n g et h em a t r i xt ob i d i a g o n a lm a t r i x ，t h e n f i n dt h es v do ft h eb i d i a g o n a lm a t r i xb yt h eg o l u b - k a h a ns v da l g o r i t h m 9 ，b y d o i n gt h i sw es o l v et h ep r o b l e m l b df i r s td e v i s e db yc o l u b 1 0 ，l a t e rs o r e n s e n i n t r o d u c e st h ei m p l i c i t l yr e s t a r t e da r n o l dm e h t h o d 2 5 t h ei m p l i c i t l yr e s t a r t e d t e c h n i q u e sw a sg e n e r a l i z e dt ol b db yb j 诎。t h e nj i a 1 3 p r e s e n t s r e f i n e da l g o - r i t h m b ya u g m e n t i n g “r e f i n e dv e c t o r r e c e n t l yk o k i o p o u l o ug e n e r a l i z e dl b d b ya u g m e n t i n gr e f i n e dv e c o rf o rf i n d i n gt h es m a l l e s ts i n g u l a rt r i p l e t s 1 5 o u rm a i nw o r ki st h ed e v e l o p m e n to fa ni m p l i c i t l yr e s t a r t d el a n c s o sb i d i g - o n a l i z a t i o nm e t h o dw i t hr e f i n e dv e c t o rf o rc o m p u t i n gs m a l l e s ts i n g u l a rt r i p l e t s o u rm a i np o i n ti st h es m a l l e s ts i n g u l a rt r i p l e t s s ow eu s et h eh a m o n r i zr i t z 嬲t h e s h i f tw h i c hc a ng e ts i n g u l a rt r i p l e tm o r ee f f i c i e n t l y i nt h ea l g o r i t h mw ea l s ou 8 e r e f i n e ds i n g u l a rv e c t o r ，t h er e f i n e dr e s i d u a l ，r e f i n e dr a y l e i g he q u o t i e n ta n dt h e d e f l a t i o nt e c h n i q u e o u rn u m e r i c a le x p e r i m e n ti l h s t r a s t et h a to u ra l g o r i t h mc a n s o l v et h el a r g e - s c a l es i n g u l a rt r i p l e t sp r o b l e me f f i c i e n t l ya n di tc o n v e r g e sr a t h e r f a s t i nc h a p t e rl ，w eg i v ea no v e r v i e wo ft h eb a z k g r o u n do ft h es i n g u l a rp r o b l e m s a n dt h er e s e a r c hh i s t o r yo fs o l v i n gs i n g u l a rp r o b l e m i nc h a p t e r2 ，啊，ef l r s t l yg i v e t h el b da l g o r i t h m t h ei m p l i c i t l yl b dt h e nt h es h i f ts e l e c t i o ns t a t e g i e s w e p r o v es o m ec o n c l u s i o no fo d t a n dp r e s e n tt h er e f i n e dv e c t o r ，r e f i n e dr e s i d u a l a n dt h er e f i n e dr a y l e i g he q u o t i e n ti nc h a p t e r3 i nc h a p t e r4 ，w ed e v e l o pt h e w h o l ea l g o r i t h m ，l a s t l yw ei l l u s t r a t et h en u m e r i c a le x p e r i m e n t s k e yw o r d ：i m p l i c i t l yr e s t a r t e dl b d ；r i t zv a l u e ；h a r m o n i cr i t zv a l u e ；d e f l a t i o n ； r e f i n e ds i n g u l a rv e t o r ；r e f i n e dr e s i d u a l 厦门大学学位论文原创性声明 兹呈交的学位论文，是本人在导师指导下独立完成的研 究成果。本人在论文写作中参考的其它个人或集体的研究成 果，均在文中以明确方式标明。本人依法享有和承担由此论 文而产生的权利和责任 责任人( 签名) ：耘障同 姗7 年j 月庐日7 厦门大学学位论文著作权使用声明 本人完全了解厦门大学有关保留、使用学位论文的规定 厦门大学有权保留并向国家主管部门或其指定机构送交论文 的纸质版和电子版，有权将学位论文用于非赢利目的的少量 复制并允许论文进入学校图书馆被查阅，有权将学位论文的 内容编入有关数据库进行检索，有权将学位论文的标题和摘 要汇编出版保密的学位论文在解密后适用本规定 本学位论文属于 1 、保密() ，在年解密后适用本授权书 2 、不保密( ) ( 请在以上相应括号内打”一) 日期：哆年功日 醐。勿7 年6 凡媚 第一章引言 l 第一章引言 在科学工程计算中经常需要计算大规模矩阵的少数最大或最小奇异 值及其对应的左( 右 奇异向量。例如：大规模矩阵的低阶逼近，图象处 理中图象的分辨率就是矩阵端部奇异值之比，诸如此类的问题还存在于总 体最小二乘问题、控制理论、信号跟踪( 见 4 1 ) 等等中 1 1 矩阵的奇异值与特征值分解 设a c m 煳是个 大二观模矩阵，不妨设m n ，否则我们用代 替令 盯i a 之够扩g 0( 1 1 ) 为a 的奇异值，相应的左奇异向量为让尹c m ，右奇异向量t ，5 ) c n ， 1 j i n ，那么我们有： a v ；a ) 二东篇，1 娜髓 ( 1 2 ) ia 垆= 垆一一” 卜一 与 a = 巧脚( t j 由) 。 其中刚= 阻( a ) ，t t ，】，诺a ) = 雠m ，】都为标准列正交 阵记与a 的最小奇异值相对应的奇异三元组为 砖m ，妒，v j a ， 我们先给出一些重要的关系式( 见【9 ，s e c t i o n8 6 】) ：以的奇异值分解 e s ，与m 舭和乳0 。：) 甜) x ( 一) 的敞蝴的关 系： u 。( a a ) 矽= 击叼 仃；，以，q ，眵 ( 1 3 ) 、_ o - - _ ， 旷( a ) v = d i a g a ；，)( 1 4 ) 分块u = 队，巩】，其中巩有n 列，有m n 列令 y = 击( ：以0 巩三) 第一章引言 2 那么y 是矩阵c 的个标准正交基，而且 y 甜= 出口g ( 一c r l ，一，曼o l ，。，) ( 1 5 ) 这样，我们可以把求解大规模矩阵的奇异值问题归结为相应的h e r m i t i a n 矩阵的特征值问题而后者已有很多专家学者进行广泛深亥8 地研究并有相 应的软件包如j a c o b i 对角化方法、幂法、子空间迭代法、l a n c z v s 方 法，a x n o l d i 方法、以及j a c o b i - d a v i d s o n 方法( 见【1 1 】) 及相应的软件包 然而，我们的目标是一些小的奇异值，这就给我们求解问题带了很多麻 烦( 见f 2 3 ，2 4 】) 事实上，当我们用小a 或a 来计算a 的奇异值时，由式( 1 3 ) ，( 1 4 ) 知需求砰，虽然这样更好地分离了大的奇异值( 1 ) ，但同时也使小的奇 异值( 1 ) 更加聚合在一起，给求解小的奇异问题带来麻烦，而且当a 的条件数k ( a ) = ii a ii li a _ 1l l 很大( 坏) 时( i l l - c o n d i t i o n d e d ) ，a + a 的条件数 k ( a a ) = k ( a a * ) = 胪( a ) 变得更坏，更增加了计算少数几个最小的奇异值 的难度( 见f 2 7 】) 如果用c 来计算a 的奇异值时，由式( 1 5 ) 知需求c 的 内部特征值，这样不仅会出现无规律的收敛现象，而且必须求对特征值 士吼，这就使计算迭代步增倍，同时要存储m + n 维奇异向量，存储空间增 大( 见【2 4 】) 1 2 奇异值及奇异向量的逼近 为了解决e 述问题，本文给出一种计算最小奇异三元组的隐式重新启 动精化l a n c z o s 双对角化方法，这是一种投影方法，即把a 投影到指定的 低维空间，再通过求解低维问题来求出a 的奇异三元组的方法先讨论 l a n c z o e 双对角化方法( 见【1 0 】) ； 对于矩阵a ，以初始酉向量t ，。c n ，k 步l a n c z o s 双对角化后满足如下 关系： a r k = u k + 1 鼠( 1 6 ) a 仉+ l = k 磁+ q 七+ 1 t + l e ：+ l = k 垅+ 弦e ：+ l ( 1 7 ) 第一幸引言 其中k c m 砧，嵋v k = k ，v k e l = i 1 ，呱l r j 七+ 1 = 厶，住= 鲰+ 1 + 1 鼠：= q 1 0 岛毗 风铂 o 风+ l r c 知+ 1 ) ( 1 8 ) 取是下双对角阵，k r 一椭为单位阵，e j 表示单位基向量， 向量我们称式( 1 6 ) 一( 1 7 ) 是a 的个部分l a n c z o s 双对角化 为残 由g o l u b - k a h a ns v d 算法( 见2 2 ) 可以很容易求出取的奇异三元 组 西凤，谚巩，弓仇 ，1 j i 七，其中 矿哆叮0 我们用它求出a 的奇异三元组的个逼近值 砖m ，妒，a ) 事实上，我们有： 黧b k ( b 仇k ) ，二然，l 娜七 ( 1 9 ) i 磺护= l 取t l 风一j 一” 、一 由左、右奇异向量构成的酉阵记为 v ：= p ，钍尹m ，u 】 巩) ：= 妒，妒，v ( b k 】 令 f 砖舢：= 蠢仇 矽：= 酞+ 1 铲，l sjs 老 i 移) ：= k 泸 ( 1 1 0 ) ( 1 1 1 ) ( 1 1 2 ) 由式( 1 。6 ) ，( 1 。7 ) ，( 1 9 ) 得 k - ( a 罢嚣懋+ r k e l u 5 刈吲s 七 n ia = 毒脚矛) 风工二j 二“ p “叫 3 第一幸引言4 当r 。, u ( ，b h ) 足够小时，可认为 孑5 m ，t - ( f a ) ，矛；枷) 是a 的奇异三元组 当m ，n 很大即a 是大规模矩阵时，部分1 正1 i i c 7 0 6 双对角化的所需的 存储空间要很大，当我们取k 较小时，虽然解决存储量的问题，但是在这 种情况下所得到的凤的奇异三元组就不定能很好地逼近a 的奇异三 元组此时，初始酉向量 0 1 的选择就变得十分重要了，最理想的情况就 是我们取初始酉向量掣，为所求的奇异值巨厅对应奇异向量的线性组合重 新启动部分l a n e z 0 8 双对角化方法有很多，他们不同在于选择不同的初始 酉向量口t 和不同的重新启动策略 最早s o r e n s e n 在文献f 2 5 】提出了种有效的隐式重新启动策略并应用 到a r n o l d i 方法和l a n c z o s 三对角化方法中，这种方法可看作是一种截断 q r 算法( c u r t a i l e dq a ) ，后来历5 疏等人把隐式重新启动策略应用到l b d 并得到相应的递推公式( 见【2 】) ，l a r s e n 在文献 16 】中对l b d 中的基进行 了部分重新正交化过程( p 哦i a lr e o r t h o g o n a l i z a t i o n ) ，接着贾仲孝通过添加 “精化向量 提出了精化算法( 见【1 3 1 ) ，最近k o k i o p o u l o u 等人应用l b d 的递推公式并添加“精化向量”来求解些奇异值及对应的向量 本文的主要工作是给出一种计算少数几个最小奇异三元组的隐式重新 启动精化l a n c z o s 双对角方法，与上述方法区别在于：我们求的是小的奇 异值，因此我们用调和r i t z 值作为位移，能有效地逼近大规模矩阵的最 小奇异值的奇异三元组在算法中我们还用精化残量，精化奇异向量和精 化r a y l e i g h 商，同时采取压缩技术压缩掉已经求出的小的奇异三元组最 后，我们用数值实验表明，我们的算法可以成功地求解大规模矩阵的小的 奇异三元组，收敛速度也较快 本文的结构和内容如下：在第二章，我们首先介绍l b d ，然后给出隐 式重启l b d 及讨论位移的选择问题在第三章中我们证明正交压缩变换 ( o d t ) 的几个结论，给出精化向量，精化残量及精化r a y l e i g h 商的求解 第四章，我们给出完整的算法，最后一章是数值实验 g ：- 章l b d 与隐式重新启动l b d5 第二章l b d 与隐式重新启动l b d 2 1l a n c z o s 双对角化( l b d ) 在算法l b d 中，生成的基向量l = f 让l ，t 2 ，让詹+ l 】，磙= f 1 ) 11 ) 2 ，1 ) k 】 会失去正交性，所以需对基向量重新正交l a r s e n 在文献f 7 ，8 1 中给出了 种对l b d 部分重新正交化方法算法如下，其中r 相当于v k ；仉相 当于巩；迭代步k 比仇，靠小很多；为e u c l i d e a n 范数 算法1 l a n c z o s 双对角化( 部分重新正交) 输入：a c m 炳或关于a 与的向量乘积函数； 烈c n 。初始酉向量；k ：迭代步 输出：最= 胁，沈，p k 】c n x k ：列正交矩阵； 瓴= q l ，他，硝c m x 七：列正交矩阵 反r k + 1 x 知：元素为与岛为下双对角阵；c n ：残向量 1 只：= p l ，吼：- - a p l ； 2 o h l ：= i i q l l ；q l ：= q l a l ；0 , 1 ：= q l ； 3 f o rj i = 1 ：k ( 3 1 ) r j ：= a t 劬一a # p j ； 重晰正交： 吩：= 巧一弓( 乎吩) ； ( 3 2 ) i f 歹 k ，t h e n 岛：= i i r j l l ；1 0 # + l ：= , j 岛；o p j + l ：- - 【弓，乃+ 1 】； 彩+ l ：= 锄+ l 一伤田； 重新正交：q j + l ：= q j + l q j ( q j t + l 留+ 1 ) ； a o + l = i i q j + , l l ；劬+ l ：- - - 毋+ l 哟+ l ；劬+ l ：= 【q ，劬+ 1 】； e n di f 4 e n df o r 第二章l b d 与隐式重新启动l b d 6 2 2 隐式重新启动l b d s o r e n s e n 在a r n o l d i 与l a n c z o s 迭代中提出了一种隐式重新启动方法一 a r p a c k ，这是种求解大规模矩阵特征值问题的有效的方法，后来c u l - h m 在文献【6 】研究了把重新启动策略应用到l b d 由于对a 进行l b d 等 价于对从进行l a n c z o s 迭代( 见【9 ，s e c t i o n9 3 2 】) ，因此z = 七+ p 步的 l b d 后，我们选取p 个位移，对三对角阵t t = 鼠犀作p 次的隐式q r 迭 代为了稳定性，我们直接对双对角阵鼠进行g o l u b - k a h a ns v d 算法即 b u l g e c h a s i n g 策略，算法如下：第个g i v e n s 旋转变换( 第4 行) 产生一 个“b u l g e ( 即次对角的非零元) ，为了保持匕双对角结构，再次用g i v e n s 旋 转变换“c h a s e 该冒泡( b u l g e ) ，这样我们得到满足 的双对角阵砑，其中q l r 1 ) 竹n ，q r r t 都是由一系列g i v e n s 旋 转变换得到的正交( 酉) 阵这样我们有 a k = 巩+ - 玩一a 盼= 奶矗。磁 其中时= v t q r ( 1 ：七，：) ，嘻1 = u t + i q l ( i ：k + 1 ，：) + 这种更新的l b d 分 解也可以由初始向量为 缸i = ( a a 。一p 2 j ) 牡1 作k 步l b d 得到，其中p 位移 算法2 ( b u l g e c h a s i n g ) g o l u b - k a h a ns v d 算法( 见【9 ，s e c t i o n8 6 2 1 ) 输入s 三对角矩阵乃= 玩磁；隐式位移步p 输出；新的匕双对角阵砑 1 令y = t l 。l p ，z = t l 。2 2 f o ri = 1 ：一l ： 3 求c = 础尚刮既使之栽引ic s 三) 邓棚 第二章l b d 与隐式重新启动l b d7 4 对玩作右g i v e n s 旋转变换：岛= b t g ( 1 ，i + l ，0 ) 5 令掣= 6 l 名= 6 + l 。i 6 拈删一柙臌龊( - - c $ 三t 咖，= 7 对玩作左g i v e n s 黼变换： b t = g ( i ，i + 1 ，鳓t & 8 i f i 一1 ，t h e n 矽= 玩，件l ，z = 良，件2 以地，肛3 ，脚为位移进行p 一1 步的g o l u b - k a h a ns v d 后，更新的 l b d 可以由初始向量为 p u ：- = ( 舭一p ；咖- = l 七步的l b d 得到 考虑隐式重新启动l b d 的位移选择问题：( 1 ) 精确的r i t z 值；( 2 ) 精 确的调和r i t z 值；( 3 ) c h e b y s h e v 多项式的零点与l e j a 点( l e j a 位移) 本文 只考虑前两者，由于求解的是小奇异值问题，所以每次隐式重启我们取大 的r i t z 值及调和r i t z 值作为位移步 第二章l b d 与隐式重新启动l b d 8 2 3 1r i t z 值的计算 定义1 ：如果是c ，i 上的个线性空间，称( o k ，钆詹) 是矩阵a 的r i t z 对，如果 a u 锹k 一v 以k , ：三 e 2 j ， 由定义得如果子空间k 的基为k = p - ，t ) 2 ，】，则( 以，u 鼍) 为t = w a k 的特征对 现给出求a 的r i t z 值的方法，进行z = 岛+ p 的l b d 迭代后，由 式( 1 6 ) ，( 1 7 ) 得：、 考虑等式两边前列有 其中 鼠= 叼a k 毋+ q 件l w a 现+ l e + l u ；u t + i 鼠磁+ 叭1 叼a 饥l 咯l ( 五，o ) 玩磁+ o r t + l 磁a t l e 知l w a 阮= & 毋， q 1 0 侥a 2 风q 3 角一l 啦 r t 妣 因此，伊2 ( 度) 是a 小的r i t z 值，我们可以用盯( 反) 来逼近盯( a ) r i t z 策略：取前p 个值盯 ( 岛) ，砖( 度) ，( 反) 作为位移步作p 次的 “b u l g e c h a s i n g ，其中砰( 反) 酲( 岛) 砰( 反) 0 需注意的是：为了 达到更好的精度，应该直接计算岛的奇异值，而避免求岛廓的特征值， 这样计算时不仅可以避免条件数变坏，而且也可以避免小的奇异值( 1 ) 更加聚合 第二幸l b d 与隐式重新启动l b d9 2 3 2 调和r i t z 值的计算 r i t z 值提供了个很好的位移策略，但是当a 的最小奇异值非常密集 ( 聚合) 时，用r i t z 值将会大大减慢隐式重新启动l a l i c z 0 6 双对角的收敛 速度为了解决这个问题，p a i g e 与p a r l e t t 等人在文献f 2 1 】引进了调和 r i t z 值的概念，后来s l e i j p e n 与v a nd e rv o r s t 在文献【2 2 】用调和r i t z 值作 为位移步，提出了一种新的位移策略：一种通过计算( a ) 以的最大r i t z 值来逼近a 的最小奇异值的方法当然在实际计算时，需要避免直接计 算( a ) 以的先给出调和r i t z 值的概念，以下都认为a 是满秩的 定义2 ：在线性空间m 。上，称瓦是矩阵g c m m 的调和r i t z 值， 如果霹1 是g 。的r i t z 值 定理1 ：令k 维线性子空间u 的基为t ，l ，v 2 ，在子空间m = a 上，( 瓯，诹) 是a 的调和r i t z 对等价于 k 二艺兰u 仁2 ， 特别的，当w 。，w 2 ，w k 张成k 维线性子空间眦= a 魄 令 l k = v l ，t j 2 ，叫； w k = l ，w 2 ，饥】； l 凰= ( 眦v k ) - 1 峨a u ； 则( 瓦，诹) 是a 的调和r i t z 对等价于 f s c k , 8 o 反s = 瓦s ( 2 3 ) l - 讥： 证明：由式2 1 在空间a v k 上，( 露1 ，舭七) 是a 一1 的r i t z 对铮( a 一一 嚣1 j ) a 诹上a 其中t 七，讥0 ， 由( a 一1 一露1 ，) a 豇七= 一万一1 ( a 雹七一百一- 豆奄) 定理第1 郭分得证 第二幸l b d 与隐式重新启动l b d1 0 等价于 等价于 即 是a 的调和r i t z 对等价于 讯0 u a 云七一以霞七上a 8 0 a k s 一9 k s 上t w ；a v i , s o k ( 畎v k ) s = 0 a s = 以s ( 2 4 ) 本文我们用斜投影方法和对应的p e t r o v - g a l e r k i n 条件求似的调和 r i t z 值( 见f 2 2 】) 此时搜索空间我们记为u 4 l ，验证空间记为吼+ l = a a + 阮+ l 时，则相应的p e t r o v - g a l e r k i n 条件是 a a + f i t + , 一巩+ l 矗“l 上a a + u t + l 其中纵1 为a 的调和r i t z 值，而且如果u 4 。，w t + 1 分别张成子空间 己“l ，m + l ，那么a a 的调和r i t z 值就是矩阵 或+ l = ( 咳1 v t + 1 ) q 瞻1 a a 。v t + l 的特征值下面定理表明只要多进行步的l b d 迭代计算，就可以直接 求出从+ 调和r i t z 值 定理2 ：在子空间m + l 上，矩阵a 小的调和r i t z 值莎及调和r i t z 向 量雪对应于特征值问题 b h l b t + l y = p y ，3 = b h l 痧 其中玩+ 1 r 1 ) 1 ) 为下双对角阵 证明：由式( 1 6 ) ，( 1 7 ) 得 从+ l = a y t + l 或+ l = 阢+ 2 岛+ l 反+ l ，il、，【，i-li、 第二章l b d 与隐式重新启动l b d 1 1 其中目“。酞2 ) 1 ) 为下双对角阵，岛+ l 也为- f x 叉x 寸角阵且可由玩+ l 去掉最后行得到即 鼠+ 1 筝 口l 0 岛q 2 岛钧 0 忍a t + l 酞1 ) 由文献【2 2 】，调和r i t z 对可由如下广义特征值问题求出： 珥1 珥1 雪= 口啄1 痧 其中丸。可由7 h 。的去掉最后行得到即 l 珥l 上辑l 上狐l 雪“l 雪= 舀反+ l 雪c + 1 雪 设岛+ 。非奇异，令箩= 直m 雪有： l 鼠+ l 可= 而 因此，求调和r i t z 值等价于求b 瓠。b t + t 的特征值，这可由反+ - 的奇 异求得而后者只需要进行z + 1 步的l b d 即可得到 和r i t z 值策略样，调和r i t z 值策略就是每次重新启动我们选择p 个大的调和r i t z 值作为位移来进行“b u l g e c h a s i n g 后面的数值实验表明 在求小的奇异值，用调和r i t z 值策略要优于r i t z 值策略 m 伊 “磁 = 巩 义定 第三章压缩变换与精化向量的逼近 1 2 第三章压缩变换与精化向量的逼近 3 1 正交压缩变换o d t 在隐式重新启动a r n o l d i 算法中，在计算重或聚合特征值时，为了提 高算法的收敛性与稳定性，常常采用压缩技术，这技术已经有效地应用到 块方法中( b l o c km e t h o d s ) 最近，结合块方法，b a g l a m a 等人提出了计算 某些特征对或奇异三元组的一种隐式重新启动a m o l d i 方法( 见【3 】) 本节 我们直接对玩进行正交压缩变换( o r t h o g o n a ld e f l a t i n gt r a n s f o r m a t i o n ) ( 见 【1 4 ， 2 6 1 ) 这样，我们可以在隐式重新启动l b d 中压缩掉所求的鼠奇异 值及对应的奇异向量，这种压缩技术最重要的步是构造正交阵q 引理1 ；令向量旷一h ，7 2 ，啦) ，l l y l l = 1 ，记y ；= h ，啦，协) ，令 ，m 、 乃= l i 殇i l ，2 jsn 定义= i 乃 l ，其中协量享，p j 三警如果 o q 三 q l ，q 2 ，训，q l = y ，那么：q + q = i ，q e l = y ，e 蚕q = 陬，一l e 。t l 】 证明：首先，有q ；= q i y r i ，i 歹，q ：y = y l y i 。l 依+ 碾风= 哿( 乏备) + 依警= 0 其中2 i n 所以= 0 ，而茁= 1 因为噶q = 霄难1 y 一l + 霹= 霄( 豪警) 2 + - 嘻w 1 - , = 竺笋= l ，其中2 曼i 他得证 由引理得到的q 有如下性质： q = r + y e ；，其中r e l = 0 ，r + y = 0 且r 是e 三角阵， q = l + y g + ，且r e l = 0 ，r 。y = e l g ，其中l 下三角g 。兰e + e ；r 去 定理3 ：进行k 步l b d 迭代后，可求出t y x 寸角阵b 的奇异三元组 ( 0 ，y l ，y r ) ，该三元组可作为a c ：，i 黼的奇异三元组的个逼近分别对向 量y l ，y r 作o d t 变换得到酉阵钆= 钆( 扎) c ( 七+ 1 ) ( 詹+ n ，锄= q r ( 期) 秽詹那么直= 砚b q r = ( 三台0 ) ，其中口为逼近的奇异值，雪，言均下 双对角阵 证明：要证言= q * 工b q r 为- v y x 寸角阵，只需证言既是上h e s s e n b e r g 第三幸压缩变换与精化向量的逼近1 3 阵，又是下三角阵即可 事实上有b = q * l b q = q b ( r r + y r e ：) ；q * l b r r + q 2 ( b ) e ：霉 q 2 b r r + q * l o y t e ；= q * l b r r + 8 e t e ；= ( l l + g l y * l ) b r a + o e l e ：一l 2 b r r + g l ( y * l b ) r r + o e l e ；= l * l b r r + 吼( y 夤冗冗) + o e l e ：= l l b r r + 口e l e ：其中三2 上三角阵，b 下双对角阵，上三角阵，所以啻为上h e s s e n b e r g 阵 现证明直为下三角阵事实上直窖q * l b q r q * l b ( l r + y 触盖) = q 2 b l r + q 2 ( 上珂r ) 鳋一q 2 b l r + p q 2 y l 鳙拳q * l b l a + 口e 1 踮一( 晚+ y l e ：) + b l a + o e l g 盖= ( 吃+ e l y 2 ) b l r 十p e l 鲧一r * l b l r + e t ( y * l b ) l s + o e , g ；z = r 2 b l r + o e t y n l r + o e x 鳋= r * l b l a + o e x ( y * r l n + 鳙) = r * l b l a + o e t e ：，其 中磁，k 均为下三角，b 下双对角阵，所以直为下三角阵得证 3 2 精化奇异向量的求解 在很多情况下，r i t z 向量收敛比r i t z 值收敛要困难得多，特别是矩 阵非h e r m i t i a a 时，为了解决这个问题，贾仲孝在文献 1 2 ，1 3 】提出了精化 方法，这是一种复合正交投影类方法 设( 凡，协) 为a 的特征对，求得a 的r i t z 对为( 以，地) ，对于某个r i t z 值巩，寻找单位长向量啦e ，使得它满足最优性质 i i ( a - 训2 爵码阻i i ( a o , z ) 札l l u e el u l 。( 3 1 ) ，ii = l 并用之逼近忱 称啦为a 在e 中对应于凡的精化近似特征向量显然 i i ( a 一侠，) 魄l l i i ( a o d ) 啦l l 所以魄至少同啦样好，诹可视为忱关于巩和e u c l i d 范数在e 上的 最佳逼近 由式( 2 1 ) ，( 3 1 ) 得e 上如下的精化正交投影类方法： f t i t ，哦印 ( a 地一岛啦) 上e ( 3 2 ) 【i l ( a 一巩啦t i l _ u 高愀一以啦i l 第三章压缩变换与精化向量的逼近 1 4 它们用以，诹逼近九， 这里我仅仅给出精化向量的求解方法我们知道l b d 分解等价予对 c = ( 二 ) 进行l 一迭代c 她s 砒以3 硼，因此可以用。 来计算精化残量( i 句j i t ) 具体过程如下：设矛为a 的最小奇异值的个 逼近，通过求下列极值问题得到精化奇异向量矗= v l + 。s 肠+ ，( 俨a ) 与 矛= v , t k ；t ( a t a ) s 忒n 叫u l + 1 暑) l | 2 1 1 6 ；t 1 1 2 = 1 矗忒倒“ ( 乏+ 。? 卜( 1 涠1 2 = 1 洲：) 1 1 2 气蠹竹蚍赢0 履) o o t 2 s i 烈* 等1 ) ( ：卜 其中 r 凹+ 2 = ( 毗曼。芋) 一矛( 如0 件1 ) 当 s ，t 】取最小奇异值n ( 勘+ z ) 所对应的右奇异向量时， gs 章压缩变换与精化向量的逼近 1 5 享1 一子( 如0 件1 ) n 2 取得最小，称这样的残量为精化残量，记为r 3 3 精化向量与精化r a y l e i g h 商 由上面方法我们可求出豇= 阮+ l s 杨+ l ( a t a ) 与雷一坳 c t ( a t a ) ， 般情况下，精化r i t z 向量比r i t z 向量能更好地逼近实际的特征( 奇异) 向量，因此为了更有效地得到最小奇异值，数值实验中我们还用精化r i t z 向量作为重新启动向量 而且由于p 比矛更好的逼近真实的奇异值( 见 2 7 ，s e c t i o n4 3 】) 因此 每步新启动时，如果0 p 子嘲( a ) 用p = 矿勘代替子翻。似) 数值实 验结果表明对于某些矩阵，这样可以更有效地得到最小奇异值 o 曰粥啦 ，fiiiiiiii一 第四章隐式重启精化l a n c r 沁s 双对角方法1 6 第四章隐式重启精化l a n c o z s 双对角方法 这章我们给出最终算法一算法3 ，其中= 七+ p 为l a n c o z s 双对角化 迭代步数；p 为隐式位移q r 迭代步数；e i g n u m 为所求的最小奇异值的 个数；t 0 1 为收敛准则t o l e r a n c e 算法3 中第步我们直接用l a r s e n 的p r o p a c k 1 8 ，1 9 里的l a n b p r o 得到l b d 分解，其中l a n b p r o 是l a n c z o s 双对角化的种部分重新正交化 方法接着我们考虑位移的选择问题：当我们用r i t z 值时，为了稳定性，我 们直接计算玩的奇异值而避免计算尻毋的特征值，而用调和r i t z 值时， 我们只需增加步的l b d 迭代求出国h 。，再求它的奇异值即可然后我 们计算出口椭。及相应的精化残量的参范数州1 2 ，当州l z t 0 1 n o r m e s t ( a ) 时，仃赫为求得的个奇异值逼近值，其中n o r m e s t ( a ) 由第步的鼠的 最大奇异值来估计此时求出盯嘲相应的左( 右) 精化奇异向量，用o d t 由左( 右) 精化奇异向量产生正交阵进行压缩( 见3 1 ) ，最后重新正交并 用l b d 把最变为鼠，如此循环，直到求出e i g n u m 个奇异值算法如下 算法3 计算最小奇异三元组的隐式重启精化l a n c z o s 双对角方法 输入：矩阵a 沪黼，初始向量u l ；七；p ；e i g n u m ；t 0 1 ；z = 七+ p 输出：e i g n u m 个最小奇异三元组 1 用l b d 计算基1 ，k 及双对角阵段 2 r e p e a t 3 位移步的选择： ( 3 1 ) i f ( s h i f t s = = r i t z ) 计算岛奇异值：吼，i = 1 ，z ( 3 2 ) e l s ei f ( s h i f t s = = h a r m o n i c ) 计算风l 及玩+ 1 的奇异值以，i = 1 ，f + 1 ( 3 3 ) e n di f 4 隐式重新启动l b d ： 第四章隐式重启精化l a l c j o s 双对角方法1 7 ( 4 1 ) 由3 求出p 个最大的( 调和) 础t z 值以，并求砰，i 毒l ，2 ，p ( 4 2 ) 以砰为位移步作p 次隐式q r 迭代于最得到a 昭= 吃。w 5 计算小的奇异值及相应的精化残量： ( 5 1 ) 纠算( a ) = m i n 以) 及对应的p = f i , a 矛 ( 5 2 ) 计算n ( a ) = m i n ( s r m i n ( a ) ，p ) 的精化残量- 6 i f l t o 幸n o r 胱s t ( a ) t h e n ( 6 1 ) 计算讥的左，右的精化奇异向量 ( 6 2 ) 计算吼，q r ，应用o d t 进行压缩 ( 6 3 ) 压缩掉阢柏k 的第列，岛的第行，第列 ( 6 4 )令七= 后一i ，e i g n u m = e i g n u m - 1 ， 7 重新正交化并作扩充： ( 7 1 ) 重新正交化让玉1 ，硭 ( 7 2 ) 把a 时= 嘻1 彤扩充到长为z = 七+ p 的l b d 分解， 8 直到e i g n u m 个奇异值都收敛 第五章数值实验 1 8 第五章数值实验 在一台a m ds e m p r o n ( t m )