+template <int dim>
+void
+MGTools::count_dofs_per_block (
+ const MGDoFHandler<dim> &dof_handler,
+ std::vector<std::vector<unsigned int> > &result,
+ std::vector<unsigned int> target_block)
+{
+ const FiniteElement<dim>& fe = dof_handler.get_fe();
+ const unsigned int n_blocks = fe.n_blocks();
+ const unsigned int nlevels = dof_handler.get_tria().n_levels();
+
+ Assert (result.size() == nlevels,
+ ExcDimensionMismatch(result.size(), nlevels));
+
+ if (target_block.size() == 0)
+ {
+ target_block.resize(n_blocks);
+ for (unsigned int i=0;i<n_blocks;++i)
+ target_block[i] = i;
+ }
+
+ Assert(target_block.size() == n_blocks,
+ ExcDimensionMismatch(target_block.size(), n_blocks));
+
+ for (unsigned int l=0;l<nlevels;++l)
+ {
+ result[l].resize (n_blocks);
+ std::fill (result[l].begin(),result[l].end(), 0U);
+
+ // special case for only one
+ // block. treat this first
+ // since it does not require any
+ // computations
+ if (n_blocks == 1)
+ result[l][0] = dof_handler.n_dofs(l);
+ else
+ {
+ // otherwise determine the number
+ // of dofs in each block
+ // separately. do so in parallel
+ std::vector<std::vector<bool> >
+ dofs_in_block (n_blocks,
+ std::vector<bool>(dof_handler.n_dofs(l), false));
+ std::vector<std::vector<bool> >
+ block_select (n_blocks, std::vector<bool>(n_blocks, false));
+ Threads::ThreadGroup<> threads;
+ for (unsigned int i=0; i<n_blocks; ++i)
+ {
+ void (*fun_ptr) (const unsigned int level,
+ const MGDoFHandler<dim> &,
+ const std::vector<bool> &,
+ std::vector<bool> &,
+ bool)
+ = &DoFTools::template extract_level_dofs<dim>;
+ block_select[i][i] = true;
+ threads += Threads::spawn (fun_ptr)(l, dof_handler, block_select[i],
+ dofs_in_block[i], true);
+ }
+ threads.join_all();
+
+ // next count what we got
+ for (unsigned int block=0;block<fe.n_blocks();++block)
+ result[l][target_block[block]]
+ += std::count(dofs_in_block[block].begin(),
+ dofs_in_block[block].end(),
+ true);
+ }
+ }
+}
+
+
+
template <int dim>
void
MGTools::count_dofs_per_component (