MatrixNnzStructure.hpp

#pragma once

#include <algorithm>
#include <vector>

#include <amdis/common/Index.hpp>
#include <amdis/linearalgebra/SymmetryStructure.hpp>
#include <amdis/linearalgebra/petsc/IndexDistribution.hpp>

#if HAVE_MPI
#include <amdis/common/parallel/Communicator.hpp>
#endif

namespace AMDiS
{
  class MatrixNnzStructure
  {
  public:
    template <class RowBasis, class ColBasis>
    MatrixNnzStructure(RowBasis const& rowBasis, ColBasis const& colBasis,
                       SymmetryStructure symmetry = SymmetryStructure::unknown)
      : symmetry_(symmetry)
    {
      init(rowBasis, rowBasis.indexDistribution(), colBasis, colBasis.indexDistribution());
    }

    // Return Number of nonzeros in the diagonal part (owner part)
    std::vector<PetscInt> const& d_nnz() const
    {
      return dnnz_;
    }

    std::vector<PetscInt> const& o_nnz() const
    {
      return onnz_;
    }

    SymmetryStructure symmetry() const
    {
      return symmetry_;
    }

  protected:
    // Create local pattern from bases
    template <class RowBasis, class ColBasis, class LI>
    void init(RowBasis const& rowBasis, PetscSequentialIndexDistribution<LI> const& rowDofMap,
              ColBasis const& colBasis, PetscSequentialIndexDistribution<LI> const& colDofMap);

    // Create distributed pattern from bases
    template <class RowBasis, class ColBasis, class GID, class LI>
    void init(RowBasis const& rowBasis, PetscParallelIndexDistribution<GID,LI> const& rowDofMap,
              ColBasis const& colBasis, PetscParallelIndexDistribution<GID,LI> const& colDofMap);
  private:
    std::vector<PetscInt> dnnz_; //< number of nonzeros in the diagonal part (owner part)
    std::vector<PetscInt> onnz_; //< number of nonzeros in the off-diagonal part (overlap part)
    SymmetryStructure symmetry_;

#if HAVE_MPI
    const Mpi::Tag tag_{916821}; //< Communication tag used internally
#endif
  };

} // end namespace AMDiS

#include <amdis/linearalgebra/petsc/MatrixNnzStructure.inc.hpp>