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_block_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_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
29
30# include <cmath>
31
33
34
35
36namespace PETScWrappers
37{
38 namespace MPI
39 {
67 class BlockSparseMatrix : public BlockMatrixBase<SparseMatrix>
68 {
69 public:
74
79
91
104
112 explicit BlockSparseMatrix(const Mat &);
113
117 template <size_t block_rows, size_t block_columns>
118 explicit BlockSparseMatrix(
119 const std::array<std::array<Mat, block_columns>, block_rows> &);
120
124 ~BlockSparseMatrix() override;
125
132
143 operator=(const double d);
144
158 void
159 reinit(const size_type n_block_rows, const size_type n_block_columns);
160
161
171 void
172 reinit(const std::vector<IndexSet> & rows,
173 const std::vector<IndexSet> & cols,
174 const BlockDynamicSparsityPattern &bdsp,
175 const MPI_Comm com);
176
177
181 void
182 reinit(const std::vector<IndexSet> & sizes,
183 const BlockDynamicSparsityPattern &bdsp,
184 const MPI_Comm com);
185
186
192 void
193 reinit(Mat A);
194
195
200 void
201 vmult(BlockVector &dst, const BlockVector &src) const;
202
207 void
208 vmult(BlockVector &dst, const Vector &src) const;
209
214 void
215 vmult(Vector &dst, const BlockVector &src) const;
216
221 void
222 vmult(Vector &dst, const Vector &src) const;
223
229 void
230 Tvmult(BlockVector &dst, const BlockVector &src) const;
231
236 void
237 Tvmult(BlockVector &dst, const Vector &src) const;
238
243 void
244 Tvmult(Vector &dst, const BlockVector &src) const;
245
250 void
251 Tvmult(Vector &dst, const Vector &src) const;
252
260 void
262
271 void
273
278 std::vector<IndexSet>
280
286 std::vector<IndexSet>
288
294 std::uint64_t
295 n_nonzero_elements() const;
296
301 get_mpi_communicator() const;
302
308
316 operator const Mat &() const;
317
328 Mat &
329 petsc_matrix();
330
331 private:
339
343 void
345
351 void
353 };
354
355
356
359 // ------------- inline and template functions -----------------
361 : BaseClass()
362 , petsc_nest_matrix(nullptr)
363 {}
364
365
366
369 {
370 this->reinit(A);
371 }
372
373
374
375 template <size_t block_rows, size_t block_columns>
377 const std::array<std::array<Mat, block_columns>, block_rows> &arrayA)
379 {
380 this->reinit(block_rows, block_columns);
381 this->sub_objects.reinit(block_rows, block_columns);
382 for (auto r = 0; r < block_rows; ++r)
383 for (auto c = 0; c < block_columns; ++c)
384 {
385 if (arrayA[r][c])
386 this->sub_objects[r][c] = new BlockType(arrayA[r][c]);
387 else
388 this->sub_objects[r][c] = nullptr;
389 }
390 this->collect_sizes();
391 }
392
393
394
395 inline BlockSparseMatrix &
397 {
399
400 for (size_type r = 0; r < this->n_block_rows(); ++r)
401 for (size_type c = 0; c < this->n_block_cols(); ++c)
402 this->block(r, c) = d;
403
404 return *this;
405 }
406
407
408
409 inline void
411 {
413 }
414
415
416
417 inline void
419 {
421 }
422
423
424
425 inline void
427 {
429 }
430
431
432
433 inline void
434 BlockSparseMatrix::vmult(Vector &dst, const Vector &src) const
435 {
437 }
438
439
440 inline void
442 {
444 }
445
446
447
448 inline void
450 {
452 }
453
454
455
456 inline void
458 {
460 }
461
462
463
464 inline void
466 {
468 }
469
470 } // namespace MPI
471
472} // namespace PETScWrappers
473
474
476
477
478#endif // DEAL_II_WITH_PETSC
479
480#endif // dealii_petsc_block_sparse_matrix_h
void Tvmult_nonblock_nonblock(VectorType &dst, const VectorType &src) const
void vmult_nonblock_nonblock(VectorType &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
unsigned int n_block_rows() const
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, true > > const_iterator
void vmult_block_block(BlockVectorType &dst, const BlockVectorType &src) const
types::global_dof_index size_type
Table< 2, SmartPointer< BlockType, BlockMatrixBase< SparseMatrix > > > sub_objects
typename BlockType::value_type value_type
void vmult_nonblock_block(VectorType &dst, const BlockVectorType &src) const
unsigned int n_block_cols() const
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, false > > iterator
BlockType & block(const unsigned int row, const unsigned int column)
void vmult_block_nonblock(BlockVectorType &dst, const VectorType &src) const
void Tvmult_block_nonblock(BlockVectorType &dst, const VectorType &src) const
void Tvmult(BlockVector &dst, const BlockVector &src) const
void vmult(BlockVector &dst, const BlockVector &src) const
void reinit(const size_type n_block_rows, const size_type n_block_columns)
BlockSparseMatrix & operator=(const BlockSparseMatrix &)
void compress(VectorOperation::values operation)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcScalarAssignmentOnlyForZeroValue()
#define Assert(cond, exc)