static void get_projection_matrix(const FiniteElement<dim> &fe1,
const FiniteElement<dim> &fe2,
FullMatrix<number> &matrix);
-
+
+ /**
+ * Compute the matrix of nodal
+ * values of a finite element
+ * applied to all its shape
+ * functions.
+ *
+ * This function is supposed to
+ * help building finite elements
+ * from polynomial spaces and
+ * should be called inside the
+ * constructor of an
+ * element. Applied to a
+ * completely initialized finite
+ * element, the result should be
+ * the unit matrix by definition
+ * of the node values.
+ *
+ * Using this matrix allows the
+ * construction of the basis of
+ * shape functions in two steps.
+ * <ol>
+ *
+ * <li>Define the space of shape
+ * functions using an arbitrary
+ * basis <i>w<sub>j</sub></i> and
+ * compute the matrix <i>M</i> of
+ * node functionals
+ * <i>N<sub>i</sub></i> applied
+ * to these basis functions.
+ *
+ * <li>Compute the basis
+ * <i>v<sub>j</sub></i> of the
+ * finite element shape function
+ * space by applying
+ * <i>M<sup>-1</sup></i> to the
+ * basis <i>w<sub>j</sub></i>.
+ * </ol>
+ *
+ * @note The FiniteElement must
+ * provide generalized support
+ * points and and interpolation
+ * functions.
+ */
+ template <int dim>
+ static void compute_node_matrix(FullMatrix<double>& M,
+ const FiniteElement<dim>& fe);
+
/**
* Compute the embedding matrices
* from a coarse cell to
}
+template<int dim>
+void
+FETools::compute_node_matrix(
+ FullMatrix<double>& N,
+ const FiniteElement<dim>& fe)
+{
+ const unsigned int n_dofs = fe.dofs_per_cell;
+ Assert (fe.has_generalized_support_points(), ExcNotInitialized());
+ Assert (N.n()==n_dofs, ExcDimensionMismatch(N.n(), n_dofs));
+ Assert (N.m()==n_dofs, ExcDimensionMismatch(N.m(), n_dofs));
+
+ const std::vector<Point<dim> >& points = fe.get_generalized_support_points();
+
+ // We need the values of the
+ // polynomials in all generalized
+ // support points.
+ std::vector<std::vector<double> >
+ values (dim, std::vector<double>(points.size()));
+
+ // In this vector, we store the
+ // result of the interpolation
+ std::vector<double> local_dofs(n_dofs);
+
+ // One row per shape
+ // function. Remember that these
+ // are the 'raw' shape functions
+ // where the inverse node matrix is
+ // empty. Otherwise, this would
+ // yield identity.
+ for (unsigned int i=0;i<n_dofs;++i)
+ {
+ for (unsigned int k=0;k<values[0].size();++k)
+ for (unsigned int d=0;d<dim;++d)
+ values[d][k] = fe.shape_value_component(i,points[k],d);
+ fe.interpolate(local_dofs, values);
+ // Enter the interpolated dofs
+ // into the matrix
+ for (unsigned int j=0;j<n_dofs;++j)
+ N(j,i) = local_dofs[j];
+ }
+}
+
+
// This function is tested by tests/fe/internals, since it produces the matrices printed there
template<int dim, typename number>
void
/*-------------- Explicit Instantiations -------------------------------*/
+template
+void FETools::compute_node_matrix(
+ FullMatrix<double>&,
+ const FiniteElement<deal_II_dimension>&);
+
template
void FETools::compute_component_wise(
const FiniteElement<deal_II_dimension>& element,