Trying to compute the kernel of a nxn zero matrix using FullPivLU throws a floating point exception. Minimal example: #include <iostream> #include <Eigen/Dense> int main(int argc, char **argv) { // Compute kernel: this line will create a floating point exception. Eigen::MatrixXd kernel = Eigen::Matrix2d::Zero().fullPivLu().kernel(); std::cout << "kernel" << kernel << std::endl; return 0; } Expected output: the kernel of a null matrix is the entire space, the return kernel should be a basis of the space (typically: the identity 2x2 matrix) Actual result: Floating point exception (core dumped) This behavior is present at least since Eigen 3.2.8, and still present in the developement build (9584fac85b3f). This exception comes from the triangular solving done inside LU/FullPivLu.h, l.685 m.topLeftCorner(rank(), rank()) .template triangularView<Upper>().solveInPlace( m.topRightCorner(rank(), dimker) ); For a zero matrix, rank() is zero, so we are solving for an empty matrix: this line is in fact has no effect, but TriangularSolver throws an exception. To fix this, we can either: - test before calling the TriangularSolver that rank is non-zero. - inside the TriangularSolver, return immediately if trying to solve an empty problem.
Good catch! I'd prefer having TriangularSolver (in fact any solver) return immediately (after asserting that the sizes are compatible) for 0x0-matrices.
Created attachment 889 [details] Minimal patch for the current bug. I agree fixing the solver is more generic. However from what I could see, there are separate implementation of triangular solvers in both Core/SolveTriangular and SparseCore/TriangularSolver. I've made a minimal patch for the current bug: it only modifies Core/SolveTriangular solveInPlace function. This is enough to get FullPivLU.kernel() to work (unit test added). If you see a better, more generic way to patch this, feel free to offer. This is my first time patching Eigen and I didn't go to far into working out the internals, but I can look deeper.
Thank you for the patch. Looks ok for dense triangular solve. https://bitbucket.org/eigen/eigen/commits/359f78f1be47/ https://bitbucket.org/eigen/eigen/commits/7d79b21a2ed8/ (3.3) Ultimately, we should likely apply similar fixes for sparse triangular solves and all solve functions.
https://bitbucket.org/eigen/eigen/commits/99e9c99cc552/ Summary: Bug 1617: add unit tests for empty triangular solve.
-- GitLab Migration Automatic Message -- This bug has been migrated to gitlab.com's GitLab instance and has been closed from further activity. You can subscribe and participate further through the new bug through this link to our GitLab instance: https://gitlab.com/libeigen/eigen/issues/1617.