Loading [MathJax]/extensions/TeX/AMSsymbols.js
 deal.II version GIT relicensing-3075-gc235bd4825 2025-04-15 08:10:00+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\}}\)
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
petsc_block_sparse_matrix.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2004 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_petsc_block_sparse_matrix_h
16#define dealii_petsc_block_sparse_matrix_h
17
18
19#include <deal.II/base/config.h>
20
21#ifdef DEAL_II_WITH_PETSC
22
28
29# include <cmath>
30# include <cstddef>
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 <std::size_t block_rows, std::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
259 void
261
270 void
272
277 std::vector<IndexSet>
279
285 std::vector<IndexSet>
287
293 std::uint64_t
294 n_nonzero_elements() const;
295
300 get_mpi_communicator() const;
301
307
315 operator const Mat &() const;
316
327 Mat &
328 petsc_matrix();
329
330 private:
338
342 void
344
350 void
352 };
353
354
355
358 // ------------- inline and template functions -----------------
360 : BaseClass()
361 , petsc_nest_matrix(nullptr)
362 {}
363
364
365
368 {
369 this->reinit(A);
370 }
371
372
373
374 template <std::size_t block_rows, std::size_t block_columns>
376 const std::array<std::array<Mat, block_columns>, block_rows> &arrayA)
378 {
379 this->reinit(block_rows, block_columns);
380 this->sub_objects.reinit(block_rows, block_columns);
381 for (auto r = 0; r < block_rows; ++r)
382 for (auto c = 0; c < block_columns; ++c)
383 {
384 if (arrayA[r][c])
385 this->sub_objects[r][c] = new BlockType(arrayA[r][c]);
386 else
387 this->sub_objects[r][c] = nullptr;
388 }
389 this->collect_sizes();
390 }
391
392
393
394 inline BlockSparseMatrix &
396 {
398
399 for (size_type r = 0; r < this->n_block_rows(); ++r)
400 for (size_type c = 0; c < this->n_block_cols(); ++c)
401 this->block(r, c) = d;
402
403 return *this;
404 }
405
406
407
408 inline void
410 {
412 }
413
414
415
416 inline void
418 {
420 }
421
422
423
424 inline void
426 {
428 }
429
430
431
432 inline void
433 BlockSparseMatrix::vmult(Vector &dst, const Vector &src) const
434 {
436 }
437
438
439 inline void
441 {
443 }
444
445
446
447 inline void
449 {
451 }
452
453
454
455 inline void
457 {
459 }
460
461
462
463 inline void
465 {
467 }
468
469 } // namespace MPI
470
471} // namespace PETScWrappers
472
473
475
476
477#endif // DEAL_II_WITH_PETSC
478
479#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
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
Table< 2, ObserverPointer< BlockType, BlockMatrixBase< SparseMatrix > > > sub_objects
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:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
static ::ExceptionBase & ExcScalarAssignmentOnlyForZeroValue()
#define Assert(cond, exc)