From 57b1de2330d0660f353037a743fada210bde41ce Mon Sep 17 00:00:00 2001
From: Rasmus Munk Larsen <4643818-rmlarsen1@users.noreply.gitlab.com>
Date: Sun, 1 Mar 2026 23:47:35 -0800
Subject: [PATCH] Fix row-major GEMV dropping rows when n8 heuristic disables
 main loop

libeigen/eigen!2233

Co-authored-by: Rasmus Munk Larsen <rmlarsen@gmail.com>
---
 Eigen/src/Core/products/GeneralMatrixVector.h |  8 ++--
 test/product_large.cpp                        | 37 +++++++++++++++++++
 2 files changed, 40 insertions(+), 5 deletions(-)

diff --git a/Eigen/src/Core/products/GeneralMatrixVector.h b/Eigen/src/Core/products/GeneralMatrixVector.h
index 5f74e37de..e0205defb 100644
--- a/Eigen/src/Core/products/GeneralMatrixVector.h
+++ b/Eigen/src/Core/products/GeneralMatrixVector.h
@@ -423,7 +423,7 @@ general_matrix_vector_product<Index, LhsScalar, LhsMapper, RowMajor, ConjugateLh
     res[(i + 6) * resIncr] += alpha * cc6;
     res[(i + 7) * resIncr] += alpha * cc7;
   }
-  if (i < n4) {
+  for (; i < n4; i += 4) {
     ResPacket c0 = pzero(ResPacket{}), c1 = pzero(ResPacket{}), c2 = pzero(ResPacket{}), c3 = pzero(ResPacket{});
 
     for (Index j = 0; j < fullColBlockEnd; j += LhsPacketSize) {
@@ -451,9 +451,8 @@ general_matrix_vector_product<Index, LhsScalar, LhsMapper, RowMajor, ConjugateLh
     res[(i + 1) * resIncr] += alpha * cc1;
     res[(i + 2) * resIncr] += alpha * cc2;
     res[(i + 3) * resIncr] += alpha * cc3;
-    i += 4;
   }
-  if (i < n2) {
+  for (; i < n2; i += 2) {
     ResPacket c0 = pzero(ResPacket{}), c1 = pzero(ResPacket{});
 
     for (Index j = 0; j < fullColBlockEnd; j += LhsPacketSize) {
@@ -473,9 +472,8 @@ general_matrix_vector_product<Index, LhsScalar, LhsMapper, RowMajor, ConjugateLh
     }
     res[(i + 0) * resIncr] += alpha * cc0;
     res[(i + 1) * resIncr] += alpha * cc1;
-    i += 2;
   }
-  if (i < rows) {
+  for (; i < rows; ++i) {
     ResPacket c0 = pzero(ResPacket{});
     ResPacketHalf c0_h = pzero(ResPacketHalf{});
     ResPacketQuarter c0_q = pzero(ResPacketQuarter{});
diff --git a/test/product_large.cpp b/test/product_large.cpp
index b3213ecd8..4b8e52d5d 100644
--- a/test/product_large.cpp
+++ b/test/product_large.cpp
@@ -94,6 +94,42 @@ void product_large_regressions() {
   }
 }
 
+// Regression test: row-major GEMV with stride*sizeof > 32000 disables the
+// 8-row main loop (n8=0). The cleanup must use `for` loops (not `if`) to
+// process all remaining rows. Without the fix, only 7 out of `rows` results
+// are computed. This manifests as loss of orthogonality in QR of tall-skinny
+// matrices, since the Householder application uses row-major GEMV internally.
+template <int>
+void bug_gemv_rowmajor_large_stride() {
+  // Direct GEMV test: row-major A with stride (= cols) triggering n8=0.
+  // The threshold is stride * sizeof(Scalar) > 32000.
+  // For double: cols > 4000. For float: cols > 8000.
+  {
+    const int rows = 100;
+    const int cols = 5000;  // cols * sizeof(double) = 40000 > 32000
+    Matrix<double, Dynamic, Dynamic, RowMajor> A(rows, cols);
+    A.setRandom();
+    VectorXd x = VectorXd::Random(cols);
+    VectorXd y = A * x;
+    VectorXd y_ref = VectorXd::Zero(rows);
+    for (int i = 0; i < rows; ++i)
+      for (int j = 0; j < cols; ++j) y_ref(i) += A(i, j) * x(j);
+    VERIFY_IS_APPROX(y, y_ref);
+  }
+
+  // QR orthogonality test: this is the high-level symptom.
+  // HouseholderQR of a col-major (m x n) matrix with m > 4000
+  // uses row-major GEMV internally during Householder application.
+  {
+    const int m = 5000;
+    const int n = 50;
+    MatrixXd A = MatrixXd::Random(m, n);
+    MatrixXd Q = A.householderQr().householderQ() * MatrixXd::Identity(m, n);
+    MatrixXd QtQ = Q.adjoint() * Q;
+    VERIFY_IS_APPROX(QtQ, MatrixXd::Identity(n, n));
+  }
+}
+
 template <int>
 void bug_1622() {
   typedef Matrix<double, 2, -1, 0, 2, -1> Mat2X;
@@ -138,6 +174,7 @@ EIGEN_DECLARE_TEST(product_large) {
   }
 
   CALL_SUBTEST_6(product_large_regressions<0>());
+  CALL_SUBTEST_6(bug_gemv_rowmajor_large_stride<0>());
 
   // Regression test for bug 714:
 #if defined EIGEN_HAS_OPENMP