amg.hh

Go to the documentation of this file.
00001 // $Id: amg.hh 1181 2010-03-04 09:19:56Z christi $
00002 #ifndef DUNE_AMG_AMG_HH
00003 #define DUNE_AMG_AMG_HH
00004 
00005 #include<memory>
00006 #include<dune/common/exceptions.hh>
00007 #include<dune/istl/paamg/smoother.hh>
00008 #include<dune/istl/paamg/transfer.hh>
00009 #include<dune/istl/paamg/hierarchy.hh>
00010 #include<dune/istl/solvers.hh>
00011 #include<dune/istl/scalarproducts.hh>
00012 #include<dune/istl/superlu.hh>
00013 #include<dune/common/typetraits.hh>
00014 
00015 namespace Dune
00016 {
00017   namespace Amg
00018   {
00043     template<class M, class X, class S, class PI=SequentialInformation,
00044              class A=std::allocator<X> >
00045     class AMG : public Preconditioner<X,X>
00046     {
00047     public:
00049       typedef M Operator;
00056       typedef PI ParallelInformation;
00058       typedef MatrixHierarchy<M, ParallelInformation, A> OperatorHierarchy;
00060       typedef typename OperatorHierarchy::ParallelInformationHierarchy ParallelInformationHierarchy;
00061 
00063       typedef X Domain;
00065       typedef X Range;
00067       typedef InverseOperator<X,X> CoarseSolver;
00073       typedef S Smoother;
00074   
00076       typedef typename SmootherTraits<Smoother>::Arguments SmootherArgs;
00077       
00078       enum {
00080         category = S::category
00081       };
00082 
00094       AMG(const OperatorHierarchy& matrices, CoarseSolver& coarseSolver, 
00095           const SmootherArgs& smootherArgs, std::size_t gamma,
00096           std::size_t preSmoothingSteps,
00097           std::size_t postSmoothingSteps, bool additive=false);
00098 
00113       template<class C>
00114       AMG(const Operator& fineOperator, const C& criterion,
00115           const SmootherArgs& smootherArgs, std::size_t gamma=1,
00116           std::size_t preSmoothingSteps=2, std::size_t postSmoothingSteps=2,
00117           bool additive=false, const ParallelInformation& pinfo=ParallelInformation());
00118 
00119       ~AMG();
00120 
00122       void pre(Domain& x, Range& b);
00123 
00125       void apply(Domain& v, const Range& d);
00126       
00128       void post(Domain& x);
00129 
00134       template<class A1>
00135       void getCoarsestAggregateNumbers(std::vector<std::size_t,A1>& cont);
00136       
00137       std::size_t levels();
00138       
00139       std::size_t maxlevels();
00140     private:
00142       void mgc(typename Hierarchy<Smoother,A>::Iterator& smoother,
00143                typename OperatorHierarchy::ParallelMatrixHierarchy::ConstIterator& matrix,
00144                typename ParallelInformationHierarchy::Iterator& pinfo,
00145                typename OperatorHierarchy::RedistributeInfoList::const_iterator& redist,
00146                typename OperatorHierarchy::AggregatesMapList::const_iterator& aggregates,
00147                typename Hierarchy<Domain,A>::Iterator& lhs,
00148                typename Hierarchy<Domain,A>::Iterator& update,
00149                typename Hierarchy<Range,A>::Iterator& rhs);
00150 
00151       void additiveMgc();
00152       
00153       //      void setupIndices(typename Matrix::ParallelIndexSet& indices, const Matrix& matrix);
00154       
00156       const OperatorHierarchy* matrices_;
00158       SmootherArgs smootherArgs_;
00160       Hierarchy<Smoother,A> smoothers_;
00162       CoarseSolver* solver_;
00164       Hierarchy<Range,A>* rhs_;
00166       Hierarchy<Domain,A>* lhs_;
00168       Hierarchy<Domain,A>* update_;
00170       typedef Dune::ScalarProductChooser<X,PI,M::category> ScalarProductChooser;
00172       typedef typename ScalarProductChooser::ScalarProduct ScalarProduct;
00174       ScalarProduct* scalarProduct_;
00176       std::size_t gamma_;
00178       std::size_t preSteps_;
00180       std::size_t postSteps_;
00181       std::size_t level;
00182       bool buildHierarchy_;
00183       bool additive;
00184       Smoother *coarseSmoother_;
00185     };
00186 
00187     template<class M, class X, class S, class P, class A>
00188     AMG<M,X,S,P,A>::AMG(const OperatorHierarchy& matrices, CoarseSolver& coarseSolver, 
00189                         const SmootherArgs& smootherArgs,
00190                         std::size_t gamma, std::size_t preSmoothingSteps,
00191                         std::size_t postSmoothingSteps, bool additive_)
00192       : matrices_(&matrices), smootherArgs_(smootherArgs),
00193         smoothers_(), solver_(&coarseSolver), scalarProduct_(0),
00194         gamma_(gamma), preSteps_(preSmoothingSteps), postSteps_(postSmoothingSteps), buildHierarchy_(false),
00195         additive(additive_), coarseSmoother_()
00196     {
00197       assert(matrices_->isBuilt());
00198       
00199       // build the necessary smoother hierarchies
00200       matrices_->coarsenSmoother(smoothers_, smootherArgs_);
00201     }
00202 
00203     template<class M, class X, class S, class P, class A>
00204     template<class C>
00205     AMG<M,X,S,P,A>::AMG(const Operator& matrix,
00206                         const C& criterion,
00207                         const SmootherArgs& smootherArgs,
00208                         std::size_t gamma, std::size_t preSmoothingSteps,
00209                         std::size_t postSmoothingSteps,
00210                         bool additive_,
00211                         const P& pinfo)
00212       : smootherArgs_(smootherArgs),
00213         smoothers_(), solver_(), scalarProduct_(0), gamma_(gamma),
00214         preSteps_(preSmoothingSteps), postSteps_(postSmoothingSteps), buildHierarchy_(true),
00215         additive(additive_), coarseSmoother_()
00216     {
00217       dune_static_assert(static_cast<int>(M::category)==static_cast<int>(S::category),
00218                          "Matrix and Solver must match in terms of category!");
00219       dune_static_assert(static_cast<int>(P::category)==static_cast<int>(S::category),
00220                          "Matrix and Solver must match in terms of category!");
00221       Timer watch;
00222       OperatorHierarchy* matrices = new OperatorHierarchy(const_cast<Operator&>(matrix), pinfo);
00223             
00224       matrices->template build<NegateSet<typename P::OwnerSet> >(criterion);
00225       
00226       matrices_ = matrices;
00227       // build the necessary smoother hierarchies
00228       matrices_->coarsenSmoother(smoothers_, smootherArgs_);
00229 
00230       if(criterion.debugLevel()>0 && matrices_->parallelInformation().finest()->communicator().rank()==0)
00231         std::cout<<"Building Hierarchy of "<<matrices_->maxlevels()<<" levels took "<<watch.elapsed()<<" seconds."<<std::endl;
00232     }
00233     
00234     template<class M, class X, class S, class P, class A>
00235     AMG<M,X,S,P,A>::~AMG()
00236     {
00237       if(buildHierarchy_){
00238         delete matrices_;
00239       }
00240       if(scalarProduct_)
00241         delete scalarProduct_;
00242     }
00243     
00244     template<class M, class X, class S, class P, class A>
00245     void AMG<M,X,S,P,A>::pre(Domain& x, Range& b)
00246     {
00247       if(smoothers_.levels()>0)
00248         smoothers_.finest()->pre(x,b);
00249       else
00250         // No smoother to make x consistent! Do it by hand
00251         matrices_->parallelInformation().coarsest()->copyOwnerToAll(x,x);
00252       Range* copy = new Range(b);
00253       rhs_ = new Hierarchy<Range,A>(*copy);
00254       Domain* dcopy = new Domain(x);
00255       lhs_ = new Hierarchy<Domain,A>(*dcopy);
00256       dcopy = new Domain(x);
00257       update_ = new Hierarchy<Domain,A>(*dcopy);
00258       matrices_->coarsenVector(*rhs_);
00259       matrices_->coarsenVector(*lhs_);
00260       matrices_->coarsenVector(*update_);
00261       
00262       // Preprocess all smoothers
00263       typedef typename Hierarchy<Smoother,A>::Iterator Iterator;
00264       typedef typename Hierarchy<Range,A>::Iterator RIterator;
00265       typedef typename Hierarchy<Domain,A>::Iterator DIterator;
00266       Iterator coarsest = smoothers_.coarsest();
00267       Iterator smoother = smoothers_.finest();
00268       RIterator rhs = rhs_->finest();
00269       DIterator lhs = lhs_->finest();
00270       if(smoothers_.levels()>0){
00271           
00272       assert(lhs_->levels()==rhs_->levels());
00273       assert(smoothers_.levels()==lhs_->levels() || matrices_->levels()==matrices_->maxlevels());
00274       assert(smoothers_.levels()+1==lhs_->levels() || matrices_->levels()<matrices_->maxlevels());
00275       
00276       if(smoother!=coarsest)
00277         for(++smoother, ++lhs, ++rhs; smoother != coarsest; ++smoother, ++lhs, ++rhs)
00278           smoother->pre(*lhs,*rhs);
00279       smoother->pre(*lhs,*rhs);
00280       }
00281       
00282       
00283       // The preconditioner might change x and b. So we have to 
00284       // copy the changes to the original vectors.
00285       x = *lhs_->finest();
00286       b = *rhs_->finest();
00287       
00288       if(buildHierarchy_ && matrices_->levels()==matrices_->maxlevels()){
00289         // We have the carsest level. Create the coarse Solver
00290         SmootherArgs sargs(smootherArgs_);
00291         sargs.iterations = 1;
00292         
00293         typename ConstructionTraits<Smoother>::Arguments cargs;
00294         cargs.setArgs(sargs);
00295         if(matrices_->redistributeInformation().back().isSetup()){
00296           // Solve on the redistributed partitioning     
00297           cargs.setMatrix(matrices_->matrices().coarsest().getRedistributed().getmat());
00298           cargs.setComm(matrices_->parallelInformation().coarsest().getRedistributed());
00299           
00300           coarseSmoother_ = ConstructionTraits<Smoother>::construct(cargs);
00301           scalarProduct_ = ScalarProductChooser::construct(matrices_->parallelInformation().coarsest().getRedistributed());
00302         }else{    
00303           cargs.setMatrix(matrices_->matrices().coarsest()->getmat());
00304           cargs.setComm(*matrices_->parallelInformation().coarsest());
00305           
00306           coarseSmoother_ = ConstructionTraits<Smoother>::construct(cargs);
00307           scalarProduct_ = ScalarProductChooser::construct(*matrices_->parallelInformation().coarsest());
00308         }
00309 #ifdef HAVE_SUPERLU
00310       // Use superlu if we are purely sequential or with only one processor on the coarsest level.
00311         if(is_same<ParallelInformation,SequentialInformation>::value // sequential mode 
00312            || matrices_->parallelInformation().coarsest()->communicator().size()==1 //parallel mode and only one processor
00313            || (matrices_->parallelInformation().coarsest().isRedistributed() 
00314                && matrices_->parallelInformation().coarsest().getRedistributed().communicator().size()==1
00315                && matrices_->parallelInformation().coarsest().getRedistributed().communicator().size()>0)){ // redistribute and 1 proc
00316           std::cout<<"Using superlu"<<std::endl;
00317           if(matrices_->parallelInformation().coarsest().isRedistributed())
00318             {
00319               if(matrices_->matrices().coarsest().getRedistributed().getmat().N()>0)
00320                 // We are still participating on this level
00321                 solver_  = new SuperLU<typename M::matrix_type>(matrices_->matrices().coarsest().getRedistributed().getmat());
00322               else
00323                 solver_ = 0;
00324             }else
00325               solver_  = new SuperLU<typename M::matrix_type>(matrices_->matrices().coarsest()->getmat());
00326         }else
00327 #endif
00328           {
00329             if(matrices_->parallelInformation().coarsest().isRedistributed())
00330               {
00331                 if(matrices_->matrices().coarsest().getRedistributed().getmat().N()>0)
00332                   // We are still participating on this level
00333                   solver_ = new BiCGSTABSolver<X>(const_cast<M&>(matrices_->matrices().coarsest().getRedistributed()), 
00334                                                   *scalarProduct_, 
00335                                                   *coarseSmoother_, 1E-2, 10000, 0);
00336                 else
00337                   solver_ = 0;
00338               }else
00339               solver_ = new BiCGSTABSolver<X>(const_cast<M&>(*matrices_->matrices().coarsest()), 
00340                                               *scalarProduct_, 
00341                                               *coarseSmoother_, 1E-2, 10000, 0);
00342           }
00343       }
00344     }
00345     template<class M, class X, class S, class P, class A>
00346     std::size_t AMG<M,X,S,P,A>::levels()
00347     {
00348       return matrices_->levels();
00349     }
00350     template<class M, class X, class S, class P, class A>
00351     std::size_t AMG<M,X,S,P,A>::maxlevels()
00352     {
00353       return matrices_->maxlevels();
00354     }
00355 
00357     template<class M, class X, class S, class P, class A>
00358     void AMG<M,X,S,P,A>::apply(Domain& v, const Range& d)
00359     {
00360       if(additive){
00361         *(rhs_->finest())=d;
00362         additiveMgc();
00363         v=*lhs_->finest();
00364       }else{
00365         typename Hierarchy<Smoother,A>::Iterator smoother = smoothers_.finest();
00366         typename OperatorHierarchy::ParallelMatrixHierarchy::ConstIterator matrix = matrices_->matrices().finest();
00367         typename ParallelInformationHierarchy::Iterator pinfo = matrices_->parallelInformation().finest();
00368         typename OperatorHierarchy::RedistributeInfoList::const_iterator redist = 
00369           matrices_->redistributeInformation().begin();
00370         typename OperatorHierarchy::AggregatesMapList::const_iterator aggregates = matrices_->aggregatesMaps().begin();
00371         typename Hierarchy<Domain,A>::Iterator lhs = lhs_->finest();
00372         typename Hierarchy<Domain,A>::Iterator update = update_->finest();
00373         typename Hierarchy<Range,A>::Iterator rhs = rhs_->finest();
00374         
00375         *lhs = v;
00376         *rhs = d;
00377         *update=0;
00378         level=0;
00379                   
00380         mgc(smoother, matrix, pinfo, redist, aggregates, lhs, update, rhs);
00381         
00382         if(postSteps_==0||matrices_->maxlevels()==1)
00383           pinfo->copyOwnerToAll(*update, *update);
00384         
00385         v=*update;
00386       }
00387       
00388     }
00389     
00390     template<class M, class X, class S, class P, class A>
00391     void AMG<M,X,S,P,A>::mgc(typename Hierarchy<Smoother,A>::Iterator& smoother,
00392                              typename OperatorHierarchy::ParallelMatrixHierarchy::ConstIterator& matrix,
00393                              typename ParallelInformationHierarchy::Iterator& pinfo,
00394                              typename OperatorHierarchy::RedistributeInfoList::const_iterator& redist,
00395                              typename OperatorHierarchy::AggregatesMapList::const_iterator& aggregates,
00396                              typename Hierarchy<Domain,A>::Iterator& lhs,
00397                              typename Hierarchy<Domain,A>::Iterator& update,
00398                              typename Hierarchy<Range,A>::Iterator& rhs){
00399       if(matrix == matrices_->matrices().coarsest() && levels()==maxlevels()){
00400         // Solve directly
00401         InverseOperatorResult res;
00402         res.converged=true; // If we do not compute this flag will not get updated
00403         if(redist->isSetup()){
00404             redist->redistribute(*rhs, rhs.getRedistributed());
00405           if(rhs.getRedistributed().size()>0){
00406             // We are still participating in the computation
00407             pinfo.getRedistributed().copyOwnerToAll(rhs.getRedistributed(), rhs.getRedistributed());
00408             if(maxlevels()==1)
00409               // prepare for iterativr solver
00410               update.getRedistributed()=0.1;
00411             solver_->apply(update.getRedistributed(), rhs.getRedistributed(), res);
00412           }
00413           redist->redistributeBackward(*update, update.getRedistributed());
00414           pinfo->copyOwnerToAll(*update, *update);
00415         }else{
00416           pinfo->copyOwnerToAll(*rhs, *rhs);
00417           if(maxlevels()==1)
00418               // prepare for iterativr solver
00419               *update=0.1;
00420           solver_->apply(*update, *rhs, res);
00421         }
00422 
00423         if(!res.converged)
00424           DUNE_THROW(MathError, "Coarse solver did not converge");
00425       }else{
00426         // presmoothing
00427           for(std::size_t i=0; i < preSteps_; ++i){
00428             *lhs=0;
00429             SmootherApplier<S>::preSmooth(*smoother, *lhs, *rhs);
00430             // Accumulate update
00431             *update += *lhs;
00432             
00433             // update defect
00434             matrix->applyscaleadd(-1,static_cast<const Domain&>(*lhs), *rhs);
00435             pinfo->project(*rhs);
00436           }
00437 
00438 #ifndef DUNE_AMG_NO_COARSEGRIDCORRECTION
00439 
00440         bool processNextLevel=true;
00441         
00442         if(redist->isSetup()){
00443           redist->redistribute(static_cast<const Range&>(*rhs), rhs.getRedistributed());
00444           processNextLevel =rhs.getRedistributed().size()>0;
00445           if(processNextLevel){         
00446             //restrict defect to coarse level right hand side.
00447             typename Hierarchy<Range,A>::Iterator fineRhs = rhs++;
00448             ++pinfo;
00449             Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
00450               ::restrict(*(*aggregates), *rhs, static_cast<const Range&>(fineRhs.getRedistributed()), *pinfo);
00451           }
00452         }else{          
00453           //restrict defect to coarse level right hand side.
00454           typename Hierarchy<Range,A>::Iterator fineRhs = rhs++;
00455           ++pinfo;
00456           Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
00457             ::restrict(*(*aggregates), *rhs, static_cast<const Range&>(*fineRhs), *pinfo);
00458         }
00459           
00460         if(processNextLevel){
00461           // prepare coarse system
00462           ++lhs;
00463           ++update;
00464           ++matrix;
00465           ++level;
00466           ++redist;
00467 
00468           if(matrix != matrices_->matrices().coarsest() || matrices_->levels()<matrices_->maxlevels()){
00469             // next level is not the globally coarsest one
00470             ++smoother;
00471             ++aggregates;
00472           }
00473         
00474           // prepare the update on the next level
00475           *update=0;
00476         
00477           // next level
00478           for(std::size_t i=0; i<gamma_; i++)
00479             mgc(smoother, matrix, pinfo, redist, aggregates, lhs, update, rhs);
00480 
00481           if(matrix != matrices_->matrices().coarsest() || matrices_->levels()<matrices_->maxlevels()){
00482             // previous level is not the globally coarsest one
00483             --smoother;
00484             --aggregates;
00485           }
00486           --redist;
00487           --level;
00488           //prolongate and add the correction (update is in coarse left hand side)
00489           --matrix;
00490         
00491           //typename Hierarchy<Domain,A>::Iterator coarseLhs = lhs--;
00492           --lhs;  
00493           --pinfo;
00494         }
00495         
00496         if(redist->isSetup()){
00497           // Need to redistribute during prolongate
00498           lhs.getRedistributed()=0;
00499           Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
00500             ::prolongate(*(*aggregates), *update, *lhs, lhs.getRedistributed(), matrices_->getProlongationDampingFactor(), 
00501                          *pinfo, *redist);
00502           }else{
00503             *lhs=0;
00504             Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
00505               ::prolongate(*(*aggregates), *update, *lhs, 
00506                            matrices_->getProlongationDampingFactor(), *pinfo);
00507           }
00508         
00509         
00510         if(processNextLevel){
00511           --update;
00512           --rhs;
00513         }
00514 
00515         *update += *lhs;
00516         
00517 
00518 #endif  
00519         
00520         // postsmoothing
00521         for(std::size_t i=0; i < preSteps_; ++i){
00522           // update defect
00523           matrix->applyscaleadd(-1,static_cast<const Domain&>(*lhs), *rhs);
00524           *lhs=0;
00525           pinfo->project(*rhs);
00526           SmootherApplier<S>::preSmooth(*smoother, *lhs, *rhs);
00527           // Accumulate update
00528           *update += *lhs;
00529         }       
00530       }
00531     }
00532 
00533     template<class M, class X, class S, class P, class A>
00534     void AMG<M,X,S,P,A>::additiveMgc(){
00535       
00536       // restrict residual to all levels
00537       typename ParallelInformationHierarchy::Iterator pinfo=matrices_->parallelInformation().finest();
00538       typename Hierarchy<Range,A>::Iterator rhs=rhs_->finest();      
00539       typename Hierarchy<Domain,A>::Iterator lhs = lhs_->finest();
00540       typename OperatorHierarchy::AggregatesMapList::const_iterator aggregates=matrices_->aggregatesMaps().begin();
00541       
00542       for(typename Hierarchy<Range,A>::Iterator fineRhs=rhs++; fineRhs != rhs_->coarsest(); fineRhs=rhs++, ++aggregates){
00543         ++pinfo;
00544         Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
00545           ::restrict(*(*aggregates), *rhs, static_cast<const Range&>(*fineRhs), *pinfo);
00546       }
00547       
00548       // pinfo is invalid, set to coarsest level
00549       //pinfo = matrices_->parallelInformation().coarsest
00550       // calculate correction for all levels
00551       lhs = lhs_->finest();
00552       typename Hierarchy<Smoother,A>::Iterator smoother = smoothers_.finest();
00553       
00554       for(rhs=rhs_->finest(); rhs != rhs_->coarsest(); ++lhs, ++rhs, ++smoother){
00555         // presmoothing
00556         *lhs=0;
00557         smoother->apply(*lhs, *rhs);
00558       }
00559       
00560       // Coarse level solve
00561 #ifndef DUNE_AMG_NO_COARSEGRIDCORRECTION 
00562       InverseOperatorResult res;
00563       pinfo->copyOwnerToAll(*rhs, *rhs);
00564       solver_->apply(*lhs, *rhs, res);
00565       
00566       if(!res.converged)
00567         DUNE_THROW(MathError, "Coarse solver did not converge");
00568 #else
00569       *lhs=0;
00570 #endif
00571       // Prologate and add up corrections from all levels
00572       --pinfo;
00573       --aggregates;
00574       
00575       for(typename Hierarchy<Domain,A>::Iterator coarseLhs = lhs--; coarseLhs != lhs_->finest(); coarseLhs = lhs--, --aggregates, --pinfo){
00576         Transfer<typename OperatorHierarchy::AggregatesMap::AggregateDescriptor,Range,ParallelInformation>
00577           ::prolongate(*(*aggregates), *coarseLhs, *lhs, 1, *pinfo);
00578       }
00579     }
00580 
00581     
00583     template<class M, class X, class S, class P, class A>
00584     void AMG<M,X,S,P,A>::post(Domain& x)
00585     {
00586       if(buildHierarchy_){
00587         if(solver_)
00588           delete solver_;
00589         if(coarseSmoother_)
00590           ConstructionTraits<Smoother>::deconstruct(coarseSmoother_);
00591       }
00592       
00593       // Postprocess all smoothers
00594       typedef typename Hierarchy<Smoother,A>::Iterator Iterator;
00595       typedef typename Hierarchy<Range,A>::Iterator RIterator;
00596       typedef typename Hierarchy<Domain,A>::Iterator DIterator;
00597       Iterator coarsest = smoothers_.coarsest();
00598       Iterator smoother = smoothers_.finest();
00599       DIterator lhs = lhs_->finest();
00600       if(smoothers_.levels()>0){
00601         if(smoother != coarsest  || matrices_->levels()<matrices_->maxlevels())
00602           smoother->post(*lhs);
00603         if(smoother!=coarsest)
00604           for(++smoother, ++lhs; smoother != coarsest; ++smoother, ++lhs)
00605             smoother->post(*lhs);
00606         smoother->post(*lhs);
00607       }
00608 
00609       delete &(*lhs_->finest());
00610       delete lhs_;
00611       delete &(*update_->finest());
00612       delete update_;
00613       delete &(*rhs_->finest());
00614       delete rhs_;
00615     }
00616 
00617     template<class M, class X, class S, class P, class A>
00618     template<class A1>
00619     void AMG<M,X,S,P,A>::getCoarsestAggregateNumbers(std::vector<std::size_t,A1>& cont)
00620     {
00621       matrices_->getCoarsestAggregatesOnFinest(cont);
00622     }
00623     
00624   } // end namespace Amg
00625 }// end namespace Dune
00626   
00627 #endif
Generated on Sat Apr 24 11:13:45 2010 for dune-istl by  doxygen 1.6.3