From 78232168799c8eaa265f37c11cfbb7549024272f Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Tue, 5 Jan 2010 17:27:07 +0000 Subject: [PATCH] Allow to say which components fit how into blocks. git-svn-id: https://svn.dealii.org/trunk@20297 0785d39b-7218-0410-832d-ea1e28bc413d --- .../deal.II/include/multigrid/mg_transfer.h | 53 ++++++++++++++++--- .../include/multigrid/mg_transfer.templates.h | 14 ++++- 2 files changed, 60 insertions(+), 7 deletions(-) diff --git a/deal.II/deal.II/include/multigrid/mg_transfer.h b/deal.II/deal.II/include/multigrid/mg_transfer.h index 600af6db41..630cd295f7 100644 --- a/deal.II/deal.II/include/multigrid/mg_transfer.h +++ b/deal.II/deal.II/include/multigrid/mg_transfer.h @@ -61,7 +61,7 @@ template class MGDoFHandler; * @author Wolfgang Bangerth, Guido Kanschat, 1999-2004 */ template -class MGTransferPrebuilt : public MGTransferBase +class MGTransferPrebuilt : public MGTransferBase { public: /** @@ -140,12 +140,45 @@ class MGTransferPrebuilt : public MGTransferBase OutVector& dst, const MGLevelObject& src) const; + /** + * If this object operates on + * BlockVector objects, we need + * to describe how the individual + * vector components are mapped + * to the blocks of a vector. For + * example, for a Stokes system, + * we have dim+1 vector + * components for velocity and + * pressure, but we may want to + * use block vectors with only + * two blocks for all velocities + * in one block, and the pressure + * variables in the other. + * + * By default, if this function + * is not called, block vectors + * have as many blocks as the + * finite element has vector + * components. However, this can + * be changed by calling this + * function with an array that + * describes how vector + * components are to be grouped + * into blocks. The meaning of + * the argument is the same as + * the one given to the + * DoFTools::count_dofs_per_component + * function. + */ + void + set_component_to_block_map (const std::vector &map); + /** * Finite element does not * provide prolongation matrices. */ DeclException0(ExcNoProlongation); - + /** * Call @p build_matrices * function first. @@ -156,7 +189,7 @@ class MGTransferPrebuilt : public MGTransferBase * Memory used by this object. */ unsigned int memory_consumption () const; - + private: @@ -179,8 +212,8 @@ class MGTransferPrebuilt : public MGTransferBase * while row indices belong to the * child cell, i.e. the fine level. */ - std::vector > > prolongation_matrices; - + std::vector > > prolongation_matrices; + /** * Mapping for the * copy_to/from_mg-functions. @@ -190,6 +223,14 @@ class MGTransferPrebuilt : public MGTransferBase std::vector > copy_indices; + /** + * The vector that stores what + * has been given to the + * set_component_to_block_map() + * function. + */ + std::vector component_to_block_map; + /** * Fill the two boolean vectors * #dofs_on_refinement_edge and @@ -219,7 +260,7 @@ class MGTransferPrebuilt : public MGTransferBase * * @todo Clean up names */ - std::vector > dofs_on_refinement_boundary; + std::vector > dofs_on_refinement_boundary; /** * The constraints of the global * system. 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 e22a2981d7..23666ab3d6 100644 --- a/deal.II/deal.II/include/multigrid/mg_transfer.templates.h +++ b/deal.II/deal.II/include/multigrid/mg_transfer.templates.h @@ -40,7 +40,7 @@ MGTransferPrebuilt::copy_to_mg ( MGLevelObject& dst, const InVector& src) const { - internal::mg::reinit_vector(mg_dof_handler, dst); + internal::mg::reinit_vector(mg_dof_handler, component_to_block_map, dst); bool first = true; for (unsigned int level=mg_dof_handler.get_tria().n_levels();level != 0;) { @@ -112,6 +112,18 @@ MGTransferPrebuilt::copy_from_mg_add ( dst(i->first) += src[level](i->second); } + + +template +void +MGTransferPrebuilt:: +set_component_to_block_map (const std::vector &map) +{ + component_to_block_map = map; +} + + + template template void MGTransferPrebuilt::find_dofs_on_refinement_edges ( -- 2.39.5