SVD: Better default to gesvd! instead of gesdd!?

linearalgebra

#21

The difference is not related to full or not but instead if the singular vectors are requested in the bidiagonal solver. The documentation for dbdsdc explains that

The code currently calls DLASDQ if singular values only are desired.

I.e. when only the values are requested then the square root free QR algorithm is used instead of the Divide and Conquer algorithm since it is extremely fast but doesn’t compute singular vectors.

Hence there are actually three different algorithms: the square root free QR (svdvals), normal QR (LAPACK.bdsqr! or svd for small n) and DnC (svd). The two QR versions give very similar results, though, and they are very precise because they use the zero shift QR algorithm when there is a risk that a shift could cause significant cancellation. This is probably why you get better results from using the QR based solver.

As other people have also pointed out, I don’t think it’s a bug in the DnC solver. It’s just that solving secular equations have a larger error than zero shift QR.

In any case, I think having a keyword argument to select which solver to use is a very good idea.


#22

Alright. I made a PR that adds such a keyword argument.