geneo_core.hh

 
#include "conbase_fork.hh"
#include "cg_fork.hh"
#include "neumann_boundary_condition.hh"
 
#include "two_level_schwarz.hh"
 
#include "subdomainprojectedcoarsespace.hh"
 
#include "partitionofunity.hh"
#include "localoperator_ovlp_region.hh"
 
#include "geneobasis.hh"
#include "liptonbabuskabasis.hh"
#include "partitionofunitybasis.hh"
#include "zembasis.hh"
 
namespace Dune{
namespace Geneo{
 
template<class GFS, class C, class VV, class CON, template<typename> class Solver>
    class OVLP_GenEO_UMFPack_Base
      : public Dune::PDELab::OVLPScalarProductImplementation<GFS>, 
        public Dune::PDELab::LinearResultStorage
    {
    public:
 
      typedef typename GFS::template ConstraintsContainer<double>::Type CC;
 
 
      OVLP_GenEO_UMFPack_Base (const GFS& gfs_, const C& c_, const CON& constraints_,unsigned maxiter_=5000,
                                              int verbose_=1)
        : Dune::PDELab::OVLPScalarProductImplementation<GFS>(gfs_), 
          gfs(gfs_), c(c_), constraints(constraints_),
          maxiter(maxiter_), 
          verbose(verbose_)
      {
 
        // Constraint Containers
        CC cg;   cg.clear(); 
        CC cg_ext; cg_ext.clear();
 
        // Setup Preconditioner
 
        Dune::PDELab::constraints_exterior(constraints,gfs,cg_ext);
        auto cc_bnd_neu_int_dir = CC();
        Dune::PDELab::PureNeumannBoundaryCondition pnbc;
        cc_bnd_neu_int_dir.clear();
        Dune::PDELab::constraints(pnbc,gfs,cc_bnd_neu_int_dir);
      }
 
      template<class M, class V, class W>
      void apply(M& A, V& z, W& r, typename Dune::template FieldTraits<typename V::ElementType >::real_type reduction)
      {
        typedef Dune::PDELab::OverlappingOperator<C,M,V,W> POP;
        POP pop(c,A);
        typedef Dune::PDELab::OVLPScalarProduct<GFS,V> PSP;
        PSP psp(*this);
#if HAVE_SUITESPARSE_UMFPACK || DOXYGEN
        typedef Dune::PDELab::UMFPackSubdomainSolver<GFS,M,V,W> PREC;
        PREC prec(gfs,A);
        int verb=0;
        if (gfs.gridView().comm().rank()==0) verb=verbose;
        Solver<V> solver(pop,psp,prec,reduction,maxiter,verb);
        Dune::InverseOperatorResult stat;
        solver.apply(z,r,stat);
        res.converged  = stat.converged;
        res.iterations = stat.iterations;
        res.elapsed    = stat.elapsed;
        res.reduction  = stat.reduction;
        res.conv_rate  = stat.conv_rate;
#else
        std::cout << "No UMFPack support, please install and configure it." << std::endl;
#endif
      }
 
    private:
      const GFS& gfs;
      const C& c;
      const CON& constraints;
      unsigned maxiter;
      int verbose;
    };
 
  }
}