《数值分析》实验报告

1、利用幂法求解矩阵的所有特征值、问题叙述在矩阵A中，特征值入的特征向量为v现如今要求任意n阶方阵A的特征值。如求任意4阶方阵的特征值。-5271225例如A=2085631428541021075045二、问题分析在矩阵A中，特征值入的特征向量为v，v中绝对值最大的元素为max(v), 则向量u的规格化向量为u=v/max(v)。设ARnxn，其最大特征值为入屮 对 应的特征向量为X,取任意非零初始向量v0=u0,它不与俎的特征向量正交，则 有lim u = lk max( x )k t ailim max( v )二九kik T a由此得出算法：输入:v,A；输出:入x1eu:=v;While

2、不满足判停准则dov:=Au;A1：=max(v);u:=v/A1；Endxi.:=u为了加快收敛，可以在算法中加入瑞列商(Rayleighquotient )。对于规格化向量u的瑞列商R(u)满足.Au , uR u )= X + O TOC o 1-5 h z :.u, u：1其中入2为第二大特征值。求解出一个特征值A1及其特征向量x1 后,可利用Householder变换矩阵H求出余下所有特征值及其特征向量。H =I一 2ww T， w = v/v ， v = x +be 2i i其中a为特征向量X的模。因此有haht =Xi T1，再对A1加以同样的循环 0 A1-i算法就可以得到其它

3、特征值。三、实验程序程序一：没有加入瑞列商的幕法程序A二52,7,12,25;2,0,85,63;14,28,54,102;10,7,50,45;u(:,1) = 1;2;1;3;fori = 1:15 v(:,i)=A*u(:,i);a(:,i)=max(v(:,i);u(:,i+1)=v(:,i)/a(:,i);end程序二：加入了瑞列商的幕法程序A二52,7,12,25;2,0,85,63;14,28,54,102;10,7,50,45;u(:,1) = 1;2;1;3;fori = 1:9 v(:,i)=A*u(:,i);a(i)二max(v(:,i);b1(i) = u(:,i)*A

4、*u(:,i);b2(i) = u(:,i)*u(:,i);b(i) = b1(i)/b2(i);u(:,i+1)=v(:,i)/a(:,i);end程序三：由得出的最大特征值对矩阵经行收缩处理e=1;0;0;0;I=eye(4);t1=sqrt(u(1/10)A2 + u(2/10)A2 + u(3/10)A2+u(4/10)A2); t=u(:,10)+t1*e;w1二sqrt(t人2+t (2)A2+t (3)人2+上(4)人2); w二t/w1;H=I-2*w*w;B=H*A*H;formatshort程序四：求出所有特征值A二52,7,12,25;2,0,85,63;14,28,54

5、,102;10,7,50,45;u = 1;2;3;4;fork=1:3u = u(1:5-k);i=1;c(1)=0;c(2)=0;b(1) = 1;j=1;whileabs(b(i)-c(i)0.0001v=A*u;a=max(v);i=i+1;b1=u*A*u;b2=u*u;b(i)二b1/b2;u二v/a;c(i+1)=b(i);endI二eye(5-k);disp(u);ifk= = 1t1=sqrt(u (1) A2 + u(2)A2 + u(3)A2 + u 人2);u(1) = u(1)+t1;w1=sqrt(u (1) A2 + u(2)A2 + u(3)A2 + u(4)A

6、2);elseifk=2t1=sqrt(u (1) A2 + u(2)A2 + u(3)A2);u(1) = u(1)+t1;w1=sqrt(u (1) A2 + u(2)A2 + u(3)A2);elseifk=3t1=sqrt(u (1)A2 + u(2)A2);u(1) = u(1)+t1;w1=sqrt(u (1)A2 + u(2)A2);endw=u/w1;H=I-2*w*w;B二H*A*H;disp(H);disp(B);A二B(2:5-k,2:5-k);disp(A);formatshort;End四、实验数据结果及分析(1)运行程序一，可以得到幕法求解的迭代过程取初始向量u=1

7、,2,1,3卩，按程序一经行迭代，计算结果列于下表表1-1幕法的迭代计算过程k012I467fk012I467f9153.0000 276.0000 430.0000 209.000010111213141547.1465 116.翦 26 126.5 孔 2 79.92 羽 56029 125.5 站 4 149必85 紀.5 呂阴 51.箔57 120762 1406 %.155孑 51.7545 126645 144.5566 %.恥 51.602百 12509弓 144.154 呂 .56 丸 51.6241 12M5469 144.2517 6.599 亍12 乱 5?园 144.2

8、255 6.5阳$ 51.6141 12M5392 144.2316 恥.591 亍 51.61 站 12主5 昭 6 144.2299 恥.5906 51.61 站 12主5孑8了 144.2 和弓 恥.5907 51.61 站 12主5孑8了 144.2 和 2 恥.5907 51.61 站 12主5孑8了 144.2 和 2 恥.5907123.5387 144.2302 86.5907 51.6133 123.5387 144.2302 86.59074 兀.0000126.5 兀 2149.丸 514左0丸6144.5566144.1548144.2517144.2255144.2

9、站 6144.2299144.2孔弓144.2 孔 2144.2302144.2 孔 2144.2 兀 21.00002.00001.0000主00000.站 580.64191.00000.职600.刀 260.91941.00000.6貂 70.站昶 0 屈4041.00000.5阳 00.站 910 屈6041.00000.60 右0.站取 0屈5551.00000.59990.站取 0屈56 音 1.00000.6005035了 90.85651.00000.600T0.站了90.呂 5661.00000.60040.站了90.呂 5651.00000.60040.站了 90.呂565

10、1.00000.60040.站了 90.呂5651.00000.60040.站了 90.呂5651.00000.60040.站了 90.呂5651.00000.60040.35790.85651.00000.60040.35790.85651.00000.6004从结果可以看出，在每次迭代布中做的规格化操作避免了溢出问题。经过12步迭代就已 经得到了很精确的结果：主特征值为入：! = 144.2302。(2)运行程序二，可以得到加速幕法求解的迭代过程在程序二中，加入了瑞列商来逼近主特征值，通过下表的计算迭代过程可以 发现只需要9步迭代就可以得到同样精确的数值：=144.1302。表1-2瑞列商

11、加速的效果011J47百153.0000 276.0000 430.0000 209.000047.1465 116.3326 126.5302 79.923353.6029 125.5 站 4 149.昶 5 眈.5 航 3011J47百153.0000 276.0000 430.0000 209.000047.1465 116.3326 126.5302 79.923353.6029 125.5 站 4 149.昶 5 眈.5 航 351.3657 123.0762 143.0386 86.155351.了545 123.6645 144.55師 6.乃恥51.6028 123.5093

12、144.15 绸 6.563051.6241 12 王 54E9 144.2517 86.59 财51.6132 123.5368 144.2255 86.589151.6141 12王5392 144.2316 86.59131.0000 2.0000 1.0000 3.0000117.4 冏殆0.0000126.5302149.3885143.0386144.55旺144.154 呂144.2517144.2255144.23160.35580.64191.00000.48600.37260.91941.00000.63170.35880.84041.00000.59300.3591 0.

13、8604 1.0000 0.60230.35800.85551.00000.59990.3580 0.8568 1.0000 0.60050 必 579 0.85651.0000 0.60 仍0.3579 0.8566 1.0000 0.60040.3579_0.8565_1.0000_0.6004144.7019142.9 了144.5114144.1600144.2495144.2259144.2315144.2299144.2303(3)运行程序三，可以得到Householder变换矩阵H1根据Householder变换得到一个正交化的矩阵H1，即利用收缩技术，可以得到矩阵B利用收缩技术

14、，可以得到矩阵B1，B1的首个元素就是A的主特征值:-0.2401-0.5746-0.6708-0.57460.7338-0.3108-0.6708-0.31080.6371-0.4027-0.1866-0.2179i-0.4027-0.1866-0.21790.8692144.2303。B =10.0000B =10.00000.00000.00002.013929.94614.4275-2.013948.09173.0058 取A1=24.94614.599528.16894.427519.47890.1564144.230256.0819-15 .4254-73.120748.09173.00584.599528.168919 .47890.1564。因此通过运行程序四，便可以得到所有A的特征值。依次所得到的矩阵分别为-0.7003有A的特征值。依次所得到的矩阵分别为-0.7003- 0.6329- 0.3301-0.63290.7644-0.1229-0.3301-0.12290.93590.9987- 0.0508H =3- 0.0508- 0.998746.8944-18.6185-7.2575B =0.0000- 35.737016.261620.00001.6353 -4.3877-36.5634-1 4.6 2 6 3B =