3#ifndef DUNE_AMG_AMG_CPR_HH
4#define DUNE_AMG_AMG_CPR_HH
10#include <opm/common/TimingMacros.hpp>
11#include <dune/common/exceptions.hh>
12#include <dune/common/version.hh>
13#include <dune/istl/paamg/amg.hh>
14#include <dune/istl/paamg/smoother.hh>
15#include <dune/istl/paamg/transfer.hh>
16#include <dune/istl/paamg/hierarchy.hh>
17#include <dune/istl/solvers.hh>
18#include <dune/istl/scalarproducts.hh>
19#include <dune/istl/superlu.hh>
20#include <dune/istl/umfpack.hh>
21#include <dune/istl/solvertype.hh>
22#include <dune/common/typetraits.hh>
23#include <dune/common/exceptions.hh>
34 template<
class M,
class T>
40 template<
class M,
class PI>
41 typename std::enable_if<!std::is_same<PI,SequentialInformation>::value,
void>::type
43 Dune::RedistributeInformation<PI>& redistInfo)
46 info.buildGlobalLookup(mat.N());
47 redistributeMatrixEntries(mat, matRedist, info, infoRedist, redistInfo);
48 info.freeGlobalLookup();
69 template<
class M,
class X,
class S,
class P,
class K,
class A>
85 template<
class M,
class X,
class S,
class PI=SequentialInformation,
86 class A=std::allocator<X> >
89 template<
class M1,
class X1,
class S1,
class P1,
class K1,
class A1>
135 const SmootherArgs& smootherArgs,
const Parameters& parms);
196 OPM_TIMEBLOCK(recalculateHierarch);
197 auto copyFlags = NegateSet<typename PI::OwnerSet>();
198 const auto& matrices = matrices_->matrices();
199 const auto& aggregatesMapHierarchy = matrices_->aggregatesMaps();
200 const auto& infoHierarchy = matrices_->parallelInformation();
201 const auto& redistInfoHierarchy = matrices_->redistributeInformation();
202 BaseGalerkinProduct productBuilder;
203 auto aggregatesMap = aggregatesMapHierarchy.begin();
204 auto info = infoHierarchy.finest();
205 auto redistInfo = redistInfoHierarchy.begin();
206 auto matrix = matrices.finest();
207 auto coarsestMatrix = matrices.coarsest();
209 using Matrix =
typename M::matrix_type;
212 if(matrix.isRedistributed()) {
214 const_cast<Matrix&
>(matrix.getRedistributed().getmat()),
215 const_cast<PI&
>(*info),
const_cast<PI&
>(info.getRedistributed()),
216 const_cast<Dune::RedistributeInformation<PI>&
>(*redistInfo));
220 for(; matrix!=coarsestMatrix; ++aggregatesMap) {
221 const Matrix& fine = (matrix.isRedistributed() ? matrix.getRedistributed() : *matrix).getmat();
225 productBuilder.calculate(fine, *(*aggregatesMap),
const_cast<Matrix&
>(matrix->getmat()), *info, copyFlags);
227 if(matrix.isRedistributed()) {
229 const_cast<Matrix&
>(matrix.getRedistributed().getmat()),
230 const_cast<PI&
>(*info),
const_cast<PI&
>(info.getRedistributed()),
231 const_cast<Dune::RedistributeInformation<PI>&
>(*redistInfo));
262 void createHierarchies(C& criterion,
Operator& matrix,
267 std::shared_ptr< Operator > op( &matrix, [](
Operator* ) {});
268 std::shared_ptr< PI > pifo(
const_cast< PI*
> (&pinfo), []( PI * ) {});
269 createHierarchies( criterion, op, pifo);
273 void createHierarchies(C& criterion, std::shared_ptr< Operator > matrix,
274 std::shared_ptr< PI > pinfo );
276 void setupCoarseSolver();
294 typename OperatorHierarchy::ParallelMatrixHierarchy::ConstIterator
matrix;
298 typename ParallelInformationHierarchy::Iterator
pinfo;
302 typename OperatorHierarchy::RedistributeInfoList::const_iterator
redist;
306 typename OperatorHierarchy::AggregatesMapList::const_iterator
aggregates;
310 typename Hierarchy<Domain,A>::Iterator
lhs;
314 typename Hierarchy<Domain,A>::Iterator
update;
318 typename Hierarchy<Range,A>::Iterator
rhs;
330 void mgc(LevelContext& levelContext);
340 void moveToFineLevel(LevelContext& levelContext,
bool processedFineLevel);
346 bool moveToCoarseLevel(LevelContext& levelContext);
352 void initIteratorsWithFineLevel(LevelContext& levelContext);
355 std::shared_ptr<OperatorHierarchy> matrices_;
359 std::shared_ptr<Hierarchy<Smoother,A> > smoothers_;
361 std::shared_ptr<CoarseSolver> solver_;
363 std::shared_ptr< Hierarchy<Range,A> > rhs_;
365 std::shared_ptr< Hierarchy<Domain,A> > lhs_;
367 std::shared_ptr< Hierarchy<Domain,A> > update_;
369 using ScalarProduct = Dune::ScalarProduct<X>;
371 std::shared_ptr<ScalarProduct> scalarProduct_;
375 std::size_t preSteps_;
377 std::size_t postSteps_;
378 bool buildHierarchy_;
380 bool coarsesolverconverged;
381 std::shared_ptr<Smoother> coarseSmoother_;
383 SolverCategory::Category category_;
385 std::size_t verbosity_;
388 template<
class M,
class X,
class S,
class PI,
class A>
390 : matrices_(amg.matrices_), smootherArgs_(amg.smootherArgs_),
391 smoothers_(amg.smoothers_), solver_(amg.solver_),
392 rhs_(), lhs_(), update_(),
393 scalarProduct_(amg.scalarProduct_), gamma_(amg.gamma_),
394 preSteps_(amg.preSteps_), postSteps_(amg.postSteps_),
395 buildHierarchy_(amg.buildHierarchy_),
396 additive(amg.additive), coarsesolverconverged(amg.coarsesolverconverged),
397 coarseSmoother_(amg.coarseSmoother_),
398 category_(amg.category_),
399 verbosity_(amg.verbosity_)
402 rhs_.reset(
new Hierarchy<Range,A>(*amg.rhs_) );
404 lhs_.reset(
new Hierarchy<Domain,A>(*amg.lhs_) );
406 update_.reset(
new Hierarchy<Domain,A>(*amg.update_) );
409 template<
class M,
class X,
class S,
class PI,
class A>
412 const Parameters& parms)
413 : matrices_(stackobject_to_shared_ptr(matrices)), smootherArgs_(smootherArgs),
414 smoothers_(new Hierarchy<
Smoother,A>), solver_(&coarseSolver),
415 rhs_(), lhs_(), update_(), scalarProduct_(0),
416 gamma_(parms.getGamma()), preSteps_(parms.getNoPreSmoothSteps()),
417 postSteps_(parms.getNoPostSmoothSteps()), buildHierarchy_(false),
418 additive(parms.getAdditive()), coarsesolverconverged(true),
421 category_(SolverCategory::category(*smoothers_->coarsest())),
422 verbosity_(parms.debugLevel())
424 assert(matrices_->isBuilt());
427 matrices_->coarsenSmoother(*smoothers_, smootherArgs_);
430 template<
class M,
class X,
class S,
class PI,
class A>
436 : smootherArgs_(smootherArgs),
437 smoothers_(new Hierarchy<
Smoother,A>), solver_(),
438 rhs_(), lhs_(), update_(), scalarProduct_(),
439 gamma_(criterion.getGamma()), preSteps_(criterion.getNoPreSmoothSteps()),
440 postSteps_(criterion.getNoPostSmoothSteps()), buildHierarchy_(true),
441 additive(criterion.getAdditive()), coarsesolverconverged(true),
443 category_(SolverCategory::category(pinfo)),
444 verbosity_(criterion.debugLevel())
446 if(SolverCategory::category(matrix) != SolverCategory::category(pinfo))
447 DUNE_THROW(InvalidSolverCategory,
"Matrix and Communication must have the same SolverCategory!");
448 createHierarchies(criterion,
const_cast<Operator&
>(matrix), pinfo);
451 template<
class M,
class X,
class S,
class PI,
class A>
454 if(buildHierarchy_) {
458 coarseSmoother_.reset();
462 template<
class M,
class X,
class S,
class PI,
class A>
469 template<
class M,
class X,
class S,
class PI,
class A>
472 OPM_TIMEBLOCK(update);
474 smoothers_.reset(
new Hierarchy<Smoother,A>);
476 coarseSmoother_.reset();
477 scalarProduct_.reset();
478 buildHierarchy_=
true;
479 coarsesolverconverged =
true;
480 smoothers_.reset(
new Hierarchy<Smoother,A>);
481 recalculateHierarchy();
482 matrices_->coarsenSmoother(*smoothers_, smootherArgs_);
484 if (verbosity_>0 && matrices_->parallelInformation().finest()->communicator().rank()==0) {
485 std::cout <<
"Recalculating galerkin and coarse smoothers "<< matrices_->maxlevels() <<
" levels "
486 << watch.elapsed() <<
" seconds." << std::endl;
490 template<
class M,
class X,
class S,
class PI,
class A>
493 std::shared_ptr< PI > pinfo )
495 OPM_TIMEBLOCK(createHierarchies);
497 matrices_.reset(
new OperatorHierarchy(matrix, pinfo));
499 matrices_->template build<NegateSet<typename PI::OwnerSet> >(criterion);
502 matrices_->coarsenSmoother(*smoothers_, smootherArgs_);
504 if(verbosity_>0 && matrices_->parallelInformation().finest()->communicator().rank()==0)
505 std::cout<<
"Building hierarchy of "<<matrices_->maxlevels()<<
" levels "
506 <<
"(inclusive coarse solver) took "<<watch.elapsed()<<
" seconds."<<std::endl;
509 template<
class M,
class X,
class S,
class PI,
class A>
510 void AMGCPR<M,X,S,PI,A>::setupCoarseSolver()
512 OPM_TIMEBLOCK(setupCoarseSolver);
517 if(buildHierarchy_ && matrices_->levels()==matrices_->maxlevels()
518 && ( ! matrices_->redistributeInformation().back().isSetup() ||
519 matrices_->parallelInformation().coarsest().getRedistributed().communicator().size() ) )
522 SmootherArgs sargs(smootherArgs_);
523 sargs.iterations = 1;
525 typename ConstructionTraits<Smoother>::Arguments cargs;
526 cargs.setArgs(sargs);
527 if(matrices_->redistributeInformation().back().isSetup()) {
529 cargs.setMatrix(matrices_->matrices().coarsest().getRedistributed().getmat());
530 cargs.setComm(matrices_->parallelInformation().coarsest().getRedistributed());
532 cargs.setMatrix(matrices_->matrices().coarsest()->getmat());
533 cargs.setComm(*matrices_->parallelInformation().coarsest());
536 coarseSmoother_ = ConstructionTraits<Smoother>::construct(cargs);
537 scalarProduct_ = createScalarProduct<X>(cargs.getComm(),category());
539 typedef DirectSolverSelector< typename M::matrix_type, X > SolverSelector;
542 if( SolverSelector::isDirectSolver &&
543 (std::is_same<ParallelInformation,SequentialInformation>::value
544 || matrices_->parallelInformation().coarsest()->communicator().size()==1
545 || (matrices_->parallelInformation().coarsest().isRedistributed()
546 && matrices_->parallelInformation().coarsest().getRedistributed().communicator().size()==1
547 && matrices_->parallelInformation().coarsest().getRedistributed().communicator().size()>0) )
550 if(matrices_->parallelInformation().coarsest().isRedistributed())
552 if(matrices_->matrices().coarsest().getRedistributed().getmat().N()>0)
555 solver_.reset(SolverSelector::create(matrices_->matrices().coarsest().getRedistributed().getmat(),
false,
false));
562 solver_.reset(SolverSelector::create(matrices_->matrices().coarsest()->getmat(),
false,
false));
564 if(verbosity_>0 && matrices_->parallelInformation().coarsest()->communicator().rank()==0)
565 std::cout<<
"Using a direct coarse solver (" << SolverSelector::name() <<
")" << std::endl;
569 if(matrices_->parallelInformation().coarsest().isRedistributed())
571 if(matrices_->matrices().coarsest().getRedistributed().getmat().N()>0)
576 solver_.reset(
new BiCGSTABSolver<X>(
const_cast<M&
>(matrices_->matrices().coarsest().getRedistributed()),
578 reinterpret_cast<typename
579 std::conditional<std::is_same<PI,SequentialInformation>::value,
580 Dune::SeqScalarProduct<X>,
581 Dune::OverlappingSchwarzScalarProduct<X,PI>
>::type&>(*scalarProduct_),
582 *coarseSmoother_, 1E-2, 1000, 0));
587 solver_.reset(
new BiCGSTABSolver<X>(
const_cast<M&
>(*matrices_->matrices().coarsest()),
589 reinterpret_cast<typename
590 std::conditional<std::is_same<PI,SequentialInformation>::value,
591 Dune::SeqScalarProduct<X>,
592 Dune::OverlappingSchwarzScalarProduct<X,PI>
>::type&>(*scalarProduct_),
593 *coarseSmoother_, 1E-2, 1000, 0));
613 template<
class M,
class X,
class S,
class PI,
class A>
621 typedef typename M::matrix_type Matrix;
622 typedef typename Matrix::ConstRowIterator RowIter;
623 typedef typename Matrix::ConstColIterator ColIter;
624 typedef typename Matrix::block_type Block;
626 zero=
typename Matrix::field_type();
628 const Matrix& mat=matrices_->matrices().finest()->getmat();
629 for(RowIter row=mat.begin(); row!=mat.end(); ++row) {
630 bool isDirichlet =
true;
631 bool hasDiagonal =
false;
633 for(ColIter col=row->begin(); col!=row->end(); ++col) {
634 if(row.index()==col.index()) {
642 if(isDirichlet && hasDiagonal)
643 diagonal.solve(x[row.index()], b[row.index()]);
646 if(smoothers_->levels()>0)
647 smoothers_->finest()->pre(x,b);
650 matrices_->parallelInformation().coarsest()->copyOwnerToAll(x,x);
653 typedef std::shared_ptr< Range > RangePtr ;
654 typedef std::shared_ptr< Domain > DomainPtr;
657 RangePtr copy(
new Range(b) );
658 rhs_.reset(
new Hierarchy<Range,A>(copy) );
659 DomainPtr dcopy(
new Domain(x) );
660 lhs_.reset(
new Hierarchy<Domain,A>(dcopy) );
661 DomainPtr dcopy2(
new Domain(x) );
662 update_.reset(
new Hierarchy<Domain,A>(dcopy2) );
664 matrices_->coarsenVector(*rhs_);
665 matrices_->coarsenVector(*lhs_);
666 matrices_->coarsenVector(*update_);
669 typedef typename Hierarchy<Smoother,A>::Iterator Iterator;
670 typedef typename Hierarchy<Range,A>::Iterator RIterator;
671 typedef typename Hierarchy<Domain,A>::Iterator DIterator;
672 Iterator coarsest = smoothers_->coarsest();
673 Iterator smoother = smoothers_->finest();
674 RIterator rhs = rhs_->finest();
675 DIterator lhs = lhs_->finest();
676 if(smoothers_->levels()>0) {
678 assert(lhs_->levels()==rhs_->levels());
679 assert(smoothers_->levels()==lhs_->levels() || matrices_->levels()==matrices_->maxlevels());
680 assert(smoothers_->levels()+1==lhs_->levels() || matrices_->levels()<matrices_->maxlevels());
682 if(smoother!=coarsest)
683 for(++smoother, ++lhs, ++rhs; smoother != coarsest; ++smoother, ++lhs, ++rhs)
684 smoother->pre(*lhs,*rhs);
685 smoother->pre(*lhs,*rhs);
695 template<
class M,
class X,
class S,
class PI,
class A>
698 return matrices_->levels();
700 template<
class M,
class X,
class S,
class PI,
class A>
703 return matrices_->maxlevels();
707 template<
class M,
class X,
class S,
class PI,
class A>
710 OPM_TIMEBLOCK(apply);
711 LevelContext levelContext;
719 initIteratorsWithFineLevel(levelContext);
722 *levelContext.lhs = v;
723 *levelContext.rhs = d;
724 *levelContext.update=0;
725 levelContext.level=0;
729 if(postSteps_==0||matrices_->maxlevels()==1)
730 levelContext.pinfo->copyOwnerToAll(*levelContext.update, *levelContext.update);
732 v=*levelContext.update;
737 template<
class M,
class X,
class S,
class PI,
class A>
740 OPM_TIMEBLOCK(initIteratorsWithFineLevel);
741 levelContext.smoother = smoothers_->finest();
742 levelContext.matrix = matrices_->matrices().finest();
743 levelContext.pinfo = matrices_->parallelInformation().finest();
744 levelContext.redist =
745 matrices_->redistributeInformation().begin();
746 levelContext.aggregates = matrices_->aggregatesMaps().begin();
747 levelContext.lhs = lhs_->finest();
748 levelContext.update = update_->finest();
749 levelContext.rhs = rhs_->finest();
752 template<
class M,
class X,
class S,
class PI,
class A>
753 bool AMGCPR<M,X,S,PI,A>
754 ::moveToCoarseLevel(LevelContext& levelContext)
756 OPM_TIMEBLOCK(moveToCoarseLevel);
757 bool processNextLevel=
true;
759 if(levelContext.redist->isSetup()) {
760 levelContext.redist->redistribute(
static_cast<const Range&
>(*levelContext.rhs),
761 levelContext.rhs.getRedistributed());
762 processNextLevel = levelContext.rhs.getRedistributed().size()>0;
763 if(processNextLevel) {
765 typename Hierarchy<Range,A>::Iterator fineRhs = levelContext.rhs++;
766 ++levelContext.pinfo;
767 Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
768 ::restrictVector(*(*levelContext.aggregates), *levelContext.rhs,
769 static_cast<const Range&
>(fineRhs.getRedistributed()),
770 *levelContext.pinfo);
774 typename Hierarchy<Range,A>::Iterator fineRhs = levelContext.rhs++;
775 ++levelContext.pinfo;
776 Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
777 ::restrictVector(*(*levelContext.aggregates),
778 *levelContext.rhs,
static_cast<const Range&
>(*fineRhs),
779 *levelContext.pinfo);
782 if(processNextLevel) {
785 ++levelContext.update;
786 ++levelContext.matrix;
787 ++levelContext.level;
788 ++levelContext.redist;
790 if(levelContext.matrix != matrices_->matrices().coarsest() || matrices_->levels()<matrices_->maxlevels()) {
792 ++levelContext.smoother;
793 ++levelContext.aggregates;
796 *levelContext.update=0;
798 return processNextLevel;
801 template<
class M,
class X,
class S,
class PI,
class A>
802 void AMGCPR<M,X,S,PI,A>
803 ::moveToFineLevel(LevelContext& levelContext,
bool processNextLevel)
805 OPM_TIMEBLOCK(moveToFineLevel);
806 if(processNextLevel) {
807 if(levelContext.matrix != matrices_->matrices().coarsest() || matrices_->levels()<matrices_->maxlevels()) {
809 --levelContext.smoother;
810 --levelContext.aggregates;
812 --levelContext.redist;
813 --levelContext.level;
815 --levelContext.matrix;
819 --levelContext.pinfo;
821 if(levelContext.redist->isSetup()) {
823 levelContext.lhs.getRedistributed()=0;
824 Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
825 ::prolongateVector(*(*levelContext.aggregates), *levelContext.update, *levelContext.lhs,
826 levelContext.lhs.getRedistributed(),
827 matrices_->getProlongationDampingFactor(),
828 *levelContext.pinfo, *levelContext.redist);
831 Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
832 ::prolongateVector(*(*levelContext.aggregates), *levelContext.update, *levelContext.lhs,
833 matrices_->getProlongationDampingFactor(),
834 *levelContext.pinfo);
838 if(processNextLevel) {
839 --levelContext.update;
843 *levelContext.update += *levelContext.lhs;
846 template<
class M,
class X,
class S,
class PI,
class A>
849 return IsDirectSolver< CoarseSolver>::value;
852 template<
class M,
class X,
class S,
class PI,
class A>
855 if(levelContext.matrix == matrices_->matrices().coarsest() && levels()==maxlevels()) {
859 if(levelContext.redist->isSetup()) {
860 levelContext.redist->redistribute(*levelContext.rhs, levelContext.rhs.getRedistributed());
861 if(levelContext.rhs.getRedistributed().size()>0) {
863 levelContext.pinfo.getRedistributed().copyOwnerToAll(levelContext.rhs.getRedistributed(),
864 levelContext.rhs.getRedistributed());
865 solver_->apply(levelContext.update.getRedistributed(),
866 levelContext.rhs.getRedistributed(), res);
868 levelContext.redist->redistributeBackward(*levelContext.update, levelContext.update.getRedistributed());
869 levelContext.pinfo->copyOwnerToAll(*levelContext.update, *levelContext.update);
871 levelContext.pinfo->copyOwnerToAll(*levelContext.rhs, *levelContext.rhs);
872 solver_->apply(*levelContext.update, *levelContext.rhs, res);
876 coarsesolverconverged =
false;
879 presmooth(levelContext, preSteps_);
881#ifndef DUNE_AMG_NO_COARSEGRIDCORRECTION
882 bool processNextLevel = moveToCoarseLevel(levelContext);
884 if(processNextLevel) {
886 for(std::size_t i=0; i<gamma_; i++)
890 moveToFineLevel(levelContext, processNextLevel);
895 if(levelContext.matrix == matrices_->matrices().finest()) {
896 coarsesolverconverged = matrices_->parallelInformation().finest()->communicator().prod(coarsesolverconverged);
897 if(!coarsesolverconverged){
902 postsmooth(levelContext, postSteps_);
907 template<
class M,
class X,
class S,
class PI,
class A>
908 void AMGCPR<M,X,S,PI,A>::additiveMgc(){
909 OPM_TIMEBLOCK(additiveMgc);
911 typename ParallelInformationHierarchy::Iterator pinfo=matrices_->parallelInformation().finest();
912 typename Hierarchy<Range,A>::Iterator rhs=rhs_->finest();
913 typename Hierarchy<Domain,A>::Iterator lhs = lhs_->finest();
914 typename OperatorHierarchy::AggregatesMapList::const_iterator aggregates=matrices_->aggregatesMaps().begin();
916 for(
typename Hierarchy<Range,A>::Iterator fineRhs=rhs++; fineRhs != rhs_->coarsest(); fineRhs=rhs++, ++aggregates) {
918 Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
919 ::restrictVector(*(*aggregates), *rhs,
static_cast<const Range&
>(*fineRhs), *pinfo);
925 lhs = lhs_->finest();
926 typename Hierarchy<Smoother,A>::Iterator smoother = smoothers_->finest();
928 for(rhs=rhs_->finest(); rhs != rhs_->coarsest(); ++lhs, ++rhs, ++smoother) {
931 smoother->apply(*lhs, *rhs);
935#ifndef DUNE_AMG_NO_COARSEGRIDCORRECTION
937 pinfo->copyOwnerToAll(*rhs, *rhs);
938 solver_->apply(*lhs, *rhs, res);
941 DUNE_THROW(MathError,
"Coarse solver did not converge");
949 for(
typename Hierarchy<Domain,A>::Iterator coarseLhs = lhs--; coarseLhs != lhs_->finest(); coarseLhs = lhs--, --aggregates, --pinfo) {
950 Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
951 ::prolongateVector(*(*aggregates), *coarseLhs, *lhs, 1.0, *pinfo);
957 template<
class M,
class X,
class S,
class PI,
class A>
962 typedef typename Hierarchy<Smoother,A>::Iterator Iterator;
963 typedef typename Hierarchy<Domain,A>::Iterator DIterator;
964 Iterator coarsest = smoothers_->coarsest();
965 Iterator smoother = smoothers_->finest();
966 DIterator lhs = lhs_->finest();
967 if(smoothers_->levels()>0) {
968 if(smoother != coarsest || matrices_->levels()<matrices_->maxlevels())
969 smoother->post(*lhs);
970 if(smoother!=coarsest)
971 for(++smoother, ++lhs; smoother != coarsest; ++smoother, ++lhs)
972 smoother->post(*lhs);
973 smoother->post(*lhs);
983 template<
class M,
class X,
class S,
class PI,
class A>
987 matrices_->getCoarsestAggregatesOnFinest(cont);
Parallel algebraic multigrid based on agglomeration.
Definition: amgcpr.hh:88
Interface class adding the update() method to the preconditioner interface.
Definition: PreconditionerWithUpdate.hpp:32
M Operator
The matrix operator type.
Definition: amgcpr.hh:96
Hierarchy< Smoother, A >::Iterator smoother
The iterator over the smoothers.
Definition: amgcpr.hh:290
AMGCPR(const OperatorHierarchy &matrices, CoarseSolver &coarseSolver, const SmootherArgs &smootherArgs, const Parameters &parms)
Construct a new amg with a specific coarse solver.
Definition: amgcpr.hh:410
Hierarchy< Domain, A >::Iterator update
The iterator over the updates.
Definition: amgcpr.hh:314
Hierarchy< Range, A >::Iterator rhs
The iterator over the right hand sided.
Definition: amgcpr.hh:318
MatrixHierarchy< M, ParallelInformation, A > OperatorHierarchy
The operator hierarchy type.
Definition: amgcpr.hh:105
ParallelInformationHierarchy::Iterator pinfo
The iterator over the parallel information.
Definition: amgcpr.hh:298
virtual SolverCategory::Category category() const
Category of the preconditioner (see SolverCategory::Category)
Definition: amgcpr.hh:167
void apply(Domain &v, const Range &d)
Definition: amgcpr.hh:708
void pre(Domain &x, Range &b)
Definition: amgcpr.hh:614
OperatorHierarchy::ParallelMatrixHierarchy::ConstIterator matrix
The iterator over the matrices.
Definition: amgcpr.hh:294
void getCoarsestAggregateNumbers(std::vector< std::size_t, A1 > &cont)
Get the aggregate number of each unknown on the coarsest level.
Definition: amgcpr.hh:985
void post(Domain &x)
Definition: amgcpr.hh:958
OperatorHierarchy::ParallelInformationHierarchy ParallelInformationHierarchy
The parallal data distribution hierarchy type.
Definition: amgcpr.hh:107
~AMGCPR()
Definition: amgcpr.hh:452
void recalculateHierarchy()
Recalculate the matrix hierarchy.
Definition: amgcpr.hh:194
SmootherTraits< Smoother >::Arguments SmootherArgs
The argument type for the construction of the smoother.
Definition: amgcpr.hh:123
X Range
The range type.
Definition: amgcpr.hh:112
X Domain
The domain type.
Definition: amgcpr.hh:110
std::size_t maxlevels()
Definition: amgcpr.hh:701
S Smoother
The type of the smoother.
Definition: amgcpr.hh:120
void updateSolver(C &criterion, const Operator &, const PI &pinfo)
Update the coarse solver and the hierarchies.
Definition: amgcpr.hh:464
OperatorHierarchy::AggregatesMapList::const_iterator aggregates
The iterator over the aggregates maps.
Definition: amgcpr.hh:306
virtual void update()
Update the coarse solver and the hierarchies.
Definition: amgcpr.hh:470
PI ParallelInformation
The type of the parallel information. Either OwnerOverlapCommunication or another type describing the...
Definition: amgcpr.hh:103
InverseOperator< X, X > CoarseSolver
the type of the coarse solver.
Definition: amgcpr.hh:114
std::size_t level
The level index.
Definition: amgcpr.hh:322
Hierarchy< Domain, A >::Iterator lhs
The iterator over the left hand side.
Definition: amgcpr.hh:310
bool usesDirectCoarseLevelSolver() const
Check whether the coarse solver used is a direct solver.
Definition: amgcpr.hh:847
Smoother SmootherType
Definition: amgcpr.hh:286
OperatorHierarchy::RedistributeInfoList::const_iterator redist
The iterator over the redistribution information.
Definition: amgcpr.hh:302
std::size_t levels()
Definition: amgcpr.hh:696
void redistributeMatrixAmg(M &, M &, SequentialInformation &, SequentialInformation &, T &)
Definition: amgcpr.hh:35
Definition: SupportsFaceTag.hpp:27
Dune::InverseOperatorResult InverseOperatorResult
Definition: BdaBridge.hpp:32