petsc_parallel_sparse_matrix.h

// ---------------------------------------------------------------------
//
// Copyright (C) 2004 - 2015 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_sparse_matrix_h
#define dealii__petsc_parallel_sparse_matrix_h

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

#ifdef DEAL_II_WITH_PETSC

#  include <deal.II/lac/exceptions.h>
#  include <deal.II/lac/petsc_matrix_base.h>
#  include <deal.II/lac/petsc_parallel_vector.h>
#  include <vector>

DEAL_II_NAMESPACE_OPEN


// forward declaration
template <typename Matrix> class BlockMatrixBase;


namespace PETScWrappers
{
  namespace MPI
  {



    class SparseMatrix : public MatrixBase
    {
    public:
      typedef types::global_dof_index size_type;

      struct Traits
      {
        static const bool zero_addition_can_be_elided = false;
      };

      SparseMatrix ();

      ~SparseMatrix ();

      SparseMatrix (const MPI_Comm  &communicator,
                    const size_type  m,
                    const size_type  n,
                    const size_type  local_rows,
                    const size_type  local_columns,
                    const size_type  n_nonzero_per_row,
                    const bool       is_symmetric = false,
                    const size_type  n_offdiag_nonzero_per_row = 0);

      SparseMatrix (const MPI_Comm               &communicator,
                    const size_type               m,
                    const size_type               n,
                    const size_type               local_rows,
                    const size_type               local_columns,
                    const std::vector<size_type> &row_lengths,
                    const bool                    is_symmetric = false,
                    const std::vector<size_type> &offdiag_row_lengths = std::vector<size_type>());

      template <typename SparsityPatternType>
      SparseMatrix (const MPI_Comm               &communicator,
                    const SparsityPatternType    &sparsity_pattern,
                    const std::vector<size_type> &local_rows_per_process,
                    const std::vector<size_type> &local_columns_per_process,
                    const unsigned int            this_process,
                    const bool                    preset_nonzero_locations = true);

      SparseMatrix &operator = (const value_type d);


      void copy_from(const SparseMatrix &other);

      void reinit (const MPI_Comm     &communicator,
                   const size_type m,
                   const size_type n,
                   const size_type local_rows,
                   const size_type local_columns,
                   const size_type n_nonzero_per_row,
                   const bool      is_symmetric = false,
                   const size_type n_offdiag_nonzero_per_row = 0);

      void reinit (const MPI_Comm               &communicator,
                   const size_type               m,
                   const size_type               n,
                   const size_type               local_rows,
                   const size_type               local_columns,
                   const std::vector<size_type> &row_lengths,
                   const bool                    is_symmetric = false,
                   const std::vector<size_type> &offdiag_row_lengths = std::vector<size_type>());

      template <typename SparsityPatternType>
      void reinit (const MPI_Comm               &communicator,
                   const SparsityPatternType    &sparsity_pattern,
                   const std::vector<size_type> &local_rows_per_process,
                   const std::vector<size_type> &local_columns_per_process,
                   const unsigned int            this_process,
                   const bool                    preset_nonzero_locations = true);

      template <typename SparsityPatternType>
      void reinit (const IndexSet            &local_rows,
                   const IndexSet            &local_columns,
                   const SparsityPatternType &sparsity_pattern,
                   const MPI_Comm            &communicator);

      void reinit (const SparseMatrix &other);

      virtual const MPI_Comm &get_mpi_communicator () const;

      DeclException2 (ExcLocalRowsTooLarge,
                      int, int,
                      << "The number of local rows " << arg1
                      << " must be larger than the total number of rows " << arg2);

      PetscScalar matrix_norm_square (const Vector &v) const;

      PetscScalar matrix_scalar_product (const Vector &u,
                                         const Vector &v) const;

    private:

      MPI_Comm communicator;

      void do_reinit (const size_type m,
                      const size_type n,
                      const size_type local_rows,
                      const size_type local_columns,
                      const size_type n_nonzero_per_row,
                      const bool      is_symmetric = false,
                      const size_type n_offdiag_nonzero_per_row = 0);

      void do_reinit (const size_type               m,
                      const size_type               n,
                      const size_type               local_rows,
                      const size_type               local_columns,
                      const std::vector<size_type> &row_lengths,
                      const bool                    is_symmetric = false,
                      const std::vector<size_type> &offdiag_row_lengths = std::vector<size_type>());

      template <typename SparsityPatternType>
      void do_reinit (const SparsityPatternType    &sparsity_pattern,
                      const std::vector<size_type> &local_rows_per_process,
                      const std::vector<size_type> &local_columns_per_process,
                      const unsigned int            this_process,
                      const bool                    preset_nonzero_locations);

      template <typename SparsityPatternType>
      void do_reinit (const IndexSet            &local_rows,
                      const IndexSet            &local_columns,
                      const SparsityPatternType &sparsity_pattern);

      friend class BlockMatrixBase<SparseMatrix>;
    };



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

    inline
    const MPI_Comm &
    SparseMatrix::get_mpi_communicator () const
    {
      return communicator;
    }
  }
}

DEAL_II_NAMESPACE_CLOSE

#endif // DEAL_II_WITH_PETSC

/*----------------------------   petsc_parallel_sparse_matrix.h     ---------------------------*/

#endif
/*----------------------------   petsc_parallel_sparse_matrix.h     ---------------------------*/