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