2020-04-25

第四章 方阵特征值和特征向量的计算

对 n 阶方阵 A,若存在常数 \lambda(可能是复数)和 n 维非零向量 x(分量可能是复数),满足 Ax=\lambda x,则称 \lambdaA 的一个特征值,称 xA 的对应于特征值 \lambda 的特征向量。
特征值和特征向量具有如下性质:

  • 特征方程 \mid (\lambda I-A)\mid=0;特征向量不唯一。
  • \lambda_1,\lambda_2,\cdots,\lambda_nA 的特征值,p(x)x 的某一多项式,则矩阵 p(A) 的特征值为 p(\lambda_1),p(\lambda_2),\cdots,p(\lambda_n)。特别地,A^k 的特征值为 \lambda_1^k,\lambda_2^k,\cdots,\lambda_n^k;若 A 可逆,则 1/\lambda_1,1/\lambda_2,\cdots,1/\lambda_nA^{-1} 的特征值,相应的特征向量保持不变。
  • A 为实对称矩阵,则 A 的所有特征值均为实数,不同特征值对应的特征向量相互正交;且存在正交矩阵 A,使得 Q^TAQ=diag(\lambda_1,\lambda_2,\cdots,\lambda_n),其中 Q^TQ=I,且 Q 的第 j 列是 \lambda_j 所对应的特征向量。
  • \mid p\mid \ne 0,B=P^{-1}AP,则称 AB 相似,相似矩阵具有相同的特征值。

4.1 乘幂法

  1. 乘幂法

用于求解矩阵按模最大特征值及其对应的特征向量。

设 n 阶矩阵 A 具有 n 个线性无关的特征向量 x_1,x_2,\cdots,x_n,相应的特征值 \lambda_1,\lambda_2,\cdots,\lambda_n 满足 |\lambda_1|>|\lambda_2|\ge|\lambda_3|\ge\cdots\ge|\lambda_n|,则对任取的非零向量 \upsilon_0,由 \upsilon_k=A\upsilon_{k-1}=\cdots=A^k\upsilon_0,k=1,2,\cdots 可产生一个向量序列 \{\upsilon_k\}

对上述向量序列 {\upsilon_k} 有:

  • 1)\underset{k\rightarrow\infty}{lim}\frac{\upsilon_k}{\lambda_1^k}=a_1x_1 (a_1为常数);
  • 2)\underset{k\rightarrow\infty}{lim}\frac{(\upsilon_{k+1})_m}{(\upsilon_k)_m}=\lambda_1.

其中 x_1\lambda_1 所对应的特征向量,(\upsilon_k)_m 表示 \upsilon_k 的第 m 个分量。

  1. 改进的乘幂法

\upsilon 为非零向量,将其规范化得到向量 u=\frac{\upsilon}{max(\upsilon)},其中 max(\upsilon) 表示向量 \upsilon 的模为最大的分量。

取初始向量 \upsilon_0\ne 0,规范化得到 u_0=\frac{\upsilon_0}{max(\upsilon_0)},构造向量序列:
\upsilon_1=Au_0=\frac{A\upsilon_0}{max(\upsilon_0)},u_1=\frac{\upsilon_1}{max(\upsilon_1)}=\frac{A\upsilon_0}{max(A\upsilon_0)}\\ \upsilon_2=Au_1=\frac{A^2\upsilon_0}{max(A\upsilon_0)}, u_2=\frac{\upsilon_2}{max(\upsilon_2)}=\frac{A^2\upsilon_0}{max(A^2\upsilon_0)}\\ \cdots\cdots\\ \upsilon_{k}=Au_{k-1}=\frac{A^k\upsilon_0}{max(A^{k-1}\upsilon_0)}, u_k=\frac{\upsilon_k}{max(\upsilon_k)}=\frac{A^k\upsilon_0}{max(A^k\upsilon_0)}
则有:
\underset{k\rightarrow\infty}{lim}u_k=\frac{x_1}{max(x_1)}, \underset{k\rightarrow\infty}{lim}max(\upsilon_k)=\lambda_1

常取初始向量 \upsilon=(1,1,\cdots,1)^T

例子
求矩阵A=\left(\begin{array}{lcr}1&-1&2\\-2&0&5\\6&-3&6\end{array}\right)

的按模最大特征值 \lambda_1 和对应的特征向量。

\begin{array}{l|lcr|lcr|r} k&\ &\upsilon_k^T &\ &\ &\upsilon_k^T &\ &\ \\ \hline 0&1.000&1.000&1.000&1.0000&1.0000&1.0000 &\ \\ \hline 1&2.000&3.000&9.000&0.2222&0.3333&1.0000&9\\ \hline 2&1.889&4.556&6.333&0.2982&0.7193&1.0000&6.333\\ \hline \vdots&\vdots&\vdots&\vdots&\vdots&\vdots &\vdots&\vdots\\ 10&1.393&4.444&5.013&0.2779&0.8865&1.0000&5.013\\ \hline 11&1.391&4.444&5.008&\ &\ &\ &5.008\\ \hline \end{array}

\lambda_1\approx5.008,x_1\approx(0.2779,0.8865,1)^T,而精确解为\lambda_1=5,x_1=(0.2778,0.8889,1)^T

  1. 反幂法

用于求矩阵 A 的按模最小特征值及其对应的特征向量。

用乘幂法求出 A^{-1} 的按模最大特征值和对应的特征向量,就可以获得 A 的按模最小特征值和对应的特征向量。

任取初始非零向量 \upsilon_0,构造向量序列:
u_0=\frac{\upsilon_0}{max(\upsilon_0)},\upsilon_k=A^{-1}u_{k-1},u_k=\frac{\upsilon_k}{max(\upsilon_k)},k=1,2,\cdots

则有:\underset{k\rightarrow\infty}{lim}u_k=\frac{x_n}{max(x_n)},\underset{k\rightarrow\infty}{lim}max(\upsilon_k)=\frac{1}{\lambda_n}

其中,方程组 A\upsilon_k=u_{k-1} 可通过矩阵三角分解的方法进行求解,进而得到向量序列 \{\upsilon_k\}

例子
求矩阵A=\left(\begin{array}{lcr}1&-1&2\\-2&0&5\\6&-3&6\end{array}\right)

的按模最大特征值 \lambda_1 和对应的特征向量。


第一步:做 ALU 分解,得到:
L=\left(\begin{array}{lcr}1&0&0\\-2&1&0\\6&-1.5&1\end{array}\right), U=\left(\begin{array}{lcr}1&-1&2\\0&-2&9\\0&0&7.5\end{array}\right)

第二步,迭代求解:

\begin{array}{l|lcr|r} k&\ &u_{k-1}^T &\ &\lambda\\ \hline 1&1.000&1.000&1.000&-0.5556\\ \hline 2&-0.3704&-1.000&-0.037006&-1.626\\ \hline \vdots&\ &\vdots &\ &\vdots\\ \hline 9&0.5003&1.000&-0.0004924&-0.9990\\ \hline 10&-0.5000&-1.000&0.00001854&-1.000 \end{array}

即: \lambda_3\approx-1.000,x_3\approx(-0.5,-1,-0.0001854)^T
该问题的精确解为: \lambda_3=-1,x_3=(-0.5,-1,0)^T

4.2 Jacobi 方法

Jacobi 方法用于求解实对称矩阵的所有特征值和对应的特征向量。

  1. 平面旋转矩阵

对平面上的双曲线 xy=1 可作正交相似变换
\left(\begin{array}{c} x\\ y \end{array} \right) =\left(\begin{array}{lr} cos \frac{\pi}{4}&-sin \frac{\pi}{4}\\ sin \frac{\pi}{4}&cos \frac{\pi}{4} \end{array}\right) \left(\begin{array}{c}u\\ \upsilon \end{array}\right)

可将其约化为实轴和虚轴均在坐标轴上的标准形式: u^2-\upsilon^2=(\sqrt 2)^2.

对一个实对称矩阵实施正交相似变换可以将其约化为对角矩阵。

Jacobi 方法就是利用一系列特殊的正交相似变换,逐次将实对称矩阵约化为对角矩阵。

得到的对角矩阵的主对角元素就是所求的特征值,而用于正交相似变换的正交矩阵的各列就是相应的特征向量。

平面旋转矩阵:
S_k=S(p,q)= \left(\begin{array}{lccccr} 1&\ &\ &\ &\ &\ \\ \ &\ddots&\ &\ &\ &\ \\ \ &\ &cos\theta_k&\cdots&sin\theta_k&\ \\ \ &\ &\vdots&\ &\vdots&\ \\ \ &\ &-sin\theta_k&\cdots&cos\theta_k\\ \ &\ &\ &\ &\ &\ddots &\ \\ \ &\ &\ &\ &\ &\ &1 \end{array}\right)

平面旋转矩阵也称为 Gicens 矩阵,其几何意义是由其定义的线性变换,使得 n 维空间的第 p 个坐标轴和第 q 个坐标轴所构成的坐标平面旋转了 \theta_k 的角度。平面旋转矩阵是正交矩阵,并且变换 S(p,q)^TAS(p,q) 只改变了矩阵A 的第 p 行、第 p 列,和第 q 行、第 q 列的元素。

  1. 古典 Jacobi 方法

例子
用古典 Jacobi 方法求矩阵 A 的全部特征值和特征向量,误差限为 0.0005.

A=\left(\begin{array}{lcr} 1.0&1.0&0.50\\ 1.0&1.0&0.25\\ 0.5&0.25&2.0 \end{array}\right)


第一步:将矩阵 A 中的 a_{12}=a_{21} 化为 0,有

S_1=\left(\begin{array}{lcr} 0.70711&0.70711&0\\ -0.70711&0.70711&0\\ 0&0&1 \end{array}\right),\\ T_1=S_1^TAS_1=\left(\begin{array}{lcr} 0&0&0.17678\\ 0&2&0.53038\\ 0.17678&0.53038&2 \end{array}\right)

第二步,将矩阵 T_1 中的 t_{23}^{(1)}=t_{32}^{(1)}化为 0,有

S_2=\left(\begin{array}{lcr} 1&0&0\\ 0&0.70711&0.70711\\ 0&-0.70711&0.70711 \end{array}\right),\\ T_2=S_2^TT_1S_2=\left(\begin{array}{lcr} 0&-0.12500&0.12500\\ -0.12500&1.4697&0\\ 0.12500&0&2.5304 \end{array}\right)

第三步,将矩阵 T_2 中的 t_{12}^{(2)}=t_{21}^{(2)}化为 0,有

S_3=\left(\begin{array}{lcr} 0.99645&-0.84145&0\\ 0.84145&0.99645&0\\ 0&0&1 \end{array}\right),\\ T_3=S_3^TT21S_3=\left(\begin{array}{lcr} -0.010556&0&0.12456\\ 0&1.4802&-0.010518\\ 0.12456&-0.010518&2.5304 \end{array}\right)

重复上述过程有:

T_5=\left(\begin{array}{lcr} -0.016647&0.000409&0.000005\\ 0.000409&1.4801&0\\ 0.000005&0&2.5365 \end{array}\right),\\ Q=S_1S_2S_3S_4S_5=\left(\begin{array}{lcr} 0.72135&0.44404&0.53584\\ -0.68616&0.56234&0.46147\\ -0.093844&-0.69757&0.71033 \end{array}\right)

由此得到 A 的误差限为 0.05% 的三个近似特征值为:\lambda_1=-0.16647,\lambda_2=1.4801,\lambda_3=2.5365

对应的近似特征向量为:
x_1=(0.72135,-0.68616,-0.093844)^T,x_2=(0.44404,0.56234,-0.69757)^T,x_3=(0.53584,0.46147,0.71033)^T

  1. Jacobi 过关法

4.3 QR 方法

QR 方法可用于求解一般矩阵的全部特征值。

  1. Householder 变换

\upsilon 是 n 维列向量,且 \upsilon^T\upsilon=1,称矩阵 H=I-2\upsilon\upsilon^T 为 Householder 矩阵,也称为镜像变换矩阵。

  • 对 n 维向量 x,y,若\Vert x\Vert_2=\Vert y\Vert_2,则存在镜像矩阵 H,使得 y=Hx
  1. LR 分解

LR 方法是目前求解矩阵所有特征值的最有效方法。

A=A_1 作分解 A_1=F_1R_1,其中 F_1 非奇异,反序相乘有 A_2=R_1F_1。又作分解 A_2=F_2R_2,其中 F_2 非奇异,反序相乘有 A_3=R_2F_2

产生矩阵序列 \{A_k\},满足:
\begin{cases} A_1=A,\\ A_k=F_kR_k,\ \ \ \ \ k=1,2,\cdots \\ A_{k+1}=R_kF_k, \end{cases}

A_{k+1}=F_k^{-1}A_kF_k=(F_1F_2\cdots F_k)^{-1}A(F_1F_2\cdots F_k) 可知,矩阵序列 \{A_k\} 是相似矩阵序列,因此它们具有相同的特征值。

例子
用 LR 方法求对称正定矩阵的所有特征值。

A=\left(\begin{array}{lcr}7&3&1\\3&4&2\\1&2&3\end{array}\right)


对 A 使用 Cholesky 分解得到:
L_1=\left(\begin{array}{lcr}2.64571&0&0\\1.133893&1.647510&0\\0.3779645&0.9538209&1.395481\end{array}\right)

反序相乘得到:
A_2=L_1^TL_1=\left(\begin{array}{lcr} 8.428568&2.228609&0.5274423\\2.228609&3.624061&1.331038\\0.5274423&1.331038&1.947367\end{array}\right)

若干次分解并反序相乘得到:
A_12=\left(\begin{array}{lcr} 9.433488&0.0160184&0.0000184\\ 0.0160184&3.4194300&0.0062214\\ 0.0000184&0.0062214&1.1470350\end{array}\right)

  1. QR分解
©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 202,723评论 5 476
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 85,080评论 2 379
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 149,604评论 0 335
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 54,440评论 1 273
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 63,431评论 5 364
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,499评论 1 281
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,893评论 3 395
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,541评论 0 256
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,751评论 1 296
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,547评论 2 319
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,619评论 1 329
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,320评论 4 318
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,890评论 3 307
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,896评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 31,137评论 1 259
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 42,796评论 2 349
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 42,335评论 2 342

推荐阅读更多精彩内容