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