《数值分析实验报告7.doc》由会员分享，可在线阅读，更多相关《数值分析实验报告7.doc（14页珍藏版）》请在金锄头文库上搜索。
1、数值剖析实验报告7.实验七、QR算法一、实验目的1、熟习matlab编程并学习QR算法原理及计算机实现；2、学惯用matlab内置函数eig和QR算法求矩阵的特色值，并比对两者差异。二、实验题目1、课本第277页第1题107872345611144567267565111已知矩阵A,B03678,H6237.861090028975910111000106711（ 1）用MATLAB函数“eig”求矩阵所有特色值；（ 2）用基本QR算法求所有特色值（可用MATLAB函数“qr”实现矩阵的QR分解）。2、用QR算法求矩阵特色值：2345662144567(i)231(ii)03678111002
2、8900010依据QR算法原理编制求（i）及（ii）中矩阵所有特色值的程序并输出计算结果（要求偏差10-5）.三、实验原理与理论基础QR方法是一种变换方法，是计算一般矩阵（中小型矩阵）所有特色值问题的最有效方法之一。当前QR方法主要用来计算上海森伯格矩阵和对称三对角矩阵的所有特色值问题，/且QR方法拥有收敛快、算法稳固等特色。关于一般矩阵nnA （或对称矩阵），第一用豪斯霍尔德方法将A化为上海森伯格矩阵B（或对称三对角矩阵），而后再用QR方法计算B 的所有特色值。1、矩阵的QR分解设Ann非奇怪，则存在正交矩阵P，使PA=R，此中R为上三角矩阵。用Householder变换结构正交矩阵P，记A
3、(0)A，它的第一列记为a(0)，不如设a(0)0，可按公式（3.2）11（Th14,约化定理设x(x1,x2,xn)T0,则存在初等反射矩阵H使Hxe1，此中）HI1uuT,sgn(x1)x2,找到矩阵H1nn,H1I11u1u1T,使uxe1,212u2(x1).H1a1(0)1e1,e1(1,0,0)Tn.于是A(1)H1A(0)(H1a1(0),H1a2(0),H1an(1)1b(1),0(1)A此中A(1)(a1(1),a2(1),an(1)1)(n1)(n1).一般地，设A(j1)D(j1)B(j1)，此中D(j1)为(j-1)阶方阵，其对角线以下元素均0()A(j1)a1(j1)
4、，可选择(n-j+1)的Householder矩阵变换为0，A为(n-j+1)阶方阵，设其第一列为Hj(nj)(nj)，使Hja1(j1)je1,e1(1,0,0)nj1,依据Hj结构n*n阶的变换矩阵Hj为HjIj10HjA(j1)D(j)B(j).它和A(j1)0,于是有A(j)0j有近似HjA的形式，不过D(j)为j阶方阵，其对角线以下元素是0，这样经过n-1步运算获得Hn1H1AA(n1)R,此中RA(n1)为上三角矩阵，PHn1H1为正交矩阵，进而有PA=R。2、QR算法设Ann，且对A进行QR分解，即AQR，此中R为上三角矩阵，Q为正交矩阵，于是可获得一个新矩阵BRQQTAQ。明显
5、，B是由A经过正交相像变换获得，所以B与A特色值同样，再对B进行QR分解，又可得一新的矩阵，重复这一过程可获得矩阵序列：设AA1将A进行QR分解AQR1111作矩阵TRQQ1AQA21111求得Ak后将Ak进行QR分解AkQkRk形成矩阵Ak1RkQkQkTAkQkQR算法，就是利用矩阵的QR分解，按上述递推法例结构矩阵序列Ak的过程。只要A为非奇怪矩阵，则由QR算法就完整确立Ak。四、实验内容1、用matlab内置函数eig求矩阵的所有特色值；2、编写求特色值的QR算法程序，并用之求矩阵特色值；3、比较两种方法的结果差异。（1）QR算法的m文件functionqrsf(A,r)Q,R=qr(
6、A);t=A(1,1)%tempA=R*Q;fork=1:50Q,R=qr(A);t=A(1,1);A=R*Q;if(abs(A(1,1)-t)r)break;endendn=size(A,1);fori=1:nformatlonggdisp(特色值,num2str(i),end%disp();fori=1:ndisp(特色值);formatlongg,A(i,i)endformatlongg,A,Q,R=,num2str(A(i,i);（2）改良后的QR算法的m文件functionqrsf(A,r)Q,R=qr(A);%t=A(1,1)%tempt(1)=max(abs(diag(R);A=
7、R*Q;fork=2:50Q,R=qr(A);z=diag(A);t(k)=max(abs(diag(R);A=R*Q;if(abs(t(k)-t(k-1) A=10787;7565;86109;75910; B=23456;44567;03678;00289;00010; H6=hilb(6);1、eig求矩阵特色值eig(A),eig(B),eig(H6)ans=3.8580574559449530.2886853458021ans=1.59565457314994-0.390788045416488-0.929096277752298ans=1.08279948406811e-0071.
8、25707571226224e-0050.242360870575212、QR算法求矩阵特色值eig(A),qrsf(A,10-8)ans=3.8580574559449530.2886853458021t=10特色值1=30.2887特色值2=3.8581特色值3=0.84311特色值4=0.01015特色值ans=30.2886853457915特色值ans=3.85805737835431特色值ans=特色值ans=eig(B),qrsf(B,10-8)ans=1.59565457314994-0.390788045416488-0.929096277752298t=2特色值1=13.1
9、724特色值2=6.5519特色值3=1.5957特色值4=-0.9291特色值5=-0.39079特色值ans=特色值ans=特色值ans=1.59565457937031特色值ans=-0.929096283974607特色值ans=-0.390788045414554eig(H6),qrsf(H6,10-8)ans=1.08279948406811e-0071.25707571226224e-0050.24236087057521t=1特色值1=1.6189特色值2=0.24236特色值3=0.016322特色值4=0.00061575特色值5=1.2571e-005特色值6=1.0828e-007特色值ans=特色值ans=0.24236087069274特色值ans=特色值ans=特色值ans=1.25707571226506e-005特色值ans=1.08279948456401e-0073、改良后的QR算法求特色值eig(A),qrsf(A,10-8)ans=3.8580574559449530.2886853458021特色值1=30.2887特色值2=3.8581特色值3=0.84311特色值4=0.01015
