From: Wolfgang Bangerth Date: Mon, 8 Sep 2003 23:40:50 +0000 (+0000) Subject: Allow vectors for some operations on BlockSparseMatrices, if the respective number... X-Git-Tag: v8.0.0~16245 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=63441ed02de300e237618187a10e0184654590cf;p=dealii.git Allow vectors for some operations on BlockSparseMatrices, if the respective number of blocks in this direction is equal to one. git-svn-id: https://svn.dealii.org/trunk@7965 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/lac/include/lac/block_sparse_matrix.h b/deal.II/lac/include/lac/block_sparse_matrix.h index af705f8eae..23722fe858 100644 --- a/deal.II/lac/include/lac/block_sparse_matrix.h +++ b/deal.II/lac/include/lac/block_sparse_matrix.h @@ -21,7 +21,9 @@ #include -template class BlockVector; +template class Vector; +template class BlockVector; + /*! @addtogroup Matrix1 *@{ @@ -530,6 +532,39 @@ class BlockSparseMatrix : public Subscriptor template void vmult (BlockVector &dst, const BlockVector &src) const; + + /** + * Matrix-vector + * multiplication. Just like the + * previous function, but only + * applicable if the matrix has + * only one block column. + */ + template + void vmult (BlockVector &dst, + const Vector &src) const; + + /** + * Matrix-vector + * multiplication. Just like the + * previous function, but only + * applicable if the matrix has + * only one block row. + */ + template + void vmult (Vector &dst, + const BlockVector &src) const; + + /** + * Matrix-vector + * multiplication. Just like the + * previous function, but only + * applicable if the matrix has + * only one block. + */ + template + void vmult (Vector &dst, + const Vector &src) const; /** * Matrix-vector multiplication: @@ -543,6 +578,39 @@ class BlockSparseMatrix : public Subscriptor void Tvmult (BlockVector &dst, const BlockVector &src) const; + /** + * Matrix-vector + * multiplication. Just like the + * previous function, but only + * applicable if the matrix has + * only one block row. + */ + template + void Tvmult (BlockVector &dst, + const Vector &src) const; + + /** + * Matrix-vector + * multiplication. Just like the + * previous function, but only + * applicable if the matrix has + * only one block column. + */ + template + void Tvmult (Vector &dst, + const BlockVector &src) const; + + /** + * Matrix-vector + * multiplication. Just like the + * previous function, but only + * applicable if the matrix has + * only one block. + */ + template + void Tvmult (Vector &dst, + const Vector &src) const; + /** * Adding Matrix-vector * multiplication. Add $M*src$ on @@ -635,6 +703,33 @@ class BlockSparseMatrix : public Subscriptor const BlockVector &src, const number omega = 1.) const; + /** + * Apply the Jacobi + * preconditioner, which + * multiplies every element of + * the @p{src} vector by the + * inverse of the respective + * diagonal element and + * multiplies the result with the + * relaxation parameter @p{omega}. + * + * All diagonal blocks must be + * square matrices for this + * operation. + * + * In contrast to the previous + * function, input and output + * vectors are not block + * vectors. This operation is + * thus only allowed if the block + * matrix consists of only one + * block. + */ + template + void precondition_Jacobi (Vector &dst, + const Vector &src, + const number omega = 1.) const; + /* Call print functions for * the SparseMatrix blocks. */ @@ -1197,6 +1292,58 @@ BlockSparseMatrix::vmult (BlockVector &dst, +template +template +void +BlockSparseMatrix::vmult (Vector &dst, + const BlockVector &src) const +{ + Assert (rows == 1, + ExcDimensionMismatch(1, rows)); + Assert (src.n_blocks() == columns, + ExcDimensionMismatch(src.n_blocks(), columns)); + + block(0,0).vmult (dst, src.block(0)); + for (unsigned int col=1; col +template +void +BlockSparseMatrix::vmult (BlockVector &dst, + const Vector &src) const +{ + Assert (dst.n_blocks() == rows, + ExcDimensionMismatch(dst.n_blocks(), rows)); + Assert (1 == columns, + ExcDimensionMismatch(1, columns)); + + for (unsigned int row=0; row +template +void +BlockSparseMatrix::vmult (Vector &dst, + const Vector &src) const +{ + Assert (1 == rows, + ExcDimensionMismatch(1, rows)); + Assert (1 == columns, + ExcDimensionMismatch(1, columns)); + + block(0,0).vmult (dst, src); +} + + + template template void @@ -1224,8 +1371,8 @@ BlockSparseMatrix::vmult_add (BlockVector &dst, template template void -BlockSparseMatrix::Tvmult (BlockVector& dst, - const BlockVector& src) const +BlockSparseMatrix::Tvmult (BlockVector &dst, + const BlockVector &src) const { Assert (dst.n_blocks() == n_block_cols(), ExcDimensionMismatch(dst.n_blocks(), n_block_cols())); @@ -1244,6 +1391,60 @@ BlockSparseMatrix::Tvmult (BlockVector& dst, +template +template +void +BlockSparseMatrix::Tvmult (BlockVector &dst, + const Vector &src) const +{ + Assert (dst.n_blocks() == n_block_cols(), + ExcDimensionMismatch(dst.n_blocks(), n_block_cols())); + Assert (1 == n_block_rows(), + ExcDimensionMismatch(1, n_block_rows())); + + dst = 0.; + + for (unsigned int col=0; col +template +void +BlockSparseMatrix::Tvmult (Vector &dst, + const BlockVector &src) const +{ + Assert (1 == n_block_cols(), + ExcDimensionMismatch(1, n_block_cols())); + Assert (src.n_blocks() == n_block_rows(), + ExcDimensionMismatch(src.n_blocks(), n_block_rows())); + + block(0,0).Tvmult (dst, src.block(0)); + + for (unsigned int row=1; row +template +void +BlockSparseMatrix::Tvmult (Vector &dst, + const Vector &src) const +{ + Assert (1 == n_block_cols(), + ExcDimensionMismatch(1, n_block_cols())); + Assert (1 == n_block_rows(), + ExcDimensionMismatch(1, n_block_rows())); + + block(0,0).Tvmult (dst, src); +} + + + template template void @@ -1384,6 +1585,25 @@ precondition_Jacobi (BlockVector &dst, +template +template +void +BlockSparseMatrix:: +precondition_Jacobi (Vector &dst, + const Vector &src, + const number omega) const +{ + Assert (rows == columns, ExcMatrixNotBlockSquare()); + Assert (1 == rows, + ExcDimensionMismatch(1, rows)); + Assert (1 == columns, + ExcDimensionMismatch(1, columns)); + + block(0,0).precondition_Jacobi (dst, src, omega); +} + + + template inline typename BlockSparseMatrix::const_iterator