1 #ifndef VIENNACL_TRAITS_SIZE_HPP_ 
    2 #define VIENNACL_TRAITS_SIZE_HPP_ 
   32 #ifdef VIENNACL_WITH_UBLAS 
   33 #include <boost/numeric/ublas/matrix_sparse.hpp> 
   34 #include <boost/numeric/ublas/matrix.hpp> 
   37 #ifdef VIENNACL_WITH_ARMADILLO 
   41 #ifdef VIENNACL_WITH_EIGEN 
   43 #include <Eigen/Sparse> 
   46 #ifdef VIENNACL_WITH_MTL4 
   47 #include <boost/numeric/mtl/mtl.hpp> 
   62 template<
typename MatrixType>
 
   65   matrix.resize(rows, cols);
 
   69 template<
typename VectorType>
 
   76 #ifdef VIENNACL_WITH_UBLAS 
   78 template<
typename ScalarType>
 
   79 void resize(boost::numeric::ublas::compressed_matrix<ScalarType> & 
matrix,
 
   83   matrix.resize(rows, cols, 
false); 
 
   88 #ifdef VIENNACL_WITH_MTL4 
   89 template<
typename ScalarType>
 
   90 void resize(mtl::compressed2D<ScalarType> & matrix,
 
   94   matrix.change_dim(rows, cols);
 
   97 template<
typename ScalarType>
 
   98 void resize(mtl::dense_vector<ScalarType> & vec,
 
  101   vec.change_dim(new_size);
 
  105 #ifdef VIENNACL_WITH_ARMADILLO 
  106 template<
typename NumericT>
 
  107 inline void resize(arma::Mat<NumericT> & A,
 
  111   A.resize(new_rows, new_cols);
 
  114 template<
typename NumericT>
 
  115 inline void resize(arma::SpMat<NumericT> & A,
 
  119   A.set_size(new_rows, new_cols);
 
  123 #ifdef VIENNACL_WITH_EIGEN 
  124 template<
typename NumericT, 
int Options>
 
  125 inline void resize(Eigen::Matrix<NumericT, Eigen::Dynamic, Eigen::Dynamic, Options> & m,
 
  129   m.resize(new_rows, new_cols);
 
  132 template<
typename T, 
int options>
 
  133 inline void resize(Eigen::SparseMatrix<T, options> & m,
 
  137   m.resize(new_rows, new_cols);
 
  140 inline void resize(Eigen::VectorXf & v,
 
  146 inline void resize(Eigen::VectorXd & v,
 
  161 template<
typename MatrixType>
 
  163 size1(MatrixType 
const & mat) { 
return mat.size1(); }
 
  166 template<
typename RowType>
 
  168 size1(std::vector< RowType > 
const & mat) { 
return mat.size(); }
 
  170 #ifdef VIENNACL_WITH_ARMADILLO 
  171 template<
typename NumericT>
 
  172 inline vcl_size_t size1(arma::Mat<NumericT> 
const & A) { 
return A.n_rows; }
 
  173 template<
typename NumericT>
 
  174 inline vcl_size_t size1(arma::SpMat<NumericT> 
const & A) { 
return A.n_rows; }
 
  177 #ifdef VIENNACL_WITH_EIGEN 
  178 template<
typename NumericT, 
int Options>
 
  179 vcl_size_t size1(Eigen::Matrix<NumericT, Eigen::Dynamic, Eigen::Dynamic, Options> 
const & m) { 
return static_cast<vcl_size_t>(m.rows()); }
 
  180 template<
typename NumericT, 
int Options>
 
  181 vcl_size_t size1(Eigen::Map< Eigen::Matrix<NumericT, Eigen::Dynamic, Eigen::Dynamic, Options> > 
const & m) { 
return static_cast<vcl_size_t>(m.rows()); }
 
  182 template<
typename T, 
int options>
 
  186 #ifdef VIENNACL_WITH_MTL4 
  187 template<
typename NumericT, 
typename T>
 
  189 template<
typename NumericT>
 
  199 template<
typename MatrixType>
 
  200 typename result_of::size_type<MatrixType>::type
 
  201 size2(MatrixType 
const & mat) { 
return mat.size2(); }
 
  204 template<
typename RowType>
 
  206 size2(std::vector< RowType > 
const & mat) { 
return mat[0].size(); }
 
  208 #ifdef VIENNACL_WITH_ARMADILLO 
  209 template<
typename NumericT>
 
  210 inline vcl_size_t size2(arma::Mat<NumericT> 
const & A) { 
return A.n_cols; }
 
  211 template<
typename NumericT>
 
  212 inline vcl_size_t size2(arma::SpMat<NumericT> 
const & A) { 
return A.n_cols; }
 
  215 #ifdef VIENNACL_WITH_EIGEN 
  216 template<
typename NumericT, 
int Options>
 
  217 inline vcl_size_t size2(Eigen::Matrix<NumericT, Eigen::Dynamic, Eigen::Dynamic, Options> 
const & m) { 
return m.cols(); }
 
  218 template<
typename NumericT, 
int Options>
 
  219 inline vcl_size_t size2(Eigen::Map< Eigen::Matrix<NumericT, Eigen::Dynamic, Eigen::Dynamic, Options> > 
const & m) { 
return m.cols(); }
 
  220 template<
typename T, 
int options>
 
  221 inline vcl_size_t size2(Eigen::SparseMatrix<T, options> & m) { 
return m.cols(); }
 
  224 #ifdef VIENNACL_WITH_MTL4 
  225 template<
typename NumericT, 
typename T>
 
  227 template<
typename NumericT>
 
  238 template<
typename VectorType>
 
  245 template<
typename SparseMatrixType, 
typename VectorType>
 
  251 template<
typename T, 
unsigned int A, 
typename VectorType>
 
  252 vcl_size_t size(vector_expression<
const circulant_matrix<T, A>, 
const VectorType, op_prod> 
const & proxy) { 
return proxy.
lhs().size1();  }
 
  254 template<
typename T, 
unsigned int A, 
typename VectorType>
 
  255 vcl_size_t size(vector_expression<
const hankel_matrix<T, A>, 
const VectorType, op_prod> 
const & proxy) { 
return proxy.
lhs().size1();  }
 
  257 template<
typename T, 
unsigned int A, 
typename VectorType>
 
  258 vcl_size_t size(vector_expression<
const toeplitz_matrix<T, A>, 
const VectorType, op_prod> 
const & proxy) { 
return proxy.
lhs().size1();  }
 
  260 template<
typename T, 
unsigned int A, 
typename VectorType>
 
  261 vcl_size_t size(vector_expression<
const vandermonde_matrix<T, A>, 
const VectorType, op_prod> 
const & proxy) { 
return proxy.
lhs().size1();  }
 
  263 template<
typename NumericT>
 
  264 vcl_size_t size(vector_expression<
const matrix_base<NumericT>, 
const vector_base<NumericT>, op_prod> 
const & proxy)  
 
  266   return proxy.
lhs().size1();
 
  269 template<
typename NumericT, 
typename LhsT, 
typename RhsT, 
typename OpT>
 
  270 vcl_size_t size(vector_expression<
const matrix_base<NumericT>, 
const vector_expression<LhsT, RhsT, OpT>, op_prod> 
const & proxy)  
 
  272   return proxy.
lhs().size1();
 
  275 template<
typename NumericT>
 
  276 vcl_size_t size(vector_expression<
const matrix_expression<
const matrix_base<NumericT>, 
const matrix_base<NumericT>, op_trans>,
 
  277                 const vector_base<NumericT>,
 
  278                 op_prod> 
const & proxy)  
 
  280   return proxy.
lhs().lhs().size2();
 
  284 #ifdef VIENNACL_WITH_MTL4 
  285 template<
typename ScalarType>
 
  286 vcl_size_t size(mtl::dense_vector<ScalarType> 
const & vec) { 
return vec.used_memory(); }
 
  289 #ifdef VIENNACL_WITH_ARMADILLO 
  290 template<
typename NumericT>
 
  291 inline vcl_size_t size(arma::Mat<NumericT> 
const & A) { 
return A.n_elem; }
 
  294 #ifdef VIENNACL_WITH_EIGEN 
  295 inline vcl_size_t size(Eigen::VectorXf 
const & v) { 
return v.rows(); }
 
  296 inline vcl_size_t size(Eigen::VectorXd 
const & v) { 
return v.rows(); }
 
  299 template<
typename LHS, 
typename RHS, 
typename OP>
 
  302   return size(proxy.lhs());
 
  305 template<
typename LHS, 
typename RHS>
 
  306 vcl_size_t size(vector_expression<LHS, 
const vector_tuple<RHS>, op_inner_prod> 
const & proxy)
 
  308   return proxy.rhs().const_size();
 
  311 template<
typename LhsT, 
typename RhsT, 
typename OpT, 
typename VectorT>
 
  312 vcl_size_t size(vector_expression<
const matrix_expression<const LhsT, const RhsT, OpT>,
 
  314                                   op_prod> 
const & proxy)
 
  316   return size1(proxy.lhs());
 
  319 template<
typename LhsT, 
typename RhsT, 
typename OpT, 
typename NumericT>
 
  320 vcl_size_t size(vector_expression<
const matrix_expression<const LhsT, const RhsT, OpT>,
 
  321                                   const vector_base<NumericT>,
 
  322                                   op_prod> 
const & proxy)
 
  324   return size1(proxy.lhs());
 
  327 template<
typename LhsT1, 
typename RhsT1, 
typename OpT1,
 
  328          typename LhsT2, 
typename RhsT2, 
typename OpT2>
 
  329 vcl_size_t size(vector_expression<
const matrix_expression<const LhsT1, const RhsT1, OpT1>,
 
  330                                   const vector_expression<const LhsT2, const RhsT2, OpT2>,
 
  331                                   op_prod> 
const & proxy)
 
  333   return size1(proxy.lhs());
 
  336 template<
typename NumericT>
 
  338                                   const matrix_base<NumericT>,
 
  339                                   op_row_sum> 
const & proxy)
 
  341   return size1(proxy.lhs());
 
  344 template<
typename LhsT, 
typename RhsT, 
typename OpT>
 
  345 vcl_size_t size(vector_expression<
const matrix_expression<const LhsT, const RhsT, OpT>,
 
  346                                   const matrix_expression<const LhsT, const RhsT, OpT>,
 
  347                                   op_row_sum> 
const & proxy)
 
  349   return size1(proxy.lhs());
 
  352 template<
typename NumericT>
 
  354                                   const matrix_base<NumericT>,
 
  355                                   op_col_sum> 
const & proxy)
 
  357   return size2(proxy.lhs());
 
  360 template<
typename LhsT, 
typename RhsT, 
typename OpT>
 
  361 vcl_size_t size(vector_expression<
const matrix_expression<const LhsT, const RhsT, OpT>,
 
  362                                   const matrix_expression<const LhsT, const RhsT, OpT>,
 
  363                                   op_col_sum> 
const & proxy)
 
  365   return size2(proxy.lhs());
 
  374 template<
typename NumericT>
 
  385 template<
typename NumericT>
 
  393 template<
typename NumericT>
 
  397 template<
typename NumericT>
 
  405 template<
typename NumericT>
 
  413 template<
typename LHS>
 
  417   int A_size1 = 
static_cast<int>(
size1(proxy.
lhs()));
 
  418   int A_size2 = 
static_cast<int>(
size2(proxy.
lhs()));
 
  420   int row_depth = 
std::min(A_size1, A_size1 + k);
 
  421   int col_depth = 
std::min(A_size2, A_size2 - k);
 
  426 template<
typename LHS>
 
  432 template<
typename LHS>
 
vcl_size_t nld(matrix_base< NumericT > const &mat)
vcl_size_t internal_size1(matrix_base< NumericT > const &mat)
Helper routine for obtaining the internal number of entries per row of a ViennaCL matrix...
vcl_size_t size1(MatrixType const &mat)
Generic routine for obtaining the number of rows of a matrix (ViennaCL, uBLAS, etc.) 
vcl_size_t internal_size2(matrix_base< NumericT > const &mat)
Helper routine for obtaining the internal number of entries per column of a ViennaCL matrix...
size_type stride2() const 
Returns the number of columns. 
This file provides the forward declarations for the main types used within ViennaCL. 
vcl_size_t internal_size(vector_base< NumericT > const &vec)
Helper routine for obtaining the buffer length of a ViennaCL vector. 
lhs_reference_type lhs() const 
Get left hand side operand. 
An expression template class that represents a binary operation that yields a vector. 
result_of::size_type< MatrixType >::type size2(MatrixType const &mat)
Generic routine for obtaining the number of columns of a matrix (ViennaCL, uBLAS, etc...
vcl_size_t size(VectorType const &vec)
Generic routine for obtaining the size of a vector (ViennaCL, uBLAS, etc.) 
vcl_size_t ld(matrix_base< NumericT > const &mat)
Helper routine for obtaining the internal number of entries per row of a ViennaCL matrix...
void resize(MatrixType &matrix, vcl_size_t rows, vcl_size_t cols)
Generic resize routine for resizing a matrix (ViennaCL, uBLAS, etc.) to a new size/dimension. 
rhs_reference_type rhs() const 
Get right hand side operand. 
size_type stride1() const 
Returns the number of rows. 
All the predicates used within ViennaCL. Checks for expressions to be vectors, etc. 
size_type internal_size2() const 
Returns the internal number of columns. Usually required for launching OpenCL kernels only...
size_type internal_size1() const 
Returns the internal number of rows. Usually required for launching OpenCL kernels only...
T min(const T &lhs, const T &rhs)
Minimum. 
size_type internal_size() const 
Returns the internal length of the vector, which is given by size() plus the extra memory due to padd...
A collection of compile time type deductions.