Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
petsc_sparse_matrix.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2004 - 2023 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 <petscis.h>
29# include <petscistypes.h>
30
31# include <vector>
32
34// forward declaration
35# ifndef DOXYGEN
36template <typename MatrixType>
37class BlockMatrixBase;
38# endif
39
40namespace PETScWrappers
41{
54 class SparseMatrix : public MatrixBase
55 {
56 public:
64 struct Traits
65 {
70 static const bool zero_addition_can_be_elided = true;
71 };
72
77
84 explicit SparseMatrix(const Mat &);
85
101 const size_type n,
102 const size_type n_nonzero_per_row,
103 const bool is_symmetric = false);
104
122 const size_type n,
123 const std::vector<size_type> &row_lengths,
124 const bool is_symmetric = false);
125
143 template <typename SparsityPatternType>
144 explicit SparseMatrix(const SparsityPatternType &sparsity_pattern,
145 const bool preset_nonzero_locations = true);
146
157 operator=(const double d);
158
162 SparseMatrix(const SparseMatrix &) = delete;
163
168 operator=(const SparseMatrix &) = delete;
169
175 void
176 reinit(const size_type m,
177 const size_type n,
178 const size_type n_nonzero_per_row,
179 const bool is_symmetric = false);
180
186 void
187 reinit(const size_type m,
188 const size_type n,
189 const std::vector<size_type> &row_lengths,
190 const bool is_symmetric = false);
191
217 template <typename SparsityPatternType>
218 void
219 reinit(const SparsityPatternType &sparsity_pattern,
220 const bool preset_nonzero_locations = true);
221
225 size_t
226 m() const;
227
231 size_t
232 n() const;
233
240 void
242 const SparseMatrix &B,
243 const MPI::Vector & V = MPI::Vector()) const;
244
252 void
254 const SparseMatrix &B,
255 const MPI::Vector & V = MPI::Vector()) const;
256
257 private:
262 void
263 do_reinit(const size_type m,
264 const size_type n,
265 const size_type n_nonzero_per_row,
266 const bool is_symmetric = false);
267
271 void
272 do_reinit(const size_type m,
273 const size_type n,
274 const std::vector<size_type> &row_lengths,
275 const bool is_symmetric = false);
276
280 template <typename SparsityPatternType>
281 void
282 do_reinit(const SparsityPatternType &sparsity_pattern,
283 const bool preset_nonzero_locations);
284
285 // To allow calling protected prepare_add() and prepare_set().
286 friend class BlockMatrixBase<SparseMatrix>;
287 };
288
289 namespace MPI
290 {
368 {
369 public:
374
382 struct Traits
383 {
391 static const bool zero_addition_can_be_elided = false;
392 };
393
397 SparseMatrix();
398
405 explicit SparseMatrix(const Mat &);
406
410 ~SparseMatrix() override;
411
434 template <typename SparsityPatternType>
435 SparseMatrix(const MPI_Comm communicator,
436 const SparsityPatternType & sparsity_pattern,
437 const std::vector<size_type> &local_rows_per_process,
438 const std::vector<size_type> &local_columns_per_process,
439 const unsigned int this_process,
440 const bool preset_nonzero_locations = true);
441
452 operator=(const value_type d);
453
454
459 void
460 copy_from(const SparseMatrix &other);
461
481 template <typename SparsityPatternType>
482 void
483 reinit(const MPI_Comm communicator,
484 const SparsityPatternType & sparsity_pattern,
485 const std::vector<size_type> &local_rows_per_process,
486 const std::vector<size_type> &local_columns_per_process,
487 const unsigned int this_process,
488 const bool preset_nonzero_locations = true);
489
496 template <typename SparsityPatternType>
497 void
498 reinit(const IndexSet & local_partitioning,
499 const SparsityPatternType &sparsity_pattern,
500 const MPI_Comm communicator);
501
508 template <typename SparsityPatternType>
509 void
510 reinit(const IndexSet & local_rows,
511 const IndexSet & local_columns,
512 const SparsityPatternType &sparsity_pattern,
513 const MPI_Comm communicator);
514
520 void
521 reinit(const SparseMatrix &other);
522
532 template <typename SparsityPatternType>
533 void
534 reinit(const IndexSet & local_rows,
535 const IndexSet & local_active_rows,
536 const IndexSet & local_columns,
537 const IndexSet & local_active_columns,
538 const SparsityPatternType &sparsity_pattern,
539 const MPI_Comm communicator);
540
549 int,
550 int,
551 << "The number of local rows " << arg1
552 << " must be larger than the total number of rows "
553 << arg2);
571 PetscScalar
572 matrix_norm_square(const Vector &v) const;
573
582 PetscScalar
583 matrix_scalar_product(const Vector &u, const Vector &v) const;
584
591
599
606 void
608 const SparseMatrix &B,
609 const MPI::Vector & V = MPI::Vector()) const;
610
618 void
620 const SparseMatrix &B,
621 const MPI::Vector & V = MPI::Vector()) const;
622
623 private:
627 template <typename SparsityPatternType>
628 void
629 do_reinit(const MPI_Comm comm,
630 const SparsityPatternType & sparsity_pattern,
631 const std::vector<size_type> &local_rows_per_process,
632 const std::vector<size_type> &local_columns_per_process,
633 const unsigned int this_process,
634 const bool preset_nonzero_locations);
635
639 template <typename SparsityPatternType>
640 void
641 do_reinit(const MPI_Comm comm,
642 const IndexSet & local_rows,
643 const IndexSet & local_columns,
644 const SparsityPatternType &sparsity_pattern);
645
650 template <typename SparsityPatternType>
651 void
653 const IndexSet & local_rows,
654 const IndexSet & local_active_rows,
655 const IndexSet & local_columns,
656 const IndexSet & local_active_columns,
657 const SparsityPatternType &sparsity_pattern);
658
659 // To allow calling protected prepare_add() and prepare_set().
660 friend class BlockMatrixBase<SparseMatrix>;
661 };
662
663 } // namespace MPI
664} // namespace PETScWrappers
665
667
668#endif // DEAL_II_WITH_PETSC
669
670#endif
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)
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
PetscScalar matrix_norm_square(const Vector &v) const
void do_reinit(const MPI_Comm comm, 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)
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)
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:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcLocalRowsTooLarge(int arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:533
unsigned int global_dof_index
Definition types.h:82
const MPI_Comm comm