Commit d84a5e30 authored by Lucas Serrano's avatar Lucas Serrano
Browse files

Added test for inversion

parent 44889340
......@@ -44,7 +44,7 @@ Matrix<T, MatrixSize, MaxVecSize> transpose(TransposedMatrix<T, MatrixSize, MaxV
}
*/
template <typename T, int NRows, int NCols, int MaxVecSize>
std::ostream & operator<<(std::ostream &os, BaseMatrix<T, NRows, NCols, MaxVecSize> &mat) {
std::ostream & operator<<(std::ostream &os, const BaseMatrix<T, NRows, NCols, MaxVecSize> &mat) {
// Matrix printing function
for (int i=0; i<NRows; i++) {
for (int j=0; j<NCols; j++) {
......@@ -243,9 +243,9 @@ class BaseMatrix {
using matrix_t = BaseMatrix<T, NumberRows, NumberCols, MaximumVectorSize>;
using pack_t = bs::pack<T, VecSize>;
using ElementType = T;
public:
using ElementType = T;
static constexpr int NRows = NumberRows;
static constexpr int NCols = NumberCols;
......@@ -306,7 +306,7 @@ class BaseMatrix {
template <typename U, int l, int m, int n, int MVS> friend void matrix_mul_mt_m(BaseMatrix<U, m, l, MVS> &a, BaseMatrix<U, m, n, MVS> &b, BaseMatrix<U, l, n, MVS> &c);
template <typename U, int l, int m, int n, int MVS> friend void matrix_mul_m_mt(BaseMatrix<U, l, m, MVS> &a, BaseMatrix<U, n, m, MVS> &b, BaseMatrix<U, l, n, MVS> &c);
friend class Inverse<T, NRows, NCols, MaxVecSize>;
friend std::ostream & operator<<<T, NRows, NCols, MaxVecSize>(std::ostream &os, matrix_t &mat);
friend std::ostream & operator<<<T, NRows, NCols, MaxVecSize>(std::ostream &os, const matrix_t &mat);
protected:
alignas(sizeof(T)*VecSize) T array[NRows*VecSize] = {0};
......
......@@ -4,14 +4,35 @@
#include <gtest/gtest-spi.h>
#include <cstdio>
#include <iomanip>
#include <cmath>
#include "../fast_matrix.hpp"
template <int size>
bool compare_float_array(const std::array<float, size> a, const std::array<float, size> b) {
bool is_equal = false;
return is_equal;
template <typename M>
::testing::AssertionResult compare_floating_point_matrices(const char* a_expr, const char* b_expr, M &a, M &b) {
using my_type = typename M::ElementType;
double precision = 1e-13;
if (typeid(my_type) == typeid(float)) {
precision = 1e-5;
}
bool is_equal = true;
std::array<my_type, M::NRows * M::NCols> array_a = a.dump_array();
std::array<my_type, M::NRows * M::NCols> array_b = b.dump_array();
for (int i = 0; i<M::NRows*M::NCols; i++) {
my_type diff = std::abs(array_a[i]-array_b[i]);
my_type largest = (array_b[i] > array_a[i]) ? array_b[i] : array_a[i];
if(diff > precision && diff > largest*precision) {
is_equal = false;
break;
}
}
if (is_equal) {
return ::testing::AssertionSuccess();
} else {
return ::testing::AssertionFailure() << a_expr << std::endl << a
<< "and " << b_expr << std::endl << b << " are not equal" << std::endl;
}
}
TEST(Inverse, Cholesky) {
......@@ -29,9 +50,10 @@ TEST(Inverse, Cholesky) {
BaseMatrix<float, 3, 3, 8> base_easy_float_mat(&base_easy_float_array[0]);
BaseMatrix<float, 3, 3, 8> inverse_easy_float_mat;
BaseMatrix<float, 3, 3, 8> result_easy_float_mat(&result_easy_float_array[0]);
Inverse<float, 3, 3, 8>::inverse(base_easy_float_mat, inverse_easy_float_mat);
//EXPECT_PRED(compare_float_array<9>, inverse_easy_float_mat.dump_array(), result_easy_float_array);
EXPECT_PRED_FORMAT2(compare_floating_point_matrices, inverse_easy_float_mat, result_easy_float_mat);
const std::array<float, 9> base_harder_float_array = {
14, 16, 11,
......@@ -45,4 +67,50 @@ TEST(Inverse, Cholesky) {
-71./529, 36./529, 38./529,
};
BaseMatrix<float, 3, 3, 8> base_harder_float_mat(&base_harder_float_array[0]);
BaseMatrix<float, 3, 3, 8> inverse_harder_float_mat;
BaseMatrix<float, 3, 3, 8> result_harder_float_mat(&result_harder_float_array[0]);
Inverse<float, 3, 3, 8>::inverse(base_harder_float_mat, inverse_harder_float_mat);
EXPECT_PRED_FORMAT2(compare_floating_point_matrices, inverse_harder_float_mat, result_harder_float_mat);
const std::array<double, 9> base_easy_double_array = {
14, 8, 3,
8, 5, 2,
3, 2, 1,
};
const std::array<double, 9> result_easy_double_array = {
1, -2, 1,
-2, 5, -4,
1, -4, 6,
};
BaseMatrix<double, 3, 3, 8> base_easy_double_mat(&base_easy_double_array[0]);
BaseMatrix<double, 3, 3, 8> inverse_easy_double_mat;
BaseMatrix<double, 3, 3, 8> result_easy_double_mat(&result_easy_double_array[0]);
Inverse<double, 3, 3, 8>::inverse(base_easy_double_mat, inverse_easy_double_mat);
EXPECT_PRED_FORMAT2(compare_floating_point_matrices, inverse_easy_double_mat, result_easy_double_mat);
const std::array<double, 9> base_harder_double_array = {
14, 16, 11,
16, 21, 10,
11, 10, 25,
};
const std::array<double, 9> result_harder_double_array = {
425./529, -290./529, -71./529,
-290./529, 229./529, 36./529,
-71./529, 36./529, 38./529,
};
BaseMatrix<double, 3, 3, 8> base_harder_double_mat(&base_harder_double_array[0]);
BaseMatrix<double, 3, 3, 8> inverse_harder_double_mat;
BaseMatrix<double, 3, 3, 8> result_harder_double_mat(&result_harder_double_array[0]);
Inverse<double, 3, 3, 8>::inverse(base_harder_double_mat, inverse_harder_double_mat);
EXPECT_PRED_FORMAT2(compare_floating_point_matrices, inverse_harder_double_mat, result_harder_double_mat);
}
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment