map_support_points_to_dofs (const Mapping<dim> &mapping,
const DH<dim> &dof_handler,
std::map<Point<dim>, unsigned int, Comp> &point_to_index_map);
+
+ /**
+ * Map a coupling table from the
+ * user friendly organization by
+ * components to the organization
+ * by blocks. Specializations of
+ * this function for DoFHandler
+ * and hp::DoFHandler are
+ * required due to the different
+ * results of their finite
+ * element access.
+ *
+ * The return vector will be
+ * initialized to the correct
+ * length inside this function.
+ */
+ template <int dim>
+ void convert_couplings_to_blocks (const hp::DoFHandler<dim>& dof_handler,
+ const Table<2, Coupling>& table_by_component,
+ std::vector<Table<2,Coupling> >& tables_by_block);
+ /**
+ * Map a coupling table from the
+ * user friendly organization by
+ * components to the organization
+ * by blocks. Specializations of
+ * this function for DoFHandler
+ * and hp::DoFHandler are
+ * required due to the different
+ * results of their finite
+ * element access.
+ *
+ * The return vector will be
+ * initialized to the correct
+ * length inside this function.
+ */
+ template <int dim>
+ void convert_couplings_to_blocks (const DoFHandler<dim>& dof_handler,
+ const Table<2, Coupling>& table_by_component,
+ std::vector<Table<2,Coupling> >& tables_by_block);
/**
* Exception
};
+/**
+ * Operator computing the maximum coupling out of two.
+ *
+ * @relates DoFTools
+ */
+DoFTools::Coupling operator |= (const DoFTools::Coupling& c1,
+ const DoFTools::Coupling& c2)
+{
+ if (c1 == DoFTools::always || c2 == DoFTools::always)
+ return DoFTools::always;
+ if (c1 == DoFTools::nonzero || c2 == DoFTools::nonzero)
+ return DoFTools::nonzero;
+ return DoFTools::none;
+}
+
+
+/**
+ * Operator computing the maximum coupling out of two.
+ *
+ * @relates DoFTools
+ */
+DoFTools::Coupling operator | (const DoFTools::Coupling& c1,
+ const DoFTools::Coupling& c2)
+{
+ if (c1 == DoFTools::always || c2 == DoFTools::always)
+ return DoFTools::always;
+ if (c1 == DoFTools::nonzero || c2 == DoFTools::nonzero)
+ return DoFTools::nonzero;
+ return DoFTools::none;
+}
// ---------------------- inline and template functions --------------------
}
+template<int dim>
+void
+DoFTools::convert_couplings_to_blocks (
+ const DoFHandler<dim>& dof_handler,
+ const Table<2, Coupling>& table,
+ std::vector<Table<2,Coupling> >& tables_by_block)
+{
+ const FiniteElement<dim>& fe = dof_handler.get_fe();
+ const unsigned int nb = fe.n_blocks();
+
+ tables_by_block.resize(1);
+ tables_by_block[0].reinit(nb, nb);
+ tables_by_block[0].fill(none);
+
+ for (unsigned int i=0;i<fe.n_components();++i)
+ {
+ const unsigned int ib = fe.component_to_block(i);
+ for (unsigned int j=0;j<fe.n_components();++j)
+ {
+ const unsigned int jb = fe.component_to_block(j);
+ tables_by_block[0](ib,jb) |= table(i,j);
+ }
+ }
+}
+
+
+template<int dim>
+void
+DoFTools::convert_couplings_to_blocks (
+ const hp::DoFHandler<dim>& dof_handler,
+ const Table<2, Coupling>& table,
+ std::vector<Table<2,Coupling> >& tables_by_block)
+{
+ const hp::FECollection<dim>& coll = dof_handler.get_fe();
+ tables_by_block.resize(coll.n_finite_elements());
+
+ for (unsigned int f=0;f<coll.n_finite_elements();++n)
+ {
+ const FiniteElement<dim>& fe = coll.ge_fe(f);
+
+ const unsigned int nb = fe.n_blocks();
+ tables_by_block[f].reinit(nb, nb);
+ tables_by_block[f].fill(none);
+ for (unsigned int i=0;i<fe.n_components();++i)
+ {
+ const unsigned int ib = fe.component_to_block(i);
+ for (unsigned int j=0;j<fe.n_components();++j)
+ {
+ const unsigned int jb = fe.component_to_block(j);
+ tables_by_block[f](ib,jb) |= table(i,j);
+ }
+ }
+ }
+}
+
+
// explicit instantiations
template void
DoFTools::compute_row_length_vector(
template
void
DoFTools::map_dofs_to_support_points<deal_II_dimension>
-(const Mapping<deal_II_dimension> &mapping,
- const DoFHandler<deal_II_dimension> &dof_handler,
- std::vector<Point<deal_II_dimension> > &support_points);
+(const Mapping<deal_II_dimension>&,
+ const DoFHandler<deal_II_dimension>&,
+ std::vector<Point<deal_II_dimension> >&);
template
void
DoFTools::compute_dof_couplings<deal_II_dimension>(
- Table<2,Coupling>& dof_couplings,
- const Table<2,Coupling>& component_couplings,
- const FiniteElement<deal_II_dimension>& fe);
+ Table<2,Coupling>&, const Table<2,Coupling>&,
+ const FiniteElement<deal_II_dimension>&);
+
+template
+void
+DoFTools::convert_couplings_to_blocks (
+ const DoFHandler<deal_II_dimension>&, const Table<2, Coupling>&,
+ std::vector<Table<2,Coupling> >&);
+
+template
+void
+DoFTools::convert_couplings_to_blocks (
+ const hp::DoFHandler<deal_II_dimension>&, const Table<2, Coupling>&,
+ std::vector<Table<2,Coupling> >&);
+