Dune Core Modules (2.6.0)

hierarchy.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_AMGHIERARCHY_HH
4#define DUNE_AMGHIERARCHY_HH
5
6#include <list>
7#include <memory>
8#include <limits>
9#include <algorithm>
10#include <tuple>
11#include "aggregates.hh"
12#include "graph.hh"
13#include "galerkin.hh"
14#include "renumberer.hh"
15#include "graphcreator.hh"
17#include <dune/common/unused.hh>
18#include <dune/common/timer.hh>
20#include <dune/istl/bvector.hh>
31
32namespace Dune
33{
34 namespace Amg
35 {
47 enum {
55 MAX_PROCESSES = 72000
56 };
57
65 template<typename T, typename A=std::allocator<T> >
67 {
68 public:
72 typedef T MemberType;
73
74 template<typename T1, typename T2>
75 class LevelIterator;
76
77 private:
81 struct Element
82 {
83 friend class LevelIterator<Hierarchy<T,A>, T>;
84 friend class LevelIterator<const Hierarchy<T,A>, const T>;
85
87 Element* coarser_;
88
90 Element* finer_;
91
93 MemberType* element_;
94
96 MemberType* redistributed_;
97 };
98 public:
99
103 typedef typename A::template rebind<Element>::other Allocator;
104
105 typedef typename ConstructionTraits<T>::Arguments Arguments;
106
112
120
125
129 Hierarchy(const Hierarchy& other);
134 void addCoarser(Arguments& args);
135
136 void addRedistributedOnCoarsest(Arguments& args);
137
142 void addFiner(Arguments& args);
143
150 template<class C, class T1>
152 : public BidirectionalIteratorFacade<LevelIterator<C,T1>,T1,T1&>
153 {
154 friend class LevelIterator<typename std::remove_const<C>::type,
155 typename std::remove_const<T1>::type >;
156 friend class LevelIterator<const typename std::remove_const<C>::type,
157 const typename std::remove_const<T1>::type >;
158
159 public:
162 : element_(0)
163 {}
164
165 LevelIterator(Element* element)
166 : element_(element)
167 {}
168
170 LevelIterator(const LevelIterator<typename std::remove_const<C>::type,
171 typename std::remove_const<T1>::type>& other)
172 : element_(other.element_)
173 {}
174
176 LevelIterator(const LevelIterator<const typename std::remove_const<C>::type,
177 const typename std::remove_const<T1>::type>& other)
178 : element_(other.element_)
179 {}
180
184 bool equals(const LevelIterator<typename std::remove_const<C>::type,
185 typename std::remove_const<T1>::type>& other) const
186 {
187 return element_ == other.element_;
188 }
189
193 bool equals(const LevelIterator<const typename std::remove_const<C>::type,
194 const typename std::remove_const<T1>::type>& other) const
195 {
196 return element_ == other.element_;
197 }
198
200 T1& dereference() const
201 {
202 return *(element_->element_);
203 }
204
207 {
208 element_ = element_->coarser_;
209 }
210
213 {
214 element_ = element_->finer_;
215 }
216
221 bool isRedistributed() const
222 {
223 return element_->redistributed_;
224 }
225
231 {
232 assert(element_->redistributed_);
233 return *element_->redistributed_;
234 }
235 void addRedistributed(T1* t)
236 {
237 element_->redistributed_ = t;
238 }
239
240 void deleteRedistributed()
241 {
242 element_->redistributed_ = nullptr;
243 }
244
245 private:
246 Element* element_;
247 };
248
251
254
260
266
267
273
279
284 std::size_t levels() const;
285
288
289 private:
291 Element* finest_;
293 Element* coarsest_;
295 Element* nonAllocated_;
297 Allocator allocator_;
299 int levels_;
300 };
301
308 template<class M, class PI, class A=std::allocator<M> >
310 {
311 public:
313 typedef M MatrixOperator;
314
316 typedef typename MatrixOperator::matrix_type Matrix;
317
320
322 typedef A Allocator;
323
326
329
332
334 typedef typename Allocator::template rebind<AggregatesMap*>::other AAllocator;
335
337 typedef std::list<AggregatesMap*,AAllocator> AggregatesMapList;
338
340 typedef RedistributeInformation<ParallelInformation> RedistributeInfoType;
341
343 typedef typename Allocator::template rebind<RedistributeInfoType>::other RILAllocator;
344
346 typedef std::list<RedistributeInfoType,RILAllocator> RedistributeInfoList;
347
353 MatrixHierarchy(const MatrixOperator& fineMatrix,
355
356
358
364 template<typename O, typename T>
365 void build(const T& criterion);
366
374 template<class F>
375 void recalculateGalerkin(const F& copyFlags);
376
381 template<class V, class BA, class TA>
382 void coarsenVector(Hierarchy<BlockVector<V,BA>, TA>& hierarchy) const;
383
389 template<class S, class TA>
390 void coarsenSmoother(Hierarchy<S,TA>& smoothers,
391 const typename SmootherTraits<S>::Arguments& args) const;
392
397 std::size_t levels() const;
398
403 std::size_t maxlevels() const;
404
405 bool hasCoarsest() const;
406
411 bool isBuilt() const;
412
417 const ParallelMatrixHierarchy& matrices() const;
418
424
429 const AggregatesMapList& aggregatesMaps() const;
430
437
438 double getProlongationDampingFactor() const
439 {
440 return prolongDamp_;
441 }
442
453 void getCoarsestAggregatesOnFinest(std::vector<std::size_t>& data) const;
454
455 private:
456 typedef typename ConstructionTraits<MatrixOperator>::Arguments MatrixArgs;
457 typedef typename ConstructionTraits<ParallelInformation>::Arguments CommunicationArgs;
459 AggregatesMapList aggregatesMaps_;
461 RedistributeInfoList redistributes_;
463 ParallelMatrixHierarchy matrices_;
465 ParallelInformationHierarchy parallelInformation_;
466
468 bool built_;
469
471 int maxlevels_;
472
473 double prolongDamp_;
474
478 template<class Matrix, bool print>
479 struct MatrixStats
480 {
481
485 static void stats(const Matrix& matrix)
486 {
487 DUNE_UNUSED_PARAMETER(matrix);
488 }
489 };
490
491 template<class Matrix>
492 struct MatrixStats<Matrix,true>
493 {
494 struct calc
495 {
496 typedef typename Matrix::size_type size_type;
497 typedef typename Matrix::row_type matrix_row;
498
499 calc()
500 {
501 min=std::numeric_limits<size_type>::max();
502 max=0;
503 sum=0;
504 }
505
506 void operator()(const matrix_row& row)
507 {
508 min=std::min(min, row.size());
509 max=std::max(max, row.size());
510 sum += row.size();
511 }
512
513 size_type min;
514 size_type max;
515 size_type sum;
516 };
520 static void stats(const Matrix& matrix)
521 {
522 calc c= for_each(matrix.begin(), matrix.end(), calc());
523 dinfo<<"Matrix row: min="<<c.min<<" max="<<c.max
524 <<" average="<<static_cast<double>(c.sum)/matrix.N()
525 <<std::endl;
526 }
527 };
528 };
529
533 template<class T>
534 class CoarsenCriterion : public T
535 {
536 public:
542
553 CoarsenCriterion(int maxLevel=100, int coarsenTarget=1000, double minCoarsenRate=1.2,
554 double prolongDamp=1.6, AccumulationMode accumulate=successiveAccu)
555 : AggregationCriterion(Dune::Amg::Parameters(maxLevel, coarsenTarget, minCoarsenRate, prolongDamp, accumulate))
556 {}
557
559 : AggregationCriterion(parms)
560 {}
561
562 };
563
564 template<typename M, typename C1>
565 bool repartitionAndDistributeMatrix(const M& origMatrix, M& newMatrix,
566 SequentialInformation& origSequentialInformationomm,
567 SequentialInformation*& newComm,
568 RedistributeInformation<SequentialInformation>& ri,
569 int nparts, C1& criterion)
570 {
571 DUNE_UNUSED_PARAMETER(origMatrix);
572 DUNE_UNUSED_PARAMETER(newMatrix);
573 DUNE_UNUSED_PARAMETER(origSequentialInformationomm);
574 DUNE_UNUSED_PARAMETER(newComm);
576 DUNE_UNUSED_PARAMETER(nparts);
577 DUNE_UNUSED_PARAMETER(criterion);
578 DUNE_THROW(NotImplemented, "Redistribution does not make sense in sequential code!");
579 }
580
581
582 template<typename M, typename C, typename C1>
583 bool repartitionAndDistributeMatrix(const M& origMatrix, M& newMatrix, C& origComm, C*& newComm,
584 RedistributeInformation<C>& ri,
585 int nparts, C1& criterion)
586 {
587 Timer time;
588#ifdef AMG_REPART_ON_COMM_GRAPH
589 // Done not repartition the matrix graph, but a graph of the communication scheme.
590 bool existentOnRedist=Dune::commGraphRepartition(origMatrix, origComm, nparts, newComm,
591 ri.getInterface(),
592 criterion.debugLevel()>1);
593
594#else
595 typedef Dune::Amg::MatrixGraph<const M> MatrixGraph;
596 typedef Dune::Amg::PropertiesGraph<MatrixGraph,
597 VertexProperties,
598 EdgeProperties,
599 IdentityMap,
600 IdentityMap> PropertiesGraph;
601 MatrixGraph graph(origMatrix);
602 PropertiesGraph pgraph(graph);
603 buildDependency(pgraph, origMatrix, criterion, false);
604
605#ifdef DEBUG_REPART
606 if(origComm.communicator().rank()==0)
607 std::cout<<"Original matrix"<<std::endl;
608 origComm.communicator().barrier();
609 printGlobalSparseMatrix(origMatrix, origComm, std::cout);
610#endif
611 bool existentOnRedist=Dune::graphRepartition(pgraph, origComm, nparts,
612 newComm, ri.getInterface(),
613 criterion.debugLevel()>1);
614#endif // if else AMG_REPART
615
616 if(origComm.communicator().rank()==0 && criterion.debugLevel()>1)
617 std::cout<<"Repartitioning took "<<time.elapsed()<<" seconds."<<std::endl;
618
619 ri.setSetup();
620
621#ifdef DEBUG_REPART
622 ri.checkInterface(origComm.indexSet(), newComm->indexSet(), origComm.communicator());
623#endif
624
625 redistributeMatrix(const_cast<M&>(origMatrix), newMatrix, origComm, *newComm, ri);
626
627#ifdef DEBUG_REPART
628 if(origComm.communicator().rank()==0)
629 std::cout<<"Original matrix"<<std::endl;
630 origComm.communicator().barrier();
631 if(newComm->communicator().size()>0)
632 printGlobalSparseMatrix(newMatrix, *newComm, std::cout);
633 origComm.communicator().barrier();
634#endif
635
636 if(origComm.communicator().rank()==0 && criterion.debugLevel()>1)
637 std::cout<<"Redistributing matrix took "<<time.elapsed()<<" seconds."<<std::endl;
638 return existentOnRedist;
639
640 }
641
642 template<typename M>
643 bool repartitionAndDistributeMatrix(M& origMatrix, M& newMatrix,
644 SequentialInformation& origComm,
645 SequentialInformation& newComm,
646 RedistributeInformation<SequentialInformation>& ri)
647 {
648 return true;
649 }
650
651 template<class M, class IS, class A>
653 const ParallelInformation& pinfo)
654 : matrices_(const_cast<MatrixOperator&>(fineOperator)),
655 parallelInformation_(const_cast<ParallelInformation&>(pinfo))
656 {
657 if (SolverCategory::category(fineOperator) != pinfo.getSolverCategory())
658 DUNE_THROW(ISTLError, "MatrixOperator and ParallelInformation must belong to the same category!");
659 }
660
661 template<class M, class IS, class A>
662 template<typename O, typename T>
663 void MatrixHierarchy<M,IS,A>::build(const T& criterion)
664 {
665 prolongDamp_ = criterion.getProlongationDampingFactor();
666 typedef O OverlapFlags;
667 typedef typename ParallelMatrixHierarchy::Iterator MatIterator;
668 typedef typename ParallelInformationHierarchy::Iterator PInfoIterator;
669
670 static const int noints=(Dune::Amg::MAX_PROCESSES/4096>0) ? (Dune::Amg::MAX_PROCESSES/4096) : 1;
671
672 typedef bigunsignedint<sizeof(int)*8*noints> BIGINT;
673 GalerkinProduct<ParallelInformation> productBuilder;
674 MatIterator mlevel = matrices_.finest();
675 MatrixStats<typename M::matrix_type,MINIMAL_DEBUG_LEVEL<=INFO_DEBUG_LEVEL>::stats(mlevel->getmat());
676
677 PInfoIterator infoLevel = parallelInformation_.finest();
678 BIGINT finenonzeros=countNonZeros(mlevel->getmat());
679 finenonzeros = infoLevel->communicator().sum(finenonzeros);
680 BIGINT allnonzeros = finenonzeros;
681
682
683 int level = 0;
684 int rank = 0;
685
686 BIGINT unknowns = mlevel->getmat().N();
687
688 unknowns = infoLevel->communicator().sum(unknowns);
689 double dunknowns=unknowns.todouble();
690 infoLevel->buildGlobalLookup(mlevel->getmat().N());
691 redistributes_.push_back(RedistributeInfoType());
692
693 for(; level < criterion.maxLevel(); ++level, ++mlevel) {
694 assert(matrices_.levels()==redistributes_.size());
695 rank = infoLevel->communicator().rank();
696 if(rank==0 && criterion.debugLevel()>1)
697 std::cout<<"Level "<<level<<" has "<<dunknowns<<" unknowns, "<<dunknowns/infoLevel->communicator().size()
698 <<" unknowns per proc (procs="<<infoLevel->communicator().size()<<")"<<std::endl;
699
700 MatrixOperator* matrix=&(*mlevel);
701 ParallelInformation* info =&(*infoLevel);
702
703 if((
704#if HAVE_PARMETIS
705 criterion.accumulate()==successiveAccu
706#else
707 false
708#endif
709 || (criterion.accumulate()==atOnceAccu
710 && dunknowns < 30*infoLevel->communicator().size()))
711 && infoLevel->communicator().size()>1 &&
712 dunknowns/infoLevel->communicator().size() <= criterion.coarsenTarget())
713 {
714 // accumulate to fewer processors
715 Matrix* redistMat= new Matrix();
716 ParallelInformation* redistComm=0;
717 std::size_t nodomains = (std::size_t)std::ceil(dunknowns/(criterion.minAggregateSize()
718 *criterion.coarsenTarget()));
719 if( nodomains<=criterion.minAggregateSize()/2 ||
720 dunknowns <= criterion.coarsenTarget() )
721 nodomains=1;
722
723 bool existentOnNextLevel =
724 repartitionAndDistributeMatrix(mlevel->getmat(), *redistMat, *infoLevel,
725 redistComm, redistributes_.back(), nodomains,
726 criterion);
727 BIGINT unknownsRedist = redistMat->N();
728 unknownsRedist = infoLevel->communicator().sum(unknownsRedist);
729 dunknowns= unknownsRedist.todouble();
730 if(redistComm->communicator().rank()==0 && criterion.debugLevel()>1)
731 std::cout<<"Level "<<level<<" (redistributed) has "<<dunknowns<<" unknowns, "<<dunknowns/redistComm->communicator().size()
732 <<" unknowns per proc (procs="<<redistComm->communicator().size()<<")"<<std::endl;
733 MatrixArgs args(*redistMat, *redistComm);
734 mlevel.addRedistributed(ConstructionTraits<MatrixOperator>::construct(args));
735 assert(mlevel.isRedistributed());
736 infoLevel.addRedistributed(redistComm);
737 infoLevel->freeGlobalLookup();
738
739 if(!existentOnNextLevel)
740 // We do not hold any data on the redistributed partitioning
741 break;
742
743 // Work on the redistributed Matrix from now on
744 matrix = &(mlevel.getRedistributed());
745 info = &(infoLevel.getRedistributed());
746 info->buildGlobalLookup(matrix->getmat().N());
747 }
748
749 rank = info->communicator().rank();
750 if(dunknowns <= criterion.coarsenTarget())
751 // No further coarsening needed
752 break;
753
754 typedef PropertiesGraphCreator<MatrixOperator,ParallelInformation> GraphCreator;
755 typedef typename GraphCreator::PropertiesGraph PropertiesGraph;
756 typedef typename GraphCreator::GraphTuple GraphTuple;
757
758 typedef typename PropertiesGraph::VertexDescriptor Vertex;
759
760 std::vector<bool> excluded(matrix->getmat().N(), false);
761
762 GraphTuple graphs = GraphCreator::create(*matrix, excluded, *info, OverlapFlags());
763
764 AggregatesMap* aggregatesMap=new AggregatesMap(std::get<1>(graphs)->maxVertex()+1);
765
766 aggregatesMaps_.push_back(aggregatesMap);
767
768 Timer watch;
769 watch.reset();
770 int noAggregates, isoAggregates, oneAggregates, skippedAggregates;
771
772 std::tie(noAggregates, isoAggregates, oneAggregates, skippedAggregates) =
773 aggregatesMap->buildAggregates(matrix->getmat(), *(std::get<1>(graphs)), criterion, level==0);
774
775 if(rank==0 && criterion.debugLevel()>2)
776 std::cout<<" Have built "<<noAggregates<<" aggregates totally ("<<isoAggregates<<" isolated aggregates, "<<
777 oneAggregates<<" aggregates of one vertex, and skipped "<<
778 skippedAggregates<<" aggregates)."<<std::endl;
779#ifdef TEST_AGGLO
780 {
781 // calculate size of local matrix in the distributed direction
782 int start, end, overlapStart, overlapEnd;
783 int procs=info->communicator().rank();
784 int n = UNKNOWNS/procs; // number of unknowns per process
785 int bigger = UNKNOWNS%procs; // number of process with n+1 unknows
786
787 // Compute owner region
788 if(rank<bigger) {
789 start = rank*(n+1);
790 end = (rank+1)*(n+1);
791 }else{
792 start = bigger + rank * n;
793 end = bigger + (rank + 1) * n;
794 }
795
796 // Compute overlap region
797 if(start>0)
798 overlapStart = start - 1;
799 else
800 overlapStart = start;
801
802 if(end<UNKNOWNS)
803 overlapEnd = end + 1;
804 else
805 overlapEnd = end;
806
807 assert((UNKNOWNS)*(overlapEnd-overlapStart)==aggregatesMap->noVertices());
808 for(int j=0; j< UNKNOWNS; ++j)
809 for(int i=0; i < UNKNOWNS; ++i)
810 {
811 if(i>=overlapStart && i<overlapEnd)
812 {
813 int no = (j/2)*((UNKNOWNS)/2)+i/2;
814 (*aggregatesMap)[j*(overlapEnd-overlapStart)+i-overlapStart]=no;
815 }
816 }
817 }
818#endif
819 if(criterion.debugLevel()>1 && info->communicator().rank()==0)
820 std::cout<<"aggregating finished."<<std::endl;
821
822 BIGINT gnoAggregates=noAggregates;
823 gnoAggregates = info->communicator().sum(gnoAggregates);
824 double dgnoAggregates = gnoAggregates.todouble();
825#ifdef TEST_AGGLO
826 BIGINT gnoAggregates=((UNKNOWNS)/2)*((UNKNOWNS)/2);
827#endif
828
829 if(criterion.debugLevel()>2 && rank==0)
830 std::cout << "Building "<<dgnoAggregates<<" aggregates took "<<watch.elapsed()<<" seconds."<<std::endl;
831
832 if(dgnoAggregates==0 || dunknowns/dgnoAggregates<criterion.minCoarsenRate())
833 {
834 if(rank==0)
835 {
836 if(dgnoAggregates>0)
837 std::cerr << "Stopped coarsening because of rate breakdown "<<dunknowns<<"/"<<dgnoAggregates
838 <<"="<<dunknowns/dgnoAggregates<<"<"
839 <<criterion.minCoarsenRate()<<std::endl;
840 else
841 std::cerr<< "Could not build any aggregates. Probably no connected nodes."<<std::endl;
842 }
843 aggregatesMap->free();
844 delete aggregatesMap;
845 aggregatesMaps_.pop_back();
846
847 if(criterion.accumulate() && mlevel.isRedistributed() && info->communicator().size()>1) {
848 // coarse level matrix was already redistributed, but to more than 1 process
849 // Therefore need to delete the redistribution. Further down it will
850 // then be redistributed to 1 process
851 delete &(mlevel.getRedistributed().getmat());
852 mlevel.deleteRedistributed();
853 delete &(infoLevel.getRedistributed());
854 infoLevel.deleteRedistributed();
855 redistributes_.back().resetSetup();
856 }
857
858 break;
859 }
860 unknowns = noAggregates;
861 dunknowns = dgnoAggregates;
862
863 CommunicationArgs commargs(info->communicator(),info->getSolverCategory());
864 parallelInformation_.addCoarser(commargs);
865
866 ++infoLevel; // parallel information on coarse level
867
869 get(VertexVisitedTag(), *(std::get<1>(graphs)));
870
871 watch.reset();
872 int aggregates = IndicesCoarsener<ParallelInformation,OverlapFlags>
873 ::coarsen(*info,
874 *(std::get<1>(graphs)),
875 visitedMap,
876 *aggregatesMap,
877 *infoLevel,
878 noAggregates);
879 GraphCreator::free(graphs);
880
881 if(criterion.debugLevel()>2) {
882 if(rank==0)
883 std::cout<<"Coarsening of index sets took "<<watch.elapsed()<<" seconds."<<std::endl;
884 }
885
886 watch.reset();
887
888 infoLevel->buildGlobalLookup(aggregates);
889 AggregatesPublisher<Vertex,OverlapFlags,ParallelInformation>::publish(*aggregatesMap,
890 *info,
891 infoLevel->globalLookup());
892
893
894 if(criterion.debugLevel()>2) {
895 if(rank==0)
896 std::cout<<"Communicating global aggregate numbers took "<<watch.elapsed()<<" seconds."<<std::endl;
897 }
898
899 watch.reset();
900 std::vector<bool>& visited=excluded;
901
902 typedef std::vector<bool>::iterator Iterator;
904 Iterator end = visited.end();
905 for(Iterator iter= visited.begin(); iter != end; ++iter)
906 *iter=false;
907
908 VisitedMap2 visitedMap2(visited.begin(), Dune::IdentityMap());
909
910 typename MatrixOperator::matrix_type* coarseMatrix;
911
912 coarseMatrix = productBuilder.build(*(std::get<0>(graphs)), visitedMap2,
913 *info,
914 *aggregatesMap,
915 aggregates,
916 OverlapFlags());
917 dverb<<"Building of sparsity pattern took "<<watch.elapsed()<<std::endl;
918 watch.reset();
919 info->freeGlobalLookup();
920
921 delete std::get<0>(graphs);
922 productBuilder.calculate(matrix->getmat(), *aggregatesMap, *coarseMatrix, *infoLevel, OverlapFlags());
923
924 if(criterion.debugLevel()>2) {
925 if(rank==0)
926 std::cout<<"Calculation entries of Galerkin product took "<<watch.elapsed()<<" seconds."<<std::endl;
927 }
928
929 BIGINT nonzeros = countNonZeros(*coarseMatrix);
930 allnonzeros = allnonzeros + infoLevel->communicator().sum(nonzeros);
931 MatrixArgs args(*coarseMatrix, *infoLevel);
932
933 matrices_.addCoarser(args);
934 redistributes_.push_back(RedistributeInfoType());
935 } // end level loop
936
937
938 infoLevel->freeGlobalLookup();
939
940 built_=true;
941 AggregatesMap* aggregatesMap=new AggregatesMap(0);
942 aggregatesMaps_.push_back(aggregatesMap);
943
944 if(criterion.debugLevel()>0) {
945 if(level==criterion.maxLevel()) {
946 BIGINT unknownsLevel = mlevel->getmat().N();
947 unknownsLevel = infoLevel->communicator().sum(unknownsLevel);
948 double dunknownsLevel = unknownsLevel.todouble();
949 if(rank==0 && criterion.debugLevel()>1) {
950 std::cout<<"Level "<<level<<" has "<<dunknownsLevel<<" unknowns, "<<dunknownsLevel/infoLevel->communicator().size()
951 <<" unknowns per proc (procs="<<infoLevel->communicator().size()<<")"<<std::endl;
952 }
953 }
954 }
955
956 if(criterion.accumulate() && !redistributes_.back().isSetup() &&
957 infoLevel->communicator().size()>1) {
958#if HAVE_MPI && !HAVE_PARMETIS
959 if(criterion.accumulate()==successiveAccu &&
960 infoLevel->communicator().rank()==0)
961 std::cerr<<"Successive accumulation of data on coarse levels only works with ParMETIS installed."
962 <<" Fell back to accumulation to one domain on coarsest level"<<std::endl;
963#endif
964
965 // accumulate to fewer processors
966 Matrix* redistMat= new Matrix();
967 ParallelInformation* redistComm=0;
968 int nodomains = 1;
969
970 repartitionAndDistributeMatrix(mlevel->getmat(), *redistMat, *infoLevel,
971 redistComm, redistributes_.back(), nodomains,criterion);
972 MatrixArgs args(*redistMat, *redistComm);
973 BIGINT unknownsRedist = redistMat->N();
974 unknownsRedist = infoLevel->communicator().sum(unknownsRedist);
975
976 if(redistComm->communicator().rank()==0 && criterion.debugLevel()>1) {
977 double dunknownsRedist = unknownsRedist.todouble();
978 std::cout<<"Level "<<level<<" redistributed has "<<dunknownsRedist<<" unknowns, "<<dunknownsRedist/redistComm->communicator().size()
979 <<" unknowns per proc (procs="<<redistComm->communicator().size()<<")"<<std::endl;
980 }
981 mlevel.addRedistributed(ConstructionTraits<MatrixOperator>::construct(args));
982 infoLevel.addRedistributed(redistComm);
983 infoLevel->freeGlobalLookup();
984 }
985
986 int levels = matrices_.levels();
987 maxlevels_ = parallelInformation_.finest()->communicator().max(levels);
988 assert(matrices_.levels()==redistributes_.size());
989 if(hasCoarsest() && rank==0 && criterion.debugLevel()>1)
990 std::cout<<"operator complexity: "<<allnonzeros.todouble()/finenonzeros.todouble()<<std::endl;
991
992 }
993
994 template<class M, class IS, class A>
997 {
998 return matrices_;
999 }
1000
1001 template<class M, class IS, class A>
1004 {
1005 return parallelInformation_;
1006 }
1007
1008 template<class M, class IS, class A>
1009 void MatrixHierarchy<M,IS,A>::getCoarsestAggregatesOnFinest(std::vector<std::size_t>& data) const
1010 {
1011 int levels=aggregatesMaps().size();
1012 int maxlevels=parallelInformation_.finest()->communicator().max(levels);
1013 std::size_t size=(*(aggregatesMaps().begin()))->noVertices();
1014 // We need an auxiliary vector for the consecutive prolongation.
1015 std::vector<std::size_t> tmp;
1016 std::vector<std::size_t> *coarse, *fine;
1017
1018 // make sure the allocated space suffices.
1019 tmp.reserve(size);
1020 data.reserve(size);
1021
1022 // Correctly assign coarse and fine for the first prolongation such that
1023 // we end up in data in the end.
1024 if(levels%2==0) {
1025 coarse=&tmp;
1026 fine=&data;
1027 }else{
1028 coarse=&data;
1029 fine=&tmp;
1030 }
1031
1032 // Number the unknowns on the coarsest level consecutively for each process.
1033 if(levels==maxlevels) {
1034 const AggregatesMap& map = *(*(++aggregatesMaps().rbegin()));
1035 std::size_t m=0;
1036
1037 for(typename AggregatesMap::const_iterator iter = map.begin(); iter != map.end(); ++iter)
1038 if(*iter< AggregatesMap::ISOLATED)
1039 m=std::max(*iter,m);
1040
1041 coarse->resize(m+1);
1042 std::size_t i=0;
1043 srand((unsigned)std::clock());
1044 std::set<size_t> used;
1045 for(typename std::vector<std::size_t>::iterator iter=coarse->begin(); iter != coarse->end();
1046 ++iter, ++i)
1047 {
1048 std::pair<std::set<std::size_t>::iterator,bool> ibpair
1049 = used.insert(static_cast<std::size_t>((((double)rand())/(RAND_MAX+1.0)))*coarse->size());
1050
1051 while(!ibpair.second)
1052 ibpair = used.insert(static_cast<std::size_t>((((double)rand())/(RAND_MAX+1.0))*coarse->size()));
1053 *iter=*(ibpair.first);
1054 }
1055 }
1056
1057 typename ParallelInformationHierarchy::Iterator pinfo = parallelInformation().coarsest();
1058 --pinfo;
1059
1060 // Now consecutively project the numbers to the finest level.
1061 for(typename AggregatesMapList::const_reverse_iterator aggregates=++aggregatesMaps().rbegin();
1062 aggregates != aggregatesMaps().rend(); ++aggregates,--levels) {
1063
1064 fine->resize((*aggregates)->noVertices());
1065 fine->assign(fine->size(), 0);
1066 Transfer<typename AggregatesMap::AggregateDescriptor, std::vector<std::size_t>, ParallelInformation>
1067 ::prolongateVector(*(*aggregates), *coarse, *fine, static_cast<std::size_t>(1), *pinfo);
1068 --pinfo;
1069 std::swap(coarse, fine);
1070 }
1071
1072 // Assertion to check that we really projected to data on the last step.
1073 assert(coarse==&data);
1074 }
1075
1076 template<class M, class IS, class A>
1079 {
1080 return aggregatesMaps_;
1081 }
1082 template<class M, class IS, class A>
1085 {
1086 return redistributes_;
1087 }
1088
1089 template<class M, class IS, class A>
1091 {
1092 typedef typename AggregatesMapList::reverse_iterator AggregatesMapIterator;
1093 typedef typename ParallelMatrixHierarchy::Iterator Iterator;
1094 typedef typename ParallelInformationHierarchy::Iterator InfoIterator;
1095
1096 AggregatesMapIterator amap = aggregatesMaps_.rbegin();
1097 InfoIterator info = parallelInformation_.coarsest();
1098 for(Iterator level=matrices_.coarsest(), finest=matrices_.finest(); level != finest; --level, --info, ++amap) {
1099 (*amap)->free();
1100 delete *amap;
1101 delete &level->getmat();
1102 if(level.isRedistributed())
1103 delete &(level.getRedistributed().getmat());
1104 }
1105 delete *amap;
1106 }
1107
1108 template<class M, class IS, class A>
1109 template<class V, class BA, class TA>
1111 {
1112 assert(hierarchy.levels()==1);
1113 typedef typename ParallelMatrixHierarchy::ConstIterator Iterator;
1114 typedef typename RedistributeInfoList::const_iterator RIter;
1115 RIter redist = redistributes_.begin();
1116
1117 Iterator matrix = matrices_.finest(), coarsest = matrices_.coarsest();
1118 int level=0;
1119 if(redist->isSetup())
1120 hierarchy.addRedistributedOnCoarsest(matrix.getRedistributed().getmat().N());
1121 Dune::dvverb<<"Level "<<level<<" has "<<matrices_.finest()->getmat().N()<<" unknowns!"<<std::endl;
1122
1123 while(matrix != coarsest) {
1124 ++matrix; ++level; ++redist;
1125 Dune::dvverb<<"Level "<<level<<" has "<<matrix->getmat().N()<<" unknowns!"<<std::endl;
1126
1127 hierarchy.addCoarser(matrix->getmat().N());
1128 if(redist->isSetup())
1129 hierarchy.addRedistributedOnCoarsest(matrix.getRedistributed().getmat().N());
1130
1131 }
1132
1133 }
1134
1135 template<class M, class IS, class A>
1136 template<class S, class TA>
1138 const typename SmootherTraits<S>::Arguments& sargs) const
1139 {
1140 assert(smoothers.levels()==0);
1141 typedef typename ParallelMatrixHierarchy::ConstIterator MatrixIterator;
1142 typedef typename ParallelInformationHierarchy::ConstIterator PinfoIterator;
1143 typedef typename AggregatesMapList::const_iterator AggregatesIterator;
1144
1145 typename ConstructionTraits<S>::Arguments cargs;
1146 cargs.setArgs(sargs);
1147 PinfoIterator pinfo = parallelInformation_.finest();
1148 AggregatesIterator aggregates = aggregatesMaps_.begin();
1149 int level=0;
1150 for(MatrixIterator matrix = matrices_.finest(), coarsest = matrices_.coarsest();
1151 matrix != coarsest; ++matrix, ++pinfo, ++aggregates, ++level) {
1152 cargs.setMatrix(matrix->getmat(), **aggregates);
1153 cargs.setComm(*pinfo);
1154 smoothers.addCoarser(cargs);
1155 }
1156 if(maxlevels()>levels()) {
1157 // This is not the globally coarsest level and therefore smoothing is needed
1158 cargs.setMatrix(matrices_.coarsest()->getmat(), **aggregates);
1159 cargs.setComm(*pinfo);
1160 smoothers.addCoarser(cargs);
1161 ++level;
1162 }
1163 }
1164
1165 template<class M, class IS, class A>
1166 template<class F>
1168 {
1169 typedef typename AggregatesMapList::iterator AggregatesMapIterator;
1170 typedef typename ParallelMatrixHierarchy::Iterator Iterator;
1171 typedef typename ParallelInformationHierarchy::Iterator InfoIterator;
1172
1173 AggregatesMapIterator amap = aggregatesMaps_.begin();
1174 BaseGalerkinProduct productBuilder;
1175 InfoIterator info = parallelInformation_.finest();
1176 typename RedistributeInfoList::iterator riIter = redistributes_.begin();
1177 Iterator level = matrices_.finest(), coarsest=matrices_.coarsest();
1178 if(level.isRedistributed()) {
1179 info->buildGlobalLookup(level->getmat().N());
1180 redistributeMatrixEntries(const_cast<Matrix&>(level->getmat()),
1181 const_cast<Matrix&>(level.getRedistributed().getmat()),
1182 *info,info.getRedistributed(), *riIter);
1183 info->freeGlobalLookup();
1184 }
1185
1186 for(; level!=coarsest; ++amap) {
1187 const Matrix& fine = (level.isRedistributed() ? level.getRedistributed() : *level).getmat();
1188 ++level;
1189 ++info;
1190 ++riIter;
1191 productBuilder.calculate(fine, *(*amap), const_cast<Matrix&>(level->getmat()), *info, copyFlags);
1192 if(level.isRedistributed()) {
1193 info->buildGlobalLookup(level->getmat().N());
1194 redistributeMatrixEntries(const_cast<Matrix&>(level->getmat()),
1195 const_cast<Matrix&>(level.getRedistributed().getmat()), *info,
1196 info.getRedistributed(), *riIter);
1197 info->freeGlobalLookup();
1198 }
1199 }
1200 }
1201
1202 template<class M, class IS, class A>
1204 {
1205 return matrices_.levels();
1206 }
1207
1208 template<class M, class IS, class A>
1210 {
1211 return maxlevels_;
1212 }
1213
1214 template<class M, class IS, class A>
1216 {
1217 return levels()==maxlevels() &&
1218 (!matrices_.coarsest().isRedistributed() ||matrices_.coarsest()->getmat().N()>0);
1219 }
1220
1221 template<class M, class IS, class A>
1223 {
1224 return built_;
1225 }
1226
1227 template<class T, class A>
1229 : finest_(0), coarsest_(0), nonAllocated_(0), allocator_(), levels_(0)
1230 {}
1231
1232 template<class T, class A>
1234 : allocator_()
1235 {
1236 finest_ = allocator_.allocate(1,0);
1237 finest_->element_ = &first;
1238 finest_->redistributed_ = nullptr;
1239 nonAllocated_ = finest_;
1240 coarsest_ = finest_;
1241 coarsest_->coarser_ = coarsest_->finer_ = nullptr;
1242 levels_ = 1;
1243 }
1244
1245 template<class T, class A>
1247 : allocator_()
1248 {
1249 finest_ = allocator_.allocate(1,0);
1250 finest_->element_ = first;
1251 finest_->redistributed_ = nullptr;
1252 nonAllocated_ = nullptr;
1253 coarsest_ = finest_;
1254 coarsest_->coarser_ = coarsest_->finer_ = nullptr;
1255 levels_ = 1;
1256 }
1257 template<class T, class A>
1259 {
1260 while(coarsest_) {
1261 Element* current = coarsest_;
1262 coarsest_ = coarsest_->finer_;
1263 if(current != nonAllocated_) {
1264 if(current->redistributed_)
1265 ConstructionTraits<T>::deconstruct(current->redistributed_);
1266 ConstructionTraits<T>::deconstruct(current->element_);
1267 }
1268 allocator_.deallocate(current, 1);
1269 current=nullptr;
1270 //coarsest_->coarser_ = nullptr;
1271 }
1272 }
1273
1274 template<class T, class A>
1276 : nonAllocated_(), allocator_(other.allocator_),
1277 levels_(other.levels_)
1278 {
1279 if(!other.finest_)
1280 {
1281 finest_=coarsest_=nonAllocated_=nullptr;
1282 return;
1283 }
1284 finest_=allocator_.allocate(1,0);
1285 Element* finer_ = nullptr;
1286 Element* current_ = finest_;
1287 Element* otherCurrent_ = other.finest_;
1288
1289 while(otherCurrent_)
1290 {
1291 T* t=new T(*(otherCurrent_->element_));
1292 current_->element_=t;
1293 current_->finer_=finer_;
1294 if(otherCurrent_->redistributed_)
1295 current_->redistributed_ = new T(*otherCurrent_->redistributed_);
1296 else
1297 current_->redistributed_= nullptr;
1298 finer_=current_;
1299 if(otherCurrent_->coarser_)
1300 {
1301 current_->coarser_=allocator_.allocate(1,0);
1302 current_=current_->coarser_;
1303 }else
1304 current_->coarser_=nullptr;
1305 otherCurrent_=otherCurrent_->coarser_;
1306 }
1307 coarsest_=current_;
1308 }
1309
1310 template<class T, class A>
1311 std::size_t Hierarchy<T,A>::levels() const
1312 {
1313 return levels_;
1314 }
1315
1316 template<class T, class A>
1317 void Hierarchy<T,A>::addRedistributedOnCoarsest(Arguments& args)
1318 {
1319 coarsest_->redistributed_ = ConstructionTraits<MemberType>::construct(args);
1320 }
1321
1322 template<class T, class A>
1323 void Hierarchy<T,A>::addCoarser(Arguments& args)
1324 {
1325 if(!coarsest_) {
1326 assert(!finest_);
1327 coarsest_ = allocator_.allocate(1,0);
1328 coarsest_->element_ = ConstructionTraits<MemberType>::construct(args);
1329 finest_ = coarsest_;
1330 coarsest_->finer_ = nullptr;
1331 }else{
1332 coarsest_->coarser_ = allocator_.allocate(1,0);
1333 coarsest_->coarser_->finer_ = coarsest_;
1334 coarsest_ = coarsest_->coarser_;
1335 coarsest_->element_ = ConstructionTraits<MemberType>::construct(args);
1336 }
1337 coarsest_->redistributed_ = nullptr;
1338 coarsest_->coarser_=nullptr;
1339 ++levels_;
1340 }
1341
1342
1343 template<class T, class A>
1344 void Hierarchy<T,A>::addFiner(Arguments& args)
1345 {
1346 if(!finest_) {
1347 assert(!coarsest_);
1348 finest_ = allocator_.allocate(1,0);
1349 finest_->element = ConstructionTraits<T>::construct(args);
1350 coarsest_ = finest_;
1351 coarsest_->coarser_ = coarsest_->finer_ = nullptr;
1352 }else{
1353 finest_->finer_ = allocator_.allocate(1,0);
1354 finest_->finer_->coarser_ = finest_;
1355 finest_ = finest_->finer_;
1356 finest_->finer = nullptr;
1357 finest_->element = ConstructionTraits<T>::construct(args);
1358 }
1359 ++levels_;
1360 }
1361
1362 template<class T, class A>
1364 {
1365 return Iterator(finest_);
1366 }
1367
1368 template<class T, class A>
1370 {
1371 return Iterator(coarsest_);
1372 }
1373
1374 template<class T, class A>
1376 {
1377 return ConstIterator(finest_);
1378 }
1379
1380 template<class T, class A>
1382 {
1383 return ConstIterator(coarsest_);
1384 }
1386 } // namespace Amg
1387} // namespace Dune
1388
1389#endif
Provides classes for the Coloring process of AMG.
Portable very large unsigned integers.
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Class providing information about the mapping of the vertices onto aggregates.
Definition: aggregates.hh:543
Base class of all aggregation criterions.
Definition: aggregates.hh:47
The criterion describing the stop criteria for the coarsening process.
Definition: hierarchy.hh:535
T AggregationCriterion
The criterion for tagging connections as strong and nodes as isolated. This might be e....
Definition: hierarchy.hh:541
CoarsenCriterion(int maxLevel=100, int coarsenTarget=1000, double minCoarsenRate=1.2, double prolongDamp=1.6, AccumulationMode accumulate=successiveAccu)
Constructor.
Definition: hierarchy.hh:553
Traits class for generically constructing non default constructable types.
Definition: construction.hh:38
Iterator over the levels in the hierarchy.
Definition: hierarchy.hh:153
LevelIterator(const LevelIterator< typename std::remove_const< C >::type, typename std::remove_const< T1 >::type > &other)
Copy constructor.
Definition: hierarchy.hh:170
T1 & dereference() const
Dereference the iterator.
Definition: hierarchy.hh:200
bool equals(const LevelIterator< typename std::remove_const< C >::type, typename std::remove_const< T1 >::type > &other) const
Equality check.
Definition: hierarchy.hh:184
bool isRedistributed() const
Check whether there was a redistribution at the current level.
Definition: hierarchy.hh:221
bool equals(const LevelIterator< const typename std::remove_const< C >::type, const typename std::remove_const< T1 >::type > &other) const
Equality check.
Definition: hierarchy.hh:193
void increment()
Move to the next coarser level.
Definition: hierarchy.hh:206
LevelIterator(const LevelIterator< const typename std::remove_const< C >::type, const typename std::remove_const< T1 >::type > &other)
Copy constructor.
Definition: hierarchy.hh:176
void decrement()
Move to the next fine level.
Definition: hierarchy.hh:212
T1 & getRedistributed() const
Get the redistributed container.
Definition: hierarchy.hh:230
LevelIterator()
Constructor.
Definition: hierarchy.hh:161
A hierarchy of coantainers (e.g. matrices or vectors)
Definition: hierarchy.hh:67
T MemberType
The type of the container we store.
Definition: hierarchy.hh:72
LevelIterator< Hierarchy< T, A >, T > Iterator
Type of the mutable iterator.
Definition: hierarchy.hh:250
LevelIterator< const Hierarchy< T, A >, const T > ConstIterator
Type of the const iterator.
Definition: hierarchy.hh:253
A::template rebind< Element >::other Allocator
The allocator to use for the list elements.
Definition: hierarchy.hh:103
The hierarchies build by the coarsening process.
Definition: hierarchy.hh:310
Dune::Amg::Hierarchy< ParallelInformation, Allocator > ParallelInformationHierarchy
The type of the parallel informarion hierarchy.
Definition: hierarchy.hh:331
std::list< AggregatesMap *, AAllocator > AggregatesMapList
The type of the aggregates maps list.
Definition: hierarchy.hh:337
PI ParallelInformation
The type of the index set.
Definition: hierarchy.hh:319
Dune::Amg::Hierarchy< MatrixOperator, Allocator > ParallelMatrixHierarchy
The type of the parallel matrix hierarchy.
Definition: hierarchy.hh:328
A Allocator
The allocator to use.
Definition: hierarchy.hh:322
RedistributeInformation< ParallelInformation > RedistributeInfoType
The type of the redistribute information.
Definition: hierarchy.hh:340
std::list< RedistributeInfoType, RILAllocator > RedistributeInfoList
The type of the list of redistribute information.
Definition: hierarchy.hh:346
Dune::Amg::AggregatesMap< typename MatrixGraph< Matrix >::VertexDescriptor > AggregatesMap
The type of the aggregates map we use.
Definition: hierarchy.hh:325
Allocator::template rebind< AggregatesMap * >::other AAllocator
Allocator for pointers.
Definition: hierarchy.hh:334
MatrixOperator::matrix_type Matrix
The type of the matrix.
Definition: hierarchy.hh:316
Allocator::template rebind< RedistributeInfoType >::other RILAllocator
Allocator for RedistributeInfoType.
Definition: hierarchy.hh:343
M MatrixOperator
The type of the matrix operator.
Definition: hierarchy.hh:313
All parameters for AMG.
Definition: parameters.hh:391
Attaches properties to the edges and vertices of a graph.
Definition: graph.hh:976
Graph::VertexDescriptor VertexDescriptor
The vertex descriptor.
Definition: graph.hh:986
Facade class for stl conformant bidirectional iterators.
Definition: iteratorfacades.hh:275
A vector of blocks with memory management.
Definition: bvector.hh:317
derive error class from the base class in common
Definition: istlexception.hh:16
Adapter to turn a random access iterator into a property map.
Definition: propertymap.hh:106
A generic dynamic dense matrix.
Definition: matrix.hh:555
A::size_type size_type
Type for indices and sizes.
Definition: matrix.hh:571
MatrixImp::DenseMatrixBase< T, A >::window_type row_type
The type implementing a matrix row.
Definition: matrix.hh:568
A simple stop watch.
Definition: timer.hh:41
void reset() noexcept
Reset timer while keeping the running/stopped state.
Definition: timer.hh:55
double elapsed() const noexcept
Get elapsed user-time from last reset until now/last stop in seconds.
Definition: timer.hh:75
Portable very large unsigned integers.
Definition: bigunsignedint.hh:70
Helper classes for the construction of classes without empty constructor.
Provides classes for initializing the link attributes of a matrix graph.
Provides a class for building the galerkin product based on a aggregation scheme.
Provdes class for identifying aggregates globally.
Provides classes for building the matrix graph.
#define DUNE_UNUSED_PARAMETER(parm)
A macro to mark intentionally unused function parameters with.
Definition: unused.hh:25
#define DUNE_THROW(E, m)
Definition: exceptions.hh:216
T accumulate(Range &&range, T value, F &&f)
Accumulate values.
Definition: hybridutilities.hh:331
const AggregatesMapList & aggregatesMaps() const
Get the hierarchy of the mappings of the nodes onto aggregates.
Definition: hierarchy.hh:1078
bool isBuilt() const
Whether the hierarchy was built.
Definition: hierarchy.hh:1222
Hierarchy(const Hierarchy &other)
Copy constructor.
Definition: hierarchy.hh:1275
std::size_t levels() const
Get the number of levels in the hierarchy.
Definition: hierarchy.hh:1311
std::size_t levels() const
Get the number of levels in the hierarchy.
Definition: hierarchy.hh:1203
ConstIterator coarsest() const
Get an iterator positioned at the coarsest level.
Definition: hierarchy.hh:1381
void addCoarser(Arguments &args)
Add an element on a coarser level.
Definition: hierarchy.hh:1323
const RedistributeInfoList & redistributeInformation() const
Get the hierarchy of the information about redistributions,.
Definition: hierarchy.hh:1084
const ParallelInformationHierarchy & parallelInformation() const
Get the hierarchy of the parallel data distribution information.
Definition: hierarchy.hh:1003
void addFiner(Arguments &args)
Add an element on a finer level.
Definition: hierarchy.hh:1344
const ParallelMatrixHierarchy & matrices() const
Get the matrix hierarchy.
Definition: hierarchy.hh:996
std::size_t maxlevels() const
Get the max number of levels in the hierarchy of processors.
Definition: hierarchy.hh:1209
static const V ISOLATED
Identifier of isolated vertices.
Definition: aggregates.hh:554
void recalculateGalerkin(const F &copyFlags)
Recalculate the galerkin products.
Definition: hierarchy.hh:1167
const void * Arguments
A type holding all the arguments needed to call the constructor.
Definition: construction.hh:44
std::size_t noVertices() const
Get the number of vertices.
Iterator coarsest()
Get an iterator positioned at the coarsest level.
Definition: hierarchy.hh:1369
Hierarchy(MemberType &first)
Construct a new hierarchy.
Definition: hierarchy.hh:1233
ConstIterator finest() const
Get an iterator positioned at the finest level.
Definition: hierarchy.hh:1375
void coarsenVector(Hierarchy< BlockVector< V, BA >, TA > &hierarchy) const
Coarsen the vector hierarchy according to the matrix hierarchy.
Definition: hierarchy.hh:1110
Hierarchy(MemberType *first)
Construct a new hierarchy.
Definition: hierarchy.hh:1246
Hierarchy()
Construct a new empty hierarchy.
Definition: hierarchy.hh:1228
~Hierarchy()
Destructor.
Definition: hierarchy.hh:1258
AccumulationMode
Identifiers for the different accumulation modes.
Definition: parameters.hh:230
static T * construct(Arguments &args)
Construct an object with the specified arguments.
Definition: construction.hh:52
Iterator finest()
Get an iterator positioned at the finest level.
Definition: hierarchy.hh:1363
void build(const T &criterion)
Build the matrix hierarchy using aggregation.
Definition: hierarchy.hh:663
MatrixHierarchy(const MatrixOperator &fineMatrix, const ParallelInformation &pinfo=ParallelInformation())
Constructor.
Definition: hierarchy.hh:652
void free()
Free the allocated memory.
static void deconstruct(T *t)
Destroys an object.
Definition: construction.hh:61
void coarsenSmoother(Hierarchy< S, TA > &smoothers, const typename SmootherTraits< S >::Arguments &args) const
Coarsen the smoother hierarchy according to the matrix hierarchy.
Definition: hierarchy.hh:1137
void buildDependency(G &graph, const typename C::Matrix &matrix, C criterion, bool finestLevel)
Build the dependency of the matrix graph.
std::tuple< int, int, int, int > buildAggregates(const M &matrix, G &graph, const C &criterion, bool finestLevel)
Build the aggregates.
void getCoarsestAggregatesOnFinest(std::vector< std::size_t > &data) const
Get the mapping of fine level unknowns to coarse level aggregates.
Definition: hierarchy.hh:1009
@ MAX_PROCESSES
Hard limit for the number of processes allowed.
Definition: hierarchy.hh:55
@ atOnceAccu
Accumulate data to on process at once.
Definition: parameters.hh:242
@ successiveAccu
Successively accumulate to fewer processes.
Definition: parameters.hh:246
int countNonZeros(const M &matrix)
Get the number of nonzero fields in the matrix.
Definition: matrixutils.hh:153
DVVerbType dvverb(std::cout)
stream for very verbose output.
Definition: stdstreams.hh:93
DInfoType dinfo(std::cout)
Stream for informative output.
Definition: stdstreams.hh:138
DVerbType dverb(std::cout)
Singleton of verbose debug stream.
Definition: stdstreams.hh:114
Provides a map between global and local indices.
Provides a class for building the index set and remote indices on the coarse level.
Functionality for redistributing a sparse matrix.
Some handy generic functions for ISTL matrices.
Dune namespace.
Definition: alignedallocator.hh:10
bool graphRepartition(const G &graph, Dune::OwnerOverlapCopyCommunication< T1, T2 > &oocomm, Metis::idx_t nparts, Dune::OwnerOverlapCopyCommunication< T1, T2 > *&outcomm, RedistributeInterface &redistInf, bool verbose=false)
execute a graph repartition for a giving graph and indexset.
Definition: repartition.hh:1265
void redistributeMatrix(M &origMatrix, M &newMatrix, C &origComm, C &newComm, RedistributeInformation< C > &ri)
Redistribute a matrix according to given domain decompositions.
Definition: matrixredistribute.hh:834
STL namespace.
Classes for the generic construction and application of the smoothers.
Standard Dune debug streams.
Traits class for getting the attribute class of a smoother.
Definition: smoother.hh:64
Tag idnetifying the visited property of a vertex.
Definition: properties.hh:27
A property map that applies the identity function to integers.
Definition: propertymap.hh:291
Selector for the property map type.
Definition: propertymap.hh:318
static Category category(const OP &op, decltype(op.category()) *=nullptr)
Helperfunction to extract the solver category either from an enum, or from the newly introduced virtu...
Definition: solvercategory.hh:32
A simple timing class.
Prolongation and restriction for amg.
Definition of the DUNE_UNUSED macro for the case that config.h is not available.
Creative Commons License   |  Legal Statements / Impressum  |  Hosted by TU Dresden  |  generated with Hugo v0.111.3 (Nov 24, 23:30, 2024)