1 #ifndef VIENNACL_LINALG_HOST_BASED_AMG_OPERATIONS_HPP 
    2 #define VIENNACL_LINALG_HOST_BASED_AMG_OPERATIONS_HPP 
   32 #ifdef VIENNACL_WITH_OPENMP 
   49 template<
typename NumericT>
 
   56   unsigned int const * A_row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle1());
 
   57   unsigned int const * A_col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle2());
 
   59   unsigned int *influences_row_ptr = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_jumper_.
handle());
 
   60   unsigned int *influences_id_ptr  = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_ids_.
handle());
 
   61   unsigned int *influences_values_ptr  = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_values_.
handle());
 
   63 #ifdef VIENNACL_WITH_OPENMP 
   64   #pragma omp parallel for 
   66   for (
long i2=0; i2<static_cast<long>(A.
size1()); ++i2)
 
   69     influences_row_ptr[i] = A_row_buffer[i];
 
   70     influences_values_ptr[i] = A_row_buffer[i+1] - A_row_buffer[i];
 
   72   influences_row_ptr[A.
size1()] = A_row_buffer[A.
size1()];
 
   74 #ifdef VIENNACL_WITH_OPENMP 
   75   #pragma omp parallel for 
   77   for (
long i=0; i<long(A.
nnz()); ++i)
 
   78     influences_id_ptr[i] = A_col_buffer[i];
 
   83 template<
typename NumericT>
 
   88   NumericT     const * A_elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(A.
handle());
 
   89   unsigned int const * A_row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle1());
 
   90   unsigned int const * A_col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle2());
 
   92   unsigned int *influences_row_ptr = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_jumper_.
handle());
 
   93   unsigned int *influences_id_ptr  = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_ids_.
handle());
 
   98 #ifdef VIENNACL_WITH_OPENMP 
   99   #pragma omp parallel for 
  101   for (
long i2=0; i2<static_cast<long>(A.
size1()); ++i2)
 
  104     unsigned int row_start = A_row_buffer[i];
 
  105     unsigned int row_stop  = A_row_buffer[i+1];
 
  109     unsigned int num_influences = 0;
 
  112     for (
unsigned int nnz_index = row_start; nnz_index < row_stop; ++nnz_index)
 
  114       unsigned int col = A_col_buffer[nnz_index];
 
  115       NumericT value   = A_elements[nnz_index];
 
  119       else if (value > largest_positive)
 
  120         largest_positive = value;
 
  121       else if (value < largest_negative)
 
  122         largest_negative = value;
 
  125     if (largest_positive <= 0 && largest_negative >= 0) 
 
  127       influences_row_ptr[i] = 0;
 
  133     for (
unsigned int nnz_index = row_start; nnz_index < row_stop; ++nnz_index)
 
  135       unsigned int col = A_col_buffer[nnz_index];
 
  140       NumericT value   = A_elements[nnz_index];
 
  149     influences_row_ptr[i] = num_influences;
 
  155   unsigned int current_entry = 0;
 
  156   for (std::size_t i=0; i<A.
size1(); ++i)
 
  158     unsigned int tmp = influences_row_ptr[i];
 
  159     influences_row_ptr[i] = current_entry;
 
  160     current_entry += tmp;
 
  162   influences_row_ptr[A.
size1()] = current_entry;
 
  168 #ifdef VIENNACL_WITH_OPENMP 
  169   #pragma omp parallel for 
  171   for (
long i2=0; i2<static_cast<long>(A.
size1()); ++i2)
 
  173     unsigned int i = 
static_cast<unsigned int>(i2);
 
  174     unsigned int row_start = A_row_buffer[i];
 
  175     unsigned int row_stop  = A_row_buffer[i+1];
 
  181     for (
unsigned int nnz_index = row_start; nnz_index < row_stop; ++nnz_index)
 
  183       unsigned int col = A_col_buffer[nnz_index];
 
  184       NumericT value   = A_elements[nnz_index];
 
  188       else if (value > largest_positive)
 
  189         largest_positive = value;
 
  190       else if (value < largest_negative)
 
  191         largest_negative = value;
 
  194     if (largest_positive <= 0 && largest_negative >= 0) 
 
  199     unsigned int *influences_id_write_ptr = influences_id_ptr + influences_row_ptr[i];
 
  200     for (
unsigned int nnz_index = row_start; nnz_index < row_stop; ++nnz_index)
 
  202       unsigned int col = A_col_buffer[nnz_index];
 
  207       NumericT value   = A_elements[nnz_index];
 
  213         *influences_id_write_ptr = col;
 
  214         ++influences_id_write_ptr;
 
  223 template<
typename NumericT>
 
  237   unsigned int *point_types_ptr  = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
point_types_.
handle());
 
  238   unsigned int *coarse_id_ptr    = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
coarse_id_.
handle());
 
  240   unsigned int coarse_id = 0;
 
  246       coarse_id_ptr[i] = coarse_id++;
 
  283 template<
typename NumericT>
 
  288   unsigned int *point_types_ptr       = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
point_types_.
handle());
 
  289   unsigned int *influences_row_ptr    = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_jumper_.
handle());
 
  290   unsigned int *influences_id_ptr     = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_ids_.
handle());
 
  291   unsigned int *influences_values_ptr = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_values_.
handle());
 
  293   std::set<amg_id_influence, std::greater<amg_id_influence> > points_by_influences;
 
  297   for (std::size_t i=0; i<A.
size1(); ++i)
 
  302   while (!points_by_influences.empty())
 
  307     points_by_influences.erase(points_by_influences.begin());
 
  320     unsigned int j_stop = influences_row_ptr[point.
id + 1];
 
  321     for (
unsigned int j = influences_row_ptr[point.
id]; j < j_stop; ++j)
 
  323       unsigned int influenced_point_id = influences_id_ptr[j];
 
  333       unsigned int k_stop = influences_row_ptr[influenced_point_id + 1];
 
  334       for (
unsigned int k = influences_row_ptr[influenced_point_id]; k < k_stop; ++k)
 
  336         unsigned int influenced_influenced_point_id = influences_id_ptr[k];
 
  340           amg_id_influence point_to_find(influenced_influenced_point_id, influences_values_ptr[influenced_influenced_point_id]);
 
  341           points_by_influences.erase(point_to_find);
 
  344           influences_values_ptr[influenced_influenced_point_id] += 1; 
 
  346           points_by_influences.insert(point_to_find);
 
  366 template<
typename NumericT>
 
  372   unsigned int *point_types_ptr       = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
point_types_.
handle());
 
  373   unsigned int *influences_row_ptr    = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_jumper_.
handle());
 
  374   unsigned int *influences_id_ptr     = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_ids_.
handle());
 
  376   for (
unsigned int i=0; i<static_cast<unsigned int>(A.
size1()); ++i)
 
  379     bool is_new_coarse_node = 
true;
 
  382     unsigned int j_stop = influences_row_ptr[i + 1];
 
  383     for (
unsigned int j = influences_row_ptr[i]; j < j_stop; ++j)
 
  385       unsigned int influenced_point_id = influences_id_ptr[j];
 
  388         is_new_coarse_node = 
false;
 
  393     if (is_new_coarse_node)
 
  396       for (
unsigned int j = influences_row_ptr[i]; j < j_stop; ++j)
 
  398         unsigned int influenced_point_id = influences_id_ptr[j];
 
  417 template<
typename NumericT>
 
  423   unsigned int  *point_types_ptr       = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
point_types_.
handle());
 
  424   unsigned int *influences_row_ptr    = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_jumper_.
handle());
 
  425   unsigned int *influences_id_ptr     = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_ids_.
handle());
 
  427   std::vector<unsigned int> random_weights(A.
size1());
 
  428   for (std::size_t i=0; i<random_weights.size(); ++i)
 
  429     random_weights[i] = static_cast<unsigned int>(rand()) % static_cast<unsigned int>(A.
size1());
 
  431   std::size_t num_threads = 1;
 
  432 #ifdef VIENNACL_WITH_OPENMP 
  433   num_threads = omp_get_max_threads();
 
  440   unsigned int *work_state_ptr     = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(work_state.handle());
 
  441   unsigned int *work_random_ptr    = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(work_random.handle());
 
  442   unsigned int *work_index_ptr     = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(work_index.handle());
 
  448   unsigned int *work_state2_ptr     = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(work_state2.handle());
 
  449   unsigned int *work_random2_ptr    = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(work_random2.handle());
 
  450   unsigned int *work_index2_ptr     = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(work_index2.handle());
 
  453   unsigned int num_undecided = 
static_cast<unsigned int>(A.
size1());
 
  454   unsigned int pmis_iters = 0;
 
  455   while (num_undecided > 0)
 
  462 #ifdef VIENNACL_WITH_OPENMP 
  463   #pragma omp parallel for 
  465     for (
long i2=0; i2<static_cast<long>(A.
size1()); ++i2)
 
  467       unsigned int i = 
static_cast<unsigned int>(i2);
 
  468       switch (point_types_ptr[i])
 
  474         throw std::runtime_error(
"Unexpected state encountered in MIS2 setup for AMG.");
 
  477       work_random_ptr[i] = random_weights[i];
 
  478       work_index_ptr[i]  = i;
 
  485     for (
unsigned int r = 0; r < 2; ++r)
 
  488 #ifdef VIENNACL_WITH_OPENMP 
  489       #pragma omp parallel for 
  491       for (
long i2=0; i2<static_cast<long>(A.
size1()); ++i2)
 
  493         unsigned int i = 
static_cast<unsigned int>(i2);
 
  495         unsigned int state  = work_state_ptr[i];
 
  496         unsigned int random = work_random_ptr[i];
 
  497         unsigned int index  = work_index_ptr[i];
 
  500         unsigned int j_stop = influences_row_ptr[i + 1];
 
  501         for (
unsigned int j = influences_row_ptr[i]; j < j_stop; ++j)
 
  503           unsigned int influenced_point_id = influences_id_ptr[j];
 
  506           if (state < work_state_ptr[influenced_point_id])
 
  508             state  = work_state_ptr[influenced_point_id];
 
  509             random = work_random_ptr[influenced_point_id];
 
  510             index  = work_index_ptr[influenced_point_id];
 
  512           else if (state == work_state_ptr[influenced_point_id])
 
  514             if (random < work_random_ptr[influenced_point_id])
 
  516               state  = work_state_ptr[influenced_point_id];
 
  517               random = work_random_ptr[influenced_point_id];
 
  518               index  = work_index_ptr[influenced_point_id];
 
  520             else if (random == work_random_ptr[influenced_point_id])
 
  522               if (index < work_index_ptr[influenced_point_id])
 
  524                 state  = work_state_ptr[influenced_point_id];
 
  525                 random = work_random_ptr[influenced_point_id];
 
  526                 index  = work_index_ptr[influenced_point_id];
 
  533         work_state2_ptr[i]  = state;
 
  534         work_random2_ptr[i] = random;
 
  535         work_index2_ptr[i]  = index;
 
  539 #ifdef VIENNACL_WITH_OPENMP 
  540       #pragma omp parallel for 
  542       for (
long i2=0; i2<static_cast<long>(A.
size1()); ++i2)
 
  544         unsigned int i = 
static_cast<unsigned int>(i2);
 
  545         work_state_ptr[i]  = work_state2_ptr[i];
 
  546         work_random_ptr[i] = work_random2_ptr[i];
 
  547         work_index_ptr[i]  = work_index2_ptr[i];
 
  554     std::vector<unsigned int> thread_buffer(num_threads);
 
  556 #ifdef VIENNACL_WITH_OPENMP 
  557     #pragma omp parallel for 
  559     for (
long i2=0; i2<static_cast<long>(A.
size1()); ++i2)
 
  561       unsigned int i = 
static_cast<unsigned int>(i2);
 
  562       unsigned int max_state  = work_state_ptr[i];
 
  563       unsigned int max_index  = work_index_ptr[i];
 
  569         else if (max_state == 2) 
 
  572 #ifdef VIENNACL_WITH_OPENMP 
  573           thread_buffer[omp_get_thread_num()] += 1;
 
  575           thread_buffer[0] += 1;
 
  581     for (std::size_t i=0; i<thread_buffer.size(); ++i)
 
  582       num_undecided += thread_buffer[i];
 
  586 #ifdef VIENNACL_WITH_OPENMP 
  587   #pragma omp parallel for 
  589   for (
long i=0; i<static_cast<long>(A.
size1()); ++i)
 
  603 template<
typename NumericT>
 
  608   unsigned int *point_types_ptr       = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
point_types_.
handle());
 
  609   unsigned int *influences_row_ptr    = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_jumper_.
handle());
 
  610   unsigned int *influences_id_ptr     = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_ids_.
handle());
 
  611   unsigned int *coarse_id_ptr         = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
coarse_id_.
handle());
 
  626 #ifdef VIENNACL_WITH_OPENMP 
  627   #pragma omp parallel for 
  629   for (
long i2=0; i2<static_cast<long>(A.
size1()); ++i2)
 
  631     unsigned int i = 
static_cast<unsigned int>(i2);
 
  634       unsigned int coarse_index = coarse_id_ptr[i];
 
  636       unsigned int j_stop = influences_row_ptr[i + 1];
 
  637       for (
unsigned int j = influences_row_ptr[i]; j < j_stop; ++j)
 
  639         unsigned int influenced_point_id = influences_id_ptr[j];
 
  640         coarse_id_ptr[influenced_point_id] = coarse_index; 
 
  642         if (influenced_point_id != i) 
 
  652 #ifdef VIENNACL_WITH_OPENMP 
  653   #pragma omp parallel for 
  655   for (
long i2=0; i2<static_cast<long>(A.
size1()); ++i2)
 
  657     unsigned int i = 
static_cast<unsigned int>(i2);
 
  660       unsigned int j_stop = influences_row_ptr[i + 1];
 
  661       for (
unsigned int j = influences_row_ptr[i]; j < j_stop; ++j)
 
  663         unsigned int influenced_point_id = influences_id_ptr[j];
 
  667           coarse_id_ptr[i] = coarse_id_ptr[influenced_point_id];
 
  678 #ifdef VIENNACL_WITH_OPENMP 
  679   #pragma omp parallel for 
  681   for (
long i=0; i<static_cast<long>(A.
size1()); ++i)
 
  696 template<
typename MatrixT>
 
  723 template<
typename NumericT>
 
  729   NumericT     const * A_elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(A.
handle());
 
  730   unsigned int const * A_row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle1());
 
  731   unsigned int const * A_col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle2());
 
  733   unsigned int *point_types_ptr       = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
point_types_.
handle());
 
  734   unsigned int *influences_row_ptr    = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_jumper_.
handle());
 
  735   unsigned int *influences_id_ptr     = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
influence_ids_.
handle());
 
  736   unsigned int *coarse_id_ptr         = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
coarse_id_.
handle());
 
  740   std::vector<std::map<unsigned int, NumericT> > P_setup(A.
size1());
 
  745 #ifdef VIENNACL_WITH_OPENMP 
  746   #pragma omp parallel for 
  748   for (
long row2=0; row2<static_cast<long>(A.
size1()); ++row2)
 
  750     unsigned int row = 
static_cast<unsigned int>(row2);
 
  751     std::map<unsigned int, NumericT> & P_setup_row = P_setup[
row];
 
  767       unsigned int row_A_start = A_row_buffer[
row];
 
  768       unsigned int row_A_end   = A_row_buffer[row + 1];
 
  769       unsigned int const *influence_iter = influences_id_ptr + influences_row_ptr[
row];
 
  770       unsigned int const *influence_end  = influences_id_ptr + influences_row_ptr[row + 1];
 
  771       for (
unsigned int index = row_A_start; index < row_A_end; ++index)
 
  773         unsigned int col = A_col_buffer[index];
 
  784           while (influence_iter != influence_end && *influence_iter < col)
 
  787           if (influence_iter != influence_end && *influence_iter == col)
 
  788             row_coarse_sum += value;
 
  797       if (std::fabs(temp_res) > 1e-2 * std::fabs(diag))
 
  800         influence_iter = influences_id_ptr + influences_row_ptr[
row];
 
  801         for (
unsigned int index = row_A_start; index < row_A_end; ++index)
 
  803           unsigned int col = A_col_buffer[index];
 
  809           while (influence_iter != influence_end && *influence_iter < col)
 
  812           if (influence_iter != influence_end && *influence_iter == col)
 
  815             P_setup_row[coarse_id_ptr[col]] = temp_res * value;
 
  824       throw std::runtime_error(
"Logic error in direct interpolation: Point is neither coarse-point nor fine-point!");
 
  840 template<
typename NumericT>
 
  849   NumericT     * P_elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(P.
handle());
 
  850   unsigned int * P_row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(P.
handle1());
 
  851   unsigned int * P_col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(P.
handle2());
 
  853   unsigned int *coarse_id_ptr = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(amg_context.
coarse_id_.
handle());
 
  856 #ifdef VIENNACL_WITH_OPENMP 
  857   #pragma omp parallel for 
  859   for (
long row2 = 0; row2 < long(A.
size1()); ++row2)
 
  861     unsigned int row = 
static_cast<unsigned int>(row2);
 
  864     P_col_buffer[
row] = coarse_id_ptr[
row];
 
  866   P_row_buffer[A.
size1()] = 
static_cast<unsigned int>(A.
size1()); 
 
  879 template<
typename NumericT>
 
  891   unsigned int const * A_row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle1());
 
  892   unsigned int const * A_col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle2());
 
  893   NumericT     const * A_elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(A.
handle());
 
  896   unsigned int * Jacobi_row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(Jacobi.handle1());
 
  897   unsigned int * Jacobi_col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(Jacobi.handle2());
 
  898   NumericT     * Jacobi_elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(Jacobi.handle());
 
  902 #ifdef VIENNACL_WITH_OPENMP 
  903   #pragma omp parallel for 
  905   for (
long row2=0; row2<static_cast<long>(A.
size1()); ++row2)
 
  907     unsigned int row = 
static_cast<unsigned int>(row2);
 
  908     unsigned int row_begin = A_row_buffer[
row];
 
  909     unsigned int row_end   = A_row_buffer[row+1];
 
  911     Jacobi_row_buffer[
row] = row_begin;
 
  915     for (
unsigned int j = row_begin; j < row_end; ++j)
 
  917       if (A_col_buffer[j] == row)
 
  919         diag = A_elements[j];
 
  925     for (
unsigned int j = row_begin; j < row_end; ++j)
 
  927       unsigned int col_index = A_col_buffer[j];
 
  928       Jacobi_col_buffer[j] = col_index;
 
  930       if (col_index == row)
 
  936   Jacobi_row_buffer[A.
size1()] = 
static_cast<unsigned int>(Jacobi.nnz()); 
 
  951 template<
typename MatrixT>
 
  962   default: 
throw std::runtime_error(
"Not implemented yet!");
 
  971 template<
typename NumericT>
 
  975   NumericT     const * A_elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(A.
handle());
 
  976   unsigned int const * A_row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle1());
 
  977   unsigned int const * A_col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle2());
 
  982   NumericT     * B_elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(B.handle());
 
  983   unsigned int * B_row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(B.handle1());
 
  984   unsigned int * B_col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(B.handle2());
 
  987   for (std::size_t i = 0; i < B.size1(); ++i)
 
  995     unsigned int row_start = A_row_buffer[
row];
 
  996     unsigned int row_stop  = A_row_buffer[
row+1];
 
  998     for (
unsigned int nnz_index = row_start; nnz_index < row_stop; ++nnz_index)
 
  999       B_row_buffer[A_col_buffer[nnz_index]] += 1;
 
 1003   unsigned int offset = B_row_buffer[0];
 
 1004   B_row_buffer[0] = 0;
 
 1005   for (std::size_t 
row = 1; 
row < B.size1(); ++
row)
 
 1007     unsigned int tmp = B_row_buffer[
row];
 
 1008     B_row_buffer[
row] = offset;
 
 1011   B_row_buffer[B.size1()] = offset;
 
 1017   std::vector<unsigned int> B_row_offsets(B.size1()); 
 
 1022     unsigned int row_start = A_row_buffer[
row];
 
 1023     unsigned int row_stop  = A_row_buffer[
row+1];
 
 1025     for (
unsigned int nnz_index = row_start; nnz_index < row_stop; ++nnz_index)
 
 1027       unsigned int col_in_A = A_col_buffer[nnz_index];
 
 1028       unsigned int B_nnz_index = B_row_buffer[col_in_A] + B_row_offsets[col_in_A];
 
 1029       B_col_buffer[B_nnz_index] = 
static_cast<unsigned int>(
row);
 
 1030       B_elements[B_nnz_index] = A_elements[nnz_index];
 
 1031       ++B_row_offsets[col_in_A];
 
 1041 template<
typename NumericT, 
unsigned int AlignmentV>
 
 1045   NumericT     const * A_elements   = detail::extract_raw_pointer<NumericT>(A.
handle());
 
 1046   unsigned int const * A_row_buffer = detail::extract_raw_pointer<unsigned int>(A.
handle1());
 
 1047   unsigned int const * A_col_buffer = detail::extract_raw_pointer<unsigned int>(A.
handle2());
 
 1049   NumericT           * B_elements   = detail::extract_raw_pointer<NumericT>(B.
handle());
 
 1051 #ifdef VIENNACL_WITH_OPENMP 
 1052   #pragma omp parallel for 
 1054   for (
long row = 0; row < static_cast<long>(A.
size1()); ++
row)
 
 1056     unsigned int row_stop  = A_row_buffer[
row+1];
 
 1058     for (
unsigned int nnz_index = A_row_buffer[
row]; nnz_index < row_stop; ++nnz_index)
 
 1059       B_elements[static_cast<unsigned int>(
row) * 
static_cast<unsigned int>(B.
internal_size2()) + A_col_buffer[nnz_index]] = A_elements[nnz_index];
 
 1073 template<
typename NumericT>
 
 1082   NumericT     const * A_elements   = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(A.
handle());
 
 1083   unsigned int const * A_row_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle1());
 
 1084   unsigned int const * A_col_buffer = viennacl::linalg::host_based::detail::extract_raw_pointer<unsigned int>(A.
handle2());
 
 1085   NumericT     const * rhs_elements = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(rhs_smooth.
handle());
 
 1087   NumericT           * x_elements     = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(x.
handle());
 
 1088   NumericT     const * x_old_elements = viennacl::linalg::host_based::detail::extract_raw_pointer<NumericT>(x_backup.
handle());
 
 1090   for (
unsigned int i=0; i<iterations; ++i)
 
 1094     #ifdef VIENNACL_WITH_OPENMP 
 1095     #pragma omp parallel for 
 1097     for (
long row2 = 0; row2 < static_cast<long>(A.
size1()); ++row2)
 
 1099       unsigned int row = 
static_cast<unsigned int>(row2);
 
 1100       unsigned int col_end   = A_row_buffer[row+1];
 
 1104       for (
unsigned int index = A_row_buffer[row]; index != col_end; ++index)
 
 1106         unsigned int col = A_col_buffer[index];
 
 1108           diag = A_elements[index];
 
 1110           sum += A_elements[index] * x_old_elements[col];
 
 1113       x_elements[
row] = weight * (rhs_elements[
row] - 
sum) / diag + (
NumericT(1) - weight) * x_old_elements[row];
 
const vcl_size_t & size2() const 
Returns the number of columns. 
Helper struct for sequential classical one-pass coarsening. 
amg_coarsening_method get_coarsening_method() const 
Returns the current coarsening strategy. 
viennacl::scalar_expression< const viennacl::vector_base< NumericT >, const viennacl::vector_base< NumericT >, viennacl::op_sum > sum(viennacl::vector_base< NumericT > const &x)
User interface function for computing the sum of all elements of a vector. 
const vcl_size_t & size1() const 
Returns the number of rows. 
viennacl::vector_expression< const viennacl::matrix_base< NumericT >, const viennacl::matrix_base< NumericT >, viennacl::op_row_sum > row_sum(viennacl::matrix_base< NumericT > const &A)
User interface function for computing the sum of all elements of each row of a matrix. 
void amg_coarse_ag(compressed_matrix< NumericT > const &A, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
AG (aggregation based) coarsening. Partially single-threaded version (VIENNACL_AMG_COARSE_AG) ...
double get_jacobi_weight() const 
Returns the Jacobi smoother weight (damping). 
void amg_interpol_ag(compressed_matrix< NumericT > const &A, compressed_matrix< NumericT > &P, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
AG (aggregation based) interpolation. Multi-Threaded! (VIENNACL_INTERPOL_AG) 
void amg_influence_advanced(compressed_matrix< NumericT > const &A, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
Routine for extracting strongly connected points considering a user-provided threshold value...
void amg_coarse_classic_onepass(compressed_matrix< NumericT > const &A, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
Classical (RS) one-pass coarsening. Single-Threaded! (VIENNACL_AMG_COARSE_CLASSIC_ONEPASS) ...
amg_id_influence(std::size_t id2, std::size_t influences2)
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. 
const vcl_size_t & nnz() const 
Returns the number of nonzero entries. 
void amg_coarse_ag_stage1_sequential(compressed_matrix< NumericT > const &A, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
AG (aggregation based) coarsening, single-threaded version of stage 1. 
VectorT prod(std::vector< std::vector< T, A1 >, A2 > const &matrix, VectorT const &vector)
void generate_row_block_information()
Builds the row block information needed for fast sparse matrix-vector multiplications. 
void amg_transpose(compressed_matrix< NumericT > const &A, compressed_matrix< NumericT > &B)
Computes B = trans(A). 
bool operator>(amg_id_influence const &a, amg_id_influence const &b)
void assign_to_dense(viennacl::compressed_matrix< NumericT, AlignmentV > const &A, viennacl::matrix_base< NumericT > &B)
double get_strong_connection_threshold() const 
Returns the strong connection threshold parameter. 
const handle_type & handle2() const 
Returns the OpenCL handle to the column index array. 
void enumerate_coarse_points(viennacl::linalg::detail::amg::amg_level_context &amg_context)
Assign IDs to coarse points. 
handle_type & handle()
Returns the OpenCL handle, non-const-version. 
vector_expression< const matrix_base< NumericT >, const int, op_matrix_diag > diag(const matrix_base< NumericT > &A, int k=0)
void smooth_jacobi(unsigned int iterations, compressed_matrix< NumericT > const &A, vector< NumericT > &x, vector< NumericT > &x_backup, vector< NumericT > const &rhs_smooth, NumericT weight)
Damped Jacobi Smoother (CUDA version) 
vector_expression< const matrix_base< NumericT, F >, const unsigned int, op_row > row(const matrix_base< NumericT, F > &A, unsigned int i)
void amg_coarse_ag_stage1_mis2(compressed_matrix< NumericT > const &A, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
AG (aggregation based) coarsening, multi-threaded version of stage 1 using parallel maximum independe...
void amg_interpol_sa(compressed_matrix< NumericT > const &A, compressed_matrix< NumericT > &P, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
Smoothed aggregation interpolation. Multi-Threaded! (VIENNACL_INTERPOL_SA) 
A tag for algebraic multigrid (AMG). Used to transport information from the user to the implementatio...
void amg_influence(compressed_matrix< NumericT > const &A, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
Dispatcher for influence processing. 
viennacl::context context(T const &t)
Returns an ID for the currently active memory domain of an object. 
void amg_influence_trivial(compressed_matrix< NumericT > const &A, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
Routine for taking all connections in the matrix as strong. 
void copy(std::vector< NumericT > &cpu_vec, circulant_matrix< NumericT, AlignmentV > &gpu_mat)
Copies a circulant matrix from the std::vector to the OpenCL device (either GPU or multi-core CPU) ...
void amg_coarse(MatrixT &A, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
Entry point and dispatcher for coarsening procedures. 
size_type size() const 
Returns the length of the vector (cf. std::vector) 
viennacl::vector< unsigned int > point_types_
amg_interpolation_method get_interpolation_method() const 
Returns the current interpolation method. 
size_type internal_size2() const 
Returns the internal number of columns. Usually required for launching OpenCL kernels only...
Helper classes and functions for the AMG preconditioner. Experimental. 
void amg_interpol(MatrixT const &A, MatrixT &P, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
Dispatcher for building the interpolation matrix. 
void amg_interpol_direct(compressed_matrix< NumericT > const &A, compressed_matrix< NumericT > &P, viennacl::linalg::detail::amg::amg_level_context &amg_context, viennacl::linalg::amg_tag &tag)
Direct interpolation. Multi-threaded! (VIENNACL_AMG_INTERPOL_DIRECT) 
viennacl::vector< unsigned int > coarse_id_
void resize(vcl_size_t new_size1, vcl_size_t new_size2, bool preserve=true)
Resize the matrix. 
const handle_type & handle() const 
Returns the memory handle. 
viennacl::vector< unsigned int > influence_values_
viennacl::vector< unsigned int > influence_ids_
viennacl::vector< unsigned int > influence_jumper_