（信号与信息处理专业论文）特征值问题的连续时间域解法.pdf

特征值问题的连续时间域解法 摘要 一直以来，在数值线性代数中，线性特征值问题都是一个重要的 研究领域。特征值问题固有的非线性特性引入了很多计算方面的问 题。除了非常特别的情形，一般来说，我们不会选择通过显式求解特 征方程来求得特征值。原因在于，我们没有一种数值计算上稳定的方 法来获得特征方程的系数。进一步，即使可以精确地得到特征方程， 我们依然无法在有限精度的限制下去计算方程的根。摄动理论显示， 对于特征方程系数的小幅度扰动即可导致根的大幅度变化。同时，数 值求解特征向量也是一个困难的事情，尤其对于那些夹角很小的特征 向量。一个比较清晰的结论是，当前的所有方法都具有迭代的特性。 在特征值问题的传统解法之外，一些连续时间域的方法开始出 现。最近，gh g o l u b 与l z l i a o 提出了一种连续方法用来求解极 端特征值问题与内侧特征值问题。然而，这一连续方法收敛的很慢， 有时会停滞。另一个问题在于，这一方法在某些稠密、坏条件的情形 下无法正常工作。在本文中，我们提出了一个新的连续方法，名为梯 度r a y l e i g h 商方法( g r q m ) 。这一方法克服了这些缺陷。g r q m 是 两方面的有机结合体，首先，是一个具有快速收敛特性的动力系统 ( g r q ) ，其次，是一个基于线性化的隐式欧拉方法、并结合了信赖 域方法来控制时间步长的常微分方程算法( t r l i e m ) 。g r q m 的主 要特点在于使用r a y l e i g h 商来更新每一步中所得到的拉格朗日乘子。 数值算例进一步证明了该方法的快速收敛性和更高的鲁棒性。 关键词：特征值特征向量动力系统r a y l e i g h 商微分方程 信赖域方法 北京邮i 【! 人学硕i j 学位论文 g r q m ：ar o b u s ta n df a s t c o n t i n u o u sm e t h o df o re x t r e m e a n di n t e r i o re i g e n v a l u ep r o b l e m s a b s t r a c t l i n e a re i g e n v a l u ep r o b l e m sc o n t i n u et ob ea ni m p o r t a n ta n dr e l e v a n t a r e ao fr e s e a r c hi nn u m e r i c a ll i n e a ra l g e b r a t h ei n h e r e n t l yn o n l i n e a r p r o p e r t yl e a d st om a n yc o m p u t a t i o n a lp r o b l e m s e x c e p tf o rp a r t i c u l a r c a s e s ，c o m p u t a t i o no ft h ee i g e n v a l u e st h r o u g ht h ee x p l i c i tc o n s t r u c t i o no f t h ec h a r a c t e r i s t i ce q u a t i o ni sn o ta no p t i o ns i n c et h ec o e f f i c i e n t so ft h e e q u a t i o n c a nn o tb ec o m p u t e df r o md e t e r m i n a n te v a l u a t i o n si na n u m e r i c a l l ys t a b l ew a y f u r t h e r m o r e e v e ni ft h ee x a c tc h a r a c t e r i s t i c e q u a t i o nc a nb eo b t a i n e d ，w ec o u l dn o tc o m p u t et h er o o t so ft h ee q u a t i o n i nt h e d e m a n d i n gp r e c i s i o n p e r t u r b a t i o nt h e o r yi n d i c a t e st h a ts m a l l p e r t u r b a t i o no ft h ec o e 伍c i e n t sc o u l d1 e a dt ot h eb i gp e r t u r b a t i o no ft h e r o o t a tt h es a m et i m e n u m e r i c a l l ys o l v i n gt h ee i g e n v e c t o r si sa l s oa d i 衔c u l tt a s k e s p e c i a l l yf o rt h o s ew i t hs m a l la n g l eb e t w e e nt h e m ac l e a r c o n c l u s i o ni st h a ta l lo ft h em o d e mm e t h o d sa r eo fa ni t e r a t i v en a t u r e r e c e n t l yac o n t i n u o u sm e t h o dh a sb e e np r o p o s e db vgh g o l u b a n dl z l i a oa sa na l t e r n a t i v ew a yt os o l v eb o t he x t r e m ea n di n t e r i o r e i g e n v a l u ep r o b l e m s h o w e v e r , t h i sc o n t i n u o u sm e t h o dc o n v e r g e ss l o w l y , s o m e t i m e ss t a g n a t e s a n o t h e rp r o b l e mi st h a ti tf a i l st ow o r ki ns o m e d e n s ea n di 1 1 c o n d i t i o n e ds i t u a t i o n s i nt h i sp a p e r , an e wc o n t i n u o u s m e t h o dn a m e dg r a d i e n tr a y l e i g hq u o t i e n tm e t h o d ( g r q m ) i sp r o p o s e d ， w h i c ho v e r c o m e st h e s ed e f i c i e n c i e s g r q mi sac o m b i n a t i o no far a p i d l y c o n v e r g e n td y n a m i c a ls y s t e m ( g r q ) a n d al i n e a r l yi m p l i c i te u l e rt y p e o d es o l v e rm i x e dw i t ht r u s t r e g i o nt i m es t e pc o n t r 0 1( t r l i e m ) t h e m a i nf e a t u r eo fg r q mi st ou s et h er a y l e i g hq u o t i e n tt ou p d a t et h e l a g r a n g em u l t i p l i e rd e r i v e df r o mas e q u e n c eo fc o n s t r a i n e do p t i m i z a t i o n s u b p r o b l e m s n u m e r i c a le x p e r i m e n t si n d i c a t ea l s ot h a tt h en e wm e t h o di s 北京邮i u 人! 二倾i ：学位论文 u s u a l l yf a s t e l a n dm o r er o b u s t k e yw o r d s ：e i g e n v a l u e e i g e n v e c t o rd y n a m i c a ls y s t e mr a y l e i g h q u o t i e n td if f e r e n t i a le q u a t i o nt r u s t r e g i o nt e c h n i q u e 北京邮l 【1 人学颂i 。学位论文 北京邮电大学学位论文独创性声明 本人声明所，? 交的沦文是本人在导师指导下进行的研究工作及取得的研究 成果。尽我所知，除了文中特别加以标注和致谢中所罗列的内容以外，论文中不 包含其他人已经发表或撰写过的研究成果，也不包含为获得北京邮电大学或其他 教育机构的学位或证书而使用过的材料。与我一同工作的同志对本研究所做的任 何贡献均已在论文i l ，作了明确的说明并表示了谢意。 申请学位论文j 资料若有不实之处，本人承担一切相关责任。 本人签名：；主4 刍鎏舡一同期：- 二盔翌旱立正l jj 关于论文使用授权的说明 学位论文作者完全了解北京邮电大学有关保留和使用学位论文的规定，即： 研究生在校攻读学位期间论文工作的知识产权单位属北京邮电大学。学校有权保 留并向国家有关部门或机构送交论文的复印件和磁盘，允许学位论文被查阅和借 阅：学校可以公布学位论文的全部或部分内容，可以允许采用影印、缩印或其它 复制手段保存、汇编学位论文。( 保密的学位论文在解密后遵守此规定) 保密论文注释：本学位论文属于保密在一年解密后适用本授权书。非保密论 文注释：本学位 本人签名： 导师签名： 适用本授权书。 同期：互竺罩：! i 羔 f i 期： 2 竺仝：乏! 羔 儿 ：i ：i l ij 、。，锄i j 学位论文 1 1 特征值问题 第一章绪论 一直以来，在数值线性代数中，线性特征值问题都是一个重要的研究领域。 大约1 5 0 年前，雅可比写下了那篇著名的论文，讨论如何求解问题出：厶，其中 彳：彳，。之后，其它的一些方法丌始流行，例如幂迭代法。幂迭代法基于如下的 思想，即当一个归一化的向量不断的通过一个矩阵进行变换，最终，它的方向将 与具有最大绝对值的特征值所对应的特征向量吻合。幂迭代法的收敛性取决于第 二大特征值( 绝对值) 与最大特征值( 绝对值) 的比率，而事实上，大多数应用 都无法接受幂迭代法的低收敛率。当前，幂迭代法仍然被人们所使用，不过大多 是隐式地应用于某些更有效率的特征值算法之中，例如k r y l o v 方法，q r 方法等。 再后来，在标准特征值问题瓜：厶之外，开始出现广义特征值问题a x ：2 b x ， 以及高阶多项式特征值问题： ( 彳o + 刎i + + 五“a 。) x = 0 当然，还有奇异值分解。 特征值问题固有的非线性特性引入了很多计算方面的问题。除了非常特别的 情形，一般来说，我们不会选择通过显式求解特征方程 d e t ( a 一2 1 = 0 来求得特征值a 。原因在于，我们没有一种数值计算上稳定的方法来获得特征方 程的系数。进一步，即使可以精确的得到特征方程，我们依然无法在有限精度的 限制下去计算方程的根。摄动理论显示，对于特征方程系数的小幅度扰动即可导 致根的大幅度变化。同时，数值求解特征向量也是一个困难的事情，尤其对于那 些夹角很小的特征向量。所以，如何精确有效的计算特征值与相应的特征向量一 直是研究的重点。 1 2 特征值问题的数值线性代数解法 一个比较清晰的结论是，当前的所有方法都具有迭代的特性。事实上，迭代 已经成为一种必须。原因在于，试图在有限的计算次数中，去求n 阶矩阵的特征 值，这本身便是与a b e l r u f f i n i 定理( 即g a l o i s 理论) 所矛盾的。该定理指出， 可以计算大于四阶的多项式的根的算法是不存在的。所以，矩阵特征值问题的算 l j j i i i l l _ f i u 人学烦i j 学位论文 法必然是迭代性质的，mj j j 题的重点是寻找更快、更精确的达1 弋弹浊。 起仞，l a n c z o s 方f = ：破认为是解决特征值问题的万能药。1 i 过，乡夸过多次研 究发现，由于舍入误差的影l l m 陔方法所得到的结果异常。直至j ：纪7 0 年代， p a i g e 发现可以通过纯迭代的方法来利用l a n c z o s 方法，从而得到i l i 确的特征值 信息【1 1 。p a i g e 所得结论中的关键点在于，揭示了在l a n c z o s 过程中产生的正交性 损失与特征对( 特征向量与相对应的特征值) 之间的关系。这为后来人的研究奠 定了基础，并最终促使将l a n c z o s 方法发展成求解大型稀疏对称矩阵的强有力方 法。 对于非对称矩阵解法，同样将经历了一个曲折、漫长的发展。这一过程中， 涌现出三个具有影响力的方法。首先，通过引入所谓的l o o k a h e a d 策略，双侧 的l a r l c z o s 过程使得鲁棒性更强。其次，隐式重启方法的发现，使得我们在a m o l d i 过程中可以满足内存的限制。再次，是j o c a b i d a v i d s o n 方法。这三个方法均为 子空间迭代方法，核心思想是通过将大型问题降至小型问题来处理。 1 3 特征值问题的连续时间域解法 在特征值问题的传统解法之外，一些连续时间域的方法开始出现，gh g o l u b 与l - z l i a o 的文章c o n t i n u o u sm e t h o d sf o ，e x t r e m ea n di n t e r i o l 。e i g e m ，a l u e p r o b l e m s 是近年来连续时间域解法中的著作1 2 j 。之后，他们的方法被成功地用于 求解对称广义特征值问题p 】。g o l u b 与l i a o 的思路在于，首先，将特征值问题转 化为等价的优化问题，而他们所建立的优化问题是一个具有单位球约束的严格凹 的函数，一般来说，此类函数在最优化理论中是易于求解的。紧接着，对于优化 问题，g o l u b 与l i a o 构建了一个常微分方程o d e ( 动力系统) 和一个价值函数。 这一动力系统与文献 4 中的线性投影方程( l p e ) 相类似，有效的利用了价值函 数中的微分信息。这样，就将特征值问题转化到了连续时间域中，用o d e 算法 可以进行求解。 然而，g o l u b 与l i a o 所提出的方法在很多情形下显现出非常慢的收敛性【2 3 】。 进一步，在某些稠密坏条件问题中，该方法无法运行，具体可以参见本文3 5 中的数值算例。 这便激发我们去研究更快更稳定的算法。本课题的研究从两个大的方向着 手：一是最优化问题，二是常微分方程组( o d e ) 的数值解法。一方面，基于文 献【5 ，6 ，我们探究了序列无约束优化问题的本质，这引导我们设计了个更适 合的动力系统，即“梯度r a y l e i g h 商”( g r q ) 动力系统。另一方面，文献 2 】中 的方法从本质上反映了最优化问题与常微分方程组( o d e ) 数值解法之间的紧密 北京i l l i ；i u ，、i 川j 论文 联系1 7 8 。o 1 2 ，13 1 。文献【9 ，1 0 ，1 1 启发我仃j i ： i j 一芒1 信赖域方法相结合的线性化隐式 i 杉： t 公式，从而得以为g r q 动力系统i 奠汁一个更有效、鲁棒性更强的o d e 7 = 法。 本文将按如下的方式进行陈述。在第一二章中，我们将重点讨论数值分析领域 中的两个大的研究方向，即优化问题与常微分方程组( o d e ) 的数值解法，并进 一步分析这二者之间的紧密联系。这将是最终形成我们的算法的理论基础。在第 三章中，我们将详细陈述研究内容。首先，是原始问题的转化策略，将最优化理 论作为桥梁，把原本的特征值问题转化为g r q 动力系统。同时，我们给出了该 动力系统的收敛性分析。在接下来的一节中，我们将会详细说明一个全新的解法 梯度r a y l e i g h 商方法( g r q m ) 。它是两方面结合的产物，一方面，是一个 快速收敛的动力系统，即g r q ；另一方面，是线性化的隐式欧拉形式的o d e 解 法( t r l i e m ) ，其中运用了信赖域的方法去自适应地控制时间步长。梯度r a y l e i g h 商方法( g r q m ) 的主要特色是，在每一步迭代中运用r a y l e i g h 商去更新当前 约束优化子问题l a g r a n g e 乘子。在3 3 中我们将会详细地讨论该算法的主要特 性。在3 4 中，我们将详细说明g o l u b 与l i a o 的处理方法，该方法与我们的不 同之处有两点：一，该方法使用的动力系统( 本文称之为g l 动力系统) 与我们 的不同：二，g o l u b 与l i a o 并没有设计一个o d e 的算法来求解g l 动力系统， 而是使用了m a t l a b 中o d e 算法库。最后，我们通过典型的数值算例来检验 g r q m 方法的性能。 l j ：tl ：l l 人学坝i ：学位论文 第二章基本理论与方法分析 g r q m 方法的构建基于对数值分析领域中的两个分支的研究：一是优化问 题，二是常微分方程组( o d e ) 数值算法。在这一章中，我们将分别讨论这两个 分支中的关键理论与方法，并说明这两个数值分析领域之间的紧密联系。 2 1 优化问题的理论与方法 首先，我们将详细叙述无约束优化与约束优化的相关理论。其次，下文将分 别讨论两大类最优化方法，即线搜索类方法与信赖域类方法。最后，我们将着重 叙述在非线性约束优化中如何构造最优化子问题序列，这是第三章中所设计的算 法的基本思想。 2 1 1 无约束优化理论 在无约束优化中【1 4 , 1 5 ，我们试图去最小化如下的目标函数， r a i n ( x ) l e 月” 其中，厂：r 一一r 是一个平滑的函数。 无约束优化理论中所运用的核心数学工具便是t a y l o r 定理。我们首先将它陈 述如下。 定理2 1t a y l o r 定理假设f ：r 一_ r 是连续可微的，且p 彤。那么 f ( x4 - p ) = f ( x ) + v f ( x + t p ) 7p ， 其中t ( o ，1 ) 。进一步，如果是二次连续可微的，那么 v f ( x + p ) = v f ( x ) + 上v 2 厂( x + t p ) p d t ， 且 厂( 工+ p ) = 厂( x ) + v 厂( x ) 7 p + i 1p 7 v2 ( x + t p ) p ， z 其中f ( o ，1 ) 。口 假设j 是一个局部最优解，无约束最优化的必要条件均是围绕v f ( x ) 和 v 二( j ) 展开的。 定理2 2 无约束优化一阶必要条件如果j 是一个局部最优解，且厂在了的 一个开邻区中是连续可微的，那么v f ( x ) ：0 。 4 北京邮i u 人学坝| 0i 0 - 证f j ： j i j ) 2i i l ：法，t r j 设v f c ) 0 。定义向量p = 一v ( 、) 汁意p t v f ( x ) = i 耵( x ) 0 ： 0 。由 二j ：确：、。i ；f , lj 王厂是连续可微的，所以存在一个卡，j “i7 使社j p r v f ( x + t p ) o ，t 【0 ，t 】 那么，基于t a y l o r 定理( 2 1 ) ，对于任意的f ( o ，刀，我们有 f ( x + 劢) = f ( x + ) + 劢7 v 厂( x + p ) ，f ( o ，f ) 所以，对于所有的i ( o ，r l ，f f x + ) 厂( z ) 。从而，我们找到了一个从工出发的 方向，可以使，减小，所以x 不是一个局部最优解，与题设矛盾。1 3 定理2 3 无约束优化二阶必要条件如果x 是一个局部最优解，且厂和v ：，在 x 的一个开邻区中是连续的，那么夥( 工) ：ol e iv ：( t ) 半j 下定的。 证明： 由定理( 2 2 ) ，我们知道夥( 工) ：0 。用反证法，假设v ：f ( x ) 不是半正定的。那么， 我们可以定义一个向量p 使得p r 夥( x ) p o 。由于v ! v 2 厂在x 的一个开邻区中是 连续的，所以存在一个标量7 _ 使得 p 7 v2 f ( x + t p ) p 0 ，t 【0 7t 】 那么，基于t a y l o r 定理( 2 1 ) ，对于任意的i 1 0 刀，我们有 f ( x + 劢) = ( 工) + 劢7 v f ( x ) + 去f2 p 7 v2 厂( x + + t p ) p ， f ( o ，f ) z 所以，对于所有的；( o ，r l ，( x + 劢) o ，使得对于丌球体d = zz - - x * l i ( ) ，工 是厂的一个严格局部最优解。口 北j 川一 j i 也j 、知；r 川。；f t 论文 定理2 5 当厂是凸的，那么仃意j 0 i i i j 2 1 ：优解x 同时也是全局最优解。进步， 如果厂是可微的，那么任意不动点、。( t t l j 丫( v ) = 0 ) 是厂的全局最优解。 证明： 对于第一部分，假设工是个局部最优解而非全局最优解。那么，我们可以 找到一个点z 使得( z ) ( x ) 。考虑点z 与点工。之间的连线，表示如下 工= 允+ ( 1 一a ) 工， 五( o ，l 】 由于函数厂具有凸的性质，所以 f ( x ) s 矽( z ) + ( 1 2 ) f ( x ) f ( x ) 由于工的任意邻域n 必然包含于点z 与点j 之间的连线，所以总存在点n 使得 ( 工) ( x ) ，所以，膏不是局部最优解，这与题设矛盾。 对于第二部分，同样假设x 不是一个全局最优解，选择点z 使得厂( z ) o 换一个表达方式，即 f v c , ( x ) 7 = 0 ， i e l ( 0e 五( z ) v c , ( x ) 7 = o ， i 彳( j ) n ，w i t h 乃 0 l v c i ( x * ) 7 缈0 ，i a ( x ) n 1w i t h 五= 0 i 定理2 7 约束优化二阶必要条件假设x 是问题( 2 1 ) 的一个局部解，且在 工。处线性独立约束限制满足。令z 为一个l a g r a n g e 乘子向量，且满足k k t 条件 ( 2 3 ) ，那么 7 v 。( x ，z ) o ，缈丘( 彳) 式( 2 4 ) 定理2 8 约束优化二阶充分条件假设对于某个可行点x r ”，存在一个 l a g r a n g e 乘子向量彳满足k k t 条件( 2 3 ) 。同时，假设如下条件满足 7 v ( x ，z ) 印0 ，山丘( 彳) ，0 式( 2 5 ) 北京邮l u 人掌卜j t 卫 圳么。是问题( 2 1 ) 的一个严格局部角华。口 2 1 3 线搜索类方法 2 1 3 1 综述 对于线搜索方法而言，算法为当前点选择一个前进方向s 。，沿着这个方向从 当b ，j 迭代点x k 前行到新的迭代点，同时使得目标函数值下降。搜索的步长取决于 求解如下的一个1 维最优化问题，它的最优解便是步长： m i n f ( x + 女) 口 式( 2 - 6 ) 事实上，我们没有必要精确求解这个问题，那样的计算代价太大。通常，线搜索 算法产生有限长度的步长序列，直至寻找到一个逼近这个问题的解的步长。线搜 索方法的标准形式如下： 以+ l = 以+ s 式( 2 7 ) 其中，吼是一个正的标量，代表步长，j 。是搜索的方向，一般是归一化的向量。 事实上，线搜索方法可以使用任意可以使闷标函数，- 下降的搜索方向。通常， 任意的下降方向可以被诠释为与一v 以德夹角严格小于万2 的方向。我们可以用 t a y l o r 定理2 1 来说明，假设搜索方向为。且步长参数为足够小的一个数s ，那 么 f ( x i + 岱) = f ( x t ) + 岱i 。v l + o ( e 2 ) 式( 2 - 8 ) 当“是一个下降方向时，s k 与的夹角q 满足c o s b 0 ，使得 s 。，= 恢洲以i i c o s 啡 0 进而可得，对于所有足够d 、的正数占，厂( _ + k t s 。) s b 0 式( 2 3 0 ) 另一方面，当较小时，由于搜索步长被限制，所以模型函数埘( ，) 中的二次项部 分对最终的解的影响很小，此时，我们可以用线性模型下求得的解来近似模型函 数m ( 。) 的解，即 s ( ) 一高，w h e n ai ss m a l l 式( 2 3 1 ) 蚓l 上面考虑了两种极端的情况，那么对于中间状态的信赖域半径来晚，解j 。( ) 通 常是一条有弧度的轨迹，如图2 2 中所示最优搜索轨迹。 夙 r ， s 一2 +s 丁， g + 厂 = d 怿州 l ，：i ：i uj j f i j j ! i ：学位论文 在d o g l e g 方法中，我们j i jh2 2 - 一所示d o g l e g 路径来近似最优搜索4 ) t j j ， 它分为两段。首先，我们找到c a u c h y 点s 。，然后从c a u c h y 点出发，阳厶：个步 长。迈进。利朋式( 2 2 7 ) j 。( 2 2 8 ) ，我们可以得到 ，一麓goo 第二段从c a u c h y 点指向全步长s e ，与信赖域边界的交点a 。我们可以将整体的 表达式写成如下形式 亨c 盯，= i “+ 。k 一。，。s 口一s u ，0 。 r t c o 是惩罚系数。通过迫使向零逼近，我们对于违反约束情形的惩罚也 越来越重。考虑一个序列。) j i l i r a 。，肌= 0 ，我们可以运用无约束优化的方法去 寻找q ( j ；以) 的一个近似解。所以，适当地选择序列饥) 以及起始的搜索点，对 于每一个肌，我们或许可以仅通过若干次步就可以得到一个无约束最小值解。 北京邮l u 人学颂l 学化论文 算法2 4 二次制函数法 注意，惩罚系数序列l 的选择是需要技巧的，这主要基于对本次迭代中所形成 的无约束优化问题(