Which algorithm does SVD function take?

Question

Bin 2024-1-29

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2075426-which-algorithm-does-svd-function-take

编辑： Matt J 2024-1-30

Hi Matlab supporter,

Could I you tell me which algorithm (algorithm name and basic description) does Matlab SVD function take for svd decomposion? (In Matlab 2020b). When I compared Matlab decomposition results with Scilab decomposition results, I find they are different, all meet SVD decomposion definition, but they have different U,V. ([U S V]=svd(A), U*S*V'=A)

Hope to get help from you.

Thank you,

Regards,

Bin

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Matt J 2024-1-29

2
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2075426-which-algorithm-does-svd-function-take#answer_1398751

编辑：Matt J 2024-1-29

The SVD algorithm is not disclosed by MathWorks, but even if you could be sure the same algorithms were used by both Matlab and SciLab, it would not help you avoid or understand differences in the results. Since you have a matrix with non-unique SVDs, the SVD computation is unstable, and can result in large differences from one implementation to another, merely from the fact that the code implementations are different.

13 个评论
显示 11更早的评论隐藏 11更早的评论

Bruno Luong 2024-1-29

编辑：Bruno Luong 2024-1-29

在 MATLAB Online 中打开

To prove the point I compute the singular vectors by 2 methods, svd and eig, of a random 1000 x 1000 matrix. They match at 11 digits, not 0.xxx. IMO if you find the difference of 1 digit between skilab and matlab the comparision is not correctly carried out. You need to sort the eigen vectors and flip the sign before comparison.

One can obserbe the error is large for smaller singular value, this is expected, but still accurate of many digits (no less than 10).

n = 1000;

A=randn(n);

[U,S,V]=svd(A);

s = diag(S); % by defaut s from SCD is descensingly sorted

B = A'*A;

B = 1/2*(B + B');

[W,D]=eig(B);

s1 = sqrt(diag(D));

[s1,is] = sort(s1,'descend');

W = W(:,is);

sf = sign(sum(W.*V,1));

W = W.*sf;

errV = vecnorm(V-W, Inf, 1);

norm(errV, inf)

ans = 1.0861e-12

B = A*A';

B = 1/2*(B + B');

[X,D]=eig(B);

s2 = sqrt(diag(D));

[s2,is] = sort(s2,'descend');

X = X(:,is);

sf = sign(sum(X.*U,1));

X = X.*sf;

errU = vecnorm(U-X, Inf, 1);

norm(errU, inf)

ans = 5.6631e-13

% Three different SVDs

norm(A-U*S*V')

ans = 4.8837e-13

norm(A-X*diag(s1)*W')

ans = 5.2702e-11

norm(A-X*diag(s2)*W')

ans = 5.2709e-11

close all

figure

subplot(2,1,1)

semilogy(errV)

xlabel('singular vector j');

ylabel('inf norm (V(:,j)-W(:,j))')

subplot(2,1,2)

semilogy(errU)

xlabel('singular vector j');

ylabel('inf norm (U(:,j)-X(:,j))')

Bin 2024-1-30

编辑：Walter Roberson 2024-1-30

在 MATLAB Online 中打开

temp1.mat

@Bruno Luong It is not allowed to upload attachment larger than 5MB, so I create 100x100 random matrix aa, calculate and save all the Matrix from Matlab and Scilab to temp1.mat file. u1, s1, v1 are Scilab SVD decomposition results, u, s, v are Matlab SVD decomposition results. You can download and run following codes in Matlab:

clear
load temp1.mat u s v u1 v1 s1 aa
disp('Scilab decomposition test');
[max(max(u1*s1*v1'-aa)) min(min(u1*s1*v1'-aa))]
disp('Matlab decomposition test');
[max(max(u*s*v'-aa)) min(min(u*s*v'-aa))]
disp('Max difference in u1 and u, Min difference in u1 and u');
[max(max(u1-u)) min(min(u1-u))]
disp('Max difference in v1 and v, Min difference in v1 and v');
[max(max(v1-v)) min(min(v1-v))]
disp('Test whether every frist element in every row in matrix u has the same sign');
[max(sign(u1(:,1))-sign(u(:,1))) min(sign(u1(:,1))-sign(u(:,1)))]
disp('Test whether every frist element in every row in matrix v has the same sign');
[max(sign(v1(:,1))-sign(v(:,1))) min(sign(v1(:,1))-sign(v(:,1)))]

The result is:

Scilab decomposition test

ans =

1.0e-14 *

0.7105 -0.4607

Matlab decomposition test

ans =

1.0e-14 *

0.5662 -0.5995

Max difference in u1 and u, Min difference in u1 and u

ans =

0.6751 -0.7299

Max difference in v1 and v, Min difference in v1 and v

ans =

0.6810 -0.7368

Test whether every frist element in every row in matrix u has the same sign

ans =

0 0

Test whether every frist element in every row in matrix v has the same sign

ans =

0 0

Bruno Luong 2024-1-30

编辑：Bruno Luong 2024-1-30

Matt's answer did not reveal the problem you had, Catalic's answer is right on.

Bruno Luong 2024-1-30

编辑：Bruno Luong 2024-1-30

@Bin " I would like to get 'same' U, V in Scilab like that in Matlab to make our calculation results identical."

You waste your time. The signs of eigen vectors are unpreditable. You could however force the sign to your own specification (*), but they come out of the algorithm almost arbitray.

(*) for example pick the largest element in absolute value of an singular vector, then set the sign so that this element is positive (assuming there no draw with opposite sign). Also you cannot change the signs of the odd number of eigen vectors, since U and V must be unitary matrices (det = 1) by convention of SVD.

请先登录，再进行评论。

Answer 2

Catalytic 2024-1-29

2
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2075426-which-algorithm-does-svd-function-take#answer_1398756

Comparing results showed that max difference in S and S1 reached 1e-14, but for U and U1 or for V and V1, max difference reached 0.xxx.

We would have to see S and S1, so you might wish to attach them for us in a .mat file. If all the S(i) are distinct, as I would expect for a random matrix A, the singular vectors are unique up to a difference in sign. So, when you calculate the max difference, you should make sure you accounted for a sign flip.

If some singular values S(i) are repeated, then the singular vectors are unstable and, as @Matt J mentions, you can expect to see all kinds of numerical differences for all kinds of reasons.