Error in svd in r2018a

4 views (last 30 days)
S
S on 17 Jul 2022
I only have access to r2018a, where there is an error in svd. I tried the recommended fix
[U,s,V] = svd(A + 1e-14*randn(size(A)));
but I still getting the same error every now and then. How do you recommend getting around this error without increasing the size of the perturbation? Could I increase the number of iterations to test convergence, if so how? Or something else?
  1 Comment
Walter Roberson
Walter Roberson on 17 Jul 2022
The problem would continue if eps(A) > 1e-14 which would be the case if abs(A) > 100 or so

Sign in to comment.

Accepted Answer

Bruno Luong
Bruno Luong on 17 Jul 2022
Edited: Bruno Luong on 17 Jul 2022
As Walter has pointed out, the workwaround given by TMW might be still flawed because it ignores the scale of A.
A better solution would be
A = randi(100,4,6)
A = 4×6
82 8 13 18 98 82 77 8 10 61 29 85 28 48 21 23 86 96 70 94 98 63 60 33
[U,s,V] = svd(A + norm(A)*1e-14*randn(size(A)))
U = 4×4
-0.5147 0.4562 -0.1289 -0.7144 -0.4295 0.2837 0.7830 0.3493 -0.4856 0.2094 -0.6069 0.5932 -0.5610 -0.8170 0.0443 -0.1255
s = 4×6
275.2436 0 0 0 0 0 0 110.2492 0 0 0 0 0 0 66.4140 0 0 0 0 0 0 42.8050 0 0
V = 6×6
-0.4656 0.0719 0.5394 -0.5573 -0.0481 -0.4174 -0.3037 -0.5517 -0.2972 0.3213 0.2364 -0.5966 -0.2767 -0.6068 -0.0339 -0.1317 -0.6095 0.4065 -0.2978 -0.1917 0.5160 0.3314 0.5339 0.4619 -0.5026 0.1989 -0.5943 -0.3829 0.3450 0.2984 -0.5226 0.4959 -0.0125 0.5589 -0.4078 -0.0469
Or might be this cumbersome code that uses EIG (hopefully not buggy in your version) instead of SVD.
[m,n] = size(A);
AAc = A*A'; [U,s2u] = eig(1/2*(AAc+AAc'));
AcA = A'*A; [V,s2v] = eig(1/2*(AcA+AcA'));
[~, isu] = sort(sqrt(max(diag(s2u),0)),'descend');
[s, isv] = sort(sqrt(max(diag(s2v),0)),'descend');
U = U(:,isu);
S = diag(s);
if n < m
S(m,1) = 0;
isv = isv(1:n);
else
S = S(1:m,:);
end
V = V(:,isv);
U
U = 4×4
0.5147 0.4562 0.1289 -0.7144 0.4295 0.2837 -0.7830 0.3493 0.4856 0.2094 0.6069 0.5932 0.5610 -0.8170 -0.0443 -0.1255
S
S = 4×6
275.2436 0 0 0 0 0 0 110.2492 0 0 0 0 0 0 66.4140 0 0 0 0 0 0 42.8050 0 0
V % Note the two last columns are arbitrary provided they span the same subspace and are orthogonal
V = 6×6
0.4656 0.0719 -0.5394 -0.5573 -0.4065 -0.1064 0.3037 -0.5517 0.2972 0.3213 -0.6239 0.1500 0.2767 -0.6068 0.0339 -0.1317 0.4883 -0.5462 0.2978 -0.1917 -0.5160 0.3314 0.3822 0.5936 0.5026 0.1989 0.5943 -0.3829 0.2469 0.3836 0.5226 0.4959 0.0125 0.5589 0.0110 -0.4103
  1 Comment
Christine Tobler
Christine Tobler on 5 Aug 2022
Thank you for pointing this out, I have included scaling by norm(A) into the proposed workaround.

Sign in to comment.

More Answers (1)

Steven Lord
Steven Lord on 17 Jul 2022
Are you using Update 3 of release R2018a or a later Update, or are you using release R2018a (no Update) or Updates 1 or 2?
If you have not installed Update 3 or a later Update, please do.
  1 Comment
S
S on 17 Jul 2022
Edited: S on 17 Jul 2022
Thanks @Steven Lord. It's the version R2018a (9.4.0.813654) 64-bit (glnxa64) from Feb. 23, 2018. It's on a cluster where I don't have admin access (I'm not able to install a new version).

Sign in to comment.

Categories

Find more on Loops and Conditional Statements in Help Center and File Exchange

Tags

Products


Release

R2018a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!