数值分析实验报告.doc

实验七、QR算法一、实验目的1、熟悉matlab编程并学习QR算法原理及计算机实现；2、学习用matlab内置函数eig和QR算法求矩阵的特征值，并比对二者差异。二、实验题目1、课本第277页第1题已知矩阵（1）用MATLAB函数“eig”求矩阵全部特征值；（2）用基本QR算法求全部特征值（可用MATLAB函数“qr”实现矩阵的QR分解）。2、用QR算法求矩阵特征值： 根据QR算法原理编制求（i）及（ii）中矩阵全部特征值的程序并输出计算结果（要求误差10 -5）.三、实验原理与理论基础QR方法是一种变换方法，是计算一般矩阵（中小型矩阵）全部特征值问题的最有效方法之一。目前QR方法主要用来计算上海森伯格矩阵和对称三对角矩阵的全部特征值问题，且QR方法具有收敛快、算法稳定等特点。对于一般矩阵（或对称矩阵），首先用豪斯霍尔德方法将A化为上海森伯格矩阵B（或对称三对角矩阵），然后再用QR方法计算B的全部特征值。1、矩阵的QR分解设非奇异，则存在正交矩阵P，使PA=R，其中R为上三角矩阵。用Householder变换构造正交矩阵P，记，它的第一列记为，不妨设，可按公式（3.2）（Th14,约化定理 设则存在初等反射矩阵H使，其中） 找到矩阵使于是其中一般地，设，其中为(j-1)阶方阵，其对角线以下元素均为0，为(n-j+1)阶方阵，设其第一列为，可选择(n-j+1)的Householder矩阵变换，使根据构造n*n阶的变换矩阵为于是有它和有类似的形式，只是为j阶方阵，其对角线以下元素是0，这样经过n-1步运算得到其中为上三角矩阵，为正交矩阵，从而有PA=R。2、QR算法设，且对A进行QR分解，即，其中R为上三角矩阵，Q为正交矩阵，于是可得到一个新矩阵。显然，B是由A经过正交相似变换得到，因此B与A特征值相同，再对B进行QR分解，又可得一新的矩阵，重复这一过程可得到矩阵序列：设将进行QR分解作矩阵求得后将进行QR分解形成矩阵QR算法，就是利用矩阵的QR分解，按上述递推法则构造矩阵序列的过程。只要A为非奇异矩阵，则由QR算法就完全确定。四、实验内容1、用matlab内置函数eig求矩阵的全部特征值；2、编写求特征值的QR算法程序，并用之求矩阵特征值；3、比较两种方法的结果差异。（1）QR算法的m文件function qrsf(A,r)Q,R=qr(A);t=A(1,1) %tempA=R*Q;for k=1:50 Q,R=qr(A); t=A(1,1); A=R*Q; if( abs( A(1,1)-t )r ) break; endendn=size(A,1);for i=1:n format long g disp( 特征值,num2str(i),=,num2str( A(i,i) ) );end%disp();for i=1:n disp(特征值); format long g,A(i,i)endformat long g,A,Q,R（2）改进后的QR算法的m文件function qrsf(A,r)Q,R=qr(A);%t=A(1,1) %tempt(1)=max(abs(diag(R);A=R*Q;for k=2:50 Q,R=qr(A); z=diag(A); t(k)=max(abs(diag(R); A=R*Q; if( abs( t(k) - t(k-1) ) A=10 7 8 7;7 5 6 5;8 6 10 9;7 5 9 10; B=2 3 4 5 6;4 4 5 6 7;0 3 6 7 8;0 0 2 8 9;0 0 0 1 0; H6=hilb(6);1、eig求矩阵特征值 eig(A),eig(B),eig(H6)ans = 0.0101500483978924 0.843107149855032 3.85805745594495 30.2886853458021ans = 13.1723513981032 6.55187835191566 1.59565457314994 -0.390788045416488 -0.929096277752298ans = 1.08279948406811e-007 1.25707571226224e-005 0.000615748354182652 0.0163215213198758 0.24236087057521 1.618899858924342、QR算法求矩阵特征值 eig(A),qrsf(A,10-8)ans = 0.0101500483978924 0.843107149855032 3.85805745594495 30.2886853458021t = 10特征值1=30.2887特征值2=3.8581特征值3=0.84311特征值4=0.01015特征值ans = 30.2886853457915特征值ans = 3.85805737835431特征值ans = 0.843107227456257特征值ans = 0.0101500483978911 eig(B),qrsf(B,10-8)ans = 13.1723513981032 6.55187835191566 1.59565457314994 -0.390788045416488 -0.929096277752298t = 2特征值1=13.1724特征值2=6.5519特征值3=1.5957特征值4=-0.9291特征值5=-0.39079特征值ans = 13.1723513891479特征值ans = 6.55187836087093特征值ans = 1.59565457937031特征值ans = -0.929096283974607特征值ans = -0.390788045414554 eig(H6),qrsf(H6,10-8)ans = 1.08279948406811e-007 1.25707571226224e-005 0.000615748354182652 0.0163215213198758 0.24236087057521 1.61889985892434t = 1特征值1=1.6189特征值2=0.24236特征值3=0.016322特征值4=0.00061575特征值5=1.2571e-005特征值6=1.0828e-007特征值ans = 1.6188998588068特征值ans = 0.24236087069274特征值ans = 0.01632152131988特征值ans = 0.000615748354182639特征值ans = 1.25707571226506e-005特征值ans = 1.08279948456401e-0073、改进后的QR算法求特征值 eig(A),qrsf(A,10-8)ans = 0.0101500483978924 0.843107149855032 3.85805745594495 30.2886853458021特征值1=30.2887特征值2=3.8581特征值3=0.84311特征值4=0.01015特征值ans = 30.2886853458019特征值ans = 3.85805745223919特征值ans = 0.843107153560957特征值ans = 0.0101500483978911 eig(B),qrsf(B,10-8)ans = 13.1723513981032 6.55187835191566 1.59565457314994 -0.390788045416488 -0.929096277752298特征值1=13.1724特征值2=6.5519特征值3=1.5957特征值4=-0.9291特征值5=-0.39079特征值ans = 13.1723513936489特征值ans = 6.55187835636998特征值ans = 1.59565456952802特征值ans = -0.929096274131193特征值ans = -0.390788045415675 eig(H6),qrsf(H6,10-8)ans = 1.08279948406811e-007 1.25707571226224e-005 0.000615748354182652 0.0163215213198758 0.24236087057521 1.61889985892434特征值1=1.6189特征值2=0.24236特征值3=0.016322特征值4=0.00061575特征值5=1.2571e-005特征值6=1.0828e-007特征值ans = 1.61889985892171特征值ans = 0.242360870577844特征值ans = 0.0163215213198758特征值ans = 0.000615748354182638特征值ans = 1.25707571226506e-005特征值ans = 1.08279948456401e-007六、实验结果分析与小结 从实验结果可以看出，用MATLAB内置函数eig求矩阵特征值与用QR算法求矩阵特征值的结果基本一致，数据只有微小差别。且单就QR算法而言，精度不同，计算出来特征值也存在一定的差异。不过这些微小差别对于计算来说影响不是很明显，除了有特殊需要要求更精确的数值外，这些结果已经能够满足计算结果的要求。七、附录3、eig与QR算法求矩阵特征值的结果比较 eig(A),qrsf(A,10-8)ans = 0.0101500483978924 0.843107149855032 3.85805745594495 30.2886853458021t = 10特征值1=30.2887特征值2=3.8581特征值3=0.84311特征值4=0.01015特征值ans = 30.2886853457915特征值ans = 3.85805737835431特征值ans = 0.843107227456257特征值ans = 0.0101500483978911A = Columns 1 through 3 30.2886853457915 1.67263719509357e-005 2.4658715813052e-009 1.67263719456478e-005 3.85805737835431 0.000483698079436215 2.46587341708519e-009 0.000483698079436138 0.843107227456257 -5.08101348050736e-023 5.50511453327879e-018 -1.40098262801473e-013 Column 4 -1.07002288191827e-015 8.90438591778815e-016 -1.38231041042597e-013 0.0101500483978911Q = Columns 1 through 3 -0.999999999990602 4.33543882155126e-006 4.37465097621799e-010 -4.33543785708671e-006 -0.999999835419719 0.000573708757424641 -2.92474424482639e-009 -0.000573708757417352 -0.999999835429117 -5.00590074187535e-021 5.42373229907221e-016 -1.38027186974381e-011 Column 4 -1.31911901611722e-021 8.46111373361065e-015 -1.38027161147482e-011 1R = Columns 1 through 3 -30.2886853454344 -0.000148041080872945 -1.00648595392028e-007 0 -3.85805646596483 -0.00269709930446301 0 0 -0.84310736620718 0 0 0 Column 4 -1.07019946786611e-015 -3.69340289178865e-015 -1.17754026710534e-011 0.0101500483978911 eig(B),qrsf(B,10-8)ans = 13.1723513981032 6.55187835191566 1.59565457314994 -0.390788045416488 -0.929096277752298t = 2特征值1=13.1724特征值2=6.5519特征值3=1.5957特征值4=-0.9291特征值5=-0.39079特征值ans = 13.1723513891479特征值ans = 6.55187836087093特征值ans = 1.59565457937031特征值ans = -0.929096283974607特征值ans = -0.390788045414554A = Columns 1 through 3 13.1723513891479 -11.2224332964075 1.38327278837818 -5.28299810482253e-009 6.55187836087093 -1.45097262375031 0 -4.57824183023911e-018 1.59565457937031 0 0 4.55829734055894e-008 0 0 0 Columns 4 through 5 12.2876570415913 2.18439169839953 -5.46614907884632 -0.424112905513863 -0.344534349380365 2.00394783492873 -0.929096283974607 -0.294914932293193 3.5295461400519e-012 -0.390788045414554Q = Columns 1 through 3 -1 -8.06333361799512e-010 -2.31352650689045e-027 8.06333361799512e-010 -1 -2.86919358232593e-018 0 2.86919358232594e-018 -0.999999999999999 0 0 4.90616249270291e-008 0 0 0 Columns 4 through 5 -1.13505369739799e-034 -1.02516554521335e-045 -1.40767299379114e-025 -1.27139170196985e-036 -4.90616249270291e-008 -4.43118132496026e-019 -0.999999999999999 -9.03186825049292e-012 -9.03186825049293e-012 1R = Columns 1 through 3 -13.1723513800989 11.2224333070288 -1.3832733912306 0 -6.55187836087093 1.45097289192847 0 0 -1.5956545624669 0 0 0 0 0 0 Columns 4 through 5 -12.2876569737454 2.18439169828855 5.46614900766307 -0.424112905464493 0.344534427647672 2.00394783493184 0.929096283977272 -0.294914932284802 0 -0.390788045414554 eig(H6),qrsf(H6,10-8)ans = 1.08279948406811e-007 1.25707571226224e-005 0.000615748354182652 0.0163215213198758 0.24236087057521 1.61889985892434t = 1特征值1=1.6189特征值2=0.24236特征值3=0.016322特征值4=0.00061575特征值5=1.2571e-005特征值6=1.0828e-007特征值ans = 1.6188998588068特征值ans = 0.24236087069274特征值ans = 0.01632152131988特征值ans = 0.000615748354182639特征值ans = 1.25707571226506e-005特征值ans = 1.08279948456401e-007A = Columns 1 through 3 1.6188998588068 1.27197394498537e-005 1.19768529656603e-012 1.2719739449838e-005 0.24236087069274 3.08869246409526e-008 1.19788968013046e-012 3.08869246092437e-008 0.01632152131988 -3.6564854489204e-021 -1.05134059527256e-016 -8.3901223728756e-011 -2.51874460802941e-031 -7.62887310800373e-027 -7.24698049380747e-021 7.33769428846597e-044 2.28755693138363e-039 2.38462685027147e-033 Columns 4 through 6 -4.65299860612185e-017 -7.85630875258594e-017 -4.99030380731949e-017 -4.69661140762161e-017 5.64100067650242e-017 -2.24039311154526e-017 -8.39012342417141e-011 -1.90672523361959e-017 9.70102034419473e-018 0.000615748354182639 8.55921938191856e-014 -5.5822503075125e-018 8.560823435926e-014 1.25707571226506e-005 -1.79372765726032e-017 -3.50196992607807e-026 -8.60412731528914e-018 1.08279948456401e-007Q = Columns 1 through 3 -0.999999998622786 5.24826445620704e-005 -2.59251450345072e-011 -5.24826445620255e-005 -0.999999998620995 1.89240475959173e-006 -7.33932613664668e-011 -1.8924047583461e-006 -0.9999999999982 5.93827888305781e-018 1.70741925354901e-013 1.36258949226309e-007 2.00365386384803e-026 6.06874592641492e-022 5.76495148470373e-016 6.77659566066424e-037 2.11263208377376e-032 2.20227926247271e-026 Columns 4 through 6 -5.09824115576308e-019 3.4008272077694e-028 1.97571398705825e-039 8.71151586217316e-014 -1.0917621219594e-022 -1.05044049606577e-033 -1.36258949226388e-007 3.51443236235137e-016 5.88051420890885e-027 -0.999999999999991 6.81010964765256e-009 2.17725918000353e-019 -6.81010964765255e-009 -1 -7.9461871176948e-011 -3.23418137522306e-019 -7.9461871176948e-011 1R = Columns 1 through 3 -1.61889985590967 -9.76838853035816e-005 -1.44084921219927e-010 0 -0.242360871026031 -4.89531790510286e-007 0 0 -0.0163215213199091 0 0 0 0 0 0 0 0 0 Columns 4 through 6 5.84784525152272e-017 7.8563087887546e-017 -4.99030380669521e-017 4.56367477791708e-014 -5.64098415545824e-017 -2.2403931119935e-017 2.30785457906863e-009 2.90479067967698e-017 9.70102034609643e-018 -0.000615748354182645 -4.2789060011646e-012 -5.58245625301427e-018 0 -1.25707571226506e-005 -1.01683315964937e-015 0 0 1.08279948456401e-007 4、eig与改进后的QR算法求矩阵特征值的结果比较 eig(A),qrsf(A,10-8)ans = 0.0101500483978924 0.843107149855032 3.85805745594495 30.2886853458021特征值1=30.2887特征值2=3.8581特征值3=0.84311特征值4=0.01015特征值ans = 30.2886853458019特征值ans = 3.85805745223919特征值ans = 0.843107153560957特征值ans = 0.0101500483978911A = Columns 1 through 3 30.2886853458019 2.13054162245146e-006 6.86375089207106e-011 2.13054161716379e-006 3.85805745223919 0.00010570327782691 6.86393453636055e-011 0.000105703277826833 0.843107153560957 1.70269961038915e-026 -1.44832419609733e-020 1.68662327113459e-015 Column 4 1.07002234228498e-015 -8.85182644520922e-016 -1.80487525890045e-016 0.0101500483978911Q = Columns 1 through 3 -0.999999999999848 5.52231697493189e-007 1.21771530170591e-011 -5.52231691626367e-007 -0.999999992140593 0.000125373478981475 -8.14123620399204e-011 -0.000125373478981449 -0.999999992140745 1.6775285630588e-024 -1.42691358634137e-018 1.66168988069558e-013 Column 4 4.42050280607066e-025 -2.22600977082374e-017 1.66168986584696e-013 1R = Columns 1 through 3 -30.2886853457962 -1.88569135427523e-005 -2.80162434078845e-009 0 -3.85805740866594 -0.000589401361893868 0 0 -0.843107160187151 0 0 0 Column 4 1.07002240146044e-015 -8.73123152404271e-016 1.3991777486471e-013 0.0101500483978911 eig(B),qrsf(B,10-8)ans = 13.1723513981032 6.55187835191566 1.59565457314994 -0.390788045416488 -0.929096277752298特征值1=13.1724特征值2=6.5519特征值3=1.5957特征值4=-0.9291特征值5=-0.39079特征值ans = 13.1723513936489特征值ans = 6.55187835636998特征值ans = 1.59565456952802特征值ans = -0.929096274131193特征值ans = -0.390788045415675A = Columns 1 through 3 13.1723513936489 -11.2224332937522 1.38327313998092 -2.62774351496931e-009 6.55187835636998 -1.4509727793467 0 -1.11499208998337e-018 1.59565456952802 0 0 -2.6541440298262e-008 0 0 0 Columns 4 through 5 12.2876570042594 -2.18439169861631 -5.46614903246668 0.424112904658541 -0.344534421512392 -2.0039478265026 -0.929096274131193 0.294914989541902 -1.48456567492261e-012 -0.390788045415675Q = Columns 1 through 3 -1 -4.01067201196511e-010 -2.80252858398934e-028 4.01067201196511e-010 -1 -6.98767831333128e-019 0 6.98767831333129e-019 -0.999999999999999 0 0 -2.85669429931243e-008 0 0 0 Columns 4 through 5 8.00596742954247e-036 -3.04138895236799e-047 1.99616608031226e-026 -7.58324027318754e-038 2.85669429931243e-008 -1.08523030585424e-019 -1 3.79890248009891e-012 3.79890248009891e-012 1R = Columns 1 through 3 -13.1723513891479 11.2224332990352 -1.38327278896012 0 -6.55187835636998 1.45097262319553 0 0 -1.59565457937031 0 0 0 0 0 0 Columns 4 through 5 -12.2876570437836 -2.18439169856963 5.46614907391814 0.424112904637775 0.344534375921806 -2.00394782650391 0.929096274132314 0.294914989538372 0 -0.390788045415675 eig(H6),qrsf(H6,10-8)ans = 1.08279948406811e-007 1.25707571226224e-005 0.000615748354182652 0.0163215213198758 0.24236087057521 1.61889985892434特征值1=1.6189特征值2=0.24236特征值3=0.016322特征值4=0.00061575特征值5=1.2571e-005特征值6=1.0828e-007特征值ans = 1.61889985892171特征值ans = 0.242360870577844特征值ans = 0.0163215213198758特征值ans = 0.000615748354182638特征值ans = 1.25707571226506e-005特征值ans = 1.08279948456401e-007A = Columns 1 through 3 1.61889985892171 1.90423583635592e-006 1.18725724292547e-014 1.90423583634025e-006 0.242360870577844 2.080045367086e-009 1.207695574630