From 15bdf87defcc7e3ee9707e4badc749d372a519ac Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?St=C3=A9phane=20Del=20Pino?= <stephane.delpino44@gmail.com>
Date: Fri, 10 Mar 2023 19:12:49 +0100
Subject: [PATCH] Fix copy for SmallMatrix and add unit tests
---
src/algebra/SmallMatrix.hpp | 11 +-
tests/CMakeLists.txt | 1 +
tests/test_SmallMatrix.cpp | 416 ++++++++++++++++++++++++++++++++++++
3 files changed, 426 insertions(+), 2 deletions(-)
create mode 100644 tests/test_SmallMatrix.cpp
diff --git a/src/algebra/SmallMatrix.hpp b/src/algebra/SmallMatrix.hpp
index d0d6ee329..b8d3b81b5 100644
--- a/src/algebra/SmallMatrix.hpp
+++ b/src/algebra/SmallMatrix.hpp
@@ -38,7 +38,7 @@ class [[nodiscard]] SmallMatrix // LCOV_EXCL_LINE
friend PUGS_INLINE SmallMatrix<std::remove_const_t<DataType>> copy(const SmallMatrix& A) noexcept
{
- return {A.m_nb_rows, A.m_nb_columns, copy(A.m_values)};
+ return SmallMatrix<std::remove_const_t<DataType>>{A.m_nb_rows, A.m_nb_columns, copy(A.m_values)};
}
friend PUGS_INLINE SmallMatrix<std::remove_const_t<DataType>> transpose(const SmallMatrix& A)
@@ -130,7 +130,8 @@ class [[nodiscard]] SmallMatrix // LCOV_EXCL_LINE
{
static_assert(std::is_same_v<std::remove_const_t<DataType>, std::remove_const_t<DataType2>>,
"incompatible data types");
- Assert((m_nb_rows == B.numberOfRows()) and (m_nb_columns == B.numberOfColumns()), "incompatible matrix sizes");
+ Assert((m_nb_rows == B.numberOfRows()) and (m_nb_columns == B.numberOfColumns()),
+ "cannot add matrix: incompatible sizes");
parallel_for(
m_values.size(), PUGS_LAMBDA(index_type i) { m_values[i] += B.m_values[i]; });
@@ -260,6 +261,12 @@ class [[nodiscard]] SmallMatrix // LCOV_EXCL_LINE
SmallMatrix(SmallMatrix &&) = default;
+ explicit SmallMatrix(size_t nb_rows, size_t nb_columns, const SmallArray<DataType>& values)
+ : m_nb_rows{nb_rows}, m_nb_columns{nb_columns}, m_values{values}
+ {
+ Assert(m_values.size() == m_nb_columns * m_nb_rows, "incompatible array size and matrix dimensions")
+ }
+
explicit SmallMatrix(size_t nb_rows, size_t nb_columns) noexcept
: m_nb_rows{nb_rows}, m_nb_columns{nb_columns}, m_values{nb_rows * nb_columns}
{}
diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt
index 0e585acac..b5301f36a 100644
--- a/tests/CMakeLists.txt
+++ b/tests/CMakeLists.txt
@@ -131,6 +131,7 @@ add_executable (unit_tests
test_RefItemList.cpp
test_RevisionInfo.cpp
test_SmallArray.cpp
+ test_SmallMatrix.cpp
test_SmallVector.cpp
test_Socket.cpp
test_SocketModule.cpp
diff --git a/tests/test_SmallMatrix.cpp b/tests/test_SmallMatrix.cpp
new file mode 100644
index 000000000..9712e13e3
--- /dev/null
+++ b/tests/test_SmallMatrix.cpp
@@ -0,0 +1,416 @@
+#include <catch2/catch_approx.hpp>
+#include <catch2/catch_test_macros.hpp>
+#include <catch2/matchers/catch_matchers_all.hpp>
+
+#include <Kokkos_Core.hpp>
+
+#include <utils/PugsAssert.hpp>
+#include <utils/Types.hpp>
+
+#include <algebra/SmallMatrix.hpp>
+
+#include <sstream>
+
+// clazy:excludeall=non-pod-global-static
+
+TEST_CASE("SmallMatrix", "[algebra]")
+{
+ SmallMatrix<int> A{3, 4};
+
+ REQUIRE(A.numberOfRows() == 3);
+ REQUIRE(A.numberOfColumns() == 4);
+
+ for (size_t i = 0; i < 3; ++i) {
+ for (size_t j = 0; j < 4; ++j) {
+ A(i, j) = 1 + 4 * i + j;
+ }
+ }
+
+ REQUIRE(not A.isSquare());
+
+ REQUIRE(((A(0, 0) == 1) and (A(0, 1) == 2) and (A(0, 2) == 3) and (A(0, 3) == 4) and //
+ (A(1, 0) == 5) and (A(1, 1) == 6) and (A(1, 2) == 7) and (A(1, 3) == 8) and //
+ (A(2, 0) == 9) and (A(2, 1) == 10) and (A(2, 2) == 11) and (A(2, 3) == 12)));
+
+ SmallMatrix<int> B{3, 4};
+ B(0, 0) = -6;
+ B(0, 1) = 5;
+ B(0, 2) = 3;
+ B(0, 3) = 8;
+
+ B(1, 0) = -4;
+ B(1, 1) = 6;
+ B(1, 2) = 5;
+ B(1, 3) = 3;
+
+ B(2, 0) = 7;
+ B(2, 1) = 1;
+ B(2, 2) = 3;
+ B(2, 3) = 6;
+
+ SECTION("copy and shallow copy")
+ {
+ SmallMatrix<const int> affected_A;
+
+ REQUIRE(affected_A.numberOfRows() == 0);
+ REQUIRE(affected_A.numberOfColumns() == 0);
+
+ affected_A = A;
+ REQUIRE(&A(0, 0) == &affected_A(0, 0));
+
+ REQUIRE(affected_A.numberOfRows() == A.numberOfRows());
+ REQUIRE(affected_A.numberOfColumns() == A.numberOfColumns());
+
+ const SmallMatrix<int> copy_A = copy(A);
+
+ REQUIRE(copy_A.numberOfRows() == A.numberOfRows());
+ REQUIRE(copy_A.numberOfColumns() == A.numberOfColumns());
+
+ bool is_same = true;
+ for (size_t i = 0; i < A.numberOfRows(); ++i) {
+ for (size_t j = 0; j < A.numberOfColumns(); ++j) {
+ is_same &= (copy_A(i, j) == A(i, j));
+ }
+ }
+ REQUIRE(is_same);
+
+ REQUIRE(&A(0, 0) != ©_A(0, 0));
+ }
+
+ SECTION("fill and special affectations")
+ {
+ SmallMatrix<int> M(5);
+ REQUIRE(M.isSquare());
+ M.fill(1);
+ for (size_t i = 0; i < M.numberOfRows(); ++i) {
+ for (size_t j = 0; j < M.numberOfColumns(); ++j) {
+ REQUIRE(M(i, j) == 1);
+ }
+ }
+
+ M.fill(-4);
+ for (size_t i = 0; i < M.numberOfRows(); ++i) {
+ for (size_t j = 0; j < M.numberOfColumns(); ++j) {
+ REQUIRE(M(i, j) == -4);
+ }
+ }
+
+ M = zero;
+ for (size_t i = 0; i < M.numberOfRows(); ++i) {
+ for (size_t j = 0; j < M.numberOfColumns(); ++j) {
+ REQUIRE(M(i, j) == 0);
+ }
+ }
+
+ M = identity;
+ for (size_t i = 0; i < M.numberOfRows(); ++i) {
+ for (size_t j = 0; j < M.numberOfColumns(); ++j) {
+ REQUIRE(M(i, j) == (i == j));
+ }
+ }
+ }
+
+ SECTION("build from TinyMatrix")
+ {
+ TinyMatrix<3, 4> T{1, 2, 3, 4, 0.5, 1, 1.5, 2, 0.3, 0.6, 0.9, 1.2};
+ SmallMatrix S{T};
+
+ REQUIRE(S.numberOfRows() == T.numberOfRows());
+ REQUIRE(S.numberOfColumns() == T.numberOfColumns());
+
+ for (size_t i = 0; i < S.numberOfRows(); ++i) {
+ for (size_t j = 0; j < S.numberOfColumns(); ++j) {
+ REQUIRE(S(i, j) == T(i, j));
+ }
+ }
+ }
+
+ SmallMatrix AT = transpose(A);
+
+ SECTION("transposed")
+ {
+ REQUIRE(AT.numberOfColumns() == A.numberOfRows());
+ REQUIRE(AT.numberOfRows() == A.numberOfColumns());
+
+ bool is_same = true;
+ for (size_t i = 0; i < A.numberOfRows(); ++i) {
+ for (size_t j = 0; j < A.numberOfColumns(); ++j) {
+ is_same &= (AT(j, i) == A(i, j));
+ }
+ }
+ REQUIRE(is_same);
+ }
+
+ SECTION("sum")
+ {
+ SmallMatrix AaddB = copy(A);
+ AaddB += B;
+
+ REQUIRE(AaddB(0, 0) == -5);
+ REQUIRE(AaddB(0, 1) == 7);
+ REQUIRE(AaddB(0, 2) == 6);
+ REQUIRE(AaddB(0, 3) == 12);
+
+ REQUIRE(AaddB(1, 0) == 1);
+ REQUIRE(AaddB(1, 1) == 12);
+ REQUIRE(AaddB(1, 2) == 12);
+ REQUIRE(AaddB(1, 3) == 11);
+
+ REQUIRE(AaddB(2, 0) == 16);
+ REQUIRE(AaddB(2, 1) == 11);
+ REQUIRE(AaddB(2, 2) == 14);
+ REQUIRE(AaddB(2, 3) == 18);
+
+ SmallMatrix ApB = A + B;
+ for (size_t i = 0; i < ApB.numberOfRows(); ++i) {
+ for (size_t j = 0; j < ApB.numberOfColumns(); ++j) {
+ REQUIRE(AaddB(i, j) == ApB(i, j));
+ }
+ }
+ }
+
+ SECTION("difference")
+ {
+ SmallMatrix AsubsB = copy(A);
+ AsubsB -= B;
+
+ REQUIRE(AsubsB(0, 0) == 7);
+ REQUIRE(AsubsB(0, 1) == -3);
+ REQUIRE(AsubsB(0, 2) == 0);
+ REQUIRE(AsubsB(0, 3) == -4);
+
+ REQUIRE(AsubsB(1, 0) == 9);
+ REQUIRE(AsubsB(1, 1) == 0);
+ REQUIRE(AsubsB(1, 2) == 2);
+ REQUIRE(AsubsB(1, 3) == 5);
+
+ REQUIRE(AsubsB(2, 0) == 2);
+ REQUIRE(AsubsB(2, 1) == 9);
+ REQUIRE(AsubsB(2, 2) == 8);
+ REQUIRE(AsubsB(2, 3) == 6);
+
+ SmallMatrix AmB = A - B;
+ for (size_t i = 0; i < AmB.numberOfRows(); ++i) {
+ for (size_t j = 0; j < AmB.numberOfColumns(); ++j) {
+ REQUIRE(AsubsB(i, j) == AmB(i, j));
+ }
+ }
+ }
+
+ SECTION("left product by a scalar")
+ {
+ int a = 2;
+
+ SmallMatrix aA = a * A;
+
+ {
+ bool is_same = true;
+ for (size_t i = 0; i < A.numberOfRows(); ++i) {
+ for (size_t j = 0; j < A.numberOfColumns(); ++j) {
+ is_same &= (aA(i, j) == a * A(i, j));
+ }
+ }
+ REQUIRE(is_same);
+ }
+
+ SmallMatrix copy_A = copy(A);
+ copy_A *= a;
+
+ {
+ bool is_same = true;
+ for (size_t i = 0; i < A.numberOfRows(); ++i) {
+ for (size_t j = 0; j < A.numberOfColumns(); ++j) {
+ is_same &= (aA(i, j) == copy_A(i, j));
+ }
+ }
+ REQUIRE(is_same);
+ }
+ }
+
+ SECTION("divide by a scalar")
+ {
+ SmallMatrix<double> M(3, 5);
+ M(0, 0) = 2;
+ M(0, 1) = 4;
+ M(0, 2) = 6;
+ M(0, 3) = 8;
+ M(0, 4) = 10;
+
+ M(1, 0) = 12;
+ M(1, 1) = 14;
+ M(1, 2) = 16;
+ M(1, 3) = 18;
+ M(1, 4) = 20;
+
+ M(2, 0) = 22;
+ M(2, 1) = 24;
+ M(2, 2) = 26;
+ M(2, 3) = 28;
+ M(2, 4) = 30;
+
+ M /= 2;
+
+ for (size_t i = 0; i < M.numberOfRows(); ++i) {
+ for (size_t j = 0; j < M.numberOfColumns(); ++j) {
+ REQUIRE(M(i, j) == 1 + i * M.numberOfColumns() + j);
+ }
+ }
+ }
+
+ SECTION("matrix-vector product")
+ {
+ SmallVector<int> u(4);
+
+ u[0] = 1;
+ u[1] = 3;
+ u[2] = -2;
+ u[3] = -1;
+
+ SmallVector v = B * u;
+ REQUIRE(((v[0] == -5) and (v[1] == 1) and (v[2] == -2)));
+ }
+
+ SECTION("matrix-matrix product")
+ {
+ SmallMatrix<int> BAT = B * AT;
+ REQUIRE(BAT.isSquare());
+ REQUIRE(BAT.numberOfRows() == 3);
+ REQUIRE(BAT.numberOfColumns() == 3);
+
+ REQUIRE(BAT(0, 0) == 45);
+ REQUIRE(BAT(0, 1) == 85);
+ REQUIRE(BAT(0, 2) == 125);
+
+ REQUIRE(BAT(1, 0) == 35);
+ REQUIRE(BAT(1, 1) == 75);
+ REQUIRE(BAT(1, 2) == 115);
+
+ REQUIRE(BAT(2, 0) == 42);
+ REQUIRE(BAT(2, 1) == 110);
+ REQUIRE(BAT(2, 2) == 178);
+
+ SmallMatrix<int> ATB = AT * B;
+ REQUIRE(ATB.isSquare());
+ REQUIRE(ATB.numberOfRows() == 4);
+ REQUIRE(ATB.numberOfColumns() == 4);
+
+ REQUIRE(ATB(0, 0) == 37);
+ REQUIRE(ATB(0, 1) == 44);
+ REQUIRE(ATB(0, 2) == 55);
+ REQUIRE(ATB(0, 3) == 77);
+
+ REQUIRE(ATB(1, 0) == 34);
+ REQUIRE(ATB(1, 1) == 56);
+ REQUIRE(ATB(1, 2) == 66);
+ REQUIRE(ATB(1, 3) == 94);
+
+ REQUIRE(ATB(2, 0) == 31);
+ REQUIRE(ATB(2, 1) == 68);
+ REQUIRE(ATB(2, 2) == 77);
+ REQUIRE(ATB(2, 3) == 111);
+
+ REQUIRE(ATB(3, 0) == 28);
+ REQUIRE(ATB(3, 1) == 80);
+ REQUIRE(ATB(3, 2) == 88);
+ REQUIRE(ATB(3, 3) == 128);
+ }
+
+ SECTION("output")
+ {
+ std::ostringstream out;
+ out << A;
+
+ std::ostringstream expected;
+ expected << "0| 0:1 1:2 2:3 3:4\n";
+ expected << "1| 0:5 1:6 2:7 3:8\n";
+ expected << "2| 0:9 1:10 2:11 3:12\n";
+
+ REQUIRE(out.str() == expected.str());
+ }
+
+#ifndef NDEBUG
+ SECTION("matrix-vector compatibility")
+ {
+ SmallVector<int> v(5);
+ REQUIRE_THROWS_WITH(A * v, "cannot compute matrix-vector product: incompatible sizes");
+ }
+
+ SECTION("matrix-matrix compatibility")
+ {
+ SmallMatrix C = copy(A);
+ REQUIRE_THROWS_WITH(C * A, "cannot compute matrix product: incompatible sizes");
+ REQUIRE_THROWS_WITH(C -= AT, "cannot substract matrix: incompatible sizes");
+ REQUIRE_THROWS_WITH(C += AT, "cannot add matrix: incompatible sizes");
+ REQUIRE_THROWS_WITH(C + AT, "cannot compute matrix sum: incompatible sizes");
+ REQUIRE_THROWS_WITH(C - AT, "cannot compute matrix difference: incompatible sizes");
+ }
+
+ SECTION("invalid element")
+ {
+ REQUIRE_THROWS_WITH(A(A.numberOfRows(), 0), "cannot access element: invalid indices");
+ REQUIRE_THROWS_WITH(A(0, A.numberOfColumns()), "cannot access element: invalid indices");
+ }
+
+ SECTION("invalid identity affectation")
+ {
+ SmallMatrix<double> C(3, 2);
+ REQUIRE_THROWS_WITH((C = identity), "identity must be a square matrix");
+ }
+
+ SECTION("invalid constructor")
+ {
+ REQUIRE_THROWS_WITH((SmallMatrix{2, 3, SmallArray<double>{7}}), "incompatible array size and matrix dimensions");
+ }
+
+ SECTION("output with signaling NaN")
+ {
+ SmallMatrix<double> A(2, 3);
+ A(0, 0) = 1;
+ A(0, 2) = 3;
+ A(1, 0) = 2;
+ std::ostringstream A_ost;
+ A_ost << A;
+
+ std::ostringstream ref_ost;
+ ref_ost << "0| 0:1 1:nan 2:3\n";
+ ref_ost << "1| 0:2 1:nan 2:nan\n";
+
+ REQUIRE(A_ost.str() == ref_ost.str());
+ }
+
+ // SECTION("checking for bounds violation")
+ // {
+ // REQUIRE_THROWS_AS(A(3, 0), AssertError);
+ // REQUIRE_THROWS_AS(A(0, 4), AssertError);
+
+ // REQUIRE_THROWS_AS(getMinor(A, 3, 0), AssertError);
+ // REQUIRE_THROWS_AS(getMinor(A, 0, 4), AssertError);
+
+ // const TinyMatrix<3, 4, int>& constA = A;
+ // REQUIRE_THROWS_AS(constA(3, 0), AssertError);
+ // REQUIRE_THROWS_AS(constA(0, 4), AssertError);
+ // }
+
+ // SECTION("checking for nan initialization")
+ // {
+ // TinyMatrix<3, 4, double> B;
+
+ // for (size_t i = 0; i < B.numberOfRows(); ++i) {
+ // for (size_t j = 0; j < B.numberOfColumns(); ++j) {
+ // REQUIRE(std::isnan(B(i, j)));
+ // }
+ // }
+ // }
+
+ // SECTION("checking for bad initialization")
+ // {
+ // TinyMatrix<3, 4, int> B;
+
+ // for (size_t i = 0; i < B.numberOfRows(); ++i) {
+ // for (size_t j = 0; j < B.numberOfColumns(); ++j) {
+ // REQUIRE(B(i, j) == std::numeric_limits<int>::max() / 2);
+ // }
+ // }
+ // }
+#endif // NDEBUG
+}
--
GitLab