]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add FETools::convert_generalized_support_point_values_to_dof_values
authorMatthias Maier <tamiko@43-1.org>
Thu, 7 Sep 2017 22:14:48 +0000 (17:14 -0500)
committerMatthias Maier <tamiko@43-1.org>
Sun, 17 Sep 2017 20:26:40 +0000 (15:26 -0500)
include/deal.II/fe/fe_tools.h
include/deal.II/fe/fe_tools.templates.h
source/fe/fe_tools.inst.in

index 202fec1088bab360f5da760fdee80eebe3e0a177..1d28e058bcee77ea2fc4c8fa17605f5082f1c11f 100644 (file)
@@ -557,13 +557,31 @@ namespace FETools
    * information.
    */
   template <int dim, int spacedim>
-  void
-  compute_projection_from_face_quadrature_points_matrix (const FiniteElement<dim, spacedim> &fe,
-                                                         const Quadrature<dim-1>    &lhs_quadrature,
-                                                         const Quadrature<dim-1>    &rhs_quadrature,
-                                                         const typename DoFHandler<dim, spacedim>::active_cell_iterator &cell,
-                                                         const unsigned int          face,
-                                                         FullMatrix<double>         &X);
+  void compute_projection_from_face_quadrature_points_matrix(
+    const FiniteElement<dim, spacedim> &fe,
+    const Quadrature<dim - 1> &lhs_quadrature,
+    const Quadrature<dim - 1> &rhs_quadrature,
+    const typename DoFHandler<dim, spacedim>::active_cell_iterator &cell,
+    const unsigned int face,
+    FullMatrix<double> &X);
+
+
+
+  /**
+   * Wrapper around
+   * FiniteElement::convert_generalized_support_point_values_to_dof_values()
+   * that works with arbitrary number types.
+   *
+   * If the number type is double this function simply calls the
+   * FiniteElement::convert_generalized_support_point_values_to_dof_values()
+   * directly. Otherwise, temporary (per-thread storage) is used to call to
+   * above function.
+   */
+  template <int dim, int spacedim, typename number>
+  void convert_generalized_support_point_values_to_dof_values(
+    const FiniteElement<dim, spacedim> &finite_element,
+    const std::vector<Vector<number> > &support_point_values,
+    std::vector<number>                &nodal_values);
 
 
 
index b10f0efd72cc78fdefe854fad9e6e1b88cec236e..cd064e4b859a7f265bff03f71c56cfbd3fb2a6fb 100644 (file)
@@ -2910,6 +2910,79 @@ namespace FETools
 
 
 
+  namespace
+  {
+    // internal helper class to partially specialize the number type
+
+    template <int dim, int spacedim, typename number>
+    struct ConvertHelper
+    {
+      static void convert_generalized_support_point_values_to_dof_values(
+        const FiniteElement<dim, spacedim> &finite_element,
+        const std::vector<Vector<number> > &support_point_values,
+        std::vector<number>                &dof_values)
+      {
+        static Threads::ThreadLocalStorage<std::vector<Vector<double> > >
+        double_support_point_values;
+        static Threads::ThreadLocalStorage<std::vector<double> >
+        double_dof_values;
+
+        double_support_point_values.get().resize(support_point_values.size());
+        double_dof_values.get().resize(dof_values.size());
+
+        for (unsigned int i = 0; i < support_point_values.size(); ++i)
+          {
+            double_support_point_values.get()[i].reinit(
+              finite_element.n_components(), false);
+            std::copy(std::begin(support_point_values[i]),
+                      std::end(support_point_values[i]),
+                      std::begin(double_support_point_values.get()[i]));
+          }
+
+        finite_element.convert_generalized_support_point_values_to_dof_values(
+          double_support_point_values.get(), double_dof_values.get());
+
+        std::copy(std::begin(double_dof_values.get()),
+                  std::end(double_dof_values.get()),
+                  std::begin(dof_values));
+      }
+    };
+
+    // In case of double simply call the finite element directly
+    template <int dim, int spacedim>
+    struct ConvertHelper<dim, spacedim, double>
+    {
+      static void convert_generalized_support_point_values_to_dof_values(
+        const FiniteElement<dim, spacedim> &finite_element,
+        const std::vector<Vector<double> > &support_point_values,
+        std::vector<double>                &dof_values)
+      {
+        finite_element.convert_generalized_support_point_values_to_dof_values(
+          support_point_values, dof_values);
+      }
+    };
+
+  } /* anonymous namespace */
+
+
+
+  template <int dim, int  spacedim, typename number>
+  void convert_generalized_support_point_values_to_dof_values(
+    const FiniteElement<dim, spacedim> &finite_element,
+    const std::vector<Vector<number> > &support_point_values,
+    std::vector<number>                &dof_values)
+  {
+    AssertDimension(support_point_values.size(),
+                    finite_element.get_generalized_support_points().size());
+    AssertDimension(dof_values.size(), finite_element.dofs_per_cell);
+
+    ConvertHelper<dim, spacedim, number>::
+    convert_generalized_support_point_values_to_dof_values(
+      finite_element, support_point_values, dof_values);
+  }
+
+
+
   template <int dim>
   void
   hierarchic_to_lexicographic_numbering (unsigned int degree, std::vector<unsigned int> &h2l)
index 53815d8f45a9edc418331ba94ee4f338707c29c4..53b7632fec7dcac664d11515fd526ac9985f5e9b 100644 (file)
@@ -258,3 +258,21 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
 #endif
     \}
 }
+
+
+for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS; number : REAL_SCALARS)
+{
+#if deal_II_dimension <= deal_II_space_dimension
+    namespace FETools
+    \{
+
+    template
+    void
+    convert_generalized_support_point_values_to_dof_values<deal_II_dimension, deal_II_space_dimension, number> (
+        const FiniteElement<deal_II_dimension, deal_II_space_dimension> &,
+        const std::vector<Vector<number> > &,
+        std::vector<number>                &);
+
+    \}
+#endif
+}

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.