ViennaCL - The Vienna Computing Library  1.5.2
amg_coarse.hpp
Go to the documentation of this file.
1 #ifndef VIENNACL_LINALG_DETAIL_AMG_AMG_COARSE_HPP
2 #define VIENNACL_LINALG_DETAIL_AMG_AMG_COARSE_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 <cmath>
27 
28 #include <map>
29 #ifdef VIENNACL_WITH_OPENMP
30 #include <omp.h>
31 #endif
32 
34 
35 namespace viennacl
36 {
37  namespace linalg
38  {
39  namespace detail
40  {
41  namespace amg
42  {
43 
51  template <typename InternalType1, typename InternalType2, typename InternalType3>
52  void amg_coarse(unsigned int level, InternalType1 & A, InternalType2 & Pointvector, InternalType3 & Slicing, amg_tag & tag)
53  {
54  switch (tag.get_coarse())
55  {
56  case VIENNACL_AMG_COARSE_RS: amg_coarse_classic (level, A, Pointvector, tag); break;
57  case VIENNACL_AMG_COARSE_ONEPASS: amg_coarse_classic_onepass (level, A, Pointvector, tag); break;
58  case VIENNACL_AMG_COARSE_RS0: amg_coarse_rs0 (level, A, Pointvector, Slicing, tag); break;
59  case VIENNACL_AMG_COARSE_RS3: amg_coarse_rs3 (level, A, Pointvector, Slicing, tag); break;
60  case VIENNACL_AMG_COARSE_AG: amg_coarse_ag (level, A, Pointvector, tag); break;
61  }
62  }
63 
70  template <typename InternalType1, typename InternalType2>
71  void amg_influence(unsigned int level, InternalType1 const & A, InternalType2 & Pointvector, amg_tag & tag)
72  {
73  typedef typename InternalType1::value_type SparseMatrixType;
74  typedef typename InternalType2::value_type PointVectorType;
75  typedef typename SparseMatrixType::value_type ScalarType;
76  typedef typename SparseMatrixType::value_type ScalarType;
77  typedef typename SparseMatrixType::const_iterator1 ConstRowIterator;
78  typedef typename SparseMatrixType::const_iterator2 ConstColIterator;
79 
80  ScalarType max;
81  int diag_sign;
82  //unsigned int i;
83 
84 #ifdef VIENNACL_WITH_OPENMP
85  #pragma omp parallel for private (max,diag_sign)
86 #endif
87  for (long i=0; i<static_cast<long>(A[level].size1()); ++i)
88  {
89  diag_sign = 1;
90  if (A[level](i,i) < 0)
91  diag_sign = -1;
92 
93  ConstRowIterator row_iter = A[level].begin1();
94  row_iter += i;
95  // Find greatest non-diagonal negative value (positive if diagonal is negative) in row
96  max = 0;
97  for (ConstColIterator col_iter = row_iter.begin(); col_iter != row_iter.end(); ++col_iter)
98  {
99  if (i == (unsigned int) col_iter.index2()) continue;
100  if (diag_sign == 1)
101  if (max > *col_iter) max = *col_iter;
102  if (diag_sign == -1)
103  if (max < *col_iter) max = *col_iter;
104  }
105 
106  // If maximum is 0 then the row is independent of the others
107  if (max == 0)
108  continue;
109 
110  // Find all points that strongly influence current point (Yang, p.5)
111  for (ConstColIterator col_iter = row_iter.begin(); col_iter != row_iter.end(); ++col_iter)
112  {
113  unsigned int j = static_cast<unsigned int>(col_iter.index2());
114  if (i == j) continue;
115  if (diag_sign * (-*col_iter) >= tag.get_threshold() * (diag_sign * (-max)))
116  {
117  // Strong influence from j to i found, save information
118  Pointvector[level][i]->add_influencing_point(Pointvector[level][j]);
119  }
120  }
121  }
122 
123  #ifdef VIENNACL_AMG_DEBUG
124  std::cout << "Influence Matrix: " << std::endl;
125  boost::numeric::ublas::matrix<bool> mat;
126  Pointvector[level].get_influence_matrix(mat);
127  printmatrix (mat);
128  #endif
129 
130  // Save influenced points
131  for (typename PointVectorType::iterator iter = Pointvector[level].begin(); iter != Pointvector[level].end(); ++iter)
132  {
133  for (typename amg_point::iterator iter2 = (*iter)->begin_influencing(); iter2 != (*iter)->end_influencing(); ++iter2)
134  {
135  (*iter2)->add_influenced_point(*iter);
136  }
137  }
138 
139  #ifdef VIENNACL_AMG_DEBUG
140  std::cout << "Influence Measures: " << std::endl;
141  boost::numeric::ublas::vector<unsigned int> temp;
142  Pointvector[level].get_influence(temp);
143  printvector (temp);
144  std::cout << "Point Sorting: " << std::endl;
145  Pointvector[level].get_sorting(temp);
146  printvector (temp);
147  #endif
148  }
149 
156  template <typename InternalType1, typename InternalType2>
157  void amg_coarse_classic_onepass(unsigned int level, InternalType1 & A, InternalType2 & Pointvector, amg_tag & tag)
158  {
159  amg_point* c_point, *point1, *point2;
160 
161  // Check and save all strong influences
162  amg_influence (level, A, Pointvector, tag);
163 
164  // Traverse through points and calculate initial influence measure
165  long i;
166 #ifdef VIENNACL_WITH_OPENMP
167  #pragma omp parallel for private (i)
168 #endif
169  for (i=0; i<static_cast<long>(Pointvector[level].size()); ++i)
170  Pointvector[level][i]->calc_influence();
171 
172  // Do initial sorting
173  Pointvector[level].sort();
174 
175  // Get undecided point with highest influence measure
176  while ((c_point = Pointvector[level].get_nextpoint()) != NULL)
177  {
178  // Make this point C point
179  Pointvector[level].make_cpoint(c_point);
180 
181  // All strongly influenced points become F points
182  for (typename amg_point::iterator iter = c_point->begin_influenced(); iter != c_point->end_influenced(); ++iter)
183  {
184  point1 = *iter;
185  // Found strong influence from C point (c_point influences point1), check whether point is still undecided, otherwise skip
186  if (!point1->is_undecided()) continue;
187  // Make this point F point if it is still undecided point
188  Pointvector[level].make_fpoint(point1);
189 
190  // Add +1 to influence measure for all undecided points that strongly influence new F point
191  for (typename amg_point::iterator iter2 = point1->begin_influencing(); iter2 != point1->end_influencing(); ++iter2)
192  {
193  point2 = *iter2;
194  // Found strong influence to F point (point2 influences point1)
195  if (point2->is_undecided())
196  Pointvector[level].add_influence(point2,1);
197  }
198  }
199  }
200 
201  // If a point is neither C nor F point but is nevertheless influenced by other points make it F point
202  // (this situation can happen when this point does not influence other points and the points that influence this point became F points already)
203  /*#pragma omp parallel for private (i,point1)
204  for (long i=0; i<static_cast<long>(Pointvector[level].size()); ++i)
205  {
206  point1 = Pointvector[level][i];
207  if (point1->is_undecided())
208  {
209  // Undecided point found. Check whether it is influenced by other point and if so: Make it F point.
210  if (point1->number_influencing() > 0)
211  {
212  #pragma omp critical
213  Pointvector[level].make_fpoint(point1);
214  }
215  }
216  }*/
217 
218  #if defined (VIENNACL_AMG_DEBUG)// or defined (VIENNACL_AMG_DEBUGBENCH)
219  unsigned int c_points = Pointvector[level].get_cpoints();
220  unsigned int f_points = Pointvector[level].get_fpoints();
221  std::cout << "1st pass: Level " << level << ": ";
222  std::cout << "No of C points = " << c_points << ", ";
223  std::cout << "No of F points = " << f_points << std::endl;
224  #endif
225 
226  #ifdef VIENNACL_AMG_DEBUG
227  std::cout << "Coarse Points:" << std::endl;
228  boost::numeric::ublas::vector<bool> C;
229  Pointvector[level].get_C(C);
230  printvector (C);
231  std::cout << "Fine Points:" << std::endl;
232  boost::numeric::ublas::vector<bool> F;
233  Pointvector[level].get_F(F);
234  printvector (F);
235  #endif
236  }
237 
244  template <typename InternalType1, typename InternalType2>
245  void amg_coarse_classic(unsigned int level, InternalType1 & A, InternalType2 & Pointvector, amg_tag & tag)
246  {
247  typedef typename InternalType2::value_type PointVectorType;
248 
249  bool add_C;
250  amg_point *c_point, *point1, *point2;
251 
252  // Use one-pass-coarsening as first pass.
253  amg_coarse_classic_onepass(level, A, Pointvector, tag);
254 
255  // 2nd pass: Add more C points if F-F connection does not have a common C point.
256  for (typename PointVectorType::iterator iter = Pointvector[level].begin(); iter != Pointvector[level].end(); ++iter)
257  {
258  point1 = *iter;
259  // If point is F point, check for strong connections.
260  if (point1->is_fpoint())
261  {
262  // Check for strong connections from influencing and influenced points.
263  amg_point::iterator iter2 = point1->begin_influencing();
264  amg_point::iterator iter3 = point1->begin_influenced();
265 
266  // Iterate over both lists at once. This makes sure that points are no checked twice when influence relation is symmetric (which is often the case).
267  // Note: Only works because influencing and influenced lists are sorted by point-index.
268  while(iter2 != point1->end_influencing() || iter3 != point1->end_influenced())
269  {
270  if (iter2 == point1->end_influencing())
271  {
272  point2 = *iter3;
273  ++iter3;
274  }
275  else if (iter3 == point1->end_influenced())
276  {
277  point2 = *iter2;
278  ++iter2;
279  }
280  else
281  {
282  if ((*iter2)->get_index() == (*iter3)->get_index())
283  {
284  point2 = *iter2;
285  ++iter2;
286  ++iter3;
287  }
288  else if ((*iter2)->get_index() < (*iter3)->get_index())
289  {
290  point2 = *iter2;
291  ++iter2;
292  }
293  else
294  {
295  point2 = *iter3;
296  ++iter3;
297  }
298  }
299  // Only check points with higher index as points with lower index have been checked already.
300  if (point2->get_index() < point1->get_index())
301  continue;
302 
303  // If there is a strong connection then it has to either be a C point or a F point with common C point.
304  // C point? Then skip as everything is ok.
305  if (point2->is_cpoint())
306  continue;
307  // F point? Then check whether F points point1 and point2 have a common C point.
308  if (point2->is_fpoint())
309  {
310  add_C = true;
311  // C point is common for two F points if they are both strongly influenced by that C point.
312  // Compare strong influences for point1 and point2.
313  for (amg_point::iterator iter3 = point1->begin_influencing(); iter3 != point1 -> end_influencing(); ++iter3)
314  {
315  c_point = *iter3;
316  // Stop search when strong common influence is found via c_point.
317  if (c_point->is_cpoint())
318  {
319  if (point2->is_influencing(c_point))
320  {
321  add_C = false;
322  break;
323  }
324  }
325  }
326  // No common C point found? Then make second F point to C point.
327  if (add_C == true)
328  Pointvector[level].switch_ftoc(point2);
329  }
330  }
331  }
332  }
333 
334  #ifdef VIENNACL_AMG_DEBUG
335  std::cout << "After 2nd pass:" << std::endl;
336  std::cout << "Coarse Points:" << std::endl;
337  boost::numeric::ublas::vector<bool> C;
338  Pointvector[level].get_C(C);
339  printvector (C);
340  std::cout << "Fine Points:" << std::endl;
341  boost::numeric::ublas::vector<bool> F;
342  Pointvector[level].get_F(F);
343  printvector (F);
344  #endif
345 
346  #ifdef VIENNACL_AMG_DEBUG
347 #ifdef VIENNACL_WITH_OPENMP
348  #pragma omp critical
349 #endif
350  {
351  std::cout << "No C and no F point: ";
352  for (typename PointVectorType::iterator iter = Pointvector[level].begin(); iter != Pointvector[level].end(); ++iter)
353  if ((*iter)->is_undecided())
354  std::cout << (*iter)->get_index() << " ";
355  std::cout << std::endl;
356  }
357  #endif
358  }
359 
367  template <typename InternalType1, typename InternalType2, typename InternalType3>
368  void amg_coarse_rs0(unsigned int level, InternalType1 & A, InternalType2 & Pointvector, InternalType3 & Slicing, amg_tag & tag)
369  {
370  unsigned int total_points;
371 
372  // Slice matrix into parts such that points are distributed among threads
373  Slicing.slice(level, A, Pointvector);
374 
375  // Run classical coarsening in parallel
376  total_points = 0;
377 #ifdef VIENNACL_WITH_OPENMP
378  #pragma omp parallel for
379 #endif
380  for (long i=0; i<static_cast<long>(Slicing.threads_); ++i)
381  {
382  amg_coarse_classic(level,Slicing.A_slice[i],Slicing.Pointvector_slice[i],tag);
383 
384  // Save C points (using Slicing.Offset on the next level as temporary memory)
385  // Note: Number of C points for point i is saved in i+1!! (makes it easier later to compute offset)
386  Slicing.Offset[i+1][level+1] = Slicing.Pointvector_slice[i][level].get_cpoints();
387  #ifdef VIENNACL_WITH_OPENMP
388  #pragma omp critical
389  #endif
390  total_points += Slicing.Pointvector_slice[i][level].get_cpoints();
391  }
392 
393  // If no coarser level can be found on any level then resume and coarsening will stop in amg_coarse()
394  if (total_points != 0)
395  {
396  #ifdef VIENNACL_WITH_OPENMP
397  #pragma omp parallel for
398  #endif
399  for (long i=0; i<static_cast<long>(Slicing.threads_); ++i)
400  {
401  // If no higher coarse level can be found on slice i (saved in Slicing.Offset[i+1][level+1]) then pull C point(s) to the next level
402  if (Slicing.Offset[i+1][level+1] == 0)
403  {
404  // All points become C points
405  for (unsigned int j=0; j<Slicing.A_slice[i][level].size1(); ++j)
406  Slicing.Pointvector_slice[i][level].make_cpoint(Slicing.Pointvector_slice[i][level][j]);
407  Slicing.Offset[i+1][level+1] = static_cast<unsigned int>(Slicing.A_slice[i][level].size1());
408  }
409  }
410 
411  // Build slicing offset from number of C points (offset = total sum of C points on threads with lower number)
412  for (unsigned int i=2; i<=Slicing.threads_; ++i)
413  Slicing.Offset[i][level+1] += Slicing.Offset[i-1][level+1];
414 
415  // Join C and F points
416  Slicing.join(level, Pointvector);
417  }
418 
419  // Calculate global influence measures for interpolation and/or RS3.
420  amg_influence(level, A, Pointvector, tag);
421 
422  #if defined(VIENNACL_AMG_DEBUG)// or defined (VIENNACL_AMG_DEBUGBENCH)
423  for (unsigned int i=0; i<Slicing._threads; ++i)
424  {
425  unsigned int c_points = Slicing.Pointvector_slice[i][level].get_cpoints();
426  unsigned int f_points = Slicing.Pointvector_slice[i][level].get_fpoints();
427  std::cout << "Thread " << i << ": ";
428  std::cout << "No of C points = " << c_points << ", ";
429  std::cout << "No of F points = " << f_points << std::endl;
430  }
431  #endif
432  }
433 
441  template <typename InternalType1, typename InternalType2, typename InternalType3>
442  void amg_coarse_rs3(unsigned int level, InternalType1 & A, InternalType2 & Pointvector, InternalType3 & Slicing, amg_tag & tag)
443  {
444  amg_point *c_point, *point1, *point2;
445  bool add_C;
446  unsigned int i, j;
447 
448  // Run RS0 first (parallel).
449  amg_coarse_rs0(level, A, Pointvector, Slicing, tag);
450 
451  // Save slicing offset
452  boost::numeric::ublas::vector<unsigned int> Offset = boost::numeric::ublas::vector<unsigned int> (Slicing.Offset.size());
453  for (i=0; i<Slicing.Offset.size(); ++i)
454  Offset[i] = Slicing.Offset[i][level];
455 
456  // Correct the coarsening with a third pass: Don't allow strong F-F connections without common C point
457  for (i=0; i<Slicing.threads_; ++i)
458  {
459  //for (j=Slicing.Offset[i][level]; j<Slicing.Offset[i+1][level]; ++j)
460  for (j=Offset[i]; j<Offset[i+1]; ++j)
461  {
462  point1 = Pointvector[level][j];
463  // If point is F point, check for strong connections.
464  if (point1->is_fpoint())
465  {
466  // Check for strong connections from influencing and influenced points.
467  amg_point::iterator iter2 = point1->begin_influencing();
468  amg_point::iterator iter3 = point1->begin_influenced();
469 
470  // Iterate over both lists at once. This makes sure that points are no checked twice when influence relation is symmetric (which is often the case).
471  // Note: Only works because influencing and influenced lists are sorted by point-index.
472  while(iter2 != point1->end_influencing() || iter3 != point1->end_influenced())
473  {
474  if (iter2 == point1->end_influencing())
475  {
476  point2 = *iter3;
477  ++iter3;
478  }
479  else if (iter3 == point1->end_influenced())
480  {
481  point2 = *iter2;
482  ++iter2;
483  }
484  else
485  {
486  if ((*iter2)->get_index() == (*iter3)->get_index())
487  {
488  point2 = *iter2;
489  ++iter2;
490  ++iter3;
491  }
492  else if ((*iter2)->get_index() < (*iter3)->get_index())
493  {
494  point2 = *iter2;
495  ++iter2;
496  }
497  else
498  {
499  point2 = *iter3;
500  ++iter3;
501  }
502  }
503 
504  // Only check points with higher index as points with lower index have been checked already.
505  if (point2->get_index() < point1->get_index())
506  continue;
507 
508  // Only check points that are outside the slicing boundaries (interior F-F connections have already been checked in second pass)
509  //if (point2->get_index() >= Slicing.Offset[i][level] || point2->get_index() < Slicing.Offset[i+1][level])
510  if (point2->get_index() >= Offset[i] && point2->get_index() < Offset[i+1])
511  continue;
512 
513  // If there is a strong connection then it has to either be a C point or a F point with common C point.
514  // C point? Then skip as everything is ok.
515  if (point2->is_cpoint())
516  continue;
517  // F point? Then check whether F points point1 and point2 have a common C point.
518  if (point2->is_fpoint())
519  {
520  add_C = true;
521  // C point is common for two F points if they are both strongly influenced by that C point.
522  // Compare strong influences for point1 and point2.
523  for (amg_point::iterator iter3 = point1->begin_influencing(); iter3 != point1 -> end_influencing(); ++iter3)
524  {
525  c_point = *iter3;
526  // Stop search when strong common influence is found via c_point.
527  if (c_point->is_cpoint())
528  {
529  if (point2->is_influencing(c_point))
530  {
531  add_C = false;
532  break;
533  }
534  }
535  }
536  // No common C point found? Then make second F point to C point.
537  if (add_C == true)
538  {
539  Pointvector[level].switch_ftoc(point2);
540  // Add +1 to offsets as one C point has been added.
541  for (unsigned int j=i+1; j<=Slicing.threads_; ++j)
542  Slicing.Offset[j][level+1]++;
543  }
544  }
545  }
546  }
547  }
548  }
549 
550  #ifdef VIENNACL_AMG_DEBUG
551  std::cout << "After 3rd pass:" << std::endl;
552  std::cout << "Coarse Points:" << std::endl;
553  boost::numeric::ublas::vector<bool> C;
554  Pointvector[level].get_C(C);
555  printvector (C);
556  std::cout << "Fine Points:" << std::endl;
557  boost::numeric::ublas::vector<bool> F;
558  Pointvector[level].get_F(F);
559  printvector (F);
560  #endif
561  }
562 
570  template <typename InternalType1, typename InternalType2>
571  void amg_coarse_ag(unsigned int level, InternalType1 & A, InternalType2 & Pointvector, amg_tag & tag)
572  {
573  typedef typename InternalType1::value_type SparseMatrixType;
574  typedef typename InternalType2::value_type PointVectorType;
575  typedef typename SparseMatrixType::value_type ScalarType;
576 
577  typedef typename SparseMatrixType::iterator1 InternalRowIterator;
578  typedef typename SparseMatrixType::iterator2 InternalColIterator;
579 
580  long x,y;
581  ScalarType diag;
582  amg_point *pointx, *pointy;
583 
584  // Cannot determine aggregates if size == 1 as then a new aggregate would always consist of this point (infinite loop)
585  if (A[level].size1() == 1) return;
586 
587  // SA algorithm (Vanek et al. p.6)
588  // Build neighborhoods
589  #ifdef VIENNACL_WITH_OPENMP
590  #pragma omp parallel for private (x,y,diag)
591  #endif
592  for (x=0; x<static_cast<long>(A[level].size1()); ++x)
593  {
594  InternalRowIterator row_iter = A[level].begin1();
595  row_iter += x;
596  diag = A[level](x,x);
597  for (InternalColIterator col_iter = row_iter.begin(); col_iter != row_iter.end(); ++col_iter)
598  {
599  y = static_cast<long>(col_iter.index2());
600  if (y == x || (std::fabs(*col_iter) >= tag.get_threshold()*pow(0.5, static_cast<double>(level-1)) * std::sqrt(std::fabs(diag*A[level](y,y)))))
601  {
602  // Neighborhood x includes point y
603  Pointvector[level][x]->add_influencing_point(Pointvector[level][y]);
604  }
605  }
606  }
607 
608  #ifdef VIENNACL_AMG_DEBUG
609  std::cout << "Neighborhoods:" << std::endl;
610  boost::numeric::ublas::matrix<bool> mat;
611  Pointvector[level].get_influence_matrix(mat);
612  printmatrix (mat);
613  #endif
614 
615  // Build aggregates from neighborhoods
616  for (typename PointVectorType::iterator iter = Pointvector[level].begin(); iter != Pointvector[level].end(); ++iter)
617  {
618  pointx = (*iter);
619 
620  if (pointx->is_undecided())
621  {
622  // Make center of aggregate to C point and include it to aggregate x.
623  Pointvector[level].make_cpoint(pointx);
624  pointx->set_aggregate (pointx->get_index());
625  for (amg_point::iterator iter2 = pointx->begin_influencing(); iter2 != pointx->end_influencing(); ++iter2)
626  {
627  pointy = (*iter2);
628 
629  if (pointy->is_undecided())
630  {
631  // Make neighbor y to F point and include it to aggregate x.
632  Pointvector[level].make_fpoint(pointy);
633  pointy->set_aggregate (pointx->get_index());
634  }
635  }
636  }
637  }
638 
639  #ifdef VIENNACL_AMG_DEBUG
640  std::cout << "After aggregation:" << std::endl;
641  std::cout << "Coarse Points:" << std::endl;
642  boost::numeric::ublas::vector<bool> C;
643  Pointvector[level].get_C(C);
644  printvector (C);
645  std::cout << "Fine Points:" << std::endl;
646  boost::numeric::ublas::vector<bool> F;
647  Pointvector[level].get_F(F);
648  printvector (F);
649  std::cout << "Aggregates:" << std::endl;
650  printvector (Aggregates[level]);
651  #endif
652  }
653  } //namespace amg
654  }
655  }
656 }
657 
658 #endif
Debug functionality for AMG. To be removed.
void add_influencing_point(amg_point *point)
Definition: amg_base.hpp:867
A class for the AMG points. Saves point index and influence measure Holds information whether point i...
Definition: amg_base.hpp:816
#define VIENNACL_AMG_COARSE_RS
Definition: amg_base.hpp:41
void amg_coarse_classic(unsigned int level, InternalType1 &A, InternalType2 &Pointvector, amg_tag &tag)
Classical (RS) two-pass coarsening. Single-Threaded! (VIENNACL_AMG_COARSE_CLASSIC) ...
Definition: amg_coarse.hpp:245
void amg_coarse_ag(unsigned int level, InternalType1 &A, InternalType2 &Pointvector, amg_tag &tag)
AG (aggregation based) coarsening. Single-Threaded! (VIENNACL_AMG_COARSE_SA)
Definition: amg_coarse.hpp:571
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
#define VIENNACL_AMG_COARSE_RS3
Definition: amg_base.hpp:44
unsigned int get_coarse() const
Definition: amg_base.hpp:90
bool is_undecided() const
Definition: amg_base.hpp:860
void printmatrix(MatrixType &, int)
Definition: amg_debug.hpp:77
void set_aggregate(unsigned int aggregate)
Definition: amg_base.hpp:855
#define VIENNACL_AMG_COARSE_RS0
Definition: amg_base.hpp:43
#define VIENNACL_AMG_COARSE_AG
Definition: amg_base.hpp:45
void make_cpoint()
Definition: amg_base.hpp:890
void amg_coarse_classic_onepass(unsigned int level, InternalType1 &A, InternalType2 &Pointvector, amg_tag &tag)
Classical (RS) one-pass coarsening. Single-Threaded! (VIENNACL_AMG_COARSE_CLASSIC_ONEPASS) ...
Definition: amg_coarse.hpp:157
A tag for algebraic multigrid (AMG). Used to transport information from the user to the implementatio...
Definition: amg_base.hpp:61
iterator begin_influencing()
Definition: amg_base.hpp:907
bool is_cpoint() const
Definition: amg_base.hpp:858
iterator begin_influenced()
Definition: amg_base.hpp:911
void amg_influence(unsigned int level, InternalType1 const &A, InternalType2 &Pointvector, amg_tag &tag)
Determines strong influences in system matrix, classical approach (RS). Multithreaded! ...
Definition: amg_coarse.hpp:71
void make_fpoint()
Definition: amg_base.hpp:897
double get_threshold() const
Definition: amg_base.hpp:96
void printvector(VectorType const &)
Definition: amg_debug.hpp:80
vector_expression< const matrix_base< NumericT, F >, const int, op_matrix_diag > diag(const matrix_base< NumericT, F > &A, int k=0)
Definition: matrix.hpp:895
unsigned int get_index() const
Definition: amg_base.hpp:853
#define VIENNACL_AMG_COARSE_ONEPASS
Definition: amg_base.hpp:42
bool is_influencing(amg_point *point) const
Definition: amg_base.hpp:865
void amg_coarse_rs3(unsigned int level, InternalType1 &A, InternalType2 &Pointvector, InternalType3 &Slicing, amg_tag &tag)
RS3 coarsening. Single-Threaded! (VIENNACL_AMG_COARSE_RS3)
Definition: amg_coarse.hpp:442
Defines an iterator for the sparse vector type.
Definition: amg_base.hpp:203
iterator end_influenced()
Definition: amg_base.hpp:912
Helper classes and functions for the AMG preconditioner. Experimental.
void amg_coarse(unsigned int level, InternalType1 &A, InternalType2 &Pointvector, InternalType3 &Slicing, amg_tag &tag)
Calls the right coarsening procedure.
Definition: amg_coarse.hpp:52
bool is_fpoint() const
Definition: amg_base.hpp:859
iterator end_influencing()
Definition: amg_base.hpp:908
void amg_coarse_rs0(unsigned int level, InternalType1 &A, InternalType2 &Pointvector, InternalType3 &Slicing, amg_tag &tag)
Parallel classical RS0 coarsening. Multi-Threaded! (VIENNACL_AMG_COARSE_RS0 || VIENNACL_AMG_COARSE_RS...
Definition: amg_coarse.hpp:368