JacobiSVD fails on small rank-deficient float matrices (U is not unitary)
Submitted by Charles Karney
Assigned to Nobody
Link to original bugzilla bug (#790)
Version: 3.2
Operating system: Linux
Description
(Found during testing of the new alignPoints function. This bug will
also cause umeyama to fail.)
JacobiSVD can fail spectacularly with small rank-deficient float
matrices. Here's an example (Eigen 3.2.1, Linux 64-bit, g++ (GCC)
4.7.2 20121109, optimize or debug):
Matrix<float,4,4> M;
M <<
0.030851028,0.275531113,-0.451450854,-0.126374155,
-0.015783134,-0.140959471,0.230958566,0.064651995,
0.051904849,0.463563203,-0.759536683,-0.212616309,
0.037682864,0.336546361,-0.551422834,-0.154359206;
JacobiSVD< Matrix<float,4,4> > svd(M, ComputeFullU | ComputeFullV);
std::cout << "U\n" << svd.matrixU() << "\n\n";
std::cout << "U.determinant() " << svd.matrixU().determinant() << "\n\n";
This prints out
U
-0.423173 -0.740408 -1.07744e-11 0
0.216492 0.378787 5.51207e-12 0
-0.711961 0.555262 8.08013e-12 0
-0.516883 0 0 0
U.determinant() 0
JacobiSVD does OK on M.transpose() and if float is promoted to double.
I realize that there are limits to what I can expect from floats. But
it would be nice if the basic properties of the SVD were more accurately
met (specifically U should be nearly unitary).