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_EIGEN
39 #include <Eigen/Sparse>
42 #ifdef VIENNACL_WITH_MTL4
43 #include <boost/numeric/mtl/mtl.hpp>
58 template<
typename MatrixType>
61 matrix.resize(rows, cols);
65 template<
typename VectorType>
72 #ifdef VIENNACL_WITH_UBLAS
74 template<
typename ScalarType>
75 void resize(boost::numeric::ublas::compressed_matrix<ScalarType> &
matrix,
79 matrix.resize(rows, cols,
false);
84 #ifdef VIENNACL_WITH_MTL4
85 template<
typename ScalarType>
86 void resize(mtl::compressed2D<ScalarType> & matrix,
90 matrix.change_dim(rows, cols);
93 template<
typename ScalarType>
94 void resize(mtl::dense_vector<ScalarType> & vec,
97 vec.change_dim(new_size);
101 #ifdef VIENNACL_WITH_EIGEN
102 inline void resize(Eigen::MatrixXf & m,
106 m.resize(new_rows, new_cols);
109 inline void resize(Eigen::MatrixXd & m,
113 m.resize(new_rows, new_cols);
116 template<
typename T,
int options>
117 inline void resize(Eigen::SparseMatrix<T, options> & m,
121 m.resize(new_rows, new_cols);
124 inline void resize(Eigen::VectorXf & v,
130 inline void resize(Eigen::VectorXd & v,
143 template<
typename VectorType>
150 template<
typename SparseMatrixType,
typename VectorType>
155 return proxy.
lhs().size1();
158 template<
typename T,
unsigned int A,
typename VectorType>
159 vcl_size_t size(vector_expression<
const circulant_matrix<T, A>,
const VectorType, op_prod>
const & proxy) {
return proxy.
lhs().size1(); }
161 template<
typename T,
unsigned int A,
typename VectorType>
162 vcl_size_t size(vector_expression<
const hankel_matrix<T, A>,
const VectorType, op_prod>
const & proxy) {
return proxy.
lhs().size1(); }
164 template<
typename T,
unsigned int A,
typename VectorType>
165 vcl_size_t size(vector_expression<
const toeplitz_matrix<T, A>,
const VectorType, op_prod>
const & proxy) {
return proxy.
lhs().size1(); }
167 template<
typename T,
unsigned int A,
typename VectorType>
168 vcl_size_t size(vector_expression<
const vandermonde_matrix<T, A>,
const VectorType, op_prod>
const & proxy) {
return proxy.
lhs().size1(); }
170 template<
typename NumericT>
171 vcl_size_t size(vector_expression<
const matrix_base<NumericT>,
const vector_base<NumericT>, op_prod>
const & proxy)
173 return proxy.
lhs().size1();
176 template<
typename NumericT>
177 vcl_size_t size(vector_expression<
const matrix_expression<
const matrix_base<NumericT>,
const matrix_base<NumericT>, op_trans>,
178 const vector_base<NumericT>,
179 op_prod>
const & proxy)
181 return proxy.
lhs().lhs().size2();
185 #ifdef VIENNACL_WITH_MTL4
186 template<
typename ScalarType>
187 vcl_size_t size(mtl::dense_vector<ScalarType>
const & vec) {
return vec.used_memory(); }
190 #ifdef VIENNACL_WITH_EIGEN
191 inline vcl_size_t size(Eigen::VectorXf
const & v) {
return v.rows(); }
192 inline vcl_size_t size(Eigen::VectorXd
const & v) {
return v.rows(); }
195 template<
typename LHS,
typename RHS,
typename OP>
198 return size(proxy.lhs());
201 template<
typename LHS,
typename RHS>
202 vcl_size_t size(vector_expression<LHS,
const vector_tuple<RHS>, op_inner_prod>
const & proxy)
204 return proxy.rhs().const_size();
214 template<
typename MatrixType>
216 size1(MatrixType
const & mat) {
return mat.size1(); }
219 template<
typename RowType>
221 size1(std::vector< RowType >
const & mat) {
return mat.size(); }
223 #ifdef VIENNACL_WITH_EIGEN
226 template<
typename T,
int options>
230 #ifdef VIENNACL_WITH_MTL4
231 template<
typename NumericT,
typename T>
233 template<
typename NumericT>
243 template<
typename MatrixType>
244 typename result_of::size_type<MatrixType>::type
245 size2(MatrixType
const & mat) {
return mat.size2(); }
248 #ifdef VIENNACL_WITH_EIGEN
249 inline vcl_size_t size2(Eigen::MatrixXf
const & m) {
return m.cols(); }
250 inline vcl_size_t size2(Eigen::MatrixXd
const & m) {
return m.cols(); }
251 template<
typename T,
int options>
252 inline vcl_size_t size2(Eigen::SparseMatrix<T, options> & m) {
return m.cols(); }
255 #ifdef VIENNACL_WITH_MTL4
256 template<
typename NumericT,
typename T>
258 template<
typename NumericT>
267 template<
typename NumericT>
278 template<
typename NumericT>
286 template<
typename NumericT>
290 template<
typename NumericT>
298 template<
typename NumericT>
306 template<
typename LHS>
310 int A_size1 =
static_cast<int>(
size1(proxy.
lhs()));
311 int A_size2 =
static_cast<int>(
size2(proxy.
lhs()));
313 int row_depth =
std::min(A_size1, A_size1 + k);
314 int col_depth =
std::min(A_size2, A_size2 - k);
319 template<
typename LHS>
325 template<
typename LHS>
Simple enable-if variant that uses the SFINAE pattern.
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.