* $id-I_h$ that is needed for evaluating $(id-I_h)z$ for e.g. the
* dual solution $z$.
*
- * @author Ralf Hartmann, 2000
+ * @author Ralf Hartmann, Guido Kanschat, 2000, 2004
*/
class FETools
{
public:
+ /**
+ * Compute the vector required to
+ * renumber the dofs of a cell by
+ * component. Furthermore,
+ * compute the vector storing the
+ * start indices of each
+ * component in the local block
+ * vector.
+ *
+ * The second vector is organized
+ * such that there is a vector
+ * for each base element
+ * containing the start index for
+ * each component served by this
+ * base element.
+ *
+ * While the first vector is
+ * checked to have the correct
+ * size, the second one is
+ * reinitialized for convenience.
+ */
+ template<int dim>
+ static void compute_component_wise(
+ const FiniteElement<dim>& fe,
+ std::vector<unsigned int>& renumbering,
+ std::vector<std::vector<unsigned int> >& start_indices);
+
/**
* Gives the interpolation matrix
* that interpolates a @p fe1-
}
}
+template<int dim>
+void FETools::compute_component_wise(
+ const FiniteElement<dim>& element,
+ std::vector<unsigned int>& renumbering,
+ std::vector<std::vector<unsigned int> >& comp_start)
+{
+ Assert(renumbering.size() == element.dofs_per_cell,
+ ExcDimensionMismatch(renumbering.size(),
+ element.dofs_per_cell));
+
+ comp_start.resize(element.n_base_elements());
+
+ unsigned int k=0;
+ for (unsigned int i=0;i<comp_start.size();++i)
+ {
+ comp_start[i].resize(element.element_multiplicity(i));
+ const unsigned int increment
+ = element.base_element(i).dofs_per_cell;
+
+ for (unsigned int j=0;j<comp_start[i].size();++j)
+ {
+ comp_start[i][j] = k;
+ k += increment;
+ }
+ }
+
+ // For each index i of the
+ // unstructured cellwise
+ // numbering, renumbering
+ // contains the index of the
+ // cell-block numbering
+ for (unsigned int i=0;i<element.dofs_per_cell;++i)
+ {
+ std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
+ indices = element.system_to_base_index(i);
+ renumbering[i] = comp_start[indices.first.first][indices.first.second]
+ +indices.second;
+ }
+}
+
+
template <int dim, typename number>
void FETools::get_interpolation_matrix (const FiniteElement<dim> &fe1,
/*-------------- Explicit Instantiations -------------------------------*/
+
+template
+void FETools::compute_component_wise(
+ const FiniteElement<deal_II_dimension>& element,
+ std::vector<unsigned int>&, std::vector<std::vector<unsigned int> >&);
template
void FETools::get_interpolation_matrix<deal_II_dimension>
(const FiniteElement<deal_II_dimension> &,