std::vector<bool> user_flags;
tria->save_user_flags(user_flags);
const_cast<Triangulation<dim> &>(*tria).clear_user_flags ();
-
- unsigned int next_free_dof = 0;
- active_cell_iterator cell = begin_active(),
- endc = end();
- for (; cell != endc; ++cell)
- next_free_dof = distribute_dofs_on_cell (cell, next_free_dof);
- used_dofs = next_free_dof;
+ /////////////////////////////////
+
+ // Step 1: distribute DoFs on all
+ // active entities
+ {
+ unsigned int next_free_dof = 0;
+ active_cell_iterator cell = begin_active(),
+ endc = end();
+
+ for (; cell != endc; ++cell)
+ next_free_dof = distribute_dofs_on_cell (cell, next_free_dof);
+
+ used_dofs = next_free_dof;
+ }
+
+ /////////////////////////////////
+
+ // Step 2: identify certain dofs
+ // if the finite element tells us
+ // that they should have the same
+ // value
+ std::vector<unsigned int> new_dof_indices (used_dofs);
+
+ // Step 2a: do so for vertices
+ {
+ typedef
+ std::vector<std::pair<unsigned int, unsigned int> > Identities;
+
+ Table<2,boost::shared_ptr<Identities> >
+ vertex_dof_identities;
+
+ // loop over all vertices and
+ // see which one we need to
+ // work on
+ for (unsigned int vertex_index=0; vertex_index<vertex_dofs_offsets.size();
+ ++vertex_index)
+ {
+ const unsigned int n_active_fe_indices
+ = n_active_vertex_fe_indices (vertex_index);
+ if (n_active_fe_indices > 1)
+ {
+ const unsigned int
+ first_fe_index = nth_active_vertex_fe_index (vertex_index, 0);
+
+ // loop over all the
+ // other FEs with which
+ // we want to identify
+ // the DoF indices of
+ // the first FE of
+ for (unsigned int f=1; f<n_active_fe_indices; ++f)
+ {
+ const unsigned int
+ other_fe_index = nth_active_vertex_fe_index (vertex_index, f);
+
+ // make sure the
+ // entry in the
+ // equivalence
+ // table exists
+ if (vertex_dof_identities[first_fe_index][other_fe_index] == 0)
+ vertex_dof_identities[first_fe_index][other_fe_index]
+ =
+ boost::shared_ptr<Identities>
+ (new Identities((*finite_elements)[first_fe_index]
+ .hp_vertex_dof_identities
+ ((*finite_elements)[other_fe_index])));
+
+ // then loop
+ // through the
+ // identities we
+ // have. first get
+ // the global
+ // numbers of the
+ // dofs we want to
+ // identify and
+ // make sure they
+ // are not yet
+ // constrained to
+ // anything else,
+ // except for to
+ // each other
+ Identities &identities
+ = *vertex_dof_identities[first_fe_index][other_fe_index];
+ for (unsigned int i=0; i<identities.size(); ++i)
+ {
+ const unsigned int first_fe_dof_index
+ = get_vertex_dof_index (vertex_index,
+ first_fe_index,
+ identities[i].first);
+ const unsigned int other_fe_dof_index
+ = get_vertex_dof_index (vertex_index,
+ other_fe_index,
+ identities[i].second);
+
+ Assert ((new_dof_indices[first_fe_dof_index] ==
+ deal_II_numbers::invalid_unsigned_int)
+ ||
+ (new_dof_indices[first_fe_dof_index] ==
+ other_fe_index),
+ ExcInternalError());
+
+ new_dof_indices[first_fe_dof_index] = other_fe_dof_index;
+ }
+ }
+ }
+
+ }
+ }
+
+
// finally restore the user flags
const_cast<Triangulation<dim> &>(*tria).load_user_flags(user_flags);
}