From: guido Date: Wed, 6 Dec 2000 16:41:04 +0000 (+0000) Subject: BlockDiagonalMatrix X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b96d414854587943ce870eeeb7c2caa803ea726b;p=dealii-svn.git BlockDiagonalMatrix git-svn-id: https://svn.dealii.org/trunk@3524 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/block_matrix.h b/deal.II/lac/include/lac/block_matrix.h new file mode 100644 index 0000000000..ac4c7fc99f --- /dev/null +++ b/deal.II/lac/include/lac/block_matrix.h @@ -0,0 +1,114 @@ +//------------------------------------------------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 1998, 1999, 2000 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//------------------------------------------------------------------- +#ifndef __deal2__block_matrix_h +#define __deal2__block_matrix_h + + +#include +#include +#include + +/** + * A matrix with several copies of the same block on the diagonal. + * + * This matrix implements an @p{m} by @p{m} block matrix. Each + * diagonal block consists of the same (non-block) matrix, while + * off-diagonal blocks are void. + * + * One special application is a one by one block matrix, allowing to + * apply the @p{vmult} of the original matrix (or preconditioner) to a + * block vector. + * + * @author Guido Kanschat, 2000 + */ +template +class BlockDiagonalMatrix : public Subscriptor +{ +public: + /** + * Constructor for an @p{n_blocks} + * by @p{n_blocks} matrix with + * diagonal blocks @p{M}. + */ + BlockDiagonalMatrix (const MATRIX& M, + unsigned int n_blocks); + + /** + * Matrix-vector-multiplication. + */ + template + void vmult (BlockVector& dst, + const BlockVector& src) const; + + /** + * Transposed matrix-vector-multiplication. + */ + template + void Tvmult (BlockVector& dst, + const BlockVector& src) const; +private: + /** + * Number of blocks. + */ + unsigned int num_blocks; + /** + * Diagonal entry. + */ + SmartPointer matrix; +}; + + +//----------------------------------------------------------------------// + +template +BlockDiagonalMatrix::BlockDiagonalMatrix (const MATRIX& M, + unsigned int num_blocks) + : num_blocks (num_blocks), + matrix(&M) +{} + + +template +template +void +BlockDiagonalMatrix::vmult (BlockVector& dst, + const BlockVector& src) const +{ + Assert (dst.n_blocks()==num_blocks, + ExcDimensionMismatch(dst.n_blocks(),num_blocks)); + Assert (src.n_blocks()==num_blocks, + ExcDimensionMismatch(src.n_blocks(),num_blocks)); + + for (unsigned int i=0;ivmult (dst.block(i), src.block(i)); +} + + +template +template +void +BlockDiagonalMatrix::Tvmult (BlockVector& dst, + const BlockVector& src) const +{ + Assert (dst.n_blocks()==num_blocks, + ExcDimensionMismatch(dst.n_blocks(),num_blocks)); + Assert (src.n_blocks()==num_blocks, + ExcDimensionMismatch(src.n_blocks(),num_blocks)); + + for (unsigned int i=0;iTvmult (dst.block(i), src.block(i)); +} + + +#endif +