deal.II version GIT relicensing-2233-g9bd70d8389 2024-12-08 16: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\}}\)
Loading...
Searching...
No Matches
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) 2000 - 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_block_sparse_matrix_h
16#define dealii_block_sparse_matrix_h
17
18
19#include <deal.II/base/config.h>
20
26
27#include <cmath>
28
30
31
48template <typename number>
49class BlockSparseMatrix : public BlockMatrixBase<SparseMatrix<number>>
50{
51public:
56
61
66 using pointer = typename BaseClass::pointer;
71 using iterator = typename BaseClass::iterator;
73
89 BlockSparseMatrix() = default;
90
104
108 virtual ~BlockSparseMatrix() override;
109
110
111
118
128 operator=(const double d);
129
138 void
140
155 virtual void
156 reinit(const BlockSparsityPattern &sparsity);
167 bool
168 empty() const;
169
174 get_row_length(const size_type row) const;
175
183
190 n_actually_nonzero_elements(const double threshold = 0.0) const;
191
202
207 std::size_t
219 template <typename block_number>
220 void
222 const BlockVector<block_number> &src) const;
223
228 template <typename block_number, typename nonblock_number>
229 void
231 const Vector<nonblock_number> &src) const;
232
237 template <typename block_number, typename nonblock_number>
238 void
240 const BlockVector<block_number> &src) const;
241
246 template <typename nonblock_number>
247 void
249
255 template <typename block_number>
256 void
258 const BlockVector<block_number> &src) const;
259
264 template <typename block_number, typename nonblock_number>
265 void
267 const Vector<nonblock_number> &src) const;
268
273 template <typename block_number, typename nonblock_number>
274 void
276 const BlockVector<block_number> &src) const;
277
282 template <typename nonblock_number>
283 void
285 const Vector<nonblock_number> &src) const;
299 template <typename BlockVectorType>
300 void
301 precondition_Jacobi(BlockVectorType &dst,
302 const BlockVectorType &src,
303 const number omega = 1.) const;
304
310 template <typename number2>
311 void
313 const Vector<number2> &src,
314 const number omega = 1.) const;
343 void
344 print_formatted(std::ostream &out,
345 const unsigned int precision = 3,
346 const bool scientific = true,
347 const unsigned int width = 0,
348 const char *zero_string = " ",
349 const double denominator = 1.,
350 const char *separator = " ") const;
363private:
371};
372
373
374
376/* ------------------------- Template functions ---------------------- */
377
378
379
380template <typename number>
383{
385
386 for (size_type r = 0; r < this->n_block_rows(); ++r)
387 for (size_type c = 0; c < this->n_block_cols(); ++c)
388 this->block(r, c) = d;
389
390 return *this;
391}
392
393
394
395template <typename number>
396template <typename block_number>
397inline void
399 const BlockVector<block_number> &src) const
400{
401 BaseClass::vmult_block_block(dst, src);
402}
403
404
405
406template <typename number>
407template <typename block_number, typename nonblock_number>
408inline void
410 const Vector<nonblock_number> &src) const
411{
412 BaseClass::vmult_block_nonblock(dst, src);
413}
414
415
416
417template <typename number>
418template <typename block_number, typename nonblock_number>
419inline void
421 const BlockVector<block_number> &src) const
422{
423 BaseClass::vmult_nonblock_block(dst, src);
424}
425
426
427
428template <typename number>
429template <typename nonblock_number>
430inline void
432 const Vector<nonblock_number> &src) const
433{
434 BaseClass::vmult_nonblock_nonblock(dst, src);
435}
436
437
438
439template <typename number>
440template <typename block_number>
441inline void
443 const BlockVector<block_number> &src) const
444{
445 BaseClass::Tvmult_block_block(dst, src);
446}
447
448
449
450template <typename number>
451template <typename block_number, typename nonblock_number>
452inline void
454 const Vector<nonblock_number> &src) const
455{
456 BaseClass::Tvmult_block_nonblock(dst, src);
457}
458
459
460
461template <typename number>
462template <typename block_number, typename nonblock_number>
463inline void
465 const BlockVector<block_number> &src) const
466{
467 BaseClass::Tvmult_nonblock_block(dst, src);
468}
469
470
471
472template <typename number>
473template <typename nonblock_number>
474inline void
476 const Vector<nonblock_number> &src) const
477{
478 BaseClass::Tvmult_nonblock_nonblock(dst, src);
479}
480
481
482
483template <typename number>
484template <typename BlockVectorType>
485inline void
487 const BlockVectorType &src,
488 const number omega) const
489{
490 Assert(this->n_block_rows() == this->n_block_cols(), ExcNotQuadratic());
491 Assert(dst.n_blocks() == this->n_block_rows(),
492 ExcDimensionMismatch(dst.n_blocks(), this->n_block_rows()));
493 Assert(src.n_blocks() == this->n_block_cols(),
494 ExcDimensionMismatch(src.n_blocks(), this->n_block_cols()));
495
496 // do a diagonal preconditioning. uses only
497 // the diagonal blocks of the matrix
498 for (size_type i = 0; i < this->n_block_rows(); ++i)
499 this->block(i, i).precondition_Jacobi(dst.block(i), src.block(i), omega);
500}
501
502
503
504template <typename number>
505template <typename number2>
506inline void
508 const Vector<number2> &src,
509 const number omega) const
510{
511 // check number of blocks. the sizes of the
512 // single block is checked in the function
513 // we call
514 Assert(this->n_block_cols() == 1,
515 ExcMessage("This function only works if the matrix has "
516 "a single block"));
517 Assert(this->n_block_rows() == 1,
518 ExcMessage("This function only works if the matrix has "
519 "a single block"));
520
521 // do a diagonal preconditioning. uses only
522 // the diagonal blocks of the matrix
523 this->block(0, 0).precondition_Jacobi(dst, src, omega);
524}
525
526
528
529#endif // dealii_block_sparse_matrix_h
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, true > > const_iterator
typename BlockType::value_type value_type
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, false > > iterator
BlockSparseMatrix(const BlockSparsityPattern &sparsity)
typename BaseClass::value_type value_type
typename BaseClass::iterator iterator
typename BaseClass::const_pointer const_pointer
BlockSparseMatrix()=default
virtual void reinit(const BlockSparsityPattern &sparsity)
size_type n_nonzero_elements() const
void precondition_Jacobi(Vector< number2 > &dst, const Vector< number2 > &src, const number omega=1.) const
std::size_t memory_consumption() const
void vmult(BlockVector< block_number > &dst, const Vector< nonblock_number > &src) const
bool empty() const
void Tvmult(Vector< nonblock_number > &dst, const Vector< nonblock_number > &src) const
virtual ~BlockSparseMatrix() override
typename BaseClass::pointer pointer
BlockSparseMatrix & operator=(const BlockSparseMatrix &)
void print_formatted(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const unsigned int width=0, const char *zero_string=" ", const double denominator=1., const char *separator=" ") const
typename BaseClass::const_reference const_reference
size_type n_actually_nonzero_elements(const double threshold=0.0) const
void Tvmult(Vector< nonblock_number > &dst, const BlockVector< block_number > &src) const
typename BaseClass::size_type size_type
void vmult(Vector< nonblock_number > &dst, const BlockVector< block_number > &src) const
void vmult(Vector< nonblock_number > &dst, const Vector< nonblock_number > &src) const
ObserverPointer< const BlockSparsityPattern, BlockSparseMatrix< number > > sparsity_pattern
typename BaseClass::const_iterator const_iterator
void Tvmult(BlockVector< block_number > &dst, const BlockVector< block_number > &src) const
void Tvmult(BlockVector< block_number > &dst, const Vector< nonblock_number > &src) const
void vmult(BlockVector< block_number > &dst, const BlockVector< block_number > &src) const
typename BaseClass::BlockType BlockType
typename BaseClass::reference reference
size_type get_row_length(const size_type row) const
void precondition_Jacobi(BlockVectorType &dst, const BlockVectorType &src, const number omega=1.) const
BlockSparseMatrix & operator=(const double d)
const BlockSparsityPattern & get_sparsity_pattern() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
#define DeclException0(Exception0)
Definition exceptions.h:466
static ::ExceptionBase & ExcScalarAssignmentOnlyForZeroValue()
#define Assert(cond, exc)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcBlockDimensionMismatch()
static ::ExceptionBase & ExcNotQuadratic()
static ::ExceptionBase & ExcMessage(std::string arg1)