]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new function for computing node matrix
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 23 Sep 2005 08:10:06 +0000 (08:10 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Fri, 23 Sep 2005 08:10:06 +0000 (08:10 +0000)
git-svn-id: https://svn.dealii.org/trunk@11510 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe_tools.h
deal.II/deal.II/source/fe/fe_tools.cc

index f11781950cb4cface8bed72ed2814fcc6321f898..2a7e69f5149f363d74b8664fa7bd1932e09c485a 100644 (file)
@@ -164,7 +164,54 @@ class FETools
     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
index 961873332f84aec5ce8c1df05e6f8dce7275f3f8..588c9ca2522bc05937d52b2283b380bf1bb2b742 100644 (file)
@@ -474,6 +474,49 @@ void FETools::get_projection_matrix (const FiniteElement<dim> &fe1,
 }
 
 
+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
@@ -1675,6 +1718,11 @@ compute_interpolation_to_quadrature_points_matrix (const FiniteElement<dim> &fe,
 /*-------------- 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,

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.