N阶实对称矩阵的三对角化
N阶实对称矩阵的三对角化
信计31 施昕 2130602018 实验原理:
用householder矩阵对n阶实矩阵上Hessenberg化,Q’AQ=T,其中Q是householder矩阵,T是上Hessenberg矩阵,因为A对称,则变换后的矩阵T是三对角矩阵。通过n个householder矩阵相乘,Hk使得每个Ak的第k行和第k列的后n-k分量变为ae(n-k)。最后得到三对角矩阵。
实验结果:
源程序:
%Householder变换:
function [v,b]=House(x) n=length(x);
x=x/norm(x,'inf'); c=x(2:n)'*x(2:n); v(2:n,1)=x(2:n); if c==0
b=0;
else
a=sqrt((x(1))^2+c); if x(1)
v(1)=x(1)-a;
else
v(1)=-c/(x(1)+a); end
b=2*(v(1))^2/(c+(v(1))^2); v=v/v(1);
end
%三对角化:
A=input('输入矩阵A:'); n=length(A);
Q=eye(n);
for k=1:n-2
[v,b]=House(A(k+1:n,k)); u=b*A(k+1:n,k+1:n)*v; w=u-(b*u'*v/2)*v; Hk=eye(n-k)-b*(v*v'); H=blkdiag(eye(k),Hk); Q=Q*H;
A(k+1,k)=norm(A(k+1:n,k)); A(k,k+1)=A(k+1,k);
A(k,k+2:n)=0;
A(k+2:n,k)=0;
A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-v*w'-w*v'; end
T=A
Q
disp(Q*T*Q');