Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
petsc_block_sparse_matrix.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2018 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_block_sparse_matrix_h
17 #define dealii_petsc_block_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/block_matrix_base.h>
25 # include <deal.II/lac/block_sparsity_pattern.h>
26 # include <deal.II/lac/exceptions.h>
27 # include <deal.II/lac/petsc_block_vector.h>
28 # include <deal.II/lac/petsc_sparse_matrix.h>
29 
30 # include <cmath>
31 
32 DEAL_II_NAMESPACE_OPEN
33 
34 
35 
36 namespace PETScWrappers
37 {
38  namespace MPI
39  {
67  class BlockSparseMatrix : public BlockMatrixBase<SparseMatrix>
68  {
69  public:
74 
79 
84  using pointer = BaseClass::pointer;
85  using const_pointer = BaseClass::const_pointer;
86  using reference = BaseClass::reference;
87  using const_reference = BaseClass::const_reference;
88  using size_type = BaseClass::size_type;
91 
103  BlockSparseMatrix() = default;
104 
108  ~BlockSparseMatrix() override = default;
109 
115  operator=(const BlockSparseMatrix &);
116 
127  operator=(const double d);
128 
142  void
143  reinit(const size_type n_block_rows, const size_type n_block_columns);
144 
145 
155  void
156  reinit(const std::vector<IndexSet> & rows,
157  const std::vector<IndexSet> & cols,
158  const BlockDynamicSparsityPattern &bdsp,
159  const MPI_Comm & com);
160 
161 
165  void
166  reinit(const std::vector<IndexSet> & sizes,
167  const BlockDynamicSparsityPattern &bdsp,
168  const MPI_Comm & com);
169 
170 
171 
176  void
177  vmult(BlockVector &dst, const BlockVector &src) const;
178 
183  void
184  vmult(BlockVector &dst, const Vector &src) const;
185 
190  void
191  vmult(Vector &dst, const BlockVector &src) const;
192 
197  void
198  vmult(Vector &dst, const Vector &src) const;
199 
205  void
206  Tvmult(BlockVector &dst, const BlockVector &src) const;
207 
212  void
213  Tvmult(BlockVector &dst, const Vector &src) const;
214 
219  void
220  Tvmult(Vector &dst, const BlockVector &src) const;
221 
226  void
227  Tvmult(Vector &dst, const Vector &src) const;
228 
236  void
237  collect_sizes();
238 
243  std::vector<IndexSet>
245 
251  std::vector<IndexSet>
253 
258  const MPI_Comm &
259  get_mpi_communicator() const;
260 
266  };
267 
268 
269 
272  // ------------- inline and template functions -----------------
273 
274  inline BlockSparseMatrix &
276  {
278 
279  for (size_type r = 0; r < this->n_block_rows(); ++r)
280  for (size_type c = 0; c < this->n_block_cols(); ++c)
281  this->block(r, c) = d;
282 
283  return *this;
284  }
285 
286 
287 
288  inline void
290  {
292  }
293 
294 
295 
296  inline void
298  {
300  }
301 
302 
303 
304  inline void
306  {
308  }
309 
310 
311 
312  inline void
313  BlockSparseMatrix::vmult(Vector &dst, const Vector &src) const
314  {
316  }
317 
318 
319  inline void
321  {
323  }
324 
325 
326 
327  inline void
329  {
331  }
332 
333 
334 
335  inline void
337  {
339  }
340 
341 
342 
343  inline void
344  BlockSparseMatrix::Tvmult(Vector &dst, const Vector &src) const
345  {
347  }
348 
349  } // namespace MPI
350 
351 } // namespace PETScWrappers
352 
353 
354 DEAL_II_NAMESPACE_CLOSE
355 
356 
357 #endif // DEAL_II_WITH_PETSC
358 
359 #endif // dealii_petsc_block_sparse_matrix_h
void vmult_block_nonblock(BlockVectorType &dst, const VectorType &src) const
std::vector< IndexSet > locally_owned_domain_indices() const
static ::ExceptionBase & ExcScalarAssignmentOnlyForZeroValue()
void vmult_nonblock_block(VectorType &dst, const BlockVectorType &src) const
void Tvmult_block_nonblock(BlockVectorType &dst, const VectorType &src) const
void Tvmult_block_block(BlockVectorType &dst, const BlockVectorType &src) const
void Tvmult_nonblock_block(VectorType &dst, const BlockVectorType &src) const
void vmult_block_block(BlockVectorType &dst, const BlockVectorType &src) const
typename BlockType::value_type value_type
void Tvmult_nonblock_nonblock(VectorType &dst, const VectorType &src) const
unsigned int n_block_cols() const
void Tvmult(BlockVector &dst, const BlockVector &src) const
#define Assert(cond, exc)
Definition: exceptions.h:1407
BlockSparseMatrix & operator=(const BlockSparseMatrix &)
void reinit(const size_type n_block_rows, const size_type n_block_columns)
void vmult_nonblock_nonblock(VectorType &dst, const VectorType &src) const
BlockType & block(const unsigned int row, const unsigned int column)
void vmult(BlockVector &dst, const BlockVector &src) const
unsigned int n_block_rows() const