dune-istl  2.8.0
smoother.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_AMGSMOOTHER_HH
4 #define DUNE_AMGSMOOTHER_HH
5 
9 #include <dune/istl/schwarz.hh>
11 #include <dune/common/propertymap.hh>
12 #include <dune/common/ftraits.hh>
13 
14 namespace Dune
15 {
16  namespace Amg
17  {
18 
34  template<class T>
36  {
40  typedef typename FieldTraits<T>::real_type RelaxationFactor;
41 
50 
55  : iterations(1), relaxationFactor(1.0)
56  {}
57  };
58 
62  template<class T>
64  {
66 
67  };
68 
69  template<class X, class Y>
71  {
73 
74  };
75 
76  template<class X, class Y, class C, class T>
78  : public SmootherTraits<T>
79  {};
80 
81  template<class C, class T>
83  : public SmootherTraits<T>
84  {};
85 
89  template<class T>
91  {
92  typedef typename T::matrix_type Matrix;
93 
95 
97 
98  public:
100  {}
101 
102  void setMatrix(const Matrix& matrix)
103  {
104  matrix_=&matrix;
105  }
106  virtual void setMatrix(const Matrix& matrix, [[maybe_unused]] const AggregatesMap& amap)
107  {
108  setMatrix(matrix);
109  }
110 
111 
112  const Matrix& getMatrix() const
113  {
114  return *matrix_;
115  }
116 
117  void setArgs(const SmootherArgs& args)
118  {
119  args_=&args;
120  }
121 
122  template<class T1>
123  void setComm([[maybe_unused]] T1& comm)
124  {}
125 
127  {
128  return comm_;
129  }
130 
131  const SmootherArgs getArgs() const
132  {
133  return *args_;
134  }
135 
136  protected:
137  const Matrix* matrix_;
138  private:
139  const SmootherArgs* args_;
140  SequentialInformation comm_;
141  };
142 
143  template<class T>
145  : public DefaultConstructionArgs<T>
146  {};
147 
148  template<class T, class C=SequentialInformation>
150  : public ConstructionArgs<T>
151  {
152  public:
154  {}
155 
156  void setComm(const C& comm)
157  {
158  comm_ = &comm;
159  }
160 
161  const C& getComm() const
162  {
163  return *comm_;
164  }
165  private:
166  const C* comm_;
167  };
168 
169 
170  template<class X, class Y>
172  {
173  typedef Richardson<X,Y> T;
174 
176 
177  public:
179  {}
180 
181  template <class... Args>
182  void setMatrix(const Args&...)
183  {}
184 
185  void setArgs(const SmootherArgs& args)
186  {
187  args_=&args;
188  }
189 
190  template<class T1>
191  void setComm([[maybe_unused]] T1& comm)
192  {}
193 
195  {
196  return comm_;
197  }
198 
199  const SmootherArgs getArgs() const
200  {
201  return *args_;
202  }
203 
204  private:
205  const SmootherArgs* args_;
206  SequentialInformation comm_;
207  };
208 
209 
210 
211  template<class T>
212  struct ConstructionTraits;
213 
217  template<class M, class X, class Y, int l>
218  struct ConstructionTraits<SeqSSOR<M,X,Y,l> >
219  {
221 
222  static inline std::shared_ptr<SeqSSOR<M,X,Y,l>> construct(Arguments& args)
223  {
224  return std::make_shared<SeqSSOR<M,X,Y,l>>
225  (args.getMatrix(), args.getArgs().iterations, args.getArgs().relaxationFactor);
226  }
227  };
228 
229 
233  template<class M, class X, class Y, int l>
234  struct ConstructionTraits<SeqSOR<M,X,Y,l> >
235  {
237 
238  static inline std::shared_ptr<SeqSOR<M,X,Y,l>> construct(Arguments& args)
239  {
240  return std::make_shared<SeqSOR<M,X,Y,l>>
241  (args.getMatrix(), args.getArgs().iterations, args.getArgs().relaxationFactor);
242  }
243  };
244 
245 
249  template<class M, class X, class Y, int l>
250  struct ConstructionTraits<SeqJac<M,X,Y,l> >
251  {
253 
254  static inline std::shared_ptr<SeqJac<M,X,Y,l>> construct(Arguments& args)
255  {
256  return std::make_shared<SeqJac<M,X,Y,l>>
257  (args.getMatrix(), args.getArgs().iterations, args.getArgs().relaxationFactor);
258  }
259  };
260 
264  template<class X, class Y>
266  {
268 
269  static inline std::shared_ptr<Richardson<X,Y>> construct(Arguments& args)
270  {
271  return std::make_shared<Richardson<X,Y>>
272  (args.getArgs().relaxationFactor);
273  }
274  };
275 
276 
277  template<class M, class X, class Y>
278  class ConstructionArgs<SeqILU<M,X,Y> >
279  : public DefaultConstructionArgs<SeqILU<M,X,Y> >
280  {
281  public:
283  : n_(n)
284  {}
285 
286  void setN(int n)
287  {
288  n_ = n;
289  }
290 
291  int getN()
292  {
293  return n_;
294  }
295 
296  private:
297  int n_;
298  };
299 
300 
304  template<class M, class X, class Y>
305  struct ConstructionTraits<SeqILU<M,X,Y> >
306  {
308 
309  static inline std::shared_ptr<SeqILU<M,X,Y>> construct(Arguments& args)
310  {
311  return std::make_shared<SeqILU<M,X,Y>>
312  (args.getMatrix(), args.getN(), args.getArgs().relaxationFactor);
313  }
314  };
315 
319  template<class M, class X, class Y, class C>
320  struct ConstructionTraits<ParSSOR<M,X,Y,C> >
321  {
323 
324  static inline std::shared_ptr<ParSSOR<M,X,Y,C>> construct(Arguments& args)
325  {
326  return std::make_shared<ParSSOR<M,X,Y,C>>
327  (args.getMatrix(), args.getArgs().iterations,
328  args.getArgs().relaxationFactor, args.getComm());
329  }
330  };
331 
332  template<class X, class Y, class C, class T>
334  {
337  static inline std::shared_ptr<BlockPreconditioner<X,Y,C,T>> construct(Arguments& args)
338  {
339  auto seqPrec = SeqConstructionTraits::construct(args);
340  return std::make_shared<BlockPreconditioner<X,Y,C,T>> (seqPrec, args.getComm());
341  }
342  };
343 
344  template<class C, class T>
346  {
349  static inline std::shared_ptr<NonoverlappingBlockPreconditioner<C,T>> construct(Arguments& args)
350  {
351  auto seqPrec = SeqConstructionTraits::construct(args);
352  return std::make_shared<NonoverlappingBlockPreconditioner<C,T>> (seqPrec, args.getComm());
353  }
354  };
355 
366  template<class T>
368  {
369  typedef T Smoother;
370  typedef typename Smoother::range_type Range;
371  typedef typename Smoother::domain_type Domain;
372 
380  static void preSmooth(Smoother& smoother, Domain& v, const Range& d)
381  {
382  smoother.apply(v,d);
383  }
384 
392  static void postSmooth(Smoother& smoother, Domain& v, const Range& d)
393  {
394  smoother.apply(v,d);
395  }
396  };
397 
403  template<typename LevelContext>
404  void presmooth(LevelContext& levelContext, size_t steps)
405  {
406  for(std::size_t i=0; i < steps; ++i) {
407  *levelContext.lhs=0;
409  ::preSmooth(*levelContext.smoother, *levelContext.lhs,
410  *levelContext.rhs);
411  // Accumulate update
412  *levelContext.update += *levelContext.lhs;
413 
414  // update defect
415  levelContext.matrix->applyscaleadd(-1, *levelContext.lhs, *levelContext.rhs);
416  levelContext.pinfo->project(*levelContext.rhs);
417  }
418  }
419 
425  template<typename LevelContext>
426  void postsmooth(LevelContext& levelContext, size_t steps)
427  {
428  for(std::size_t i=0; i < steps; ++i) {
429  // update defect
430  levelContext.matrix->applyscaleadd(-1, *levelContext.lhs,
431  *levelContext.rhs);
432  *levelContext.lhs=0;
433  levelContext.pinfo->project(*levelContext.rhs);
435  ::postSmooth(*levelContext.smoother, *levelContext.lhs, *levelContext.rhs);
436  // Accumulate update
437  *levelContext.update += *levelContext.lhs;
438  }
439  }
440 
441  template<class M, class X, class Y, int l>
442  struct SmootherApplier<SeqSOR<M,X,Y,l> >
443  {
445  typedef typename Smoother::range_type Range;
446  typedef typename Smoother::domain_type Domain;
447 
448  static void preSmooth(Smoother& smoother, Domain& v, Range& d)
449  {
450  smoother.template apply<true>(v,d);
451  }
452 
453 
454  static void postSmooth(Smoother& smoother, Domain& v, Range& d)
455  {
456  smoother.template apply<false>(v,d);
457  }
458  };
459 
460  template<class M, class X, class Y, class C, int l>
461  struct SmootherApplier<BlockPreconditioner<X,Y,C,SeqSOR<M,X,Y,l> > >
462  {
464  typedef typename Smoother::range_type Range;
465  typedef typename Smoother::domain_type Domain;
466 
467  static void preSmooth(Smoother& smoother, Domain& v, Range& d)
468  {
469  smoother.template apply<true>(v,d);
470  }
471 
472 
473  static void postSmooth(Smoother& smoother, Domain& v, Range& d)
474  {
475  smoother.template apply<false>(v,d);
476  }
477  };
478 
479  template<class M, class X, class Y, class C, int l>
481  {
483  typedef typename Smoother::range_type Range;
484  typedef typename Smoother::domain_type Domain;
485 
486  static void preSmooth(Smoother& smoother, Domain& v, Range& d)
487  {
488  smoother.template apply<true>(v,d);
489  }
490 
491 
492  static void postSmooth(Smoother& smoother, Domain& v, Range& d)
493  {
494  smoother.template apply<false>(v,d);
495  }
496  };
497 
498  } // end namespace Amg
499 
500  // forward declarations
501  template<class M, class X, class MO, class MS, class A>
502  class SeqOverlappingSchwarz;
503 
504  struct MultiplicativeSchwarzMode;
505 
506  namespace Amg
507  {
508  template<class M, class X, class MS, class TA>
510  MS,TA> >
511  {
513  typedef typename Smoother::range_type Range;
514  typedef typename Smoother::domain_type Domain;
515 
516  static void preSmooth(Smoother& smoother, Domain& v, const Range& d)
517  {
518  smoother.template apply<true>(v,d);
519  }
520 
521 
522  static void postSmooth(Smoother& smoother, Domain& v, const Range& d)
523  {
524  smoother.template apply<false>(v,d);
525 
526  }
527  };
528 
529  // template<class M, class X, class TM, class TA>
530  // class SeqOverlappingSchwarz;
531 
532  template<class T>
534  : public DefaultSmootherArgs<T>
535  {
537 
539  bool onthefly;
540 
542  bool onthefly_=false)
543  : overlap(overlap_), onthefly(onthefly_)
544  {}
545  };
546 
547  template<class M, class X, class TM, class TS, class TA>
548  struct SmootherTraits<SeqOverlappingSchwarz<M,X,TM,TS,TA> >
549  {
551  };
552 
553  template<class M, class X, class TM, class TS, class TA>
555  : public DefaultConstructionArgs<SeqOverlappingSchwarz<M,X,TM,TS,TA> >
556  {
558 
559  public:
564  typedef typename Vector::value_type Subdomain;
565 
566  virtual void setMatrix(const M& matrix, const AggregatesMap& amap)
567  {
568  Father::setMatrix(matrix);
569 
570  std::vector<bool> visited(amap.noVertices(), false);
571  typedef IteratorPropertyMap<std::vector<bool>::iterator,IdentityMap> VisitedMapType;
572  VisitedMapType visitedMap(visited.begin());
573 
574  MatrixGraph<const M> graph(matrix);
575 
577 
578  switch(Father::getArgs().overlap) {
579  case SmootherArgs::vertex :
580  {
581  VertexAdder visitor(subdomains, amap);
582  createSubdomains(matrix, graph, amap, visitor, visitedMap);
583  }
584  break;
585  case SmootherArgs::pairwise :
586  {
587  createPairDomains(graph);
588  }
589  break;
590  case SmootherArgs::aggregate :
591  {
592  AggregateAdder<VisitedMapType> visitor(subdomains, amap, graph, visitedMap);
593  createSubdomains(matrix, graph, amap, visitor, visitedMap);
594  }
595  break;
596  case SmootherArgs::none :
597  NoneAdder visitor;
598  createSubdomains(matrix, graph, amap, visitor, visitedMap);
599  break;
600  default :
601  DUNE_THROW(NotImplemented, "This overlapping scheme is not supported!");
602  }
603  }
604  void setMatrix(const M& matrix)
605  {
606  Father::setMatrix(matrix);
607 
608  /* Create aggregates map where each aggregate is just one vertex. */
609  AggregatesMap amap(matrix.N());
610  VertexDescriptor v=0;
611  for(typename AggregatesMap::iterator iter=amap.begin();
612  iter!=amap.end(); ++iter)
613  *iter=v++;
614 
615  std::vector<bool> visited(amap.noVertices(), false);
616  typedef IteratorPropertyMap<std::vector<bool>::iterator,IdentityMap> VisitedMapType;
617  VisitedMapType visitedMap(visited.begin());
618 
619  MatrixGraph<const M> graph(matrix);
620 
622 
623  switch(Father::getArgs().overlap) {
624  case SmootherArgs::vertex :
625  {
626  VertexAdder visitor(subdomains, amap);
627  createSubdomains(matrix, graph, amap, visitor, visitedMap);
628  }
629  break;
630  case SmootherArgs::aggregate :
631  {
632  DUNE_THROW(NotImplemented, "Aggregate overlap is not supported yet");
633  /*
634  AggregateAdder<VisitedMapType> visitor(subdomains, amap, graph, visitedMap);
635  createSubdomains(matrix, graph, amap, visitor, visitedMap);
636  */
637  }
638  break;
639  case SmootherArgs::pairwise :
640  {
641  createPairDomains(graph);
642  }
643  break;
644  case SmootherArgs::none :
645  NoneAdder visitor;
646  createSubdomains(matrix, graph, amap, visitor, visitedMap);
647 
648  }
649  }
650 
652  {
653  return subdomains;
654  }
655 
656  private:
657  struct VertexAdder
658  {
659  VertexAdder(Vector& subdomains_, const AggregatesMap& aggregates_)
660  : subdomains(subdomains_), max(-1), subdomain(-1), aggregates(aggregates_)
661  {}
662  template<class T>
663  void operator()(const T& edge)
664  {
665  if(aggregates[edge.target()]!=AggregatesMap::ISOLATED)
666  subdomains[subdomain].insert(edge.target());
667  }
668  int setAggregate(const AggregateDescriptor& aggregate_)
669  {
670  subdomain=aggregate_;
671  max = std::max(subdomain, aggregate_);
672  return subdomain;
673  }
674  int noSubdomains() const
675  {
676  return max+1;
677  }
678  private:
679  Vector& subdomains;
681  AggregateDescriptor subdomain;
682  const AggregatesMap& aggregates;
683  };
684  struct NoneAdder
685  {
686  template<class T>
687  void operator()(const T& edge)
688  {}
689  int setAggregate(const AggregateDescriptor& aggregate_)
690  {
691  return -1;
692  }
693  int noSubdomains() const
694  {
695  return -1;
696  }
697  };
698 
699  template<class VM>
700  struct AggregateAdder
701  {
702  AggregateAdder(Vector& subdomains_, const AggregatesMap& aggregates_,
703  const MatrixGraph<const M>& graph_, VM& visitedMap_)
704  : subdomains(subdomains_), subdomain(-1), aggregates(aggregates_),
705  adder(subdomains_, aggregates_), graph(graph_), visitedMap(visitedMap_)
706  {}
707  template<class T>
708  void operator()(const T& edge)
709  {
710  subdomains[subdomain].insert(edge.target());
711  // If we (the neighbouring vertex of the aggregate)
712  // are not isolated, add the aggregate we belong to
713  // to the same subdomain using the OneOverlapAdder
714  if(aggregates[edge.target()]!=AggregatesMap::ISOLATED) {
715  assert(aggregates[edge.target()]!=aggregate);
716  typename AggregatesMap::VertexList vlist;
717  aggregates.template breadthFirstSearch<true,false>(edge.target(), aggregate,
718  graph, vlist, adder, adder,
719  visitedMap);
720  }
721  }
722 
723  int setAggregate(const AggregateDescriptor& aggregate_)
724  {
725  adder.setAggregate(aggregate_);
726  aggregate=aggregate_;
727  return ++subdomain;
728  }
729  int noSubdomains() const
730  {
731  return subdomain+1;
732  }
733 
734  private:
735  AggregateDescriptor aggregate;
736  Vector& subdomains;
737  int subdomain;
738  const AggregatesMap& aggregates;
739  VertexAdder adder;
740  const MatrixGraph<const M>& graph;
741  VM& visitedMap;
742  };
743 
744  void createPairDomains(const MatrixGraph<const M>& graph)
745  {
746  typedef typename MatrixGraph<const M>::ConstVertexIterator VIter;
747  typedef typename MatrixGraph<const M>::ConstEdgeIterator EIter;
748  typedef typename M::size_type size_type;
749 
750  std::set<std::pair<size_type,size_type> > pairs;
751  int total=0;
752  for(VIter v=graph.begin(), ve=graph.end(); ve != v; ++v)
753  for(EIter e = v.begin(), ee=v.end(); ee!=e; ++e)
754  {
755  ++total;
756  if(e.source()<e.target())
757  pairs.insert(std::make_pair(e.source(),e.target()));
758  else
759  pairs.insert(std::make_pair(e.target(),e.source()));
760  }
761 
762 
763  subdomains.resize(pairs.size());
764  Dune::dinfo <<std::endl<< "Created "<<pairs.size()<<" ("<<total<<") pair domains"<<std::endl<<std::endl;
765  typedef typename std::set<std::pair<size_type,size_type> >::const_iterator SIter;
766  typename Vector::iterator subdomain=subdomains.begin();
767 
768  for(SIter s=pairs.begin(), se =pairs.end(); se!=s; ++s)
769  {
770  subdomain->insert(s->first);
771  subdomain->insert(s->second);
772  ++subdomain;
773  }
774  std::size_t minsize=10000;
775  std::size_t maxsize=0;
776  int sum=0;
777  for(typename Vector::size_type i=0; i < subdomains.size(); ++i) {
778  sum+=subdomains[i].size();
779  minsize=std::min(minsize, subdomains[i].size());
780  maxsize=std::max(maxsize, subdomains[i].size());
781  }
782  Dune::dinfo<<"Subdomain size: min="<<minsize<<" max="<<maxsize<<" avg="<<(sum/subdomains.size())
783  <<" no="<<subdomains.size()<<std::endl;
784  }
785 
786  template<class Visitor>
787  void createSubdomains(const M& matrix, const MatrixGraph<const M>& graph,
788  const AggregatesMap& amap, Visitor& overlapVisitor,
789  IteratorPropertyMap<std::vector<bool>::iterator,IdentityMap>& visitedMap )
790  {
791  // count number ag aggregates. We assume that the
792  // aggregates are numbered consecutively from 0 except
793  // for the isolated ones. All isolated vertices form
794  // one aggregate, here.
795  int isolated=0;
796  AggregateDescriptor maxAggregate=0;
797 
798  for(std::size_t i=0; i < amap.noVertices(); ++i)
799  if(amap[i]==AggregatesMap::ISOLATED)
800  isolated++;
801  else
802  maxAggregate = std::max(maxAggregate, amap[i]);
803 
804  subdomains.resize(maxAggregate+1+isolated);
805 
806  // reset the subdomains
807  for(typename Vector::size_type i=0; i < subdomains.size(); ++i)
808  subdomains[i].clear();
809 
810  // Create the subdomains from the aggregates mapping.
811  // For each aggregate we mark all entries and the
812  // neighbouring vertices as belonging to the same subdomain
813  VertexAdder aggregateVisitor(subdomains, amap);
814 
815  for(VertexDescriptor i=0; i < amap.noVertices(); ++i)
816  if(!get(visitedMap, i)) {
817  AggregateDescriptor aggregate=amap[i];
818 
819  if(amap[i]==AggregatesMap::ISOLATED) {
820  // isolated vertex gets its own aggregate
821  subdomains.push_back(Subdomain());
822  aggregate=subdomains.size()-1;
823  }
824  overlapVisitor.setAggregate(aggregate);
825  aggregateVisitor.setAggregate(aggregate);
826  subdomains[aggregate].insert(i);
827  typename AggregatesMap::VertexList vlist;
828  amap.template breadthFirstSearch<false,false>(i, aggregate, graph, vlist, aggregateVisitor,
829  overlapVisitor, visitedMap);
830  }
831 
832  std::size_t minsize=10000;
833  std::size_t maxsize=0;
834  int sum=0;
835  for(typename Vector::size_type i=0; i < subdomains.size(); ++i) {
836  sum+=subdomains[i].size();
837  minsize=std::min(minsize, subdomains[i].size());
838  maxsize=std::max(maxsize, subdomains[i].size());
839  }
840  Dune::dinfo<<"Subdomain size: min="<<minsize<<" max="<<maxsize<<" avg="<<(sum/subdomains.size())
841  <<" no="<<subdomains.size()<<" isolated="<<isolated<<std::endl;
842 
843 
844 
845  }
846  Vector subdomains;
847  };
848 
849 
850  template<class M, class X, class TM, class TS, class TA>
852  {
854 
855  static inline std::shared_ptr<SeqOverlappingSchwarz<M,X,TM,TS,TA>> construct(Arguments& args)
856  {
857  return std::make_shared<SeqOverlappingSchwarz<M,X,TM,TS,TA>>
858  (args.getMatrix(),
859  args.getSubDomains(),
860  args.getArgs().relaxationFactor,
861  args.getArgs().onthefly);
862  }
863  };
864 
865 
866  } // namespace Amg
867 } // namespace Dune
868 
869 
870 
871 #endif
Provides classes for the Coloring process of AMG.
Helper classes for the construction of classes without empty constructor.
Define general preconditioner interface.
const Matrix & getMatrix() const
Definition: smoother.hh:112
DefaultSmootherArgs< typename X::field_type > Arguments
Definition: smoother.hh:72
static std::shared_ptr< SeqILU< M, X, Y > > construct(Arguments &args)
Definition: smoother.hh:309
static void postSmooth(Smoother &smoother, Domain &v, Range &d)
Definition: smoother.hh:492
const SequentialInformation & getComm()
Definition: smoother.hh:194
static void preSmooth(Smoother &smoother, Domain &v, const Range &d)
Definition: smoother.hh:516
ConstructionArgs< SeqILU< M, X, Y > > Arguments
Definition: smoother.hh:307
int setAggregate(const AggregateDescriptor &aggregate_)
Definition: smoother.hh:689
DefaultSmootherArgs< typename T::matrix_type::field_type > Arguments
Definition: smoother.hh:65
int setAggregate(const AggregateDescriptor &aggregate_)
Definition: smoother.hh:668
ConstructionTraits< T > SeqConstructionTraits
Definition: smoother.hh:348
void setArgs(const SmootherArgs &args)
Definition: smoother.hh:117
ConstructionArgs< SeqOverlappingSchwarz< M, X, TM, TS, TA > > Arguments
Definition: smoother.hh:853
NonoverlappingBlockPreconditioner< C, SeqSOR< M, X, Y, l > > Smoother
Definition: smoother.hh:482
static void postSmooth(Smoother &smoother, Domain &v, Range &d)
Definition: smoother.hh:454
AggregateAdder(Vector &subdomains_, const AggregatesMap &aggregates_, const MatrixGraph< const M > &graph_, VM &visitedMap_)
Definition: smoother.hh:702
void setComm(const C &comm)
Definition: smoother.hh:156
DefaultConstructionArgs< Richardson< X, Y > > Arguments
Definition: smoother.hh:267
virtual ~DefaultConstructionArgs()
Definition: smoother.hh:178
SeqOverlappingSchwarzSmootherArgs< typename M::field_type > Arguments
Definition: smoother.hh:550
static void preSmooth(Smoother &smoother, Domain &v, Range &d)
Definition: smoother.hh:486
int getN()
Definition: smoother.hh:291
SeqSOR< M, X, Y, l > Smoother
Definition: smoother.hh:444
static std::shared_ptr< SeqSSOR< M, X, Y, l > > construct(Arguments &args)
Definition: smoother.hh:222
AggregatesMap::AggregateDescriptor AggregateDescriptor
Definition: smoother.hh:562
static std::shared_ptr< Richardson< X, Y > > construct(Arguments &args)
Definition: smoother.hh:269
bool onthefly
Definition: smoother.hh:539
DefaultConstructionArgs< SeqSOR< M, X, Y, l > > Arguments
Definition: smoother.hh:236
virtual void setMatrix(const M &matrix, const AggregatesMap &amap)
Definition: smoother.hh:566
void setMatrix(const Args &...)
Definition: smoother.hh:182
static std::shared_ptr< SeqJac< M, X, Y, l > > construct(Arguments &args)
Definition: smoother.hh:254
DefaultConstructionArgs< SeqJac< M, X, Y, l > > Arguments
Definition: smoother.hh:252
DefaultConstructionArgs< SeqSSOR< M, X, Y, l > > Arguments
Definition: smoother.hh:220
static std::shared_ptr< NonoverlappingBlockPreconditioner< C, T > > construct(Arguments &args)
Definition: smoother.hh:349
const SmootherArgs getArgs() const
Definition: smoother.hh:131
VertexAdder(Vector &subdomains_, const AggregatesMap &aggregates_)
Definition: smoother.hh:659
static std::shared_ptr< SeqSOR< M, X, Y, l > > construct(Arguments &args)
Definition: smoother.hh:238
SeqOverlappingSchwarz< M, X, TM, TS, TA >::subdomain_vector Vector
Definition: smoother.hh:563
static void preSmooth(Smoother &smoother, Domain &v, Range &d)
Definition: smoother.hh:467
const_iterator begin() const
Definition: aggregates.hh:723
void setMatrix(const Matrix &matrix)
Definition: smoother.hh:102
T Smoother
Definition: smoother.hh:369
static void preSmooth(Smoother &smoother, Domain &v, Range &d)
Definition: smoother.hh:448
BlockPreconditioner< X, Y, C, SeqSOR< M, X, Y, l > > Smoother
Definition: smoother.hh:463
static std::shared_ptr< BlockPreconditioner< X, Y, C, T > > construct(Arguments &args)
Definition: smoother.hh:337
AggregateDescriptor * iterator
Definition: aggregates.hh:733
void presmooth(LevelContext &levelContext, size_t steps)
Apply pre smoothing on the current level.
Definition: smoother.hh:404
SeqOverlappingSchwarzSmootherArgs(Overlap overlap_=vertex, bool onthefly_=false)
Definition: smoother.hh:541
DefaultParallelConstructionArgs< T, C > Arguments
Definition: smoother.hh:335
DefaultParallelConstructionArgs< T, C > Arguments
Definition: smoother.hh:347
const Matrix * matrix_
Definition: smoother.hh:137
const_iterator end() const
Definition: aggregates.hh:728
V AggregateDescriptor
The aggregate descriptor type.
Definition: aggregates.hh:578
static const V ISOLATED
Identifier of isolated vertices.
Definition: aggregates.hh:569
const SequentialInformation & getComm()
Definition: smoother.hh:126
DefaultSmootherArgs()
Default constructor.
Definition: smoother.hh:54
void setMatrix(const M &matrix)
Definition: smoother.hh:604
void setComm([[maybe_unused]] T1 &comm)
Definition: smoother.hh:191
std::size_t noVertices() const
Get the number of vertices.
static void postSmooth(Smoother &smoother, Domain &v, const Range &d)
apply post smoothing in forward direction
Definition: smoother.hh:392
FieldTraits< T >::real_type RelaxationFactor
The type of the relaxation factor.
Definition: smoother.hh:40
static std::shared_ptr< ParSSOR< M, X, Y, C > > construct(Arguments &args)
Definition: smoother.hh:324
virtual void setMatrix(const Matrix &matrix, [[maybe_unused]] const AggregatesMap &amap)
Definition: smoother.hh:106
Smoother::domain_type Domain
Definition: smoother.hh:371
void setN(int n)
Definition: smoother.hh:286
SLList< VertexDescriptor, Allocator > VertexList
The type of a single linked list of vertex descriptors.
Definition: aggregates.hh:590
SeqOverlappingSchwarz< M, X, MultiplicativeSchwarzMode, MS, TA > Smoother
Definition: smoother.hh:512
MatrixGraph< M >::VertexDescriptor VertexDescriptor
Definition: smoother.hh:560
const C & getComm() const
Definition: smoother.hh:161
void postsmooth(LevelContext &levelContext, size_t steps)
Apply post smoothing on the current level.
Definition: smoother.hh:426
Dune::Amg::AggregatesMap< VertexDescriptor > AggregatesMap
Definition: smoother.hh:561
const SmootherArgs getArgs() const
Definition: smoother.hh:199
static void postSmooth(Smoother &smoother, Domain &v, const Range &d)
Definition: smoother.hh:522
DefaultParallelConstructionArgs< M, C > Arguments
Definition: smoother.hh:322
RelaxationFactor relaxationFactor
The relaxation factor to use.
Definition: smoother.hh:49
virtual ~DefaultParallelConstructionArgs()
Definition: smoother.hh:153
static std::shared_ptr< SeqOverlappingSchwarz< M, X, TM, TS, TA > > construct(Arguments &args)
Definition: smoother.hh:855
Smoother::domain_type Domain
Definition: smoother.hh:446
int setAggregate(const AggregateDescriptor &aggregate_)
Definition: smoother.hh:723
Smoother::range_type Range
Definition: smoother.hh:445
virtual ~DefaultConstructionArgs()
Definition: smoother.hh:99
static void preSmooth(Smoother &smoother, Domain &v, const Range &d)
apply pre smoothing in forward direction
Definition: smoother.hh:380
static void postSmooth(Smoother &smoother, Domain &v, Range &d)
Definition: smoother.hh:473
ConstructionTraits< T > SeqConstructionTraits
Definition: smoother.hh:336
ConstructionArgs(int n=0)
Definition: smoother.hh:282
void setComm([[maybe_unused]] T1 &comm)
Definition: smoother.hh:123
void setArgs(const SmootherArgs &args)
Definition: smoother.hh:185
int iterations
The numbe of iterations to perform.
Definition: smoother.hh:45
Smoother::range_type Range
Definition: smoother.hh:370
Overlap overlap
Definition: smoother.hh:538
@ aggregate
Definition: smoother.hh:536
@ pairwise
Definition: smoother.hh:536
@ vertex
Definition: smoother.hh:536
Definition: allocator.hh:9
PropertyMapTypeSelector< Amg::VertexVisitedTag, Amg::PropertiesGraph< G, Amg::VertexProperties, EP, VM, EM > >::Type get([[maybe_unused]] const Amg::VertexVisitedTag &tag, Amg::PropertiesGraph< G, Amg::VertexProperties, EP, VM, EM > &graph)
Definition: dependency.hh:291
Sequential overlapping Schwarz preconditioner.
Definition: overlappingschwarz.hh:753
X range_type
The range type of the preconditioner.
Definition: overlappingschwarz.hh:768
X domain_type
The domain type of the preconditioner.
Definition: overlappingschwarz.hh:763
Traits class for generically constructing non default constructable types.
Definition: construction.hh:37
Nonoverlapping parallel preconditioner.
Definition: novlpschwarz.hh:277
P::range_type range_type
The range type of the preconditioner.
Definition: novlpschwarz.hh:285
P::domain_type domain_type
The domain type of the preconditioner.
Definition: novlpschwarz.hh:283
Tag that tells the Schwarz method to be multiplicative.
Definition: overlappingschwarz.hh:124
Class providing information about the mapping of the vertices onto aggregates.
Definition: aggregates.hh:558
VertexIterator end()
Get an iterator over the vertices.
M::size_type VertexDescriptor
The vertex descriptor.
Definition: graph.hh:71
VertexIterator begin()
Get an iterator over the vertices.
Definition: pinfo.hh:26
The default class for the smoother arguments.
Definition: smoother.hh:36
Traits class for getting the attribute class of a smoother.
Definition: smoother.hh:64
Construction Arguments for the default smoothers.
Definition: smoother.hh:91
Definition: smoother.hh:146
Helper class for applying the smoothers.
Definition: smoother.hh:368
Sequential SSOR preconditioner.
Definition: preconditioners.hh:139
Sequential SOR preconditioner.
Definition: preconditioners.hh:259
X domain_type
The domain type of the preconditioner.
Definition: preconditioners.hh:264
Y range_type
The range type of the preconditioner.
Definition: preconditioners.hh:266
The sequential jacobian preconditioner.
Definition: preconditioners.hh:410
Sequential ILU preconditioner.
Definition: preconditioners.hh:530
Richardson preconditioner.
Definition: preconditioners.hh:711
A parallel SSOR preconditioner.
Definition: schwarz.hh:176
Block parallel preconditioner.
Definition: schwarz.hh:279
X domain_type
The domain type of the preconditioner.
Definition: schwarz.hh:286
Y range_type
The range type of the preconditioner.
Definition: schwarz.hh:291