From: guido Date: Fri, 8 Aug 2003 19:41:54 +0000 (+0000) Subject: selection with block grouping X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=a34911c06ab0d6e8e6e6c53a0cae54b35e6d9996;p=dealii-svn.git selection with block grouping git-svn-id: https://svn.dealii.org/trunk@7907 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/multigrid/mg_base.h b/deal.II/deal.II/include/multigrid/mg_base.h index 2be97591b5..bc0859e92f 100644 --- a/deal.II/deal.II/include/multigrid/mg_base.h +++ b/deal.II/deal.II/include/multigrid/mg_base.h @@ -110,6 +110,30 @@ class MGCoarseGridBase : public Subscriptor * context. This class is an abstract one and has no implementations of * possible algorithms for these operations. * + * There are several derived classes, reflecting the fact that vector + * types and numbering of the fine-grid discretization and of the + * multi-level implementation are independent. + * + * If you use multigrid for your complete system of equations, you + * will use @ref{MGTransferPrebuilt} together with Multigrid for + * Vector objects. The fine-grid vectors may still be of type + * BlockVector. + * + * For mixed systems, it may be required to do multigrid only for a + * single component or for some components. The classes + * @ref{MGTransferSelect} and @ref{MGTransferBlock} handle these cases. + * + * MGTransferSelect is used if you use mutligrid (on Vector oblects) + * for a single component, possibly grouped using + * @p{mg_target_component}. + * + * The class MGTransferBlock handles the case where your multigrid + * method operates on BlockVector objects. These can contain all or a + * consecutive set of the blocks of the complete system. Since most + * smoothers cannot operate on block structures, it is not clear + * whether this case is really useful. Therefore, a tested + * implementation of this case will be supplied when needed. + * * @author Wolfgang Bangerth, Guido Kanschat, 1999, 2002 */ template diff --git a/deal.II/deal.II/include/multigrid/mg_dof_tools.h b/deal.II/deal.II/include/multigrid/mg_dof_tools.h index 99ce9dcb26..3976262247 100644 --- a/deal.II/deal.II/include/multigrid/mg_dof_tools.h +++ b/deal.II/deal.II/include/multigrid/mg_dof_tools.h @@ -180,6 +180,38 @@ class MGTools MGLevelObject >& v, const std::vector& selected, const std::vector& target_component); + /** + * Adjust vectors on all levels + * to correct size. Count the + * numbers of degrees of freedom + * on each level component-wise + * in a single component. Then, + * assign @p{vector} the + * corresponding size. + * + * The boolean field @p{selected} + * may be nonzero in a single + * component, indicating the + * block of a block vector the + * argument @p{v} corresponds to. + * + * Degrees of freedom must be + * sorted by component in order + * to obtain reasonable results + * from this function. + * + * The argument + * @p{target_component} allows to + * re-sort and groupt components + * as in + * @p{DoFRenumbering::component_wise}. + */ + template + static void + reinit_vector (const MGDoFHandler& mg_dof, + MGLevelObject >& v, + const std::vector& selected, + const std::vector& target_component); }; diff --git a/deal.II/deal.II/include/multigrid/mg_transfer.h b/deal.II/deal.II/include/multigrid/mg_transfer.h index 3d405e64c3..8a5388648c 100644 --- a/deal.II/deal.II/include/multigrid/mg_transfer.h +++ b/deal.II/deal.II/include/multigrid/mg_transfer.h @@ -41,6 +41,9 @@ template class MGDoFHandler; * once by looping over all cells and storing the result in a matrix for * each level, but requires additional memory. * + * See @ref{MGTransferBase} to find out which of the transfer classes + * is best for your needs. + * * @author Wolfgang Bangerth, Guido Kanschat, 1999-2003 */ template @@ -163,11 +166,6 @@ class MGTransferPrebuilt : public MGTransferBase > private: - /** - * Selected component. - */ - unsigned int selected; - /** * Sizes of the multi-level vectors. */ @@ -253,7 +251,7 @@ class MGTransferPrebuilt : public MGTransferBase > /** * Implementation of matrix generation for @ref{MGTransferBlock} and - * @p{MGTransferSelected}. + * @p{MGTransferSelect}. * * @author Guido Kanschat, 2001-2003 */ @@ -291,6 +289,13 @@ class MGTransferBlockBase */ std::vector target_component; + /** + * Target component if + * renumbering of level vectors + * is required. + */ + std::vector mg_target_component; + /** * Sizes of the multi-level vectors. */ @@ -343,11 +348,18 @@ class MGTransferBlockBase #endif }; - +//TODO:[GK] Update this class /** * Implementation of the @p{MGTransferBase} interface for block - * matrices and block vectors. In addition to the functionality of + * matrices and block vectors. + * + * Warning! Due to additional requirements on MGTransferSelect, the + * implementation in the base class has changed. Therefore, this class + * is left in an untested state. If you use it and you encounter + * problems, please contact Guido Kanschat. + * + * In addition to the functionality of * @ref{MGTransferPrebuilt}, the operation may be restricted to * certain blocks of the vector. * @@ -355,6 +367,9 @@ class MGTransferBlockBase * transfer routines may only have as many components as there are * @p{true}s in the selected-field. * + * See @ref{MGTransferBase} to find out which of the transfer classes + * is best for your needs. + * * @author Guido Kanschat, 2001, 2002 */ template @@ -391,6 +406,8 @@ class MGTransferBlock : public MGTransferBase >, void build_matrices (const MGDoFHandler &mg_dof, std::vector selected = std::vector(), const std::vector& target_component + = std::vector(), + const std::vector& mg_target_component =std::vector()); /** @@ -534,13 +551,17 @@ class MGTransferBlock : public MGTransferBase >, }; +//TODO:[GK] Update documentation for copy_* functions /** * Implementation of the @p{MGTransferBase} interface for block - * matrices and simple vectors. This class uses @ref{MGTransferBlock} + * matrices and simple vectors. This class uses @ref{MGTransferBlockBase} * selecting a single component or grouping several components into a * single block. The transfer operators themselves are implemented for - * simple vectors again. + * Vector and BlockVector objects. + * + * See @ref{MGTransferBase} to find out which of the transfer classes + * is best for your needs. * * @author Guido Kanschat, 2001, 2002, 2003 */ @@ -558,8 +579,13 @@ class MGTransferSelect : public MGTransferBase >, * Actually build the prolongation * matrices for grouped components. * - * @p{selected} is the number of - * the component for which the + * @p{selected} tells the copy + * functions operating on single + * vectors, which component this + * vector belongs to. + * + * @p{mg_selected} is the number + * of the component for which the * transfer matrices should be * built. * @@ -580,9 +606,18 @@ class MGTransferSelect : public MGTransferBase >, template void build_matrices (const MGDoFHandler &mg_dof, unsigned int selected, + unsigned int mg_selected, const std::vector& target_component + = std::vector(), + const std::vector& mg_target_component = std::vector()); + /** + * Change selected + * component. Handle with care! + */ + void select (const unsigned int component); + /** * Prolongate a vector from level * @p{to_level-1} to level @@ -648,11 +683,11 @@ class MGTransferSelect : public MGTransferBase >, * on each of the levels * separately, i.a. an @p{MGVector}. */ - template + template void copy_to_mg (const MGDoFHandler &mg_dof, MGLevelObject > &dst, - const InVector &src) const; + const Vector &src) const; /** * Transfer from multi-level vector to @@ -666,10 +701,10 @@ class MGTransferSelect : public MGTransferBase >, * constrained degrees of freedom * are set to zero. */ - template + template void copy_from_mg (const MGDoFHandler &mg_dof, - OutVector &dst, + Vector &dst, const MGLevelObject > &src) const; /** @@ -680,18 +715,90 @@ class MGTransferSelect : public MGTransferBase >, * function, but probably not for * continuous elements. */ - template + template void copy_from_mg_add (const MGDoFHandler &mg_dof, - OutVector &dst, + Vector &dst, + const MGLevelObject > &src) const; + + /** + * Transfer from a vector on the + * global grid to vectors defined + * on each of the levels + * separately, i.a. an @p{MGVector}. + */ + template + void + copy_to_mg (const MGDoFHandler &mg_dof, + MGLevelObject > &dst, + const BlockVector &src) const; + + /** + * Transfer from multi-level vector to + * normal vector. + * + * Copies data from active + * portions of an MGVector into + * the respective positions of a + * @p{Vector}. In order to + * keep the result consistent, + * constrained degrees of freedom + * are set to zero. + */ + template + void + copy_from_mg (const MGDoFHandler &mg_dof, + BlockVector &dst, + const MGLevelObject > &src) const; + + /** + * Add a multi-level vector to a + * normal vector. + * + * Works as the previous + * function, but probably not for + * continuous elements. + */ + template + void + copy_from_mg_add (const MGDoFHandler &mg_dof, + BlockVector &dst, const MGLevelObject > &src) const; private: - /** - * Selected component. - */ - unsigned int selected; + /** + * Transfer from multi-level vector to + * normal vector. + * + * Copies data from active + * portions of an MGVector into + * the respective positions of a + * @p{Vector}. In order to + * keep the result consistent, + * constrained degrees of freedom + * are set to zero. + */ + template + void + do_copy_from_mg (const MGDoFHandler &mg_dof, + OutVector &dst, + const MGLevelObject > &src, + const unsigned int offset) const; + /** + * Add a multi-level vector to a + * normal vector. + * + * Works as the previous + * function, but probably not for + * continuous elements. + */ + template + void + do_copy_from_mg_add (const MGDoFHandler &mg_dof, + OutVector &dst, + const MGLevelObject > &src, + const unsigned int offset) const; /** * Implementation of the @@ -712,10 +819,11 @@ class MGTransferSelect : public MGTransferBase >, */ template void - copy_to_mg (const MGDoFHandler &mg_dof, - MGLevelObject > &dst, - const InVector &src, - const is_1d &) const; + do_copy_to_mg (const MGDoFHandler &mg_dof, + MGLevelObject > &dst, + const InVector &src, + const unsigned int offset, + const is_1d &) const; /** * Same for all other space @@ -723,12 +831,25 @@ class MGTransferSelect : public MGTransferBase >, */ template void - copy_to_mg (const MGDoFHandler &mg_dof, - MGLevelObject > &dst, - const InVector &src, - const is_1d &) const; + do_copy_to_mg (const MGDoFHandler &mg_dof, + MGLevelObject > &dst, + const InVector &src, + const unsigned int offset, + const is_1d &) const; + + /** + * Selected component. + */ + unsigned int selected_component; }; +//----------------------------------------------------------------------// +template +inline void +MGTransferSelect::select(const unsigned int component) +{ + selected_component = component; +} #endif diff --git a/deal.II/deal.II/include/multigrid/mg_transfer.templates.h b/deal.II/deal.II/include/multigrid/mg_transfer.templates.h index 173ad81027..53c3b513b5 100644 --- a/deal.II/deal.II/include/multigrid/mg_transfer.templates.h +++ b/deal.II/deal.II/include/multigrid/mg_transfer.templates.h @@ -34,9 +34,10 @@ template template void -MGTransferPrebuilt::copy_to_mg (const MGDoFHandler &mg_dof_handler, - MGLevelObject > &dst, - const InVector &src) const +MGTransferPrebuilt::copy_to_mg ( + const MGDoFHandler &mg_dof_handler, + MGLevelObject > &dst, + const InVector &src) const { // forward to the correct // specialization @@ -47,10 +48,11 @@ MGTransferPrebuilt::copy_to_mg (const MGDoFHandler &mg_dof_h template template void -MGTransferPrebuilt::copy_to_mg (const MGDoFHandler &, - MGLevelObject > &, - const InVector &, - const is_1d &) const +MGTransferPrebuilt::copy_to_mg ( + const MGDoFHandler &, + MGLevelObject > &, + const InVector &, + const is_1d &) const { Assert(false, ExcNotImplemented()); } @@ -59,10 +61,11 @@ MGTransferPrebuilt::copy_to_mg (const MGDoFHandler &, template template void -MGTransferPrebuilt::copy_to_mg (const MGDoFHandler &mg_dof_handler, - MGLevelObject > &dst, - const InVector &src, - const is_1d &) const +MGTransferPrebuilt::copy_to_mg ( + const MGDoFHandler &mg_dof_handler, + MGLevelObject > &dst, + const InVector &src, + const is_1d &) const { // Make src a real finite element function // InVector src = osrc; @@ -159,10 +162,10 @@ MGTransferPrebuilt::copy_to_mg (const MGDoFHandler &mg_dof_h template template void -MGTransferPrebuilt:: -copy_from_mg(const MGDoFHandler &mg_dof_handler, - OutVector &dst, - const MGLevelObject > &src) const +MGTransferPrebuilt::copy_from_mg( + const MGDoFHandler &mg_dof_handler, + OutVector &dst, + const MGLevelObject > &src) const { const unsigned int dofs_per_cell = mg_dof_handler.get_fe().dofs_per_cell; @@ -207,10 +210,10 @@ copy_from_mg(const MGDoFHandler &mg_dof_handler, template template void -MGTransferPrebuilt:: -copy_from_mg_add (const MGDoFHandler &mg_dof_handler, - OutVector &dst, - const MGLevelObject > &src) const +MGTransferPrebuilt::copy_from_mg_add ( + const MGDoFHandler &mg_dof_handler, + OutVector &dst, + const MGLevelObject > &src) const { const unsigned int dofs_per_cell = mg_dof_handler.get_fe().dofs_per_cell; @@ -259,16 +262,33 @@ copy_from_mg_add (const MGDoFHandler &mg_dof_handler, template -template +template void -MGTransferSelect:: -copy_to_mg (const MGDoFHandler &mg_dof_handler, - MGLevelObject > &dst, - const InVector &src) const +MGTransferSelect::copy_to_mg ( + const MGDoFHandler &mg_dof_handler, + MGLevelObject > &dst, + const BlockVector &src) const { // forward to the correct // specialization - copy_to_mg (mg_dof_handler, dst, src, is_1d<(dim==1)>()); + do_copy_to_mg (mg_dof_handler, dst, src, 0, is_1d<(dim==1)>()); +} + + + +template +template +void +MGTransferSelect::copy_to_mg ( + const MGDoFHandler &mg_dof_handler, + MGLevelObject > &dst, + const Vector &src) const +{ + // forward to the correct + // specialization + do_copy_to_mg (mg_dof_handler, dst, src, + component_start[selected_component], + is_1d<(dim==1)>()); } @@ -276,11 +296,12 @@ copy_to_mg (const MGDoFHandler &mg_dof_handler, template template void -MGTransferSelect:: -copy_to_mg (const MGDoFHandler &, - MGLevelObject > &, - const InVector &, - const is_1d &) const +MGTransferSelect::do_copy_to_mg ( + const MGDoFHandler&, + MGLevelObject >&, + const InVector&, + const unsigned int, + const is_1d&) const { Assert(false, ExcNotImplemented()); } @@ -290,14 +311,15 @@ copy_to_mg (const MGDoFHandler &, template template void -MGTransferSelect:: -copy_to_mg (const MGDoFHandler &mg_dof_handler, - MGLevelObject > &dst, - const InVector &osrc, - const is_1d &) const +MGTransferSelect::do_copy_to_mg ( + const MGDoFHandler& mg_dof_handler, + MGLevelObject >& dst, + const InVector& src, + const unsigned int offset, + const is_1d&) const { // Make src a real finite element function - InVector src = osrc; +// InVector src = osrc; // constraints->distribute(src); const FiniteElement& fe = mg_dof_handler.get_fe(); @@ -313,17 +335,13 @@ copy_to_mg (const MGDoFHandler &mg_dof_handler, Assert(sizes.size()==mg_dof_handler.get_tria().n_levels(), ExcMatricesNotBuilt()); - for (unsigned int l=minlevel;l<=maxlevel;++l) - dst[l].reinit(sizes[l][selected]); + + MGTools::reinit_vector(mg_dof_handler, dst, selected, mg_target_component); std::vector global_dof_indices (dofs_per_cell); std::vector level_dof_indices (dofs_per_cell); std::vector level_face_indices (dofs_per_face); - // Start of input vector inside a - // block vector. - const unsigned int start = component_start[selected]; - // traverse the grid top-down // (i.e. starting with the most // refined grid). this way, we can @@ -336,10 +354,6 @@ copy_to_mg (const MGDoFHandler &mg_dof_handler, // already have built. for (int level=maxlevel; level>=static_cast(minlevel); --level) { - // Start of treated component - // for this level - const unsigned int level_start = mg_component_start[level][selected]; - typename MGDoFHandler::active_cell_iterator level_cell = mg_dof_handler.begin_active(level); const typename MGDoFHandler::active_cell_iterator @@ -366,9 +380,16 @@ copy_to_mg (const MGDoFHandler &mg_dof_handler, // into the level-wise one for (unsigned int i=0; i::faces_per_cell; ++face_n) @@ -402,13 +423,68 @@ copy_to_mg (const MGDoFHandler &mg_dof_handler, +template +template +void +MGTransferSelect::copy_from_mg ( + const MGDoFHandler& mg_dof_handler, + BlockVector& dst, + const MGLevelObject >& src) const +{ + do_copy_from_mg (mg_dof_handler, dst, src, 0); +} + + + +template +template +void +MGTransferSelect::copy_from_mg ( + const MGDoFHandler& mg_dof_handler, + Vector& dst, + const MGLevelObject >& src) const +{ + do_copy_from_mg (mg_dof_handler, dst, src, + component_start[selected_component]); +} + + + +template +template +void +MGTransferSelect::copy_from_mg_add ( + const MGDoFHandler& mg_dof_handler, + BlockVector& dst, + const MGLevelObject >& src) const +{ + do_copy_from_mg_add (mg_dof_handler, dst, src, 0); +} + + + +template +template +void +MGTransferSelect::copy_from_mg_add ( + const MGDoFHandler& mg_dof_handler, + Vector& dst, + const MGLevelObject >& src) const +{ + do_copy_from_mg_add (mg_dof_handler, dst, src, + component_start[selected_component]); +} + + + template template void -MGTransferSelect:: -copy_from_mg (const MGDoFHandler &mg_dof_handler, - OutVector &dst, - const MGLevelObject > &src) const +MGTransferSelect::do_copy_from_mg ( + const MGDoFHandler &mg_dof_handler, + OutVector &dst, + const MGLevelObject > &src, + const unsigned int offset) const { const FiniteElement& fe = mg_dof_handler.get_fe(); @@ -422,10 +498,6 @@ copy_from_mg (const MGDoFHandler &mg_dof_handler, const typename MGDoFHandler::active_cell_iterator endc = mg_dof_handler.end(); - // Start of input vector inside a - // block vector. - const unsigned int start = component_start[selected]; - // traverse all cells and copy the // data appropriately to the output // vector @@ -437,10 +509,6 @@ copy_from_mg (const MGDoFHandler &mg_dof_handler, DoFObjectAccessor& global_cell = *level_cell; const unsigned int level = level_cell->level(); - // Start of treated component - // for this level - const unsigned int level_start = mg_component_start[level][selected]; - // get the dof numbers of // this cell for the global // and the level-wise @@ -451,9 +519,17 @@ copy_from_mg (const MGDoFHandler &mg_dof_handler, // copy level-wise data to // global vector for (unsigned int i=0; i &mg_dof_handler, template template void -MGTransferSelect:: -copy_from_mg_add (const MGDoFHandler &mg_dof_handler, - OutVector &dst, - const MGLevelObject > &src) const +MGTransferSelect::do_copy_from_mg_add ( + const MGDoFHandler &mg_dof_handler, + OutVector &dst, + const MGLevelObject > &src, + const unsigned int offset) const { const FiniteElement& fe = mg_dof_handler.get_fe(); @@ -477,10 +554,6 @@ copy_from_mg_add (const MGDoFHandler &mg_dof_handler, std::vector global_dof_indices (dofs_per_cell); std::vector level_dof_indices (dofs_per_cell); - // Start of input vector inside a - // block vector. - const unsigned int start = component_start[selected]; - typename MGDoFHandler::active_cell_iterator level_cell = mg_dof_handler.begin_active(); const typename MGDoFHandler::active_cell_iterator @@ -497,10 +570,6 @@ copy_from_mg_add (const MGDoFHandler &mg_dof_handler, DoFObjectAccessor& global_cell = *level_cell; const unsigned int level = level_cell->level(); - // Start of treated component - // for this level - const unsigned int level_start = mg_component_start[level][selected]; - // get the dof numbers of // this cell for the global // and the level-wise @@ -508,18 +577,25 @@ copy_from_mg_add (const MGDoFHandler &mg_dof_handler, global_cell.get_dof_indices (global_dof_indices); level_cell->get_mg_dof_indices(level_dof_indices); -//TODO:[GK] Probably wrong for continuous elements +//TODO:[GK+WB] Probably wrong for continuous elements // copy level-wise data to // global vector for (unsigned int i=0; iset_zero(dst); +//TODO:[GK+WB] constraints->set_zero(dst); } @@ -531,10 +607,10 @@ copy_from_mg_add (const MGDoFHandler &mg_dof_handler, template template void -MGTransferBlock:: -copy_to_mg (const MGDoFHandler &mg_dof_handler, - MGLevelObject > &dst, - const InVector &src) const +MGTransferBlock::copy_to_mg ( + const MGDoFHandler &mg_dof_handler, + MGLevelObject > &dst, + const InVector &src) const { // forward to the correct // specialization @@ -546,11 +622,11 @@ copy_to_mg (const MGDoFHandler &mg_dof_handler, template template void -MGTransferBlock:: -copy_to_mg (const MGDoFHandler &, - MGLevelObject > &, - const InVector &, - const is_1d &) const +MGTransferBlock::copy_to_mg ( + const MGDoFHandler &, + MGLevelObject > &, + const InVector &, + const is_1d &) const { Assert(false, ExcNotImplemented()); } @@ -560,11 +636,11 @@ copy_to_mg (const MGDoFHandler &, template template void -MGTransferBlock:: -copy_to_mg (const MGDoFHandler &mg_dof_handler, - MGLevelObject > &dst, - const InVector &src, - const is_1d &) const +MGTransferBlock::copy_to_mg ( + const MGDoFHandler &mg_dof_handler, + MGLevelObject > &dst, + const InVector &src, + const is_1d &) const { // Make src a real finite element // function @@ -582,7 +658,7 @@ copy_to_mg (const MGDoFHandler &mg_dof_handler, dst.clear(); - MGTools::reinit_vector(mg_dof_handler, dst, selected, target_component); + MGTools::reinit_vector(mg_dof_handler, dst, selected, mg_target_component); std::vector global_dof_indices (dofs_per_cell); std::vector level_dof_indices (dofs_per_cell); @@ -646,10 +722,10 @@ copy_to_mg (const MGDoFHandler &mg_dof_handler, template template void -MGTransferBlock:: -copy_from_mg (const MGDoFHandler &mg_dof_handler, - OutVector &dst, - const MGLevelObject > &src) const +MGTransferBlock::copy_from_mg ( + const MGDoFHandler &mg_dof_handler, + OutVector &dst, + const MGLevelObject > &src) const { const FiniteElement& fe = mg_dof_handler.get_fe(); @@ -693,10 +769,10 @@ copy_from_mg (const MGDoFHandler &mg_dof_handler, template template void -MGTransferBlock:: -copy_from_mg_add (const MGDoFHandler &mg_dof_handler, - OutVector &dst, - const MGLevelObject > &src) const +MGTransferBlock::copy_from_mg_add ( + const MGDoFHandler &mg_dof_handler, + OutVector &dst, + const MGLevelObject > &src) const { const FiniteElement& fe = mg_dof_handler.get_fe(); const unsigned int dofs_per_cell = fe.dofs_per_cell; diff --git a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc index 42bc9452af..d04abdde0c 100644 --- a/deal.II/deal.II/source/multigrid/mg_dof_tools.cc +++ b/deal.II/deal.II/source/multigrid/mg_dof_tools.cc @@ -388,10 +388,6 @@ MGTools::reinit_vector (const MGDoFHandler& mg_dof, } -//TODO:[GK] Check if this function can be simplified -// It took quite a lot of try and error to make it work with -// target_component, so it is more grown than designed - template void MGTools::reinit_vector (const MGDoFHandler& mg_dof, @@ -418,14 +414,12 @@ MGTools::reinit_vector (const MGDoFHandler& mg_dof, Assert (selected.size() == target_component.size(), ExcDimensionMismatch(selected.size(), target_component.size())); - - std::vector target_selected(selected.size(), false); - for (unsigned int i=0;i > ndofs(mg_dof.get_tria().n_levels(), @@ -450,6 +444,46 @@ MGTools::reinit_vector (const MGDoFHandler& mg_dof, } +template +void +MGTools::reinit_vector (const MGDoFHandler& mg_dof, + MGLevelObject >& v, + const std::vector& selected, + const std::vector& target_component) +{ + const unsigned int ncomp = mg_dof.get_fe().n_components(); + + Assert (selected.size() == target_component.size(), + ExcDimensionMismatch(selected.size(), target_component.size())); + + // Compute the number of blocks needed + if (DEBUG) + { + const unsigned int n_selected + = std::accumulate(selected.begin(), + selected.end(), + 0U); + Assert(n_selected == 1, ExcDimensionMismatch(n_selected, 1)); + } + + unsigned int selected_block = 0; + while (!selected[selected_block]) + ++selected_block; + + std::vector > + ndofs(mg_dof.get_tria().n_levels(), + std::vector(target_component.size())); + + count_dofs_per_component (mg_dof, ndofs, target_component); + + for (unsigned int level=v.get_minlevel(); + level<=v.get_maxlevel();++level) + { + v[level].reinit(ndofs[level][selected_block]); + } +} + + @@ -562,6 +596,17 @@ template void MGTools::reinit_vector ( const std::vector&, const std::vector&); +template void MGTools::reinit_vector ( + const MGDoFHandler&, + MGLevelObject >&, + const std::vector&, + const std::vector&); +template void MGTools::reinit_vector ( + const MGDoFHandler&, + MGLevelObject >&, + const std::vector&, + const std::vector&); + template void MGTools::count_dofs_per_component ( const MGDoFHandler&, diff --git a/deal.II/deal.II/source/multigrid/mg_transfer_block.cc b/deal.II/deal.II/source/multigrid/mg_transfer_block.cc index b10ff85066..0562a86dbe 100644 --- a/deal.II/deal.II/source/multigrid/mg_transfer_block.cc +++ b/deal.II/deal.II/source/multigrid/mg_transfer_block.cc @@ -34,12 +34,16 @@ void MGTransferBlockBase::build_matrices ( const MGDoFHandler& mg_dof, const std::vector& select) { + // Fill target component with + // standard values (identity) if it + // is empty if (target_component.size() == 0) { target_component.resize(mg_dof.get_fe().n_components()); for (unsigned int i=0;i& fe = mg_dof.get_fe(); + // Effective number of components + // is the maximum entry in + // mg_target_component. This + // assumes that the values in that + // vector don't have holes. const unsigned int n_components = - *std::max_element(target_component.begin(), target_component.end()) + 1; + *std::max_element(mg_target_component.begin(), mg_target_component.end()) + 1; const unsigned int dofs_per_cell = fe.dofs_per_cell; const unsigned int n_levels = mg_dof.get_tria().n_levels(); + // Selected refers to the component + // in mg_target_component. selected = select; Assert (selected.size() == fe.n_components(), - ExcDimensionMismatch(selected.size(), fe.n_components())) - + ExcDimensionMismatch(selected.size(), fe.n_components())); + + // Compute the lengths of all blocks sizes.resize(n_levels); - MGTools::count_dofs_per_component(mg_dof, sizes, target_component); + MGTools::count_dofs_per_component(mg_dof, sizes, mg_target_component); + + // Fill some index vectors + // for later use. + mg_component_start = sizes; + // Compute start indices from sizes + for (unsigned int l=0;ladd(dof_indices_child[i], dof_indices_parent[j]); }; @@ -217,7 +276,7 @@ void MGTransferBlockBase::build_matrices ( { const unsigned int icomp = fe.system_to_component_index(i).first; const unsigned int jcomp = fe.system_to_component_index(j).first; - if ((icomp==jcomp) && selected[target_component[icomp]]) + if ((icomp==jcomp) && selected[mg_target_component[icomp]]) prolongation_matrices[level]->set(dof_indices_child[i], dof_indices_parent[j], prolongation(i,j)); @@ -225,31 +284,6 @@ void MGTransferBlockBase::build_matrices ( }; }; }; - // Finally, fill some index vectors - // for later use. - mg_component_start = sizes; - // Compute start indices from sizes - for (unsigned int l=0;l void MGTransferBlock::build_matrices ( const MGDoFHandler &mg_dof, std::vector select, - const std::vector& t_component) + const std::vector& t_component, + const std::vector& mg_t_component) { if (select.size() == 0) select = std::vector (mg_dof.get_fe().n_components(), true); target_component = t_component; + mg_target_component = mg_t_component; MGTransferBlockBase::build_matrices (mg_dof, select); } @@ -287,15 +323,19 @@ template void MGTransferSelect::build_matrices ( const MGDoFHandler &mg_dof, unsigned int select, - const std::vector& t_component) + unsigned int mg_select, + const std::vector& t_component, + const std::vector& mg_t_component) { unsigned int ncomp = mg_dof.get_fe().n_components(); - selected = select; + selected_component = select; std::vector s(ncomp, false); - s[select] = true; - + s[mg_select] = true; +//TODO:[GK] Assert that select and mg_select access the same component +// or figure out how to use different components target_component = t_component; + mg_target_component = mg_t_component; MGTransferBlockBase::build_matrices (mg_dof, s); } @@ -312,7 +352,7 @@ void MGTransferBlock::build_matrices template void MGTransferSelect::build_matrices (const MGDoFHandler &mg_dof, - unsigned int, const std::vector&); + unsigned int, unsigned int, const std::vector&); template void MGTransferBlock::build_matrices @@ -322,7 +362,7 @@ void MGTransferBlock::build_matrices template void MGTransferSelect::build_matrices (const MGDoFHandler &mg_dof, - unsigned int, const std::vector&); + unsigned int, unsigned int, const std::vector&); template void diff --git a/deal.II/deal.II/source/multigrid/multigrid.all_dimensions.cc b/deal.II/deal.II/source/multigrid/multigrid.all_dimensions.cc index 093099414a..70ba223c79 100644 --- a/deal.II/deal.II/source/multigrid/multigrid.all_dimensions.cc +++ b/deal.II/deal.II/source/multigrid/multigrid.all_dimensions.cc @@ -185,7 +185,8 @@ void MGTransferSelect::prolongate ( Assert ((to_level >= 1) && (to_level<=prolongation_matrices.size()), ExcIndexRange (to_level, 1, prolongation_matrices.size()+1)); - prolongation_matrices[to_level-1]->block(selected, selected) + prolongation_matrices[to_level-1]->block(selected_component, + selected_component) .vmult (dst, src); } @@ -199,7 +200,8 @@ void MGTransferSelect::restrict_and_add ( Assert ((from_level >= 1) && (from_level<=prolongation_matrices.size()), ExcIndexRange (from_level, 1, prolongation_matrices.size()+1)); - prolongation_matrices[from_level-1]->block(selected, selected) + prolongation_matrices[from_level-1]->block(selected_component, + selected_component) .Tvmult_add (dst, src); } diff --git a/tests/multigrid/transfer.cc b/tests/multigrid/transfer.cc index 77353a49e2..6eba0bc4ff 100644 --- a/tests/multigrid/transfer.cc +++ b/tests/multigrid/transfer.cc @@ -18,6 +18,7 @@ #include #include #include +#include #include #include #include @@ -25,9 +26,11 @@ #include #include #include +#include #include #include +#include #include @@ -67,6 +70,7 @@ void check_simple(const FiniteElement& fe) template void check_select(const FiniteElement& fe, + unsigned int selected, std::vector target_component, std::vector mg_target_component) { @@ -82,24 +86,24 @@ void check_select(const FiniteElement& fe, for (unsigned int l=0;l > ndofs(mgdof.get_tria().n_levels()); - MGTools::count_dofs_per_component(mgdof, ndofs, mg_target_component); - for (unsigned int l=0;l > mg_ndofs(mgdof.get_tria().n_levels()); + MGTools::count_dofs_per_component(mgdof, mg_ndofs, mg_target_component); + for (unsigned int l=0;l transfer; - transfer.build_matrices(mgdof, 0, - mg_target_component); + transfer.build_matrices(mgdof, selected, selected, + target_component, mg_target_component); - Vector u2(ndofs[2][0]); - Vector u1(ndofs[1][0]); - Vector u0(ndofs[0][0]); + Vector u2(mg_ndofs[2][0]); + Vector u1(mg_ndofs[1][0]); + Vector u0(mg_ndofs[0][0]); u0 = 1; transfer.prolongate(1,u1,u0); @@ -114,6 +118,23 @@ void check_select(const FiniteElement& fe, deallog << "\tu1 " << u1.l2_norm() << "\tu0 " << u0.l2_norm() << std::endl; + + // Check copy to mg and back + std::vector ndofs (target_component.size()); + DoFTools::count_dofs_per_component(mgdof, ndofs, target_component); + BlockVector u; + u.reinit (ndofs); + for (unsigned int i=0;i > v; + v.resize(2,2); + v[2].reinit(mg_ndofs[2][selected]); + + transfer.copy_to_mg(mgdof, v, u); + for (unsigned int i=0; i(0)); // check_simple (FE_DGP<2>(1)); @@ -140,7 +161,14 @@ int main() v3[i] = i/2; } - check_select (FESystem<2>(FE_DGQ<2>(1), 4), v1, v1); - check_select (FESystem<2>(FE_DGQ<2>(1), 4), v2, v2); - check_select (FESystem<2>(FE_DGQ<2>(1), 4), v3, v3); + check_select (FESystem<2>(FE_DGQ<2>(1), 4), 0, v1, v1); + check_select (FESystem<2>(FE_DGQ<2>(1), 4), 1, v1, v1); + check_select (FESystem<2>(FE_DGQ<2>(1), 4), 0, v2, v2); + check_select (FESystem<2>(FE_DGQ<2>(1), 4), 0, v3, v3); + check_select (FESystem<2>(FE_DGQ<2>(1), 4), 1, v3, v3); + check_select (FESystem<2>(FE_DGQ<2>(1), 4), 0, v1, v3); + check_select (FESystem<2>(FE_DGQ<2>(1), 4), 1, v1, v3); } + + + diff --git a/tests/results/i686-pc-linux-gnu+gcc2.95/multigrid/transfer.output b/tests/results/i686-pc-linux-gnu+gcc2.95/multigrid/transfer.output index c9b7473583..76aa1fd83a 100644 --- a/tests/results/i686-pc-linux-gnu+gcc2.95/multigrid/transfer.output +++ b/tests/results/i686-pc-linux-gnu+gcc2.95/multigrid/transfer.output @@ -5,14 +5,76 @@ DEAL::FESystem<2>[FE_DGQ<2>(1)^4] DEAL::Level 0 dofs: 4 4 4 4 DEAL::Level 1 dofs: 16 16 16 16 DEAL::Level 2 dofs: 64 64 64 64 +DEAL:Transfer::Selected 1 0 0 0 +DEAL:Transfer::com_start 0 64 128 192 +DEAL:Transfer::com_start[0] 0 4 8 12 +DEAL:Transfer::com_start[1] 0 16 32 48 +DEAL:Transfer::com_start[2] 0 64 128 192 DEAL::u0 2.00 u1 4.00 u2 8.00 u1 16.0 u0 32.0 +DEAL:: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 +DEAL::FESystem<2>[FE_DGQ<2>(1)^4] +DEAL::Level 0 dofs: 4 4 4 4 +DEAL::Level 1 dofs: 16 16 16 16 +DEAL::Level 2 dofs: 64 64 64 64 +DEAL:Transfer::Selected 0 1 0 0 +DEAL:Transfer::com_start 0 64 128 192 +DEAL:Transfer::com_start[0] 0 4 8 12 +DEAL:Transfer::com_start[1] 0 16 32 48 +DEAL:Transfer::com_start[2] 0 64 128 192 +DEAL::u0 2.00 u1 4.00 u2 8.00 u1 16.0 u0 32.0 +DEAL:: 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 DEAL::FESystem<2>[FE_DGQ<2>(1)^4] DEAL::Level 0 dofs: 16 0 0 0 DEAL::Level 1 dofs: 64 0 0 0 DEAL::Level 2 dofs: 256 0 0 0 +DEAL:Transfer::Selected 1 0 0 0 +DEAL:Transfer::com_start 0 256 256 256 +DEAL:Transfer::com_start[0] 0 16 16 16 +DEAL:Transfer::com_start[1] 0 64 64 64 +DEAL:Transfer::com_start[2] 0 256 256 256 DEAL::u0 4.00 u1 8.00 u2 16.0 u1 32.0 u0 64.0 +DEAL:: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 +DEAL::FESystem<2>[FE_DGQ<2>(1)^4] +DEAL::Level 0 dofs: 8 8 0 0 +DEAL::Level 1 dofs: 32 32 0 0 +DEAL::Level 2 dofs: 128 128 0 0 +DEAL:Transfer::Selected 1 0 0 0 +DEAL:Transfer::com_start 0 128 256 256 +DEAL:Transfer::com_start[0] 0 8 16 16 +DEAL:Transfer::com_start[1] 0 32 64 64 +DEAL:Transfer::com_start[2] 0 128 256 256 +DEAL::u0 2.83 u1 5.66 u2 11.3 u1 22.6 u0 45.3 +DEAL:: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 +DEAL::FESystem<2>[FE_DGQ<2>(1)^4] +DEAL::Level 0 dofs: 8 8 0 0 +DEAL::Level 1 dofs: 32 32 0 0 +DEAL::Level 2 dofs: 128 128 0 0 +DEAL:Transfer::Selected 0 1 0 0 +DEAL:Transfer::com_start 0 128 256 256 +DEAL:Transfer::com_start[0] 0 8 16 16 +DEAL:Transfer::com_start[1] 0 32 64 64 +DEAL:Transfer::com_start[2] 0 128 256 256 +DEAL::u0 2.83 u1 5.66 u2 11.3 u1 22.6 u0 45.3 +DEAL:: 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 +DEAL::FESystem<2>[FE_DGQ<2>(1)^4] +DEAL::Level 0 dofs: 8 8 0 0 +DEAL::Level 1 dofs: 32 32 0 0 +DEAL::Level 2 dofs: 128 128 0 0 +DEAL:Transfer::Selected 1 0 0 0 +DEAL:Transfer::com_start 0 64 128 192 +DEAL:Transfer::com_start[0] 0 8 16 16 +DEAL:Transfer::com_start[1] 0 32 64 64 +DEAL:Transfer::com_start[2] 0 128 256 256 +DEAL::u0 2.83 u1 5.66 u2 11.3 u1 22.6 u0 45.3 +DEAL:: 0 1 2 3 64 65 66 67 4 5 6 7 68 69 70 71 8 9 10 11 72 73 74 75 12 13 14 15 76 77 78 79 16 17 18 19 80 81 82 83 20 21 22 23 84 85 86 87 24 25 26 27 88 89 90 91 28 29 30 31 92 93 94 95 32 33 34 35 96 97 98 99 36 37 38 39 100 101 102 103 40 41 42 43 104 105 106 107 44 45 46 47 108 109 110 111 48 49 50 51 112 113 114 115 52 53 54 55 116 117 118 119 56 57 58 59 120 121 122 123 60 61 62 63 124 125 126 127 DEAL::FESystem<2>[FE_DGQ<2>(1)^4] DEAL::Level 0 dofs: 8 8 0 0 DEAL::Level 1 dofs: 32 32 0 0 DEAL::Level 2 dofs: 128 128 0 0 +DEAL:Transfer::Selected 0 1 0 0 +DEAL:Transfer::com_start 0 64 128 192 +DEAL:Transfer::com_start[0] 0 8 16 16 +DEAL:Transfer::com_start[1] 0 32 64 64 +DEAL:Transfer::com_start[2] 0 128 256 256 DEAL::u0 2.83 u1 5.66 u2 11.3 u1 22.6 u0 45.3 +DEAL:: 128 129 130 131 192 193 194 195 132 133 134 135 196 197 198 199 136 137 138 139 200 201 202 203 140 141 142 143 204 205 206 207 144 145 146 147 208 209 210 211 148 149 150 151 212 213 214 215 152 153 154 155 216 217 218 219 156 157 158 159 220 221 222 223 160 161 162 163 224 225 226 227 164 165 166 167 228 229 230 231 168 169 170 171 232 233 234 235 172 173 174 175 236 237 238 239 176 177 178 179 240 241 242 243 180 181 182 183 244 245 246 247 184 185 186 187 248 249 250 251 188 189 190 191 252 253 254 255