]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added Scalar::get_function_laplacians_from_local_dof_values
authorJean-Paul Pelteret <jppelteret@gmail.com>
Thu, 10 Aug 2017 22:53:56 +0000 (16:53 -0600)
committerJean-Paul Pelteret <jppelteret@gmail.com>
Sat, 12 Aug 2017 17:09:30 +0000 (11:09 -0600)
include/deal.II/fe/fe_values.h
source/fe/fe_values.cc
source/fe/fe_values.inst.in

index 3b6b591c9ac93157108b46691b7ee986439d43e0..e17caaf5cc4f5c1b05fd54d240ac3380504a025e 100644 (file)
@@ -432,6 +432,14 @@ namespace FEValuesViews
     void get_function_laplacians (const InputVector &fe_function,
                                   std::vector<typename ProductType<value_type,typename InputVector::value_type>::type> &laplacians) const;
 
+    /**
+     * @copydoc FEValuesViews::Scalar::get_function_values_from_local_dof_values()
+     */
+    template <class InputVector>
+    void get_function_laplacians_from_local_dof_values (const InputVector &dof_values,
+                                                        std::vector<typename OutputType<typename InputVector::value_type>::laplacian_type> &laplacians) const;
+
+
     /**
      * Return the third derivatives of the selected scalar component of the
      * finite element function characterized by <tt>fe_function</tt> at the
index 496042749c73f1fba91d6b35bc4015cc8e442807..3f67cd9e3dd406dad8594a093143501e00f300c9 100644 (file)
@@ -527,14 +527,15 @@ namespace FEValuesViews
     do_function_laplacians (const ArrayView<Number> &dof_values,
                             const Table<2,dealii::Tensor<2,spacedim> > &shape_hessians,
                             const std::vector<typename Scalar<dim,spacedim>::ShapeFunctionData> &shape_function_data,
-                            std::vector<typename ProductType<Number,double>::type>           &laplacians)
+                            std::vector<typename Scalar<dim,spacedim>::template OutputType<Number>::laplacian_type> &laplacians)
     {
       const unsigned int dofs_per_cell = dof_values.size();
       const unsigned int n_quadrature_points = dofs_per_cell > 0 ?
                                                shape_hessians[0].size() : laplacians.size();
       AssertDimension (laplacians.size(), n_quadrature_points);
 
-      std::fill (laplacians.begin(), laplacians.end(), typename ProductType<Number,double>::type());
+      std::fill (laplacians.begin(), laplacians.end(),
+                 typename Scalar<dim,spacedim>::template OutputType<Number>::laplacian_type());
 
       for (unsigned int shape_function=0;
            shape_function<dofs_per_cell; ++shape_function)
@@ -1420,6 +1421,26 @@ namespace FEValuesViews
 
 
 
+  template <int dim, int spacedim>
+  template <class InputVector>
+  void
+  Scalar<dim,spacedim>::
+  get_function_laplacians_from_local_dof_values(const InputVector &dof_values,
+                                                std::vector<typename OutputType<typename InputVector::value_type>::laplacian_type> &laplacians) const
+  {
+    Assert (fe_values->update_flags & update_hessians,
+            (typename FEValuesBase<dim,spacedim>::ExcAccessToUninitializedField("update_hessians")));
+    Assert (fe_values->present_cell.get() != nullptr,
+            ExcMessage ("FEValues object is not reinit'ed to any cell"));
+    AssertDimension (dof_values.size(), fe_values->dofs_per_cell);
+
+    internal::do_function_laplacians<dim,spacedim>
+    (make_array_view(dof_values.begin(), dof_values.end()),
+     fe_values->finite_element_output.shape_hessians, shape_function_data, laplacians);
+  }
+
+
+
   template <int dim, int spacedim>
   template <class InputVector>
   void
index cc403b24d4f4ac12467871d1b1b7c77a7e81130a..067db65566d3655aa9f380115653783a40452ffc 100644 (file)
@@ -159,6 +159,11 @@ for (VEC : GENERAL_CONTAINER_TYPES;
     void FEValuesViews::Scalar<deal_II_dimension, deal_II_space_dimension>
     ::get_function_hessians_from_local_dof_values<VEC<Number>>
             (const VEC<Number>&, std::vector<typename OutputType<Number>::hessian_type>&) const;
+
+    template
+    void FEValuesViews::Scalar<deal_II_dimension, deal_II_space_dimension>
+    ::get_function_laplacians_from_local_dof_values<VEC<Number>>
+            (const VEC<Number>&, std::vector<typename OutputType<Number>::value_type>&) const;
 #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.