SolverWrapper.hpp

#pragma once

#include <memory>
#include <utility>

#include <dune/common/shared_ptr.hh>
#include <dune/istl/operators.hh>
#include <dune/istl/preconditioner.hh>
#include <dune/istl/scalarproducts.hh>
#include <dune/istl/solver.hh>
#include <dune/istl/solvercategory.hh>

#include <amdis/common/TypeTraits.hpp>

namespace AMDiS
{
  template <class S>
  class IterativeSolverWrapper
      : public Dune::InverseOperator<typename S::domain_type, typename S::range_type>
  {
    using Solver = S;

  public:
    using domain_type = typename S::domain_type;
    using range_type = typename S::range_type;

    using LinOp = Dune::LinearOperator<domain_type,range_type>;
    using Prec = Dune::Preconditioner<domain_type,range_type>;
    using ScalProd = Dune::ScalarProduct<domain_type>;

  public:
    template <class... Args>
    IterativeSolverWrapper(LinOp& linOp, Prec& prec, Args&&... args)
      : linOp_(Dune::stackobject_to_shared_ptr(linOp))
      , prec_(Dune::stackobject_to_shared_ptr(prec))
      , solver_(*linOp_, *prec_, FWD(args)...)
    {}

    template <class... Args>
    IterativeSolverWrapper(std::shared_ptr<LinOp> linOp, std::shared_ptr<Prec> prec, Args&&... args)
      : linOp_(std::move(linOp))
      , prec_(std::move(prec))
      , solver_(*linOp_, *prec_, FWD(args)...)
    {}

    template <class... Args>
    IterativeSolverWrapper(LinOp& linOp, ScalProd& sp, Prec& prec, Args&&... args)
      : linOp_(Dune::stackobject_to_shared_ptr(linOp))
      , sp_(Dune::stackobject_to_shared_ptr(sp))
      , prec_(Dune::stackobject_to_shared_ptr(prec))
      , solver_(*linOp_, *sp_, *prec_, FWD(args)...)
    {}

    template <class... Args>
    IterativeSolverWrapper(std::shared_ptr<LinOp> linOp, std::shared_ptr<ScalProd> sp, std::shared_ptr<Prec> prec, Args&&... args)
      : linOp_(std::move(linOp))
      , sp_(std::move(sp))
      , prec_(std::move(prec))
      , solver_(*linOp_, *sp_, *prec_, FWD(args)...)
    {}

    void apply(domain_type& x, range_type& b, Dune::InverseOperatorResult& res) override
    {
      solver_.apply(x, b, res);
    }

    void apply(domain_type& x, range_type& b, double reduction, Dune::InverseOperatorResult& res) override
    {
      solver_.apply(x, b, reduction, res);
    }

    Dune::SolverCategory::Category category() const override
    {
      return solver_.category();
    }

  private:
    std::shared_ptr<LinOp> linOp_;
    std::shared_ptr<ScalProd> sp_;
    std::shared_ptr<Prec> prec_;

    Solver solver_;
  };

} // end namespace AMDiS