Select Git revision
IBoundaryConditionDescriptor.hpp
-
Stéphane Del Pino authored
Now BC types are - Dirichlet - Fourier - Neumann - Symmetry One can set tags for Dirichlet, Neumann and Fourier in order to distinguish variables concerned by these conditions. In the future, it is planned to allow freefem like conditions such as `u = g on XMIN` or `alpha * u + dnu (u) = g on 3`,... However this will only make sense when "algorithms" will be set in pugs' language.
Stéphane Del Pino authoredNow BC types are - Dirichlet - Fourier - Neumann - Symmetry One can set tags for Dirichlet, Neumann and Fourier in order to distinguish variables concerned by these conditions. In the future, it is planned to allow freefem like conditions such as `u = g on XMIN` or `alpha * u + dnu (u) = g on 3`,... However this will only make sense when "algorithms" will be set in pugs' language.
test_TinyMatrix.cpp 12.14 KiB
#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/TinyMatrix.hpp>
#include <sstream>
// Instantiate to ensure full coverage is performed
template class TinyMatrix<1, 1, int>;
template class TinyMatrix<2, 2, int>;
template class TinyMatrix<3, 4, int>;
template class TinyMatrix<4, 4, double>;
// clazy:excludeall=non-pod-global-static
TEST_CASE("TinyMatrix", "[algebra]")
{
REQUIRE(TinyMatrix<1, 1, int>::Dimension == 1);
REQUIRE(TinyMatrix<1, 1, int>::NumberOfRows == 1);
REQUIRE(TinyMatrix<1, 1, int>::NumberOfColumns == 1);
REQUIRE(TinyMatrix<2, 3, int>::Dimension == 6);
REQUIRE(TinyMatrix<2, 3, int>::NumberOfRows == 2);
REQUIRE(TinyMatrix<2, 3, int>::NumberOfColumns == 3);
REQUIRE(TinyMatrix<5, 4, int>::Dimension == 20);
REQUIRE(TinyMatrix<5, 4, int>::NumberOfRows == 5);
REQUIRE(TinyMatrix<5, 4, int>::NumberOfColumns == 4);
TinyMatrix<3, 4, int> A(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12);
TinyVector<3, int> x1(1, 5, 9);
TinyVector<3, int> x2(2, 6, 10);
TinyVector<3, int> x3(3, 7, 11);
TinyVector<3, int> x4(4, 8, 12);
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)));
REQUIRE(A == TinyMatrix<3, 4, int>(x1, x2, x3, x4));
TinyMatrix<3, 4, int> B(6, 5, 3, 8, 34, 6, 35, 6, 7, 1, 3, 6);
SECTION("checking for opposed matrix")
{
const TinyMatrix minus_A = -A;
REQUIRE(
((minus_A(0, 0) == -1) and (minus_A(0, 1) == -2) and (minus_A(0, 2) == -3) and (minus_A(0, 3) == -4) and //
(minus_A(1, 0) == -5) and (minus_A(1, 1) == -6) and (minus_A(1, 2) == -7) and (minus_A(1, 3) == -8) and //
(minus_A(2, 0) == -9) and (minus_A(2, 1) == -10) and (minus_A(2, 2) == -11) and (minus_A(2, 3) == -12)));
}
SECTION("checking for equality and difference tests")
{
const TinyMatrix copy_A = A;
REQUIRE(((copy_A(0, 0) == 1) and (copy_A(0, 1) == 2) and (copy_A(0, 2) == 3) and (copy_A(0, 3) == 4) and //
(copy_A(1, 0) == 5) and (copy_A(1, 1) == 6) and (copy_A(1, 2) == 7) and (copy_A(1, 3) == 8) and //
(copy_A(2, 0) == 9) and (copy_A(2, 1) == 10) and (copy_A(2, 2) == 11) and (copy_A(2, 3) == 12)));
REQUIRE(copy_A == A);
REQUIRE_FALSE(copy_A != A);
TinyMatrix<3, 4, int> affected_A;
affected_A = A;
REQUIRE(affected_A == A);
REQUIRE_FALSE(affected_A != A);
REQUIRE(A != B);
REQUIRE_FALSE(A == B);
}
SECTION("checking for scalar left product")
{
const int a = 2;
const TinyMatrix<3, 4, int> aA = a * A;
REQUIRE(aA == TinyMatrix<3, 4, int>(2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22, 24));
}
SECTION("checking for scalar seft product")
{
const int a = 2;
TinyMatrix copy_A = A;
REQUIRE((copy_A *= a) == TinyMatrix<3, 4, int>(2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22, 24));
}
SECTION("checking for null matrix management")
{
TinyMatrix<4, 3, int> Z = zero;
REQUIRE(Z == TinyMatrix<4, 3, int>(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0));
TinyMatrix<4, 3, int> affected_Z;
affected_Z = zero;
REQUIRE(affected_Z == Z);
}
SECTION("checking for identity management")
{
TinyMatrix<3, 3, int> I = identity;
REQUIRE(I == TinyMatrix<3, 3, int>(1, 0, 0, 0, 1, 0, 0, 0, 1));
TinyMatrix<3, 3, int> affected_I;
affected_I = identity;
REQUIRE(affected_I == I);
}
SECTION("checking for matrices sum")
{
REQUIRE(A + B == TinyMatrix<3, 4, int>(7, 7, 6, 12, 39, 12, 42, 14, 16, 11, 14, 18));
TinyMatrix<3, 4, int> ApB = A;
ApB += B;
REQUIRE(ApB == A + B);
TinyMatrix<3, 4, int> Ap2B = A + 2 * B;
REQUIRE(Ap2B == ApB + B);
}
SECTION("checking for matrices difference ")
{
REQUIRE(A - B == TinyMatrix<3, 4, int>(-5, -3, 0, -4, -29, 0, -28, 2, 2, 9, 8, 6));
TinyMatrix<3, 4, int> AmB = A;
AmB -= B;
REQUIRE(AmB == A - B);
TinyMatrix<3, 4, int> Am2B = A - 2 * B;
REQUIRE(Am2B == AmB - B);
}
SECTION("checking for matrices product")
{
TinyMatrix<4, 2, int> C{3, -2, 2, 6, -2, 5, 7, 2};
REQUIRE(A * C == TinyMatrix<3, 2, int>(29, 33, 69, 77, 109, 121));
TinyMatrix<2, 3, int> D{-3, 2, 3, 2, 5, -7};
REQUIRE(D * A == TinyMatrix<2, 4, int>(34, 36, 38, 40, -36, -36, -36, -36));
}
SECTION("checking for matrix-vector product")
{
REQUIRE(A * TinyVector<4, int>(2, -3, 5, 2) == TinyVector<3, int>(19, 43, 67));
}
SECTION("checking for tensor product")
{
const TinyVector<4, int> u(1, 3, 7, 5);
const TinyVector<3, int> v(6, 2, -3);
REQUIRE(tensorProduct(u, v) == TinyMatrix<4, 3, int>(6, 2, -3, 18, 6, -9, 42, 14, -21, 30, 10, -15));
}
SECTION("checking for minor calculation")
{
REQUIRE(getMinor(A, 0, 0) == TinyMatrix<2, 3, int>(6, 7, 8, 10, 11, 12));
REQUIRE(getMinor(A, 1, 0) == TinyMatrix<2, 3, int>(2, 3, 4, 10, 11, 12));
REQUIRE(getMinor(A, 2, 0) == TinyMatrix<2, 3, int>(2, 3, 4, 6, 7, 8));
REQUIRE(getMinor(A, 0, 1) == TinyMatrix<2, 3, int>(5, 7, 8, 9, 11, 12));
REQUIRE(getMinor(A, 1, 1) == TinyMatrix<2, 3, int>(1, 3, 4, 9, 11, 12));
REQUIRE(getMinor(A, 2, 1) == TinyMatrix<2, 3, int>(1, 3, 4, 5, 7, 8));
REQUIRE(getMinor(A, 0, 2) == TinyMatrix<2, 3, int>(5, 6, 8, 9, 10, 12));
REQUIRE(getMinor(A, 1, 2) == TinyMatrix<2, 3, int>(1, 2, 4, 9, 10, 12));
REQUIRE(getMinor(A, 2, 2) == TinyMatrix<2, 3, int>(1, 2, 4, 5, 6, 8));
REQUIRE(getMinor(A, 0, 3) == TinyMatrix<2, 3, int>(5, 6, 7, 9, 10, 11));
REQUIRE(getMinor(A, 1, 3) == TinyMatrix<2, 3, int>(1, 2, 3, 9, 10, 11));
REQUIRE(getMinor(A, 2, 3) == TinyMatrix<2, 3, int>(1, 2, 3, 5, 6, 7));
}
SECTION("checking for cofactors")
{
TinyMatrix<3, 3, int> A(1, 2, 3, 4, 5, 6, 7, 8, 9);
REQUIRE(cofactor(A, 0, 0) == (5 * 9 - 8 * 6));
REQUIRE(cofactor(A, 1, 0) == -(2 * 9 - 8 * 3));
REQUIRE(cofactor(A, 2, 0) == (2 * 6 - 5 * 3));
REQUIRE(cofactor(A, 0, 1) == -(4 * 9 - 7 * 6));
REQUIRE(cofactor(A, 1, 1) == (1 * 9 - 7 * 3));
REQUIRE(cofactor(A, 2, 1) == -(1 * 6 - 4 * 3));
REQUIRE(cofactor(A, 0, 2) == (4 * 8 - 5 * 7));
REQUIRE(cofactor(A, 1, 2) == -(1 * 8 - 7 * 2));
REQUIRE(cofactor(A, 2, 2) == (1 * 5 - 4 * 2));
}
SECTION("checking for determinant calculations")
{
REQUIRE(det(TinyMatrix<1, 1, int>(6)) == 6);
REQUIRE(det(TinyMatrix<2, 2, int>(3, 1, -3, 6)) == 21);
REQUIRE(det(TinyMatrix<3, 3, int>(1, 1, 1, 1, 2, 1, 2, 1, 3)) == 1);
REQUIRE(det(TinyMatrix<3, 3, int>(6, 5, 3, 8, 34, 6, 35, 6, 7)) == -1444);
REQUIRE(det(TinyMatrix<4>(1, 2.3, 7, -6.2, 3, 4, 9, 1, 4.1, 5, 2, -3, 2, 27, 3, 17.5)) ==
Catch::Approx(6661.455).epsilon(1E-14));
REQUIRE(det(TinyMatrix<4>(1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 2, 0, 0, 2, 2)) == 0);
}
SECTION("checking for trace calculations")
{
REQUIRE(trace(TinyMatrix<1, 1, int>(6)) == 6);
REQUIRE(trace(TinyMatrix<2, 2, int>(5, 1, -3, 6)) == 5 + 6);
REQUIRE(trace(TinyMatrix<3, 3, int>(1, 1, 1, 1, 2, 1, 2, 1, 3)) == 1 + 2 + 3);
REQUIRE(trace(TinyMatrix<3, 3, int>(6, 5, 3, 8, 34, 6, 35, 6, 7)) == 6 + 34 + 7);
REQUIRE(trace(TinyMatrix<4>(1, 2.3, 7, -6.2, 3, 4, 9, 1, 4.1, 5, 2, -3, 2, 27, 3, 17.5)) ==
Catch::Approx(1 + 4 + 2 + 17.5).epsilon(1E-14));
REQUIRE(trace(TinyMatrix<4>(1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 1, 2, 0, 0, 2, 2)) == 1 + 0 + 1 + 2);
}
SECTION("checking for inverse calculations")
{
{
const TinyMatrix<1> A1(2);
REQUIRE(inverse(A1)(0, 0) == Catch::Approx(0.5).epsilon(1E-14));
}
{
const TinyMatrix<2> A2(2, 3, 4, 1);
const TinyMatrix<2> I = inverse(A2) * A2;
REQUIRE(I(0, 0) == Catch::Approx(1).epsilon(1E-14));
REQUIRE(I(0, 1) == Catch::Approx(0).margin(1E-14));
REQUIRE(I(1, 0) == Catch::Approx(0).margin(1E-14));
REQUIRE(I(1, 1) == Catch::Approx(1).epsilon(1E-14));
}
{
const TinyMatrix<3> A3(2, 3, 1, 4, -1, 5, -2, 3, 4);
const TinyMatrix<3> I = inverse(A3) * A3;
REQUIRE(I(0, 0) == Catch::Approx(1).epsilon(1E-14));
REQUIRE(I(0, 1) == Catch::Approx(0).margin(1E-14));
REQUIRE(I(0, 2) == Catch::Approx(0).margin(1E-14));
REQUIRE(I(1, 0) == Catch::Approx(0).margin(1E-14));
REQUIRE(I(1, 1) == Catch::Approx(1).epsilon(1E-14));
REQUIRE(I(1, 2) == Catch::Approx(0).margin(1E-14));
REQUIRE(I(2, 0) == Catch::Approx(0).margin(1E-14));
REQUIRE(I(2, 1) == Catch::Approx(0).margin(1E-14));
REQUIRE(I(2, 2) == Catch::Approx(1).epsilon(1E-14));
}
}
SECTION("checking for sizes")
{
REQUIRE(TinyMatrix<1>{}.numberOfRows() == 1);
REQUIRE(TinyMatrix<1>{}.numberOfColumns() == 1);
REQUIRE(TinyMatrix<1>{}.dimension() == 1);
REQUIRE(TinyMatrix<1>{}.numberOfValues() == 1);
REQUIRE(TinyMatrix<2>{}.numberOfRows() == 2);
REQUIRE(TinyMatrix<2>{}.numberOfColumns() == 2);
REQUIRE(TinyMatrix<2>{}.numberOfValues() == 4);
REQUIRE(TinyMatrix<3>{}.numberOfRows() == 3);
REQUIRE(TinyMatrix<3>{}.numberOfColumns() == 3);
REQUIRE(TinyMatrix<3>{}.dimension() == 9);
REQUIRE(TinyMatrix<3>{}.numberOfValues() == 9);
REQUIRE(TinyMatrix<3, 4>{}.numberOfRows() == 3);
REQUIRE(TinyMatrix<3, 4>{}.numberOfColumns() == 4);
REQUIRE(TinyMatrix<3, 4>{}.dimension() == 12);
REQUIRE(TinyMatrix<3, 4>{}.numberOfValues() == 12);
}
SECTION("is square matrix")
{
REQUIRE(TinyMatrix<1>{}.isSquare());
REQUIRE(TinyMatrix<2>{}.isSquare());
REQUIRE(TinyMatrix<3>{}.isSquare());
REQUIRE_FALSE(TinyMatrix<3, 4>{}.isSquare());
}
SECTION("transpose")
{
TinyMatrix tA = transpose(A);
REQUIRE(((tA(0, 0) == 1) and (tA(1, 0) == 2) and (tA(2, 0) == 3) and (tA(3, 0) == 4) and //
(tA(0, 1) == 5) and (tA(1, 1) == 6) and (tA(2, 1) == 7) and (tA(3, 1) == 8) and //
(tA(0, 2) == 9) and (tA(1, 2) == 10) and (tA(2, 2) == 11) and (tA(3, 2) == 12)));
TinyMatrix ttA = transpose(tA);
REQUIRE(ttA == A);
}
SECTION("checking for matrices output")
{
REQUIRE(Catch::Detail::stringify(A) == "[[1,2,3,4],[5,6,7,8],[9,10,11,12]]");
REQUIRE(Catch::Detail::stringify(TinyMatrix<1, 1, int>(7)) == "[[7]]");
}
SECTION("checking scalarProduct")
{
TinyMatrix<3, 4, int> B(0, -1, 1, -1, 1, -1, 1, -1, 1, -1, 1, -1);
// TinyMatrix<3, 4, int> A(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12);
REQUIRE(scalarProduct(A, B) == -7);
}
SECTION("checking norm")
{
TinyMatrix<3, 4, int> B(0, 0, 0, -1, 1, -1, 1, -1, 1, -1, 1, -1);
// TinyMatrix<3, 4, int> A(1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12);
REQUIRE(norm(B) == 3);
}
#ifndef NDEBUG
SECTION("output with signaling NaN")
{
TinyMatrix<2, 3> A;
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 << "[[1,nan,3],[2,nan,nan]]";
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);
}
}
}
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
}