不多废话,直接贴代码
function [A_, T_, T] = my_qr_givens( A ) %利用givens旋转进行qr分解 %输出 %A_ 每次变换后的A矩阵 %T_ 对应于A_的变换矩阵 A_ = sym([]); T_ = sym([]); A = sym(A); n = size(A,1); T = sym(eye(n)); sum = 1; for j = 1:n-1 %从第一列到第n-1列,全部变换到第j列第j个元素上 for i = j+1:n %从第j+1行到第n行 if A(i,j) ~= 0 r = (A(j,j)*A(j,j)'+ A(i,j)*A(i,j)')^(1/2); %r = simplify(sym(real(A(j,j))^2 + imag(A(j,j))^2 + real(A(i,j))^2 + imag(A(i,j))^2)^(1/2)); T_(:,:,sum) = eye(n); T_(j,j,sum) = A(j,j)/r; T_(i,i,sum) = A(j,j)/r; T_(j,i,sum) = A(i,j)/r; T_(i,j,sum) = -A(i,j)/r; A_(:,:,sum) = T_(:,:,sum)*A; T = simplify(T_(:,:,sum)*T); A = A_(:,:,sum); sum = sum + 1; end end end end