Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.4.1
\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}} \newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=} \newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]} \newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
petsc_sparse_matrix.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2004 - 2021 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
27
28# include <vector>
29
31// forward declaration
32# ifndef DOXYGEN
33template <typename MatrixType>
34class BlockMatrixBase;
35# endif
36
37namespace PETScWrappers
38{
51 class SparseMatrix : public MatrixBase
52 {
53 public:
61 struct Traits
62 {
67 static const bool zero_addition_can_be_elided = true;
68 };
69
74
90 const size_type n,
91 const size_type n_nonzero_per_row,
92 const bool is_symmetric = false);
93
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
146 operator=(const double d);
147
151 SparseMatrix(const SparseMatrix &) = delete;
152
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
239 const SparseMatrix &B,
240 const MPI::Vector & V = MPI::Vector()) const;
241
249 void
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
283 // To allow calling protected prepare_add() and prepare_set().
284 friend class BlockMatrixBase<SparseMatrix>;
285 };
286
287 namespace MPI
288 {
366 {
367 public:
372
380 struct Traits
381 {
389 static const bool zero_addition_can_be_elided = false;
390 };
391
395 SparseMatrix();
396
400 ~SparseMatrix() override;
401
424 template <typename SparsityPatternType>
426 const SparsityPatternType & sparsity_pattern,
427 const std::vector<size_type> &local_rows_per_process,
428 const std::vector<size_type> &local_columns_per_process,
429 const unsigned int this_process,
430 const bool preset_nonzero_locations = true);
431
442 operator=(const value_type d);
443
444
449 void
450 copy_from(const SparseMatrix &other);
451
471 template <typename SparsityPatternType>
472 void
474 const SparsityPatternType & sparsity_pattern,
475 const std::vector<size_type> &local_rows_per_process,
476 const std::vector<size_type> &local_columns_per_process,
477 const unsigned int this_process,
478 const bool preset_nonzero_locations = true);
479
486 template <typename SparsityPatternType>
487 void
488 reinit(const IndexSet & local_rows,
489 const IndexSet & local_columns,
490 const SparsityPatternType &sparsity_pattern,
491 const MPI_Comm & communicator);
492
498 void
499 reinit(const SparseMatrix &other);
500
505 virtual const MPI_Comm &
506 get_mpi_communicator() const override;
507
516 int,
517 int,
518 << "The number of local rows " << arg1
519 << " must be larger than the total number of rows "
520 << arg2);
522
538 PetscScalar
539 matrix_norm_square(const Vector &v) const;
540
549 PetscScalar
550 matrix_scalar_product(const Vector &u, const Vector &v) const;
551
558
566
573 void
575 const SparseMatrix &B,
576 const MPI::Vector & V = MPI::Vector()) const;
577
585 void
587 const SparseMatrix &B,
588 const MPI::Vector & V = MPI::Vector()) const;
589
590 private:
595
599 template <typename SparsityPatternType>
600 void
601 do_reinit(const SparsityPatternType & sparsity_pattern,
602 const std::vector<size_type> &local_rows_per_process,
603 const std::vector<size_type> &local_columns_per_process,
604 const unsigned int this_process,
605 const bool preset_nonzero_locations);
606
610 template <typename SparsityPatternType>
611 void
612 do_reinit(const IndexSet & local_rows,
613 const IndexSet & local_columns,
614 const SparsityPatternType &sparsity_pattern);
615
616 // To allow calling protected prepare_add() and prepare_set().
617 friend class BlockMatrixBase<SparseMatrix>;
618 };
619
620
621
622 // -------- template and inline functions ----------
623
624 inline const MPI_Comm &
626 {
627 return communicator;
628 }
629 } // namespace MPI
630} // namespace PETScWrappers
631
633
634# endif // DEAL_II_WITH_PETSC
635
636#endif
637/*--------------------------- petsc_sparse_matrix.h -------------------------*/
SparseMatrix & operator=(const value_type d)
void copy_from(const SparseMatrix &other)
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)
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
void do_reinit(const size_type m, const size_type n, const size_type n_nonzero_per_row, const bool is_symmetric=false)
virtual const MPI_Comm & get_mpi_communicator() const override
void Tmmult(SparseMatrix &C, const SparseMatrix &B, const MPI::Vector &V=MPI::Vector()) const
void reinit(const size_type m, const size_type n, const size_type n_nonzero_per_row, const bool is_symmetric=false)
SparseMatrix(const SparseMatrix &)=delete
SparseMatrix & operator=(const SparseMatrix &)=delete
SparseMatrix & operator=(const double d)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
void mmult(SparseMatrix &C, const SparseMatrix &B, const MPI::Vector &V=MPI::Vector()) const
PetscScalar matrix_scalar_product(const Vector &u, const Vector &v) const
void Tmmult(SparseMatrix &C, const SparseMatrix &B, const MPI::Vector &V=MPI::Vector()) const
static ::ExceptionBase & ExcLocalRowsTooLarge(int arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:532
PetscScalar matrix_norm_square(const Vector &v) const
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)
unsigned int global_dof_index
Definition: types.h:76