MatrixFacade.hpp

#pragma once
 
#include <type_traits>
 
#include <dune/common/shared_ptr.hh>
 
#include <amdis/common/Concepts.hpp>
#include <amdis/common/ConceptsBase.hpp>
#include <amdis/common/TypeTraits.hpp>
#include <amdis/functions/NodeIndices.hpp>
#include <amdis/typetree/MultiIndex.hpp>
 
namespace AMDiS
{
  template <class T, template <class> class MatrixImpl>
  class MatrixFacade
  {
    using Self = MatrixFacade;
    using Impl = MatrixImpl<T>;
 
  public:
    template <class RowBasis, class ColBasis>
    MatrixFacade(RowBasis const& rowBasis, ColBasis const& colBasis)
      : impl_(rowBasis.indexDistribution(), colBasis.indexDistribution())
    {}
 
    Impl const& impl() const { return impl_; }
    Impl&       impl()       { return impl_; }
 
    template <class SparsityPattern>
    void init(SparsityPattern const& pattern)
    {
      impl_.init(pattern);
    }
 
    void init()
    {
      impl_.init();
    }
 
    void finish()
    {
      impl_.finish();
    }
 
    template <class RowIndex, class ColIndex,
      REQUIRES(Concepts::MultiIndex<RowIndex>),
      REQUIRES(Concepts::MultiIndex<ColIndex>)>
    void insert(RowIndex const& row, ColIndex const& col, typename Impl::value_type const& value)
    {
      impl_.insert(row, col, value);
    }
 
    template <class RowLocalView, class ColLocalView, class LocalMatrix,
      REQUIRES(Concepts::LocalView<RowLocalView>),
      REQUIRES(Concepts::LocalView<ColLocalView>)>
    void scatter(RowLocalView const& r, ColLocalView const& c, LocalMatrix const& localMatrix)
    {
      assert(r.size() * c.size() == localMatrix.size());
      assert(r.size() == localMatrix.rows());
      assert(c.size() == localMatrix.cols());
 
      const bool optimized = std::is_same_v<RowLocalView,ColLocalView>
        && std::uintptr_t(&r) == std::uintptr_t(&c);
 
      if (optimized)
        impl_.scatter(nodeIndices(r), localMatrix);
      else
        impl_.scatter(nodeIndices(r), nodeIndices(c), localMatrix);
    }
 
    template <class Idx,
      REQUIRES(Concepts::MultiIndex<Idx>)>
    void zeroRows(std::vector<Idx> const& ind, bool diag)
    {
      impl_.zeroRows(ind, diag);
    }
 
    template <class RowIdx, class ColIdx,
      REQUIRES(Concepts::MultiIndex<RowIdx>),
      REQUIRES(Concepts::MultiIndex<ColIdx>)>
    void zeroRows(std::vector<RowIdx> const& rowInd, std::vector<ColIdx> const& colInd, bool diag)
    {
      assert(rowInd.size() == colInd.size());
      impl_.zeroRows(rowInd, colInd, diag);
    }
 
    template <class Idx, class VecX, class VecB,
      REQUIRES(Concepts::MultiIndex<Idx>)>
    void zeroRowsColumns(std::vector<Idx> const& ind, bool diag, VecX const& x, VecB& b)
    {
      impl_.zeroRowsColumns(ind, diag, x.impl(), b.impl());
    }
 
    template <class Idx,
      REQUIRES(Concepts::MultiIndex<Idx>)>
    void zeroRowsColumns(std::vector<Idx> const& ind, bool diag)
    {
      impl_.zeroRowsColumns(ind, diag);
    }
 
    template <class RowIdx, class ColIdx, class VecX, class VecB,
      REQUIRES(Concepts::MultiIndex<RowIdx>),
      REQUIRES(Concepts::MultiIndex<ColIdx>)>
    void zeroRowsColumns(std::vector<RowIdx> const& rowInd, std::vector<ColIdx> const& colInd, bool diag, VecX const& x, VecB& b)
    {
      assert(rowInd.size() == colInd.size());
      impl_.zeroRowsColumns(rowInd, colInd, diag, x.impl(), b.impl());
    }
 
    template <class RowIdx, class ColIdx,
      REQUIRES(Concepts::MultiIndex<RowIdx>),
      REQUIRES(Concepts::MultiIndex<ColIdx>)>
    void zeroRowsColumns(std::vector<RowIdx> const& rowInd, std::vector<ColIdx> const& colInd, bool diag)
    {
      assert(rowInd.size() == colInd.size());
      impl_.zeroRowsColumns(rowInd, colInd, diag);
    }
 
    std::size_t nnz() const
    {
      return impl_.nnz();
    }
 
  protected:
    Impl impl_;
  };
 
} // end namespace AMDiS