jacobitezhengzhieight4.m

来自「数值分析方法,主要里面涵盖了计算方法方面matlab程序！」· M 代码 · 共 64 行

clear;
clc;
eps=1e-5;
N=0;
A=[4,-2,7,3,-1,8;
   -2,5,1,1,4,7;
   7,1,7,2,3,5;
   3,1,2,6,5,1;
   -1,4,3,5,3,2;
   8,7,5,1,2,4];
[m,n]=size(A);
while (1)
   N=N+1;
   B=abs(A);
   max=B(1,2);
   t=1;k=2;
   for i=1:m-1
       for j=i+1:n
           if (B(i,j)>max)
               max=B(i,j);
               t=i;
               k=j;
           end
       end
   end
   if max<eps
       break;
   end
   max=A(t,k);
   fprintf('max=%f,i=%d,j=%d\n',max,t,k);
   if A(t,t)~=A(k,k)
       jz=1/2*(atan(2*A(t,k)/(A(t,t)-A(k,k))));
   elseif A(t,t)==A(k,k)&A(t,k)>0
       jz=pi/4;
   elseif A(t,t)==A(k,k)&A(t,k)<0
       jz=-pi/4;
   end
   fprintf('第%d次旋转角度为：%f\n',N,jz*180/pi);
   A1=A;
   for l=1:n
       A1(t,l)=A(t,l)*cos(jz)+A(k,l)*sin(jz);
   end
   for l=1:n
       A1(l,t)=A(t,l)*cos(jz)+A(k,l)*sin(jz);
   end
   for l=1:n
       A1(k,l)=-A(t,l)*sin(jz)+A(k,l)*cos(jz);
   end
   for l=1:n
       A1(l,k)=-A(t,l)*sin(jz)+A(k,l)*cos(jz);
   end
   A1(t,t)=A(t,t)*(cos(jz)^2)+A(k,k)*(sin(jz)^2)+A(t,k)*sin(2*jz);
   A1(k,k)=A(t,t)*(sin(jz)^2)+A(k,k)*(cos(jz)^2)-A(t,k)*sin(2*jz);
   A1(t,k)=0;
   A1(k,t)=0;
   A=A1;
end
for i=1:n
    Y(i)=A(i,i);
end 
%输出特征值
Y'