ViennaCL - The Vienna Computing Library  1.6.0
Free open-source GPU-accelerated linear algebra and solver library.
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
size.hpp
Go to the documentation of this file.
1 #ifndef VIENNACL_TRAITS_SIZE_HPP_
2 #define VIENNACL_TRAITS_SIZE_HPP_
3 
4 /* =========================================================================
5  Copyright (c) 2010-2014, Institute for Microelectronics,
6  Institute for Analysis and Scientific Computing,
7  TU Wien.
8  Portions of this software are copyright by UChicago Argonne, LLC.
9 
10  -----------------
11  ViennaCL - The Vienna Computing Library
12  -----------------
13 
14  Project Head: Karl Rupp rupp@iue.tuwien.ac.at
15 
16  (A list of authors and contributors can be found in the PDF manual)
17 
18  License: MIT (X11), see file LICENSE in the base directory
19 ============================================================================= */
20 
25 #include <string>
26 #include <fstream>
27 #include <sstream>
28 #include "viennacl/forwards.h"
31 
32 #ifdef VIENNACL_WITH_UBLAS
33 #include <boost/numeric/ublas/matrix_sparse.hpp>
34 #include <boost/numeric/ublas/matrix.hpp>
35 #endif
36 
37 #ifdef VIENNACL_WITH_EIGEN
38 #include <Eigen/Core>
39 #include <Eigen/Sparse>
40 #endif
41 
42 #ifdef VIENNACL_WITH_MTL4
43 #include <boost/numeric/mtl/mtl.hpp>
44 #endif
45 
46 #include <vector>
47 #include <map>
48 
49 namespace viennacl
50 {
51 namespace traits
52 {
53 
54 //
55 // Resize: Change the size of vectors and matrices
56 //
58 template<typename MatrixType>
59 void resize(MatrixType & matrix, vcl_size_t rows, vcl_size_t cols)
60 {
61  matrix.resize(rows, cols);
62 }
63 
65 template<typename VectorType>
66 void resize(VectorType & vec, vcl_size_t new_size)
67 {
68  vec.resize(new_size);
69 }
70 
72 #ifdef VIENNACL_WITH_UBLAS
73 //ublas needs separate treatment:
74 template<typename ScalarType>
75 void resize(boost::numeric::ublas::compressed_matrix<ScalarType> & matrix,
76  vcl_size_t rows,
77  vcl_size_t cols)
78 {
79  matrix.resize(rows, cols, false); //Note: omitting third parameter leads to compile time error (not implemented in ublas <= 1.42)
80 }
81 #endif
82 
83 
84 #ifdef VIENNACL_WITH_MTL4
85 template<typename ScalarType>
86 void resize(mtl::compressed2D<ScalarType> & matrix,
87  vcl_size_t rows,
88  vcl_size_t cols)
89 {
90  matrix.change_dim(rows, cols);
91 }
92 
93 template<typename ScalarType>
94 void resize(mtl::dense_vector<ScalarType> & vec,
95  vcl_size_t new_size)
96 {
97  vec.change_dim(new_size);
98 }
99 #endif
100 
101 #ifdef VIENNACL_WITH_EIGEN
102 inline void resize(Eigen::MatrixXf & m,
103  vcl_size_t new_rows,
104  vcl_size_t new_cols)
105 {
106  m.resize(new_rows, new_cols);
107 }
108 
109 inline void resize(Eigen::MatrixXd & m,
110  vcl_size_t new_rows,
111  vcl_size_t new_cols)
112 {
113  m.resize(new_rows, new_cols);
114 }
115 
116 template<typename T, int options>
117 inline void resize(Eigen::SparseMatrix<T, options> & m,
118  vcl_size_t new_rows,
119  vcl_size_t new_cols)
120 {
121  m.resize(new_rows, new_cols);
122 }
123 
124 inline void resize(Eigen::VectorXf & v,
125  vcl_size_t new_size)
126 {
127  v.resize(new_size);
128 }
129 
130 inline void resize(Eigen::VectorXd & v,
131  vcl_size_t new_size)
132 {
133  v.resize(new_size);
134 }
135 #endif
136 
139 //
140 // size: Returns the length of vectors
141 //
143 template<typename VectorType>
144 vcl_size_t size(VectorType const & vec)
145 {
146  return vec.size();
147 }
148 
150 template<typename SparseMatrixType, typename VectorType>
152 vcl_size_t >::type
154 {
155  return proxy.lhs().size1();
156 }
157 
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(); }
160 
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(); }
163 
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(); }
166 
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(); }
169 
170 template<typename NumericT>
171 vcl_size_t size(vector_expression<const matrix_base<NumericT>, const vector_base<NumericT>, op_prod> const & proxy) //matrix-vector product
172 {
173  return proxy.lhs().size1();
174 }
175 
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) //transposed matrix-vector product
180 {
181  return proxy.lhs().lhs().size2();
182 }
183 
184 
185 #ifdef VIENNACL_WITH_MTL4
186 template<typename ScalarType>
187 vcl_size_t size(mtl::dense_vector<ScalarType> const & vec) { return vec.used_memory(); }
188 #endif
189 
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(); }
193 #endif
194 
195 template<typename LHS, typename RHS, typename OP>
196 vcl_size_t size(vector_expression<LHS, RHS, OP> const & proxy)
197 {
198  return size(proxy.lhs());
199 }
200 
201 template<typename LHS, typename RHS>
202 vcl_size_t size(vector_expression<LHS, const vector_tuple<RHS>, op_inner_prod> const & proxy)
203 {
204  return proxy.rhs().const_size();
205 }
206 
210 //
211 // size1: No. of rows for matrices
212 //
214 template<typename MatrixType>
216 size1(MatrixType const & mat) { return mat.size1(); }
217 
219 template<typename RowType>
221 size1(std::vector< RowType > const & mat) { return mat.size(); }
222 
223 #ifdef VIENNACL_WITH_EIGEN
224 inline vcl_size_t size1(Eigen::MatrixXf const & m) { return static_cast<vcl_size_t>(m.rows()); }
225 inline vcl_size_t size1(Eigen::MatrixXd const & m) { return static_cast<vcl_size_t>(m.rows()); }
226 template<typename T, int options>
227 inline vcl_size_t size1(Eigen::SparseMatrix<T, options> & m) { return static_cast<vcl_size_t>(m.rows()); }
228 #endif
229 
230 #ifdef VIENNACL_WITH_MTL4
231 template<typename NumericT, typename T>
232 vcl_size_t size1(mtl::dense2D<NumericT, T> const & m) { return static_cast<vcl_size_t>(m.num_rows()); }
233 template<typename NumericT>
234 vcl_size_t size1(mtl::compressed2D<NumericT> const & m) { return static_cast<vcl_size_t>(m.num_rows()); }
235 #endif
236 
239 //
240 // size2: No. of columns for matrices
241 //
243 template<typename MatrixType>
244 typename result_of::size_type<MatrixType>::type
245 size2(MatrixType const & mat) { return mat.size2(); }
246 
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(); }
253 #endif
254 
255 #ifdef VIENNACL_WITH_MTL4
256 template<typename NumericT, typename T>
257 vcl_size_t size2(mtl::dense2D<NumericT, T> const & m) { return static_cast<vcl_size_t>(m.num_cols()); }
258 template<typename NumericT>
259 vcl_size_t size2(mtl::compressed2D<NumericT> const & m) { return static_cast<vcl_size_t>(m.num_cols()); }
260 #endif
261 
263 //
264 // internal_size: Returns the internal (padded) length of vectors
265 //
267 template<typename NumericT>
269 {
270  return vec.internal_size();
271 }
272 
273 
274 //
275 // internal_size1: No. of internal (padded) rows for matrices
276 //
278 template<typename NumericT>
280 
281 
282 //
283 // internal_size2: No. of internal (padded) columns for matrices
284 //
286 template<typename NumericT>
288 
290 template<typename NumericT>
292 {
293  if (mat.row_major())
294  return mat.internal_size2();
295  return mat.internal_size1();
296 }
297 
298 template<typename NumericT>
300 {
301  if (mat.row_major())
302  return mat.stride2();
303  return mat.stride1();
304 }
305 
306 template<typename LHS>
308 {
309  int k = proxy.rhs();
310  int A_size1 = static_cast<int>(size1(proxy.lhs()));
311  int A_size2 = static_cast<int>(size2(proxy.lhs()));
312 
313  int row_depth = std::min(A_size1, A_size1 + k);
314  int col_depth = std::min(A_size2, A_size2 - k);
315 
316  return vcl_size_t(std::min(row_depth, col_depth));
317 }
318 
319 template<typename LHS>
321 {
322  return size2(proxy.lhs());
323 }
324 
325 template<typename LHS>
327 {
328  return size1(proxy.lhs());
329 }
330 
331 } //namespace traits
332 } //namespace viennacl
333 
334 
335 #endif
Simple enable-if variant that uses the SFINAE pattern.
Definition: enable_if.hpp:30
vcl_size_t nld(matrix_base< NumericT > const &mat)
Definition: size.hpp:299
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...
Definition: size.hpp:279
vcl_size_t size1(MatrixType const &mat)
Generic routine for obtaining the number of rows of a matrix (ViennaCL, uBLAS, etc.)
Definition: size.hpp:216
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...
Definition: size.hpp:287
size_type stride2() const
Returns the number of columns.
Definition: matrix_def.hpp:225
This file provides the forward declarations for the main types used within ViennaCL.
A dense matrix class.
Definition: forwards.h:374
vcl_size_t internal_size(vector_base< NumericT > const &vec)
Helper routine for obtaining the buffer length of a ViennaCL vector.
Definition: size.hpp:268
lhs_reference_type lhs() const
Get left hand side operand.
Definition: vector.hpp:74
An expression template class that represents a binary operation that yields a vector.
Definition: forwards.h:238
result_of::size_type< MatrixType >::type size2(MatrixType const &mat)
Generic routine for obtaining the number of columns of a matrix (ViennaCL, uBLAS, etc...
Definition: size.hpp:245
vcl_size_t size(VectorType const &vec)
Generic routine for obtaining the size of a vector (ViennaCL, uBLAS, etc.)
Definition: size.hpp:144
vcl_size_t ld(matrix_base< NumericT > const &mat)
Helper routine for obtaining the internal number of entries per row of a ViennaCL matrix...
Definition: size.hpp:291
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.
Definition: size.hpp:59
rhs_reference_type rhs() const
Get right hand side operand.
Definition: vector.hpp:77
size_type stride1() const
Returns the number of rows.
Definition: matrix_def.hpp:223
std::size_t vcl_size_t
Definition: forwards.h:74
All the predicates used within ViennaCL. Checks for expressions to be vectors, etc.
bool row_major() const
Definition: matrix_def.hpp:239
size_type internal_size2() const
Returns the internal number of columns. Usually required for launching OpenCL kernels only...
Definition: matrix_def.hpp:231
size_type internal_size1() const
Returns the internal number of rows. Usually required for launching OpenCL kernels only...
Definition: matrix_def.hpp:229
T min(const T &lhs, const T &rhs)
Minimum.
Definition: util.hpp:45
size_type internal_size() const
Returns the internal length of the vector, which is given by size() plus the extra memory due to padd...
Definition: vector_def.hpp:120
A collection of compile time type deductions.