ViennaCL - The Vienna Computing Library  1.7.0
Free open-source GPU-accelerated linear algebra and solver library.
prod.hpp
Go to the documentation of this file.
1 #ifndef VIENNACL_LINALG_PROD_HPP_
2 #define VIENNACL_LINALG_PROD_HPP_
3 
4 /* =========================================================================
5  Copyright (c) 2010-2015, 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 manual)
17 
18  License: MIT (X11), see file LICENSE in the base directory
19 ============================================================================= */
20 
27 #include "viennacl/forwards.h"
28 #include "viennacl/tools/tools.hpp"
29 #include "viennacl/meta/enable_if.hpp"
30 #include "viennacl/meta/tag_of.hpp"
31 #include <vector>
32 #include <map>
33 
34 namespace viennacl
35 {
36  //
37  // generic prod function
38  // uses tag dispatch to identify which algorithm
39  // should be called
40  //
41  namespace linalg
42  {
43  #ifdef VIENNACL_WITH_MTL4
44  // ----------------------------------------------------
45  // mtl4
46  //
47  template< typename MatrixT, typename VectorT >
48  typename viennacl::enable_if< viennacl::is_mtl4< typename viennacl::traits::tag_of< MatrixT >::type >::value,
49  VectorT>::type
50  prod(MatrixT const& matrix, VectorT const& vector)
51  {
52  return VectorT(matrix * vector);
53  }
54  #endif
55 
56  #ifdef VIENNACL_WITH_ARMADILLO
57  // ----------------------------------------------------
58  // Armadillo
59  //
60  template<typename NumericT, typename VectorT>
61  VectorT prod(arma::SpMat<NumericT> const& A, VectorT const& vector)
62  {
63  return A * vector;
64  }
65  #endif
66 
67  #ifdef VIENNACL_WITH_EIGEN
68  // ----------------------------------------------------
69  // Eigen
70  //
71  template< typename MatrixT, typename VectorT >
72  typename viennacl::enable_if< viennacl::is_eigen< typename viennacl::traits::tag_of< MatrixT >::type >::value,
73  VectorT>::type
74  prod(MatrixT const& matrix, VectorT const& vector)
75  {
76  return matrix * vector;
77  }
78  #endif
79 
80  #ifdef VIENNACL_WITH_UBLAS
81  // ----------------------------------------------------
82  // UBLAS
83  //
84  template< typename MatrixT, typename VectorT >
85  typename viennacl::enable_if< viennacl::is_ublas< typename viennacl::traits::tag_of< MatrixT >::type >::value,
86  VectorT>::type
87  prod(MatrixT const& matrix, VectorT const& vector)
88  {
89  // std::cout << "ublas .. " << std::endl;
90  return boost::numeric::ublas::prod(matrix, vector);
91  }
92  #endif
93 
94 
95  // ----------------------------------------------------
96  // STL type
97  //
98 
99  // dense matrix-vector product:
100  template< typename T, typename A1, typename A2, typename VectorT >
101  VectorT
102  prod(std::vector< std::vector<T, A1>, A2 > const & matrix, VectorT const& vector)
103  {
104  VectorT result(matrix.size());
105  for (typename std::vector<T, A1>::size_type i=0; i<matrix.size(); ++i)
106  {
107  result[i] = 0; //we will not assume that VectorT is initialized to zero
108  for (typename std::vector<T, A1>::size_type j=0; j<matrix[i].size(); ++j)
109  result[i] += matrix[i][j] * vector[j];
110  }
111  return result;
112  }
113 
114  // sparse matrix-vector product:
115  template< typename KEY, typename DATA, typename COMPARE, typename AMAP, typename AVEC, typename VectorT >
116  VectorT
117  prod(std::vector< std::map<KEY, DATA, COMPARE, AMAP>, AVEC > const& matrix, VectorT const& vector)
118  {
119  typedef std::vector< std::map<KEY, DATA, COMPARE, AMAP>, AVEC > MatrixType;
120 
121  VectorT result(matrix.size());
122  for (typename MatrixType::size_type i=0; i<matrix.size(); ++i)
123  {
124  result[i] = 0; //we will not assume that VectorT is initialized to zero
125  for (typename std::map<KEY, DATA, COMPARE, AMAP>::const_iterator row_entries = matrix[i].begin();
126  row_entries != matrix[i].end();
127  ++row_entries)
128  result[i] += row_entries->second * vector[row_entries->first];
129  }
130  return result;
131  }
132 
133 
134  /*template< typename MatrixT, typename VectorT >
135  VectorT
136  prod(MatrixT const& matrix, VectorT const& vector,
137  typename viennacl::enable_if< viennacl::is_stl< typename viennacl::traits::tag_of< MatrixT >::type >::value
138  >::type* dummy = 0)
139  {
140  // std::cout << "std .. " << std::endl;
141  return prod_impl(matrix, vector);
142  }*/
143 
144  // ----------------------------------------------------
145  // VIENNACL
146  //
147 
148  // standard product:
149  template<typename NumericT>
150  viennacl::matrix_expression< const viennacl::matrix_base<NumericT>,
151  const viennacl::matrix_base<NumericT>,
152  viennacl::op_mat_mat_prod >
153  prod(viennacl::matrix_base<NumericT> const & A,
154  viennacl::matrix_base<NumericT> const & B)
155  {
156  return viennacl::matrix_expression< const viennacl::matrix_base<NumericT>,
157  const viennacl::matrix_base<NumericT>,
158  viennacl::op_mat_mat_prod >(A, B);
159  }
160 
161  // right factor is a matrix expression:
162  template<typename NumericT, typename LhsT, typename RhsT, typename OpT>
163  viennacl::matrix_expression< const viennacl::matrix_base<NumericT>,
164  const viennacl::matrix_expression<const LhsT, const RhsT, OpT>,
165  viennacl::op_mat_mat_prod >
166  prod(viennacl::matrix_base<NumericT> const & A,
167  viennacl::matrix_expression<const LhsT, const RhsT, OpT> const & B)
168  {
169  return viennacl::matrix_expression< const viennacl::matrix_base<NumericT>,
170  const viennacl::matrix_expression<const LhsT, const RhsT, OpT>,
171  viennacl::op_mat_mat_prod >(A, B);
172  }
173 
174  // left factor is a matrix expression:
175  template<typename LhsT, typename RhsT, typename OpT, typename NumericT>
176  viennacl::matrix_expression< const viennacl::matrix_expression<const LhsT, const RhsT, OpT>,
177  const viennacl::matrix_base<NumericT>,
178  viennacl::op_mat_mat_prod >
179  prod(viennacl::matrix_expression<const LhsT, const RhsT, OpT> const & A,
180  viennacl::matrix_base<NumericT> const & B)
181  {
182  return viennacl::matrix_expression< const viennacl::matrix_expression<const LhsT, const RhsT, OpT>,
183  const viennacl::matrix_base<NumericT>,
184  viennacl::op_mat_mat_prod >(A, B);
185  }
186 
187 
188  // both factors transposed:
189  template<typename LhsT1, typename RhsT1, typename OpT1,
190  typename LhsT2, typename RhsT2, typename OpT2>
191  viennacl::matrix_expression< const viennacl::matrix_expression<const LhsT1, const RhsT1, OpT1>,
192  const viennacl::matrix_expression<const LhsT2, const RhsT2, OpT2>,
193  viennacl::op_mat_mat_prod >
194  prod(viennacl::matrix_expression<const LhsT1, const RhsT1, OpT1> const & A,
195  viennacl::matrix_expression<const LhsT2, const RhsT2, OpT2> const & B)
196  {
197  return viennacl::matrix_expression< const viennacl::matrix_expression<const LhsT1, const RhsT1, OpT1>,
198  const viennacl::matrix_expression<const LhsT2, const RhsT2, OpT2>,
199  viennacl::op_mat_mat_prod >(A, B);
200  }
201 
202 
203 
204  // matrix-vector product
205  template< typename NumericT>
206  viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
207  const viennacl::vector_base<NumericT>,
208  viennacl::op_prod >
209  prod(viennacl::matrix_base<NumericT> const & A,
210  viennacl::vector_base<NumericT> const & x)
211  {
212  return viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
213  const viennacl::vector_base<NumericT>,
214  viennacl::op_prod >(A, x);
215  }
216 
217  // matrix-vector product (resolve ambiguity)
218  template<typename NumericT, typename F>
219  viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
220  const viennacl::vector_base<NumericT>,
221  viennacl::op_prod >
222  prod(viennacl::matrix<NumericT, F> const & A,
223  viennacl::vector_base<NumericT> const & x)
224  {
225  return viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
226  const viennacl::vector_base<NumericT>,
227  viennacl::op_prod >(A, x);
228  }
229 
230  // matrix-vector product (resolve ambiguity)
231  template<typename MatrixT, typename NumericT>
232  viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
233  const viennacl::vector_base<NumericT>,
234  viennacl::op_prod >
235  prod(viennacl::matrix_range<MatrixT> const & A,
236  viennacl::vector_base<NumericT> const & x)
237  {
238  return viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
239  const viennacl::vector_base<NumericT>,
240  viennacl::op_prod >(A, x);
241  }
242 
243  // matrix-vector product (resolve ambiguity)
244  template<typename MatrixT, typename NumericT>
245  viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
246  const viennacl::vector_base<NumericT>,
247  viennacl::op_prod >
248  prod(viennacl::matrix_slice<MatrixT> const & A,
249  viennacl::vector_base<NumericT> const & x)
250  {
251  return viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
252  const viennacl::vector_base<NumericT>,
253  viennacl::op_prod >(A, x);
254  }
255 
256  // matrix-vector product with matrix expression (including transpose)
257  template< typename NumericT, typename LhsT, typename RhsT, typename OpT>
258  viennacl::vector_expression< const viennacl::matrix_expression<const LhsT, const RhsT, OpT>,
259  const viennacl::vector_base<NumericT>,
260  viennacl::op_prod >
261  prod(viennacl::matrix_expression<const LhsT, const RhsT, OpT> const & A,
262  viennacl::vector_base<NumericT> const & x)
263  {
264  return viennacl::vector_expression< const viennacl::matrix_expression<const LhsT, const RhsT, OpT>,
265  const viennacl::vector_base<NumericT>,
266  viennacl::op_prod >(A, x);
267  }
268 
269 
270  // matrix-vector product with vector expression
271  template< typename NumericT, typename LhsT, typename RhsT, typename OpT>
272  viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
273  const viennacl::vector_expression<const LhsT, const RhsT, OpT>,
274  viennacl::op_prod >
275  prod(viennacl::matrix_base<NumericT> const & A,
276  viennacl::vector_expression<const LhsT, const RhsT, OpT> const & x)
277  {
278  return viennacl::vector_expression< const viennacl::matrix_base<NumericT>,
279  const viennacl::vector_expression<const LhsT, const RhsT, OpT>,
280  viennacl::op_prod >(A, x);
281  }
282 
283 
284  // matrix-vector product with matrix expression (including transpose) and vector expression
285  template<typename LhsT1, typename RhsT1, typename OpT1,
286  typename LhsT2, typename RhsT2, typename OpT2>
287  viennacl::vector_expression< const viennacl::matrix_expression<const LhsT1, const RhsT1, OpT1>,
288  const viennacl::vector_expression<const LhsT2, const RhsT2, OpT2>,
289  viennacl::op_prod >
290  prod(viennacl::matrix_expression<const LhsT1, const RhsT1, OpT1> const & A,
291  viennacl::vector_expression<const LhsT2, const RhsT2, OpT2> const & x)
292  {
293  return viennacl::vector_expression< const viennacl::matrix_expression<const LhsT1, const RhsT1, OpT1>,
294  const viennacl::vector_expression<const LhsT2, const RhsT2, OpT2>,
295  viennacl::op_prod >(A, x);
296  }
297 
298 
299 
300 
301  template< typename SparseMatrixType, typename SCALARTYPE>
302  typename viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value,
303  viennacl::matrix_expression<const SparseMatrixType,
304  const matrix_base <SCALARTYPE>,
305  op_prod >
306  >::type
307  prod(const SparseMatrixType & sp_mat,
308  const viennacl::matrix_base<SCALARTYPE> & d_mat)
309  {
310  return viennacl::matrix_expression<const SparseMatrixType,
311  const viennacl::matrix_base<SCALARTYPE>,
312  op_prod >(sp_mat, d_mat);
313  }
314 
315  // right factor is transposed
316  template< typename SparseMatrixType, typename SCALARTYPE>
317  typename viennacl::enable_if< viennacl::is_any_sparse_matrix<SparseMatrixType>::value,
318  viennacl::matrix_expression< const SparseMatrixType,
319  const viennacl::matrix_expression<const viennacl::matrix_base<SCALARTYPE>,
320  const viennacl::matrix_base<SCALARTYPE>,
321  op_trans>,
322  viennacl::op_prod >
323  >::type
324  prod(const SparseMatrixType & A,
325  viennacl::matrix_expression<const viennacl::matrix_base<SCALARTYPE>,
326  const viennacl::matrix_base<SCALARTYPE>,
327  op_trans> const & B)
328  {
329  return viennacl::matrix_expression< const SparseMatrixType,
330  const viennacl::matrix_expression<const viennacl::matrix_base<SCALARTYPE>,
331  const viennacl::matrix_base<SCALARTYPE>,
332  op_trans>,
333  viennacl::op_prod >(A, B);
334  }
335 
336 
338  template<typename NumericT>
339  viennacl::matrix_expression<const compressed_matrix<NumericT>,
340  const compressed_matrix<NumericT>,
341  op_prod >
342  prod(compressed_matrix<NumericT> const & A,
343  compressed_matrix<NumericT> const & B)
344  {
345  return viennacl::matrix_expression<const compressed_matrix<NumericT>,
346  const compressed_matrix<NumericT>,
347  op_prod >(A, B);
348  }
349 
351  template<typename SparseMatrixType, typename NumericT>
352  vector_expression<const SparseMatrixType,
353  const vector_base<NumericT>,
354  op_prod >
355  prod(const SparseMatrixType & A,
356  const vector_base<NumericT> & x)
357  {
358  return vector_expression<const SparseMatrixType,
359  const vector_base<NumericT>,
360  op_prod >(A, x);
361  }
362 
363  } // end namespace linalg
364 } // end namespace viennacl
365 #endif
366 
367 
368 
369 
370 
viennacl::enable_if
Simple enable-if variant that uses the SFINAE pattern.
Definition: enable_if.hpp:30
tag_of.hpp
Dispatch facility for distinguishing between ublas, STL and ViennaCL types.
viennacl::matrix_slice
Class for representing strided submatrices of a bigger matrix A.
Definition: forwards.h:443
tools.hpp
Various little tools used here and there in ViennaCL.
viennacl::matrix_base< NumericT >
viennacl::matrix_expression
Expression template class for representing a tree of expressions which ultimately result in a matrix...
Definition: forwards.h:341
forwards.h
This file provides the forward declarations for the main types used within ViennaCL.
viennacl::matrix
A dense matrix class.
Definition: forwards.h:375
viennacl::vector_expression
An expression template class that represents a binary operation that yields a vector.
Definition: forwards.h:239
viennacl
Main namespace in ViennaCL. Holds all the basic types such as vector, matrix, etc. and defines operations upon them.
Definition: cpu_ram.hpp:34
viennacl::linalg::prod
VectorT prod(std::vector< std::vector< T, A1 >, A2 > const &matrix, VectorT const &vector)
Definition: prod.hpp:102
viennacl::vector_base< NumericT >
viennacl::op_mat_mat_prod
A tag class representing matrix-matrix products.
Definition: forwards.h:96
viennacl::op_prod
A tag class representing matrix-vector products and element-wise multiplications. ...
Definition: forwards.h:94
prod
void prod(std::vector< std::map< IndexT, NumericT > > const &stl_A, std::vector< std::map< IndexT, NumericT > > const &stl_B, std::vector< std::map< IndexT, NumericT > > &stl_C)
Definition: sparse_prod.cpp:114
viennacl::op_trans
A tag class representing transposed matrices.
Definition: forwards.h:220
viennacl::compressed_matrix< NumericT >
viennacl::matrix_range
Class for representing non-strided submatrices of a bigger matrix A.
Definition: forwards.h:440
enable_if.hpp
Simple enable-if variant that uses the SFINAE pattern.