Page 1 of 1

SVD error for M != N

PostPosted: Tue Aug 25, 2009 8:36 am
by jrk1982
The function 'sgesvd' for M>N and M<N, seems to return incorrect U and V. U and V are not orthogonal(especially U). The singular values seem to be computed correctly. I use Matlab's output as the reference.
==============
Input Matrix:
==============
A =

0.1209 0.4949
0.6432 0.6757
0.8746 0.4142
0.3369 0.8696
0.5964 0.1947

CULA output (u):
=================
u*u'

ans =

0.9543 0.1110 -0.3099 0.6114 0.2013
0.1110 0.6765 0.6049 0.6424 0.1549
-0.3099 0.6049 1.1470 -0.0406 0.6843
0.6114 0.6424 -0.0406 1.2464 -0.3316
0.2013 0.1549 0.6843 -0.3316 0.9758

u =

-0.2545 -0.4181 -0.0937 -0.4804 -0.6894
-0.5448 -0.0350 -0.4989 -0.2217 0.2837
-0.5322 0.5167 -0.6783 0.3002 0.2157
-0.4986 -0.5950 -0.2613 -0.7310 0.2029
-0.3266 0.4505 -0.4625 0.3093 -0.5972

v =

-0.7064 0.7078
0.5936 0.5925

v*v'

ans =

1.0000 -0.0000
-0.0000 0.7033

========================
u*s*v'

ans =

0.1209 -0.4150
0.6432 -0.5667
0.8746 -0.3474
0.3369 -0.7293
0.5964 -0.1633
====================
MATLAB outputs(u1, s1, v1)
u1*s1*v1'

ans =

0.1209 0.4949
0.6432 0.6757
0.8746 0.4142
0.3369 0.8696
0.5964 0.1947
=====================
u1*u1'

ans =

1.0000 0.0000 0.0000 0.0000 -0.0000
0.0000 1.0000 0.0000 0.0000 0
0.0000 0.0000 1.0000 -0.0000 0.0000
0.0000 0.0000 -0.0000 1.0000 -0.0000
-0.0000 0 0.0000 -0.0000 1.0000

v1*v1'

ans =

1 0
0 1
======================

Re:SVD error for M != N

PostPosted: Tue Aug 25, 2009 1:47 pm
by kyle
Hi Jrk1982,

Thanks for the bug report. Your issue is something we have experienced internally and we have been working to resolve the matter for the final release.

Specially, we've noticed that the reconstruction (A = U*SIGMA*V) is always correct but for M>N caess U isn't orthogonal (it's close) and that for N>M cases that V isn't orthogonal (also close). However square cases are always correct.

Some other quick points:

CULA & LAPACK both create VT (aka V') while MATLAB creates V. Be sure to consider this when checking for orthogonality.

Also, a 2x5 matrix is very small (and not a good candidate for GPU acceleration). If you look at some larger matrices it much easier to a see a pattern to the errors.

For example, this is what U*U' looks like for a larger matrix. You can see that identity is mostly visable:

Image

We'll try to keep you updated on the bug. Let us know if find any other beta quirks!

Re:SVD error for M != N

PostPosted: Tue Aug 25, 2009 2:14 pm
by jrk1982
Thanks for the response.
I'm not able to see a 'near orthogonal' U though. V is near orthogonal as you mention(for the case M>N). Things are the other way around when M<N. I'm surprised it happens only for the non-square case!

Re:SVD error for M != N

PostPosted: Tue Aug 25, 2009 2:24 pm
by jrk1982
Just curious: could you reproduce what I observed for the matrix I mentioned above?

Re:SVD error for M != N

PostPosted: Tue Aug 25, 2009 7:06 pm
by john
Thank you for the bug report - we have reproduced it and we're analyzing it. We will update when we have a solution.

Re:SVD error for M != N

PostPosted: Thu Aug 27, 2009 1:55 pm
by john
JRK, we have applied a fix in Beta 2. Can you try it against some data to determine if the problem has been alleviated?

Thank you,
John

Re:SVD error for M != N

PostPosted: Fri Aug 28, 2009 1:10 pm
by jrk1982
Indeed, the u and v are orthonormal now! Thanks a lot!

And in terms of Norm, match Matlab's orthogonality pretty well.

But the V returned is not orthonormal, to single precision. Here is the example output I ran...

A = grand(5,2)

A =

0.0620 0.1531
0.6746 0.2626
0.7100 0.6719
0.7672 0.9992
0.0207 0.5136

>> [u s v] = svd(A)

u =

-0.0876 0.1196 -0.5456 -0.7591 -0.3226
-0.3673 -0.6640 -0.3386 -0.0542 0.5538
-0.5515 -0.1448 0.7018 -0.4001 -0.1495
-0.7099 0.2260 -0.3081 0.4837 -0.3407
-0.2222 0.6876 0.0170 -0.1635 0.6714


s =

1.7679 0
0 0.4822
0 0
0 0
0 0


v =

-0.6754 -0.7375
-0.5047 0.4622

>> v*v'

ans =

1.0000 0.0000
0.0000 0.4683

>> u*u'

ans =

1.0000 0.0000 0.0000 0.0000 -0.0000
0.0000 1.0000 -0.0000 -0.0000 0.0000
0.0000 -0.0000 1.0000 0.0000 -0.0000
0.0000 -0.0000 0.0000 1.0000 -0.0000
-0.0000 0.0000 -0.0000 -0.0000 1.0000