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