From: Guido Kanschat Date: Sun, 29 Jan 2006 20:21:16 +0000 (+0000) Subject: counting by blocks X-Git-Tag: v8.0.0~12473 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=33652465e7ac0127275a44bcdf0c99883bf0e81d;p=dealii.git counting by blocks git-svn-id: https://svn.dealii.org/trunk@12205 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index b4a3c7f86a..dc7c48027a 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -751,27 +751,29 @@ class DoFTools * Extract the indices of the * degrees of freedom belonging * to certain vector components - * of a vector-valued finite - * element. The bit vector - * @p select defines, which - * components of an - * FESystem are to be - * extracted from the - * DoFHandler @p dof. The - * entries in @p selected_dofs - * corresponding to degrees of - * freedom belonging to these - * components are then flagged - * @p true, while all others are - * set to @p false. + * or blocks (if the last + * argument is true) of + * a vector-valued finite + * element. The bit vector @p + * select defines, which + * components or blocks of an + * FESystem are to be extracted + * from the DoFHandler @p + * dof. The entries in @p + * selected_dofs corresponding to + * degrees of freedom belonging + * to these components are then + * flagged @p true, while all + * others are set to @p false. * - * The size of - * @p component_select shall + * The size of @p select must * equal the number of components - * in the finite element used by - * @p dof. The size of - * @p selected_dofs shall equal - * dof_handler.n_dofs(). Previous + * or blocks in the FiniteElement + * used by @p dof, depending on + * the argument + * blocks. The size of + * @p selected_dofs must equal + * DoFHandler::n_dofs(). Previous * contents of this array are * overwritten. * @@ -781,8 +783,8 @@ class DoFTools * of its shape functions are * non-zero in more than one * vector component (which holds, - * for example, for Nedelec or - * Raviart-Thomas elements), then + * for example, for FE_Nedelec or + * FE_RaviartThomas elements), then * shape functions cannot be * associated with a single * vector component. In this @@ -799,12 +801,13 @@ class DoFTools template static void extract_dofs (const DoFHandler &dof_handler, - const std::vector &component_select, - std::vector &selected_dofs); + const std::vector &select, + std::vector &selected_dofs, + bool blocks = false); /** * Do the same thing as - * @p extract_dofs for one level + * extract_dofs() for one level * of a multi-grid DoF numbering. */ template @@ -812,7 +815,8 @@ class DoFTools extract_level_dofs (const unsigned int level, const MGDoFHandler &dof, const std::vector &select, - std::vector &selected_dofs); + std::vector &selected_dofs, + bool blocks = false); /** * Extract all degrees of freedom diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 0c345d6bd0..c944344a1f 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -1732,12 +1732,25 @@ void DoFTools::distribute_cell_to_dof_vector ( template void -DoFTools::extract_dofs (const DoFHandler &dof, - const std::vector &component_select, - std::vector &selected_dofs) +DoFTools::extract_dofs ( + const DoFHandler &dof, + const std::vector &component_select, + std::vector &selected_dofs, + bool count_by_blocks) { - Assert(component_select.size() == n_components(dof), - ExcDimensionMismatch(component_select.size(), n_components(dof))); + const FiniteElement &fe = dof.get_fe(); + + if (count_by_blocks) + { + Assert(component_select.size() == fe.n_blocks(), + ExcDimensionMismatch(component_select.size(), fe.n_blocks())); + } + else + { + Assert(component_select.size() == n_components(dof), + ExcDimensionMismatch(component_select.size(), n_components(dof))); + } + Assert(selected_dofs.size() == dof.n_dofs(), ExcDimensionMismatch(selected_dofs.size(), dof.n_dofs())); @@ -1759,8 +1772,6 @@ DoFTools::extract_dofs (const DoFHandler &dof, }; - const FiniteElement &fe = dof.get_fe(); - // preset all values by false std::fill_n (selected_dofs.begin(), dof.n_dofs(), false); @@ -1770,50 +1781,54 @@ DoFTools::extract_dofs (const DoFHandler &dof, // something interesting or not std::vector local_selected_dofs (fe.dofs_per_cell, false); for (unsigned int i=0; i &dof, template void -DoFTools::extract_level_dofs(const unsigned int level, - const MGDoFHandler &dof, - const std::vector &component_select, - std::vector &selected_dofs) +DoFTools::extract_level_dofs( + const unsigned int level, + const MGDoFHandler &dof, + const std::vector &component_select, + std::vector &selected_dofs, + bool count_by_blocks) { const FiniteElement& fe = dof.get_fe(); - Assert(component_select.size() == fe.n_components(), - ExcDimensionMismatch(component_select.size(), fe.n_components())); + + if (count_by_blocks) + { + Assert(component_select.size() == fe.n_blocks(), + ExcDimensionMismatch(component_select.size(), fe.n_blocks())); + } + else + { + Assert(component_select.size() == fe.n_components(), + ExcDimensionMismatch(component_select.size(), fe.n_components())); + } + Assert(selected_dofs.size() == dof.n_dofs(level), ExcDimensionMismatch(selected_dofs.size(), dof.n_dofs(level))); @@ -1867,51 +1894,55 @@ DoFTools::extract_level_dofs(const unsigned int level, // something interesting or not std::vector local_selected_dofs (fe.dofs_per_cell, false); for (unsigned int i=0; i indices(fe.dofs_per_cell); @@ -2465,11 +2496,12 @@ DoFTools::count_dofs_per_component ( { void (*fun_ptr) (const DoFHandler &, const std::vector &, - std::vector &) + std::vector &, + bool) = &DoFTools::template extract_dofs; component_select[i][i] = true; threads += Threads::spawn (fun_ptr)(dof_handler, component_select[i], - dofs_in_component[i]); + dofs_in_component[i], false); }; threads.join_all (); @@ -4007,15 +4039,12 @@ DoFTools::distribute_cell_to_dof_vector template void DoFTools::extract_dofs -(const DoFHandler& dof, - const std::vector& component_select, - std::vector& selected_dofs); +(const DoFHandler&, + const std::vector&, std::vector&, bool); template void DoFTools::extract_level_dofs -(const unsigned int level, - const MGDoFHandler& dof, - const std::vector& component_select, - std::vector& selected_dofs); +(const unsigned int level, const MGDoFHandler&, + const std::vector&, std::vector&, bool); #if deal_II_dimension != 1 template