From becefd59e2186fa340c8014c2e060400bb7a6c0e Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Antonio=20S=C3=A1nchez?= <cantonios@google.com>
Date: Wed, 12 Feb 2025 07:09:20 +0000
Subject: [PATCH] Returns condition number of zero if matrix is not invertible.

---
 Eigen/src/LU/FullPivLU.h |  5 ++++-
 test/lu.cpp              | 29 +++++++++++++++++++++++++++--
 2 files changed, 31 insertions(+), 3 deletions(-)

diff --git a/Eigen/src/LU/FullPivLU.h b/Eigen/src/LU/FullPivLU.h
index 466834ada..3e57764ac 100644
--- a/Eigen/src/LU/FullPivLU.h
+++ b/Eigen/src/LU/FullPivLU.h
@@ -243,7 +243,10 @@ class FullPivLU : public SolverBase<FullPivLU<MatrixType_, PermutationIndex_> >
       the LU decomposition.
     */
   inline RealScalar rcond() const {
-    eigen_assert(m_isInitialized && "PartialPivLU is not initialized.");
+    eigen_assert(m_isInitialized && "FullPivLU is not initialized.");
+    if (!isInvertible()) {
+      return RealScalar(0);
+    }
     return internal::rcond_estimate_helper(m_l1_norm, *this);
   }
 
diff --git a/test/lu.cpp b/test/lu.cpp
index 1792c2bcb..b20bcfc80 100644
--- a/test/lu.cpp
+++ b/test/lu.cpp
@@ -195,6 +195,29 @@ void lu_verify_assert() {
   VERIFY_RAISES_ASSERT(plu.inverse())
 }
 
+// Rank-deficient matrix returns 0.
+// https://gitlab.com/libeigen/eigen/-/issues/2889
+void test_2889() {
+  Eigen::MatrixXd A =
+      Eigen::MatrixXd({{0.0000000000000000, 0.0000000000000000, 1.0000000000000000, 0.0000000000000000,
+                        0.34149999916553497, 0.0000000000000000, 0.79877008515664061},
+                       {0.0000000000000000, 1.0000000000000000, 0.0000000000000000, 0.29200000315904617,
+                        0.0000000000000000, -0.37149999849498272, -0.16425902650844920},
+                       {0.0000000000000000, 0.0000000000000000, 1.0000000000000000, 0.0000000000000000,
+                        0.34149999916553497, 0.0000000000000000, 0.79877008515664061},
+                       {0.0000000000000000, 1.0000000000000000, 0.0000000000000000, 0.040500000119209290,
+                        0.0000000000000000, -0.30099999904632568, -0.081170580429391403},
+                       {1.0000000000000000, 0.0000000000000000, 0.0000000000000000, 0.0000000000000000,
+                        0.0000000000000000, 0.0000000000000000, -0.0000000000000000},
+                       {0.0000000000000000, 0.70710689672598170, 0.70710666564709435, 0.027000000700354562,
+                        0.025455838867477515, -0.025455847186317101, -0.0068972271572272821},
+                       {1.0000000000000000, 0.0000000000000000, 0.0000000000000000, 0.0000000000000000,
+                        0.0000000000000000, 0.0000000000000000, -0.0000000000000000}});
+  Eigen::FullPivLU<Eigen::MatrixXd> lu_factorization(A);
+  double rcond = lu_factorization.rcond();
+  VERIFY_IS_EQUAL(rcond, 0.0);
+}
+
 EIGEN_DECLARE_TEST(lu) {
   for (int i = 0; i < g_repeat; i++) {
     CALL_SUBTEST_1(lu_non_invertible<Matrix3f>());
@@ -229,7 +252,9 @@ EIGEN_DECLARE_TEST(lu) {
     CALL_SUBTEST_7((lu_non_invertible<Matrix<float, Dynamic, 16> >()));
 
     // Test problem size constructors
-    CALL_SUBTEST_9(PartialPivLU<MatrixXf>(10));
-    CALL_SUBTEST_9(FullPivLU<MatrixXf>(10, 20););
+    CALL_SUBTEST_8(PartialPivLU<MatrixXf>(10));
+    CALL_SUBTEST_8(FullPivLU<MatrixXf>(10, 20););
+
+    CALL_SUBTEST_9(test_2889());
   }
 }