ViennaCL - The Vienna Computing Library  1.2.0
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-2011, Institute for Microelectronics,
6  Institute for Analysis and Scientific Computing,
7  TU Wien.
8 
9  -----------------
10  ViennaCL - The Vienna Computing Library
11  -----------------
12 
13  Project Head: Karl Rupp rupp@iue.tuwien.ac.at
14 
15  (A list of authors and contributors can be found in the PDF manual)
16 
17  License: MIT (X11), see file LICENSE in the base directory
18 ============================================================================= */
19 
26 #include "viennacl/forwards.h"
27 #include "viennacl/tools/tools.hpp"
29 #include "viennacl/meta/tag_of.hpp"
30 #include <vector>
31 #include <map>
32 
33 namespace viennacl
34 {
35  //
36  // generic prod function
37  // uses tag dispatch to identify which algorithm
38  // should be called
39  //
40  namespace linalg
41  {
42  #ifdef VIENNACL_HAVE_MTL4
43  // ----------------------------------------------------
44  // mtl4
45  //
46  template< typename MatrixT, typename VectorT >
47  VectorT
48  prod(MatrixT const& matrix, VectorT const& vector,
50  >::type* dummy = 0)
51  {
52  // std::cout << "mtl4 .. " << std::endl;
53  return VectorT(matrix * vector);
54  }
55  #endif
56 
57  #ifdef VIENNACL_HAVE_EIGEN
58  // ----------------------------------------------------
59  // Eigen
60  //
61  template< typename MatrixT, typename VectorT >
62  VectorT
63  prod(MatrixT const& matrix, VectorT const& vector,
65  >::type* dummy = 0)
66  {
67  // std::cout << "ublas .. " << std::endl;
68  return matrix * vector;
69  }
70  #endif
71 
72  #ifdef VIENNACL_HAVE_UBLAS
73  // ----------------------------------------------------
74  // UBLAS
75  //
76  template< typename MatrixT, typename VectorT >
77  VectorT
78  prod(MatrixT const& matrix, VectorT const& vector,
80  >::type* dummy = 0)
81  {
82  // std::cout << "ublas .. " << std::endl;
83  return boost::numeric::ublas::prod(matrix, vector);
84  }
85  #endif
86 
87 
88  // ----------------------------------------------------
89  // STL type
90  //
91 
92  // dense matrix-vector product:
93  template< typename T, typename A1, typename A2, typename VectorT >
94  VectorT
95  prod_impl(std::vector< std::vector<T, A1>, A2 > const& matrix, VectorT const& vector)
96  {
97  VectorT result(matrix.size());
98  for (typename std::vector<T, A1>::size_type i=0; i<matrix.size(); ++i)
99  {
100  result[i] = 0; //we will not assume that VectorT is initialized to zero
101  for (typename std::vector<T, A1>::size_type j=0; j<matrix[i].size(); ++j)
102  result[i] += matrix[i][j] * vector[j];
103  }
104  return result;
105  }
106 
107  // sparse matrix-vector product:
108  template< typename KEY, typename DATA, typename COMPARE, typename AMAP, typename AVEC, typename VectorT >
109  VectorT
110  prod_impl(std::vector< std::map<KEY, DATA, COMPARE, AMAP>, AVEC > const& matrix, VectorT const& vector)
111  {
112  typedef std::vector< std::map<KEY, DATA, COMPARE, AMAP>, AVEC > MatrixType;
113 
114  VectorT result(matrix.size());
115  for (typename MatrixType::size_type i=0; i<matrix.size(); ++i)
116  {
117  result[i] = 0; //we will not assume that VectorT is initialized to zero
118  for (typename std::map<KEY, DATA, COMPARE, AMAP>::const_iterator row_entries = matrix[i].begin();
119  row_entries != matrix[i].end();
120  ++row_entries)
121  result[i] += row_entries->second * vector[row_entries->first];
122  }
123  return result;
124  }
125 
126 
127  template< typename MatrixT, typename VectorT >
128  VectorT
129  prod(MatrixT const& matrix, VectorT const& vector,
131  >::type* dummy = 0)
132  {
133  // std::cout << "std .. " << std::endl;
134  return prod_impl(matrix, vector);
135  }
136 
137  // ----------------------------------------------------
138  // VIENNACL
139  //
140  template< typename MatrixT1, typename MatrixT2 >
141  viennacl::matrix_expression< const MatrixT1,
143  viennacl::op_prod >
144  prod(MatrixT1 const& A,
147  >::type* dummy = 0)
148  {
149  // std::cout << "viennacl .. " << std::endl;
150  return viennacl::matrix_expression< const MatrixT1,
152  viennacl::op_prod >(A, B);
153  }
154 
155 
156  template< typename MatrixT1, typename MatrixT2 >
157  viennacl::matrix_expression< const MatrixT1,
160  op_trans>,
161  viennacl::op_prod >
162  prod(MatrixT1 const & A,
165  op_trans> const & B,
167  >::type* dummy = 0)
168  {
169  // std::cout << "viennacl .. " << std::endl;
170  return viennacl::matrix_expression< const MatrixT1,
173  op_trans>,
174  viennacl::op_prod >(A, B);
175  }
176 
177 
178 
179 
180 
181 
182 
183  template< typename MatrixT, typename NumericT, unsigned int ALIGNMENT >
184  viennacl::vector_expression< const MatrixT,
186  viennacl::op_prod >
187  prod(MatrixT const& matrix,
190  >::type* dummy = 0)
191  {
192  // std::cout << "viennacl .. " << std::endl;
193  return viennacl::linalg::prod_impl(matrix, vector);
194  }
195 
196  template< typename MatrixT, typename NumericT, typename F, unsigned int ALIGNMENT >
197  viennacl::matrix_expression< const MatrixT,
199  viennacl::op_prod >
200  prod(MatrixT const& matrix_A,
203  >::type* dummy = 0)
204  {
205  // std::cout << "viennacl .. " << std::endl;
206  return viennacl::matrix_expression< const MatrixT,
208  viennacl::op_prod >(matrix_A, matrix_B);
209  }
210 
211  template< typename MatrixT, typename NumericT, typename F, unsigned int ALIGNMENT >
212  viennacl::matrix_expression< const MatrixT,
215  viennacl::op_trans >,
216  viennacl::op_prod >
217  prod(MatrixT const& matrix_A,
220  viennacl::op_trans > & matrix_B,
222  >::type* dummy = 0)
223  {
224  // std::cout << "viennacl .. " << std::endl;
225  return viennacl::matrix_expression< const MatrixT,
228  viennacl::op_trans >,
229  viennacl::op_prod >(matrix_A, matrix_B);
230  //return viennacl::linalg::prod_impl(matrix_A, matrix_B);
231  }
232 
233  } // end namespace linalg
234 } // end namespace viennacl
235 #endif
236 
237 
238 
239 
240