Reference documentation for deal.II version GIT d8dacc551e 2022-08-19 06:50:03+00:00
\(\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\}}\)
petsc_matrix_free.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2012 - 2020 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_matrix_free_h
17 # define dealii_petsc_matrix_free_h
18 
19 
20 # include <deal.II/base/config.h>
21 
22 # ifdef DEAL_II_WITH_PETSC
23 # include <deal.II/lac/exceptions.h>
27 
28 
29 
30 namespace PETScWrappers
31 {
59  class MatrixFree : public MatrixBase
60  {
61  public:
65  MatrixFree();
66 
81  MatrixFree(const MPI_Comm & communicator,
82  const unsigned int m,
83  const unsigned int n,
84  const unsigned int local_rows,
85  const unsigned int local_columns);
86 
98  MatrixFree(const MPI_Comm & communicator,
99  const unsigned int m,
100  const unsigned int n,
101  const std::vector<unsigned int> &local_rows_per_process,
102  const std::vector<unsigned int> &local_columns_per_process,
103  const unsigned int this_process);
104 
110  MatrixFree(const unsigned int m,
111  const unsigned int n,
112  const unsigned int local_rows,
113  const unsigned int local_columns);
114 
120  MatrixFree(const unsigned int m,
121  const unsigned int n,
122  const std::vector<unsigned int> &local_rows_per_process,
123  const std::vector<unsigned int> &local_columns_per_process,
124  const unsigned int this_process);
125 
131  void
132  reinit(const MPI_Comm & communicator,
133  const unsigned int m,
134  const unsigned int n,
135  const unsigned int local_rows,
136  const unsigned int local_columns);
137 
143  void
144  reinit(const MPI_Comm & communicator,
145  const unsigned int m,
146  const unsigned int n,
147  const std::vector<unsigned int> &local_rows_per_process,
148  const std::vector<unsigned int> &local_columns_per_process,
149  const unsigned int this_process);
150 
155  void
156  reinit(const unsigned int m,
157  const unsigned int n,
158  const unsigned int local_rows,
159  const unsigned int local_columns);
160 
165  void
166  reinit(const unsigned int m,
167  const unsigned int n,
168  const std::vector<unsigned int> &local_rows_per_process,
169  const std::vector<unsigned int> &local_columns_per_process,
170  const unsigned int this_process);
171 
176  void
177  clear();
178 
183  const MPI_Comm &
184  get_mpi_communicator() const override;
185 
197  virtual void
198  vmult(VectorBase &dst, const VectorBase &src) const = 0;
199 
212  virtual void
213  Tvmult(VectorBase &dst, const VectorBase &src) const = 0;
214 
226  virtual void
227  vmult_add(VectorBase &dst, const VectorBase &src) const = 0;
228 
241  virtual void
242  Tvmult_add(VectorBase &dst, const VectorBase &src) const = 0;
243 
252  virtual void
253  vmult(Vec &dst, const Vec &src) const;
254 
255  private:
260  MPI_Comm communicator;
261 
273  static int
274  matrix_free_mult(Mat A, Vec src, Vec dst);
275 
281  void
282  do_reinit(const unsigned int m,
283  const unsigned int n,
284  const unsigned int local_rows,
285  const unsigned int local_columns);
286  };
287 
288 
289 
290  // -------- template and inline functions ----------
291 
292  inline const MPI_Comm &
294  {
295  return communicator;
296  }
297 } // namespace PETScWrappers
298 
300 
301 # endif // DEAL_II_WITH_PETSC
302 
303 #endif
304 /*---------------------------- petsc_matrix_free.h --------------------------*/
const MPI_Comm & get_mpi_communicator() const override
virtual void Tvmult(VectorBase &dst, const VectorBase &src) const =0
void do_reinit(const unsigned int m, const unsigned int n, const unsigned int local_rows, const unsigned int local_columns)
void reinit(const MPI_Comm &communicator, const unsigned int m, const unsigned int n, const unsigned int local_rows, const unsigned int local_columns)
virtual void vmult_add(VectorBase &dst, const VectorBase &src) const =0
static int matrix_free_mult(Mat A, Vec src, Vec dst)
virtual void Tvmult_add(VectorBase &dst, const VectorBase &src) const =0
virtual void vmult(VectorBase &dst, const VectorBase &src) const =0
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
static const char A