ViennaCL - The Vienna Computing Library  1.6.2
Free open-source GPU-accelerated linear algebra and solver library.
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
matrix_operations.hpp
Go to the documentation of this file.
1 #ifndef VIENNACL_LINALG_MATRIX_OPERATIONS_HPP_
2 #define VIENNACL_LINALG_MATRIX_OPERATIONS_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 "viennacl/forwards.h"
26 #include "viennacl/scalar.hpp"
27 #include "viennacl/vector.hpp"
29 #include "viennacl/tools/tools.hpp"
33 #include "viennacl/traits/size.hpp"
37 #include "viennacl/vector.hpp"
39 
40 #ifdef VIENNACL_WITH_OPENCL
42 #endif
43 
44 #ifdef VIENNACL_WITH_CUDA
46 #endif
47 
48 namespace viennacl
49 {
50  namespace linalg
51  {
52 
53 
54  template<typename NumericT,
55  typename SizeT, typename DistanceT>
57  matrix_base<NumericT> & temp_trans)
58  {
59  switch (viennacl::traits::handle(proxy).get_active_handle_id())
60  {
62  viennacl::linalg::host_based::trans(proxy, temp_trans);
63  break;
64 #ifdef VIENNACL_WITH_OPENCL
66  viennacl::linalg::opencl::trans(proxy,temp_trans);
67  break;
68 #endif
69 #ifdef VIENNACL_WITH_CUDA
71  viennacl::linalg::cuda::trans(proxy,temp_trans);
72  break;
73 #endif
75  throw memory_exception("not initialised!");
76  default:
77  throw memory_exception("not implemented");
78  }
79  }
80 
81 
82  template<typename NumericT,
83  typename ScalarType1>
85  matrix_base<NumericT> const & mat2, ScalarType1 const & alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha)
86  {
87  switch (viennacl::traits::handle(mat1).get_active_handle_id())
88  {
90  viennacl::linalg::host_based::am(mat1, mat2, alpha, len_alpha, reciprocal_alpha, flip_sign_alpha);
91  break;
92 #ifdef VIENNACL_WITH_OPENCL
94  viennacl::linalg::opencl::am(mat1, mat2, alpha, len_alpha, reciprocal_alpha, flip_sign_alpha);
95  break;
96 #endif
97 #ifdef VIENNACL_WITH_CUDA
99  viennacl::linalg::cuda::am(mat1, mat2, alpha, len_alpha, reciprocal_alpha, flip_sign_alpha);
100  break;
101 #endif
103  throw memory_exception("not initialised!");
104  default:
105  throw memory_exception("not implemented");
106  }
107  }
108 
109 
110  template<typename NumericT,
111  typename ScalarType1, typename ScalarType2>
113  matrix_base<NumericT> const & mat2, ScalarType1 const & alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha,
114  matrix_base<NumericT> const & mat3, ScalarType2 const & beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
115  {
116  switch (viennacl::traits::handle(mat1).get_active_handle_id())
117  {
120  mat2, alpha, len_alpha, reciprocal_alpha, flip_sign_alpha,
121  mat3, beta, len_beta, reciprocal_beta, flip_sign_beta);
122  break;
123 #ifdef VIENNACL_WITH_OPENCL
126  mat2, alpha, len_alpha, reciprocal_alpha, flip_sign_alpha,
127  mat3, beta, len_beta, reciprocal_beta, flip_sign_beta);
128  break;
129 #endif
130 #ifdef VIENNACL_WITH_CUDA
133  mat2, alpha, len_alpha, reciprocal_alpha, flip_sign_alpha,
134  mat3, beta, len_beta, reciprocal_beta, flip_sign_beta);
135  break;
136 #endif
138  throw memory_exception("not initialised!");
139  default:
140  throw memory_exception("not implemented");
141  }
142  }
143 
144 
145  template<typename NumericT,
146  typename ScalarType1, typename ScalarType2>
148  matrix_base<NumericT> const & mat2, ScalarType1 const & alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha,
149  matrix_base<NumericT> const & mat3, ScalarType2 const & beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
150  {
151  switch (viennacl::traits::handle(mat1).get_active_handle_id())
152  {
155  mat2, alpha, len_alpha, reciprocal_alpha, flip_sign_alpha,
156  mat3, beta, len_beta, reciprocal_beta, flip_sign_beta);
157  break;
158 #ifdef VIENNACL_WITH_OPENCL
161  mat2, alpha, len_alpha, reciprocal_alpha, flip_sign_alpha,
162  mat3, beta, len_beta, reciprocal_beta, flip_sign_beta);
163  break;
164 #endif
165 #ifdef VIENNACL_WITH_CUDA
168  mat2, alpha, len_alpha, reciprocal_alpha, flip_sign_alpha,
169  mat3, beta, len_beta, reciprocal_beta, flip_sign_beta);
170  break;
171 #endif
173  throw memory_exception("not initialised!");
174  default:
175  throw memory_exception("not implemented");
176  }
177  }
178 
179 
180  template<typename NumericT>
181  void matrix_assign(matrix_base<NumericT> & mat, NumericT s, bool clear = false)
182  {
183  switch (viennacl::traits::handle(mat).get_active_handle_id())
184  {
187  break;
188 #ifdef VIENNACL_WITH_OPENCL
191  break;
192 #endif
193 #ifdef VIENNACL_WITH_CUDA
196  break;
197 #endif
199  throw memory_exception("not initialised!");
200  default:
201  throw memory_exception("not implemented");
202  }
203  }
204 
205 
206  template<typename NumericT>
208  {
209  switch (viennacl::traits::handle(mat).get_active_handle_id())
210  {
213  break;
214 #ifdef VIENNACL_WITH_OPENCL
217  break;
218 #endif
219 #ifdef VIENNACL_WITH_CUDA
222  break;
223 #endif
225  throw memory_exception("not initialised!");
226  default:
227  throw memory_exception("not implemented");
228  }
229  }
230 
231 
233  template<typename NumericT>
235  {
236  switch (viennacl::traits::handle(v).get_active_handle_id())
237  {
240  break;
241 #ifdef VIENNACL_WITH_OPENCL
244  break;
245 #endif
246 #ifdef VIENNACL_WITH_CUDA
249  break;
250 #endif
252  throw memory_exception("not initialised!");
253  default:
254  throw memory_exception("not implemented");
255  }
256  }
257 
259  template<typename NumericT>
261  {
262  switch (viennacl::traits::handle(A).get_active_handle_id())
263  {
266  break;
267 #ifdef VIENNACL_WITH_OPENCL
270  break;
271 #endif
272 #ifdef VIENNACL_WITH_CUDA
275  break;
276 #endif
278  throw memory_exception("not initialised!");
279  default:
280  throw memory_exception("not implemented");
281  }
282  }
283 
284  template<typename NumericT>
285  void matrix_row(const matrix_base<NumericT> & A, unsigned int i, vector_base<NumericT> & v)
286  {
287  switch (viennacl::traits::handle(A).get_active_handle_id())
288  {
291  break;
292 #ifdef VIENNACL_WITH_OPENCL
295  break;
296 #endif
297 #ifdef VIENNACL_WITH_CUDA
300  break;
301 #endif
303  throw memory_exception("not initialised!");
304  default:
305  throw memory_exception("not implemented");
306  }
307  }
308 
309  template<typename NumericT>
310  void matrix_column(const matrix_base<NumericT> & A, unsigned int j, vector_base<NumericT> & v)
311  {
312  switch (viennacl::traits::handle(A).get_active_handle_id())
313  {
316  break;
317 #ifdef VIENNACL_WITH_OPENCL
320  break;
321 #endif
322 #ifdef VIENNACL_WITH_CUDA
325  break;
326 #endif
328  throw memory_exception("not initialised!");
329  default:
330  throw memory_exception("not implemented");
331  }
332  }
333 
341  template<typename T>
343  scalar<T> & result)
344  {
345  typedef typename matrix_base<T>::handle_type HandleType;
346 
347  if ((A.start1() > 0) || (A.start2() > 0) || (A.stride1() > 1) || (A.stride2() > 1)) {
348  if (A.row_major()) {
350  viennacl::vector_base<T> temp(const_cast<HandleType &>(temp_A.handle()), temp_A.internal_size(), 0, 1);
351  norm_2_impl(temp, result);
352  } else {
354  viennacl::vector_base<T> temp(const_cast<HandleType &>(temp_A.handle()), temp_A.internal_size(), 0, 1);
355  norm_2_impl(temp, result);
356  }
357  } else {
358  viennacl::vector_base<T> temp(const_cast<HandleType &>(A.handle()), A.internal_size(), 0, 1);
359  norm_2_impl(temp, result);
360  }
361 
362  }
363 
371  template<typename T>
373  T & result)
374  {
375  typedef typename matrix_base<T>::handle_type HandleType;
376 
377  if ((A.start1() > 0) || (A.start2() > 0) || (A.stride1() > 1) || (A.stride2() > 1)) {
378  if (A.row_major()) {
380  viennacl::vector_base<T> temp(const_cast<HandleType &>(temp_A.handle()), temp_A.internal_size(), 0, 1);
381  norm_2_cpu(temp, result);
382  } else {
384  viennacl::vector_base<T> temp(const_cast<HandleType &>(temp_A.handle()), temp_A.internal_size(), 0, 1);
385  norm_2_cpu(temp, result);
386  }
387  } else {
388  viennacl::vector_base<T> temp(const_cast<HandleType &>(A.handle()), A.internal_size(), 0, 1);
389  norm_2_cpu(temp, result);
390  }
391 
392  }
393 
394  //
396  //
397 
398 
399 
400  // A * x
401 
410  template<typename NumericT>
412  const vector_base<NumericT> & vec,
413  vector_base<NumericT> & result)
414  {
415  assert( (viennacl::traits::size1(mat) == viennacl::traits::size(result)) && bool("Size check failed at v1 = prod(A, v2): size1(A) != size(v1)"));
416  assert( (viennacl::traits::size2(mat) == viennacl::traits::size(vec)) && bool("Size check failed at v1 = prod(A, v2): size2(A) != size(v2)"));
417 
418  switch (viennacl::traits::handle(mat).get_active_handle_id())
419  {
421  viennacl::linalg::host_based::prod_impl(mat, false, vec, result);
422  break;
423 #ifdef VIENNACL_WITH_OPENCL
425  viennacl::linalg::opencl::prod_impl(mat, false, vec, result);
426  break;
427 #endif
428 #ifdef VIENNACL_WITH_CUDA
430  viennacl::linalg::cuda::prod_impl(mat, false, vec, result);
431  break;
432 #endif
434  throw memory_exception("not initialised!");
435  default:
436  throw memory_exception("not implemented");
437  }
438  }
439 
440 
441  // trans(A) * x
442 
451  template<typename NumericT>
453  const vector_base<NumericT> & vec,
454  vector_base<NumericT> & result)
455  {
456  assert( (viennacl::traits::size1(mat_trans.lhs()) == viennacl::traits::size(vec)) && bool("Size check failed at v1 = trans(A) * v2: size1(A) != size(v2)"));
457  assert( (viennacl::traits::size2(mat_trans.lhs()) == viennacl::traits::size(result)) && bool("Size check failed at v1 = trans(A) * v2: size2(A) != size(v1)"));
458 
459  switch (viennacl::traits::handle(mat_trans.lhs()).get_active_handle_id())
460  {
462  viennacl::linalg::host_based::prod_impl(mat_trans.lhs(), true, vec, result);
463  break;
464 #ifdef VIENNACL_WITH_OPENCL
466  viennacl::linalg::opencl::prod_impl(mat_trans.lhs(), true, vec, result);
467  break;
468 #endif
469 #ifdef VIENNACL_WITH_CUDA
471  viennacl::linalg::cuda::prod_impl(mat_trans.lhs(), true, vec, result);
472  break;
473 #endif
475  throw memory_exception("not initialised!");
476  default:
477  throw memory_exception("not implemented");
478  }
479  }
480 
481 
482  //
484  //
485 
491  template<typename NumericT, typename ScalarType >
493  const matrix_base<NumericT> & B,
495  ScalarType alpha,
496  ScalarType beta)
497  {
498  assert( (viennacl::traits::size1(A) == viennacl::traits::size1(C)) && bool("Size check failed at C = prod(A, B): size1(A) != size1(C)"));
499  assert( (viennacl::traits::size2(A) == viennacl::traits::size1(B)) && bool("Size check failed at C = prod(A, B): size2(A) != size1(B)"));
500  assert( (viennacl::traits::size2(B) == viennacl::traits::size2(C)) && bool("Size check failed at C = prod(A, B): size2(B) != size2(C)"));
501 
502 
503  switch (viennacl::traits::handle(A).get_active_handle_id())
504  {
506  viennacl::linalg::host_based::prod_impl(A, false, B, false, C, alpha, beta);
507  break;
508 #ifdef VIENNACL_WITH_OPENCL
510  viennacl::linalg::opencl::prod_impl(A, false, B, false, C, alpha, beta);
511  break;
512 #endif
513 #ifdef VIENNACL_WITH_CUDA
515  viennacl::linalg::cuda::prod_impl(A, false, B, false, C, alpha, beta);
516  break;
517 #endif
519  throw memory_exception("not initialised!");
520  default:
521  throw memory_exception("not implemented");
522  }
523  }
524 
525 
526 
532  template<typename NumericT, typename ScalarType >
534  const matrix_base<NumericT>,
535  op_trans> & A,
536  const matrix_base<NumericT> & B,
538  ScalarType alpha,
539  ScalarType beta)
540  {
541  assert(viennacl::traits::size2(A.lhs()) == viennacl::traits::size1(C) && bool("Size check failed at C = prod(trans(A), B): size2(A) != size1(C)"));
542  assert(viennacl::traits::size1(A.lhs()) == viennacl::traits::size1(B) && bool("Size check failed at C = prod(trans(A), B): size1(A) != size1(B)"));
543  assert(viennacl::traits::size2(B) == viennacl::traits::size2(C) && bool("Size check failed at C = prod(trans(A), B): size2(B) != size2(C)"));
544 
545  switch (viennacl::traits::handle(A.lhs()).get_active_handle_id())
546  {
548  viennacl::linalg::host_based::prod_impl(A.lhs(), true, B, false, C, alpha, beta);
549  break;
550 #ifdef VIENNACL_WITH_OPENCL
552  viennacl::linalg::opencl::prod_impl(A.lhs(), true, B, false, C, alpha, beta);
553  break;
554 #endif
555 #ifdef VIENNACL_WITH_CUDA
557  viennacl::linalg::cuda::prod_impl(A.lhs(), true, B, false, C, alpha, beta);
558  break;
559 #endif
561  throw memory_exception("not initialised!");
562  default:
563  throw memory_exception("not implemented");
564  }
565  }
566 
567 
568 
569 
575  template<typename NumericT, typename ScalarType >
579  ScalarType alpha,
580  ScalarType beta)
581  {
582  assert(viennacl::traits::size1(A) == viennacl::traits::size1(C) && bool("Size check failed at C = prod(A, trans(B)): size1(A) != size1(C)"));
583  assert(viennacl::traits::size2(A) == viennacl::traits::size2(B.lhs()) && bool("Size check failed at C = prod(A, trans(B)): size2(A) != size2(B)"));
584  assert(viennacl::traits::size1(B.lhs()) == viennacl::traits::size2(C) && bool("Size check failed at C = prod(A, trans(B)): size1(B) != size2(C)"));
585 
587  {
589  viennacl::linalg::host_based::prod_impl(A, false, B.lhs(), true, C, alpha, beta);
590  break;
591 #ifdef VIENNACL_WITH_OPENCL
593  viennacl::linalg::opencl::prod_impl(A, false, B.lhs(), true, C, alpha, beta);
594  break;
595 #endif
596 #ifdef VIENNACL_WITH_CUDA
598  viennacl::linalg::cuda::prod_impl(A, false, B.lhs(), true, C, alpha, beta);
599  break;
600 #endif
602  throw memory_exception("not initialised!");
603  default:
604  throw memory_exception("not implemented");
605  }
606  }
607 
608 
609 
615  template<typename NumericT, typename ScalarType >
619  ScalarType alpha,
620  ScalarType beta)
621  {
622  assert(viennacl::traits::size2(A.lhs()) == viennacl::traits::size1(C) && bool("Size check failed at C = prod(trans(A), trans(B)): size2(A) != size1(C)"));
623  assert(viennacl::traits::size1(A.lhs()) == viennacl::traits::size2(B.lhs()) && bool("Size check failed at C = prod(trans(A), trans(B)): size1(A) != size2(B)"));
624  assert(viennacl::traits::size1(B.lhs()) == viennacl::traits::size2(C) && bool("Size check failed at C = prod(trans(A), trans(B)): size1(B) != size2(C)"));
625 
626  switch (viennacl::traits::handle(A.lhs()).get_active_handle_id())
627  {
629  viennacl::linalg::host_based::prod_impl(A.lhs(), true, B.lhs(), true, C, alpha, beta);
630  break;
631 #ifdef VIENNACL_WITH_OPENCL
633  viennacl::linalg::opencl::prod_impl(A.lhs(), true, B.lhs(), true, C, alpha, beta);
634  break;
635 #endif
636 #ifdef VIENNACL_WITH_CUDA
638  viennacl::linalg::cuda::prod_impl(A.lhs(), true, B.lhs(), true, C, alpha, beta);
639  break;
640 #endif
642  throw memory_exception("not initialised!");
643  default:
644  throw memory_exception("not implemented");
645  }
646  }
647 
648 
650 
651 
652 
658  template<typename T, typename OP>
660  matrix_expression<const matrix_base<T>, const matrix_base<T>, OP> const & proxy)
661  {
662  assert( (viennacl::traits::size1(A) == viennacl::traits::size1(proxy)) && bool("Size check failed at A = element_op(B): size1(A) != size1(B)"));
663  assert( (viennacl::traits::size2(A) == viennacl::traits::size2(proxy)) && bool("Size check failed at A = element_op(B): size2(A) != size2(B)"));
664 
665  switch (viennacl::traits::handle(A).get_active_handle_id())
666  {
669  break;
670 #ifdef VIENNACL_WITH_OPENCL
673  break;
674 #endif
675 #ifdef VIENNACL_WITH_CUDA
678  break;
679 #endif
681  throw memory_exception("not initialised!");
682  default:
683  throw memory_exception("not implemented");
684  }
685  }
686 
687 
688 #define VIENNACL_MAKE_BINARY_OP(OPNAME)\
689  template<typename T>\
690  viennacl::matrix_expression<const matrix_base<T>, const matrix_base<T>, op_element_binary<op_##OPNAME> >\
691  element_##OPNAME(matrix_base<T> const & A, matrix_base<T> const & B)\
692  {\
693  return viennacl::matrix_expression<const matrix_base<T>, const matrix_base<T>, op_element_binary<op_##OPNAME> >(A, B);\
694  }\
695 \
696  template<typename M1, typename M2, typename OP, typename T>\
697  viennacl::matrix_expression<const matrix_expression<const M1, const M2, OP>,\
698  const matrix_base<T>,\
699  op_element_binary<op_##OPNAME> >\
700  element_##OPNAME(matrix_expression<const M1, const M2, OP> const & proxy, matrix_base<T> const & B)\
701  {\
702  return viennacl::matrix_expression<const matrix_expression<const M1, const M2, OP>,\
703  const matrix_base<T>,\
704  op_element_binary<op_##OPNAME> >(proxy, B);\
705  }\
706 \
707  template<typename T, typename M2, typename M3, typename OP>\
708  viennacl::matrix_expression<const matrix_base<T>,\
709  const matrix_expression<const M2, const M3, OP>,\
710  op_element_binary<op_##OPNAME> >\
711  element_##OPNAME(matrix_base<T> const & A, matrix_expression<const M2, const M3, OP> const & proxy)\
712  {\
713  return viennacl::matrix_expression<const matrix_base<T>,\
714  const matrix_expression<const M2, const M3, OP>,\
715  op_element_binary<op_##OPNAME> >(A, proxy);\
716  }\
717 \
718  template<typename M1, typename M2, typename OP1,\
719  typename M3, typename M4, typename OP2>\
720  viennacl::matrix_expression<const matrix_expression<const M1, const M2, OP1>,\
721  const matrix_expression<const M3, const M4, OP2>,\
722  op_element_binary<op_##OPNAME> >\
723  element_##OPNAME(matrix_expression<const M1, const M2, OP1> const & proxy1,\
724  matrix_expression<const M3, const M4, OP2> const & proxy2)\
725  {\
726  return viennacl::matrix_expression<const matrix_expression<const M1, const M2, OP1>,\
727  const matrix_expression<const M3, const M4, OP2>,\
728  op_element_binary<op_##OPNAME> >(proxy1, proxy2);\
729  }
730 
734 
737  VIENNACL_MAKE_BINARY_OP(greater)
741 
742 #undef VIENNACL_GENERATE_BINARY_OP_OVERLOADS
743 
744 
745 
746 #define VIENNACL_MAKE_UNARY_ELEMENT_OP(funcname) \
747  template<typename T> \
748  viennacl::matrix_expression<const matrix_base<T>, const matrix_base<T>, op_element_unary<op_##funcname> > \
749  element_##funcname(matrix_base<T> const & A) \
750  { \
751  return viennacl::matrix_expression<const matrix_base<T>, const matrix_base<T>, op_element_unary<op_##funcname> >(A, A); \
752  } \
753  template<typename LHS, typename RHS, typename OP> \
754  viennacl::matrix_expression<const matrix_expression<const LHS, const RHS, OP>, \
755  const matrix_expression<const LHS, const RHS, OP>, \
756  op_element_unary<op_##funcname> > \
757  element_##funcname(matrix_expression<const LHS, const RHS, OP> const & proxy) \
758  { \
759  return viennacl::matrix_expression<const matrix_expression<const LHS, const RHS, OP>, \
760  const matrix_expression<const LHS, const RHS, OP>, \
761  op_element_unary<op_##funcname> >(proxy, proxy); \
762  } \
763 
781 
782 #undef VIENNACL_MAKE_UNARY_ELEMENT_OP
783 
784 
785  //
787  //
788 
789 
795  template<typename NumericT>
796  viennacl::matrix_expression<const vector_base<NumericT>, const vector_base<NumericT>, op_prod>
798  {
800  }
801 
802 
815  template<typename NumericT, typename S1>
817  S1 const & alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha,
818  const vector_base<NumericT> & vec1,
819  const vector_base<NumericT> & vec2)
820  {
821  switch (viennacl::traits::handle(mat1).get_active_handle_id())
822  {
825  alpha, len_alpha, reciprocal_alpha, flip_sign_alpha,
826  vec1, vec2);
827  break;
828 #ifdef VIENNACL_WITH_OPENCL
831  alpha, len_alpha, reciprocal_alpha, flip_sign_alpha,
832  vec1, vec2);
833  break;
834 #endif
835 #ifdef VIENNACL_WITH_CUDA
838  alpha, len_alpha, reciprocal_alpha, flip_sign_alpha,
839  vec1, vec2);
840  break;
841 #endif
843  throw memory_exception("not initialised!");
844  default:
845  throw memory_exception("not implemented");
846  }
847  }
848 
857  template <typename NumericT, typename VectorType>
859  VectorType & dh,
860  VectorType & sh
861  )
862  {
863  switch (viennacl::traits::handle(A).get_active_handle_id())
864  {
867  break;
868 #ifdef VIENNACL_WITH_OPENCL
871  break;
872 #endif
873 
874 #ifdef VIENNACL_WITH_CUDA
877  break;
878 #endif
879 
881  throw memory_exception("not initialised!");
882  default:
883  throw memory_exception("not implemented");
884  }
885 
886 
887  }
898  template <typename SCALARTYPE>
901  vcl_size_t row_start,
902  vcl_size_t col_start,
903  bool copy_col
904  )
905  {
906  switch (viennacl::traits::handle(A).get_active_handle_id())
907  {
909  viennacl::linalg::host_based::copy_vec(A, V, row_start, col_start, copy_col);
910  break;
911 #ifdef VIENNACL_WITH_OPENCL
913  viennacl::linalg::opencl::copy_vec(A, V, row_start, col_start, copy_col);
914  break;
915 #endif
916 
917 #ifdef VIENNACL_WITH_CUDA
919  viennacl::linalg::cuda::copy_vec(A, V, row_start, col_start, copy_col);
920  break;
921 #endif
922 
924  throw memory_exception("not initialised!");
925  default:
926  throw memory_exception("not implemented");
927  }
928 
929  }
930 
937  template <typename NumericT>
941  {
942  switch (viennacl::traits::handle(A).get_active_handle_id())
943  {
946  break;
947 #ifdef VIENNACL_WITH_OPENCL
950  break;
951 #endif
952 
953 #ifdef VIENNACL_WITH_CUDA
956  break;
957 #endif
958 
960  throw memory_exception("not initialised!");
961  default:
962  throw memory_exception("not implemented");
963  }
964  }
965 
966 
974  template <typename NumericT>
977  {
978  switch (viennacl::traits::handle(A).get_active_handle_id())
979  {
982  break;
983 #ifdef VIENNACL_WITH_OPENCL
986  break;
987 #endif
988 
989 #ifdef VIENNACL_WITH_CUDA
992  break;
993 #endif
994 
996  throw memory_exception("not initialised!");
997  default:
998  throw memory_exception("not implemented");
999  }
1000  }
1001 
1009  template <typename NumericT>
1012  vcl_size_t A_size1)
1013  {
1014  switch (viennacl::traits::handle(Q).get_active_handle_id())
1015  {
1016  case viennacl::MAIN_MEMORY:
1018  break;
1019 #ifdef VIENNACL_WITH_OPENCL
1022  break;
1023 #endif
1024 
1025 #ifdef VIENNACL_WITH_CUDA
1026  case viennacl::CUDA_MEMORY:
1028  break;
1029 #endif
1030 
1032  throw memory_exception("not initialised!");
1033  default:
1034  throw memory_exception("not implemented");
1035  }
1036  }
1037 
1038 
1048  template<typename NumericT>
1050  vector_base<NumericT> & tmp1,
1051  vector_base<NumericT> & tmp2,
1052  int l,
1053  int m
1054  )
1055  {
1056  switch (viennacl::traits::handle(Q).get_active_handle_id())
1057  {
1058  case viennacl::MAIN_MEMORY:
1059  viennacl::linalg::host_based::givens_next(Q, tmp1, tmp2, l, m);
1060  break;
1061 #ifdef VIENNACL_WITH_OPENCL
1063  viennacl::linalg::opencl::givens_next(Q, tmp1, tmp2, l, m);
1064  break;
1065 #endif
1066 
1067 #ifdef VIENNACL_WITH_CUDA
1068  case viennacl::CUDA_MEMORY:
1069  viennacl::linalg::cuda::givens_next(Q, tmp1, tmp2, l, m);
1070  break;
1071 #endif
1072 
1074  throw memory_exception("not initialised!");
1075  default:
1076  throw memory_exception("not implemented");
1077  }
1078  }
1079 
1086  template<typename NumericT>
1088  vector_base<NumericT> & vec1,
1089  vector_base<NumericT> & vec2
1090  )
1091  {
1092  switch (viennacl::traits::handle(vec1).get_active_handle_id())
1093  {
1094  case viennacl::MAIN_MEMORY:
1096  break;
1097 #ifdef VIENNACL_WITH_OPENCL
1100  break;
1101 #endif
1102 
1103 #ifdef VIENNACL_WITH_CUDA
1104  case viennacl::CUDA_MEMORY:
1106  break;
1107 #endif
1108 
1110  throw memory_exception("not initialised!");
1111  default:
1112  throw memory_exception("not implemented");
1113  }
1114  }
1115 
1122  template<typename NumericT>
1124  vector_base<NumericT> & vec1,
1125  vector_base<NumericT> & vec2
1126  )
1127  {
1128  switch (viennacl::traits::handle(vec1).get_active_handle_id())
1129  {
1130  case viennacl::MAIN_MEMORY:
1132  break;
1133 #ifdef VIENNACL_WITH_OPENCL
1136  break;
1137 #endif
1138 
1139 #ifdef VIENNACL_WITH_CUDA
1140  case viennacl::CUDA_MEMORY:
1142  break;
1143 #endif
1144 
1146  throw memory_exception("not initialised!");
1147  default:
1148  throw memory_exception("not implemented");
1149  }
1150  }
1151 
1152  } //namespace linalg
1153 
1154 
1155 
1156 
1157  //
1159  //
1160 
1161 
1162  //v += A * x
1168  template<typename NumericT>
1169  vector<NumericT>
1172  {
1173  assert(viennacl::traits::size1(proxy.lhs()) == v1.size() && bool("Size check failed for v1 += A * v2: size1(A) != size(v1)"));
1174 
1175  vector<NumericT> result(viennacl::traits::size1(proxy.lhs()));
1176  viennacl::linalg::prod_impl(proxy.lhs(), proxy.rhs(), result);
1177  v1 += result;
1178  return v1;
1179  }
1180 
1186  template<typename NumericT>
1187  vector<NumericT>
1190  {
1191  assert(viennacl::traits::size1(proxy.lhs()) == v1.size() && bool("Size check failed for v1 -= A * v2: size1(A) != size(v1)"));
1192 
1193  vector<NumericT> result(viennacl::traits::size1(proxy.lhs()));
1194  viennacl::linalg::prod_impl(proxy.lhs(), proxy.rhs(), result);
1195  v1 -= result;
1196  return v1;
1197  }
1198 
1199 
1200 
1201 
1202 
1203  //free functions:
1209  template<typename NumericT>
1213  {
1214  assert(viennacl::traits::size1(proxy.lhs()) == viennacl::traits::size(v1) && bool("Size check failed for v1 + A * v2: size1(A) != size(v1)"));
1215 
1217  viennacl::linalg::prod_impl(proxy.lhs(), proxy.rhs(), result);
1218  result += v1;
1219  return result;
1220  }
1221 
1227  template<typename NumericT>
1231  {
1232  assert(viennacl::traits::size1(proxy.lhs()) == viennacl::traits::size(v1) && bool("Size check failed for v1 - A * v2: size1(A) != size(v1)"));
1233 
1235  viennacl::linalg::prod_impl(proxy.lhs(), proxy.rhs(), result);
1236  result = v1 - result;
1237  return result;
1238  }
1239 
1240 
1242 
1243 
1244  //v += A^T * x
1250  template<typename NumericT>
1251  vector<NumericT>
1254  const vector_base<NumericT>,
1255  op_prod> & proxy)
1256  {
1257  assert(viennacl::traits::size2(proxy.lhs()) == v1.size() && bool("Size check failed in v1 += trans(A) * v2: size2(A) != size(v1)"));
1258 
1259  vector<NumericT> result(viennacl::traits::size2(proxy.lhs()));
1260  viennacl::linalg::prod_impl(proxy.lhs(), proxy.rhs(), result);
1261  v1 += result;
1262  return v1;
1263  }
1264 
1265  //v -= A^T * x
1271  template<typename NumericT>
1272  vector<NumericT>
1275  const vector_base<NumericT>,
1276  op_prod> & proxy)
1277  {
1278  assert(viennacl::traits::size2(proxy.lhs()) == v1.size() && bool("Size check failed in v1 += trans(A) * v2: size2(A) != size(v1)"));
1279 
1280  vector<NumericT> result(viennacl::traits::size2(proxy.lhs()));
1281  viennacl::linalg::prod_impl(proxy.lhs(), proxy.rhs(), result);
1282  v1 -= result;
1283  return v1;
1284  }
1285 
1286 
1287  //free functions:
1293  template<typename NumericT>
1294  vector<NumericT>
1297  const vector_base<NumericT>,
1298  op_prod> & proxy)
1299  {
1300  assert(viennacl::traits::size2(proxy.lhs()) == viennacl::traits::size(v1) && bool("Size check failed in v1 + trans(A) * v2: size2(A) != size(v1)"));
1301 
1303  viennacl::linalg::prod_impl(proxy.lhs(), proxy.rhs(), result);
1304  result += v1;
1305  return result;
1306  }
1307 
1313  template<typename NumericT>
1314  vector<NumericT>
1317  const vector_base<NumericT>,
1318  op_prod> & proxy)
1319  {
1320  assert(viennacl::traits::size2(proxy.lhs()) == viennacl::traits::size(v1) && bool("Size check failed in v1 - trans(A) * v2: size2(A) != size(v1)"));
1321 
1323  viennacl::linalg::prod_impl(proxy.lhs(), proxy.rhs(), result);
1324  result = v1 - result;
1325  return result;
1326  }
1327 
1328 
1329 } //namespace viennacl
1330 
1331 
1332 #endif
void norm_frobenius_cpu(matrix_base< T > const &vec, T &result)
Computes the Frobenius norm of a vector with final reduction on the CPU.
Implementations of dense matrix related operations, including matrix-vector products, using a plain single-threaded or OpenMP-enabled execution on CPU.
void matrix_diag_from_vector(const vector_base< NumericT > &vec, int k, matrix_base< NumericT > &A)
void trans(const matrix_expression< const matrix_base< NumericT, SizeT, DistanceT >, const matrix_base< NumericT, SizeT, DistanceT >, op_trans > &proxy, matrix_base< NumericT > &temp_trans)
void matrix_diag_to_vector(const matrix_base< NumericT > &A, int k, vector_base< NumericT > &v)
Dispatcher interface for v = diag(A, k)
void house_update_QL(matrix_base< NumericT > &Q, vector_base< NumericT > &D, vcl_size_t A_size1)
This function updates the matrix Q, which is needed for the computation of the eigenvectors.
void am(matrix_base< NumericT > &A, matrix_base< NumericT > const &B, ScalarT1 const &alpha, vcl_size_t, bool reciprocal_alpha, bool flip_sign_alpha)
void matrix_column(const matrix_base< NumericT > &A, unsigned int j, vector_base< NumericT > &vec)
This class represents a single scalar value on the GPU and behaves mostly like a built-in scalar type...
Definition: forwards.h:226
Implementations of dense matrix related operations, including matrix-vector products, using OpenCL.
void house_update_A_right(matrix_base< NumericT > &A, vector_base< NumericT > &D)
This function applies a householder transformation to a matrix: A <- A * P with a householder reflect...
vector< NumericT > operator-=(vector_base< NumericT > &v1, const viennacl::vector_expression< const matrix_base< NumericT >, const vector_base< NumericT >, viennacl::op_prod > &proxy)
Implementation of the operation v1 -= A * v2, where A is a matrix.
void matrix_diag_to_vector(const matrix_base< NumericT > &mat, int k, vector_base< NumericT > &vec)
size_type internal_size() const
Returns the total amount of allocated memory in multiples of sizeof(NumericT)
Definition: matrix_def.hpp:233
Exception class in case of memory errors.
Definition: forwards.h:571
void exclusive_scan(vector_base< NumericT > &vec1, vector_base< NumericT > &vec2)
This function implements an exclusive scan.
void copy_vec(matrix_base< SCALARTYPE > &A, vector_base< SCALARTYPE > &V, vcl_size_t row_start, vcl_size_t col_start, bool copy_col)
This function copies a row or a column from a matrix to a vector.
Generic size and resize functionality for different vector and matrix types.
void trans(matrix_expression< const matrix_base< NumericT, SizeT, DistanceT >, const matrix_base< NumericT, SizeT, DistanceT >, op_trans > const &proxy, matrix_base< NumericT > &temp_trans)
void exclusive_scan(vector_base< NumericT > &vec1, vector_base< NumericT > &vec2)
Extracts the underlying OpenCL start index handle from a vector, a matrix, an expression etc...
void matrix_assign(matrix_base< NumericT > &mat, NumericT s, bool clear=false)
Various little tools used here and there in ViennaCL.
void house_update_A_right(matrix_base< NumericT > &A, vector_base< NumericT > &D)
This function applies a householder transformation to a matrix: A <- A * P with a householder reflect...
vector< NumericT > operator+=(vector_base< NumericT > &v1, const viennacl::vector_expression< const matrix_base< NumericT >, const vector_base< NumericT >, viennacl::op_prod > &proxy)
Implementation of the operation v1 += A * v2, where A is a matrix.
void matrix_assign(matrix_base< NumericT > &mat, NumericT s, bool clear=false)
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
void norm_2_cpu(vector_base< T > const &vec, T &result)
Computes the l^2-norm of a vector with final reduction on the CPU - dispatcher interface.
void house_update_QL(matrix_base< NumericT > &Q, vector_base< NumericT > &D, vcl_size_t A_size1)
This function updates the matrix Q, which is needed for the computation of the eigenvectors.
void am(matrix_base< NumericT > &mat1, matrix_base< NumericT > const &mat2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha)
Expression template class for representing a tree of expressions which ultimately result in a matrix...
Definition: forwards.h:340
void ambm(matrix_base< NumericT > &mat1, matrix_base< NumericT > const &mat2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT > const &mat3, ScalarType2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
size_type stride2() const
Returns the number of columns.
Definition: matrix_def.hpp:225
void clear(VectorType &vec)
Generic routine for setting all entries of a vector to zero. This is the version for non-ViennaCL obj...
Definition: clear.hpp:57
This file provides the forward declarations for the main types used within ViennaCL.
A dense matrix class.
Definition: forwards.h:374
void ambm(matrix_base< NumericT > &mat1, matrix_base< NumericT > const &mat2, ScalarT1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT > const &mat3, ScalarT2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
Determines row and column increments for matrices and matrix proxies.
void bidiag_pack(matrix_base< NumericT > &A, viennacl::vector< NumericT > &dh, viennacl::vector< NumericT > &sh)
void prod_impl(const matrix_base< NumericT > &mat, bool mat_transpose, const vector_base< NumericT > &vec, vector_base< NumericT > &result)
Carries out matrix-vector multiplication.
An expression template class that represents a binary operation that yields a vector.
Definition: forwards.h:238
void element_op(matrix_base< T > &A, matrix_expression< const matrix_base< T >, const matrix_base< T >, OP > const &proxy)
Implementation of the element-wise operation A = B .* C and A = B ./ C for matrices (using MATLAB syn...
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
void scaled_rank_1_update(matrix_base< NumericT > &A, ScalarT1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, const vector_base< NumericT > &vec1, const vector_base< NumericT > &vec2)
The implementation of the operation mat += alpha * vec1 * vec2^T, i.e. a scaled rank 1 update...
viennacl::vector< NumericT > operator-(const vector_base< NumericT > &v1, const vector_expression< const matrix_base< NumericT >, const vector_base< NumericT >, op_prod > &proxy)
Implementation of the operation 'result = v1 - A * v2', where A is a matrix.
void house_update_A_right(matrix_base< NumericT > &A, vector_base< NumericT > &D)
This function applies a householder transformation to a matrix: A <- A * P with a householder reflect...
void copy_vec(matrix_base< NumericT > &A, vector_base< NumericT > &V, vcl_size_t row_start, vcl_size_t col_start, bool copy_col)
void house_update_A_right(matrix_base< NumericT > &A, vector_base< NumericT > &D)
void element_op(matrix_base< NumericT > &A, matrix_expression< const matrix_base< NumericT >, const matrix_base< NumericT >, op_element_binary< OpT > > const &proxy)
Implementation of binary element-wise operations A = OP(B,C)
void norm_2_impl(vector_base< T > const &vec, scalar< T > &result)
Computes the l^2-norm of a vector - dispatcher interface.
void ambm(matrix_base< NumericT > &A, matrix_base< NumericT > const &B, ScalarT1 const &alpha, vcl_size_t, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT > const &C, ScalarT2 const &beta, vcl_size_t, bool reciprocal_beta, bool flip_sign_beta)
void norm_frobenius_impl(matrix_base< T > const &vec, scalar< T > &result)
Computes the Frobenius norm of a matrix - dispatcher interface.
void house_update_QL(matrix_base< NumericT > &Q, vector_base< NumericT > &D, vcl_size_t A_size1)
viennacl::vector< float > v1
VectorT prod(std::vector< std::vector< T, A1 >, A2 > const &matrix, VectorT const &vector)
Definition: prod.hpp:91
vcl_size_t size(VectorType const &vec)
Generic routine for obtaining the size of a vector (ViennaCL, uBLAS, etc.)
Definition: size.hpp:144
void matrix_row(const matrix_base< NumericT > &A, unsigned int i, vector_base< NumericT > &v)
void copy_vec(matrix_base< NumericT > &A, vector_base< S1 > &V, vcl_size_t row_start, vcl_size_t col_start, bool copy_col)
This function copies a row or a column from a matrix to a vector.
void prod_impl(const matrix_base< NumericT > &mat, bool trans, const vector_base< NumericT > &vec, vector_base< NumericT > &result)
Carries out matrix-vector multiplication.
void ambm_m(matrix_base< NumericT > &mat1, matrix_base< NumericT > const &mat2, ScalarT1 const &alpha, vcl_size_t, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT > const &mat3, ScalarT2 const &beta, vcl_size_t, bool reciprocal_beta, bool flip_sign_beta)
void am(matrix_base< NumericT > &mat1, matrix_base< NumericT > const &mat2, ScalarT const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha)
void matrix_diag_to_vector(matrix_base< NumericT > const &mat, int k, vector_base< NumericT > &vec)
#define VIENNACL_MAKE_UNARY_ELEMENT_OP(funcname)
void matrix_column(const matrix_base< NumericT > &mat, unsigned int j, vector_base< NumericT > &vec)
void bidiag_pack(matrix_base< NumericT > &A, VectorType &dh, VectorType &sh)
This function stores the diagonal and the superdiagonal of a matrix in two vectors.
void house_update_A_left(matrix_base< NumericT > &A, vector_base< NumericT > &D, vcl_size_t start)
This function applies a householder transformation to a matrix. A <- P * A with a householder reflect...
void matrix_diagonal_assign(matrix_base< NumericT > &mat, NumericT s)
void element_op(matrix_base< NumericT > &A, matrix_expression< const matrix_base< NumericT >, const matrix_base< NumericT >, op_element_binary< OpT > > const &proxy)
Implementation of the element-wise operations A = B .* C and A = B ./ C (using MATLAB syntax) ...
result_of::size_type< T >::type start(T const &obj)
Definition: start.hpp:44
void house_update_QL(matrix_base< NumericT > &Q, vector_base< NumericT > &D, vcl_size_t A_size1)
This function updates the matrix Q, which is needed for the computation of the eigenvectors.
void inclusive_scan(vector_base< NumericT > &vec1, vector_base< NumericT > &vec2)
This function implements an inclusive scan.
void scaled_rank_1_update(matrix_base< NumericT > &mat1, ScalarT const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, const vector_base< NumericT > &vec1, const vector_base< NumericT > &vec2)
The implementation of the operation mat += alpha * vec1 * vec2^T, i.e. a scaled rank 1 update...
void ambm_m(matrix_base< NumericT > &mat1, matrix_base< NumericT > const &mat2, ScalarT1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT > const &mat3, ScalarT2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
size_type stride1() const
Returns the number of rows.
Definition: matrix_def.hpp:223
void house_update_A_left(matrix_base< NumericT > &A, vector_base< NumericT > &D, vcl_size_t start)
This function applies a householder transformation to a matrix. A <- P * A with a householder reflect...
#define VIENNACL_MAKE_BINARY_OP(OPNAME)
void scaled_rank_1_update(matrix_base< NumericT > &mat1, S1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, const vector_base< NumericT > &vec1, const vector_base< NumericT > &vec2)
The implementation of the operation mat += alpha * vec1 * vec2^T, i.e. a scaled rank 1 update...
void matrix_diag_from_vector(const vector_base< NumericT > &vec, int k, matrix_base< NumericT > &mat)
void matrix_diagonal_assign(matrix_base< NumericT > &mat, NumericT s)
void ambm(matrix_base< NumericT > &mat1, matrix_base< NumericT > const &mat2, ScalarT1 const &alpha, vcl_size_t, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT > const &mat3, ScalarT2 const &beta, vcl_size_t, bool reciprocal_beta, bool flip_sign_beta)
std::size_t vcl_size_t
Definition: forwards.h:74
handle_type & handle()
Returns the OpenCL handle, non-const-version.
Definition: matrix_def.hpp:235
void trans(const matrix_expression< const matrix_base< NumericT, SizeT, DistanceT >, const matrix_base< NumericT, SizeT, DistanceT >, op_trans > &proxy, matrix_base< NumericT > &temp_trans)
void house_update_A_left(matrix_base< NumericT > &A, vector_base< NumericT > &D, vcl_size_t start)
void matrix_column(const matrix_base< NumericT > &A, unsigned int j, vector_base< NumericT > &v)
void inclusive_scan(vector_base< NumericT > &vec1, vector_base< NumericT > &vec2)
void inclusive_scan(vector_base< NumericT > &vec1, vector_base< NumericT > &vec2)
This function implements an inclusive scan.
Proxy classes for vectors.
void scaled_rank_1_update(matrix_base< NumericT > &mat1, ScalarT const &alpha, vcl_size_t, bool reciprocal_alpha, bool flip_sign_alpha, const vector_base< NumericT > &vec1, const vector_base< NumericT > &vec2)
The implementation of the operation mat += alpha * vec1 * vec2^T, i.e. a scaled rank 1 update...
All the predicates used within ViennaCL. Checks for expressions to be vectors, etc.
void matrix_row(const matrix_base< NumericT > &mat, unsigned int i, vector_base< NumericT > &vec)
viennacl::vector< NumericT > operator+(const vector_base< NumericT > &v1, const vector_expression< const matrix_base< NumericT >, const vector_base< NumericT >, op_prod > &proxy)
Implementation of the operation 'result = v1 + A * v2', where A is a matrix.
void am(matrix_base< NumericT > &mat1, matrix_base< NumericT > const &mat2, ScalarT1 const &alpha, vcl_size_t, bool reciprocal_alpha, bool flip_sign_alpha)
void matrix_column(const matrix_base< NumericT > &mat, unsigned int j, vector_base< NumericT > &vec)
viennacl::matrix_expression< const vector_base< NumericT >, const vector_base< NumericT >, op_prod > outer_prod(const vector_base< NumericT > &vec1, const vector_base< NumericT > &vec2)
Returns a proxy class for the operation mat += vec1 * vec2^T, i.e. a rank 1 update.
A tag class representing matrix-vector products and element-wise multiplications. ...
Definition: forwards.h:93
void matrix_diag_to_vector(const matrix_base< NumericT > &A, int k, vector_base< NumericT > &vec)
Implementations of dense matrix related operations, including matrix-vector products, using CUDA.
void element_op(matrix_base< NumericT, SizeT > &A, matrix_expression< const matrix_base< NumericT, SizeT >, const matrix_base< NumericT, SizeT >, op_element_binary< OpT > > const &proxy)
void inclusive_scan(vector_base< NumericT > &vec1, vector_base< NumericT > &vec2)
This function implements an inclusive scan.
void matrix_row(matrix_base< NumericT > const &mat, unsigned int i, vector_base< NumericT > &vec)
The vector type with operator-overloads and proxy classes is defined here. Linear algebra operations ...
bool row_major() const
Definition: matrix_def.hpp:239
void matrix_assign(matrix_base< NumericT > &A, NumericT s, bool up_to_internal_size=false)
void givens_next(matrix_base< NumericT > &Q, vector_base< NumericT > &tmp1, vector_base< NumericT > &tmp2, int l, int m)
This function updates the matrix Q. It is part of the tql2 algorithm.
void trans(const matrix_expression< const matrix_base< NumericT, SizeT, DistanceT >, const matrix_base< NumericT, SizeT, DistanceT >, op_trans > &proxy, matrix_base< NumericT > &temp_trans)
size_type size() const
Returns the length of the vector (cf. std::vector)
Definition: vector_def.hpp:118
float ScalarType
Definition: fft_1d.cpp:42
void bidiag_pack(matrix_base< NumericT > &A, VectorType &dh, VectorType &sh)
This function stores the diagonal and the superdiagonal of a matrix in two vectors.
void bidiag_pack(matrix_base< NumericT > &A, VectorType &dh, VectorType &sh)
Main abstraction class for multiple memory domains. Represents a buffer in either main RAM...
Definition: mem_handle.hpp:89
void givens_next(matrix_base< NumericT > &matrix, vector_base< NumericT > &tmp1, vector_base< NumericT > &tmp2, int l, int m)
A tag class representing transposed matrices.
Definition: forwards.h:219
size_type start2() const
Returns the number of columns.
Definition: matrix_def.hpp:221
void house_update_A_left(matrix_base< NumericT > &A, vector_base< NumericT > &D, vcl_size_t start)
This function applies a householder transformation to a matrix. A <- P * A with a householder reflect...
void exclusive_scan(vector_base< NumericT > &vec1, vector_base< NumericT > &vec2)
This function implements an exclusive scan.
void matrix_diagonal_assign(matrix_base< NumericT > &mat, NumericT s)
void matrix_row(const matrix_base< NumericT > &A, unsigned int i, vector_base< NumericT > &vec)
void givens_next(matrix_base< NumericT > &Q, vector_base< NumericT > &tmp1, vector_base< NumericT > &tmp2, int l, int m)
This function updates the matrix Q. It is part of the tql2 algorithm.
Extracts the underlying OpenCL handle from a vector, a matrix, an expression etc. ...
void prod_impl(const matrix_base< NumericT > &mat, const vector_base< NumericT > &vec, vector_base< NumericT > &result)
Carries out matrix-vector multiplication.
viennacl::backend::mem_handle & handle(T &obj)
Returns the generic memory handle of an object. Non-const version.
Definition: handle.hpp:41
void matrix_diag_from_vector(const vector_base< NumericT > &v, int k, matrix_base< NumericT > &A)
Dispatcher interface for A = diag(v, k)
void matrix_diag_from_vector(const vector_base< NumericT > &vec, int k, matrix_base< NumericT > &mat)
void ambm_m(matrix_base< NumericT > &mat1, matrix_base< NumericT > const &mat2, ScalarType1 const &alpha, vcl_size_t len_alpha, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT > const &mat3, ScalarType2 const &beta, vcl_size_t len_beta, bool reciprocal_beta, bool flip_sign_beta)
Implementation of the ViennaCL scalar class.
A collection of compile time type deductions.
void matrix_assign(matrix_base< NumericT > &mat, NumericT s, bool clear=false)
void prod_impl(const matrix_base< NumericT > &A, bool trans_A, const vector_base< NumericT > &vec, vector_base< NumericT > &result)
Carries out matrix-vector multiplication.
void ambm_m(matrix_base< NumericT > &A, matrix_base< NumericT > const &B, ScalarT1 const &alpha, vcl_size_t, bool reciprocal_alpha, bool flip_sign_alpha, matrix_base< NumericT > const &C, ScalarT2 const &beta, vcl_size_t, bool reciprocal_beta, bool flip_sign_beta)
void givens_next(matrix_base< NumericT > &Q, vector_base< NumericT > &tmp1, vector_base< NumericT > &tmp2, int l, int m)
This function updates the matrix Q. It is part of the tql2 algorithm.
void matrix_diagonal_assign(matrix_base< NumericT > &A, NumericT s)
void copy_vec(matrix_base< NumericT > &A, vector_base< NumericT > &V, vcl_size_t row_start, vcl_size_t col_start, bool copy_col)
This function copies a row or a column from a matrix to a vector.
Simple enable-if variant that uses the SFINAE pattern.
size_type start1() const
Returns the number of rows.
Definition: matrix_def.hpp:219
void exclusive_scan(vector_base< NumericT, F > &vec1, vector_base< NumericT, F > &vec2)
This function implements an exclusive scan.
memory_types get_active_handle_id() const
Returns an ID for the currently active memory buffer. Other memory buffers might contain old or no da...
Definition: mem_handle.hpp:118