% Jacobi方法
% 注释的地方需要修改
function Jacobimethod;
A = [3,1,2;1,3,4;2,4,6]; %所求矩阵
[V,D]=jacobi(A);
output(A,V,D);
function [p,q,convergence]=pq(A)
U=triu(A)-diag(diag(A));
[M,I] = max(abs(U));
[N,J] = max(M);
p=I(J); q=J;
s = abs(A(p,q));
if (p==q) || (s <= eps),
convergence=true;
else
convergence=false;
end
function R= givens(p,q,A)
b = (A(p,p) - A(q,q))/(2 * A(p,q));
t = 1/(abs(b) + sqrt(1 + b^2));
if b < 0 t = -t; end
c = 1 /sqrt(1 + t^2);
s = t * c;
R = eye(size(A));
R(p,p) = c; R(p,q) = -s;
R(q,p) = s; R(q,q) = c;
function [V,D]=jacobi(A)
V=eye(size(A));
convergence=false;
while ~convergence
[p,q,convergence]= pq(A);
if convergence, break; end;
R= givens(p,q,A);
A = R' * A * R;
V = V * R;
end
D = diag(A);
[D,IX] = sort(D);
V=V(:,IX);
function output(A,V,D)
n=length(A);
for i = 1:n
fprintf('lambda%1d = %.8f\nv%1d = ',i,D(i),i); % 结果有效数字位数为.和f之间的数
for j = 1:n
fprintf('%10.8f ',V(j,i)); % 结果有效数字位数为.和f之间的数
end
fprintf('\n\n');
end