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
block_sparse_matrix.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2000 - 2022 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_block_sparse_matrix_h
17#define dealii_block_sparse_matrix_h
18
19
20#include <deal.II/base/config.h>
21
27
28#include <cmath>
29
31
32
49template <typename number>
50class BlockSparseMatrix : public BlockMatrixBase<SparseMatrix<number>>
51{
52public:
57
62
67 using pointer = typename BaseClass::pointer;
72 using iterator = typename BaseClass::iterator;
74
90 BlockSparseMatrix() = default;
91
105
109 virtual ~BlockSparseMatrix() override;
110
111
112
119
129 operator=(const double d);
130
139 void
141
156 virtual void
157 reinit(const BlockSparsityPattern &sparsity);
168 bool
169 empty() const;
170
175 get_row_length(const size_type row) const;
176
184
191 n_actually_nonzero_elements(const double threshold = 0.0) const;
192
203
208 std::size_t
220 template <typename block_number>
221 void
223 const BlockVector<block_number> &src) const;
224
229 template <typename block_number, typename nonblock_number>
230 void
232 const Vector<nonblock_number> &src) const;
233
238 template <typename block_number, typename nonblock_number>
239 void
241 const BlockVector<block_number> &src) const;
242
247 template <typename nonblock_number>
248 void
250
256 template <typename block_number>
257 void
259 const BlockVector<block_number> &src) const;
260
265 template <typename block_number, typename nonblock_number>
266 void
268 const Vector<nonblock_number> &src) const;
269
274 template <typename block_number, typename nonblock_number>
275 void
277 const BlockVector<block_number> &src) const;
278
283 template <typename nonblock_number>
284 void
286 const Vector<nonblock_number> &src) const;
300 template <class BlockVectorType>
301 void
302 precondition_Jacobi(BlockVectorType & dst,
303 const BlockVectorType &src,
304 const number omega = 1.) const;
305
311 template <typename number2>
312 void
314 const Vector<number2> &src,
315 const number omega = 1.) const;
342 void
343 print_formatted(std::ostream & out,
344 const unsigned int precision = 3,
345 const bool scientific = true,
346 const unsigned int width = 0,
347 const char * zero_string = " ",
348 const double denominator = 1.) const;
361private:
369};
370
371
372
374/* ------------------------- Template functions ---------------------- */
375
376
377
378template <typename number>
381{
383
384 for (size_type r = 0; r < this->n_block_rows(); ++r)
385 for (size_type c = 0; c < this->n_block_cols(); ++c)
386 this->block(r, c) = d;
387
388 return *this;
389}
390
391
392
393template <typename number>
394template <typename block_number>
395inline void
397 const BlockVector<block_number> &src) const
398{
399 BaseClass::vmult_block_block(dst, src);
400}
401
402
403
404template <typename number>
405template <typename block_number, typename nonblock_number>
406inline void
408 const Vector<nonblock_number> &src) const
409{
410 BaseClass::vmult_block_nonblock(dst, src);
411}
412
413
414
415template <typename number>
416template <typename block_number, typename nonblock_number>
417inline void
419 const BlockVector<block_number> &src) const
420{
421 BaseClass::vmult_nonblock_block(dst, src);
422}
423
424
425
426template <typename number>
427template <typename nonblock_number>
428inline void
430 const Vector<nonblock_number> &src) const
431{
432 BaseClass::vmult_nonblock_nonblock(dst, src);
433}
434
435
436
437template <typename number>
438template <typename block_number>
439inline void
441 const BlockVector<block_number> &src) const
442{
443 BaseClass::Tvmult_block_block(dst, src);
444}
445
446
447
448template <typename number>
449template <typename block_number, typename nonblock_number>
450inline void
452 const Vector<nonblock_number> &src) const
453{
454 BaseClass::Tvmult_block_nonblock(dst, src);
455}
456
457
458
459template <typename number>
460template <typename block_number, typename nonblock_number>
461inline void
463 const BlockVector<block_number> &src) const
464{
465 BaseClass::Tvmult_nonblock_block(dst, src);
466}
467
468
469
470template <typename number>
471template <typename nonblock_number>
472inline void
474 const Vector<nonblock_number> &src) const
475{
476 BaseClass::Tvmult_nonblock_nonblock(dst, src);
477}
478
479
480
481template <typename number>
482template <class BlockVectorType>
483inline void
485 const BlockVectorType &src,
486 const number omega) const
487{
488 Assert(this->n_block_rows() == this->n_block_cols(), ExcNotQuadratic());
489 Assert(dst.n_blocks() == this->n_block_rows(),
490 ExcDimensionMismatch(dst.n_blocks(), this->n_block_rows()));
491 Assert(src.n_blocks() == this->n_block_cols(),
492 ExcDimensionMismatch(src.n_blocks(), this->n_block_cols()));
493
494 // do a diagonal preconditioning. uses only
495 // the diagonal blocks of the matrix
496 for (size_type i = 0; i < this->n_block_rows(); ++i)
497 this->block(i, i).precondition_Jacobi(dst.block(i), src.block(i), omega);
498}
499
500
501
502template <typename number>
503template <typename number2>
504inline void
506 const Vector<number2> &src,
507 const number omega) const
508{
509 // check number of blocks. the sizes of the
510 // single block is checked in the function
511 // we call
512 Assert(this->n_block_cols() == 1,
513 ExcMessage("This function only works if the matrix has "
514 "a single block"));
515 Assert(this->n_block_rows() == 1,
516 ExcMessage("This function only works if the matrix has "
517 "a single block"));
518
519 // do a diagonal preconditioning. uses only
520 // the diagonal blocks of the matrix
521 this->block(0, 0).precondition_Jacobi(dst, src, omega);
522}
523
524
526
527#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 &)
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
SmartPointer< 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 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
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:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define DeclException0(Exception0)
Definition exceptions.h:465
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)