From 079d09206279db4e274dd4a70de721563b119de9 Mon Sep 17 00:00:00 2001 From: bangerth Date: Sat, 18 Oct 2008 03:10:44 +0000 Subject: [PATCH] Provide an alternative count_dofs_with_subdomain_association function that splits degrees of freedom among the various vector components. git-svn-id: https://svn.dealii.org/trunk@17271 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_tools.h | 22 ++++++++ deal.II/deal.II/source/dofs/dof_tools.cc | 72 +++++++++++++++++++++++- deal.II/doc/news/changes.h | 10 ++++ 3 files changed, 101 insertions(+), 3 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index 2634881624..e2a3962092 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -1160,6 +1160,28 @@ class DoFTools count_dofs_with_subdomain_association (const DH &dof_handler, const unsigned int subdomain); + /** + * Count how many degrees of freedom are + * uniquely associated with the given + * @p subdomain index. + * + * This function does what the previous + * one does except that it splits the + * result among the vector components of + * the finite element in use by the + * DoFHandler object. The last argument + * (which must have a length equal to the + * number of vector components) will + * therefore store how many degrees of + * freedom of each vector component are + * associated with the given subdomain. + */ + template + static void + count_dofs_with_subdomain_association (const DH &dof_handler, + const unsigned int subdomain, + std::vector &n_dofs_on_subdomain); + /** * Count how many degrees of * freedom out of the total diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index fbbac19936..559c480adc 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -3738,9 +3738,8 @@ DoFTools::get_subdomain_association (const DH &dof_handler, template unsigned int -DoFTools::count_dofs_with_subdomain_association ( - const DH &dof_handler, - const unsigned int subdomain) +DoFTools::count_dofs_with_subdomain_association (const DH &dof_handler, + const unsigned int subdomain) { // in debug mode, make sure that there are // some cells at least with this subdomain @@ -3771,6 +3770,55 @@ DoFTools::count_dofs_with_subdomain_association ( +template +void +DoFTools::count_dofs_with_subdomain_association (const DH &dof_handler, + const unsigned int subdomain, + std::vector &n_dofs_on_subdomain) +{ + Assert (n_dofs_on_subdomain.size() == dof_handler.get_fe().n_components(), + ExcDimensionMismatch (n_dofs_on_subdomain.size(), + dof_handler.get_fe().n_components())); + std::fill (n_dofs_on_subdomain.begin(), n_dofs_on_subdomain.end(), 0); + + // in debug mode, make sure that there are + // some cells at least with this subdomain + // id +#ifdef DEBUG + { + bool found = false; + for (typename Triangulation::active_cell_iterator + cell=dof_handler.get_tria().begin_active(); + cell!=dof_handler.get_tria().end(); ++cell) + if (cell->subdomain_id() == subdomain) + { + found = true; + break; + } + Assert (found == true, + ExcMessage ("There are no cells for the given subdomain!")); + } +#endif + + std::vector subdomain_association (dof_handler.n_dofs()); + extract_subdomain_dofs (dof_handler, subdomain, subdomain_association); + + std::vector component_association (dof_handler.n_dofs()); + for (unsigned int c=0; c component_mask (dof_handler.get_fe().n_components(), false); + component_mask[c] = true; + extract_dofs (dof_handler, component_mask, component_association); + + for (unsigned int i=0; i void DoFTools::count_dofs_per_component ( @@ -5788,15 +5836,33 @@ DoFTools::count_dofs_with_subdomain_association > (const DoFHandler &, const unsigned int); template +void +DoFTools::count_dofs_with_subdomain_association > +(const DoFHandler &, + const unsigned int, + std::vector &); +template unsigned int DoFTools::count_dofs_with_subdomain_association > (const hp::DoFHandler &, const unsigned int); template +void +DoFTools::count_dofs_with_subdomain_association > +(const hp::DoFHandler &, + const unsigned int, + std::vector &); +template unsigned int DoFTools::count_dofs_with_subdomain_association > (const MGDoFHandler &, const unsigned int); +template +void +DoFTools::count_dofs_with_subdomain_association > +(const MGDoFHandler &, + const unsigned int, + std::vector &); template diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 00a7296b5f..94005cd216 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -353,6 +353,16 @@ inconvenience this causes.

deal.II

    +
  1. +

    + New: There is now a second DoFTools::count_dofs_with_subdomain_association function that + calculates the number of degrees of freedom associated with a certain subdomain and + splits the result up according to the vector component of each degree of freedom. This + function is needed when splitting block matrices in parallel computations. +
    + (WB 2008/10/07) +

    +
  2. Fixed: The GridOut::write_gnuplot function had a bug that made it output only the -- 2.39.5