3 #ifndef DUNE_AMGHIERARCHY_HH 4 #define DUNE_AMGHIERARCHY_HH 16 #include <dune/common/stdstreams.hh> 17 #include <dune/common/unused.hh> 18 #include <dune/common/timer.hh> 19 #include <dune/common/bigunsignedint.hh> 21 #include <dune/common/parallel/indexset.hh> 65 template<
typename T,
typename A=std::allocator<T> >
74 template<
typename T1,
typename T2>
96 MemberType* redistributed_;
103 typedef typename A::template rebind<Element>::other
Allocator;
150 template<
class C,
class T1>
152 :
public BidirectionalIteratorFacade<LevelIterator<C,T1>,T1,T1&>
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 >;
171 typename std::remove_const<T1>::type>& other)
172 : element_(other.element_)
177 const typename std::remove_const<T1>::type>& other)
178 : element_(other.element_)
185 typename std::remove_const<T1>::type>& other)
const 187 return element_ == other.element_;
194 const typename std::remove_const<T1>::type>& other)
const 196 return element_ == other.element_;
202 return *(element_->element_);
208 element_ = element_->coarser_;
214 element_ = element_->finer_;
223 return element_->redistributed_;
232 assert(element_->redistributed_);
233 return *element_->redistributed_;
237 element_->redistributed_ = t;
242 element_->redistributed_ =
nullptr;
272 ConstIterator
finest()
const;
284 std::size_t
levels()
const;
295 Element* nonAllocated_;
297 Allocator allocator_;
308 template<
class M,
class PI,
class A=std::allocator<M> >
316 typedef typename MatrixOperator::matrix_type
Matrix;
334 typedef typename Allocator::template rebind<AggregatesMap*>::other
AAllocator;
343 typedef typename Allocator::template rebind<RedistributeInfoType>::other
RILAllocator;
354 const ParallelInformation& pinfo=ParallelInformation());
364 template<
typename O,
typename T>
365 void build(
const T& criterion);
375 void recalculateGalerkin(
const F& copyFlags);
381 template<
class V,
class BA,
class TA>
389 template<
class S,
class TA>
397 std::size_t
levels()
const;
403 std::size_t maxlevels()
const;
405 bool hasCoarsest()
const;
411 bool isBuilt()
const;
417 const ParallelMatrixHierarchy& matrices()
const;
423 const ParallelInformationHierarchy& parallelInformation()
const;
429 const AggregatesMapList& aggregatesMaps()
const;
436 const RedistributeInfoList& redistributeInformation()
const;
453 void getCoarsestAggregatesOnFinest(std::vector<std::size_t>& data)
const;
459 AggregatesMapList aggregatesMaps_;
461 RedistributeInfoList redistributes_;
463 ParallelMatrixHierarchy matrices_;
465 ParallelInformationHierarchy parallelInformation_;
478 template<
class Matrix,
bool pr
int>
485 static void stats(
const Matrix& matrix)
487 DUNE_UNUSED_PARAMETER(matrix);
491 template<
class Matrix>
492 struct MatrixStats<Matrix,true>
501 min=std::numeric_limits<size_type>::max();
508 min=std::min(min, row.size());
509 max=std::max(max, row.size());
520 static void stats(
const Matrix& matrix)
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()
555 : AggregationCriterion(
Dune::Amg::
Parameters(maxLevel, coarsenTarget, minCoarsenRate, prolongDamp, accumulate))
559 : AggregationCriterion(parms)
564 template<
typename M,
typename C1>
569 int nparts, C1& criterion)
571 DUNE_UNUSED_PARAMETER(origMatrix);
572 DUNE_UNUSED_PARAMETER(newMatrix);
573 DUNE_UNUSED_PARAMETER(origSequentialInformationomm);
574 DUNE_UNUSED_PARAMETER(newComm);
575 DUNE_UNUSED_PARAMETER(ri);
576 DUNE_UNUSED_PARAMETER(nparts);
577 DUNE_UNUSED_PARAMETER(criterion);
578 DUNE_THROW(NotImplemented,
"Redistribution does not make sense in sequential code!");
582 template<
typename M,
typename C,
typename C1>
585 int nparts, C1& criterion)
588 #ifdef AMG_REPART_ON_COMM_GRAPH 592 criterion.debugLevel()>1);
601 MatrixGraph graph(origMatrix);
602 PropertiesGraph pgraph(graph);
606 if(origComm.communicator().rank()==0)
607 std::cout<<
"Original matrix"<<std::endl;
608 origComm.communicator().barrier();
612 newComm, ri.getInterface(),
613 criterion.debugLevel()>1);
614 #endif // if else AMG_REPART 616 if(origComm.communicator().rank()==0 && criterion.debugLevel()>1)
617 std::cout<<
"Repartitioning took "<<time.elapsed()<<
" seconds."<<std::endl;
622 ri.checkInterface(origComm.indexSet(), newComm->indexSet(), origComm.communicator());
628 if(origComm.communicator().rank()==0)
629 std::cout<<
"Original matrix"<<std::endl;
630 origComm.communicator().barrier();
631 if(newComm->communicator().size()>0)
633 origComm.communicator().barrier();
636 if(origComm.communicator().rank()==0 && criterion.debugLevel()>1)
637 std::cout<<
"Redistributing matrix took "<<time.elapsed()<<
" seconds."<<std::endl;
638 return existentOnRedist;
651 template<
class M,
class IS,
class A>
653 const ParallelInformation& pinfo)
654 : matrices_(const_cast<MatrixOperator&>(fineOperator)),
655 parallelInformation_(const_cast<ParallelInformation&>(pinfo))
658 DUNE_THROW(
ISTLError,
"MatrixOperator and ParallelInformation must belong to the same category!");
661 template<
class M,
class IS,
class A>
662 template<
typename O,
typename T>
665 prolongDamp_ = criterion.getProlongationDampingFactor();
666 typedef O OverlapFlags;
672 typedef bigunsignedint<sizeof(int)*8*noints> BIGINT;
674 MatIterator mlevel = matrices_.
finest();
675 MatrixStats<typename M::matrix_type,MINIMAL_DEBUG_LEVEL<=INFO_DEBUG_LEVEL>::stats(mlevel->getmat());
677 PInfoIterator infoLevel = parallelInformation_.
finest();
679 finenonzeros = infoLevel->communicator().sum(finenonzeros);
680 BIGINT allnonzeros = finenonzeros;
686 BIGINT unknowns = mlevel->getmat().N();
688 unknowns = infoLevel->communicator().sum(unknowns);
689 double dunknowns=unknowns.todouble();
690 infoLevel->buildGlobalLookup(mlevel->getmat().N());
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;
700 MatrixOperator* matrix=&(*mlevel);
701 ParallelInformation* info =&(*infoLevel);
710 && dunknowns < 30*infoLevel->communicator().size()))
711 && infoLevel->communicator().size()>1 &&
712 dunknowns/infoLevel->communicator().size() <= criterion.coarsenTarget())
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() )
723 bool existentOnNextLevel =
725 redistComm, redistributes_.back(), nodomains,
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);
735 assert(mlevel.isRedistributed());
736 infoLevel.addRedistributed(redistComm);
737 infoLevel->freeGlobalLookup();
739 if(!existentOnNextLevel)
744 matrix = &(mlevel.getRedistributed());
745 info = &(infoLevel.getRedistributed());
746 info->buildGlobalLookup(matrix->getmat().N());
749 rank = info->communicator().rank();
750 if(dunknowns <= criterion.coarsenTarget())
756 typedef typename GraphCreator::GraphTuple GraphTuple;
760 std::vector<bool> excluded(matrix->getmat().N(),
false);
762 GraphTuple graphs = GraphCreator::create(*matrix, excluded, *info, OverlapFlags());
764 AggregatesMap* aggregatesMap=
new AggregatesMap(std::get<1>(graphs)->maxVertex()+1);
766 aggregatesMaps_.push_back(aggregatesMap);
770 int noAggregates, isoAggregates, oneAggregates, skippedAggregates;
772 std::tie(noAggregates, isoAggregates, oneAggregates, skippedAggregates) =
773 aggregatesMap->
buildAggregates(matrix->getmat(), *(std::get<1>(graphs)), criterion, level==0);
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;
782 int start, end, overlapStart, overlapEnd;
783 int procs=info->communicator().rank();
784 int n = UNKNOWNS/procs;
785 int bigger = UNKNOWNS%procs;
790 end = (rank+1)*(n+1);
792 start = bigger + rank * n;
793 end = bigger + (rank + 1) * n;
798 overlapStart = start - 1;
800 overlapStart = start;
803 overlapEnd = end + 1;
807 assert((UNKNOWNS)*(overlapEnd-overlapStart)==aggregatesMap->
noVertices());
808 for(
int j=0; j< UNKNOWNS; ++j)
809 for(
int i=0; i < UNKNOWNS; ++i)
811 if(i>=overlapStart && i<overlapEnd)
813 int no = (j/2)*((UNKNOWNS)/2)+i/2;
814 (*aggregatesMap)[j*(overlapEnd-overlapStart)+i-overlapStart]=no;
819 if(criterion.debugLevel()>1 && info->communicator().rank()==0)
820 std::cout<<
"aggregating finished."<<std::endl;
822 BIGINT gnoAggregates=noAggregates;
823 gnoAggregates = info->communicator().sum(gnoAggregates);
824 double dgnoAggregates = gnoAggregates.todouble();
826 BIGINT gnoAggregates=((UNKNOWNS)/2)*((UNKNOWNS)/2);
829 if(criterion.debugLevel()>2 && rank==0)
830 std::cout <<
"Building "<<dgnoAggregates<<
" aggregates took "<<watch.elapsed()<<
" seconds."<<std::endl;
832 if(dgnoAggregates==0 || dunknowns/dgnoAggregates<criterion.minCoarsenRate())
837 std::cerr <<
"Stopped coarsening because of rate breakdown "<<dunknowns<<
"/"<<dgnoAggregates
838 <<
"="<<dunknowns/dgnoAggregates<<
"<" 839 <<criterion.minCoarsenRate()<<std::endl;
841 std::cerr<<
"Could not build any aggregates. Probably no connected nodes."<<std::endl;
843 aggregatesMap->
free();
844 delete aggregatesMap;
845 aggregatesMaps_.pop_back();
847 if(criterion.accumulate() && mlevel.isRedistributed() && info->communicator().size()>1) {
851 delete &(mlevel.getRedistributed().getmat());
852 mlevel.deleteRedistributed();
853 delete &(infoLevel.getRedistributed());
854 infoLevel.deleteRedistributed();
855 redistributes_.back().resetSetup();
860 unknowns = noAggregates;
861 dunknowns = dgnoAggregates;
863 CommunicationArgs commargs(info->communicator(),info->getSolverCategory());
868 typename PropertyMapTypeSelector<VertexVisitedTag,PropertiesGraph>::Type visitedMap =
874 *(std::get<1>(graphs)),
879 GraphCreator::free(graphs);
881 if(criterion.debugLevel()>2) {
883 std::cout<<
"Coarsening of index sets took "<<watch.elapsed()<<
" seconds."<<std::endl;
888 infoLevel->buildGlobalLookup(aggregates);
891 infoLevel->globalLookup());
894 if(criterion.debugLevel()>2) {
896 std::cout<<
"Communicating global aggregate numbers took "<<watch.elapsed()<<
" seconds."<<std::endl;
900 std::vector<bool>& visited=excluded;
902 typedef std::vector<bool>::iterator Iterator;
903 typedef IteratorPropertyMap<Iterator, IdentityMap> VisitedMap2;
904 Iterator end = visited.end();
905 for(Iterator iter= visited.begin(); iter != end; ++iter)
908 VisitedMap2 visitedMap2(visited.begin(), Dune::IdentityMap());
910 typename MatrixOperator::matrix_type* coarseMatrix;
912 coarseMatrix = productBuilder.
build(*(std::get<0>(graphs)), visitedMap2,
917 dverb<<
"Building of sparsity pattern took "<<watch.elapsed()<<std::endl;
919 info->freeGlobalLookup();
921 delete std::get<0>(graphs);
922 productBuilder.
calculate(matrix->getmat(), *aggregatesMap, *coarseMatrix, *infoLevel, OverlapFlags());
924 if(criterion.debugLevel()>2) {
926 std::cout<<
"Calculation entries of Galerkin product took "<<watch.elapsed()<<
" seconds."<<std::endl;
930 allnonzeros = allnonzeros + infoLevel->communicator().sum(nonzeros);
931 MatrixArgs args(*coarseMatrix, *infoLevel);
938 infoLevel->freeGlobalLookup();
942 aggregatesMaps_.push_back(aggregatesMap);
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;
956 if(criterion.accumulate() && !redistributes_.back().isSetup() &&
957 infoLevel->communicator().size()>1) {
958 #if HAVE_MPI && !HAVE_PARMETIS 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;
966 Matrix* redistMat=
new Matrix();
967 ParallelInformation* redistComm=0;
971 redistComm, redistributes_.back(), nodomains,criterion);
972 MatrixArgs args(*redistMat, *redistComm);
973 BIGINT unknownsRedist = redistMat->N();
974 unknownsRedist = infoLevel->communicator().sum(unknownsRedist);
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;
982 infoLevel.addRedistributed(redistComm);
983 infoLevel->freeGlobalLookup();
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;
994 template<
class M,
class IS,
class A>
1001 template<
class M,
class IS,
class A>
1005 return parallelInformation_;
1008 template<
class M,
class IS,
class A>
1012 int maxlevels=parallelInformation_.
finest()->communicator().max(levels);
1015 std::vector<std::size_t> tmp;
1016 std::vector<std::size_t> *coarse, *fine;
1033 if(levels==maxlevels) {
1039 m=std::max(*iter,m);
1041 coarse->resize(m+1);
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();
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());
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);
1061 for(
typename AggregatesMapList::const_reverse_iterator aggregates=++
aggregatesMaps().rbegin();
1064 fine->resize((*aggregates)->noVertices());
1065 fine->assign(fine->size(), 0);
1067 ::prolongateVector(*(*aggregates), *coarse, *fine, static_cast<std::size_t>(1), *pinfo);
1069 std::swap(coarse, fine);
1073 assert(coarse==&data);
1076 template<
class M,
class IS,
class A>
1080 return aggregatesMaps_;
1082 template<
class M,
class IS,
class A>
1086 return redistributes_;
1089 template<
class M,
class IS,
class A>
1092 typedef typename AggregatesMapList::reverse_iterator AggregatesMapIterator;
1096 AggregatesMapIterator amap = aggregatesMaps_.rbegin();
1097 InfoIterator info = parallelInformation_.
coarsest();
1098 for(Iterator level=matrices_.
coarsest(), finest=matrices_.
finest(); level != finest; --level, --info, ++amap) {
1101 delete &level->getmat();
1102 if(level.isRedistributed())
1103 delete &(level.getRedistributed().getmat());
1108 template<
class M,
class IS,
class A>
1109 template<
class V,
class BA,
class TA>
1112 assert(hierarchy.levels()==1);
1114 typedef typename RedistributeInfoList::const_iterator RIter;
1115 RIter redist = redistributes_.begin();
1117 Iterator matrix = matrices_.
finest(), coarsest = matrices_.
coarsest();
1119 if(redist->isSetup())
1120 hierarchy.addRedistributedOnCoarsest(matrix.getRedistributed().getmat().N());
1121 Dune::dvverb<<
"Level "<<level<<
" has "<<matrices_.
finest()->getmat().N()<<
" unknowns!"<<std::endl;
1123 while(matrix != coarsest) {
1124 ++matrix; ++level; ++redist;
1125 Dune::dvverb<<
"Level "<<level<<
" has "<<matrix->getmat().N()<<
" unknowns!"<<std::endl;
1127 hierarchy.addCoarser(matrix->getmat().N());
1128 if(redist->isSetup())
1129 hierarchy.addRedistributedOnCoarsest(matrix.getRedistributed().getmat().N());
1135 template<
class M,
class IS,
class A>
1136 template<
class S,
class TA>
1140 assert(smoothers.
levels()==0);
1143 typedef typename AggregatesMapList::const_iterator AggregatesIterator;
1146 cargs.setArgs(sargs);
1147 PinfoIterator pinfo = parallelInformation_.
finest();
1148 AggregatesIterator aggregates = aggregatesMaps_.begin();
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);
1158 cargs.setMatrix(matrices_.
coarsest()->getmat(), **aggregates);
1159 cargs.setComm(*pinfo);
1165 template<
class M,
class IS,
class A>
1169 typedef typename AggregatesMapList::iterator AggregatesMapIterator;
1173 AggregatesMapIterator amap = aggregatesMaps_.begin();
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());
1181 const_cast<Matrix&>(level.getRedistributed().getmat()),
1182 *info,info.getRedistributed(), *riIter);
1183 info->freeGlobalLookup();
1186 for(; level!=coarsest; ++amap) {
1187 const Matrix& fine = (level.isRedistributed() ? level.getRedistributed() : *level).getmat();
1191 productBuilder.
calculate(fine, *(*amap), const_cast<Matrix&>(level->getmat()), *info, copyFlags);
1192 if(level.isRedistributed()) {
1193 info->buildGlobalLookup(level->getmat().N());
1195 const_cast<Matrix&>(level.getRedistributed().getmat()), *info,
1196 info.getRedistributed(), *riIter);
1197 info->freeGlobalLookup();
1202 template<
class M,
class IS,
class A>
1205 return matrices_.
levels();
1208 template<
class M,
class IS,
class A>
1214 template<
class M,
class IS,
class A>
1218 (!matrices_.
coarsest().isRedistributed() ||matrices_.
coarsest()->getmat().N()>0);
1221 template<
class M,
class IS,
class A>
1227 template<
class T,
class A>
1229 : finest_(0), coarsest_(0), nonAllocated_(0), allocator_(), levels_(0)
1232 template<
class T,
class A>
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;
1245 template<
class T,
class A>
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;
1257 template<
class T,
class A>
1261 Element* current = coarsest_;
1262 coarsest_ = coarsest_->finer_;
1263 if(current != nonAllocated_) {
1264 if(current->redistributed_)
1268 allocator_.deallocate(current, 1);
1274 template<
class T,
class A>
1276 : nonAllocated_(), allocator_(other.allocator_),
1277 levels_(other.levels_)
1281 finest_=coarsest_=nonAllocated_=
nullptr;
1284 finest_=allocator_.allocate(1,0);
1285 Element* finer_ =
nullptr;
1286 Element* current_ = finest_;
1287 Element* otherCurrent_ = other.finest_;
1289 while(otherCurrent_)
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_);
1297 current_->redistributed_=
nullptr;
1299 if(otherCurrent_->coarser_)
1301 current_->coarser_=allocator_.allocate(1,0);
1302 current_=current_->coarser_;
1304 current_->coarser_=
nullptr;
1305 otherCurrent_=otherCurrent_->coarser_;
1310 template<
class T,
class A>
1316 template<
class T,
class A>
1322 template<
class T,
class A>
1327 coarsest_ = allocator_.allocate(1,0);
1329 finest_ = coarsest_;
1330 coarsest_->finer_ =
nullptr;
1332 coarsest_->coarser_ = allocator_.allocate(1,0);
1333 coarsest_->coarser_->finer_ = coarsest_;
1334 coarsest_ = coarsest_->coarser_;
1337 coarsest_->redistributed_ =
nullptr;
1338 coarsest_->coarser_=
nullptr;
1343 template<
class T,
class A>
1348 finest_ = allocator_.allocate(1,0);
1350 coarsest_ = finest_;
1351 coarsest_->coarser_ = coarsest_->finer_ =
nullptr;
1353 finest_->finer_ = allocator_.allocate(1,0);
1354 finest_->finer_->coarser_ = finest_;
1355 finest_ = finest_->finer_;
1356 finest_->finer =
nullptr;
1362 template<
class T,
class A>
1368 template<
class T,
class A>
1374 template<
class T,
class A>
1380 template<
class T,
class A>
T MemberType
The type of the container we store.
Definition: hierarchy.hh:72
Functionality for redistributing a sparse matrix.
void printGlobalSparseMatrix(const M &mat, C &ooc, std::ostream &os)
Definition: matrixutils.hh:175
const void * Arguments
A type holding all the arguments needed to call the constructor.
Definition: construction.hh:44
const_iterator begin() const
Definition: aggregates.hh:712
Definition: allocator.hh:7
Dune::Amg::Hierarchy< MatrixOperator, Allocator > ParallelMatrixHierarchy
The type of the parallel matrix hierarchy.
Definition: hierarchy.hh:328
std::size_t noVertices() const
Get the number of vertices.
size_type max
Definition: hierarchy.hh:514
bool hasCoarsest() const
Definition: hierarchy.hh:1215
T1 & dereference() const
Dereference the iterator.
Definition: hierarchy.hh:200
const AggregateDescriptor * const_iterator
Definition: aggregates.hh:710
const RedistributeInfoList & redistributeInformation() const
Get the hierarchy of the information about redistributions,.
Definition: hierarchy.hh:1084
LevelIterator(const LevelIterator< typename std::remove_const< C >::type, typename std::remove_const< T1 >::type > &other)
Copy constructor.
Definition: hierarchy.hh:170
Helper classes for the construction of classes without empty constructor.
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
Traits class for generically constructing non default constructable types.
Definition: novlpschwarz.hh:247
PI ParallelInformation
The type of the index set.
Definition: hierarchy.hh:319
MatrixHierarchy(const MatrixOperator &fineMatrix, const ParallelInformation &pinfo=ParallelInformation())
Constructor.
Definition: hierarchy.hh:652
size_type sum
Definition: hierarchy.hh:515
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
LevelIterator< Hierarchy< T, A >, T > Iterator
Type of the mutable iterator.
Definition: hierarchy.hh:250
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 calculate(const M &fine, const AggregatesMap< V > &aggregates, M &coarse, const I &pinfo, const O ©)
Calculate the galerkin product.
All parameters for AMG.
Definition: parameters.hh:390
void build(const T &criterion)
Build the matrix hierarchy using aggregation.
Definition: hierarchy.hh:663
This file implements a vector space as a tensor product of a given vector space. The number of compon...
Matrix::row_type matrix_row
Definition: hierarchy.hh:497
std::list< AggregatesMap *, AAllocator > AggregatesMapList
The type of the aggregates maps list.
Definition: hierarchy.hh:337
void operator()(const matrix_row &row)
Definition: hierarchy.hh:506
RedistributeInformation< ParallelInformation > RedistributeInfoType
The type of the redistribute information.
Definition: hierarchy.hh:340
MatrixOperator::matrix_type Matrix
The type of the matrix.
Definition: hierarchy.hh:316
bool repartitionAndDistributeMatrix(const M &origMatrix, M &newMatrix, SequentialInformation &origSequentialInformationomm, SequentialInformation *&newComm, RedistributeInformation< SequentialInformation > &ri, int nparts, C1 &criterion)
Definition: hierarchy.hh:565
T1 & getRedistributed() const
Get the redistributed container.
Definition: hierarchy.hh:230
T AggregationCriterion
The criterion for tagging connections as strong and nodes as isolated. This might be e...
Definition: hierarchy.hh:541
void deleteRedistributed()
Definition: hierarchy.hh:240
std::list< RedistributeInfoType, RILAllocator > RedistributeInfoList
The type of the list of redistribute information.
Definition: hierarchy.hh:346
void free()
Free the allocated memory.
Hierarchy()
Construct a new empty hierarchy.
Definition: hierarchy.hh:1228
void increment()
Move to the next coarser level.
Definition: hierarchy.hh:206
Accumulate data to on process at once.
Definition: parameters.hh:242
calc()
Definition: hierarchy.hh:499
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
Iterator over the levels in the hierarchy.
Definition: hierarchy.hh:75
Allocator::template rebind< AggregatesMap * >::other AAllocator
Allocator for pointers.
Definition: hierarchy.hh:334
int countNonZeros(const M &matrix)
Get the number of nonzero fields in the matrix.
Definition: matrixutils.hh:153
Definition: transfer.hh:30
bool equals(const LevelIterator< typename std::remove_const< C >::type, typename std::remove_const< T1 >::type > &other) const
Equality check.
Definition: hierarchy.hh:184
void decrement()
Move to the next fine level.
Definition: hierarchy.hh:212
A::size_type size_type
Type for indices and sizes.
Definition: matrix.hh:571
AccumulationMode
Identifiers for the different accumulation modes.
Definition: parameters.hh:230
LevelIterator(const LevelIterator< const typename std::remove_const< C >::type, const typename std::remove_const< T1 >::type > &other)
Copy constructor.
Definition: hierarchy.hh:176
Definition: galerkin.hh:115
Attaches properties to the edges and vertices of a graph.
Definition: graph.hh:975
Tag idnetifying the visited property of a vertex.
Definition: properties.hh:26
Iterator coarsest()
Get an iterator positioned at the coarsest level.
Definition: hierarchy.hh:1369
~MatrixHierarchy()
Definition: hierarchy.hh:1090
Matrix::size_type size_type
Definition: hierarchy.hh:496
A::template rebind< Element >::other Allocator
The allocator to use for the list elements.
Definition: hierarchy.hh:103
double getProlongationDampingFactor() const
Definition: hierarchy.hh:438
MatrixImp::DenseMatrixBase< T, A >::window_type row_type
The type implementing a matrix row.
Definition: matrix.hh:568
void buildDependency(G &graph, const typename C::Matrix &matrix, C criterion, bool finestLevel)
Build the dependency of the matrix graph.
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:1266
const ParallelInformationHierarchy & parallelInformation() const
Get the hierarchy of the parallel data distribution information.
Definition: hierarchy.hh:1003
Classes for the generic construction and application of the smoothers.
std::size_t levels() const
Get the number of levels in the hierarchy.
Definition: hierarchy.hh:1203
const_iterator end() const
Definition: aggregates.hh:717
Class providing information about the mapping of the vertices onto aggregates.
Definition: aggregates.hh:542
derive error class from the base class in common
Definition: istlexception.hh:16
Definition: graphcreator.hh:20
void getCoarsestAggregatesOnFinest(std::vector< std::size_t > &data) const
Get the mapping of fine level unknowns to coarse level aggregates.
Definition: hierarchy.hh:1009
Hard limit for the number of processes allowed.
Definition: hierarchy.hh:55
ConstructionTraits< T >::Arguments Arguments
Definition: hierarchy.hh:105
Definition: galerkin.hh:97
void addRedistributed(T1 *t)
Definition: hierarchy.hh:235
Prolongation and restriction for amg.
bool isRedistributed() const
Check whether there was a redistribution at the current level.
Definition: hierarchy.hh:221
Dune::Amg::Hierarchy< ParallelInformation, Allocator > ParallelInformationHierarchy
The type of the parallel informarion hierarchy.
Definition: hierarchy.hh:331
static const V ISOLATED
Identifier of isolated vertices.
Definition: aggregates.hh:554
LevelIterator< const Hierarchy< T, A >, const T > ConstIterator
Type of the const iterator.
Definition: hierarchy.hh:253
The hierarchies build by the coarsening process.
Definition: hierarchy.hh:309
Allocator::template rebind< RedistributeInfoType >::other RILAllocator
Allocator for RedistributeInfoType.
Definition: hierarchy.hh:343
~Hierarchy()
Destructor.
Definition: hierarchy.hh:1258
Provides a class for building the galerkin product based on a aggregation scheme. ...
A hierarchy of coantainers (e.g. matrices or vectors)
Definition: hierarchy.hh:66
void addRedistributedOnCoarsest(Arguments &args)
Definition: hierarchy.hh:1317
Traits class for getting the attribute class of a smoother.
Definition: smoother.hh:63
CoarsenCriterion(const Dune::Amg::Parameters &parms)
Definition: hierarchy.hh:558
bool commGraphRepartition(const M &mat, Dune::OwnerOverlapCopyCommunication< T1, T2 > &oocomm, Metis::idx_t nparts, Dune::OwnerOverlapCopyCommunication< T1, T2 > *&outcomm, RedistributeInterface &redistInf, bool verbose=false)
Definition: repartition.hh:856
const AggregatesMapList & aggregatesMaps() const
Get the hierarchy of the mappings of the nodes onto aggregates.
Definition: hierarchy.hh:1078
LevelIterator(Element *element)
Definition: hierarchy.hh:165
Class representing a node in the matrix graph.
Definition: dependency.hh:125
A vector of blocks with memory management.
Definition: bvector.hh:316
Definition: indicescoarsener.hh:34
Definition: globalaggregates.hh:129
G::MutableMatrix * build(G &fineGraph, V &visitedMap, const ParallelInformation &pinfo, AggregatesMap< typename G::VertexDescriptor > &aggregates, const typename G::Matrix::size_type &size, const Set ©)
Calculates the coarse matrix via a Galerkin product.
Definition: galerkin.hh:564
void coarsenVector(Hierarchy< BlockVector< V, BA >, TA > &hierarchy) const
Coarsen the vector hierarchy according to the matrix hierarchy.
Definition: hierarchy.hh:1110
M MatrixOperator
The type of the matrix operator.
Definition: hierarchy.hh:313
void recalculateGalerkin(const F ©Flags)
Recalculate the galerkin products.
Definition: hierarchy.hh:1167
std::size_t maxlevels() const
Get the max number of levels in the hierarchy of processors.
Definition: hierarchy.hh:1209
Graph::VertexDescriptor VertexDescriptor
The vertex descriptor.
Definition: graph.hh:986
Provides classes for the Coloring process of AMG.
std::tuple< int, int, int, int > buildAggregates(const M &matrix, G &graph, const C &criterion, bool finestLevel)
Build the aggregates.
std::size_t levels() const
Get the number of levels in the hierarchy.
Definition: hierarchy.hh:1311
Provides a class for building the index set and remote indices on the coarse level.
Definition: matrixredistribute.hh:20
Class representing the properties of an ede in the matrix graph.
Definition: dependency.hh:37
The criterion describing the stop criteria for the coarsening process.
Definition: hierarchy.hh:534
void addFiner(Arguments &args)
Add an element on a finer level.
Definition: hierarchy.hh:1344
Provdes class for identifying aggregates globally.
Successively accumulate to fewer processes.
Definition: parameters.hh:246
Dune::Amg::AggregatesMap< typename MatrixGraph< Matrix >::VertexDescriptor > AggregatesMap
The type of the aggregates map we use.
Definition: hierarchy.hh:325
Provides classes for initializing the link attributes of a matrix graph.
size_type min
Definition: hierarchy.hh:513
void addCoarser(Arguments &args)
Add an element on a coarser level.
Definition: hierarchy.hh:1323
const ParallelMatrixHierarchy & matrices() const
Get the matrix hierarchy.
Definition: hierarchy.hh:996
CoarsenCriterion(int maxLevel=100, int coarsenTarget=1000, double minCoarsenRate=1.2, double prolongDamp=1.6, AccumulationMode accumulate=successiveAccu)
Constructor.
Definition: hierarchy.hh:553
void redistributeMatrixEntries(M &origMatrix, M &newMatrix, C &origComm, C &newComm, RedistributeInformation< C > &ri)
Definition: matrixredistribute.hh:771
bool isBuilt() const
Whether the hierarchy was built.
Definition: hierarchy.hh:1222
Iterator finest()
Get an iterator positioned at the finest level.
Definition: hierarchy.hh:1363
Provides classes for building the matrix graph.
Some handy generic functions for ISTL matrices.