China-pub.coM 下载 第8章 特征值和特征向量 MATLAB中的命令计算特征值和特征向量很方便,可以得到不同的子结果和分解,这在 线性代数教学时很有用。注意,本章中的命令只能对二维矩阵操作。 8.1特征值和特征向量的计算 假设A是一个m×n的矩阵,A的特征值问题就是找到方程组的解: AX=入x 其中入是一个标量,x是一个长度为的列向量。标量入是A的特征值,x是相对应的特征向 量。对于实数矩阵A来说,特征值和特征向量可能是复数。一个n×n的矩阵有n个特征值,表 示为入,2,,入m0 MATLAB中用命令eig来确定矩阵A的特征值和特征向量。特征向量的规格化,就是每个 特征向量的欧几里得范数为1:参见7.6节。 命令eig自动完成对矩阵A的平衡化。这就要求MATLAB找出一个相似变换矩阵Q,满足 条件A。求A=QAQ的特征值比求A的特征值条件更好些。万一A有一个和机器错误大小一 样的元素,平衡化对于计算过程是没有好处的。带有参数nobalance的命令eig可用来计算 没有这个变换矩阵的特征值和特征向量。 命令集79 特征值和特征向量 eig(A) 求包含矩阵A的特征值的向量。 [X,D]=eig(A) 产生一个矩阵A的特征值在对角线上的对角矩阵D和矩阵 X,它们的列是相应的特征向量,满足AX=XD。为了得到 有更好条件特征值的矩阵要进行相似变换。 [X,D]= 不经过平衡处理求得矩阵A的特征值和特征向量,也就是 eig(A,'nobalance') 不进行平衡相似变换。 balance(A) 求平衡矩阵。 [T,B]=balance(A) 找到一个相似变换矩阵T和矩阵B,使得它们满足B=T-AT。 B是用命令balance求得的平衡矩阵。 eigs(A) 返回一个由矩阵A的部分特征值组成的向量,和命令eig 样,但是不返回全部的特征值。如果不带有参量,则计 算出最大的特征值。当计算所有特征值时,如果矩阵A的 秩不小于6,则计算出6个特征值来。 eigs(f,n) 求出矩阵A的部分特征值。在使用一个矩阵列的线性运算 符时,字符串中包含的是M文件的文件名,n指定问题的 阶次。用这种方法来求特征值比开始就用运算符来求要快
下载 第8章 特征值和特征向量 M AT L A B中的命令计算特征值和特征向量很方便,可以得到不同的子结果和分解,这在 线性代数教学时很有用。注意,本章中的命令只能对二维矩阵操作。 8.1 特征值和特征向量的计算 假设A是一个m×n的矩阵,A的特征值问题就是找到方程组的解: 其中l是一个标量,x是一个长度为n的列向量。标量l是A的特征值,x是相对应的特征向 量。对于实数矩阵 A来说,特征值和特征向量可能是复数。一个 n×n的矩阵有n个特征值,表 示为l1,l2,. . .,ln。 M AT L A B中用命令e i g来确定矩阵A的特征值和特征向量。特征向量的规格化,就是每个 特征向量的欧几里得范数为 1;参见7 . 6节。 命令e i g自动完成对矩阵 A的平衡化。这就要求 M AT L A B找出一个相似变换矩阵 Q,满足 条件 。求 的特征值比求A的特征值条件更好些。万一 A有一个和机器错误大小一 样的元素,平衡化对于计算过程是没有好处的。带有参数 n o b a l a n c e的命令e i g可用来计算 没有这个变换矩阵的特征值和特征向量。 命令集7 9 特征值和特征向量 e i g ( A ) 求包含矩阵A的特征值的向量。 [ X , D ] = e i g ( A ) 产生一个矩阵 A的特征值在对角线上的对角矩阵 D和矩阵 X,它们的列是相应的特征向量,满足 A X=X D。为了得到 有更好条件特征值的矩阵要进行相似变换。 [ X , D ] = 不经过平衡处理求得矩阵A的特征值和特征向量,也就是 e i g ( A , ’ n o b a l a n c e ’ ) 不进行平衡相似变换。 b a l a n c e ( A ) 求平衡矩阵。 [ T , B ] = b a l a n c e ( A ) 找到一个相似变换矩阵T和矩阵B,使得它们满足B=T-1AT。 B是用命令b a l a n c e求得的平衡矩阵。 e i g s ( A ) 返回一个由矩阵 A的部分特征值组成的向量,和命令 e i g 一样,但是不返回全部的特征值。如果不带有参量,则计 算出最大的特征值。当计算所有特征值时,如果矩阵 A的 秩不小于6,则计算出6个特征值来。 e i g s ( f , n ) 求出矩阵A的部分特征值。在使用一个矩阵列的线性运算 符时,字符串f中包含的是M文件的文件名, n指定问题的 阶次。用这种方法来求特征值比开始就用运算符来求要快
114 China-pub.com MATLAB5手册 下载 eigs (A,B,k,sigma) 求矩阵A的部分特征值,矩阵B的大小和A相同:如果没有 给出B=eye(size(A),那么k就是要计算的特征值的个数: 如果k没有给出,就用小于6的数或者A的秩:变量sigma是 一个实数或者复数的移位参数,或者下列文本字符串中的 一个,文本字符串指明需要的是哪种特征值: “1m'最大的特征值(缺省) “sm’最小的特征值 ‘r’最大的实数部分 ‘sr'最小的实数部分 “be’同时求得最大和最小的实数部分 condeig(A) 返回一个由矩阵A的特征值条件数组成的向量。 [v,D,s]=condeig(A)返回[V,D]=eig(A)和s=condeig(A)。 如果A是实数矩阵,MATLAB在计算中用QR因式分解:否则用QZ因式分解。 左特征向量是满足下面条件的非零行向量y: yA=hy 如果用命令eig对A'作用,也可以计算出左特征向量,因为: A'y'=iy 这里的撇号'代表矩阵的转置和共轭复数(见3.4节),入上的短杠表示共轭复数。矩阵特征值 的集合称为矩阵的谱,谱半径p(A)定义为max(abs(eig(A)))。矩阵A的特征值的乘积等于 det(A),和等于trace(A),这是矩阵A主对角线上元素的和。 如果X一个列向量为A的特征向量的矩阵,并且它的秩为,那么特征向量线性无关。如 果不是这样,则称矩阵为缺陷阵。如果XX=I,则特征向量正交,这对于对称矩阵是成立的。 ■例8.1 矩阵A定义为: 9 -3 -16 13 16 3 10 (a)运行命令[Evect,Evalue]=eig(A),得到结果为: Evect -0.7071 -0.5774 -0.5774 0.7071 0.5774 -0.5774 0.0000 0.5774 0.5774 Evalue -6.0000 0 0 0 10.0000 0 0 0 4.0000 可知特征值都是非零数,矩阵是满秩的,可以用therank=rank(Evect)来确认: therank=
e i g s ( A ,B , k , s i g m a) 求矩阵A的部分特征值,矩阵B的大小和A相同;如果没有 给出B= e y e ( s i z e ( A ) ),那么k就是要计算的特征值的个数; 如果k没有给出,就用小于6的数或者A的秩;变量s i g m a是 一个实数或者复数的移位参数,或者下列文本字符串中的 一个,文本字符串指明需要的是哪种特征值: ‘l m’最大的特征值(缺省) ‘s m’最小的特征值 ‘l r’最大的实数部分 ‘s r’最小的实数部分 ‘b e’同时求得最大和最小的实数部分 c o n d e i g ( A ) 返回一个由矩阵A的特征值条件数组成的向量。 [ V , D , s ] = c o n d e i g ( A ) 返回[ V , D ] = e i g ( A )和s = c o n d e i g ( A )。 如果A是实数矩阵,M AT L A B在计算中用Q R因式分解;否则用Q Z因式分解。 左特征向量是满足下面条件的非零行向量 y: y A =ly 如果用命令e i g对A¢作用,也可以计算出左特征向量,因为: 这里的撇号'代表矩阵的转置和共轭复数 (见3 . 4节), 上的短杠表示共轭复数。矩阵特征值 的集合称为矩阵的谱,谱半径 ( A )定义为m a x ( a b s ( e i g ( A ) ) )。矩阵A的特征值的乘积等于 d e t ( A ),和等于t r a c e ( A ),这是矩阵A主对角线上元素的和。 如果X一个列向量为 A的特征向量的矩阵,并且它的秩为 n,那么特征向量线性无关。如 果不是这样,则称矩阵为缺陷阵。如果 X ' X=I,则特征向量正交,这对于对称矩阵是成立的。 ■ 例8 . 1 矩阵A定义为: (a) 运行命令[ E v e c t , E v a l u e ] = e i g ( A ),得到结果为: 可知特征值都是非零数,矩阵是满秩的,可以用 t h e r a n k = r a n k ( E v e c t )来确认: t h e r a n k = 3 1 1 4 M ATLAB 5 手册 下载
China-pub.com 第8章特征值和特征向量 115 下载 令M=Evect'*Evect得: M= 1.0000 0.8165 -0.0000 0.8165 1.0000 0.3333 -0.0000 0.3333 1.0000 可知特征向量没有相互正交。 (b)determinant prod(diag(Evalue)),... determinant2 det(A) 给出结果为: determinant -240.0000 determinant2 -240 可知行列式等于特征值的积。 (c)theTrace =trace(A),theTrace2 =sum(diag(Evalue)) 结果为: theTrace 8 theTrace2= 8.0000 ■ 可知矩阵的迹等于特征值的和。 如果矩阵A是实数矩阵,但是有复数特征值,那么这些特征值是以共轭复数的形式出现的。 如果[X,D]=eig(A),可以用命令cdf2rdf将矩阵D转换为一个实数块对角矩阵。在对角线 上用一个2×2实数块代替共轭复数对。 命令集80 复对角矩阵变成实对角矩阵 [Y,E]=cdf2rdf(X,D)将复对角矩0变成实对角矩阵氧,Y的列不是A的特征向量。 ■例8.2 假设矩阵A为: 01 0 A= -100 003 则运行[X,D]=eig(A)可得: X= 0.7071 0.7071 0 0+0.7071i 0-0.70711 0 0 1.0000 D= 0+1.0000i 0 0-1.00001 0 0 0 3.0000
令M = E v e c t '*E v e c t得: 可知特征向量没有相互正交。 给出结果为: 可知行列式等于特征值的积。 结果为: 可知矩阵的迹等于特征值的和。 如果矩阵A是实数矩阵,但是有复数特征值,那么这些特征值是以共轭复数的形式出现的。 如果[ X , D ] = e i g ( A ),可以用命令c d f 2 r d f将矩阵D转换为一个实数块对角矩阵。在对角线 上用一个2×2实数块代替共轭复数对。 命令集8 0 复对角矩阵变成实对角矩阵 [ Y , E ] = c d f 2 r d f ( X , D ) 将复对角矩阵D变成实对角矩阵E,Y的列不是A的特征向量。 ■ 例8 . 2 假设矩阵A为: 则运行[ X , D ] = e i g ( A )可得: 第8章 特征值和特征向量 1 1 5 下载 ■
116 China-pub.coM MATLAB5手册 下载 X和D都是复数矩阵,运行命令:[Y,E]=cdf2rdf(X,D),结果为: Y= 0.7071 0 0 00.7071 0 01.0000 E= 0 1 0 -1 0 0 0 3 所得的矩阵E正好是A矩阵。 ■ 注意特征向量是特征多项式det()I一A)=0的根,其中I是单位矩阵。用命令poly来求 特征多项式,参见11.1节。命令roots也可求得特征值,但是用命令eig求得的特征 值更准确,精度更高。 广义特征值问题就是找到方程组Ax=入Bx的重要解,其中B也是一个n×n的矩阵。入值和向 量x分别称为广义特征值和广义特征向量。 如果B是一个奇异矩阵,则用QZ算法来求解。 标准和广义特征值问题都属于矩阵多项式特征问题,都可以用命令polyeig来求它们的解。 命令集81 广义特征值和广义特征向量 eig(A,B) 返回一个含有广义特征值的向量,A和B都是方阵。 [X,D]=eig(A,B) 返回一个对角线上为广义特征值的对角矩阵D和矩阵X, X的列是相对应的特征向量,因此有AX=BXD。 [X,v]= 给出度为k的特征问题(A+A+入A+..+2A)x=0的特征值和 polyeig(A0,A1,,.·,特征向量。向量v的长度为nk,包含有特征值:nXk的矩阵 AK) X的列是特征向量。如果有A=A和A,=一I,那么这就是标 准特征值问题。 为了检查特征值的条件或者它的敏感性,可以计算出条件数cod(X)=KX-‖,矩阵X 的列是A的特征向量。条件数大表示坏条件,也就是对扰动很敏感。 为了检查特征向量的条件或者它的敏感性可以查看特征值,多个重复的特征值或者特征 值彼此相差很小就表示是坏条件问题。 ■例8.3 假设: 3.75 -0.5 -0.375 0.495 -1.37\ 0.25 2.5 0.375 -0.495 -0.63 A= 1.25 -0.5 2.875 0.495 -2.12 0.25 -0.5 -0.625 2.505 0.37 0.25 -0.5 -0.625 0.495 2.38 运行命令[XX,DD]=eig(A)结果为:
X和D都是复数矩阵,运行命令: [ Y , E ] = c d f 2 r d f ( X , D ),结果为: 所得的矩阵E正好是A矩阵。 注意 特征向量是特征多项式d e t (lI-A) = 0的根,其中I是单位矩阵。用命令p o l y来求 特征多项式,参见11 . 1节。命令r o o t s也可求得特征值,但是用命令 e i g求得的特征 值更准确,精度更高。 广义特征值问题就是找到方程组 A x=lB x的重要解,其中B也是一个n×n的矩阵。l值和向 量x分别称为广义特征值和广义特征向量。 如果B是一个奇异矩阵,则用Q Z算法来求解。 标准和广义特征值问题都属于矩阵多项式特征问题,都可以用命令p o l y e i g来求它们的解。 命令集8 1 广义特征值和广义特征向量 e i g ( A , B ) 返回一个含有广义特征值的向量, A和B都是方阵。 [ X , D ] = e i g ( A , B ) 返回一个对角线上为广义特征值的对角矩阵 D和矩阵 X, X的列是相对应的特征向量,因此有 A X=B X D。 [ X , v ] = 给出度为k的特征问题(A0+lA1+lA2+ . . . +l kAk)x = 0的特征值和 p o l y e i g ( A 0 , A 1 , . . . , 特征向量。向量v的长度为n k,包含有特征值;n×nk的矩阵 A K ) X的列是特征向量。如果有 A0=A和A1=-I,那么这就是标 准特征值问题。 为了检查特征值的条件或者它的敏感性,可以计算出条件数 c o n d (X) = | |X|| ||X-1 | |,矩阵X 的列是A的特征向量。条件数大表示坏条件,也就是对扰动很敏感。 为了检查特征向量的条件或者它的敏感性可以查看特征值,多个重复的特征值或者特征 值彼此相差很小就表示是坏条件问题。 ■ 例8 . 3 假设: 运行命令[XX, DD]=eig(A),结果为: 1 1 6 M ATLAB 5 手册 下载 ■
China-bub.coM 第8章特征值和特征向量 117 下载 XX -0.0000 0.0000-0.44721 0.0000+0.4472i-0.4472 0.4472 0.5000 -0.0000+0.4472i -0.0000-0.4472i -0.4472-0.4472 0.5000 0.0000-0.4472i 0.0000+0.44721 -0.4472 0.4472 -0.5000 0.0000-0.44721 0.0000+0.4472i -0.4472-0.4472 -0.5000 0.0000-0.44721 0.0000+0.44721 -0.44720.4472 显然两个特征向量,列2和列3是复数。 DD 4.0000 0 0 0 0 03.0000+0.0000i 0 0 0 0 0 3.0000-0.00001 0 0 0 0 0 2.0000 0 0 0 0 02.0100 从上可以看出特征值为2,2.01,3,3和4,可以知道这个特征值问题是一个坏条件问题。 输入badMatrix=cond(xX)可求得条件数: badMatrix 5.0156e+07 将这个数和例8.2中矩阵的特征值条件数niceMatrix=cond(X)比较: niceMatrix 1.0000 它们是不同的。 ■ 8.2上海森伯形式、QR和QZ因式分解 如果只求特征值和特征向量,推荐用上一节中提到的方法。然而,有时要求更详细了解 计算过程,可用在这一节和下一节中定义的命令来满足这样的要求。 如果矩陶的第一子对角线下元素都是零,则它是一个上海森伯但essenberg)矩阵。如果矩阵是 对称矩阵,则它的海森伯形式是对角三角阵。MATLAB可以通过相似变换将矩阵变换成这种形式。 命令集82 上海森伯形式 hess(A) 返回矩阵A的上海森伯形式。 [P,H]=hess(A)返回一个酉矩阵P和上海森伯矩阵H,使A=PHP和PP'=I。 在MATLAB中,QR算法是计算矩阵所有特征值的一种有效的数学方法,也可以用命令 e1g来求。在用这种方法时,建议将矩阵转换成相似的上海森伯形式,参见例8.4。 QR算法是基于QR因式分解的一种算法,每个m×n的矩阵A可以表示成: A=QR 其中Q是一个m×m的酉矩阵,R是一个m×的上三角矩阵。如果A是一个方阵,R也还是 这样的一个矩阵。当用命令qr时,会返回矩阵Q和R,也可参见例7.7。 命令集83 QR因式分解 [Q,R]=gr(A) 产生一个m×m的酉矩阵Q和一个m×n的上三角矩阵
显然两个特征向量,列2和列3是复数。 从上可以看出特征值为 2,2 . 0 1,3,3和4,可以知道这个特征值问题是一个坏条件问题。 输入b a d M a t r i x = c o n d ( X X )可求得条件数: 将这个数和例8 . 2中矩阵的特征值条件数n i c e M a t r i x = c o n d ( X )比较: 它们是不同的。 8.2 上海森伯形式、Q R和Q Z因式分解 如果只求特征值和特征向量,推荐用上一节中提到的方法。然而,有时要求更详细了解 计算过程,可用在这一节和下一节中定义的命令来满足这样的要求。 如果矩阵H的第一子对角线下元素都是零,则它是一个上海森伯(H e s s e n b e r g)矩阵。如果矩阵是 对称矩阵,则它的海森伯形式是对角三角阵。M AT L A B可以通过相似变换将矩阵变换成这种形式。 命令集8 2 上海森伯形式 h e s s ( A ) 返回矩阵A的上海森伯形式。 [ P , H ] = h e s s ( A ) 返回一个酉矩阵P和上海森伯矩阵H,使A=P H P ´和P P´ =I。 在M AT L A B中,Q R算法是计算矩阵所有特征值的一种有效的数学方法,也可以用命令 e i g来求。在用这种方法时,建议将矩阵转换成相似的上海森伯形式,参见例 8 . 4。 Q R算法是基于Q R因式分解的一种算法,每个 m×n的矩阵A可以表示成: A = Q R 其中Q是一个m×m的酉矩阵,R是一个m×n的上三角矩阵。如果A是一个方阵,R也还是 这样的一个矩阵。当用命令 q r时,会返回矩阵Q和R,也可参见例7 . 7。 命令集8 3 Q R因式分解 [ Q , R ] = q r ( A ) 产生一个 m×m的酉矩阵 Q和一个 m×n的上三角矩阵 第8章 特征值和特征向量 1 1 7 下载 ■
118 MATLAB5手册 China-bub.com 下载 R,使得A=QR。 [Q,R,P]=qr (A) 产生一个大小为m×m、列正交的酉矩阵Q,一个对角 线元素递减的m×n的上三角矩阵R和一个置换矩阵P, 使得AP=QR。 [Q,R]=qrinsert 由于在矩阵A的列后插入一个额外的列b而得到新的 (Q,R,j,b) QR因式分解,Q和R是对矩阵A进行QR因式分解得到 的矩阵。如果=n+1,那么b就插入在矩阵A的最后一 列。 [Q,R]= 由于去掉矩阵A的第列而得到新的QR分解,Q和R是 qrdelete(Q,R,j) 对矩阵A进行QR分解得到的矩阵。 [Q1,R1]= 给出A+xy'的QR分解,也就是用秩为1的矩阵改变A的 grupdate(Q,R,x,y) QR分解。 如果A是上海森伯矩阵,则Q也是一样。对于QR算法,下面给出一些简短的描述: QR算法: 1)令A=A,k=0: 2)找到A的分解:A=QRe: 3)迭代计算下一个矩阵:A+=QR,令=k+1: 4)返回到2。 这种方法也称为不移位的Q方法,就是在某种约定下逼近于上三角矩阵。因为所有的矩 阵A和A。=A相似,所以有和原始矩阵相同的特征值,即最后的上三角矩阵的对角线元素就是 A的特征值。 如果矩阵一开始就转换成有接近一半元素是零的上海森伯形式,就可以减少可观的计算 步骤。QR方法作为MATLAB的一个内建函数,为了加快逼近速度也可进行移位。 ■例8.4 用不移位的QR因式分解算法,计算例8.1中矩阵A的特征值 -3-16 3 716 3 10 正确的特征值为入=10,入2=4和,=一6。 第1步: AO hess(A)[QO,RO]qr(AO);A1 RO*QO 返回得到: A1= 1.7992 26.8770-12.6126 2.3625 4.5085 -0.1434 0 4.9518 1.6923 第2步:[Q1,R1]=qr(A1):A2=R1*Q1,得到:
R,使得A=Q R。 [ Q , R , P ] = q r ( A ) 产生一个大小为 m×m、列正交的酉矩阵 Q,一个对角 线元素递减的 m×n的上三角矩阵 R和一个置换矩阵 P, 使得A P=Q R。 [ Q , R ] = q r i n s e r t 由于在矩阵A的j列后插入一个额外的列b而得到新的 ( Q , R , j , b ) Q R因式分解,Q和R是对矩阵A进行Q R因式分解得到 的矩阵。如果 j=n+ 1,那么b就插入在矩阵 A的最后一 列。 [ Q , R ] = 由于去掉矩阵A的第j列而得到新的Q R分解,Q和R是 q r d e l e t e ( Q , R , j ) 对矩阵A进行Q R分解得到的矩阵。 [ Q 1 , R 1 ]= 给出A+x y´的Q R分解,也就是用秩为1的矩阵改变A的 q r u p d a t e ( Q , R , x , y ) Q R分解。 如果A是上海森伯矩阵,则Q也是一样。对于Q R算法,下面给出一些简短的描述: Q R算法: 1) 令A0=A,k= 0; 2) 找到Ak的分解:Ak=QkRk; 3) 迭代计算下一个矩阵:Ak + 1=QkRk,令k=k+ 1; 4) 返回到2。 这种方法也称为不移位的 Q R方法,就是在某种约定下逼近于上三角矩阵。因为所有的矩 阵Ak和A0=A相似,所以有和原始矩阵相同的特征值,即最后的上三角矩阵的对角线元素就是 A的特征值。 如果矩阵一开始就转换成有接近一半元素是零的上海森伯形式,就可以减少可观的计算 步骤。Q R方法作为M AT L A B的一个内建函数,为了加快逼近速度也可进行移位。 ■ 例8 . 4 用不移位的Q R因式分解算法,计算例8 . 1中矩阵A的特征值 正确的特征值为l1= 1 0,l2= 4和l3=-6。 第1步: 返回得到: 第2步:[Q1, R1]=qr(A1) ;A 2 = R 1*Q 1,得到: 1 1 8 M ATLAB 5 手册 下载
China-pub.com 第8章特征值和特征向量 119 下载 A2= 17.607711.3432 5.0128 -15.3516-13.6557 -5.8721 0 1.0748 4.0480 在计算开始时,看不出要逼近什么样的矩阵。但是在计算了10步以后,就可看出主对角 线以下的元素较小。 [Q9,R9]=qr(A9);A10=R9*Q9 结果为: A10= 10.1297 22.6238 15.3505 -0.0924 -6.1616 -5.8036 0 0.0562 4.0319 注意,在整个计算过程中一直保留着上海森伯形式。 上例中的迭代过程可用MATLAB的内建编程语言简明地写出。在12.2节中有相关的例子。 QZ算法是用来计算复数矩阵的复数特征向量对和广义特征值的。在MATLAB中,依据下 面的命令集84来调用命令qz。 命令集84 QZ算法 [C,D,Q,Z,V]= 得到对角线元素是广义特征值的上三角矩阵C、D和广义特征向量矩阵V。 qz (A,B) 矩阵Q和Z是变换矩阵,使得QAZ=C和QBZ=D。 QZ方法是基于QZ分解的方法。 8.3 舒尔分解和奇异值分解 如果A是一个方阵,则有一个这样的酉矩阵U,使得: U-AU=U'AU=T 其中T是上三角矩阵。这是一个相似变换,因此矩阵A和T有相同的特征值。因为T是一个 对角矩阵,因此其对角线元素就是它的特征值。 如果A是实数矩阵且对称,那么T有对角形式,U的列就是A的特征向量。 如果A是实数矩阵,但是有复数特征值,那么T就是一个复数矩阵。为了避免复杂的计算, 用2×2的实数矩阵来代表每一对共轭复数特征值,参见例8.2。这样T就是一个块三角实数矩 阵。MATLAB中用命令schur来进行实数和复数矩阵A的舒尔分解。 如果A是实数矩阵,那么schur(A)返回实数的舒尔形式,但是如果A是复数矩阵,则给 出复数形式。它们的差别在于实数形式中在对角线上用2×2的矩阵块来表示共轭复数特征值 对,而复数形式给出的对角线元素是复数。函数rsf2csf可以将实数形式转换成复数形式。 命令集85 舒尔分解 schur(A) 给出矩阵A的舒尔分解,也就是如上所述的矩阵T。 [U,T]=schur(A) 给出矩阵A的舒尔分解和一个酉矩阵U,使得A=UTU'。 [V,S]=rsf2csf(U,T)将实数舒尔形式矩阵U和T转变成复数舒尔形式矩阵V和S
在计算开始时,看不出要逼近什么样的矩阵。但是在计算了 1 0步以后,就可看出主对角 线以下的元素较小。 结果为: 注意,在整个计算过程中一直保留着上海森伯形式。 上例中的迭代过程可用M AT L A B的内建编程语言简明地写出。在 1 2 . 2节中有相关的例子。 Q Z算法是用来计算复数矩阵的复数特征向量对和广义特征值的。在 M AT L A B中,依据下 面的命令集8 4来调用命令q z。 命令集8 4 Q Z算法 [ C , D , Q , Z , V ]= 得到对角线元素是广义特征值的上三角矩阵C、D和广义特征向量矩阵V。 q z ( A , B ) 矩阵Q和Z是变换矩阵,使得Q A Z=C和Q B Z=D。 Q Z方法是基于Q Z分解的方法。 8.3 舒尔分解和奇异值分解 如果A是一个方阵,则有一个这样的酉矩阵 U,使得: 其中T是上三角矩阵。这是一个相似变换,因此矩阵 A和T有相同的特征值。因为T是一个 对角矩阵,因此其对角线元素就是它的特征值。 如果A是实数矩阵且对称,那么T有对角形式,U的列就是A的特征向量。 如果A是实数矩阵,但是有复数特征值,那么 T就是一个复数矩阵。为了避免复杂的计算, 用2×2的实数矩阵来代表每一对共轭复数特征值,参见例 8 . 2。这样T就是一个块三角实数矩 阵。M AT L A B中用命令s c h u r来进行实数和复数矩阵A的舒尔分解。 如果A是实数矩阵,那么s c h u r ( A )返回实数的舒尔形式,但是如果 A是复数矩阵,则给 出复数形式。它们的差别在于实数形式中在对角线上用 2×2的矩阵块来表示共轭复数特征值 对,而复数形式给出的对角线元素是复数。函数 r s f 2 c s f可以将实数形式转换成复数形式。 命令集8 5 舒尔分解 s c h u r ( A ) 给出矩阵A的舒尔分解,也就是如上所述的矩阵 T。 [ U , T ] = s c h u r ( A ) 给出矩阵A的舒尔分解和一个酉矩阵U,使得A=U T U´。 [ V , S ] = r s f 2 c s f ( U , T ) 将实数舒尔形式矩阵U和T转变成复数舒尔形式矩阵V和S。 第8章 特征值和特征向量 1 1 9 下载 ■
120 MATLAB5手册 China-bub.com 下载 ■例8.5 将A1,A2,A3定义为: A1=(?)2=(62)A3=(-98) 运行的命令为: Sch1 schur(A1),Sch2 schur(A2)},... [U,Sch3]=schur(A3) 结果为: Sch1 1 0 0 3 Sch2 3 2 U= 1 0 0 1 Sch3 0 1 -1 0 因为有复数特征值,所以可以看出矩阵Sc3不是一个上三角矩阵。为了验证一下可以输 入:[V,S]=rsf2csf(U,Sch3),得到: V= 0+0.7071i0.7071 -0.7071 0-0.7071i S= 0+1.0000i 0 0 0-1.00001 这里的特征值是i和一i。 ■ MATLAB还可以计算奇异值分解,即SVD和矩阵的奇异值。这些数都是非负数,在某种 特定情况下,它们和矩阵的特征值相同。命令svds和命令eigs相似,也返回一些奇异值。 命令集86 SVD分解 svd(A) 返回一个包含矩阵A奇异值的向量。 [U,S,V]=svd(A) 返回一个对角矩阵S和大小分别为m×m和n×n的酉矩 阵U和V,矩阵S的大小和A相同,也是m×n的,而且 奇异值在矩阵的对角线上。奇异值是非负数且按降序
■ 例8 . 5 将A 1,A 2,A 3定义为: 运行的命令为: 结果为: 因为有复数特征值,所以可以看出矩阵 S c h 3不是一个上三角矩阵。为了验证一下可以输 入:[ V , S ] = r s f 2 c s f ( U , S c h 3 ),得到: 这里的特征值是i和-i。 M AT L A B还可以计算奇异值分解,即 S V D和矩阵的奇异值。这些数都是非负数,在某种 特定情况下,它们和矩阵的特征值相同。命令 s v d s和命令e i g s相似,也返回一些奇异值。 命令集8 6 S V D分解 s v d ( A ) 返回一个包含矩阵A奇异值的向量。 [ U , S , V ] = s v d ( A ) 返回一个对角矩阵 S和大小分别为 m×m和n×n的酉矩 阵U和V,矩阵S的大小和 A相同,也是 m×n的,而且 奇异值在矩阵的对角线上。奇异值是非负数且按降序 1 2 0 M ATLAB 5 手册 下载 ■
China-pub.coM 第8章特征值和特征向量 121 下载 排列。这些矩阵满足A=USV'和U'AV=S。 [U,S,V]=svd(A,0) 得到一个“有效大小”的分解,只计算出矩阵U的前n 列。矩阵S的大小为n义n。 svds (A,k,0) 计算出k个最大的奇异值和相应矩阵A的向量。如果k没 有给出,则缺省值为5。如果0作为最后一个参量值, 则计算最小值:否则计算最大值。 gsvd(A) 给出广义奇异分解,参见gsvd的帮助可得更多相关信息。 在MATLAB中,矩阵A的伪逆可以用命令piv(A)来求,也可用SVD分解来求矩阵的伪逆, 参见7.1节。 如果s,是奇异值,那么IA2=maxs=s1,IA-1l2=(mins)1=s1和cond(A)=s,/s,其中sn是最 小奇异值。对于非奇异矩阵和满秩的长方形矩阵(m>),最后一个表达式都成立。 ■例8.6 假设矩阵A,B: a-()s-(G8 (a)求奇异值分解: [Ua,Sa,Va]svd(A),[Ub,Sb,Vb]svd(B) Ua 0.3231 -0.8538 0.4082 0.5475 -0.1832 -0.8165 0.7719 0.4873 0.4082 Sa= 4.0791 0 0 0.6005 0 0 Va 0.4027 -0.9153 0.9153 0.4027 Ub -0.1641 0.9668 0.1961 0.5653 0.2551 -0.7845 0.8084 0.0178 0.5883 Sb 3.9116 0 0 0 1.3036 0 0 0 0
排列。这些矩阵满足A=U S V′和U′AV=S。 [ U , S , V ] = s v d ( A , 0 ) 得到一个“有效大小”的分解,只计算出矩阵 U的前n 列。矩阵S的大小为n×n。 s v d s ( A , k ,0) 计算出k个最大的奇异值和相应矩阵 A的向量。如果k没 有给出,则缺省值为 5。如果 0作为最后一个参量值, 则计算最小值;否则计算最大值。 g s v d ( A ) 给出广义奇异分解,参见g s v d的帮助可得更多相关信息。 在M AT L A B中,矩阵A的伪逆可以用命令p i n v ( A )来求,也可用S V D分解来求矩阵的伪逆, 参见7 . 1节。 如果si是奇异值,那么 和c o n d ( A ) =s1/sn,其中s n是最 小奇异值。对于非奇异矩阵和满秩的长方形矩阵 (m>n),最后一个表达式都成立。 ■ 例8 . 6 假设矩阵A,B: (a) 求奇异值分解: 第8章 特征值和特征向量 1 2 1 下载
122 China-bub.coM MATLAB5手册 下载 Vb= 0.3092 0.9510 0 0.9510 -0.3092 0 0 0 1.0000 可以看出矩阵A和B有两个非零的奇异值,也就是它们的秩为2。 (b)这些矩阵的逆是没有办法求的,但是伪逆可以用PseudoA=pinv(A)和 PseudoB=pinv(B)来求: PseudoA 1.3333 0.3333 -0.6667 -0.5000 0.0000 0.5000 PseudoB 0.6923 0.2308 0.0769 -0.2692 0.0769 0.1923 0 0 0 (c)矩阵A和B的范数和条件数为: normA norm(A),normB norm(B),... condA cond(A),condB =cond(B) 卫ormA= 4.0791 normB 3.9116 condA 6.7930 condB Inf 可以看出欧几里得范数等于最大奇异值,条件数等于s/s。 ■
可以看出矩阵A和B有两个非零的奇异值,也就是它们的秩为 2。 (b) 这 些 矩 阵 的 逆 是 没 有 办 法 求 的 , 但 是 伪 逆 可 以 用 P s e u d o A = p i n v ( A )和 P s e u d o B = p i n v ( B )来求: (c) 矩阵A和B的范数和条件数为: 可以看出欧几里得范数等于最大奇异值,条件数等于 s1/sn。 1 2 2 M ATLAB 5 手册 下载 ■