33 #ifndef LINALG_BACKEND_EIGEN_H__ 34 #define LINALG_BACKEND_EIGEN_H__ 48 #define DEFINE_FOR_ALL_PTYPE(METHODNAME, Container) \ 49 METHODNAME(bool, Container); \ 50 METHODNAME(char, Container); \ 51 METHODNAME(int8_t, Container); \ 52 METHODNAME(uint8_t, Container); \ 53 METHODNAME(int16_t, Container); \ 54 METHODNAME(uint16_t, Container); \ 55 METHODNAME(int32_t, Container); \ 56 METHODNAME(uint32_t, Container); \ 57 METHODNAME(int64_t, Container); \ 58 METHODNAME(uint64_t, Container); \ 59 METHODNAME(float32_t, Container); \ 60 METHODNAME(float64_t, Container); \ 61 METHODNAME(floatmax_t, Container); \ 62 METHODNAME(complex128_t, Container); \ 64 #define DEFINE_FOR_REAL_PTYPE(METHODNAME, Container) \ 65 METHODNAME(bool, Container); \ 66 METHODNAME(char, Container); \ 67 METHODNAME(int8_t, Container); \ 68 METHODNAME(uint8_t, Container); \ 69 METHODNAME(int16_t, Container); \ 70 METHODNAME(uint16_t, Container); \ 71 METHODNAME(int32_t, Container); \ 72 METHODNAME(uint32_t, Container); \ 73 METHODNAME(int64_t, Container); \ 74 METHODNAME(uint64_t, Container); \ 75 METHODNAME(float32_t, Container); \ 76 METHODNAME(float64_t, Container); \ 77 METHODNAME(floatmax_t, Container); 79 #define DEFINE_FOR_NON_INTEGER_PTYPE(METHODNAME, Container) \ 80 METHODNAME(float32_t, Container); \ 81 METHODNAME(float64_t, Container); \ 82 METHODNAME(floatmax_t, Container); \ 83 METHODNAME(complex128_t, Container); 85 #define DEFINE_FOR_NUMERIC_PTYPE(METHODNAME, Container) \ 86 METHODNAME(char, Container); \ 87 METHODNAME(int8_t, Container); \ 88 METHODNAME(uint8_t, Container); \ 89 METHODNAME(int16_t, Container); \ 90 METHODNAME(uint16_t, Container); \ 91 METHODNAME(int32_t, Container); \ 92 METHODNAME(uint32_t, Container); \ 93 METHODNAME(int64_t, Container); \ 94 METHODNAME(uint64_t, Container); \ 95 METHODNAME(float32_t, Container); \ 96 METHODNAME(float64_t, Container); \ 97 METHODNAME(floatmax_t, Container); 100 #define BACKEND_GENERIC_IN_PLACE_ADD(Type, Container) \ 101 virtual void add(Container<Type>& a, Container<Type>& b, Type alpha, \ 102 Type beta, Container<Type>& result) const \ 104 add_impl(a, b, alpha, beta, result); \ 108 #undef BACKEND_GENERIC_IN_PLACE_ADD 111 #define BACKEND_GENERIC_CHOLESKY_FACTOR(Type, Container) \ 112 virtual Container<Type> cholesky_factor(const Container<Type>& A, \ 113 const bool lower) const \ 115 return cholesky_factor_impl(A, lower); \ 118 #undef BACKEND_GENERIC_CHOLESKY_FACTOR 121 #define BACKEND_GENERIC_CHOLESKY_SOLVER(Type, Container) \ 122 virtual SGVector<Type> cholesky_solver(const Container<Type>& L, \ 123 const SGVector<Type>& b, const bool lower) const \ 125 return cholesky_solver_impl(L, b, lower); \ 128 #undef BACKEND_GENERIC_CHOLESKY_SOLVER 131 #define BACKEND_GENERIC_DOT(Type, Container) \ 132 virtual Type dot(const Container<Type>& a, const Container<Type>& b) const \ 134 return dot_impl(a, b); \ 137 #undef BACKEND_GENERIC_DOT 140 #define BACKEND_GENERIC_IN_PLACE_ELEMENT_PROD(Type, Container) \ 141 virtual void element_prod(Container<Type>& a, Container<Type>& b,\ 142 Container<Type>& result) const \ 144 element_prod_impl(a, b, result); \ 147 #undef BACKEND_GENERIC_IN_PLACE_ELEMENT_PROD 150 #define BACKEND_GENERIC_IN_PLACE_BLOCK_ELEMENT_PROD(Type, Container) \ 151 virtual void element_prod(linalg::Block<Container<Type>>& a, \ 152 linalg::Block<Container<Type>>& b, Container<Type>& result) const \ 154 element_prod_impl(a, b, result); \ 157 #undef BACKEND_GENERIC_IN_PLACE_BLOCK_ELEMENT_PROD 160 #define BACKEND_GENERIC_LOGISTIC(Type, Container) \ 161 virtual void logistic(Container<Type>& a, Container<Type>& result) const \ 163 logistic_impl(a, result); \ 166 #undef BACKEND_GENERIC_LOGISTIC 169 #define BACKEND_GENERIC_IN_PLACE_MATRIX_PROD(Type, Container) \ 170 virtual void matrix_prod(SGMatrix<Type>& a, Container<Type>& b,\ 171 Container<Type>& result, bool transpose_A, bool transpose_B) const \ 173 matrix_prod_impl(a, b, result, transpose_A, transpose_B); \ 177 #undef BACKEND_GENERIC_IN_PLACE_MATRIX_PROD 180 #define BACKEND_GENERIC_MAX(Type, Container) \ 181 virtual Type max(const Container<Type>& a) const \ 183 return max_impl(a); \ 187 #undef BACKEND_GENERIC_MAX 190 #define BACKEND_GENERIC_REAL_MEAN(Type, Container) \ 191 virtual float64_t mean(const Container<Type>& a) const \ 193 return mean_impl(a); \ 197 #undef BACKEND_GENERIC_REAL_MEAN 200 #define BACKEND_GENERIC_COMPLEX_MEAN(Container) \ 201 virtual complex128_t mean(const Container<complex128_t>& a) const \ 203 return mean_impl(a); \ 207 #undef BACKEND_GENERIC_COMPLEX_MEAN 210 #define BACKEND_GENERIC_RANGE_FILL(Type, Container) \ 211 virtual void range_fill(Container<Type>& a, const Type start) const \ 213 range_fill_impl(a, start); \ 217 #undef BACKEND_GENERIC_RANGE_FILL 220 #define BACKEND_GENERIC_IN_PLACE_SCALE(Type, Container) \ 221 virtual void scale(Container<Type>& a, Type alpha, Container<Type>& result) const \ 223 scale_impl(a, alpha, result); \ 227 #undef BACKEND_GENERIC_IN_PLACE_SCALE 230 #define BACKEND_GENERIC_SET_CONST(Type, Container) \ 231 virtual void set_const(Container<Type>& a, const Type value) const \ 233 set_const_impl(a, value); \ 237 #undef BACKEND_GENERIC_SET_CONST 240 #define BACKEND_GENERIC_SUM(Type, Container) \ 241 virtual Type sum(const Container<Type>& a, bool no_diag) const \ 243 return sum_impl(a, no_diag); \ 247 #undef BACKEND_GENERIC_SUM 250 #define BACKEND_GENERIC_BLOCK_SUM(Type, Container) \ 251 virtual Type sum(const linalg::Block<Container<Type>>& a, bool no_diag) const \ 253 return sum_impl(a, no_diag); \ 256 #undef BACKEND_GENERIC_BLOCK_SUM 259 #define BACKEND_GENERIC_SYMMETRIC_SUM(Type, Container) \ 260 virtual Type sum_symmetric(const Container<Type>& a, bool no_diag) const \ 262 return sum_symmetric_impl(a, no_diag); \ 265 #undef BACKEND_GENERIC_SYMMETRIC_SUM 268 #define BACKEND_GENERIC_SYMMETRIC_BLOCK_SUM(Type, Container) \ 269 virtual Type sum_symmetric(const linalg::Block<Container<Type>>& a, bool no_diag) const \ 271 return sum_symmetric_impl(a, no_diag); \ 274 #undef BACKEND_GENERIC_SYMMETRIC_BLOCK_SUM 277 #define BACKEND_GENERIC_COLWISE_SUM(Type, Container) \ 278 virtual SGVector<Type> colwise_sum(const Container<Type>& a, bool no_diag) const \ 280 return colwise_sum_impl(a, no_diag); \ 283 #undef BACKEND_GENERIC_COLWISE_SUM 286 #define BACKEND_GENERIC_BLOCK_COLWISE_SUM(Type, Container) \ 287 virtual SGVector<Type> colwise_sum(const linalg::Block<Container<Type>>& a, bool no_diag) const \ 289 return colwise_sum_impl(a, no_diag); \ 292 #undef BACKEND_GENERIC_BLOCK_COLWISE_SUM 295 #define BACKEND_GENERIC_ROWWISE_SUM(Type, Container) \ 296 virtual SGVector<Type> rowwise_sum(const Container<Type>& a, bool no_diag) const \ 298 return rowwise_sum_impl(a, no_diag); \ 301 #undef BACKEND_GENERIC_ROWWISE_SUM 304 #define BACKEND_GENERIC_BLOCK_ROWWISE_SUM(Type, Container) \ 305 virtual SGVector<Type> rowwise_sum(const linalg::Block<Container<Type>>& a, bool no_diag) const \ 307 return rowwise_sum_impl(a, no_diag); \ 310 #undef BACKEND_GENERIC_BLOCK_ROWWISE_SUM 312 #undef DEFINE_FOR_ALL_PTYPE 316 template <
typename T>
317 void add_impl(SGVector<T>& a, SGVector<T>& b, T alpha, T beta, SGVector<T>& result)
const 323 result_eig = alpha * a_eig + beta * b_eig;
327 template <
typename T>
328 void add_impl(SGMatrix<T>& a, SGMatrix<T>& b, T alpha, T beta, SGMatrix<T>& result)
const 334 result_eig = alpha * a_eig + beta * b_eig;
338 template <
typename T>
339 SGMatrix<T> cholesky_factor_impl(
const SGMatrix<T>& A,
const bool lower)
const 341 SGMatrix<T> c(A.num_rows, A.num_cols);
342 set_const_impl<T>(c, 0);
346 Eigen::LLT<Eigen::Matrix<T,Eigen::Dynamic,Eigen::Dynamic> > llt(A_eig);
350 c_eig = llt.matrixU();
352 c_eig = llt.matrixL();
360 REQUIRE(llt.info()!=Eigen::NumericalIssue,
"Matrix is not Hermitian positive definite!\n");
366 template <
typename T>
367 SGVector<T> cholesky_solver_impl(
const SGMatrix<T>& L,
const SGVector<T>& b,
368 const bool lower)
const 370 SGVector<T> x(b.size());
371 set_const_impl<T>(x, 0);
378 Eigen::TriangularView<Eigen::Map<typename SGMatrix<T>::EigenMatrixXt,
381 x_eig = (tlv.transpose()).solve(tlv.solve(b_eig));
385 Eigen::TriangularView<Eigen::Map<typename SGMatrix<T>::EigenMatrixXt,
387 x_eig = (tlv.transpose()).solve(tlv.solve(b_eig));
394 template <
typename T>
395 T dot_impl(
const SGVector<T>& a,
const SGVector<T>& b)
const 401 template <
typename T>
402 void element_prod_impl(SGMatrix<T>& a, SGMatrix<T>& b, SGMatrix<T>& result)
const 408 result_eig = a_eig.array() * b_eig.array();
412 template <
typename T>
413 void logistic_impl(SGMatrix<T>& a, SGMatrix<T>& result)
const 418 result_eig = (T)1 / (1 + ((-1 * a_eig).array()).exp());
422 template <
typename T>
423 void element_prod_impl(linalg::Block<SGMatrix<T>>& a,
424 linalg::Block<SGMatrix<T>>& b, SGMatrix<T>& result)
const 430 Eigen::Block<typename SGMatrix<T>::EigenMatrixXtMap> a_block =
431 a_eig.block(a.m_row_begin, a.m_col_begin, a.m_row_size, a.m_col_size);
432 Eigen::Block<typename SGMatrix<T>::EigenMatrixXtMap> b_block =
433 b_eig.block(b.m_row_begin, b.m_col_begin, b.m_row_size, b.m_col_size);
435 result_eig = a_block.array() * b_block.array();
439 template <
typename T>
440 void matrix_prod_impl(SGMatrix<T>& a, SGVector<T>& b, SGVector<T>& result,
441 bool transpose,
bool transpose_B=
false)
const 448 result_eig = a_eig.transpose() * b_eig;
450 result_eig = a_eig * b_eig;
454 template <
typename T>
455 void matrix_prod_impl(SGMatrix<T>& a, SGMatrix<T>& b, SGMatrix<T>& result,
456 bool transpose_A,
bool transpose_B)
const 462 if (transpose_A && transpose_B)
463 result_eig = a_eig.transpose() * b_eig.transpose();
465 else if (transpose_A)
466 result_eig = a_eig.transpose() * b_eig;
468 else if (transpose_B)
469 result_eig = a_eig * b_eig.transpose();
472 result_eig = a_eig * b_eig;
476 template <
typename T>
477 T max_impl(
const SGVector<T>& vec)
const 483 template <
typename T>
484 T max_impl(
const SGMatrix<T>& mat)
const 490 template <
typename T,
template <
typename>
class Container>
491 typename std::enable_if<!std::is_same<T, complex128_t>::value,
float64_t>::type
492 mean_impl(
const Container<T>& a)
const 494 return sum_impl(a)/(
float64_t(a.size()));
498 template<
template <
typename>
class Container>
499 complex128_t mean_impl(
const Container<complex128_t>& a)
const 505 template <
typename T,
template <
typename>
class Container>
506 void range_fill_impl(Container<T>& a,
const T start)
const 508 for (
index_t i = 0; i < a.size(); ++i)
513 template <
typename T>
514 void scale_impl(SGVector<T>& a, T alpha, SGVector<T>& result)
const 519 result_eig = alpha * a_eig;
523 template <
typename T>
524 void scale_impl(SGMatrix<T>& a, T alpha, SGMatrix<T>& result)
const 529 result_eig = alpha * a_eig;
533 template <
typename T,
template <
typename>
class Container>
534 void set_const_impl(Container<T>& a, T value)
const 536 for (
index_t i = 0; i < a.size(); ++i)
541 template <
typename T>
542 T sum_impl(
const SGVector<T>& vec,
bool no_diag=
false)
const 548 template <
typename T>
549 T sum_impl(
const SGMatrix<T>& mat,
bool no_diag=
false)
const 554 result -= m.diagonal().sum();
560 template <
typename T>
561 T sum_impl(
const linalg::Block<SGMatrix<T>>& mat,
bool no_diag=
false)
const 564 Eigen::Block<typename SGMatrix<T>::EigenMatrixXtMap> m_block = m.block(
565 mat.m_row_begin, mat.m_col_begin, mat.m_row_size, mat.m_col_size);
567 T result = m_block.sum();
569 result -= m_block.diagonal().sum();
575 template <
typename T>
576 T sum_symmetric_impl(
const SGMatrix<T>& mat,
bool no_diag=
false)
const 583 m.template triangularView<Eigen::StrictlyUpper>();
584 T result = m_upper.sum();
588 result += m.diagonal().sum();
593 template <
typename T>
594 T sum_symmetric_impl(
const linalg::Block<SGMatrix<T>>& mat,
bool no_diag=
false)
const 597 Eigen::Block<typename SGMatrix<T>::EigenMatrixXtMap> m_block = m.block(
598 mat.m_row_begin, mat.m_col_begin, mat.m_row_size, mat.m_col_size);
603 m_block.template triangularView<Eigen::StrictlyUpper>();
604 T result = m_upper.sum();
608 result += m_block.diagonal().sum();
613 template <
typename T>
614 SGVector<T> colwise_sum_impl(
const SGMatrix<T>& mat,
bool no_diag)
const 616 SGVector<T> result(mat.num_cols);
621 result_eig = mat_eig.colwise().sum();
626 index_t len_major_diag = mat_eig.rows() < mat_eig.cols()
627 ? mat_eig.rows() : mat_eig.cols();
628 for (
index_t i = 0; i < len_major_diag; ++i)
629 result_eig[i] -= mat_eig(i,i);
636 template <
typename T>
637 SGVector<T> colwise_sum_impl(
const linalg::Block<SGMatrix<T>>& mat,
bool no_diag)
const 639 SGVector<T> result(mat.m_col_size);
642 Eigen::Block<typename SGMatrix<T>::EigenMatrixXtMap> m_block = m.block(
643 mat.m_row_begin, mat.m_col_begin, mat.m_row_size, mat.m_col_size);
646 result_eig = m_block.colwise().sum();
651 index_t len_major_diag = m_block.rows() < m_block.cols()
652 ? m_block.rows() : m_block.cols();
653 for (
index_t i = 0; i < len_major_diag; ++i)
654 result_eig[i] -= m_block(i,i);
661 template <
typename T>
662 SGVector<T> rowwise_sum_impl(
const SGMatrix<T>& mat,
bool no_diag)
const 664 SGVector<T> result(mat.num_rows);
669 result_eig = mat_eig.rowwise().sum();
674 index_t len_major_diag = mat_eig.rows() < mat_eig.cols()
675 ? mat_eig.rows() : mat_eig.cols();
676 for (
index_t i = 0; i < len_major_diag; ++i)
677 result_eig[i] -= mat_eig(i,i);
684 template <
typename T>
685 SGVector<T> rowwise_sum_impl(
const linalg::Block<SGMatrix<T>>& mat,
bool no_diag)
const 687 SGVector<T> result(mat.m_row_size);
690 Eigen::Block<typename SGMatrix<T>::EigenMatrixXtMap> m_block = m.block(
691 mat.m_row_begin, mat.m_col_begin, mat.m_row_size, mat.m_col_size);
694 result_eig = m_block.rowwise().sum();
699 index_t len_major_diag = m_block.rows() < m_block.cols()
700 ? m_block.rows() : m_block.cols();
701 for (
index_t i = 0; i < len_major_diag; ++i)
702 result_eig[i] -= m_block(i,i);
708 #undef DEFINE_FOR_ALL_PTYPE 709 #undef DEFINE_FOR_REAL_PTYPE 710 #undef DEFINE_FOR_NON_INTEGER_PTYPE 711 #undef DEFINE_FOR_NUMERIC_PTYPE 716 #endif //LINALG_BACKEND_EIGEN_H__ #define BACKEND_GENERIC_BLOCK_COLWISE_SUM(Type, Container)
#define BACKEND_GENERIC_CHOLESKY_FACTOR(Type, Container)
std::complex< float64_t > complex128_t
#define BACKEND_GENERIC_COMPLEX_MEAN(Container)
#define DEFINE_FOR_ALL_PTYPE(METHODNAME, Container)
Eigen::Map< EigenMatrixXt, 0, Eigen::Stride< 0, 0 > > EigenMatrixXtMap
#define BACKEND_GENERIC_ROWWISE_SUM(Type, Container)
#define BACKEND_GENERIC_IN_PLACE_BLOCK_ELEMENT_PROD(Type, Container)
Eigen::Matrix< T,-1,-1, 0,-1,-1 > EigenMatrixXt
#define BACKEND_GENERIC_DOT(Type, Container)
Eigen::Map< EigenVectorXt, 0, Eigen::Stride< 0, 0 > > EigenVectorXtMap
#define BACKEND_GENERIC_BLOCK_SUM(Type, Container)
#define BACKEND_GENERIC_REAL_MEAN(Type, Container)
#define BACKEND_GENERIC_COLWISE_SUM(Type, Container)
#define BACKEND_GENERIC_BLOCK_ROWWISE_SUM(Type, Container)
#define BACKEND_GENERIC_MAX(Type, Container)
#define BACKEND_GENERIC_IN_PLACE_MATRIX_PROD(Type, Container)
#define BACKEND_GENERIC_SUM(Type, Container)
#define BACKEND_GENERIC_CHOLESKY_SOLVER(Type, Container)
#define BACKEND_GENERIC_IN_PLACE_ELEMENT_PROD(Type, Container)
#define DEFINE_FOR_REAL_PTYPE(METHODNAME, Container)
#define BACKEND_GENERIC_LOGISTIC(Type, Container)
#define BACKEND_GENERIC_IN_PLACE_SCALE(Type, Container)
#define BACKEND_GENERIC_SYMMETRIC_SUM(Type, Container)
all of classes and functions are contained in the shogun namespace
Base interface of generic linalg methods and generic memory transfer methods.
#define BACKEND_GENERIC_IN_PLACE_ADD(Type, Container)
Linalg methods with Eigen3 backend.
#define DEFINE_FOR_NON_INTEGER_PTYPE(METHODNAME, Container)
#define BACKEND_GENERIC_SYMMETRIC_BLOCK_SUM(Type, Container)
#define BACKEND_GENERIC_SET_CONST(Type, Container)
#define DEFINE_FOR_NUMERIC_PTYPE(METHODNAME, Container)
#define BACKEND_GENERIC_RANGE_FILL(Type, Container)