吉文斯旋转 QR 分解
Givens rotation QR decomposition
我正在尝试创建一个函数来计算 Givens 旋转 QR 分解,遵循这个伪代码。
function [Q,R] = givens(A)
[m,n] = size(A);
indexI = zeros(m,n);
indexJ = zeros(m,n);
C = zeros(m,n);
S = zeros(m,n);
for i = 1:n
for j = i+1:m
c = A(i,i)/((A(i,i))^2 + (A(j,i)^2))^0.5;
s = A(j,i)/((A(i,i))^2 + (A(j,i)^2))^0.5;
A(i,:) = c*A(i,:) + s*A(j,:);
A(j,:) = -s*A(i,:) + c*A(j,:);
indexI(j,i) = i;
indexJ(j,i) = j;
C(j,i) = c;
S(j,i) = s;
end
end
R = A;
Q = eye(m);
for i = 1:n
for j= j+1:m
Q(:,i) = c*Q(:,i) + s*Q(:,j);
Q(:,j) = -s*Q(:,i) + c*Q(:,j);
end
end
但是,我得到的 R 矩阵不是上三角矩阵。我似乎找不到这里的错误。
幻灯片上有些歧义。
Givens旋转实际上是一次对两行进行矩阵乘法。
假设 [ri;rj]
是你的两行,Q
是相应的 givens 旋转矩阵。
更新是 [ri; rj] = Q*[ri; rj]
但在您的代码中,您先更新 ri
然后使用更新后的 ri
立即更新 rj
.
A(i,:) = c*A(i,:) + s*A(j,:);
A(j,:) = -s*A(i,:) + c*A(j,:);
这里是错误的修正:
tmp = c*A(i,:) + s*A(j,:);
A(j,:) = -s*A(i,:) + c*A(j,:);
A(i,:) = tmp;
我正在尝试创建一个函数来计算 Givens 旋转 QR 分解,遵循这个伪代码。
function [Q,R] = givens(A)
[m,n] = size(A);
indexI = zeros(m,n);
indexJ = zeros(m,n);
C = zeros(m,n);
S = zeros(m,n);
for i = 1:n
for j = i+1:m
c = A(i,i)/((A(i,i))^2 + (A(j,i)^2))^0.5;
s = A(j,i)/((A(i,i))^2 + (A(j,i)^2))^0.5;
A(i,:) = c*A(i,:) + s*A(j,:);
A(j,:) = -s*A(i,:) + c*A(j,:);
indexI(j,i) = i;
indexJ(j,i) = j;
C(j,i) = c;
S(j,i) = s;
end
end
R = A;
Q = eye(m);
for i = 1:n
for j= j+1:m
Q(:,i) = c*Q(:,i) + s*Q(:,j);
Q(:,j) = -s*Q(:,i) + c*Q(:,j);
end
end
但是,我得到的 R 矩阵不是上三角矩阵。我似乎找不到这里的错误。
幻灯片上有些歧义。
Givens旋转实际上是一次对两行进行矩阵乘法。
假设 [ri;rj]
是你的两行,Q
是相应的 givens 旋转矩阵。
更新是 [ri; rj] = Q*[ri; rj]
但在您的代码中,您先更新 ri
然后使用更新后的 ri
立即更新 rj
.
A(i,:) = c*A(i,:) + s*A(j,:);
A(j,:) = -s*A(i,:) + c*A(j,:);
这里是错误的修正:
tmp = c*A(i,:) + s*A(j,:);
A(j,:) = -s*A(i,:) + c*A(j,:);
A(i,:) = tmp;