From: kanschat Date: Thu, 18 Mar 2010 04:46:44 +0000 (+0000) Subject: introduce MatrixBlockVector X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=b08337705c5740fcc850544f025f276d8c331756;p=dealii-svn.git introduce MatrixBlockVector git-svn-id: https://svn.dealii.org/trunk@20842 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/numerics/mesh_worker_assembler.h b/deal.II/deal.II/include/numerics/mesh_worker_assembler.h index 36524292e1..015bfca1e9 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_assembler.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_assembler.h @@ -662,9 +662,6 @@ namespace MeshWorker class MatrixLocalBlocksToGlobalBlocks { public: - /// The object that is stored - typedef boost::shared_ptr > MatrixPtr; - /** * Constructor, initializing * the #threshold, which @@ -681,7 +678,7 @@ namespace MeshWorker * cell matrix vectors. */ void initialize(const BlockInfo* block_info, - std::vector& matrices); + MatrixBlockVector& matrices); /** * Initialize the local data * in the @@ -737,7 +734,8 @@ namespace MeshWorker * stored as a vector of * pointers. */ - std::vector matrices; + SmartPointer, + MatrixLocalBlocksToGlobalBlocks > matrices; /** * A pointer to the object containing the block structure. @@ -1400,10 +1398,12 @@ namespace MeshWorker template inline void - MatrixLocalBlocksToGlobalBlocks::initialize(const BlockInfo* b, std::vector& m) + MatrixLocalBlocksToGlobalBlocks::initialize( + const BlockInfo* b, + MatrixBlockVector& m) { block_info = b; - matrices = m; + matrices = &m; } @@ -1415,7 +1415,7 @@ namespace MeshWorker DoFInfo& info, bool interior_face) const { - info.initialize_matrices(matrices, interior_face); + info.initialize_matrices(*matrices, interior_face); } @@ -1467,14 +1467,14 @@ namespace MeshWorker MatrixLocalBlocksToGlobalBlocks::assemble( const DoFInfo& info) { - for (unsigned int i=0;isize();++i) { // Row and column index of // the block we are dealing with - const unsigned int row = matrices[i]->row; - const unsigned int col = matrices[i]->column; + const unsigned int row = matrices->block(i).row; + const unsigned int col = matrices->block(i).column; - assemble(matrices[i]->matrix, info.matrix(i,false).matrix, row, col, info.indices, info.indices); + assemble(matrices->matrix(i), info.matrix(i,false).matrix, row, col, info.indices, info.indices); } } @@ -1486,17 +1486,17 @@ namespace MeshWorker const DoFInfo& info1, const DoFInfo& info2) { - for (unsigned int i=0;isize();++i) { // Row and column index of // the block we are dealing with - const unsigned int row = matrices[i]->row; - const unsigned int col = matrices[i]->column; + const unsigned int row = matrices->block(i).row; + const unsigned int col = matrices->block(i).column; - assemble(matrices[i]->matrix, info1.matrix(i,false).matrix, row, col, info1.indices, info1.indices); - assemble(matrices[i]->matrix, info1.matrix(i,true).matrix, row, col, info1.indices, info2.indices); - assemble(matrices[i]->matrix, info2.matrix(i,false).matrix, row, col, info2.indices, info2.indices); - assemble(matrices[i]->matrix, info2.matrix(i,true).matrix, row, col, info2.indices, info1.indices); + assemble(matrices->matrix(i), info1.matrix(i,false).matrix, row, col, info1.indices, info1.indices); + assemble(matrices->matrix(i), info1.matrix(i,true).matrix, row, col, info1.indices, info2.indices); + assemble(matrices->matrix(i), info2.matrix(i,false).matrix, row, col, info2.indices, info2.indices); + assemble(matrices->matrix(i), info2.matrix(i,true).matrix, row, col, info2.indices, info1.indices); } } diff --git a/deal.II/deal.II/include/numerics/mesh_worker_info.h b/deal.II/deal.II/include/numerics/mesh_worker_info.h index 5e69d90ee3..9b13df27d8 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_info.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_info.h @@ -84,8 +84,8 @@ namespace MeshWorker * instantiation in order to * provide row and column info. */ - template - void initialize_matrices(const std::vector& matrices, + template + void initialize_matrices(const MatrixBlockVector& matrices, bool both); /** @@ -672,7 +672,8 @@ namespace MeshWorker * #values, #gradients and * #hessians. */ - void fill_local_data(const DoFInfo& info, const bool split_fevalues); + template + void fill_local_data(const DoFInfo& info, bool split_fevalues); /** * The global data vector @@ -782,10 +783,10 @@ namespace MeshWorker template - template + template inline void LocalResults::initialize_matrices( - const std::vector& matrices, + const MatrixBlockVector& matrices, bool both) { M1.resize(matrices.size()); @@ -793,8 +794,8 @@ namespace MeshWorker M2.resize(matrices.size()); for (unsigned int i=0;irow; - const unsigned int col = matrices[i]->column; + const unsigned int row = matrices.block(i).row; + const unsigned int col = matrices.block(i).column; M1[i].row = row; M1[i].column = col; @@ -1142,6 +1143,35 @@ namespace MeshWorker //----------------------------------------------------------------------// + template + template + inline void + IntegrationInfo::initialize( + const FiniteElement& el, + const Mapping& mapping, + const Quadrature& quadrature, + const UpdateFlags flags, + const BlockInfo* block_info) + { + if (block_info == 0 || block_info->local().size() == 0) + { + fevalv.resize(1); + fevalv[0] = boost::shared_ptr > ( + new FEVALUES (mapping, el, quadrature, flags)); + } + else + { + fevalv.resize(el.n_base_elements()); + for (unsigned int i=0;i > ( + new FEVALUES (mapping, el.base_element(i), quadrature, flags)); + } + } + n_components = el.n_components(); + } + + template inline const FEValuesBase& IntegrationInfo::fe_values() const diff --git a/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h b/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h index dadacaa15d..b22fbf024c 100644 --- a/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h +++ b/deal.II/deal.II/include/numerics/mesh_worker_info.templates.h @@ -119,35 +119,6 @@ namespace MeshWorker } - template - template - void - IntegrationInfo::initialize( - const FiniteElement& el, - const Mapping& mapping, - const Quadrature& quadrature, - const UpdateFlags flags, - const BlockInfo* block_info) - { - if (block_info == 0 || block_info->local().size() == 0) - { - fevalv.resize(1); - fevalv[0] = boost::shared_ptr > ( - new FEVALUES (mapping, el, quadrature, flags)); - } - else - { - fevalv.resize(el.n_base_elements()); - for (unsigned int i=0;i > ( - new FEVALUES (mapping, el.base_element(i), quadrature, flags)); - } - } - n_components = el.n_components(); - } - - template void IntegrationInfo::initialize_data( @@ -212,8 +183,9 @@ namespace MeshWorker template + template void - IntegrationInfo::fill_local_data(const DoFInfo& info, bool split_fevalues) + IntegrationInfo::fill_local_data(const DoFInfo& info, bool split_fevalues) { if (split_fevalues) { diff --git a/deal.II/deal.II/source/numerics/mesh_worker_info.cc b/deal.II/deal.II/source/numerics/mesh_worker_info.cc index c53b15956a..83c742746f 100644 --- a/deal.II/deal.II/source/numerics/mesh_worker_info.cc +++ b/deal.II/deal.II/source/numerics/mesh_worker_info.cc @@ -26,25 +26,30 @@ DEAL_II_NAMESPACE_OPEN namespace MeshWorker { + template class IntegrationInfo; + template class LocalResults; template class DoFInfo; + template void IntegrationInfo::fill_local_data( + const DoFInfo&, bool); template class LocalResults; template class DoFInfo; - template class IntegrationInfo; + template void IntegrationInfo::fill_local_data( + const DoFInfo&, bool); - template void IntegrationInfo - ::initialize >( - const FiniteElement&, const Mapping&, - const Quadrature::integral_dimension>&, const UpdateFlags, const BlockInfo*); - template void IntegrationInfo - ::initialize >( - const FiniteElement&, const Mapping&, - const Quadrature::integral_dimension>&, const UpdateFlags, const BlockInfo*); - template void IntegrationInfo - ::initialize >( - const FiniteElement&, const Mapping&, - const Quadrature::integral_dimension>&, const UpdateFlags, const BlockInfo*); +// template void IntegrationInfo +// ::initialize >( +// const FiniteElement&, const Mapping&, +// const Quadrature::integral_dimension>&, const UpdateFlags, const BlockInfo*); +// template void IntegrationInfo +// ::initialize >( +// const FiniteElement&, const Mapping&, +// const Quadrature::integral_dimension>&, const UpdateFlags, const BlockInfo*); +// template void IntegrationInfo +// ::initialize >( +// const FiniteElement&, const Mapping&, +// const Quadrature::integral_dimension>&, const UpdateFlags, const BlockInfo*); } #endif diff --git a/deal.II/lac/include/lac/matrix_block.h b/deal.II/lac/include/lac/matrix_block.h index 6a0a6db419..93d26832fa 100644 --- a/deal.II/lac/include/lac/matrix_block.h +++ b/deal.II/lac/include/lac/matrix_block.h @@ -2,7 +2,7 @@ // $Id$ // Version: $Name$ // -// Copyright (C) 2007, 2008, 2009 by the deal.II authors +// Copyright (C) 2007, 2008, 2009, 2010 by the deal.II authors // // This file is subject to QPL and may not be distributed // without copyright and license information. Please refer @@ -14,6 +14,9 @@ #define __deal2__matrix_block_h #include +#include + +#include DEAL_II_NAMESPACE_OPEN @@ -90,6 +93,41 @@ struct MatrixBlock }; +/** + * A vector of MatrixBlock, which is implemented using shared + * pointers, in order to allow for copying and rearranging. Each + * matrix block can be identified by name. + * + * @author Baerbel Janssen, Guido Kanschat, 2010 + */ +template +class MatrixBlockVector : public NamedData > > +{ + public: + /** + * Add a new matrix block at the + * position + * in the block system. + */ + void add(unsigned int row, unsigned int column, const std::string& name); + /** + * Access a constant reference to + * the block at position i. + */ + const MatrixBlock& block(unsigned int i) const; + /** + * Access the matrix at position + * i for read and write + * access. + */ + MATRIX& matrix(unsigned int i); + + +}; + + +//----------------------------------------------------------------------// + template MatrixBlock::MatrixBlock() : @@ -113,6 +151,33 @@ MatrixBlock::MatrixBlock(unsigned int i, unsigned int j) row(i), column(j) {} +//----------------------------------------------------------------------// + +template +inline void +MatrixBlockVector::add(unsigned int row, unsigned int column, const std::string& name) +{ + boost::shared_ptr > p(new MatrixBlock(row, column)); + NamedData > >::add(p, name); +} + + +template +inline const MatrixBlock& +MatrixBlockVector::block(unsigned int i) const +{ + return *this->read(i); +} + + +template +inline MATRIX& +MatrixBlockVector::matrix(unsigned int i) +{ + return (*this)(i)->matrix; +} + + DEAL_II_NAMESPACE_CLOSE