Feat/block implementation

See merge request oss/rotgen!10
This commit is contained in:
Joel Falcou 2025-06-23 15:22:11 +02:00
commit 8e545dd51a
6 changed files with 406 additions and 287 deletions

View file

@ -35,163 +35,183 @@ namespace rotgen
block(parent const& base) : parent(base) {} block(parent const& base) : parent(base) {}
/* concrete_type transpose() const
block transpose() const {
{ return concrete_type(static_cast<parent const &>(*this).transpose());
return block(static_cast<parent const&>(*this).transpose()); }
}
concrete_type conjugate() const
block conjugate() const {
{ return concrete_type(static_cast<parent const &>(*this).conjugate());
return block(static_cast<parent const&>(*this).conjugate()); }
}
concrete_type adjoint() const
block adjoint() const {
{ return concrete_type(static_cast<parent const &>(*this).adjoint());
return block(static_cast<parent const&>(*this).adjoint()); }
}
*/ void transposeInPlace() { parent::transposeInPlace(); }
friend bool operator==(block const& lhs, block const& rhs) void adjointInPlace() { parent::adjointInPlace(); }
{
return static_cast<parent const&>(lhs) == static_cast<parent const&>(rhs); friend bool operator==(block const& lhs, block const& rhs)
} {
return static_cast<parent const&>(lhs) == static_cast<parent const&>(rhs);
block& operator+=(block const& rhs) }
{
static_cast<parent&>(*this) += static_cast<parent const&>(rhs); block& operator+=(block const& rhs)
return *this; {
} static_cast<parent&>(*this) += static_cast<parent const&>(rhs);
return *this;
block& operator-=(block const& rhs) }
{
static_cast<parent&>(*this) -= static_cast<parent const&>(rhs); block& operator-=(block const& rhs)
return *this; {
} static_cast<parent&>(*this) -= static_cast<parent const&>(rhs);
return *this;
concrete_type operator-() const }
{
return concrete_type(static_cast<parent const&>(*this).operator-()); concrete_type operator-() const
} {
return concrete_type(static_cast<parent const&>(*this).operator-());
block& operator*=(block const& rhs) }
{
static_cast<parent&>(*this) *= static_cast<parent const&>(rhs); concrete_type matrix_addition(block const& rhs) const
return *this; {
} return concrete_type(static_cast<parent const&>(*this).matrix_addition(rhs));
}
block& operator*=(scalar_type rhs)
{ concrete_type matrix_subtraction(block const& rhs) const
static_cast<parent&>(*this) *= rhs; {
return *this; return concrete_type(static_cast<parent const&>(*this).matrix_subtraction(rhs));
} }
block& operator/=(scalar_type rhs) concrete_type matrix_multiplication(block const& rhs) const
{ {
static_cast<parent&>(*this) /= rhs; return concrete_type(static_cast<parent const&>(*this).matrix_multiplication(rhs));
return *this; }
}
concrete_type matrix_multiplication(scalar_type rhs) const
static concrete_type Zero() {
requires (Rows != -1 && Cols != -1) return concrete_type(static_cast<parent const&>(*this).matrix_multiplication(rhs));
{ }
return parent::Zero(Rows, Cols);
} concrete_type matrix_division(scalar_type rhs) const
{
static concrete_type Zero(int rows, int cols) return concrete_type(static_cast<parent const&>(*this).matrix_division(rhs));
{ }
if constexpr(Rows != -1) assert(rows == Rows && "Mismatched between dynamic and static row size");
if constexpr(Cols != -1) assert(cols == Cols && "Mismatched between dynamic and static column size"); block& operator*=(block const& rhs)
return parent::Zero(rows, cols); {
} static_cast<parent&>(*this) *= static_cast<parent const&>(rhs);
return *this;
static concrete_type Constant(scalar_type value) }
requires (Rows != -1 && Cols != -1)
{ block& operator*=(scalar_type rhs)
return parent::Constant(Rows, Cols, static_cast<double>(value)); {
} static_cast<parent&>(*this) *= rhs;
return *this;
static concrete_type Constant(int rows, int cols, scalar_type value) }
{
if constexpr(Rows != -1) assert(rows == Rows && "Mismatched between dynamic and static row size"); block& operator/=(scalar_type rhs)
if constexpr(Cols != -1) assert(cols == Cols && "Mismatched between dynamic and static column size"); {
return parent::Constant(rows, cols, static_cast<double>(value)); static_cast<parent&>(*this) /= rhs;
} return *this;
}
static concrete_type Random()
requires (Rows != -1 && Cols != -1) static concrete_type Zero()
{ requires (Rows != -1 && Cols != -1)
return parent::Random(Rows, Cols); {
} return parent::Zero(Rows, Cols);
}
static concrete_type Random(int rows, int cols)
{ static concrete_type Zero(int rows, int cols)
if constexpr(Rows != -1) assert(rows == Rows && "Mismatched between dynamic and static row size"); {
if constexpr(Cols != -1) assert(cols == Cols && "Mismatched between dynamic and static column size"); if constexpr(Rows != -1) assert(rows == Rows && "Mismatched between dynamic and static row size");
return parent::Random(rows, cols); if constexpr(Cols != -1) assert(cols == Cols && "Mismatched between dynamic and static column size");
} return parent::Zero(rows, cols);
}
static concrete_type Identity()
requires (Rows != -1 && Cols != -1) static concrete_type Constant(scalar_type value)
{ requires (Rows != -1 && Cols != -1)
return parent::Identity(Rows, Cols); {
} return parent::Constant(Rows, Cols, static_cast<double>(value));
}
static concrete_type Identity(int rows, int cols)
{ static concrete_type Constant(int rows, int cols, scalar_type value)
if constexpr(Rows != -1) assert(rows == Rows && "Mismatched between dynamic and static row size"); {
if constexpr(Cols != -1) assert(cols == Cols && "Mismatched between dynamic and static column size"); if constexpr(Rows != -1) assert(rows == Rows && "Mismatched between dynamic and static row size");
return parent::Identity(rows, cols); if constexpr(Cols != -1) assert(cols == Cols && "Mismatched between dynamic and static column size");
} return parent::Constant(rows, cols, static_cast<double>(value));
}
template<int P>
double lpNorm() const static concrete_type Random()
{ requires (Rows != -1 && Cols != -1)
assert(P == 1 || P == 2 || P == Infinity); {
return parent::lpNorm(P); return parent::Random(Rows, Cols);
} }
};
static concrete_type Random(int rows, int cols)
/* {
template<typename S, int R, int C, int O, int MR, int MC> if constexpr(Rows != -1) assert(rows == Rows && "Mismatched between dynamic and static row size");
block<S,R,C,O,MR,MC> operator+(block<S,R,C,O,MR,MC> const& lhs, block<S,R,C,O,MR,MC> const& rhs) if constexpr(Cols != -1) assert(cols == Cols && "Mismatched between dynamic and static column size");
{ return parent::Random(rows, cols);
block<S,R,C,O,MR,MC> that(lhs); }
return that += rhs;
} static concrete_type Identity()
requires (Rows != -1 && Cols != -1)
template<typename S, int R, int C, int O, int MR, int MC> {
block<S,R,C,O,MR,MC> operator-(block<S,R,C,O,MR,MC> const& lhs, block<S,R,C,O,MR,MC> const& rhs) return parent::Identity(Rows, Cols);
{ }
block<S,R,C,O,MR,MC> that(lhs);
return that -= rhs; static concrete_type Identity(int rows, int cols)
} {
if constexpr(Rows != -1) assert(rows == Rows && "Mismatched between dynamic and static row size");
template<typename S, int R, int C, int O, int MR, int MC> if constexpr(Cols != -1) assert(cols == Cols && "Mismatched between dynamic and static column size");
block<S,R,C,O,MR,MC> operator*(block<S,R,C,O,MR,MC> const& lhs, block<S,R,C,O,MR,MC> const& rhs) return parent::Identity(rows, cols);
{ }
block<S,R,C,O,MR,MC> that(lhs);
return that *= rhs; template<int P>
} double lpNorm() const
{
template<typename S, int R, int C, int O, int MR, int MC> assert(P == 1 || P == 2 || P == Infinity);
block<S,R,C,O,MR,MC> operator*(block<S,R,C,O,MR,MC> const& lhs, double rhs) return parent::lpNorm(P);
{ }
block<S,R,C,O,MR,MC> that(lhs); };
return that *= rhs;
} template<typename Ref, int R, int C, bool I, int F>
typename block<Ref, R, C, I, F>::concrete_type operator+(block<Ref, R, C, I, F> const& lhs, block<Ref, R, C, I, F> const& rhs)
template<typename S, int R, int C, int O, int MR, int MC> {
block<S,R,C,O,MR,MC> operator*(double lhs, block<S,R,C,O,MR,MC> const& rhs) return lhs.matrix_addition(rhs);
{ }
return rhs * lhs;
} template<typename Ref, int R, int C, bool I, int F>
typename block<Ref, R, C, I, F>::concrete_type operator-(block<Ref, R, C, I, F> const& lhs, block<Ref, R, C, I, F> const& rhs)
template<typename S, int R, int C, int O, int MR, int MC> {
block<S,R,C,O,MR,MC> operator/(block<S,R,C,O,MR,MC> const& lhs, double rhs) return lhs.matrix_subtraction(rhs);
{ }
block<S,R,C,O,MR,MC> that(lhs);
return that /= rhs; template<typename Ref, int R, int C, bool I, int F>
} typename block<Ref, R, C, I, F>::concrete_type operator*(block<Ref, R, C, I, F> const& lhs, block<Ref, R, C, I, F> const& rhs)
*/ {
return lhs.matrix_multiplication(rhs);
}
template<typename Ref, int R, int C, bool I, int F>
typename block<Ref, R, C, I, F>::concrete_type operator*(block<Ref, R, C, I, F> const& lhs, double rhs)
{
return lhs.matrix_multiplication(rhs);
}
template<typename Ref, int R, int C, bool I, int F>
typename block<Ref, R, C, I, F>::concrete_type operator*(double lhs, block<Ref, R, C, I, F> const& rhs)
{
return rhs.matrix_multiplication(lhs);
}
template<typename Ref, int R, int C, bool I, int F>
typename block<Ref, R, C, I, F>::concrete_type operator/(block<Ref, R, C, I, F> const& lhs, double rhs)
{
return lhs.matrix_division(rhs);
}
} }

View file

@ -28,12 +28,12 @@ class CLASSNAME
std::size_t cols() const; std::size_t cols() const;
std::size_t size() const; std::size_t size() const;
// CLASSNAME transpose() const; SOURCENAME transpose() const;
// CLASSNAME conjugate() const; SOURCENAME conjugate() const;
// CLASSNAME adjoint() const; SOURCENAME adjoint() const;
// void transposeInPlace(); void transposeInPlace();
// void adjointInPlace(); void adjointInPlace();
TYPE sum() const; TYPE sum() const;
TYPE prod() const; TYPE prod() const;
@ -61,6 +61,12 @@ class CLASSNAME
CLASSNAME& operator*=(TYPE d); CLASSNAME& operator*=(TYPE d);
CLASSNAME& operator/=(TYPE d); CLASSNAME& operator/=(TYPE d);
SOURCENAME matrix_addition(CLASSNAME const& rhs) const;
SOURCENAME matrix_subtraction(CLASSNAME const& rhs) const;
SOURCENAME matrix_multiplication(CLASSNAME const& rhs) const;
SOURCENAME matrix_multiplication(TYPE s) const;
SOURCENAME matrix_division(TYPE s) const;
friend std::ostream& operator<<(std::ostream&,CLASSNAME const&); friend std::ostream& operator<<(std::ostream&,CLASSNAME const&);
friend bool operator==(CLASSNAME const& lhs, CLASSNAME const& rhs); friend bool operator==(CLASSNAME const& lhs, CLASSNAME const& rhs);

View file

@ -27,6 +27,8 @@ namespace rotgen
payload(data_type&& matrix) : data(std::move(matrix)) {} payload(data_type&& matrix) : data(std::move(matrix)) {}
void assign(Eigen::Block<data_type> ref) { data = ref; } void assign(Eigen::Block<data_type> ref) { data = ref; }
void assign(data_type const& mat) { data = mat; }
void assign(data_type&& mat) { data = std::move(mat); }
}; };
struct matrix_impl64_row::payload struct matrix_impl64_row::payload
@ -39,6 +41,8 @@ namespace rotgen
payload(data_type&& matrix) : data(std::move(matrix)) {} payload(data_type&& matrix) : data(std::move(matrix)) {}
void assign(Eigen::Block<data_type> ref) { data = ref; } void assign(Eigen::Block<data_type> ref) { data = ref; }
void assign(data_type const& mat) { data = mat; }
void assign(data_type&& mat) { data = std::move(mat); }
}; };
struct matrix_impl32_col::payload struct matrix_impl32_col::payload
@ -51,6 +55,8 @@ namespace rotgen
payload(data_type&& matrix) : data(std::move(matrix)) {} payload(data_type&& matrix) : data(std::move(matrix)) {}
void assign(Eigen::Block<data_type> ref) { data = ref; } void assign(Eigen::Block<data_type> ref) { data = ref; }
void assign(data_type const& mat) { data = mat; }
void assign(data_type&& mat) { data = std::move(mat); }
}; };
struct matrix_impl32_row::payload struct matrix_impl32_row::payload
@ -63,5 +69,7 @@ namespace rotgen
payload(data_type&& matrix) : data(std::move(matrix)) {} payload(data_type&& matrix) : data(std::move(matrix)) {}
void assign(Eigen::Block<data_type> ref) { data = ref; } void assign(Eigen::Block<data_type> ref) { data = ref; }
void assign(data_type const& mat) { data = mat; }
void assign(data_type&& mat) { data = std::move(mat); }
}; };
} }

View file

@ -21,6 +21,7 @@ struct CLASSNAME::payload
using data_type = Eigen::Block<base_type>; using data_type = Eigen::Block<base_type>;
data_type data; data_type data;
payload (data_type const& o) : data(o) {} payload (data_type const& o) : data(o) {}
payload (base_type& r, std::size_t i0, std::size_t j0, std::size_t ni, std::size_t nj) payload (base_type& r, std::size_t i0, std::size_t j0, std::size_t ni, std::size_t nj)
@ -59,114 +60,144 @@ struct CLASSNAME::payload
TYPE& CLASSNAME::operator()(std::size_t i, std::size_t j) { return storage_->data(i,j); } TYPE& CLASSNAME::operator()(std::size_t i, std::size_t j) { return storage_->data(i,j); }
TYPE const& CLASSNAME::operator()(std::size_t i, std::size_t j) const { return storage_->data(i,j); } TYPE const& CLASSNAME::operator()(std::size_t i, std::size_t j) const { return storage_->data(i,j); }
/* /*
TYPE& CLASSNAME::operator()(std::size_t index) { return storage_->data(index); } TYPE& CLASSNAME::operator()(std::size_t index) { return storage_->data(index); }
TYPE const& CLASSNAME::operator()(std::size_t index) const { return storage_->data(index); } TYPE const& CLASSNAME::operator()(std::size_t index) const { return storage_->data(index); }
*/ */
const TYPE* CLASSNAME::data() const { return storage_->data.data(); } const TYPE* CLASSNAME::data() const { return storage_->data.data(); }
/* SOURCENAME CLASSNAME::transpose() const
CLASSNAME CLASSNAME::transpose() const {
{ SOURCENAME result;
CLASSNAME result(*this); result.storage()->assign(storage_->data.transpose().eval());
result.storage_->data.transposeInPlace(); return result;
return result; }
} SOURCENAME CLASSNAME::conjugate() const {
SOURCENAME result;
result.storage()->assign(storage_->data.conjugate().eval());
return result;
}
CLASSNAME CLASSNAME::conjugate() const SOURCENAME CLASSNAME::adjoint() const {
{ SOURCENAME result;
CLASSNAME result(*this); result.storage()->assign(storage_->data.adjoint().eval());
result.storage_->data = storage_->data.conjugate(); return result;
return result; }
}
CLASSNAME CLASSNAME::adjoint() const void CLASSNAME::transposeInPlace()
{ {
CLASSNAME result(*this); storage_->data.transposeInPlace();
result.storage_->data.adjointInPlace(); }
return result;
}
void CLASSNAME::transposeInPlace() void CLASSNAME::adjointInPlace()
{ {
storage_->data.transposeInPlace(); storage_->data.adjointInPlace();
} }
void CLASSNAME::adjointInPlace() TYPE CLASSNAME::sum() const { return storage_->data.sum(); }
{ TYPE CLASSNAME::prod() const { return storage_->data.prod(); }
storage_->data.adjointInPlace(); TYPE CLASSNAME::mean() const { return storage_->data.mean(); }
} TYPE CLASSNAME::trace() const { return storage_->data.trace(); }
*/
TYPE CLASSNAME::sum() const { return storage_->data.sum(); } TYPE CLASSNAME::minCoeff() const { return storage_->data.minCoeff(); }
TYPE CLASSNAME::prod() const { return storage_->data.prod(); } TYPE CLASSNAME::maxCoeff() const { return storage_->data.maxCoeff(); }
TYPE CLASSNAME::mean() const { return storage_->data.mean(); }
TYPE CLASSNAME::trace() const { return storage_->data.trace(); }
TYPE CLASSNAME::minCoeff() const { return storage_->data.minCoeff(); } TYPE CLASSNAME::minCoeff(std::ptrdiff_t* row, std::ptrdiff_t* col) const { return storage_->data.minCoeff(row, col); }
TYPE CLASSNAME::maxCoeff() const { return storage_->data.maxCoeff(); } TYPE CLASSNAME::maxCoeff(std::ptrdiff_t* row, std::ptrdiff_t* col) const { return storage_->data.maxCoeff(row, col); }
TYPE CLASSNAME::minCoeff(std::ptrdiff_t* row, std::ptrdiff_t* col) const { return storage_->data.minCoeff(row, col); } TYPE CLASSNAME::squaredNorm() const { return storage_->data.squaredNorm(); }
TYPE CLASSNAME::maxCoeff(std::ptrdiff_t* row, std::ptrdiff_t* col) const { return storage_->data.maxCoeff(row, col); } TYPE CLASSNAME::norm() const { return storage_->data.norm(); }
TYPE CLASSNAME::squaredNorm() const { return storage_->data.squaredNorm(); } TYPE CLASSNAME::lpNorm(int p) const
TYPE CLASSNAME::norm() const { return storage_->data.norm(); } {
if (p == 1) return storage_->data.lpNorm<1>();
else if (p == 2) return storage_->data.lpNorm<2>();
else return storage_->data.lpNorm<Eigen::Infinity>();
}
TYPE CLASSNAME::lpNorm(int p) const //==================================================================================================
{ // Operators
if (p == 1) return storage_->data.lpNorm<1>(); //==================================================================================================
else if (p == 2) return storage_->data.lpNorm<2>(); std::ostream& operator<<(std::ostream& os,CLASSNAME const& m)
else return storage_->data.lpNorm<Eigen::Infinity>(); {
} return os << m.storage_->data;
}
//================================================================================================== bool operator==(CLASSNAME const& lhs, CLASSNAME const& rhs)
// Operators {
//================================================================================================== return lhs.storage_->data == rhs.storage_->data;
std::ostream& operator<<(std::ostream& os,CLASSNAME const& m) }
{
return os << m.storage_->data;
}
bool operator==(CLASSNAME const& lhs, CLASSNAME const& rhs) CLASSNAME& CLASSNAME::operator+=(CLASSNAME const& rhs)
{ {
return lhs.storage_->data == rhs.storage_->data; storage_->data += rhs.storage_->data;
} return *this;
}
CLASSNAME& CLASSNAME::operator+=(CLASSNAME const& rhs) CLASSNAME& CLASSNAME::operator-=(CLASSNAME const& rhs)
{ {
storage_->data += rhs.storage_->data; storage_->data -= rhs.storage_->data;
return *this; return *this;
} }
CLASSNAME& CLASSNAME::operator-=(CLASSNAME const& rhs) SOURCENAME CLASSNAME::operator-() const
{ {
storage_->data -= rhs.storage_->data; SOURCENAME result;
return *this; result.storage()->assign(storage_->data);
} return -result;
}
SOURCENAME CLASSNAME::operator-() const CLASSNAME& CLASSNAME::operator*=(CLASSNAME const& rhs)
{ {
SOURCENAME result; storage_->data *= rhs.storage_->data;
result.storage()->assign(storage_->data); return *this;
return -result; }
}
CLASSNAME& CLASSNAME::operator*=(CLASSNAME const& rhs) CLASSNAME& CLASSNAME::operator*=(TYPE s)
{ {
storage_->data *= rhs.storage_->data; storage_->data *= s;
return *this; return *this;
} }
CLASSNAME& CLASSNAME::operator*=(TYPE s) CLASSNAME& CLASSNAME::operator/=(TYPE s)
{ {
storage_->data *= s; storage_->data /= s;
return *this; return *this;
} }
CLASSNAME& CLASSNAME::operator/=(TYPE s) SOURCENAME CLASSNAME::matrix_addition(CLASSNAME const& rhs) const
{ {
storage_->data /= s; SOURCENAME result;
return *this; result.storage()->assign(storage_->data + rhs.storage_->data);
} return result;
}
SOURCENAME CLASSNAME::matrix_subtraction(CLASSNAME const& rhs) const
{
SOURCENAME result;
result.storage()->assign(storage_->data - rhs.storage_->data);
return result;
}
SOURCENAME CLASSNAME::matrix_multiplication(CLASSNAME const& rhs) const
{
SOURCENAME result;
result.storage()->assign(storage_->data * rhs.storage_->data);
return result;
}
SOURCENAME CLASSNAME::matrix_multiplication(TYPE s) const
{
SOURCENAME result;
result.storage()->assign(storage_->data * s);
return result;
}
SOURCENAME CLASSNAME::matrix_division(TYPE s) const
{
SOURCENAME result;
result.storage()->assign(storage_->data / s);
return result;
}

View file

@ -11,6 +11,37 @@
#include <rotgen/extract.hpp> #include <rotgen/extract.hpp>
#include <Eigen/Dense> #include <Eigen/Dense>
template<typename Type1, typename Type2>
void test_comparison(const Type1& t1, const Type2& t2)
{
TTS_EQUAL(static_cast<std::ptrdiff_t>(t1.rows()), static_cast<std::ptrdiff_t>(t2.rows()));
TTS_EQUAL(static_cast<std::ptrdiff_t>(t1.cols()), static_cast<std::ptrdiff_t>(t2.cols()));
for (std::size_t r = 0; r < static_cast<std::size_t>(t1.rows()); ++r)
for (std::size_t c = 0; c < static_cast<std::size_t>(t1.cols()); ++c)
TTS_EQUAL(t1(r, c), t2(r, c));
}
template<typename Matrix1, typename Matrix2, typename Block1, typename Block2>
void test_block_unary_ops(const Matrix1& original_matrix, const Matrix2& ref_matrix,
Block1 original_block, Block2 ref_block)
{
test_comparison(original_block.transpose(), ref_block.transpose());
test_comparison(original_block.conjugate(), ref_block.conjugate());
test_comparison(original_block.adjoint(), ref_block.adjoint());
if (original_block.rows() == original_block.cols()) {
original_block.transposeInPlace();
ref_block.transposeInPlace();
test_comparison(original_block, ref_block);
test_comparison(original_matrix, ref_matrix);
original_block.adjointInPlace();
ref_block.adjointInPlace();
test_comparison(original_block, ref_block);
test_comparison(original_matrix, ref_matrix);
}
}
template<typename Block1, typename Block2> template<typename Block1, typename Block2>
void compare_reductions(const Block1& block, const Block2& ref) void compare_reductions(const Block1& block, const Block2& ref)
{ {
@ -68,8 +99,11 @@ void test_dynamic_block_reductions(rotgen::tests::matrix_block_test_case<MatrixT
ref.bottomRows(matrix_construct.ni) }, ref.bottomRows(matrix_construct.ni) },
}; };
for (const auto& [original_block, ref_block] : test_cases) for (const auto& [original_block, ref_block] : test_cases) {
compare_reductions(original_block, ref_block); compare_reductions(original_block, ref_block);
test_block_unary_ops(original, ref, original_block, ref_block);
}
} }
TTS_CASE_TPL("Test dynamic block reductions", rotgen::tests::types) TTS_CASE_TPL("Test dynamic block reductions", rotgen::tests::types)
@ -79,19 +113,19 @@ TTS_CASE_TPL("Test dynamic block reductions", rotgen::tests::types)
std::vector<rotgen::tests::matrix_block_test_case<mat_t>> test_cases = std::vector<rotgen::tests::matrix_block_test_case<mat_t>> test_cases =
{ {
{6, 5, [](auto r, auto c) {return r + c; }, 1, 2, 3, 2}, {6, 5, [](auto r, auto c) { return T(r + c); }, 1, 2, 3, 2},
{9, 11, [](auto r, auto c) {return r + c; }, 0, 1, 4, 9}, {9, 11, [](auto r, auto c) {return T(r + c); }, 0, 1, 4, 9},
{3, 3, [](auto , auto ) {return 0.0; }, 1, 1, 1, 1}, {3, 3, [](auto , auto ) {return T(0.0); }, 1, 1, 1, 1},
{1, 4, [](auto r, auto c) {return -r -c*c - 1234; }, 0, 0, 1, 1}, {1, 4, [](auto r, auto c) {return T(-r -c*c - 1234); }, 0, 0, 1, 1},
{4, 1, [](auto , auto ) {return 7.0; }, 2, 0, 2, 1}, {9, 9, [](auto r, auto c) {return T(-r + 2*c); }, 0, 1, 3, 3},
{1, 1, [](auto , auto ) {return 42.0; }, 0, 0, 1, 1}, {11, 13, [](auto r, auto c) {return T(std::tan(r+c)); }, 1, 1, 2, 2},
{12, 13, [](auto r, auto c) {return std::sin(r + c); }, 2, 3, 4, 5 }, {4, 1, [](auto , auto ) {return T(7.0); }, 2, 0, 2, 1},
{4, 9, [](auto r, auto c) {return -1.5 * r + 2.56 * c; }, 0, 1, 2, 3 }, {1, 1, [](auto , auto ) {return T(42.0); }, 0, 0, 1, 1},
{2, 5, [](auto r, auto c) {return (r == c ? 1.0 : 0.0); }, 1, 1, 1, 1}, {12, 13, [](auto r, auto c) {return T(std::sin(r + c)); }, 2, 3, 4, 5 },
{4, 9, [](auto r, auto c) {return T(-1.5 * r + 2.56 * c); }, 0, 1, 2, 3 },
{2, 5, [](auto r, auto c) {return T(r == c ? 1.0 : 0.0); }, 1, 1, 1, 1},
}; };
for (const auto& test_case : test_cases) for (const auto& test_case : test_cases)
test_dynamic_block_reductions<mat_t, T>(test_case); test_dynamic_block_reductions<mat_t, T>(test_case);
}; };

View file

@ -12,7 +12,7 @@
template <typename MatrixType, typename T> template <typename MatrixType, typename T>
void test_block_matrix_operations(rotgen::tests::matrix_block_test_case<MatrixType> const& matrix_construct, void test_block_matrix_operations(rotgen::tests::matrix_block_test_case<MatrixType> const& matrix_construct,
auto b_init_fn, auto self_ops) auto b_init_fn, auto ops, auto self_ops)
{ {
using EigenMatrix = Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic>; using EigenMatrix = Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic>;
@ -39,12 +39,12 @@ void test_block_matrix_operations(rotgen::tests::matrix_block_test_case<MatrixTy
auto ref_b_block = ref_b.block(matrix_construct.i0, matrix_construct.j0, auto ref_b_block = ref_b.block(matrix_construct.i0, matrix_construct.j0,
matrix_construct.ni, matrix_construct.nj); matrix_construct.ni, matrix_construct.nj);
/*auto result_block = ops(a_block, b_block); auto result_block = ops(a_block, b_block);
auto ref_result_block = ops(ref_a_block, ref_b_block); auto ref_result_block = ops(ref_a_block, ref_b_block);
for (std::size_t r = 0; r < matrix_construct.ni; ++r) for (std::size_t r = 0; r < matrix_construct.ni; ++r)
for (std::size_t c = 0; c < matrix_construct.nj; ++c) for (std::size_t c = 0; c < matrix_construct.nj; ++c)
TTS_EQUAL(result_block(r, c), ref_result_block(r, c));*/ TTS_EQUAL(result_block(r, c), ref_result_block(r, c));
self_ops(a_block,b_block); self_ops(a_block,b_block);
self_ops(ref_a_block, ref_b_block); self_ops(ref_a_block, ref_b_block);
@ -65,7 +65,7 @@ void test_block_matrix_operations(rotgen::tests::matrix_block_test_case<MatrixTy
template <typename MatrixType, typename T> template <typename MatrixType, typename T>
void test_block_scalar_operations(rotgen::tests::matrix_block_test_case<MatrixType> const& matrix_construct, void test_block_scalar_operations(rotgen::tests::matrix_block_test_case<MatrixType> const& matrix_construct,
auto scalar, auto self_ops) auto scalar, auto ops, auto self_ops)
{ {
using EigenMatrix = Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic>; using EigenMatrix = Eigen::Matrix<T, Eigen::Dynamic, Eigen::Dynamic>;
@ -81,7 +81,12 @@ void test_block_scalar_operations(rotgen::tests::matrix_block_test_case<MatrixTy
auto ref_a_block = ref_a.block(matrix_construct.i0, matrix_construct.j0, auto ref_a_block = ref_a.block(matrix_construct.i0, matrix_construct.j0,
matrix_construct.ni, matrix_construct.nj); matrix_construct.ni, matrix_construct.nj);
//TTS_EQUAL(ops(a, s), ref); auto result = ops(a_block, scalar);
auto ref_result = ops(ref_a_block, scalar);
for (std::size_t r = 0; r < matrix_construct.ni; ++r)
for (std::size_t c = 0; c < matrix_construct.nj; ++c)
TTS_EQUAL(result(r, c), ref_result(r, c));
self_ops(a_block,scalar); self_ops(a_block,scalar);
self_ops(ref_a_block, scalar); self_ops(ref_a_block, scalar);
@ -113,8 +118,19 @@ void test_scalar_block_multiplications(rotgen::tests::matrix_block_test_case<Mat
auto ref_a_block = ref_a.block(matrix_construct.i0, matrix_construct.j0, auto ref_a_block = ref_a.block(matrix_construct.i0, matrix_construct.j0,
matrix_construct.ni, matrix_construct.nj); matrix_construct.ni, matrix_construct.nj);
//TTS_EQUAL(a * s, ref); auto a_scalar_multiplication = a_block * scalar;
//TTS_EQUAL(s * a, ref); auto scalar_a_multiplication = scalar * a_block;
auto a_scalar_multiplication_ref = ref_a_block * scalar;
auto scalar_a_multiplication_ref = scalar * ref_a_block;
for (std::size_t r = 0; r < matrix_construct.ni; ++r)
{
for (std::size_t c = 0; c < matrix_construct.nj; ++c)
{
TTS_EQUAL(a_scalar_multiplication (r, c), a_scalar_multiplication_ref(r, c));
TTS_EQUAL(scalar_a_multiplication(r, c), scalar_a_multiplication_ref(r, c));
}
}
a_block *= scalar; a_block *= scalar;
ref_a_block *= scalar; ref_a_block *= scalar;
@ -157,7 +173,13 @@ void test_block_multiplication(rotgen::tests::matrix_block_test_case<MatrixType>
auto ref_b_block = ref_b.block(b_matrix_construct.i0, b_matrix_construct.j0, auto ref_b_block = ref_b.block(b_matrix_construct.i0, b_matrix_construct.j0,
b_matrix_construct.ni, b_matrix_construct.nj); b_matrix_construct.ni, b_matrix_construct.nj);
// a * b auto a_b_product_original = a_block * b_block;
auto a_b_product_ref = ref_a_block * ref_b_block;
for (std::size_t r = 0; r < a_matrix_construct.ni; ++r)
for (std::size_t c = 0; c < a_matrix_construct.nj; ++c)
TTS_EQUAL(a_b_product_original (r, c), a_b_product_ref(r, c));
a_block *= b_block; a_block *= b_block;
ref_a_block *= ref_b_block; ref_a_block *= ref_b_block;
@ -183,19 +205,18 @@ TTS_CASE_TPL("Check block addition", rotgen::tests::types)
<typename T, typename O>( tts::type< tts::types<T,O>> ) <typename T, typename O>( tts::type< tts::types<T,O>> )
{ {
using mat_t = rotgen::matrix<T,rotgen::Dynamic,rotgen::Dynamic,O::value>; using mat_t = rotgen::matrix<T,rotgen::Dynamic,rotgen::Dynamic,O::value>;
//auto op = [](auto a, auto b) { return a + b; }; auto op = [](auto a, auto b) { return a + b; };
auto s_op = [](auto& a, auto b) { return a += b; }; auto s_op = [](auto& a, auto b) { return a += b; };
//test_block_matrix_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, init_b, op, s_op); test_block_matrix_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, init_b, s_op); test_block_matrix_operations<mat_t, T>({13 , 15, init_a, 1, 2, 3, 4}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({13 , 15, init_a, 1, 2, 3, 4}, init_b, s_op); test_block_matrix_operations<mat_t, T>({5 , 9, init_a, 2, 2, 2, 2}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({5 , 9, init_a, 2, 2, 2, 2}, init_b, s_op); test_block_matrix_operations<mat_t, T>({15 , 15, init_a, 3, 4, 5, 5}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({15 , 15, init_a, 3, 4, 5, 5}, init_b, s_op); test_block_matrix_operations<mat_t, T>({5 , 5, init_b, 1, 0, 3, 2}, init_a, op, s_op);
test_block_matrix_operations<mat_t, T>({5 , 5, init_b, 1, 0, 3, 2}, init_a, s_op); test_block_matrix_operations<mat_t, T>({10, 1, init_a, 0, 0, 5, 1}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({10, 1, init_a, 0, 0, 5, 1}, init_b, s_op); test_block_matrix_operations<mat_t, T>({1 , 10, init_a, 0, 0, 1, 5}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({1 , 10, init_a, 0, 0, 1, 5}, init_b, s_op); test_block_matrix_operations<mat_t, T>({21 , 5, init_0, 4, 4, 10, 1}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({21 , 5, init_0, 4, 4, 10, 1}, init_b, s_op); test_block_matrix_operations<mat_t, T>({11 , 7, init_a, 2, 0, 7, 5}, init_0, op, s_op);
test_block_matrix_operations<mat_t, T>({11 , 7, init_a, 2, 0, 7, 5}, init_0, s_op);
}; };
@ -203,19 +224,19 @@ TTS_CASE_TPL("Check block subtraction", rotgen::tests::types)
<typename T, typename O>( tts::type< tts::types<T,O>> ) <typename T, typename O>( tts::type< tts::types<T,O>> )
{ {
using mat_t = rotgen::matrix<T,rotgen::Dynamic,rotgen::Dynamic,O::value>; using mat_t = rotgen::matrix<T,rotgen::Dynamic,rotgen::Dynamic,O::value>;
//auto op = [](auto a, auto b) { return a - b; }; auto op = [](auto a, auto b) { return a - b; };
auto s_op = [](auto& a, auto b) { return a -= b; }; auto s_op = [](auto& a, auto b) { return a -= b; };
//test_block_matrix_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, init_b, op, s_op); test_block_matrix_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, init_b, s_op); test_block_matrix_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({13 , 15, init_a, 1, 2, 3, 4}, init_b, s_op); test_block_matrix_operations<mat_t, T>({13 , 15, init_a, 1, 2, 3, 4}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({5 , 9, init_a, 2, 2, 2, 2}, init_b, s_op); test_block_matrix_operations<mat_t, T>({5 , 9, init_a, 2, 2, 2, 2}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({15 , 15, init_a, 3, 4, 5, 5}, init_b, s_op); test_block_matrix_operations<mat_t, T>({15 , 15, init_a, 3, 4, 5, 5}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({5 , 5, init_b, 1, 0, 3, 2}, init_a, s_op); test_block_matrix_operations<mat_t, T>({5 , 5, init_b, 1, 0, 3, 2}, init_a, op, s_op);
test_block_matrix_operations<mat_t, T>({10, 1, init_a, 0, 0, 5, 1}, init_b, s_op); test_block_matrix_operations<mat_t, T>({10, 1, init_a, 0, 0, 5, 1}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({1 , 10, init_a, 0, 0, 1, 5}, init_b,s_op); test_block_matrix_operations<mat_t, T>({1 , 10, init_a, 0, 0, 1, 5}, init_b,op, s_op);
test_block_matrix_operations<mat_t, T>({21 , 5, init_0, 4, 4, 10, 1}, init_b, s_op); test_block_matrix_operations<mat_t, T>({21 , 5, init_0, 4, 4, 10, 1}, init_b, op, s_op);
test_block_matrix_operations<mat_t, T>({11 , 7, init_a, 2, 0, 7, 5}, init_0, s_op); test_block_matrix_operations<mat_t, T>({11 , 7, init_a, 2, 0, 7, 5}, init_0, op, s_op);
}; };
@ -258,17 +279,16 @@ TTS_CASE_TPL("Check block division with scalar", rotgen::tests::types)
<typename T, typename O>( tts::type< tts::types<T,O>> ) <typename T, typename O>( tts::type< tts::types<T,O>> )
{ {
using mat_t = rotgen::matrix<T,rotgen::Dynamic,rotgen::Dynamic,O::value>; using mat_t = rotgen::matrix<T,rotgen::Dynamic,rotgen::Dynamic,O::value>;
//auto op = [](auto a, auto b) { return a / b; }; auto op = [](auto a, auto b) { return a / b; };
auto s_op = [](auto& a, auto b) { return a /= b; }; auto s_op = [](auto& a, auto b) { return a /= b; };
//test_block_scalar_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, T{ 3.5}, op, s_op); test_block_scalar_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, T{ 3.5}, op, s_op);
test_block_scalar_operations<mat_t, T>({1 , 1, init_a, 0, 0, 1, 1}, T{ 3.5}, s_op); test_block_scalar_operations<mat_t, T>({13 , 15, init_a, 1, 2, 3, 4}, T{-2.5}, op, s_op);
test_block_scalar_operations<mat_t, T>({13 , 15, init_a, 1, 2, 3, 4}, T{-2.5}, s_op); test_block_scalar_operations<mat_t, T>({5 , 9, init_a, 2, 2, 2, 2}, T{ 42. }, op, s_op);
test_block_scalar_operations<mat_t, T>({5 , 9, init_a, 2, 2, 2, 2}, T{ 42. }, s_op); test_block_scalar_operations<mat_t, T>({15 , 15, init_a, 3, 4, 5, 5}, T{-5. }, op, s_op);
test_block_scalar_operations<mat_t, T>({15 , 15, init_a, 3, 4, 5, 5}, T{-5. }, s_op); test_block_scalar_operations<mat_t, T>({5 , 5, init_b, 1, 0, 3, 2}, T{ 1. }, op, s_op);
test_block_scalar_operations<mat_t, T>({5 , 5, init_b, 1, 0, 3, 2}, T{ 1. }, s_op); test_block_scalar_operations<mat_t, T>({10, 1, init_a, 0, 0, 5, 1}, T{ 0. }, op, s_op);
test_block_scalar_operations<mat_t, T>({10, 1, init_a, 0, 0, 5, 1}, T{ 0. }, s_op); test_block_scalar_operations<mat_t, T>({1 , 10, init_a, 0, 0, 1, 5}, T{ 6. }, op, s_op);
test_block_scalar_operations<mat_t, T>({1 , 10, init_a, 0, 0, 1, 5}, T{ 6. },s_op); test_block_scalar_operations<mat_t, T>({21 , 5, init_0, 4, 4, 10, 1}, T{ 10.}, op, s_op);
test_block_scalar_operations<mat_t, T>({21 , 5, init_0, 4, 4, 10, 1}, T{ 10.}, s_op); test_block_scalar_operations<mat_t, T>({11 , 7, init_a, 2, 0, 7, 5}, T{-0.5}, op, s_op);
test_block_scalar_operations<mat_t, T>({11 , 7, init_a, 2, 0, 7, 5}, T{-0.5}, s_op);
}; };