Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
petsc_sparse_matrix.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2019 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_petsc_sparse_matrix_h
17 # define dealii_petsc_sparse_matrix_h
18 
19 
20 # include <deal.II/base/config.h>
21 
22 # ifdef DEAL_II_WITH_PETSC
23 
24 # include <deal.II/lac/exceptions.h>
25 # include <deal.II/lac/petsc_matrix_base.h>
26 # include <deal.II/lac/petsc_vector.h>
27 
28 # include <vector>
29 
30 DEAL_II_NAMESPACE_OPEN
31 // forward declaration
32 template <typename MatrixType>
33 class BlockMatrixBase;
34 
35 
36 namespace PETScWrappers
37 {
51  class SparseMatrix : public MatrixBase
52  {
53  public:
61  struct Traits
62  {
67  static const bool zero_addition_can_be_elided = true;
68  };
69 
73  SparseMatrix();
74 
89  SparseMatrix(const size_type m,
90  const size_type n,
91  const size_type n_nonzero_per_row,
92  const bool is_symmetric = false);
93 
110  SparseMatrix(const size_type m,
111  const size_type n,
112  const std::vector<size_type> &row_lengths,
113  const bool is_symmetric = false);
114 
132  template <typename SparsityPatternType>
133  explicit SparseMatrix(const SparsityPatternType &sparsity_pattern,
134  const bool preset_nonzero_locations = true);
135 
145  SparseMatrix &
146  operator=(const double d);
147 
151  SparseMatrix(const SparseMatrix &) = delete;
152 
156  SparseMatrix &
157  operator=(const SparseMatrix &) = delete;
158 
164  void
165  reinit(const size_type m,
166  const size_type n,
167  const size_type n_nonzero_per_row,
168  const bool is_symmetric = false);
169 
175  void
176  reinit(const size_type m,
177  const size_type n,
178  const std::vector<size_type> &row_lengths,
179  const bool is_symmetric = false);
180 
206  template <typename SparsityPatternType>
207  void
208  reinit(const SparsityPatternType &sparsity_pattern,
209  const bool preset_nonzero_locations = true);
210 
216  virtual const MPI_Comm &
217  get_mpi_communicator() const override;
218 
222  size_t
223  m() const;
224 
228  size_t
229  n() const;
230 
237  void
238  mmult(SparseMatrix & C,
239  const SparseMatrix &B,
240  const MPI::Vector & V = MPI::Vector()) const;
241 
249  void
250  Tmmult(SparseMatrix & C,
251  const SparseMatrix &B,
252  const MPI::Vector & V = MPI::Vector()) const;
253 
254  private:
260  void
261  do_reinit(const size_type m,
262  const size_type n,
263  const size_type n_nonzero_per_row,
264  const bool is_symmetric = false);
265 
269  void
270  do_reinit(const size_type m,
271  const size_type n,
272  const std::vector<size_type> &row_lengths,
273  const bool is_symmetric = false);
274 
278  template <typename SparsityPatternType>
279  void
280  do_reinit(const SparsityPatternType &sparsity_pattern,
281  const bool preset_nonzero_locations);
282 
287  };
288 
289  namespace MPI
290  {
368  class SparseMatrix : public MatrixBase
369  {
370  public:
375 
383  struct Traits
384  {
392  static const bool zero_addition_can_be_elided = false;
393  };
394 
398  SparseMatrix();
399 
403  ~SparseMatrix() override;
404 
427  DEAL_II_DEPRECATED
428  SparseMatrix(const MPI_Comm &communicator,
429  const size_type m,
430  const size_type n,
431  const size_type local_rows,
432  const size_type local_columns,
433  const size_type n_nonzero_per_row,
434  const bool is_symmetric = false,
435  const size_type n_offdiag_nonzero_per_row = 0);
436 
460  DEAL_II_DEPRECATED
461  SparseMatrix(const MPI_Comm & communicator,
462  const size_type m,
463  const size_type n,
464  const size_type local_rows,
465  const size_type local_columns,
466  const std::vector<size_type> &row_lengths,
467  const bool is_symmetric = false,
468  const std::vector<size_type> &offdiag_row_lengths =
469  std::vector<size_type>());
470 
493  template <typename SparsityPatternType>
494  SparseMatrix(const MPI_Comm & communicator,
495  const SparsityPatternType & sparsity_pattern,
496  const std::vector<size_type> &local_rows_per_process,
497  const std::vector<size_type> &local_columns_per_process,
498  const unsigned int this_process,
499  const bool preset_nonzero_locations = true);
500 
510  SparseMatrix &
511  operator=(const value_type d);
512 
513 
518  void
519  copy_from(const SparseMatrix &other);
520 
529  DEAL_II_DEPRECATED
530  void
531  reinit(const MPI_Comm &communicator,
532  const size_type m,
533  const size_type n,
534  const size_type local_rows,
535  const size_type local_columns,
536  const size_type n_nonzero_per_row,
537  const bool is_symmetric = false,
538  const size_type n_offdiag_nonzero_per_row = 0);
539 
548  DEAL_II_DEPRECATED
549  void
550  reinit(const MPI_Comm & communicator,
551  const size_type m,
552  const size_type n,
553  const size_type local_rows,
554  const size_type local_columns,
555  const std::vector<size_type> &row_lengths,
556  const bool is_symmetric = false,
557  const std::vector<size_type> &offdiag_row_lengths =
558  std::vector<size_type>());
559 
579  template <typename SparsityPatternType>
580  void
581  reinit(const MPI_Comm & communicator,
582  const SparsityPatternType & sparsity_pattern,
583  const std::vector<size_type> &local_rows_per_process,
584  const std::vector<size_type> &local_columns_per_process,
585  const unsigned int this_process,
586  const bool preset_nonzero_locations = true);
587 
594  template <typename SparsityPatternType>
595  void
596  reinit(const IndexSet & local_rows,
597  const IndexSet & local_columns,
598  const SparsityPatternType &sparsity_pattern,
599  const MPI_Comm & communicator);
600 
606  void
607  reinit(const SparseMatrix &other);
608 
613  virtual const MPI_Comm &
614  get_mpi_communicator() const override;
615 
624  int,
625  int,
626  << "The number of local rows " << arg1
627  << " must be larger than the total number of rows "
628  << arg2);
630 
646  PetscScalar
647  matrix_norm_square(const Vector &v) const;
648 
657  PetscScalar
658  matrix_scalar_product(const Vector &u, const Vector &v) const;
659 
664  IndexSet
666 
672  IndexSet
674 
681  void
682  mmult(SparseMatrix & C,
683  const SparseMatrix &B,
684  const MPI::Vector & V = MPI::Vector()) const;
685 
693  void
694  Tmmult(SparseMatrix & C,
695  const SparseMatrix &B,
696  const MPI::Vector & V = MPI::Vector()) const;
697 
698  private:
702  MPI_Comm communicator;
703 
712  DEAL_II_DEPRECATED
713  void
714  do_reinit(const size_type m,
715  const size_type n,
716  const size_type local_rows,
717  const size_type local_columns,
718  const size_type n_nonzero_per_row,
719  const bool is_symmetric = false,
720  const size_type n_offdiag_nonzero_per_row = 0);
721 
728  DEAL_II_DEPRECATED
729  void
730  do_reinit(const size_type m,
731  const size_type n,
732  const size_type local_rows,
733  const size_type local_columns,
734  const std::vector<size_type> &row_lengths,
735  const bool is_symmetric = false,
736  const std::vector<size_type> &offdiag_row_lengths =
737  std::vector<size_type>());
738 
742  template <typename SparsityPatternType>
743  void
744  do_reinit(const SparsityPatternType & sparsity_pattern,
745  const std::vector<size_type> &local_rows_per_process,
746  const std::vector<size_type> &local_columns_per_process,
747  const unsigned int this_process,
748  const bool preset_nonzero_locations);
749 
753  template <typename SparsityPatternType>
754  void
755  do_reinit(const IndexSet & local_rows,
756  const IndexSet & local_columns,
757  const SparsityPatternType &sparsity_pattern);
758 
763  };
764 
765 
766 
767  // -------- template and inline functions ----------
768 
769  inline const MPI_Comm &
771  {
772  return communicator;
773  }
774  } // namespace MPI
775 } // namespace PETScWrappers
776 
777 DEAL_II_NAMESPACE_CLOSE
778 
779 # endif // DEAL_II_WITH_PETSC
780 
781 #endif
782 /*--------------------------- petsc_sparse_matrix.h -------------------------*/
static ::ExceptionBase & ExcLocalRowsTooLarge(int arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:541
void copy_from(const SparseMatrix &other)
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)
virtual const MPI_Comm & get_mpi_communicator() const override
PetscBool is_symmetric(const double tolerance=1.e-12)
void mmult(SparseMatrix &C, const SparseMatrix &B, const MPI::Vector &V=MPI::Vector()) const
SparseMatrix & operator=(const double d)
void reinit(const size_type m, const size_type n, const size_type n_nonzero_per_row, const bool is_symmetric=false)
SparseMatrix & operator=(const value_type d)
PetscScalar matrix_scalar_product(const Vector &u, const Vector &v) const
virtual const MPI_Comm & get_mpi_communicator() const override
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)
unsigned int global_dof_index
Definition: types.h:89
types::global_dof_index size_type
void Tmmult(SparseMatrix &C, const SparseMatrix &B, const MPI::Vector &V=MPI::Vector()) const
void mmult(SparseMatrix &C, const SparseMatrix &B, const MPI::Vector &V=MPI::Vector()) const
PetscScalar matrix_norm_square(const Vector &v) const
void Tmmult(SparseMatrix &C, const SparseMatrix &B, const MPI::Vector &V=MPI::Vector()) const
void do_reinit(const size_type m, const size_type n, const size_type n_nonzero_per_row, const bool is_symmetric=false)