数值分析矩阵的特征值

数值分析矩阵的特征值矩阵的特征值与特征向量的计算4-1第1页，课件共101页，创作于2023年2月第九章

矩阵的特征值

与特征向量的计算§1乘幂法与反幂法§2Jacobi方法

§3QR方法

2第2页，课件共101页，创作于2023年2月第九章矩阵的特征值与特征向量的计算概述线性代数中已接触过这个概念，矩阵的特征值与特征向量能反映矩阵的性态，在理论上重要，而工程技术中的许多问题，如桥梁的振动，机械的振动，建筑物的振动及飞机机翼的颤动，这些问题的求解常归结为求矩阵的特征值问题，另外一些稳定性分析问题也会转化为求特征值与特征向量的问题。先简单提一下有关概念n阶方阵A的特征值是这样的：它使Ax=x，即有(A-I)x=0（其中特征向量x为非零向量），将其看作方程组，则是n个未知数（为n阶向量），n个方程的齐次线性方程组，它有非零解x的充分必要条件是系数矩阵所成行列式，称为关于的特征方程：3第3页，课件共101页，创作于2023年2月矩阵的特征值与特征向量的计算概述（续1）将行列式展开，得关于的n次多项式：

()称为A的特征多项式，一般()=0的n个根都是A的特征值，对应于n个特征向量，且若x为对应的特征向量，则kx也是对应的特征向量（允许相差一常数k）。当n不大时，如n4

解特征方程，可求出全部特征值（n3较难）当n较大（n>5），计算量会增大得惊人，且不可能求得准确结果，还可能出现不稳定，所以当n稍大一般不直接求解特征方程，而根据实际问题的需要，介绍相应的一些行之有效的数值解法

4第4页，课件共101页，创作于2023年2月§1

乘幂法与反幂法

像我们在范数，谱半径中知道的那样，有时只需求出矩阵的按模最大的特征值与相应的特征向量。

1.1乘幂法

乘幂法是一种迭代法；先找初始向量x(0)反复乘以给定的方阵A，依次得x(1),x(2)……，直到满足精度为止，可从中分离出绝对值最大的特征值。

设nn阶实矩阵A的特征值i

(i=1,2,…，n)满足

5第5页，课件共101页，创作于2023年2月构造向量序列假定i(i=1,2,…，n)对应的特征向量U1,U2,…,Un线性无关，这组特征向量构成Rn中的一个基底，则任一非零向量可表为Ui(i=1,2,…，n)的线性组合，即存在n个不全为o的常数i(i=1,2,…，n)使

可构造向量序列

所以乘幂法实际上是，对于给定的初始向量（零向量）由迭代法：6第6页，课件共101页，创作于2023年2月求出最大的特征值1产生迭代向量序列，可以证明，当k充分大时有

（由x的某一分量的相邻二次结果之比可得出1），而相应的特征向量为。实际上，由式（4-1）可得

：7第7页，课件共101页，创作于2023年2月求出最大的特征值1（续1）因此可看作为与对应的特征向量的近似，而由式，

取它们的任一分作比值即可得：

8第8页，课件共101页，创作于2023年2月几点注释注2：当而k充分大时，会随k的增大而无限增大或无限趋于0，这样上机计算时会产生溢出（上溢或下溢）的情况，为避免这种情形出现，实际计算时，每次迭代求得的向量x(k)要进行归一化（规范化）处理：取x(k)中绝对值最大的一个分量除x(k)，这样将x(k)的分量全部控制在[-1,1]中，而1是由相邻二次分量的比值所决定，因此不会受到影响。注1：一般有10，若恰好x(0)使1为0，也不影响上述法，因为实际计算中，由于有舍入误差的影响，迭代n次后所得到的向量x(k)在u1方向上的分量不会为0，因此，可得x(2)为初始向量。继续下去。9第9页，课件共101页，创作于2023年2月几点注释（续）

注3：因此乘幂法实际使用的公式及算法为：

10第10页，课件共101页，创作于2023年2月

算法4.1

4.计算

，，。

注4：初始向量x(0)的选取对迭代有影响，若收敛太慢可考虑重新选初值。6.若k<N，置k+1k,

，转3；否则，输出失败信息，停机。11第11页，课件共101页，创作于2023年2月例题1例1用乘幂法求按模最大的特征值和特征向量，取x(0)=(0，0，1)T，要求误差不超过10-3。解由式（4-3）

12第12页，课件共101页，创作于2023年2月例题1（续1）如此继续下去，计算结果见表4-1k000110011012200.5120.522.52.50.20.8131.22.62.82.80.42857140.9285714141.78571422.85714282.9285714

0.60975600.9756097152.19512182.9513942.9756097

0.73770480.9918619162.46727172.98372382.9918619

0.82466090.9972799172.64660182.99455982.9972799

0.88300120.9990924182.76509482.99818482.99909024

0.92197720.9996973192.84365172.99939462.9996973

13第13页，课件共101页，创作于2023年2月例题1（续2）因为 ∴ 其对应的特征向量：

可与准确的特值值；1=3，2=2，3=1及1的特征向量（1,-1,1）T相比较，再次说明：特征向量允许相差一常数。14第14页，课件共101页，创作于2023年2月两点注释

注1：幂法的收敛速度依赖于，比值越小，收敛越快，对上例，按准确的值，比值不是很小。因此对=103也迭代了9次，不算收敛快。对按乘幂法计算其最大特征值，取x(0)=(1,1,1)T。因为，所以对于

=103，只需迭代5次。当收敛慢时，还要考虑加快收敛的技术。15第15页，课件共101页，创作于2023年2月两点注释（续1）注2：前面假定，严格大于其他特征值，也有可能，这样的1有多个，如m个，那么上述幂法还行吗？讨论特殊情况如下：（1）1是m重根。即1=2=…=m，幂法仍有效，此时有，式（4-2）变成：

只要1,2,…,m不全为零，当k充分大时

16第16页，课件共101页，创作于2023年2月两点注释（续2）x(k+1)为1对应的特征向量收敛到1u1+…+mum17第17页，课件共101页，创作于2023年2月两点注释（续3）表明x(2k)是1对应的特征向量的近似。

式（4-4），（4-5）为最大的特征值1,2与对应的特征向量的计算公式。18第18页，课件共101页，创作于2023年2月1.2幂法的加速：原点移位法

当比值接近于1，则幂法收敛慢，应配以加速技术。下面介绍两种加速技术：

所谓原点移位法是：以A-0I代替A，用乘幂法求得A-0I的最大特征值i，则A的最大特征值1=i+0，其特征向量不变。而对A-0I按乘幂法（4-1）有：适当地选取0满足且这样在用幂法求A-0I的最大特征值1-0时，收敛速度比对A用幂法要快。1.原点移位法19第19页，课件共101页，创作于2023年2月0的估计

原点移位法较简单，然而0如何选取是很因难的，一般情况下，如果对特征值的分布情况有大概的了解，可粗略地估计出0。

此不等式表明，原点移位法求1，加快了收敛速度。20第20页，课件共101页，创作于2023年2月例题2

取0=2.9，用原点移位法求A的最大特征值，要求误差<10-4。按式（4-3）计算，结果见表4-2

21第21页，课件共101页，创作于2023年2月例题2（续1）表4-2k0111111117.15.11.17.110.71830980.154929523.15633722.2549290.98450713.156337210.71441320.311914433.1017852.21557330.96880863.10178510.71428970.31233943.10005682.2143260.96876613.100056810.71428560.312499453.09999842.21428460.96875013.09999841

22第22页，课件共101页，创作于2023年2月例题2（续2）所以，矩阵A的按模最大的特征值为不难求出，A的特征值为，若对A直接用幂法，则比值，而用原点移位法，则有：因此收敛速度明显加快。23第23页，课件共101页，创作于2023年2月2.幂法的埃特肯（Aitken）加速

1）首先介绍Aitken加速法，并且说明对线性收敛速度的迭代均可用此法加快收敛；2）进一步说明幂法是线性收敛的；3）Aitken加速法用于幂法，加速收敛。下面是详细内容：1）Aitken加速法

设序列{ak}线性收敛到a,记误差k=aka即有于是当k充分大，且k，k+1同号时有

：可解出：24第24页，课件共101页，创作于2023年2月埃特肯（Aitken）加速（续）

因此所有具线性的收敛速度的序列均可使用Aitken法加快收敛

。上面结果表明：分子超于0的速度比分母超于0的速度快，亦即分子是比分母更高价的无穷小，因此比{ak}快。25第25页，课件共101页，创作于2023年2月算法4.2

2）乘幂法收敛速度是线性的，即乘幂法所得向量序列

{x(k+1)}收敛到x1的速度是线性的。

3）Aitken加速技术用于乘幂法，得如下算法4.2

1.输入A=(aij)，初始向量x=(x1,x2,…xn)T，误差限，最大迭代次数N；3.

4.计， 置26第26页，课件共101页，创作于2023年2月例题3例3用Aitken加速法求例1中矩阵的按模最大特征值与对应的特征向量，取x(0)=(0,0,1)T。

解（一）2.

， ，，转327第27页，课件共101页，创作于2023年2月例题3（续1）（二）3.

（三）3。28第28页，课件共101页，创作于2023年2月例题3（续2）（四）3.计算结果见表4-3：29第29页，课件共101页，创作于2023年2月例题3（续3）k10122

20.522.52.5

31.22.62.82.83.2541.78571422.85714282.92857142.92857143.024999952.19512182.9513942.97560972.97560973.002747262.4672172.98372382.99186192.99186193.000441630第30页，课件共101页，创作于2023年2月例题3（续4）

计算结果与例1的计算结果比较，显然Aitken加速法的收敛速度比幂法快（这里加速后的收敛快）。

也可上述算法4.2是用幂法迭代一次，加速一次，修改为用幂法迭代二次或三次，加速一次。

Aitken加速序列的收敛速度是原序列收敛速度的两倍。31第31页，课件共101页，创作于2023年2月1.3反幂法

在乘幂法中，以A-1代替A，即为反幂法，用于求A的最小特征值及对应的特征向量。

因为以A-1代替A，由此式表明A-1的特征值是A的特征值的倒数，而相应特征向量不变。这样，若对A用乘幂法求最大特征值1

：则对A-1用乘幂法求最大特征值应满足：32第32页，课件共101页，创作于2023年2月反幂法（续1）若按乘幂法，以A-1代替A，有x(k+1)=A-1x(k)，为避免求逆阵A-1，实际运算时，常化为求解Ax(k+1)=x(k)，这样迭代一次的算法，转化为求解一次方程组。由于矩阵A在迭代过程中不会改变，因此可先对其进行分解A=LU，于是在每次迭代时，就只需转化为求解两个三角方程组：也就是说，对A-1用乘幂法可计算A-1的最大特征值，而实际上是A的最小特征值。

33第33页，课件共101页，创作于2023年2月反幂法（续2）

注2：反幂法通常用于：已知矩阵的某一个近似特征值时，求其对应的特征向量，并对此近似特征值加以修正，且与原点移位法合起来用。

注1：反幂法的收敛率为比值（假定求的最大：，而在A中最小）：收敛速度与 收敛率有关。34第34页，课件共101页，创作于2023年2月注3

带原点移位的反幂法设已知A的一个特征值的近似值，一般很小，即，这表明是矩阵（原点移位法中0=*）的按模最小的特征值可对用反幂法，求最小特征值。由于实际上收敛率较小所以此时收敛很快。往往只需要二，三步迭代，就可达到很高的精度。

算法4.3

1.输入A=(aij)，初始向量x=(x1,x2,…xn)T，误差限，最大迭代次数N；

35第35页，课件共101页，创作于2023年2月算法4.3（续）3.作三角分解

4.求整数r,使得；

7.若，则置，，转4，否则输出失败信息，停机。

36第36页，课件共101页，创作于2023年2月注4

注4若已知A的全部特征值的近似值，用上述反幂法可求A的全部特征向量，同时可对这些近似特征值加以修正改正。

说明：如果*=0，则算法4.3求出A的按模最小的特征值与特征向量。

例4用反幂法求矩阵

接近2.93的特征值，并求相应的特征向量，取x(0)=(0,0,1)T。

37第37页，课件共101页，创作于2023年2月例4（续1）

[解]对A2.93I作三角分解得：按算法4.3计算的数值结果见表4-438第38页，课件共101页，创作于2023年2月例4（续2）表4-4k10017.95905957.40192546.88379067.95905953.0752688210.931.864912.69231312.80385112.83758112.8375813.008787830.98868410.99737252.072443614.27843614.26762914.26626814.2784363.000095439第39页，课件共101页，创作于2023年2月例4（续3）

迭代3次，即得3.0000954，与准确值=3的误差小于10-4。相应的特征向量为：与准确值u=(1,-1,1)T比较，残量为：40第40页，课件共101页，创作于2023年2月§2Jacobi方法序Jacobi法是用于求实对称阵的全部特征值及特征向量的一种迭代法，其基本思想是：利用一系列的正交相似变换阵Q，可将n阶实对称阵A化为一对角阵，这些对角元就是A的特征值，而利用这些正交变换阵的乘积，可求出对应的特征向量。先介绍一些预备知识和相关结论：

（1）A为n阶实对称阵，则A有n个实特征值，且有n个对应的正交的线性无关的特征向量。

Q为正交阵，则：Q

QT=I，即QT=Q1，正交阵的乘积仍为正交阵。

（2）

A，B为n阶方阵，方阵R可逆，指R使：AR=RA=I。

A，B相似：存在可逆方阵R使，R1AR=B。相似矩阵的特征值相同，A，B相似，若A为对称阵，则B也是对称阵。

41第41页，课件共101页，创作于2023年2月Jacobi方法序（续1）（3）A为实对称阵，Q为正交阵，则B=QTAQ也是对称阵。（4）对于n阶实对称阵A，一定存在正交阵Q，使：

即A，D相似，由（2）知i(i=1,2,…,n)就是A的特征值，而Q的第i列qi(i=1,2,…,n)就是A的，对应于i的特征向量,

因为：：42第42页，课件共101页，创作于2023年2月Jacobi方法序（续2）（5）在正交相似变换下，矩阵元素的平方和不变。

Jacobi方法正是建立在上述结论上的：通过一次正交变换，将A中的两个非零的非对角线元素化为零，使得非对角元素的平方和减少，反复进行上述过程，使变换后的矩阵的非对角元素的平方和趋于零，从而使该矩阵近似为对角矩阵，得到全部特征值和特征向量。综上所述，Jacobi法的关键，在于找到合适的正交阵Q，为了说明这个问题，我们从最简单的情况开始。则43第43页，课件共101页，创作于2023年2月

2.1平面旋转变换

设在平面上有一条二次曲线：

可以通过坐标轴的旋转

化为标准形状：

如果把（4-7）式写成矩阵形式，即为

44第44页，课件共101页，创作于2023年2月平面旋转变换（续1）其中a12=a21，而（4-8）式即为：

把（4-10）式代入（4-9）式，便有：

45第45页，课件共101页，创作于2023年2月平面旋转变换（续2）

则上式可简写为

两个非对角线元素已化为零46第46页，课件共101页，创作于2023年2月平面旋转变换（例5）容易验证Q是一个正交矩阵，称（4-7）式是一个正交变换。正交变换Q把对称矩阵A变成为对角阵，而1与2即是矩阵的特征值。正交矩阵Q的两个列向量分别对应于1,2的两个单位特征向量，即1所对应的特征向量为，2所对应的特征向量为。为了把上述结果推广到一般情况，我们再用一个具体例子来说明。

例5椭球

与坐标平面Ox1x2的交线是：如果把Ox1,Ox2轴旋转，则可知此二次曲线是一个椭圆47第47页，课件共101页，创作于2023年2月平面旋转变换例5（续1）为此令变换：（4-11）式经过此变换以后，得新方程为 把（4-11）和（4-13）式均写成矩阵形式可知经过变换以后，使矩阵：48第48页，课件共101页，创作于2023年2月平面旋转变换例5（续2）或完整地写为：

49第49页，课件共101页，创作于2023年2月平面旋转变换例5（续3）下面我们来考察经过这个变换以后，（4-11）中矩阵元素的变化情况。1

对角线上元素的平方和由19增加到27。2

非对角线上元素的平方和由减少到，而矩阵所有元素的平方和不变。由于（4-13）式中，仍然保留着y1y3与y2y3的乘积项，如果仿上再作一次变换，例如把Oy2y3平面与（4-13）式截口化成标准形，可以作如下旋转变换：50第50页，课件共101页，创作于2023年2月平面旋转变换例5（续4）这时椭球方程化为：

写成矩阵形式，得系数矩阵：这时对角线上元素的平方和达到，而非对角线元素的平方和减少成。51第51页，课件共101页，创作于2023年2月平面旋转变换例5（续5）综合上述两次变换的结果，我们可以得出如下结论：

（1）实对称矩阵A经过正交变换以后，对角线上的元素的平方和在不断增加，非对角线上的元素的平方和在不断减少，而矩阵所有元素的平方和是不改变的。（2）同时也看到经过第二次变换后，上一次变换时已经化为零的元素又可能会变成不是零。但不管怎样，经过这样反复变换，总可达到目的：使对角线上元素的平方和不断增大而非对角线上元素平方和趋向于零。

这个例子的具体作法，就是Jacob法的基本思想。52第52页，课件共101页，创作于2023年2月Jacob法

设阶实对称矩阵A=(aij)。经过正交变换以后，得：又设A中以aij(ij)的绝对值为最大，当然aij0（否则非对角线上所有元素全为零），作正交变换：53第53页，课件共101页，创作于2023年2月Jacob法（续1）矩阵Vij()中，对角线上元素除Vii=Vjj=cos

外，其它皆为1，非对角线上元素除-Vii=Vjj=sin

外，其它皆为0，称Vij()为平面旋转阵。容易直接验证Vij()有如下性质：1

.

2

如果A是对称阵，则，所以是对称阵，就是说，对称阵经过正交变换以后，仍是对称阵。

3

矩阵A经过变换后，A1中第i行，第j行及第i列、第j列元素的变化如下：

54第54页，课件共101页，创作于2023年2月Jacob法（续2）

其它元素不变，即：55第55页，课件共101页，创作于2023年2月Jacob法（续3）

同样可以直接验证：1） 即经过正交变换后，矩阵所有素的平方和不变。56第56页，课件共101页，创作于2023年2月Jacob法（续4）

利用条件（4-18）得：

对A(1)重复上述过程，可得A(2)，如此继续下去，得到一个矩阵序列{A(k)}。可以证明，虽然这种变换不一定能使矩阵中非对角元中零元素的个数单调增加，但可以保证非对角元的平方和递减，我们以A与A(1)为例进行讨论。

由此可知，经过正交变换后矩阵A中对角线元素的平方和比原来增加了，而非对角线元素的平方和减少了。57第57页，课件共101页，创作于2023年2月Jacob法（续5）

可得：（4-19）式：

58第58页，课件共101页，创作于2023年2月Jacob法（续6）

式（4-19）表明，在上述旋转变换下，非对角元的平方和严格单调递减，因而由式（4-6）即得，对角元的平方和单调递增。正是利用这一点，导出了Jacobi方法。2.2Jacobi方法

通过一系列旋转相似变换将A变成A(k+1)，求得A的全部特征值与特征向量的方法称为Jacobi方法。如果在对A作相似变换的过程中，每一步都选绝对值最大的非对角元素，以此确定旋转矩阵，这种方法称为古典Jacobi方法。59第59页，课件共101页，创作于2023年2月古典Jacob法（3）计算旋转矩阵：

（2）求整数i,j,使得：

古典Jacobi方法。其计算过程如下：60第60页，课件共101页，创作于2023年2月古典Jacob法（续1）61第61页，课件共101页，创作于2023年2月古典Jacob法（续2）一般地，Jacobi法不能在有限步内将A化成对角阵，但有以下收敛性结果。

定理4.1

设A为n阶实对称阵，对A用古典Jacobi法得到序列{A(k)}，其中A(0)=A，则：即古典Jacobi法收敛。

[证明]由Jacobi法计算过程：另一方面，由计算A(k+1)的公式可以得出：62第62页，课件共101页，创作于2023年2月古典Jacob法（续3）

定理4.1表明，古典Jacobi法是收敛的，进一步分析还可以得出Jacobi法收敛较快。另外，这种方法对舍入误差有较强的稳定性，因而解的精度高，且所求得的特征向量正交性很好。它的不足之处是运算量大，且不能保持矩阵的特殊形状（如稀疏性），因此Jacobi法是求中小型稠密实对称矩阵的全部特征值与特征向量的较好方法。63第63页，课件共101页，创作于2023年2月两点说明

在实际计算中还可采用一些措施来提高精度和节省工作量。1

减少舍入误差的影响。从公式中可知，具体计算时只需用到sin

,cos的值，为了提高精度，舍入误差越小越好。常常利用三角函数之间的关系，写成便于计算的公式

即得 （4-21）64第64页，课件共101页，创作于2023年2月两点说明（续）

2

节省工作时间。在雅可比法中，每次变换是把非对角元绝对值最大者化为零，但在n阶矩阵中要去寻找这个最大元素要花较多的机器时间，所以一般不选最大元。改进的一种方法是设某些“关口”，如a1,a2,…,ak，先按次序用aij(ij,j=1,2,…,n)与a比较，若，则通过不加运算，若，就进行一次旋转变换，使之化为0，一遍轮流过以后，再用a2来比较，作同样处理，直至达到所需精度为止，这种方法称为雅可比过关法。例6用Jacobi方法求下列矩阵的特征值。65第65页，课件共101页，创作于2023年2月例6（续1）66第66页，课件共101页，创作于2023年2月例6（续2）下面应取i=2,

j=3，重复上述过程。如此继续下去，可得：67第67页，课件共101页，创作于2023年2月例6（续3）所以A的特征值为：

与其准确值比较，最大误差为0.0002036。68第68页，课件共101页，创作于2023年2月古典Jacobi法一种改进古典Jacobi法在计算每个旋转矩阵V(k)前都需要对个非对角元素进行比较，从中找出绝对值最大的元素。为减少运算量常用的一种改进方法是取定正{k},k0(k)，以k为限，逐行检查非对角元，若就跳过，否则以Vij()消去元aij和aji，反复进行上述过程，直到所有非对角元的绝对值均小于k，再以k+1为限，进行第k+1轮循环消元。当k充分小时，所得到的矩阵的对角元即为A的全部特征值。69第69页，课件共101页，创作于2023年2月§3QR方法3.1基本QR方法六十年代出现的QR算法是目前计算一般中小型矩阵的全部特征值与特征向量的最有效的方法。这里仅讨论实矩阵，并假定矩阵非奇异。因为否则，矩阵AI（不是A的特征值）必定是非奇异的，而由AI

的特征值与特征向量容易得到A的特征值与特征向量。因为任一非奇异实矩阵A都可以分解成一个正交矩阵Q和一个上三角矩阵R的乘积，而且当R的对角元符号取定时，分解是唯一的。基本QR方法的基本思想是利用矩阵的QR分解，通过迭代格式：70第70页，课件共101页，创作于2023年2月QR方法（续1）可将A=A(1)化成相似的上三角阵（或分块上三角阵），从而求出矩阵A的全部特征值与特征量。

即A(2)与A相似。同理可得，A(k)~A(k=2,3,…)。故它们有相同的特征值。

可以证明，在一定条件下，基本QR方法产生的矩阵序列{A(k)}

“基本”收敛于一个上三角矩阵（或分块上三角阵），即主对角线（或主对角线子块）及其以下元素均收敛，主对角线（或主对角线子块）以上元素可以不收敛。特别地，如果A是实对称阵，则{A(k)}

“基本”收敛于对角矩阵。

71第71页，课件共101页，创作于2023年2月QR方法（续2）

因为上三角阵的主对角元（或分块上三角阵中，主对角子块的特征值）即为该矩阵的特征值，故当K充分大时，A(k)的主对角元（或主对角子块的特征值）就可以作为A的特征值的近似。

基本QR方法的主要运算是对矩阵作QR分解，分解的方法有多种，下面以Schmit正交化方法为例说明。

设A为n阶非奇异实矩阵，记为

72第72页，课件共101页，创作于2023年2月QR方法（续3）一般地，取：

（4-23）

于是：这就是用Schmit正交化方法对矩阵进行QR分解的过程。73第73页，课件共101页，创作于2023年2月

例7用Schmit正交化方法将矩阵进行QR分解。

[解]

74第74页，课件共101页，创作于2023年2月例7（续1）所以

75第75页，课件共101页，创作于2023年2月例7（续2）

基本QR方法每次迭代都需作一次QR分解与矩阵乘法，计算量大，而且收敛速度慢。因此实际使用的QR方法是先用一系列相似变换将A化成拟上三角矩阵（称为上Hessenberg矩阵），然后对此矩阵用基本QR方法。因为拟上三角矩阵中具有较多的零元素，这样做可减少运算量。化A为相似的拟上三角阵的方法有多种，本节介绍其中的一种：Householder变换

76第76页，课件共101页，创作于2023年2月3.2豪斯豪尔德（Householder）变换

为Householder矩阵或反射矩阵。容易证明，Householder矩阵具有以下性质：

（1）H是实对称的正交矩阵，即；（2）；（3）H仅有两个不等的特征值±1，其中1是n－1重特征值，－1是单重特征值，为其相应的特征向量。77第77页，课件共101页，创作于2023年2月定理9-2从图4-1可以看出，x+aw与H(x+aw)关于超平面(span{w})

对称，反射变换的名称由此而来。定理9-2设x,y为Rn中任意非零向量，且，则存在Householder矩阵H，使得： （4-25）78第78页，课件共101页，创作于2023年2月定理9-2（续1）

[证明]

由性质4容易得出，对任意xRn,w与xHx平行。故要使式（9-25）成立，应取

（4-26）令H=I2wwT，于是：

79第79页，课件共101页，创作于2023年2月定理9-2（续2）将上式代入式（4-27）即得：。定理9.2表明，对任一非零向量x，都可以构造一个Householder变换，它将x变成事先给定的单位向量的倍数。特别地，若取y=ei，则x经过Householder变换后可变成只有一个分量不为零。由2－范数的定义80第80页，课件共101页，创作于2023年2月定理9-2（续3）

使得Householder矩阵H=I2wwT将x变成与ei共线的向量，即

实际计算时，若，从而在计算w时会产生较大的误差，为此我们总取：81第81页，课件共101页，创作于2023年2月9.3化一般矩阵为拟上三角阵

的矩阵为拟上三角阵，也称为上海森堡（Hessenberg）阵。如果次对角线元hii1

(i=2,3,…,n)全不为零，则称该矩阵为不可约的上Hessenberg矩阵。下面讨论如何利用Householder变换将一般矩阵A相似变换成（4-29）的形式。首先，选取Householder矩阵H1，使得经H1相似变换后的矩阵H1AH1的第1列中有尽可能多的零元素。为此，应取H1为如下形式：

称形如82第82页，课件共101页，创作于2023年2月化一般矩阵为拟上三角阵（续1）其中为n-1阶

Householder矩阵。于是有：83第83页，课件共101页，创作于2023年2月化一般矩阵为拟上三角阵（续2）由定理4.2，只要取使得：就会使得变换后的矩阵H1AH1的第1列出现n-2个零元。完全类似地，可构造如下形式的Householder矩阵：使得：84第84页，课件共101页，创作于2023年2月化一般矩阵为拟上三角阵（续3）如此进行n-2次，可以构造n-2个Householder矩阵H1,H2,…,Hn-2，使得其中H为形如（4-29）的上Hessenberg矩阵。特别地，如果A为实对称矩阵，则经过上述正交相似变换后所得的矩阵H为三对角阵。

例8用Householder变换将矩阵化成拟上三角阵。

[解]因为a1=(1,0,0)T，由式（4-30）得H1=I。记：85第85页，课件共101页，创作于2023年2月例8

（续1）为使Householder矩阵满足由式（4-28），取：

86第86页，课件共101页，创作于2023年2月例8

（续2）于是有：

87第87页，课件共101页，创作于2023年2月3.4拟上三角矩阵的QR分解

因为