From 63f9fe6fc41f2e4267cdcdefc2e552d208488e8a Mon Sep 17 00:00:00 2001 From: kanschat Date: Wed, 29 Sep 2010 17:08:51 +0000 Subject: [PATCH] we NEED different row and column indices for the edge matrices in multigrid git-svn-id: https://svn.dealii.org/trunk@22196 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/lac/include/lac/matrix_block.h | 135 +++++++++++++++++-------- 1 file changed, 92 insertions(+), 43 deletions(-) diff --git a/deal.II/lac/include/lac/matrix_block.h b/deal.II/lac/include/lac/matrix_block.h index 260d493855..b2edc2a8c5 100644 --- a/deal.II/lac/include/lac/matrix_block.h +++ b/deal.II/lac/include/lac/matrix_block.h @@ -125,7 +125,8 @@ class MatrixBlock * Reinitialize the matrix for a * new BlockSparsityPattern. This * adujusts the #matrix as well - * as the #block_indices. + * as the #row_indices and + * #column_indices. * * @note The row and column block * structure of the sparsity @@ -133,7 +134,8 @@ class MatrixBlock */ void reinit(const BlockSparsityPattern& sparsity); - operator MATRIX() const; + operator MATRIX&(); + operator const MATRIX&() const; /** * Add value to the @@ -333,7 +335,7 @@ class MatrixBlock /** * The block number computed from * an index by using - * #block_indices does not match + * BlockIndices does not match * the block coordinates stored * in this object. */ @@ -360,14 +362,22 @@ class MatrixBlock MATRIX matrix; private: /** - * The BlockIndices of the whole - * system. Using row() and - * column(), this allows us to - * find the index of the first - * row and column degrees of - * freedom for this block. + * The rwo BlockIndices of the + * whole system. Using row(), + * this allows us to find the + * index of the first row degree + * of freedom for this block. */ - BlockIndices block_indices; + BlockIndices row_indices; + /** + * The column BlockIndices of the + * whole system. Using column(), + * this allows us to find the + * index of the first column + * degree of freedom for this + * block. + */ + BlockIndices column_indices; friend void internal::reinit<>(MatrixBlock&, const BlockSparsityPattern&); }; @@ -383,9 +393,16 @@ class MatrixBlock * @author Baerbel Janssen, Guido Kanschat, 2010 */ template -class MatrixBlockVector : public NamedData > > +class MatrixBlockVector + : + private NamedData > > { public: + /** + * The type of object stored. + */ + typedef MatrixBlock value_type; + /** * Add a new matrix block at the * position (row,column) @@ -430,13 +447,13 @@ class MatrixBlockVector : public NamedDatai. */ - const MatrixBlock& block(unsigned int i) const; + const value_type& block(unsigned int i) const; /** * Access a reference to * the block at position i. */ - MatrixBlock& block(unsigned int i); + value_type& block(unsigned int i); /** * Access the matrix at position @@ -445,7 +462,9 @@ class MatrixBlockVector : public NamedData >::size; }; @@ -460,9 +479,14 @@ class MatrixBlockVector : public NamedData class MGMatrixBlockVector - : public NamedData > > > + : private NamedData > > > { public: + /** + * The type of object stored. + */ + typedef MGLevelObject > value_type; + /** * Add a new matrix block at the * position (row,column) @@ -506,18 +530,24 @@ class MGMatrixBlockVector * Access a constant reference to * the block at position i. */ - const MGLevelObject >& block(unsigned int i) const; + const value_type& block(unsigned int i) const; /** * Access a reference to * the block at position i. */ - MGLevelObject >& block(unsigned int i); + value_type& block(unsigned int i); /** * The memory used by this object. */ unsigned int memory_consumption () const; + + Subscriptor::subscribe; + Subscriptor::unsubscribe; + NamedData >::size; + + private: }; @@ -529,8 +559,8 @@ namespace internal void reinit(MatrixBlock& v, const BlockSparsityPattern& p) { - Assert(p.get_row_indices() == p.get_column_indices(), ExcNotImplemented()); - v.block_indices = p.get_row_indices(); + v.row_indices = p.get_row_indices(); + v.column_indices = p.get_column_indices(); } @@ -538,8 +568,8 @@ namespace internal void reinit(MatrixBlock >& v, const BlockSparsityPattern& p) { - Assert(p.get_row_indices() == p.get_column_indices(), ExcNotImplemented()); - v.block_indices = p.get_row_indices(); + v.row_indices = p.get_row_indices(); + v.column_indices = p.get_column_indices(); v.matrix.reinit(p.block(v.row, v.column)); } } @@ -562,7 +592,8 @@ MatrixBlock::MatrixBlock(const MatrixBlock& M) row(M.row), column(M.column), matrix(M.matrix), - block_indices(M.block_indices) + row_indices(M.row_indices), + column_indices(M.column_indices) {} @@ -585,7 +616,15 @@ MatrixBlock::reinit(const BlockSparsityPattern& sparsity) template inline -MatrixBlock::operator MATRIX() const +MatrixBlock::operator MATRIX&() +{ + return matrix; +} + + +template +inline +MatrixBlock::operator const MATRIX&() const { return matrix; } @@ -598,12 +637,13 @@ MatrixBlock::add ( const unsigned int gj, const typename MATRIX::value_type value) { - Assert(block_indices.size() != 0, ExcNotInitialized()); + Assert(row_indices.size() != 0, ExcNotInitialized()); + Assert(column_indices.size() != 0, ExcNotInitialized()); const std::pair bi - = block_indices.global_to_local(gi); + = row_indices.global_to_local(gi); const std::pair bj - = block_indices.global_to_local(gj); + = column_indices.global_to_local(gj); Assert (bi.first == row, ExcBlockIndexMismatch(bi.first, row)); Assert (bj.first == column, ExcBlockIndexMismatch(bj.first, column)); @@ -616,18 +656,19 @@ template template inline void -MatrixBlock::add (const std::vector& row_indices, - const std::vector& col_indices, +MatrixBlock::add (const std::vector& r_indices, + const std::vector& c_indices, const FullMatrix& values, const bool elide_zero_values) { - Assert(block_indices.size() != 0, ExcNotInitialized()); - - AssertDimension (row_indices.size(), values.m()); - AssertDimension (col_indices.size(), values.n()); + Assert(row_indices.size() != 0, ExcNotInitialized()); + Assert(column_indices.size() != 0, ExcNotInitialized()); + + AssertDimension (r_indices.size(), values.m()); + AssertDimension (c_indices.size(), values.n()); for (unsigned int i=0; i::add (const unsigned int b_row, const bool, const bool) { - Assert(block_indices.size() != 0, ExcNotInitialized()); + Assert(row_indices.size() != 0, ExcNotInitialized()); + Assert(column_indices.size() != 0, ExcNotInitialized()); + const std::pair bi - = block_indices.global_to_local(b_row); + = row_indices.global_to_local(b_row); // In debug mode, we check whether // all indices are in the correct @@ -661,7 +704,7 @@ MatrixBlock::add (const unsigned int b_row, for (unsigned int j=0;j bj - = block_indices.global_to_local(col_indices[j]); + = column_indices.global_to_local(col_indices[j]); Assert(bj.first == column, ExcBlockIndexMismatch(bj.first, column)); matrix.add(bi.second, bj.second, values[j]); @@ -678,8 +721,9 @@ MatrixBlock::add (const std::vector &indices, const FullMatrix &values, const bool elide_zero_values) { - Assert(block_indices.size() != 0, ExcNotInitialized()); - + Assert(row_indices.size() != 0, ExcNotInitialized()); + Assert(column_indices.size() != 0, ExcNotInitialized()); + AssertDimension (indices.size(), values.m()); Assert (values.n() == values.m(), ExcNotQuadratic()); @@ -699,7 +743,9 @@ MatrixBlock::add (const unsigned int row, const std::vector &values, const bool elide_zero_values) { - Assert(block_indices.size() != 0, ExcNotInitialized()); + Assert(row_indices.size() != 0, ExcNotInitialized()); + Assert(column_indices.size() != 0, ExcNotInitialized()); + AssertDimension (col_indices.size(), values.size()); add (row, col_indices.size(), &col_indices[0], &values[0], elide_zero_values); @@ -765,8 +811,8 @@ MatrixBlockVector::add( unsigned int row, unsigned int column, const std::string& name) { - std_cxx1x::shared_ptr > p(new MatrixBlock(row, column)); - NamedData > >::add(p, name); + std_cxx1x::shared_ptr p(new value_type(row, column)); + NamedData >::add(p, name); } @@ -832,8 +878,11 @@ MGMatrixBlockVector::add( const std::string& name) { std_cxx1x::shared_ptr > > - p(new MGLevelObject >(row, column)); - NamedData > >::add(p, name); + p(new value_type(0, 1)); + (*p)[0].row = row; + (*p)[0].column = column; + + NamedData >::add(p, name); } -- 2.39.5