petsc_parallel_block_sparse_matrix.h

// ---------------------------------------------------------------------
//
// Copyright (C) 2004 - 2017 by the deal.II authors
//
// This file is part of the deal.II library.
//
// The deal.II library is free software; you can use it, redistribute
// it, and/or modify it under the terms of the GNU Lesser General
// Public License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// The full text of the license can be found in the file LICENSE at
// the top level of the deal.II distribution.
//
// ---------------------------------------------------------------------

#ifndef dealii_petsc_parallel_block_sparse_matrix_h
#define dealii_petsc_parallel_block_sparse_matrix_h


#include <deal.II/base/config.h>

#ifdef DEAL_II_WITH_PETSC

#  include <deal.II/base/table.h>
#  include <deal.II/lac/block_matrix_base.h>
#  include <deal.II/lac/block_sparsity_pattern.h>
#  include <deal.II/lac/petsc_parallel_sparse_matrix.h>
#  include <deal.II/lac/petsc_parallel_block_vector.h>
#  include <deal.II/lac/exceptions.h>
#  include <cmath>

DEAL_II_NAMESPACE_OPEN



namespace PETScWrappers
{
  namespace MPI
  {

    class BlockSparseMatrix : public BlockMatrixBase<SparseMatrix>
    {
    public:
      typedef BlockMatrixBase<SparseMatrix> BaseClass;

      typedef BaseClass::BlockType  BlockType;

      typedef BaseClass::value_type      value_type;
      typedef BaseClass::pointer         pointer;
      typedef BaseClass::const_pointer   const_pointer;
      typedef BaseClass::reference       reference;
      typedef BaseClass::const_reference const_reference;
      typedef BaseClass::size_type       size_type;
      typedef BaseClass::iterator        iterator;
      typedef BaseClass::const_iterator  const_iterator;

      BlockSparseMatrix () = default;

      ~BlockSparseMatrix () = default;

      BlockSparseMatrix &
      operator = (const BlockSparseMatrix &);

      BlockSparseMatrix &
      operator = (const double d);

      void reinit (const size_type n_block_rows,
                   const size_type n_block_columns);


      void reinit(const std::vector<IndexSet> &rows,
                  const std::vector<IndexSet> &cols,
                  const BlockDynamicSparsityPattern &bdsp,
                  const MPI_Comm &com);


      void reinit(const std::vector<IndexSet> &sizes,
                  const BlockDynamicSparsityPattern &bdsp,
                  const MPI_Comm &com);



      void vmult (BlockVector       &dst,
                  const BlockVector &src) const;

      void vmult (BlockVector          &dst,
                  const Vector &src) const;

      void vmult (Vector    &dst,
                  const BlockVector &src) const;

      void vmult (Vector       &dst,
                  const Vector &src) const;

      void Tvmult (BlockVector       &dst,
                   const BlockVector &src) const;

      void Tvmult (BlockVector  &dst,
                   const Vector &src) const;

      void Tvmult (Vector    &dst,
                   const BlockVector &src) const;

      void Tvmult (Vector       &dst,
                   const Vector &src) const;

      void collect_sizes ();

      std::vector< IndexSet > locally_owned_domain_indices() const;

      std::vector< IndexSet > locally_owned_range_indices() const;

      const MPI_Comm &get_mpi_communicator () const;

      using BlockMatrixBase<SparseMatrix>::clear;
    };



// ------------- inline and template functions -----------------

    inline
    BlockSparseMatrix &
    BlockSparseMatrix::operator = (const double d)
    {
      Assert (d==0, ExcScalarAssignmentOnlyForZeroValue());

      for (size_type r=0; r<this->n_block_rows(); ++r)
        for (size_type c=0; c<this->n_block_cols(); ++c)
          this->block(r,c) = d;

      return *this;
    }



    inline
    void
    BlockSparseMatrix::vmult (BlockVector       &dst,
                              const BlockVector &src) const
    {
      BaseClass::vmult_block_block (dst, src);
    }



    inline
    void
    BlockSparseMatrix::vmult (BlockVector  &dst,
                              const Vector &src) const
    {
      BaseClass::vmult_block_nonblock (dst, src);
    }



    inline
    void
    BlockSparseMatrix::vmult (Vector            &dst,
                              const BlockVector &src) const
    {
      BaseClass::vmult_nonblock_block (dst, src);
    }



    inline
    void
    BlockSparseMatrix::vmult (Vector       &dst,
                              const Vector &src) const
    {
      BaseClass::vmult_nonblock_nonblock (dst, src);
    }


    inline
    void
    BlockSparseMatrix::Tvmult (BlockVector       &dst,
                               const BlockVector &src) const
    {
      BaseClass::Tvmult_block_block (dst, src);
    }



    inline
    void
    BlockSparseMatrix::Tvmult (BlockVector  &dst,
                               const Vector &src) const
    {
      BaseClass::Tvmult_block_nonblock (dst, src);
    }



    inline
    void
    BlockSparseMatrix::Tvmult (Vector            &dst,
                               const BlockVector &src) const
    {
      BaseClass::Tvmult_nonblock_block (dst, src);
    }



    inline
    void
    BlockSparseMatrix::Tvmult (Vector       &dst,
                               const Vector &src) const
    {
      BaseClass::Tvmult_nonblock_nonblock (dst, src);
    }

  }

}


DEAL_II_NAMESPACE_CLOSE


#endif    // DEAL_II_WITH_PETSC

#endif    // dealii_petsc_parallel_block_sparse_matrix_h