From: Guido Kanschat Date: Fri, 23 Sep 2005 08:10:06 +0000 (+0000) Subject: new function for computing node matrix X-Git-Tag: v8.0.0~13093 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=0486ca4c3a76123cf8a384eb90eaed28baf4f9a9;p=dealii.git new function for computing node matrix git-svn-id: https://svn.dealii.org/trunk@11510 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_tools.h b/deal.II/deal.II/include/fe/fe_tools.h index f11781950c..2a7e69f514 100644 --- a/deal.II/deal.II/include/fe/fe_tools.h +++ b/deal.II/deal.II/include/fe/fe_tools.h @@ -164,7 +164,54 @@ class FETools static void get_projection_matrix(const FiniteElement &fe1, const FiniteElement &fe2, FullMatrix &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. + *
    + * + *
  1. Define the space of shape + * functions using an arbitrary + * basis wj and + * compute the matrix M of + * node functionals + * Ni applied + * to these basis functions. + * + *
  2. Compute the basis + * vj of the + * finite element shape function + * space by applying + * M-1 to the + * basis wj. + *
+ * + * @note The FiniteElement must + * provide generalized support + * points and and interpolation + * functions. + */ + template + static void compute_node_matrix(FullMatrix& M, + const FiniteElement& fe); + /** * Compute the embedding matrices * from a coarse cell to diff --git a/deal.II/deal.II/source/fe/fe_tools.cc b/deal.II/deal.II/source/fe/fe_tools.cc index 961873332f..588c9ca252 100644 --- a/deal.II/deal.II/source/fe/fe_tools.cc +++ b/deal.II/deal.II/source/fe/fe_tools.cc @@ -474,6 +474,49 @@ void FETools::get_projection_matrix (const FiniteElement &fe1, } +template +void +FETools::compute_node_matrix( + FullMatrix& N, + const FiniteElement& 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 >& points = fe.get_generalized_support_points(); + + // We need the values of the + // polynomials in all generalized + // support points. + std::vector > + values (dim, std::vector(points.size())); + + // In this vector, we store the + // result of the interpolation + std::vector 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 void @@ -1675,6 +1718,11 @@ compute_interpolation_to_quadrature_points_matrix (const FiniteElement &fe, /*-------------- Explicit Instantiations -------------------------------*/ +template +void FETools::compute_node_matrix( + FullMatrix&, + const FiniteElement&); + template void FETools::compute_component_wise( const FiniteElement& element,