Skip to content
Snippets Groups Projects
Commit 15bdf87d authored by Stéphane Del Pino's avatar Stéphane Del Pino
Browse files

Fix copy for SmallMatrix and add unit tests

parent 6b1805ae
No related branches found
No related tags found
1 merge request!165Feature/algebra
...@@ -38,7 +38,7 @@ class [[nodiscard]] SmallMatrix // LCOV_EXCL_LINE ...@@ -38,7 +38,7 @@ class [[nodiscard]] SmallMatrix // LCOV_EXCL_LINE
friend PUGS_INLINE SmallMatrix<std::remove_const_t<DataType>> copy(const SmallMatrix& A) noexcept 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) friend PUGS_INLINE SmallMatrix<std::remove_const_t<DataType>> transpose(const SmallMatrix& A)
...@@ -130,7 +130,8 @@ class [[nodiscard]] SmallMatrix // LCOV_EXCL_LINE ...@@ -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>>, static_assert(std::is_same_v<std::remove_const_t<DataType>, std::remove_const_t<DataType2>>,
"incompatible data types"); "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( parallel_for(
m_values.size(), PUGS_LAMBDA(index_type i) { m_values[i] += B.m_values[i]; }); 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 ...@@ -260,6 +261,12 @@ class [[nodiscard]] SmallMatrix // LCOV_EXCL_LINE
SmallMatrix(SmallMatrix &&) = default; 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 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} : m_nb_rows{nb_rows}, m_nb_columns{nb_columns}, m_values{nb_rows * nb_columns}
{} {}
......
...@@ -131,6 +131,7 @@ add_executable (unit_tests ...@@ -131,6 +131,7 @@ add_executable (unit_tests
test_RefItemList.cpp test_RefItemList.cpp
test_RevisionInfo.cpp test_RevisionInfo.cpp
test_SmallArray.cpp test_SmallArray.cpp
test_SmallMatrix.cpp
test_SmallVector.cpp test_SmallVector.cpp
test_Socket.cpp test_Socket.cpp
test_SocketModule.cpp test_SocketModule.cpp
......
#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) != &copy_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
}
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment