1 #ifndef VIENNACL_LINALG_HOST_BASED_SPARSE_MATRIX_OPERATIONS_HPP_
2 #define VIENNACL_LINALG_HOST_BASED_SPARSE_MATRIX_OPERATIONS_HPP_
46 template<
typename NumericT,
unsigned int AlignmentV>
51 NumericT * result_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
52 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
53 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle1());
54 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle2());
59 unsigned int row_end = row_buffer[
row+1];
61 switch (info_selector)
64 for (
unsigned int i = row_buffer[
row]; i < row_end; ++i)
65 value = std::max<NumericT>(value, std::fabs(elements[i]));
69 for (
unsigned int i = row_buffer[
row]; i < row_end; ++i)
70 value += std::fabs(elements[i]);
74 for (
unsigned int i = row_buffer[
row]; i < row_end; ++i)
75 value += elements[i] * elements[i];
76 value = std::sqrt(value);
80 for (
unsigned int i = row_buffer[
row]; i < row_end; ++i)
82 if (col_buffer[i] ==
row)
90 result_buf[
row] = value;
104 template<
typename NumericT,
unsigned int AlignmentV>
109 NumericT * result_buf = detail::extract_raw_pointer<NumericT>(result.
handle());
110 NumericT
const * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
111 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
112 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle1());
113 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle2());
115 #ifdef VIENNACL_WITH_OPENMP
116 #pragma omp parallel for
118 for (
long row = 0; row < static_cast<long>(mat.
size1()); ++
row)
123 dot_prod += elements[i] * vec_buf[col_buffer[i] * vec.
stride() + vec.
start()];
137 template<
typename NumericT,
unsigned int AlignmentV>
142 NumericT
const * sp_mat_elements = detail::extract_raw_pointer<NumericT>(sp_mat.
handle());
143 unsigned int const * sp_mat_row_buffer = detail::extract_raw_pointer<unsigned int>(sp_mat.
handle1());
144 unsigned int const * sp_mat_col_buffer = detail::extract_raw_pointer<unsigned int>(sp_mat.
handle2());
146 NumericT
const * d_mat_data = detail::extract_raw_pointer<NumericT>(d_mat);
147 NumericT * result_data = detail::extract_raw_pointer<NumericT>(result);
164 d_mat_wrapper_row(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
166 d_mat_wrapper_col(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
169 result_wrapper_row(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
171 result_wrapper_col(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
174 #ifdef VIENNACL_WITH_OPENMP
175 #pragma omp parallel for
177 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row) {
182 for (
vcl_size_t k = row_start; k < row_end; ++k) {
183 temp += sp_mat_elements[k] * d_mat_wrapper_row(static_cast<vcl_size_t>(sp_mat_col_buffer[k]), col);
186 result_wrapper_row(
row, col) = temp;
188 result_wrapper_col(
row, col) = temp;
193 #ifdef VIENNACL_WITH_OPENMP
194 #pragma omp parallel for
196 for (
long col = 0; col < static_cast<long>(d_mat.
size2()); ++col) {
197 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row) {
201 for (
vcl_size_t k = row_start; k < row_end; ++k) {
202 temp += sp_mat_elements[k] * d_mat_wrapper_col(static_cast<vcl_size_t>(sp_mat_col_buffer[k]), static_cast<vcl_size_t>(col));
205 result_wrapper_row(
row, col) = temp;
207 result_wrapper_col(
row, col) = temp;
223 template<
typename NumericT,
unsigned int AlignmentV>
230 NumericT
const * sp_mat_elements = detail::extract_raw_pointer<NumericT>(sp_mat.
handle());
231 unsigned int const * sp_mat_row_buffer = detail::extract_raw_pointer<unsigned int>(sp_mat.
handle1());
232 unsigned int const * sp_mat_col_buffer = detail::extract_raw_pointer<unsigned int>(sp_mat.
handle2());
234 NumericT
const * d_mat_data = detail::extract_raw_pointer<NumericT>(d_mat.lhs());
235 NumericT * result_data = detail::extract_raw_pointer<NumericT>(result);
252 d_mat_wrapper_row(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
254 d_mat_wrapper_col(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
257 result_wrapper_row(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
259 result_wrapper_col(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
261 if ( d_mat.lhs().row_major() ) {
262 #ifdef VIENNACL_WITH_OPENMP
263 #pragma omp parallel for
265 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row) {
268 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col) {
270 for (
vcl_size_t k = row_start; k < row_end; ++k) {
271 temp += sp_mat_elements[k] * d_mat_wrapper_row(col, static_cast<vcl_size_t>(sp_mat_col_buffer[k]));
274 result_wrapper_row(
row, col) = temp;
276 result_wrapper_col(
row, col) = temp;
281 #ifdef VIENNACL_WITH_OPENMP
282 #pragma omp parallel for
284 for (
long col = 0; col < static_cast<long>(d_mat.size2()); ++col) {
289 for (
vcl_size_t k = row_start; k < row_end; ++k) {
290 temp += sp_mat_elements[k] * d_mat_wrapper_col(col, static_cast<vcl_size_t>(sp_mat_col_buffer[k]));
293 result_wrapper_row(
row, col) = temp;
295 result_wrapper_col(
row, col) = temp;
308 template<
typename NumericT,
typename ConstScalarArrayT,
typename ScalarArrayT,
typename IndexArrayT>
310 IndexArrayT
const & col_buffer,
311 ConstScalarArrayT
const & element_buffer,
312 ScalarArrayT & vec_buffer,
319 NumericT vec_entry = vec_buffer[
row];
321 for (
vcl_size_t i = row_begin; i < row_end; ++i)
325 vec_entry -= vec_buffer[col_index] * element_buffer[i];
327 vec_buffer[
row] = vec_entry;
332 template<
typename NumericT,
typename ConstScalarArrayT,
typename ScalarArrayT,
typename IndexArrayT>
334 IndexArrayT
const & col_buffer,
335 ConstScalarArrayT
const & element_buffer,
336 ScalarArrayT & vec_buffer,
343 NumericT vec_entry = vec_buffer[
row];
347 NumericT diagonal_entry = 0;
348 for (
vcl_size_t i = row_begin; i < row_end; ++i)
352 vec_entry -= vec_buffer[col_index] * element_buffer[i];
353 else if (col_index ==
row)
354 diagonal_entry = element_buffer[i];
357 vec_buffer[
row] = vec_entry / diagonal_entry;
363 template<
typename NumericT,
typename ConstScalarArrayT,
typename ScalarArrayT,
typename IndexArrayT>
365 IndexArrayT
const & col_buffer,
366 ConstScalarArrayT
const & element_buffer,
367 ScalarArrayT & vec_buffer,
371 for (
vcl_size_t row2 = 1; row2 < num_cols; ++row2)
374 NumericT vec_entry = vec_buffer[
row];
377 for (
vcl_size_t i = row_begin; i < row_end; ++i)
381 vec_entry -= vec_buffer[col_index] * element_buffer[i];
383 vec_buffer[
row] = vec_entry;
387 template<
typename NumericT,
typename ConstScalarArrayT,
typename ScalarArrayT,
typename IndexArrayT>
389 IndexArrayT
const & col_buffer,
390 ConstScalarArrayT
const & element_buffer,
391 ScalarArrayT & vec_buffer,
395 for (
vcl_size_t row2 = 0; row2 < num_cols; ++row2)
398 NumericT vec_entry = vec_buffer[
row];
403 NumericT diagonal_entry = 0;
404 for (
vcl_size_t i = row_begin; i < row_end; ++i)
408 vec_entry -= vec_buffer[col_index] * element_buffer[i];
409 else if (col_index == row)
410 diagonal_entry = element_buffer[i];
413 vec_buffer[
row] = vec_entry / diagonal_entry;
427 template<
typename NumericT,
unsigned int AlignmentV>
432 NumericT * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
433 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(L.
handle());
434 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(L.
handle1());
435 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(L.
handle2());
437 detail::csr_inplace_solve<NumericT>(row_buffer, col_buffer, elements, vec_buf, L.
size2(), tag);
446 template<
typename NumericT,
unsigned int AlignmentV>
451 NumericT * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
452 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(L.
handle());
453 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(L.
handle1());
454 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(L.
handle2());
456 detail::csr_inplace_solve<NumericT>(row_buffer, col_buffer, elements, vec_buf, L.
size2(), tag);
466 template<
typename NumericT,
unsigned int AlignmentV>
471 NumericT * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
472 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(U.
handle());
473 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(U.
handle1());
474 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(U.
handle2());
476 detail::csr_inplace_solve<NumericT>(row_buffer, col_buffer, elements, vec_buf, U.
size2(), tag);
485 template<
typename NumericT,
unsigned int AlignmentV>
490 NumericT * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
491 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(U.
handle());
492 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(U.
handle1());
493 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(U.
handle2());
495 detail::csr_inplace_solve<NumericT>(row_buffer, col_buffer, elements, vec_buf, U.
size2(), tag);
510 template<
typename NumericT,
typename ConstScalarArrayT,
typename ScalarArrayT,
typename IndexArrayT>
512 IndexArrayT
const & col_buffer,
513 ConstScalarArrayT
const & element_buffer,
514 ScalarArrayT & vec_buffer,
519 for (
vcl_size_t col = 0; col < num_cols; ++col)
521 NumericT vec_entry = vec_buffer[col];
523 for (
vcl_size_t i = col_begin; i < col_end; ++i)
525 unsigned int row_index = col_buffer[i];
527 vec_buffer[row_index] -= vec_entry * element_buffer[i];
533 template<
typename NumericT,
typename ConstScalarArrayT,
typename ScalarArrayT,
typename IndexArrayT>
535 IndexArrayT
const & col_buffer,
536 ConstScalarArrayT
const & element_buffer,
537 ScalarArrayT & vec_buffer,
542 for (
vcl_size_t col = 0; col < num_cols; ++col)
547 NumericT diagonal_entry = 0;
548 for (
vcl_size_t i = col_begin; i < col_end; ++i)
551 if (row_index == col)
553 diagonal_entry = element_buffer[i];
559 NumericT vec_entry = vec_buffer[col] / diagonal_entry;
560 vec_buffer[col] = vec_entry;
561 for (
vcl_size_t i = col_begin; i < col_end; ++i)
565 vec_buffer[row_index] -= vec_entry * element_buffer[i];
571 template<
typename NumericT,
typename ConstScalarArrayT,
typename ScalarArrayT,
typename IndexArrayT>
573 IndexArrayT
const & col_buffer,
574 ConstScalarArrayT
const & element_buffer,
575 ScalarArrayT & vec_buffer,
579 for (
vcl_size_t col2 = 0; col2 < num_cols; ++col2)
583 NumericT vec_entry = vec_buffer[col];
586 for (
vcl_size_t i = col_begin; i < col_end; ++i)
590 vec_buffer[row_index] -= vec_entry * element_buffer[i];
596 template<
typename NumericT,
typename ConstScalarArrayT,
typename ScalarArrayT,
typename IndexArrayT>
598 IndexArrayT
const & col_buffer,
599 ConstScalarArrayT
const & element_buffer,
600 ScalarArrayT & vec_buffer,
604 for (
vcl_size_t col2 = 0; col2 < num_cols; ++col2)
611 NumericT diagonal_entry = 0;
612 for (
vcl_size_t i = col_begin; i < col_end; ++i)
615 if (row_index == col)
617 diagonal_entry = element_buffer[i];
623 NumericT vec_entry = vec_buffer[col] / diagonal_entry;
624 vec_buffer[col] = vec_entry;
625 for (
vcl_size_t i = col_begin; i < col_end; ++i)
629 vec_buffer[row_index] -= vec_entry * element_buffer[i];
638 template<
typename NumericT,
unsigned int AlignmentV>
649 unsigned int const * row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(L.lhs().handle1());
650 unsigned int const * col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(L.lhs().handle2());
651 NumericT
const * elements = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(L.lhs().handle());
652 NumericT * vec_buffer = detail::extract_raw_pointer<NumericT>(vec.
handle());
655 for (
vcl_size_t col = 0; col < L.lhs().size1(); ++col)
657 NumericT vec_entry = vec_buffer[col];
659 for (
vcl_size_t i = col_begin; i < col_end; ++i)
661 unsigned int row_index = col_buffer[i];
663 vec_buffer[row_index] -= vec_entry * elements[i];
669 template<
typename NumericT,
unsigned int AlignmentV>
680 unsigned int const * row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(L.lhs().handle1());
681 unsigned int const * col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(L.lhs().handle2());
682 NumericT
const * elements = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(L.lhs().handle());
683 NumericT
const * diagonal_buffer = detail::extract_raw_pointer<NumericT>(L_diagonal.
handle());
684 NumericT * vec_buffer = detail::extract_raw_pointer<NumericT>(vec.
handle());
687 for (
vcl_size_t col = 0; col < L.lhs().size1(); ++col)
691 NumericT vec_entry = vec_buffer[col] / diagonal_buffer[col];
692 vec_buffer[col] = vec_entry;
693 for (
vcl_size_t i = col_begin; i < col_end; ++i)
697 vec_buffer[row_index] -= vec_entry * elements[i];
705 template<
typename NumericT,
unsigned int AlignmentV>
716 unsigned int const * row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(U.lhs().handle1());
717 unsigned int const * col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(U.lhs().handle2());
718 NumericT
const * elements = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(U.lhs().handle());
719 NumericT * vec_buffer = detail::extract_raw_pointer<NumericT>(vec.
handle());
721 for (
vcl_size_t col2 = 0; col2 < U.lhs().size1(); ++col2)
723 vcl_size_t col = (U.lhs().size1() - col2) - 1;
725 NumericT vec_entry = vec_buffer[col];
728 for (
vcl_size_t i = col_begin; i < col_end; ++i)
732 vec_buffer[row_index] -= vec_entry * elements[i];
738 template<
typename NumericT,
unsigned int AlignmentV>
749 unsigned int const * row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(U.lhs().handle1());
750 unsigned int const * col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(U.lhs().handle2());
751 NumericT
const * elements = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(U.lhs().handle());
752 NumericT
const * diagonal_buffer = detail::extract_raw_pointer<NumericT>(U_diagonal.
handle());
753 NumericT * vec_buffer = detail::extract_raw_pointer<NumericT>(vec.
handle());
755 for (
vcl_size_t col2 = 0; col2 < U.lhs().size1(); ++col2)
757 vcl_size_t col = (U.lhs().size1() - col2) - 1;
762 NumericT vec_entry = vec_buffer[col] / diagonal_buffer[col];
763 vec_buffer[col] = vec_entry;
764 for (
vcl_size_t i = col_begin; i < col_end; ++i)
768 vec_buffer[row_index] -= vec_entry * elements[i];
782 template<
typename NumericT,
unsigned int AlignmentV>
789 NumericT * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
790 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(proxy.lhs().handle());
791 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(proxy.lhs().handle1());
792 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(proxy.lhs().handle2());
794 detail::csr_trans_inplace_solve<NumericT>(row_buffer, col_buffer, elements, vec_buf, proxy.lhs().size1(), tag);
803 template<
typename NumericT,
unsigned int AlignmentV>
810 NumericT * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
811 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(proxy.lhs().handle());
812 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(proxy.lhs().handle1());
813 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(proxy.lhs().handle2());
815 detail::csr_trans_inplace_solve<NumericT>(row_buffer, col_buffer, elements, vec_buf, proxy.lhs().size1(), tag);
825 template<
typename NumericT,
unsigned int AlignmentV>
832 NumericT * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
833 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(proxy.lhs().handle());
834 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(proxy.lhs().handle1());
835 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(proxy.lhs().handle2());
837 detail::csr_trans_inplace_solve<NumericT>(row_buffer, col_buffer, elements, vec_buf, proxy.lhs().size1(), tag);
847 template<
typename NumericT,
unsigned int AlignmentV>
854 NumericT * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
855 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(proxy.lhs().handle());
856 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(proxy.lhs().handle1());
857 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(proxy.lhs().handle2());
859 detail::csr_trans_inplace_solve<NumericT>(row_buffer, col_buffer, elements, vec_buf, proxy.lhs().size1(), tag);
876 template<
typename NumericT>
881 NumericT * result_buf = detail::extract_raw_pointer<NumericT>(result.
handle());
882 NumericT
const * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
883 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
884 unsigned int const * row_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle1());
885 unsigned int const * row_indices = detail::extract_raw_pointer<unsigned int>(mat.
handle3());
886 unsigned int const * col_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle2());
890 #ifdef VIENNACL_WITH_OPENMP
891 #pragma omp parallel for
893 for (
long i = 0; i < static_cast<long>(mat.
nnz1()); ++i)
897 for (
vcl_size_t j = row_buffer[i]; j < row_end; ++j)
898 dot_prod += elements[j] * vec_buf[col_buffer[j] * vec.
stride() + vec.
start()];
912 template<
typename NumericT,
unsigned int AlignmentV>
917 NumericT * result_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
918 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
919 unsigned int const * coord_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle12());
922 unsigned int last_row = 0;
926 unsigned int current_row = coord_buffer[2*i];
928 if (current_row != last_row)
931 value = std::sqrt(value);
933 result_buf[last_row] = value;
935 last_row = current_row;
938 switch (info_selector)
941 value = std::max<NumericT>(value, std::fabs(elements[i]));
945 value += std::fabs(elements[i]);
949 value += elements[i] * elements[i];
953 if (coord_buffer[2*i+1] == current_row)
963 value = std::sqrt(value);
965 result_buf[last_row] = value;
977 template<
typename NumericT,
unsigned int AlignmentV>
982 NumericT * result_buf = detail::extract_raw_pointer<NumericT>(result.
handle());
983 NumericT
const * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
984 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
985 unsigned int const * coord_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle12());
988 result_buf[i * result.
stride() + result.
start()] = 0;
991 result_buf[coord_buffer[2*i] * result.
stride() + result.
start()]
992 += elements[i] * vec_buf[coord_buffer[2*i+1] * vec.
stride() + vec.
start()];
1003 template<
typename NumericT,
unsigned int AlignmentV>
1008 NumericT
const * sp_mat_elements = detail::extract_raw_pointer<NumericT>(sp_mat.
handle());
1009 unsigned int const * sp_mat_coords = detail::extract_raw_pointer<unsigned int>(sp_mat.
handle12());
1011 NumericT
const * d_mat_data = detail::extract_raw_pointer<NumericT>(d_mat);
1012 NumericT * result_data = detail::extract_raw_pointer<NumericT>(result);
1029 d_mat_wrapper_row(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1031 d_mat_wrapper_col(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1034 result_wrapper_row(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1036 result_wrapper_col(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1040 #ifdef VIENNACL_WITH_OPENMP
1041 #pragma omp parallel for
1043 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row)
1047 result_wrapper_row(
row, col) = (NumericT)0;
1050 result_wrapper_col(
row, col) = (NumericT)0;
1053 #ifdef VIENNACL_WITH_OPENMP
1054 #pragma omp parallel for
1056 for (
long i = 0; i < static_cast<long>(sp_mat.
nnz()); ++i) {
1057 NumericT x =
static_cast<NumericT
>(sp_mat_elements[i]);
1061 NumericT y = d_mat_wrapper_row( c, col);
1063 result_wrapper_row(r, col) += x * y;
1065 result_wrapper_col(r, col) += x * y;
1072 #ifdef VIENNACL_WITH_OPENMP
1073 #pragma omp parallel for
1075 for (
long col = 0; col < static_cast<long>(d_mat.
size2()); ++col)
1079 result_wrapper_row(
row, col) = (NumericT)0;
1082 result_wrapper_col(
row, col) = (NumericT)0;
1085 #ifdef VIENNACL_WITH_OPENMP
1086 #pragma omp parallel for
1088 for (
long col = 0; col < static_cast<long>(d_mat.
size2()); ++col) {
1092 NumericT x =
static_cast<NumericT
>(sp_mat_elements[i]);
1095 NumericT y = d_mat_wrapper_col( c, col);
1098 result_wrapper_row( r, col) += x*y;
1100 result_wrapper_col( r, col) += x*y;
1117 template<
typename NumericT,
unsigned int AlignmentV>
1124 NumericT
const * sp_mat_elements = detail::extract_raw_pointer<NumericT>(sp_mat.
handle());
1125 unsigned int const * sp_mat_coords = detail::extract_raw_pointer<unsigned int>(sp_mat.
handle12());
1127 NumericT
const * d_mat_data = detail::extract_raw_pointer<NumericT>(d_mat.lhs());
1128 NumericT * result_data = detail::extract_raw_pointer<NumericT>(result);
1145 d_mat_wrapper_row(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1147 d_mat_wrapper_col(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1150 result_wrapper_row(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1152 result_wrapper_col(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1154 if ( d_mat.lhs().row_major() )
1156 #ifdef VIENNACL_WITH_OPENMP
1157 #pragma omp parallel for
1159 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row)
1162 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col)
1163 result_wrapper_row(
row, col) = (NumericT)0;
1165 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col)
1166 result_wrapper_col(
row, col) = (NumericT)0;
1169 #ifdef VIENNACL_WITH_OPENMP
1170 #pragma omp parallel for
1172 for (
long i = 0; i < static_cast<long>(sp_mat.
nnz()); ++i) {
1173 NumericT x =
static_cast<NumericT
>(sp_mat_elements[i]);
1178 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col) {
1179 NumericT y = d_mat_wrapper_row( col, c);
1180 result_wrapper_row(r, col) += x * y;
1185 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col) {
1186 NumericT y = d_mat_wrapper_row( col, c);
1187 result_wrapper_col(r, col) += x * y;
1196 #ifdef VIENNACL_WITH_OPENMP
1197 #pragma omp parallel for
1199 for (
long col = 0; col < static_cast<long>(d_mat.size2()); ++col)
1203 result_wrapper_row(
row, col) = (NumericT)0;
1206 result_wrapper_col(
row, col) = (NumericT)0;
1209 #ifdef VIENNACL_WITH_OPENMP
1210 #pragma omp parallel for
1212 for (
long i = 0; i < static_cast<long>(sp_mat.
nnz()); ++i) {
1213 NumericT x =
static_cast<NumericT
>(sp_mat_elements[i]);
1218 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col) {
1219 NumericT y = d_mat_wrapper_col( col, c);
1220 result_wrapper_row(r, col) += x * y;
1225 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col) {
1226 NumericT y = d_mat_wrapper_col( col, c);
1227 result_wrapper_col(r, col) += x * y;
1245 template<
typename NumericT,
unsigned int AlignmentV>
1250 NumericT * result_buf = detail::extract_raw_pointer<NumericT>(result.
handle());
1251 NumericT
const * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
1252 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
1253 unsigned int const * coords = detail::extract_raw_pointer<unsigned int>(mat.
handle2());
1259 for (
unsigned int item_id = 0; item_id < mat.
internal_maxnnz(); ++item_id)
1262 NumericT val = elements[offset];
1264 if (val > 0 || val < 0)
1266 unsigned int col = coords[offset];
1267 sum += (vec_buf[col * vec.
stride() + vec.
start()] * val);
1283 template<
typename NumericT,
unsigned int AlignmentV>
1288 NumericT
const * sp_mat_elements = detail::extract_raw_pointer<NumericT>(sp_mat.
handle());
1289 unsigned int const * sp_mat_coords = detail::extract_raw_pointer<unsigned int>(sp_mat.
handle2());
1291 NumericT
const * d_mat_data = detail::extract_raw_pointer<NumericT>(d_mat);
1292 NumericT * result_data = detail::extract_raw_pointer<NumericT>(result);
1309 d_mat_wrapper_row(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1311 d_mat_wrapper_col(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1314 result_wrapper_row(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1316 result_wrapper_col(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1319 #ifdef VIENNACL_WITH_OPENMP
1320 #pragma omp parallel for
1322 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row)
1326 result_wrapper_row(
row, col) = (NumericT)0;
1329 result_wrapper_col(
row, col) = (NumericT)0;
1332 #ifdef VIENNACL_WITH_OPENMP
1333 #pragma omp parallel for
1335 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row)
1337 for (
long item_id = 0; item_id < static_cast<long>(sp_mat.
maxnnz()); ++item_id)
1340 NumericT sp_mat_val =
static_cast<NumericT
>(sp_mat_elements[offset]);
1343 if (sp_mat_val < 0 || sp_mat_val > 0)
1347 result_wrapper_row(static_cast<vcl_size_t>(
row), col) += sp_mat_val * d_mat_wrapper_row( sp_mat_col, col);
1350 result_wrapper_col(static_cast<vcl_size_t>(
row), col) += sp_mat_val * d_mat_wrapper_row( sp_mat_col, col);
1356 #ifdef VIENNACL_WITH_OPENMP
1357 #pragma omp parallel for
1359 for (
long col = 0; col < static_cast<long>(d_mat.
size2()); ++col)
1362 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row)
1363 result_wrapper_row(
row, col) = (NumericT)0;
1365 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row)
1366 result_wrapper_col(
row, col) = (NumericT)0;
1369 #ifdef VIENNACL_WITH_OPENMP
1370 #pragma omp parallel for
1372 for (
long col = 0; col < static_cast<long>(d_mat.
size2()); ++col) {
1374 for (
unsigned int item_id = 0; item_id < sp_mat.
maxnnz(); ++item_id) {
1379 NumericT sp_mat_val =
static_cast<NumericT
>(sp_mat_elements[offset]);
1382 if (sp_mat_val < 0 || sp_mat_val > 0)
1385 result_wrapper_row(
row, col) += sp_mat_val * d_mat_wrapper_col( sp_mat_col, col);
1387 result_wrapper_col(
row, col) += sp_mat_val * d_mat_wrapper_col( sp_mat_col, col);
1405 template<
typename NumericT,
unsigned int AlignmentV>
1412 NumericT
const * sp_mat_elements = detail::extract_raw_pointer<NumericT>(sp_mat.
handle());
1413 unsigned int const * sp_mat_coords = detail::extract_raw_pointer<unsigned int>(sp_mat.
handle2());
1415 NumericT
const * d_mat_data = detail::extract_raw_pointer<NumericT>(d_mat.lhs());
1416 NumericT * result_data = detail::extract_raw_pointer<NumericT>(result);
1433 d_mat_wrapper_row(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1435 d_mat_wrapper_col(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1438 result_wrapper_row(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1440 result_wrapper_col(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1442 if ( d_mat.lhs().row_major() )
1444 #ifdef VIENNACL_WITH_OPENMP
1445 #pragma omp parallel for
1447 for (
long row = 0; row < static_cast<long>(sp_mat.
size1()); ++
row)
1450 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col)
1451 result_wrapper_row(
row, col) = (NumericT)0;
1453 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col)
1454 result_wrapper_col(
row, col) = (NumericT)0;
1457 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col) {
1459 for (
unsigned int item_id = 0; item_id < sp_mat.
maxnnz(); ++item_id) {
1464 NumericT sp_mat_val =
static_cast<NumericT
>(sp_mat_elements[offset]);
1467 if (sp_mat_val < 0 || sp_mat_val > 0)
1470 result_wrapper_row(
row, col) += sp_mat_val * d_mat_wrapper_row( col, sp_mat_col);
1472 result_wrapper_col(
row, col) += sp_mat_val * d_mat_wrapper_row( col, sp_mat_col);
1480 #ifdef VIENNACL_WITH_OPENMP
1481 #pragma omp parallel for
1483 for (
long col = 0; col < static_cast<long>(d_mat.size2()); ++col)
1487 result_wrapper_row(
row, col) = (NumericT)0;
1490 result_wrapper_col(
row, col) = (NumericT)0;
1493 #ifdef VIENNACL_WITH_OPENMP
1494 #pragma omp parallel for
1496 for (
long item_id = 0; item_id < static_cast<long>(sp_mat.
maxnnz()); ++item_id) {
1501 NumericT sp_mat_val =
static_cast<NumericT
>(sp_mat_elements[offset]);
1504 if (sp_mat_val < 0 || sp_mat_val > 0)
1507 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col)
1508 result_wrapper_row(
row, col) += sp_mat_val * d_mat_wrapper_col( col, sp_mat_col);
1510 for (
vcl_size_t col = 0; col < d_mat.size2(); ++col)
1511 result_wrapper_col(
row, col) += sp_mat_val * d_mat_wrapper_col( col, sp_mat_col);
1531 template<
typename NumericT,
typename IndexT>
1536 NumericT * result_buf = detail::extract_raw_pointer<NumericT>(result.
handle());
1537 NumericT
const * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
1538 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
1539 IndexT
const * columns_per_block = detail::extract_raw_pointer<IndexT>(mat.
handle1());
1540 IndexT
const * column_indices = detail::extract_raw_pointer<IndexT>(mat.
handle2());
1541 IndexT
const * block_start = detail::extract_raw_pointer<IndexT>(mat.
handle3());
1545 #ifdef VIENNACL_WITH_OPENMP
1546 #pragma omp parallel for
1548 for (
long block_idx2 = 0; block_idx2 < static_cast<long>(num_blocks); ++block_idx2)
1551 vcl_size_t current_columns_per_block = columns_per_block[block_idx];
1555 for (IndexT column_entry_index = 0;
1556 column_entry_index < current_columns_per_block;
1557 ++column_entry_index)
1562 for (IndexT row_in_block = 0; row_in_block < mat.
rows_per_block(); ++row_in_block)
1564 NumericT val = elements[stride_start + row_in_block];
1566 result_values[row_in_block] += (val > 0 || val < 0) ? vec_buf[column_indices[stride_start + row_in_block] * vec.
stride() + vec.
start()] * val : 0;
1571 for (IndexT row_in_block = 0; row_in_block < mat.
rows_per_block(); ++row_in_block)
1573 if (first_row_in_matrix + row_in_block < result.
size())
1574 result_buf[(first_row_in_matrix + row_in_block) * result.
stride() + result.
start()] = result_values[row_in_block];
1591 template<
typename NumericT,
unsigned int AlignmentV>
1596 NumericT * result_buf = detail::extract_raw_pointer<NumericT>(result.
handle());
1597 NumericT
const * vec_buf = detail::extract_raw_pointer<NumericT>(vec.
handle());
1598 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
1599 unsigned int const * coords = detail::extract_raw_pointer<unsigned int>(mat.
handle2());
1600 NumericT
const * csr_elements = detail::extract_raw_pointer<NumericT>(mat.
handle5());
1601 unsigned int const * csr_row_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle3());
1602 unsigned int const * csr_col_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle4());
1612 for (
unsigned int item_id = 0; item_id < mat.
internal_ellnnz(); ++item_id)
1615 NumericT val = elements[offset];
1617 if (val > 0 || val < 0)
1619 unsigned int col = coords[offset];
1620 sum += (vec_buf[col * vec.
stride() + vec.
start()] * val);
1630 for (
vcl_size_t item_id = col_begin; item_id < col_end; item_id++)
1632 sum += (vec_buf[csr_col_buffer[item_id] * vec.
stride() + vec.
start()] * csr_elements[item_id]);
1651 template<
typename NumericT,
unsigned int AlignmentV>
1656 NumericT
const * d_mat_data = detail::extract_raw_pointer<NumericT>(d_mat);
1657 NumericT * result_data = detail::extract_raw_pointer<NumericT>(result);
1674 d_mat_wrapper_row(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1676 d_mat_wrapper_col(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1679 result_wrapper_row(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1681 result_wrapper_col(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1683 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
1684 unsigned int const * coords = detail::extract_raw_pointer<unsigned int>(mat.
handle2());
1685 NumericT
const * csr_elements = detail::extract_raw_pointer<NumericT>(mat.
handle5());
1686 unsigned int const * csr_row_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle3());
1687 unsigned int const * csr_col_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle4());
1690 for (
vcl_size_t result_col = 0; result_col < result.
size2(); ++result_col)
1699 for (
unsigned int item_id = 0; item_id < mat.
internal_ellnnz(); ++item_id)
1702 NumericT val = elements[offset];
1704 if (val < 0 || val > 0)
1708 sum += d_mat_wrapper_row(col, result_col) * val;
1710 sum += d_mat_wrapper_col(col, result_col) * val;
1721 for (
vcl_size_t item_id = col_begin; item_id < col_end; item_id++)
1722 sum += d_mat_wrapper_row(static_cast<vcl_size_t>(csr_col_buffer[item_id]), result_col) * csr_elements[item_id];
1724 for (
vcl_size_t item_id = col_begin; item_id < col_end; item_id++)
1725 sum += d_mat_wrapper_col(static_cast<vcl_size_t>(csr_col_buffer[item_id]), result_col) * csr_elements[item_id];
1728 result_wrapper_row(
row, result_col) =
sum;
1730 result_wrapper_col(
row, result_col) =
sum;
1744 template<
typename NumericT,
unsigned int AlignmentV>
1751 NumericT
const * d_mat_data = detail::extract_raw_pointer<NumericT>(d_mat);
1752 NumericT * result_data = detail::extract_raw_pointer<NumericT>(result);
1769 d_mat_wrapper_row(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1771 d_mat_wrapper_col(d_mat_data, d_mat_start1, d_mat_start2, d_mat_inc1, d_mat_inc2, d_mat_internal_size1, d_mat_internal_size2);
1774 result_wrapper_row(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1776 result_wrapper_col(result_data, result_start1, result_start2, result_inc1, result_inc2, result_internal_size1, result_internal_size2);
1778 NumericT
const * elements = detail::extract_raw_pointer<NumericT>(mat.
handle());
1779 unsigned int const * coords = detail::extract_raw_pointer<unsigned int>(mat.
handle2());
1780 NumericT
const * csr_elements = detail::extract_raw_pointer<NumericT>(mat.
handle5());
1781 unsigned int const * csr_row_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle3());
1782 unsigned int const * csr_col_buffer = detail::extract_raw_pointer<unsigned int>(mat.
handle4());
1785 for (
vcl_size_t result_col = 0; result_col < result.
size2(); ++result_col)
1794 for (
unsigned int item_id = 0; item_id < mat.
internal_ellnnz(); ++item_id)
1797 NumericT val = elements[offset];
1799 if (val < 0 || val > 0)
1802 if (d_mat.lhs().row_major())
1803 sum += d_mat_wrapper_row(result_col, col) * val;
1805 sum += d_mat_wrapper_col(result_col, col) * val;
1815 if (d_mat.lhs().row_major())
1816 for (
vcl_size_t item_id = col_begin; item_id < col_end; item_id++)
1817 sum += d_mat_wrapper_row(result_col, static_cast<vcl_size_t>(csr_col_buffer[item_id])) * csr_elements[item_id];
1819 for (
vcl_size_t item_id = col_begin; item_id < col_end; item_id++)
1820 sum += d_mat_wrapper_col(result_col, static_cast<vcl_size_t>(csr_col_buffer[item_id])) * csr_elements[item_id];
1823 result_wrapper_row(
row, result_col) =
sum;
1825 result_wrapper_col(
row, result_col) =
sum;
const vcl_size_t & size2() const
Returns the number of columns.
vcl_size_t internal_ellnnz() const
Sparse matrix class using a hybrid format composed of the ELL and CSR format for storing the nonzeros...
result_of::size_type< matrix_base< NumericT > >::type stride1(matrix_base< NumericT > const &s)
const handle_type & handle3() const
const vcl_size_t & size1() const
Returns the number of rows.
const handle_type & handle2() const
Returns the OpenCL handle to the column index array.
const handle_type & handle1() const
Returns the OpenCL handle to the row index array.
const handle_type & handle() const
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...
const handle_type & handle12() const
Returns the OpenCL handle to the (row, column) index array.
A tag class representing a lower triangular matrix.
vcl_size_t internal_size1() const
void inplace_solve(matrix_base< NumericT > const &A, matrix_base< NumericT > &B, SolverTagT)
Direct inplace solver for triangular systems with multiple right hand sides, i.e. A \ B (MATLAB notat...
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...
Expression template class for representing a tree of expressions which ultimately result in a matrix...
This file provides the forward declarations for the main types used within ViennaCL.
result_of::size_type< T >::type start1(T const &obj)
vcl_size_t nnz() const
Returns the number of nonzero entries.
void dot_prod(MatrixT const &A, unsigned int beg_ind, NumericT &res)
Dot prod of particular column of martix A with it's self starting at a certain index beg_ind...
const handle_type & handle4() const
vcl_size_t rows_per_block() const
void vector_assign(vector_base< NumericT > &vec1, const NumericT &alpha, bool up_to_internal_size=false)
Assign a constant value to a vector (-range/-slice)
statement sum(scalar< NumericT > const *s, vector_base< NumericT > const *x)
const handle_type & handle() const
Returns the OpenCL handle to the matrix entry array.
const handle_type & handle1() const
Returns the OpenCL handle to the row index array.
vcl_size_t internal_size1() const
size_type stride() const
Returns the stride within the buffer (in multiples of sizeof(NumericT))
const handle_type & handle2() const
const handle_type & handle() const
Returns the OpenCL handle to the matrix entry array.
result_of::size_type< T >::type start2(T const &obj)
Helper array for accessing a strided submatrix embedded in a larger matrix.
Sparse matrix class using the ELLPACK format for storing the nonzeros.
A tag class representing an upper triangular matrix.
void prod_impl(const matrix_base< NumericT > &mat, bool trans, const vector_base< NumericT > &vec, vector_base< NumericT > &result)
Carries out matrix-vector multiplication.
Sparse matrix class using the sliced ELLPACK with parameters C, .
const handle_type & handle3() const
Returns the OpenCL handle to the row index array.
A sparse square matrix in compressed sparse rows format optimized for the case that only a few rows c...
const handle_type & handle2() const
Returns the OpenCL handle to the column index array.
void block_inplace_solve(const matrix_expression< const compressed_matrix< NumericT, AlignmentV >, const compressed_matrix< NumericT, AlignmentV >, op_trans > &L, viennacl::backend::mem_handle const &, vcl_size_t, vector_base< NumericT > const &, vector_base< NumericT > &vec, viennacl::linalg::unit_lower_tag)
size_type size2() const
Returns the number of columns.
Common routines for single-threaded or OpenMP-enabled execution on CPU.
vcl_size_t maxnnz() const
result_of::size_type< matrix_base< NumericT > >::type stride2(matrix_base< NumericT > const &s)
vector_expression< const matrix_base< NumericT, F >, const unsigned int, op_row > row(const matrix_base< NumericT, F > &A, unsigned int i)
const handle_type & handle() const
Returns the OpenCL handle to the matrix entry array.
The vector type with operator-overloads and proxy classes is defined here. Linear algebra operations ...
size_type size() const
Returns the length of the vector (cf. std::vector)
const vcl_size_t & nnz1() const
Returns the number of nonzero entries.
void csr_inplace_solve(IndexArrayT const &row_buffer, IndexArrayT const &col_buffer, ConstScalarArrayT const &element_buffer, ScalarArrayT &vec_buffer, vcl_size_t num_cols, viennacl::linalg::unit_lower_tag)
A tag class representing a lower triangular matrix with unit diagonal.
Main abstraction class for multiple memory domains. Represents a buffer in either main RAM...
A tag class representing transposed matrices.
void row_info(compressed_matrix< NumericT, AlignmentV > const &mat, vector_base< NumericT > &vec, viennacl::linalg::detail::row_info_types info_selector)
A sparse square matrix in compressed sparse rows format.
const handle_type & handle5() const
size_type start() const
Returns the offset within the buffer.
vcl_size_t size1() const
Returns the number of rows.
vcl_size_t internal_maxnnz() const
Implementation of the ViennaCL scalar class.
const handle_type & handle() const
Returns the memory handle.
void csr_trans_inplace_solve(IndexArrayT const &row_buffer, IndexArrayT const &col_buffer, ConstScalarArrayT const &element_buffer, ScalarArrayT &vec_buffer, vcl_size_t num_cols, viennacl::linalg::unit_lower_tag)
A tag class representing an upper triangular matrix with unit diagonal.
A sparse square matrix, where entries are stored as triplets (i,j, val), where i and j are the row an...
Implementations of NMF operations using a plain single-threaded or OpenMP-enabled execution on CPU...