From 0f69e0de2858c694d2f3e80400ca3c9a9ea03511 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Del=20Pino?= <stephane.delpino44@gmail.com>
Date: Wed, 7 Jul 2021 07:45:03 +0200
Subject: [PATCH] Fix the case of incorrect 0 number of values prediction in
rows
---
src/algebra/CRSMatrixDescriptor.hpp | 12 ++++++++++--
tests/test_CG.cpp | 13 +++++++++++++
tests/test_CRSMatrixDescriptor.cpp | 8 ++++++++
3 files changed, 31 insertions(+), 2 deletions(-)
diff --git a/src/algebra/CRSMatrixDescriptor.hpp b/src/algebra/CRSMatrixDescriptor.hpp
index cbf82a904..8f123896d 100644
--- a/src/algebra/CRSMatrixDescriptor.hpp
+++ b/src/algebra/CRSMatrixDescriptor.hpp
@@ -35,6 +35,7 @@ class CRSMatrixDescriptor
Array<IndexType> row_map(m_nb_rows + 1);
row_map[0] = 0;
for (IndexType i = 0; i < m_nb_rows; ++i) {
+ Assert(non_zeros[i] >= 0);
row_map[i + 1] = row_map[i] + non_zeros[i];
}
@@ -83,7 +84,7 @@ class CRSMatrixDescriptor
const auto& overflow_row = A.m_overflow_per_row[i];
os << i << "|";
- if (A.m_nb_values_per_row[i] > 0) {
+ if (A.m_nb_values_per_row[i] + overflow_row.size() > 0) {
IndexType j = 0;
auto i_overflow = overflow_row.begin();
while (j < A.m_nb_values_per_row[i] or i_overflow != overflow_row.end()) {
@@ -207,7 +208,7 @@ class CRSMatrixDescriptor
for (IndexType i = 0; i < m_nb_rows; ++i) {
const auto& overflow_row = m_overflow_per_row[i];
- if (m_nb_values_per_row[i] > 0) {
+ if (m_nb_values_per_row[i] + overflow_row.size() > 0) {
IndexType j = 0;
auto i_overflow = overflow_row.begin();
while (j < m_nb_values_per_row[i] or i_overflow != overflow_row.end()) {
@@ -257,6 +258,13 @@ class CRSMatrixDescriptor
{
m_nb_values_per_row.fill(0);
m_values.fill(0);
+
+ // Diagonal is always set to fulfill PETSc's requirements
+ if (m_nb_columns == m_nb_rows) {
+ for (IndexType i = 0; i < m_nb_rows; ++i) {
+ this->operator()(i, i) = DataType{0};
+ }
+ }
}
~CRSMatrixDescriptor() = default;
diff --git a/tests/test_CG.cpp b/tests/test_CG.cpp
index b76d90a04..dbc1c3595 100644
--- a/tests/test_CG.cpp
+++ b/tests/test_CG.cpp
@@ -62,8 +62,21 @@ TEST_CASE("CG", "[algebra]")
non_zeros.fill(0);
CRSMatrixDescriptor<double> S{5, 5, non_zeros};
+ std::ostringstream Sout;
+ Sout << S;
+ REQUIRE(Sout.str() == R"(0| 0:0
+1| 1:0
+2| 2:0
+3| 3:0
+4| 4:0
+)");
+
CRSMatrix<double> A{S.getCRSMatrix()};
+ std::ostringstream Aout;
+ Aout << A;
+ REQUIRE(Sout.str() == Aout.str());
+
Vector<double> b{5};
b = 0;
diff --git a/tests/test_CRSMatrixDescriptor.cpp b/tests/test_CRSMatrixDescriptor.cpp
index 2bf59af22..2b311d84b 100644
--- a/tests/test_CRSMatrixDescriptor.cpp
+++ b/tests/test_CRSMatrixDescriptor.cpp
@@ -162,5 +162,13 @@ TEST_CASE("CRSMatrixDescriptor", "[algebra]")
REQUIRE_THROWS_AS(S(0, 5) = 1, AssertError);
}
+ SECTION("negative number of non-zeros")
+ {
+ Array<int> non_zeros(2);
+ non_zeros[0] = 1;
+ non_zeros[1] = -2;
+ REQUIRE_THROWS_AS((CRSMatrixDescriptor<int>{2, 4, non_zeros}), AssertError);
+ }
+
#endif // NDEBUG
}
--
GitLab