用Matlab实现：双边Jacobi的SVD

双边Jacobi SVD是一种计算矩阵奇异值分解（SVD）的方法，它通过不断地施加正交变换来把原矩阵变换成对角矩阵。以下是使用Matlab实现双边Jacobi SVD的示例代码：

function [U, S, V] = bidiagonal_jacobi_svd(A)
% 双边Jacobi SVD算法
% 输入：矩阵A
% 输出：左奇异向量U、奇异值S、右奇异向量V
[m, n] = size(A);
U = eye(m);
V = eye(n);
B = A;
for i = 1:min(m, n)
    % 对左下角的子矩阵进行双边Jacobi旋转
    [B, c, s] = bidiagonal_jacobi_rotate(B, i);
    % 更新左奇异向量U和右奇异向量V
    U(:, i:m) = U(:, i:m) * [c -s; s c];
    V(:, i:n) = V(:, i:n) * [c -s; s c];
end
% 提取对角线元素作为奇异值
S = diag(B);
end

function [B, c, s] = bidiagonal_jacobi_rotate(A, k)
% 对左下角的子矩阵进行双边Jacobi旋转
% 输入：矩阵A、旋转轴k
% 输出：旋转后的矩阵B、旋转矩阵的cos和sin值
[m, n] = size(A);
% 初始化旋转矩阵
c = 1;
s = 0;
for i = k+1:m
    % 计算旋转角度
    theta = atan2(A(i, k), A(k, k));
    c_old = c;
    s_old = s;
    c = cos(theta);
    s = sin(theta);
    % 应用旋转矩阵
    A(k, k) = A(k, k) * c + A(i, k) * s;
    A(i, k) = A(i, k) * c - A(k, k) * s;
    A(i, k+1:n) = c * A(i, k+1:n) - s * A(k, k+1:n);
    A(k, k+1:n) = c_old * A(k, k+1:n) + s_old * A(i, k+1:n);
    % 更新旋转矩阵
    [c, s] = givens_rotate(c_old, s_old, c, s);
end
B = A;
end

function [c, s] = givens_rotate(a, b, c, s)
% 计算Givens旋转矩阵的cos和sin值
% 输入：旋转前的cos和sin值a、b，旋转后的cos和sin值c、s
% 输出：旋转后的cos和sin值c、s
if b == 0
    c = sign(a);
    s = 0;
elseif abs(b) > abs(a)
    tau = a / b;
    s = sign(b) / sqrt(1 + tau^2);
    c = s * tau;
else
    tau = b / a;
    c = sign(a) / sqrt(1 + tau^2);
    s = c * tau;
end
end

使用该函数可以计算任意大小的矩阵的SVD，例如：

A = rand(5, 3);
[U, S, V] = bidiagonal_jacobi_svd(A);
disp(U);
disp(S);
disp(V);