17 #ifndef __deal2__petsc_matrix_free_h
18 #define __deal2__petsc_matrix_free_h
21 #include <deal.II/base/config.h>
23 #ifdef DEAL_II_WITH_PETSC
25 # include <deal.II/lac/exceptions.h>
26 # include <deal.II/lac/petsc_matrix_base.h>
27 # include <deal.II/lac/petsc_vector.h>
29 DEAL_II_NAMESPACE_OPEN
33 namespace PETScWrappers
99 const unsigned int local_rows,
100 const unsigned int local_columns);
120 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);
133 const unsigned int n,
134 const unsigned int local_rows,
135 const unsigned int local_columns);
144 const unsigned int n,
145 const std::vector<unsigned int> &local_rows_per_process,
146 const std::vector<unsigned int> &local_columns_per_process,
147 const unsigned int this_process);
157 void reinit (
const MPI_Comm &communicator,
158 const unsigned int m,
159 const unsigned int n,
160 const unsigned int local_rows,
161 const unsigned int local_columns);
171 void reinit (
const MPI_Comm &communicator,
172 const unsigned int m,
173 const unsigned int n,
174 const std::vector<unsigned int> &local_rows_per_process,
175 const std::vector<unsigned int> &local_columns_per_process,
176 const unsigned int this_process);
182 void reinit (
const unsigned int m,
183 const unsigned int n,
184 const unsigned int local_rows,
185 const unsigned int local_columns);
191 void reinit (
const unsigned int m,
192 const unsigned int n,
193 const std::vector<unsigned int> &local_rows_per_process,
194 const std::vector<unsigned int> &local_columns_per_process,
195 const unsigned int this_process);
315 void vmult (Vec &dst,
const Vec &src)
const;
354 const unsigned int n,
355 const unsigned int local_rows,
356 const unsigned int local_columns);
371 DEAL_II_NAMESPACE_CLOSE
373 #endif // DEAL_II_WITH_PETSC
virtual void vmult_add(VectorBase &dst, const VectorBase &src) const =0
virtual void Tvmult_add(VectorBase &dst, const VectorBase &src) const =0
virtual void vmult(VectorBase &dst, const VectorBase &src) const =0
virtual void Tvmult(VectorBase &dst, const VectorBase &src) const =0
const MPI_Comm & get_mpi_communicator() const
static int matrix_free_mult(Mat A, Vec src, Vec dst)
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)