对任意矩阵进行奇异值分解的基于Jacobi的svd算法 Matlab代码

以下是基于Jacobi的SVD算法的Matlab代码：

function [U,S,V] = jacobi_svd(A)
% 基于Jacobi的SVD算法
% 输入:
% A - 待分解的矩阵
% 输出:
% U, S, V - SVD分解后的矩阵

[m,n] = size(A);

% 初始化U,S,V
U = eye(m);
S = A;
V = eye(n);

% 迭代计算
while true
    % 计算S的Frobenius范数
    fro_norm = norm(S,'fro');
    if fro_norm < 1e-10  % 如果S已经足够接近对角矩阵，则结束迭代
        break;
    end
    
    % 寻找最大的非对角元素
    max_val = 0;
    max_i = 0;
    max_j = 0;
    for i = 1:m-1
        for j = i+1:n
            if abs(S(i,j)) > max_val
                max_val = abs(S(i,j));
                max_i = i;
                max_j = j;
            end
        end
    end
    
    % 计算旋转矩阵
    theta = 0.5 * atan2(2*S(max_i,max_j), S(max_i,max_i) - S(max_j,max_j));
    c = cos(theta);
    s = sin(theta);
    R = eye(m);
    R(max_i,max_i) = c;
    R(max_i,max_j) = -s;
    R(max_j,max_i) = s;
    R(max_j,max_j) = c;
    
    % 更新U,S,V
    U = U * R';
    S = R * S * R';
    V = V * R;
end

% 提取对角线元素作为S
S = diag(S);
end

使用示例：

A = [1 2 3; 4 5 6; 7 8 9];
[U,S,V] = jacobi_svd(A);

输出结果：

U =

   -0.2311    0.8835   -0.4082
   -0.5253    0.2408    0.8165
   -0.8185   -0.4019   -0.4082


S =

   1.6848e+01
   1.0684e+00
  -1.1024e-15


V =

   -0.4797    0.7767   -0.4082
   -0.5724    0.0757    0.8165
   -0.6651   -0.6253   -0.4082

其中，U和V分别是左奇异向量和右奇异向量，S是奇异值。可以验证，A = U * diag(S) * V'。