Preconditioners.hpp

#pragma once
 
#include <type_traits>
 
#include <dune/istl/preconditioners.hh>
#include <dune/istl/novlpschwarz.hh>
#include <dune/istl/schwarz.hh>
 
#include <amdis/CreatorInterface.hpp>
#include <amdis/CreatorMap.hpp>
#include <amdis/common/TypeTraits.hpp>
#include <amdis/linearalgebra/istl/ISTLPreconCreator.hpp>
 
#if AMDIS_ENABLE_PAMG
#include <amdis/linearalgebra/istl/AMGPrecon.hpp>
#endif
 
namespace AMDiS
{
 
  template <class Traits>
  class DefaultCreators<tag::preconditioner<Traits>>
  {
    using M = typename Traits::M;
    using X = typename Traits::X;
    using Y = typename Traits::Y;
 
    template <class Precon>
    using Creator = ISTLPreconCreator<Precon, Traits>;
 
    using Map = CreatorMap<tag::preconditioner<Traits>>;
    using FTraits = Dune::FieldTraits<typename M::field_type>;
 
  public:
    static void init()
    {
      auto jacobi = new Creator<Dune::SeqJac<M,X,Y>>;
      Map::addCreator("diag", jacobi);
      Map::addCreator("jacobi", jacobi);
 
      auto gs = new Creator<Dune::SeqGS<M,X,Y>>;
      Map::addCreator("gs", gs);
      Map::addCreator("gauss_seidel", gs);
 
      auto sor = new Creator<Dune::SeqSOR<M,X,Y>>;
      Map::addCreator("sor", sor);
 
      auto ssor = new Creator<Dune::SeqSSOR<M,X,Y>>;
      Map::addCreator("ssor", ssor);
 
      init_ilu(std::is_arithmetic<typename FTraits::field_type>{});
 
#if AMDIS_ENABLE_PAMG
      init_amg(std::is_same<typename FTraits::real_type, double>{});
#endif
 
      auto richardson = new Creator<Dune::Richardson<X,Y>>;
      Map::addCreator("richardson", richardson);
      Map::addCreator("default", richardson);
 
      auto solver = new Creator<tag::solver<Traits>>;
      Map::addCreator("solver", solver);
 
      init_bjacobi(Types<TYPEOF(std::declval<typename Traits::Comm>().impl())>{}, Dune::PriorityTag<10>{});
    }
 
    static void init_ilu(std::false_type)
    {
      warning("ILU preconditioners not created for the matrix with field_type = {}.",
        Dune::className<typename FTraits::field_type>());
    }
 
    static void init_ilu(std::true_type)
    {
      auto ilu = new Creator<Dune::SeqILU<M,X,Y>>;
      Map::addCreator("ilu", ilu);
      Map::addCreator("ilu0", ilu);
 
      auto ildl = new Creator<Dune::SeqILDL<M,X,Y>>;
      Map::addCreator("ildl", ildl);
    }
 
    static void init_amg(std::false_type)
    {
      warning("AMG preconditioners not created for the matrix with real_type = {}.",
        Dune::className<typename FTraits::real_type>());
    }
 
#if AMDIS_ENABLE_PAMG
    static void init_amg(std::true_type)
    {
      auto amg = new AMGPreconCreator<Dune::Amg::AMG, Traits>;
      Map::addCreator("amg", amg);
      auto fastamg = new AMGPreconCreator<Dune::Amg::FastAMG, Traits>;
      Map::addCreator("fastamg", fastamg);
      auto kamg = new AMGPreconCreator<Dune::Amg::KAMG, Traits>;
      Map::addCreator("kamg", kamg);
    }
#endif
 
    static void init_bjacobi(Types<Dune::Amg::SequentialInformation>, Dune::PriorityTag<2>) {}
 
    template <class Comm>
    static void init_bjacobi(Types<Comm>, Dune::PriorityTag<1>)
    {
      auto pssor = new Creator<Dune::ParSSOR<M,X,Y,Comm>>;
      Map::addCreator("pssor", pssor);
 
      auto bjacobi = new Creator<tag::bjacobi>;
      Map::addCreator("bjacobi", bjacobi);
    }
  };
 
} // end namespace AMDiS