（计算数学专业论文）右端依赖位移线性方程组的krylov子空间方法.pdf

2 0 0 4 年上海大学硕士学位论文 摘要 本文讨论下列右端依赖位移线性方程组的数值方法 + a l i ) x ( ) = b ( a f ) ，z = 1 ，p ( 1 ) 我们利用种子投影思想 3 3 】，提出了两种形式的g m r e s 种子投影方法。我们研 究了算法的残量估计数值实验显示了算法的有效性。 在一些实际领域中，模型( 1 ) 还可以转换成下列多右端位移方程组， ( a + a t i ) x ( ) = f ii = 1 ，p ( 2 ) 其中f c “。q 与巩无关我们运用块q m r 方法【8 】于( 2 ) ，着重讨论了结构力学 中二次参数复对称方程组的处理，得到了简化的块q m r 算法二个实际问题的数 值实验显示了算法具有很高的效益 2 0 0 4 年上海大学硕士学位论文 i i a b s t r a c t w ea r ei n t r e s t e di nt h es o l u t i o no ft h ef o l l o w i n gp a r a m e t e r - d e p e n d e n tl i n e a rs y s t e m s ( a + 印j ) z ( ) = b ( 0 1 ) ，l = 1 ，一，p ( 1 ) a c c o r d i n g t ot h ei d e ao ft h es e e dp r o j e c t i o n 3 3 ，w ep r e s e n tt w of o r m so fg m r e s s e e dp r o - j e c t i o nm e t h o d s t h e r e s i d u a le v a l u a t e so ft h em e t h o d sa r eg i v e n 。n u m e r i c a le x p e r i m e n t s a r er e p o r t e dt oi l l u s t r a t et h e e f f e c i e n c yo f t h em e t h o d s i ns o m ea p p l i c a t i o np r o b l e m s ，t h es y s t e m s ( 1 ) c o u l db et r a n s f e r e dt ot h ef o l l o w i n g s h i f t e ds y s t e m sw i t hm u l t i p l er i g h t h a n ds i d e s ， ( a + o l j ) x ( ”= f i z = l ，p ( 2 ) w h e r ef c “qi si n d e p e n d e n to n 川w e a p p l yt h eb l o c kq m r 8 m e t h o dt o ( 2 ) ，d i s c u s s h o wt oh a n d l et h ec o m p l e xs y m m e t r i cl i n e a rs y s t e m sw i t ht w o - o r d e rp a r a m e t e rf r o mt h e s t r u c t u r a ld y n a m i c sa n d p r e s e n tt h es i m p l i f i e db l o c kq m r m e t h o d t h en u m e r i c a lr e s u l t s o ft w o a p p l i c a t i o np r o b l e m ss h o wt h a tt h em e t h o d o b t a i n sm u c hb e n e f i t 原创性声明 本人声明：所里交的论文是本人在导师指导下进行的研究工作。 除了文中特别加以标注和致谢的地方外，论文中不包含其他人已发表 或撰写过的研究成果。参与同一工作的其他同志对本研究所做的任何 贡献均已在论文中作了明确的说明并表示了谢意。 签名：基垫达日期超兰：! 本论文使用授权说明 本人完全了解上海大学有关保留、使用学位论文的规定，即；学 校有j 汉保留论文及送交论文复印件，允许论文被查阅和借阅；学校可 以公布论文的全部或部分内容。 ( 保密的论文在解密后应遵守此规定) 签名：墓二起数导师签名：刁啵日期： 毒。r 2 0 0 4 年上海大学硕士学位论文 第一章绪论 1 1 问题介绍 我们考虑下列多组非对称线性方程组的求解 a ( o o ( f ) ：b ( 0 其中 ( ) 是n n 非奇矩阵，且一般a ( ) a ，6 屯，i j 由于n 很大。使用 迭代法是合适的，如k r y l o v 子空间方法。 如果各系数矩阵a ( ) 和右端向量6 ( ) 完全是任意的，一般没有希望找到有效 的方法去同时求解( 1 1 ) 此时我们不得不逐个独立地求解( i i ) 的p 个方程组幸 运地是许多实际应用问题归结出来的各系数矩阵4 ( d 和右端向量彬) 都有某种关 联，或者说在某种程度上都有一些共同的信息可分享，例如渡散射问题【1 , 8 ，1 2 ， 积分方程的数值方法【2 2 ，最小二乘问题的递归计算f 2 5 和图像恢复 2 0 】等等， 这些问题导出的矩阵或是位移矩阵a = a + 仃l ，或相差一个秩一l 矩阵，或具有相 同的特征向量等等。在这种情况下，有可能建立起有效的方法同时求解( i i ) 。 问题( 1 1 ) 可以分成下列三种情况： 多右端问题在( 1 1 ) 中，当a ( ) ea 时，即成为多右端线性方程组 a x = b e 6 ( 1 ) ，纠 ( 1 2 ) 这类问题在许多领域有重要应用，如结构力学、控制论和电磁场多方向入射点的 散射计算等等【4 ，3 3 。对于多右端问题( 1 2 ) ，目前主要有两大类方法，块迭代法 和种子投影方法。 块迭代法是求解( 1 2 ) 的自然选择，对其的研究始于1 9 8 0 年0 l e a r y 2 4 的 块c g 方法( b c g ) ；之后有块g m r e s 方法( b g m r e s ) 2 9 ，3 0 】及其变形【1 6 ，3 l 】，块 q m r 方法( b q m r ) 【8 ，块e n 方法( b e n ) j 1 7 等等这些块迭代法一般都具有较好 2 0 0 4 年上海大学硕士学位论文 2 的收敛性态。目前对块算法的研究还在继续中 3 5 】，特别是应用于矩阵方程的求 解【1 3 ，1 9 ，2 1 ，2 6 1 然而随着块尺寸的增大，块算法的计算量和存储量将迅速增大；其次块算法 还须处理线性相关性问题【8 ，1 7 ，3 4 在1 9 8 9 年，s m i t h 等【3 3 】提出了求解( 1 2 ) 的 另一类比较有效的方法，称之为种子( s e e d ) 投影方法当方程组的各右端在某种 程度上比较“接近”时，这种方法特别有效种子投影方法的基本思想是：首先 按某种方式构造( 或选取) 一个方程组作为种子方程组，对该种子方程组使用某 种k r y t o v ( 型) 子空间方法进行求解，同时形成了一个k r y l o v ( 型) 子空间，随 后把非种子方程组的当前残量正交投影到上，以获得一个更好的近似一旦种 子方程组收敛到精度时，则在余下还没有收敛的非种子方程组中，再构造一个新 的种子方程组，重复上述过程；直到获得所有方程组的解在【l9 中，针对c g 方 法，c h a r t 和w a n 分析了种子投影方法的超收敛性态，即把投影近似解作为迭代 初值时，其收敛性比通常的c g 法要快；特别是当各右端比较“接近”时。 位移方程组问题在( 1 1 ) 中，当a ( 2 ) = a + 田f ，b ；6 ( 1 ) ，即成为位移方程组 ( a 一印，) 。 】= b ，f = 1 ，p ( 1 3 ) 在利用高阶隐式方法求解p d e 问题【3 4 ，控制论 5 】，q c d 问题( q u a n t u m c h r o - m o d y n a m i c s ) 1 2 】以及结构动力学【3 2 中，都会遇到位移方程组( 1 3 ) 的求解由于 需要求解上百个位移方程组，如何建立有效的数值方法是一个感兴趣的问题。我 们知道，k r y l o v 子空间 苊m ( a ，口) = s p a n v ，a v ，一，a ”一1 ， ( 1 4 ) 对于位移具有不变性，即j c 。c a ，7 3 ) = b ：( a c r i ，口) ，因此使用k r y l o v 子空间方法 求解( 1 3 ) 是自然的选择这些方法中包括q m r 方法 7 j ，g m r e s 方法( f u l l o r r e s t a r t i n gv e r s i o n ) 1 1 】及其变形 1 5 ，1 8 ，f o m 方法 2 8 】以及b i c g s t a b 方法 1 0 等 等这些方法的基本思想类似于多右端的种子投影方法【3 3 】，但不是投影过程 取某个位移方程组作为种子方程( 不妨设其为零位移方程) ，即 a x = b ( 1 5 ) 2 0 0 4 年上海大学硕士学位论文 3 用某种k r y l o v 子空间方法求解( 1 5 ) ，同时在其所产生的k r y l o v 子空间。( a ，r 0 ) 中，用很小的计算量产生所有余下位移方程组的近似解，由于k r y l o v 子空间的位 移不变性( 只要位移方程组的初始残量而与r o 保持线性相关) 。在一定条件下， 这些近似解也收敛 众所周知，k r y l o v 子空间方法有长递推和短递推两种，前着是基于l a n c z o s 过 程产生的算法，如c g 方法，q m r 方法等等，后着是基于a r n o l d i 过程产生的算 法，如g m r e s 方法，f o m 方法等等在用短递推方法求解位移方程组( 1 3 ) 时， 会遇到存储量的限制 7 ，1 l 】而长递推方法在计算量和存储量上也受到限制，为此 一般采用重启动形式；为使重启动的初始向量r 0 仍与确保持线性相关，f r o m m e r 【1 1 】修正了重新启动的g m r e s 方法以求解位移方程组为改善其收敛性态，g u 提出了扩展的g m r e s ( 重新启动形式) 方法f 1 5 ，18 】而重新启动的f o m 方法则 自然保持重启动初始向量r o 与确的线性相关性【2 8 】，但后者的收敛性态并不一定 好目前位移方程组的预条件处理是一个感兴趣的、也是具有挑战性的问题( 6 ，3 2 右端依赖位移问题在许多应用问题中【2 0 ，2 5 ，3 2 】，( 1 1 ) 中的系数矩阵a ( ) 具 有位移形式a = a + a t i ，即 ( a + a t i ) x ( 。) = 6 ( 2 ) ! b ( a 1 )z = 1 ，p ( 1 6 ) 我们称此为右端依赖位移的线性方程组方程组( 1 6 ) 可等价地写成s y l v e s t e r 方 程， a x + x a = b 其中a = d i a g ( a l ，o p ) ，x = p ( ”，z ( p 】) b = 【6 ( 1 ) ，6 ( p 】对一般的s y l v e s t e r 方程a x + x d = c ，若对矩阵d 进行s c h u r 分解，可转化成问题( 1 6 ) 的求解。 令d = q 兄q r ，其中q 是正交阵，r 是上三角阵。于是有 a x + x q r o t = a - - - - - a x q + x q r = c o 令y = x q ，f = g q ，则a y + y r = f ，即 i i + r u i ) y = ) 一r i t y ( ”， f _ 1 ，p i = l 2 0 0 4 年上海大学硕士学位论文 4 注意到此时右端项依赖于前面z 一1 个解 目前求解一般的多组线性方程组问题( 1 | i ) 和方程组( 1 6 ) 的有效方法不是很 多，研究也不多在 2 中，针对对称正定的系数矩阵a ( f ) ，t c h a n 等建立了c g 型的投影方法以求解多组线性方程组( 1 1 ) 多右端位移方程组在一些实际问题中，方程( 1 6 ) 可转换成下列多右端位移 方程组的求解( 详见第三章) ， ( a + o j i ) x ( j ) = b ，j = l ，p ， ( 1 7 ) 其中b 是n q 阶矩阵，如结构动力学中的频域分析与计算同题【3 2 】。由于块k r y l o v 子空间 。( 4 ，v ) = s t u n v , a v , 一，a “一1 y ) ( 1 8 ) 仍具有位移不变性： 尼m ( a + 盯j ，v ) = k m ( a ，y ) 因此，用块k r y l o v 子空间方法仍是求解方程组( 1 ，7 ) 的自然选择 b g m r e s 方法，b q m r 方法等所谓块k r y l o v 子空间 j c 。( a ，y ) = s p a n v , a t , 一，a “一1 y ) 如b c g 方法， 是所有生成元k a v , ，a m - 1 y 的列向量线性组合全体，其中v r “q 1 2 背景知识 ( 1 9 ) a r n o l d i 过程首先我们介绍重要的生成k r y l o v 子空间 靠( a ，v 1 ) 一组标准正交 基，口。的a r n o l d i 过程【2 7 】，其基本形式是下列算法。 算法1 1a r n o l d i ( m 步) 过程 ( 1 ) 取v l r “，慨0 2 = 1 ( 2 ) 对j = 1 ，m ， ( 3 ) h i j = ( a v j ，q ) ，i = 1 ，- 一，j ， ( 4 ) w j = a v j 一名lh i j v i ， 2 0 0 4 年上海大学硕士学位论文5 ( 5 ) b + z z = f l w j l l 2 ， ( 6 ) 若b + 1 ，j = 0 ，停止 ( 7 ) v i + l = w j h j + 1 ，j ( 8 ) e n d ( j ) 考虑到数值计算的稳定性，在a r n o i d i 过程的实际计算中，一般都用修正的 g r a m - s c h m i d t 正交化过程，既在上述算法的第( 3 ) 行和第( 4 ) 行中，改成 ( 3 ) w j = a q ， ( 4 ) 对i = l ，j ，岣= ( 叼，v i ) ，q = 吣一h i j v i e n d ( i ) 由a r n o l d i 过程，可以得到下列重要的a r n o l d i 分解： a i 么= v m + 1 曰m = 日k + h m + l ，m m + l e 磊，( 1 9 ) 舯耻h ，川，= ( h 属r n ，) 是c 帅的仙s s e 慨蹴 k 十1 。舔， 。 。 g m r e s 方法 g m r e s 方法( g e n e r a l i z e dm i n i m u mr e s i d u a l ，广义最小残量法) 是由s a a d 和 s c h u l t z 在1 9 8 6 年( 2 8 】提出的求解般非对称线性方程组的最有效的方法之一。 考虑k r y l o v 子空间 尼m 兰咒m ( a ，t o ) = s p a n r o ，a r o ，- ，且”一1 r 0 ) ，( 1 1 0 ) 其中r o = b - a x o 初始残向量。令”1 = 啊成p = l i t o 阮由a r n o l d i 过程生成 ( a ，伯) 的一组标准正交基，v r n g m r e s 方法便是在 z o ) + g m ( a ，r o ) 中计算迭代解 。m = x o + d m ，使得其对应的残向量= b a x 。与空间a 。正交，即 ( a ) 丁。r 。= 0 ( 1 1 1 ) 注意到r o = 触= 芦十1 e l 于是正交化条件( 11 1 ) 可以写成：矗三v 量。( r 0 一 a d m ) = 矗- 。tv 。t + l + 1 e l 一+ l 凰d m ) = 矗三( 卢e l 一取) = 0 ，即d ( m ) 是最小 二乘问题的解： d m 2 ”g r a i n 。慷l 一m ( 1 1 2 ) 2 0 0 4 年上海大学硕士学位论文 6 ( 1 1 2 ) 式表明，在空间d o + c 。( a ，r o ) 中，z 。具有残量模最小的这一重要最优性 质，即 m 0 22 n r g 。缸r a 。+ i n 丘。l i b a 。1 1 2 2 o r g d m i n 。l i b 算法1 2 ：g m r e s 方法 取初值蛳 ( 1 ) 计算= b a x o ，声= 怖峙口l = r o 3 ( 2 ) 执行m 步a r n o l d i 过程( 算法1 1 ) ( 3 ) d 。= a r g m i n d e r i i 肚1 一日n d l l 2 ( 4 ) z 。= 2 ：0 + d m a ( 2 ：o + 酬2 2 n r g r a i n 。l i f o 一硎2 ( 1 1 3 ) 注记1 1 在m - 步a r n o l d i 过程中，看遇到h + i d 一0 ，则令m = j ，直接转弟( 3 ) 步，此时得到的。是a 2 ：= b 的精确解，这里因为当h j + 1 ，= 0 时，a r n o l d i 分解 ( 1 9 ) 成为a 巧= 岣岛，而d s 成为冯d = 3 e l 的解。因此b - a z 。= b - a x o - a v j d j = r o 一巧马嘭= v j ( 3 e l h j e s ) = 0 这种现象称之为”l u c k ”中断 注记1 2 在实际计算中，可以通过计算残向量的模l i r m l h 来判断迭代是否结 束我们知道，由于差k 为上h e s s e n b e r g 矩阵，求解最小二乘问题( 1 1 2 ) ，一般用 g i v e n s 旋转变换q ，把豆。化至上三角矩阵。，即 ( ? ) ， 其中p = n m n 1 记9 m + l = p ( f l e l ) = ( 7 h 一，+ 1 ) 丁则 a 。= n r 。m 。r i n 。| | ，。+ ，一( ? ) d f | z c - a ， 于是残向量 r 。= 。一a 。= v 赢+ - c 卢e - 一h m d m ) = u 订+ p t c 。+ - 一( ? ) 商。， 2 0 0 4 年上海大学硕士学位论文 7 l l ：| 鲕+ l f 1d m ) 1 1 2 ：l 川 ( 1 1 5 ) 、0 重启动的g m r e s ( m ) 方法标准的g m r e s 算法( 称之为f u l lg m r e s ) 是一种 长递推算法；当i i l 2 = l + 1 i e 不成立时，在m 一步a r n o l d i 过程的基础上，再 做第( m + 1 ) 步，即生成h m + l 和+ 2 求出如+ 1 及z m + l = x o + 十l d m 十1 显然 这种做法要存储+ 。，而h m + 1 ，”。+ 2 的生成则是一个长递推正交化过程随着m 的增大，算法的存储量和计算量迅速增大 因此在大多数情况下，g m r e s 方法往往采用重新启动( r e s t a r t i n g ) 的形式 【27 ，即当i l 1 1 2 e 不成立时，令z o ：= z 。，重新执行算法1 2 的第( 1 ) 步到第( 4 ) 步，这样算法只需存放k ，而正交过程也至多是n l 步 残量估计在本文中，e ( c ，e ，a ) 表示一个包含某矩阵所有特征值的椭圆，其中 c ，c 千e ，o 分别表示椭圆的中心，焦点和长半轴；并假设c r + ，e 2 r 。 定理1 1 f 2 7 】设a 可对角化，a = x a x ，a = d i a g ( a 1 ) r ，k ) 若包含a ，h 的椭圆e ( c ，e a ) 不包含原点，则1 i 2 步g m r e s 算法的残量满足 。姒7 ( 1 、c 啪m ( a e ) e ) | i i r 0 ( 硼( 高需) ， ( 1 1 6 ) 其中c k ( t ) 是1 t i 阶c h e b y s h e v 多项式 定理1 2 2 7 】若a 是正定矩阵，则对任意的m 1 ，g m r e s ( m ) 收敛。 块l a n c z o s 双正交过程在第三章中，我们讨论用块q m r 方法求解多右端位移 方程组( 1 6 ) ，而块q m r 方法是基于块l a n c z o s 双正交过程 2 7 ，该过程形成了块 k r y l o v 子空间b ；( a ，y ) 和咒。( a t ，w ) 的一组基。 算法1 3 快l a n c z o s 双正交过程过程 ( 1 ) 取u ，肌r “口，满足w 甲h = j ( 2 ) 芦1 = 6 1 = o 口q ，w o = v o = 瓯口 ( 3 ) 对j = 1 ，m ， ( 4 ) 嘶= w 于a 岣， 2 0 0 4 年上海大学硕士学位论文 8 ( 5 ) 巧+ l = a 巧一巧q k l 岛， ( 6 ) 嘭+ ，= 一哆一一可， ( 7 ) 作q r 分解； 嚷1 巧十l = 岛+ - 妨+ ( 8 ) i 巧+ l = 巧+ 1 岛+ 1 ， ( 9 ) 巧+ l - 岣+ l 群1 ( 1 0 )e n d ( j ) 记= m ，】，w m = w 1 ，】 性质1 3 若算法在m 步前不发生中断，则 n 景= f ( 1 t t ) 此外，w 。分别是块k r y l o v 子空间c 。( a ，v ) 和 c m ( a t ，w ) 的基；并且成立 a = + l 磊= v m + + l 碥+ 1 碍，( 1 1 8 ) 妒w m = w 。蠕十+ 1 藤+ 1 碟，( 1 t 9 ) 震a 、k = z 仇，( 1 2 0 ) 其中 = o l 函 如啦 e m = 0 ，】t r m 口。q 风 d mo 仃1 矗m g o m 口 块q m r 方法 8 】设凰c n 。是方程组( 1 2 ) 的初值。在凰+ 尼。( a ，h ) 中 找方程组( 1 2 ) 的近似解 x m = x o + z 。( 1 2 1 ) 其残量r 。= b a x 。= 矗。一a 、k z 。= 凰一+ l 露z 。取是的q r 正交 因子，即 r o = h r ( 1 2 2 ) 嘏 | l 露 2 0 0 4 年上海大学硕士学位论文 9 于是瑞产、k + 1 ( e l r 窖k z 。) 取z m 为下列最小二乘问题的解 2 哪z 毒卿。，l i e i r 一露圳f ( 1 2 3 ) 这便是用来求解多右端方程组( 1 2 ) 的b q m r 方法 1 3 我们的工作 ( 1 ) 基于s m i t h 等 3 3 提出的用以求解多右端方程组( 1 2 ) 的种子投影思想， 针对非对称右端依赖位移问题( 1 6 ) ，我们建立了g m r e s 种子投影方法，算法2 1 和算法2 2 ；特别是算法2 2 ，也可以用来求解多组线性方程组( 1 1 ) 我们研究了 算法的残量性质；并进行了数值实验我们选取不同的右端，不同的位移以及位 移的个数进行试算，结果显示了我们算法的有效性 ( 2 ) 在第三章，我们把块q m r 方法应用于多右端位移方程组( 1 7 ) 的数值求 解对于结构动力学中的二次参数方程组问题，我们把其转换成方程组( 1 7 ) 的形 式，并针对矩阵的复对称形式，得到了简化形式的块q m r 算法；并进行了两个实 际问题的数值计算，得到了很好的效果 本论文的数值实验都是在i n t e l p e n t i u m3 0 0m h z 计算机上进行，而算法的程 序则是用双精度f o r t r a n 语言编制 2 0 0 4 年上海大学硕士学位论文 1 0 第二章g m r e s 种子投影方法 在本章中，我们把s m i t h 等【3 3 1 提出的求解多右端方程组( 1 2 ) 的种子投影思 想，推广到右端依赖位移的线性方程组( 1 6 ) ， ( a + a t i ) x ( 。) = 6 ( ) eb ( a 1 )f = 1 ，一，p ， 我们提出了两种形式的g m r e s 种子投影方法，算法2 , 1 和算法2 ，2 ；特别是算法 2 2 ，也可以用来求解多组线性方程组( 1 1 ) 我们研究了算法的残量性质，并进行 了数值实验。相对于单独逐个求解个位移方程组，数值实验显示我们的投影算法 是有效的 2 1 两个算法 设r ( 。) 是种子方程组的初始残向量，由a r n o l d i 过程产生k r y l o v 子空间k 。i 芘。( 也，r c 。) ) 的一组正交基”l ，”。，其中”1 = r 和) f i r ( 3 ) 1 1 2 。记= p 1 ， 。】，于 是成立 厶= + 1 e k ( 2 1 ) 种子方程组的g m r e s 近似解是i ( 。) = z ( s ) + 扣) ，其中d ( s ) 满足 d 5 = n r g d r a i n | | 唿l r 3 一厩吼20 r g d r a r m i n 慨l 一 m d l l 2 ， ( 2 - 2 ) 卢= l i t ( 5 ) 怯 考虑非种子方程组的近似解设x 0 ) 是第j 个非种子方程组的初始近似解。我 们在空间如u ) + 肥。( 也，r ( s ) ) 中计算新的近似解：2 0 ) = x ( 2 + v m d ( j ) 由种子投影 的思想，使孟( j ) 的残量e l i ) = 6 ( j ) 一向圣o ) 0 s ) 在空间s p a n v m + l = k 。+ 1 ( a 。，r ( 5 ) ) 中的正交投影模i i + l 唿1 e ( j ) i h = l 嘿l f ( ，1 1 2 达到最小注意到a j = a + 勺j = a 。+ ( 唧一口。) j ，故成立 如吐帕，惦砜讹t ( ：) t 马坤渤 2 0 0 4 年上海大学硕士学位论文 1 1 其中岛= + c q 一以，( ：) 于是有+ ，唿。朋= “嗡。扣l 马饥 即d ( ，) 满足下列最小二乘问题， d r g r a i n 。0 嚆1 r 一岛北 ( 2 4 ) 如此我们可以得到求解右端依赖位移的位移方程组( 1 6 ) 的g m r e s 种子投影 方法。记作g s s m a ( i ) + 算法2 1 ( g s s h l f t ( i ) 算法) ： ( 1 ) 取初值z ( ) = 0 ，r ( ) = 删，i = l ，p ( 2 ) 对l = 1 ，p ，( 求解所有方程组的解) ( 3 ) 令s = f ，( 即取当前第z 个方程组作为种子方程组) ( 4 ) r ( s ) = “，) 一a 。( “，p = i i t ( 。) 1 1 2 ， ( 5 ) 对k = 1 ，( 直到求出该种子方程组的解) ( 6 ) ”1 = r ( 。) 卢， ( 7 ) i l l - 步a r n o l d i ( f l 。， 1 ) 过程，形成ta 3 = 十l 豆k ( 8 ) d ( 引= a r g m i n d 舻l | 肚l 一日n 训2 ( 9 ) 雪( 。) = z ( s ) + d ( “， ( 1 0 ) f ( $ ) = 删一a 。牙( s ) ( o ri ( s ) = r ( 5 ) 一+ l 凰d ( 5 ) ) ， ( 1 1 ) 对j = f + 1 ，p ，( 求解所有非种子方程组的投影解) ( 1 2 ) d ( j ) = a r gm i n d 舻i i v y + l r ( j ) 一马u l l 2 ( 1 3 ) 面0 ) = x 0 ) + d 0 ) ， ( 1 4 ) 如) = 删一a j 孟o ( o r 删= r ( ，) 一+ 1 岛棚) ， ( 1 5 ) 若舻( j ) i 2 t ，删除该已收敛的方程组 ( 1 6 ) e n d f i ) ( 1 7 ) 若峭3 1 1 2 e ，转( 1 9 ) ( 1 8 ) e n d ( k ) ( 1 9 ) e n d ( 1 ) 2 0 0 4 年上海大学硕士学位论文 1 2 把上述算法用于问题( 1 1 ) 的求解，非种子方程组的残量 f 0 ) = b u ) 一a ( j ) 孟u ) = r u ) 一a ( i ) v m d o ) 由于a ( d v 并没有a r n o l d i 分解，因而雪( j ) 的残量模极小问题 m 。i ni l r o ) 一a ( i ) v m d l 2 d e r m 1 。 或其正交投影模的极小问题 d r a 胛i n 1 1 唿l ( r 力一v m d ) 1 1 2 的计算都是不现实的由于种子方程组成立a r n o l d i 分解a ( 5 ) v = ，m + 1 鼠。，我们 定义非种子方程组的伪残量；r s ( j ) = 6 ( j ) 一a ( a ) 。( 力这样我们可以得到舟( j ) = r s ( i ) 一a o ) v m d u ) = r s ( j ) 一+ l 豆。d 我们仍在空间 z u ) + ( a ( “，r ( 3 ) ) 计算 非种子方程组的近似解：孟( ，) = 。+ v m d ( j ) ，陡其对应的伪残量内o ) 在空间 。+ 1 ( a ( “，r ( 。) 的正交投影模i i + 1 唿1 厅( 钏2 达到最小，从而 d k 口r g m i n 。1 1 唿1 r s 一岛讹 ( 2 5 ) 这样我们得到求解多组线性方程组( 1 1 ) 或右端依赖位移的位移方程组( 1 6 ) 的 g m r e s 种子投影方法，记作g s s h i f t ( i f ) 算法2 2 ( g s s h i f t ( i i ) 算法) ： ( 1 ) 取初值$ ( ) = 0 ，f = 1 ，p r s ( i ) = 删，i = 2 ，p ( 2 ) 对z = 1 ，p ，( 求解所有方程组的解) ( 3 ) 令s = t ，( 即取当前第f 个方程组作为种子方程组) ( 4 )r ( s ) = 6 ( s ) 一a ( 8 ) 。( “，芦= i i t ( 5 ) 1 1 2 ， ( 5 ) 对k = 1 ，( 直到求出该种子方程组的解) ( 6 ) v l = r 卢， ( 7 ) i n - 步a r n o l d i ( a ( “，u 1 ) 过程，形成：a ( s ) = + 1 i l m ( 8 ) d = a r gm i n d e r mi i 肛l z k r i l l 2 ( 9 ) f f ：c s ) = 。( s ) + v r m d ( “， ( 1 0 ) f ( s ) = 6 ( s ) 一a ( 5 ) 孟( ，) ( o rf ( 5 ) = r ( 3 ) 一l ，m + 1 砧d ( 8 ) ) ， 2 0 0 4 年上海大学硕士学位论文 1 3 ( 1 1 ) 对j = l + 1 ，p ，( 求解所有非种子方程组的投影解) ( 1 2 ) d ( j ) = a r g m i n d r m | f 嘿1 r s ( j ) 一国- m 吼 ( 1 3 ) 毒( j ) = 。仃) + v m d o ) ， ( 1 4 ) ，- 8 ( ) = 6 ( j ) 一a ( 8 ) 雪( ) ( o rv s ( j ) = r s o ) 一v m + i d 0 ) ) ， ( 1 5 ) e n d 0 ) ( 1 6 ) 若彬5 ) 1 1 2 e ，转( 1 8 ) ( 1 7 ) e n d ( k ) ( 1 8 ) e n d ( 1 ) 下面我们统计算法g s - s h i r ( i ) 每重新启动一步( 即关于) 的基本计算量 计算量按矩阵与向量的乘积a x 和n 维向量的运算( 意指内积运算d o t 和数乘 运算d a x p y ) 来统计；另外还有( m + 1 ) mh e s s e n b e r g 矩阵的最小二乘问题( 记 作l s ) 表2 1 列出用于种子方程组的基本g m r e s ( m ) 部分的计算量、用于( 一 个) 非种子方程组的投影部分( 算法2 1 中第( 1 1 ) 行到第( 1 6 ) 行) 的计算量和总 计算量 表2 1 一g s s h i a ( i ) 每迭代一步( 关于k ) 的计算量 基本g m r e s ( m ) 部分附加投影部分 合计 a xm + l1m + 2 d o t 和d a x p ym 2 + 4 m + 22 m + 1m 2 + 6 m + 3 l sl2 数据显示算法的计算量主要集中在基本g m r e s ( m ) 部分：矩阵与向量乘积是投影 部分的m + l 倍，向量运算是；m 倍丽最小二乘问题由于m 一般很小，其计算 量可以忽略不计相对而言，i t i 取得适当大，投影方法的效果可能更好些。 g s s h i f f ( i i ) 算法的计算量基本上与g s s h i n ( i ) 相同，除了在最小二乘问题的 计算量上前者有一些改善，但这些改善可以忽略不计 2 0 0 4 年上海大学硕士学位论文1 4 2 2 残量性质 首先我们分析算法2 1 的残量性质由( 2 3 ) ， 剐“。一山列) = r ( j ) 一如) = r ( ) 一y m + l 岛棚， 从而唿l f 1 = v 瓢1 r o ) 一岛d ( 力另外。我们记p = + l 嘿l 是7 c 。+ l ( 也，“a ) 上 的正交投影算子则有 f ( ) = t o ) 一a j v m d ( j ) = ( j p ) r ( j ) + ( p r o ) 一a j v m d g ) ) = ( j p ) r ( j + i ，仇+ 1 ( 唿l r u ) 岛d u ) ( 3 1 ) 于是，由( 2 4 ) 和上述表达式，我们得到下列残量的性质 性质2 1 算法g s s h m ( i ) 中的非种子方程组残量力) 满足 i 唏t 十1 蚣；珊。o 嘿，( 删一a t z ) l i 2 ( 2 7 ) 和 幅：j f ( ，一p ) r 1 ；+ 。r a 。i n | | 唿i r 一岛碓 ( 2 8 ) 式( 2 7 ) 表明，牙( j ) 是$ o + ( a 。，r ( s ) 中第j 个非种子方程组的残量在空间 j l c 。+ d a 。，r ( 8 ) ) 中的正交投影模达到最小的解 定理2 2 设。e 赶。( j 4 。，r ( s ) ) 为m 维子空间，则关于算法2 1 ，成立 对于种子方程组， 办。) 3 ) ) + 咒。，其残量f ( s ) 品+ 1 ， 且f ( 3 ) 上a 。咒。c + 1 夥对于非种子方程组， 牙( ) ( ) ) + 。，其残量f ( j ) r ( j ) 十_ i c 州一l ， 且f ( j ) 上a f 证明t1 ) 的结论是g m r e s 方法的自然结论我们证嘎2 ) 牙( ) = 。( 1 + 删 z ( ，) ) + 。也是自然的。 棚= b u ) 山铷= r o ) 一向剥= r ( ，) 一+ 1 岛和 2 0 0 4 年上海大学硕士学位论文15 - 二二： 用如作内积， ( 山) r f u = 霹嘿l ( r 一+ l 岛d ) = 霹( 嘿1 r 。) 一岛d 。) 由于d o = 。g m i r l d 聊i i 唿l r 。) 一岛硎2 是最小二乘解，其满足正规方程： 蟹岛d o = 霹唿 ，r c j ) ， 从而( a j ) 于f 0 ) = 0 o 定理2 - 3 设也有谱分解a = z a z ，a = d i a g ( a 1 ，k ) 记旗( t ) 是残量r ( ，) 在+ 1 中的投影多项式，即p r 。) = 旗( a 。) 并设其德( q 一) 0 假设含 有a j 所有特征值的椭圆e j ( c j ，勺，q ) 不包含原点，且砂g ( q 一吼) o ，则算法2 j 的近似解牙“) 成立 “fu。s“cjp，r。+足。cz，i坛)c一如，i(i考；耄蔫)”cz。， 证明： 记“u = 扣) 由( 2 6 ) 及最小性质( 2 8 ) ， 一讹i i ( - p ) r 制删一a j u u ) | | 2 = i i ( x p ) r e j ) i i m r a i n 。l i p r 。) 一山“1 1 2 由于t ，记“= 一i ( a 。) ”l ，d = q a 。，其中毋。l ( t ) 是m 一1 阶多项式。则 。r a i n 。r l p 7 u 一山j j 2 ：十一m i n 一，懒饥) 口l 一如一i ( 也) 北 。一m ，i n 一。l l 渺g ( a ，+ 。，) 一南一l ( 冯+ a 列 l i l 2 。i 姥一l ，眺p 糍) ( 。) ：，l i 螺( a a - 1 1 2 l 城( 一i 嘏，；r a 蒯i n ( 。) ：1l i p 鲫( m - 其中p 瓣( ) 2 礤栖( 旗o 一一) 一t 一 一) 现取 螂= 鬻， ( 3 5 ) 其中表示第一类m 次c h e b y s b o y 多项式。于是利用c h e b y s h o v 多项式的性质， 可得到定弹的结诊 2 0 0 4 年上海大学硕士学位论文 1 6 当qi 口。，彤时，即对于多右端问题( 1 2 ) ，此时a i = a 。，结论( 2 9 ) 便是文 【1 4 ，2 9 ，3 1 】中的结论式( 2 9 ) 右端第一项表明，非种子方程组的残量r ( 2 在投影空 间m + t ( 也，r ( 5 ) ) 中的成份多少决定了该项的大小，成份越多，该项越小 下面我们分析算法2 2 的( 伪) 残量性质。注意到伪残量仍有表达式 f j ( j ) = 6 ( j ) 一a ( 5 ) 孟( j ) = r s ( j ) 一a ( 5 ) d 0 ) = r s u ) 一l 乞+ l 歌d ( 川， 和 俺( j ) = r s ( j ) 一a ( 5 ) d ( j ) = ( ，一p ) r s ( j ) + i 厶+ 1 ( 1 儡t + l r s ( j ) 一f i n d ( j ) ) 完全类似地成立下列结论 性质2 4 算法g s - s h i _ c c ( i i ) 中的非种子方程组伪残量而( j ) 满足 蠕1 f s ( j ) f f 2 = 。溉。f f 唿1 ( 6 一a i s ) x ) 1 1 2 ( 2 t o ) 和 一嘲= 一p ) r s ( j ) 肚枷m i n 。8 嘿1 r s 一砜礁 ( 2 1 1 ) o 定理25设厄。；丘。( a ( “，r ( 5 ) ) 为m 维子空间，则关于算法2 2 成立 ) 对于种子方程组，i ( 8 ) e ( 2 ( 3 + 赶。，其残量i ( 3 ) 瓦。+ 1 ， 且f ( s ) 上a ( s ) e 。k m + 1 剀对于非种子方程组，面( ，) z ( ，) ) + 坛。，其伪残量而( ，) r s ( j ) ) + 。+ 1 ， 且f j ( ，) 上a ( 3 ) k m 瓦m + 1 定理2 6设4 ( s ) 有谱分解a ( s ) = z a z 一，妒g (