Classes | |
class | ExcMultipleDiagonal |
class | ExcNoDiagonal |
Public Member Functions | |
BlockTrianglePrecondition () | |
BlockTrianglePrecondition (unsigned int n_block_rows, VectorMemory< Vector< number > > &mem, bool backward=false) | |
void | initialize (const unsigned int n_block_rows, VectorMemory< Vector< number > > &mem, bool backward=false) |
void | reinit (const unsigned int n_block_rows) |
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< double > > &mem, const MATRIX &matrix, const unsigned int row, const unsigned int col, const double prefix=1., const bool transpose=false) |
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 |
Private Member Functions | |
void | do_row (BlockVector< number > &dst, unsigned int row_num) const |
Private Attributes | |
bool | backward |
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.
deallog << "Error " << x.l2_norm() << std::endl; deallog << "Error A-norm " << std::sqrt(matrix.matrix_norm_square(x)) << std::endl; FullMatrix<float> Ainv(4,4); Ainv.invert(A); FullMatrix<float> Cinv(2,2); Cinv.invert(C);
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.
BlockTrianglePrecondition<double> precondition(2, simple_mem); precondition.enter(Ainv,0,0,.5); precondition.enter(Cinv,1,1);
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.
cg.solve(matrix, x, y, precondition); x.add(-1., result); deallog << "Error " << x.l2_norm() << std::endl;
Now, we enter the subdiagonal block. This is the same as in matrix
.
precondition.enter(B1,1,0,-1., true);
precondition.enter(B2,1,0,1.);
Since the preconditioner is not symmetric anymore, we use the GMRES method for solving.
SolverGMRES<BlockVector<double> > gmres(control, mem); gmres.solve(matrix, x, y, precondition); x.add(-1., result); deallog << "Error " << x.l2_norm() << std::endl;
BlockTrianglePrecondition< number >::BlockTrianglePrecondition | ( | ) |
Default constructor creating a useless object. initialize() must be called before using it.
BlockTrianglePrecondition< number >::BlockTrianglePrecondition | ( | unsigned int | n_block_rows, | |
VectorMemory< Vector< number > > & | mem, | |||
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, | |||
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 unsigned int | row, | |||
const unsigned int | col, | |||
const double | prefix = 1. , |
|||
const bool | transpose = false | |||
) | [inline] |
Enter a block. This calls BlockMatrixArray::enter(). Remember that the diagonal blocks should actually be inverse matrices or preconditioners.
Reimplemented from BlockMatrixArray< number >.
void BlockTrianglePrecondition< number >::enter_aux | ( | VectorMemory< Vector< double > > & | mem, | |
const MATRIX & | matrix, | |||
const unsigned int | row, | |||
const unsigned int | col, | |||
const double | prefix = 1. , |
|||
const bool | transpose = false | |||
) | [inline] |
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.
Reimplemented from BlockMatrixArray< number >.
void BlockTrianglePrecondition< number >::vmult_add | ( | BlockVector< number > & | dst, | |
const BlockVector< number > & | src | |||
) | const |
Preconditioning adding to dst
.
Reimplemented from BlockMatrixArray< number >.
void BlockTrianglePrecondition< number >::Tvmult | ( | BlockVector< number > & | dst, | |
const BlockVector< number > & | src | |||
) | const |
Transposed preconditioning
Reimplemented from BlockMatrixArray< number >.
void BlockTrianglePrecondition< number >::Tvmult_add | ( | BlockVector< number > & | dst, | |
const BlockVector< number > & | src | |||
) | const |
Transposed preconditioning adding to dst
.
Reimplemented from BlockMatrixArray< number >.
void BlockTrianglePrecondition< number >::do_row | ( | BlockVector< number > & | dst, | |
unsigned int | row_num | |||
) | const [private] |
Add all off-diagonal contributions and return the entry of the diagonal element for one row.
bool BlockTrianglePrecondition< number >::backward [private] |
Flag for backward insertion.