![]() |
Reference documentation for deal.II version 8.1.0
|
#include <block_matrix_array.h>
Public Types | |
typedef types::global_dof_index | size_type |
Public Member Functions | |
BlockTrianglePrecondition () | |
BlockTrianglePrecondition (const unsigned int n_blocks) | |
BlockTrianglePrecondition (const unsigned int n_block_rows, VectorMemory< Vector< number > > &mem, const bool backward=false) DEAL_II_DEPRECATED | |
void | initialize (const unsigned int n_block_rows, VectorMemory< Vector< number > > &mem, const bool backward=false) DEAL_II_DEPRECATED |
void | reinit (const unsigned int n_block_rows) |
template<class MATRIX > | |
void | enter (const MATRIX &matrix, const size_type row, const size_type col, const double prefix=1., const bool transpose=false) |
template<class MATRIX > | |
void | enter_aux (VectorMemory< Vector< double > > &mem, const MATRIX &matrix, const size_type row, const size_type col, const double prefix=1., const bool transpose=false) DEAL_II_DEPRECATED |
void | vmult (BlockVector< number > &dst, const BlockVector< number > &src) const |
void | vmult_add (BlockVector< number > &dst, const BlockVector< number > &src) const |
void | Tvmult (BlockVector< number > &dst, const BlockVector< number > &src) const |
void | Tvmult_add (BlockVector< number > &dst, const BlockVector< number > &src) const |
DeclException1 (ExcNoDiagonal, size_type,<< "No diagonal entry was added for block "<< arg1) | |
DeclException1 (ExcMultipleDiagonal, size_type,<< "Inverse diagonal entries may not be added in block "<< arg1) | |
Private Member Functions | |
void | do_row (BlockVector< number > &dst, size_type row_num) const |
![]() | |
BlockMatrixArray () | |
BlockMatrixArray (const unsigned int n_block_rows, const unsigned int n_block_cols) | |
void | initialize (const unsigned int n_block_rows, const unsigned int n_block_cols) |
BlockMatrixArray (const unsigned int n_block_rows, const unsigned int n_block_cols, VectorMemory< Vector< number > > &mem) DEAL_II_DEPRECATED | |
void | initialize (const unsigned int n_block_rows, const unsigned int n_block_cols, VectorMemory< Vector< number > > &mem) DEAL_II_DEPRECATED |
void | reinit (const unsigned int n_block_rows, const unsigned int n_block_cols) |
template<class MATRIX > | |
void | enter (const MATRIX &matrix, const unsigned int row, const unsigned int col, const double prefix=1., const bool transpose=false) |
template<class MATRIX > | |
void | enter_aux (VectorMemory< Vector< number > > &mem, const MATRIX &matrix, const unsigned int row, const unsigned int col, const double prefix=1., const bool transpose=false) DEAL_II_DEPRECATED |
void | clear () |
unsigned int | n_block_rows () const |
unsigned int | n_block_cols () const |
void | vmult (BlockVector< number > &dst, const BlockVector< number > &src) const |
void | vmult_add (BlockVector< number > &dst, const BlockVector< number > &src) const |
void | Tvmult (BlockVector< number > &dst, const BlockVector< number > &src) const |
void | Tvmult_add (BlockVector< number > &dst, const BlockVector< number > &src) const |
number | matrix_scalar_product (const BlockVector< number > &u, const BlockVector< number > &v) const |
number | matrix_norm_square (const BlockVector< number > &u) const |
template<class STREAM > | |
void | print_latex (STREAM &out) const |
![]() | |
Subscriptor () | |
Subscriptor (const Subscriptor &) | |
virtual | ~Subscriptor () |
Subscriptor & | operator= (const Subscriptor &) |
void | subscribe (const char *identifier=0) const |
void | unsubscribe (const char *identifier=0) const |
unsigned int | n_subscriptions () const |
void | list_subscribers () const |
DeclException3 (ExcInUse, int, char *, std::string &,<< "Object of class "<< arg2<< " is still used by "<< arg1<< " other objects.\n"<< "(Additional information: "<< arg3<< ")\n"<< "Note the entry in the Frequently Asked Questions of "<< "deal.II (linked to from http://www.dealii.org/) for "<< "more information on what this error means.") | |
DeclException2 (ExcNoSubscriber, char *, char *,<< "No subscriber with identifier \""<< arg2<< "\" did subscribe to this object of class "<< arg1) | |
template<class Archive > | |
void | serialize (Archive &ar, const unsigned int version) |
Private Attributes | |
bool | backward |
![]() | |
std::vector< Entry > | entries |
Additional Inherited Members | |
![]() | |
typedef types::global_dof_index | size_type |
Inversion of a block-triangular matrix.
In this block matrix, the inverses of the diagonal blocks are stored together with the off-diagonal blocks of a block matrix. Then, forward or backward insertion is performed block-wise. The diagonal blocks are NOT inverted for this purpose!
Like for all preconditioners, the preconditioning operation is performed by the vmult() member function.
The implementation may be a little clumsy, but it should be sufficient as long as the block sizes are much larger than the number of blocks.
Here, we document the second part of examples/doxygen/block_matrix_array.cc
. For the beginning of this file, see BlockMatrixArray.
In order to set up the preconditioner, we have to compute the inverses of the diagonal blocks ourselves. Since we used FullMatrix objects, this is fairly easy.
After creating a 2x2 BlockTrianglePrecondition object, we only fill its diagonals. The scaling factor 1/2 used for A
is the reciprocal of the scaling factor used for the matrix
itself. Remember, this preconditioner actually multiplies with the diagonal blocks.
Now, we have a block Jacobi preconditioner, which is still symmetric, since the blocks are symmetric. Therefore, we can still use the preconditioned conjugate gradient method.
Now, we enter the subdiagonal block. This is the same as in matrix
.
Since the preconditioner is not symmetric anymore, we use the GMRES method for solving.
Definition at line 480 of file block_matrix_array.h.
typedef types::global_dof_index BlockTrianglePrecondition< number >::size_type |
Declare type for container size.
Definition at line 487 of file block_matrix_array.h.
BlockTrianglePrecondition< number >::BlockTrianglePrecondition | ( | ) |
Default constructor creating a useless object. initialize() must be called before using it.
BlockTrianglePrecondition< number >::BlockTrianglePrecondition | ( | const unsigned int | n_blocks | ) |
Constructor. This matrix must be block-quadratic, and n_blocks
is the number of blocks in each direction.
BlockTrianglePrecondition< number >::BlockTrianglePrecondition | ( | const unsigned int | n_block_rows, |
VectorMemory< Vector< number > > & | mem, | ||
const bool | backward = false |
||
) |
Constructor. This matrix must be block-quadratic. The additional parameter allows for backward insertion instead of forward.
void BlockTrianglePrecondition< number >::initialize | ( | const unsigned int | n_block_rows, |
VectorMemory< Vector< number > > & | mem, | ||
const bool | backward = false |
||
) |
Initialize object completely. This is the function to call for an object created by the default constructor.
void BlockTrianglePrecondition< number >::reinit | ( | const unsigned int | n_block_rows | ) |
Resize preconditioner to a new size and clear all blocks.
void BlockTrianglePrecondition< number >::enter | ( | const MATRIX & | matrix, |
const size_type | row, | ||
const size_type | col, | ||
const double | prefix = 1. , |
||
const bool | transpose = false |
||
) |
Enter a block. This calls BlockMatrixArray::enter(). Remember that the diagonal blocks should actually be inverse matrices or preconditioners.
void BlockTrianglePrecondition< number >::enter_aux | ( | VectorMemory< Vector< double > > & | mem, |
const MATRIX & | matrix, | ||
const size_type | row, | ||
const size_type | col, | ||
const double | prefix = 1. , |
||
const bool | transpose = false |
||
) |
Enter a block. This calls BlockMatrixArray::enter_aux(). Remember that the diagonal blocks should actually be inverse matrices or preconditioners.
void BlockTrianglePrecondition< number >::vmult | ( | BlockVector< number > & | dst, |
const BlockVector< number > & | src | ||
) | const |
Preconditioning.
void BlockTrianglePrecondition< number >::vmult_add | ( | BlockVector< number > & | dst, |
const BlockVector< number > & | src | ||
) | const |
Preconditioning adding to dst
.
void BlockTrianglePrecondition< number >::Tvmult | ( | BlockVector< number > & | dst, |
const BlockVector< number > & | src | ||
) | const |
Transposed preconditioning
void BlockTrianglePrecondition< number >::Tvmult_add | ( | BlockVector< number > & | dst, |
const BlockVector< number > & | src | ||
) | const |
Transposed preconditioning adding to dst
.
|
private |
Add all off-diagonal contributions and return the entry of the diagonal element for one row.
|
private |
Flag for backward insertion.
Definition at line 655 of file block_matrix_array.h.