#include <lac/vector.h>
#include <algorithm>
-
+#include <numeric>
// if necessary try to work around a bug in the IBM xlC compiler
#ifdef XLC_WORK_AROUND_STD_BUG
+template <int dim>
+void
+DoFTools::count_dofs_per_component (const DoFHandler<dim> &dof_handler,
+ std::vector<unsigned int> &dofs_per_component)
+{
+ const unsigned int n_components = dof_handler.get_fe().n_components();
+ dofs_per_component.resize (n_components);
+
+ // special case for only one
+ // component. treat this first
+ // since it does not require any
+ // computations
+ if (n_components == 1)
+ {
+ dofs_per_component[0] = dof_handler.n_dofs();
+ return;
+ };
+
+ // otherwise determine the number
+ // of dofs in each component
+ // separately. do so in parallel
+ std::vector<std::vector<bool> > dofs_in_component (n_components,
+ std::vector<bool>(dof_handler.n_dofs(),
+ false));
+ std::vector<std::vector<bool> > component_select (n_components,
+ std::vector<bool>(n_components, false));
+ Threads::ThreadManager thread_manager;
+ for (unsigned int i=0; i<n_components; ++i)
+ {
+ component_select[i][i] = true;
+ Threads::spawn (thread_manager,
+ Threads::encapsulate (&DoFTools::template extract_dofs<dim>)
+ .collect_args (dof_handler, component_select[i],
+ dofs_in_component[i]));
+ };
+ thread_manager.wait();
+
+ // next count what we got
+ for (unsigned int i=0; i<n_components; ++i)
+ dofs_per_component[i] = std::count(dofs_in_component[i].begin(),
+ dofs_in_component[i].end(),
+ true);
+
+ // finally sanity check
+ Assert (std::accumulate (dofs_per_component.begin(), dofs_per_component.end(), 0U)
+ ==
+ dof_handler.n_dofs(),
+ ExcInternalError());
+
+};
+
+
+
template <int dim>
void
DoFTools::compute_intergrid_constraints (const DoFHandler<dim> &coarse_grid,
std::vector<bool> &);
#endif
+
+template
+void
+DoFTools::count_dofs_per_component (const DoFHandler<deal_II_dimension> &dof_handler,
+ std::vector<unsigned int> &dofs_per_component);
+
+
template
void
DoFTools::compute_intergrid_constraints (const DoFHandler<deal_II_dimension> &,
<h3>deal.II</h3>
<ol>
+ <li> <p>
+ New: The <code class="class">DoFTools</code> class now has
+ a function <code class="member">count_dofs_per_component</code>
+ that counts the number of degrees of freedom in each of the
+ components of the finite element, i.e. how many degrees of
+ freedom there are on the global mesh for each variable (field).
+ <br>
+ (WB 2001/02/22)
+ </p>
+
<li> <p>
- New: The <code class="class">CellAccessor</code> now has a function
+ New: The <code class="class">CellAccessor</code> class now has a function
<code class="member">has_boundary_lines</code> that mostly has
the same semantics that <code class="member">at_boundary</code>
has, but also covers the case that a hexahedron in 3d may be at