]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
point_value now also works with hp::dofhandler.
authorgoll <goll@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 26 Jan 2012 09:47:37 +0000 (09:47 +0000)
committergoll <goll@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 26 Jan 2012 09:47:37 +0000 (09:47 +0000)
git-svn-id: https://svn.dealii.org/trunk@24932 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/numerics/vectors.h
deal.II/include/deal.II/numerics/vectors.templates.h

index 38b8a6cb8915210ca1767e21c56af19932d064e2..0bf46f634bb8e6b3b5c6579fbfea12d0783b10ee 100644 (file)
@@ -1874,6 +1874,18 @@ namespace VectorTools
               const Point<spacedim>      &point,
               Vector<double>        &value);
 
+  /**
+   * Same as above for hp.
+   */
+
+  template <int dim, class InVector, int spacedim>
+  void
+  point_value (const hp::DoFHandler<dim,spacedim> &dof,
+         const InVector        &fe_function,
+         const Point<spacedim>      &point,
+         Vector<double>        &value);
+
+
                                   /**
                                    * Evaluate a scalar finite
                                    * element function defined by
@@ -1898,6 +1910,16 @@ namespace VectorTools
               const InVector        &fe_function,
               const Point<spacedim>      &point);
 
+  /**
+   * Same as above for hp.
+   */
+
+  template <int dim, class InVector, int spacedim>
+  double
+  point_value (const hp::DoFHandler<dim,spacedim> &dof,
+         const InVector        &fe_function,
+         const Point<spacedim>      &point);
+
                                   /**
                                    * Evaluate a possibly
                                    * vector-valued finite element
@@ -1921,6 +1943,18 @@ namespace VectorTools
               const Point<spacedim>      &point,
               Vector<double>        &value);
 
+  /**
+   * Same as above for hp.
+   */
+
+  template <int dim, class InVector, int spacedim>
+  void
+  point_value (const hp::MappingCollection<dim, spacedim>    &mapping,
+         const hp::DoFHandler<dim,spacedim> &dof,
+         const InVector        &fe_function,
+         const Point<spacedim>      &point,
+         Vector<double>        &value);
+
                                   /**
                                    * Evaluate a scalar finite
                                    * element function defined by
@@ -1941,6 +1975,17 @@ namespace VectorTools
               const InVector        &fe_function,
               const Point<spacedim>      &point);
 
+  /**
+   * Same as above for hp.
+   */
+
+  template <int dim, class InVector, int spacedim>
+  double
+  point_value (const hp::MappingCollection<dim,spacedim>    &mapping,
+         const hp::DoFHandler<dim,spacedim> &dof,
+         const InVector        &fe_function,
+         const Point<spacedim>      &point);
+
                                   //@}
                                   /**
                                    * Mean value operations
index ff66a2d8c35457e7a03ca108005dd39b4f4abdf1..578c21b52b5a8a02242e9cebb7ef6836aa95b53a 100644 (file)
@@ -5359,6 +5359,20 @@ namespace VectorTools
   }
 
 
+  template <int dim, class InVector, int spacedim>
+  void
+  point_value (const hp::DoFHandler<dim,spacedim> &dof,
+         const InVector        &fe_function,
+         const Point<spacedim>      &point,
+         Vector<double>        &value)
+  {
+      point_value(hp::StaticMappingQ1<dim,spacedim>::mapping_collection,
+            dof,
+            fe_function,
+            point,
+            value);
+  }
+
 
   template <int dim, class InVector, int spacedim>
   double
@@ -5372,6 +5386,20 @@ namespace VectorTools
                        point);
   }
 
+
+  template <int dim, class InVector, int spacedim>
+  double
+  point_value (const hp::DoFHandler<dim,spacedim> &dof,
+         const InVector        &fe_function,
+         const Point<spacedim>      &point)
+  {
+      return point_value(hp::StaticMappingQ1<dim,spacedim>::mapping_collection,
+          dof,
+          fe_function,
+          point);
+  }
+
+
   template <int dim, class InVector, int spacedim>
   void
   point_value (const Mapping<dim, spacedim>    &mapping,
@@ -5410,6 +5438,42 @@ namespace VectorTools
   }
 
 
+  template <int dim, class InVector, int spacedim>
+   void
+   point_value (const hp::MappingCollection<dim, spacedim>    &mapping,
+       const hp::DoFHandler<dim,spacedim> &dof,
+       const InVector        &fe_function,
+       const Point<spacedim>      &point,
+       Vector<double>        &value)
+   {
+     const hp::FECollection<dim, spacedim>& fe = dof.get_fe();
+
+     Assert(value.size() == fe.n_components(),
+     ExcDimensionMismatch(value.size(), fe.n_components()));
+
+                                      // first find the cell in which this point
+                                      // is, initialize a quadrature rule with
+                                      // it, and then a FEValues object
+     const std::pair<typename hp::DoFHandler<dim,spacedim>::active_cell_iterator, Point<spacedim> >
+     cell_point =  GridTools::find_active_cell_around_point (mapping, dof, point);
+
+     Assert(GeometryInfo<dim>::distance_to_unit_cell(cell_point.second) < 1e-10,
+         ExcInternalError());
+
+     const Quadrature<dim>
+     quadrature (GeometryInfo<dim>::project_to_unit_cell(cell_point.second));
+     hp::FEValues<dim, spacedim> hp_fe_values(mapping, fe, hp::QCollection<dim>(quadrature), update_values);
+     hp_fe_values.reinit(cell_point.first);
+     const FEValues<dim, spacedim>& fe_values = hp_fe_values.get_present_fe_values();
+
+                                      // then use this to get at the values of
+                                      // the given fe_function at this point
+     std::vector<Vector<double> > u_value(1, Vector<double> (fe.n_components()));
+     fe_values.get_function_values(fe_function, u_value);
+
+     value = u_value[0];
+   }
+
 
   template <int dim, class InVector, int spacedim>
   double
@@ -5446,6 +5510,42 @@ namespace VectorTools
   }
 
 
+  template <int dim, class InVector, int spacedim>
+  double
+  point_value (const hp::MappingCollection<dim, spacedim>    &mapping,
+      const hp::DoFHandler<dim,spacedim> &dof,
+      const InVector        &fe_function,
+      const Point<spacedim>      &point)
+  {
+      const hp::FECollection<dim, spacedim>& fe = dof.get_fe();
+
+      Assert(fe.n_components() == 1,
+          ExcMessage ("Finite element is not scalar as is necessary for this function"));
+
+      // first find the cell in which this point
+      // is, initialize a quadrature rule with
+      // it, and then a FEValues object
+      const std::pair<typename hp::DoFHandler<dim,spacedim>::active_cell_iterator, Point<spacedim> >
+      cell_point =  GridTools::find_active_cell_around_point (mapping, dof, point);
+
+      Assert(GeometryInfo<dim>::distance_to_unit_cell(cell_point.second) < 1e-10,
+          ExcInternalError());
+
+      const Quadrature<dim>
+      quadrature (GeometryInfo<dim>::project_to_unit_cell(cell_point.second));
+      hp::FEValues<dim, spacedim> hp_fe_values(mapping, fe, hp::QCollection<dim>(quadrature), update_values);
+      hp_fe_values.reinit(cell_point.first);
+      const FEValues<dim, spacedim>& fe_values = hp_fe_values.get_present_fe_values();
+
+      // then use this to get at the values of
+      // the given fe_function at this point
+      std::vector<double> u_value(1);
+      fe_values.get_function_values(fe_function, u_value);
+
+      return u_value[0];
+  }
+
+
   template <int dim, class InVector, int spacedim>
   double
   compute_mean_value (const Mapping<dim, spacedim>    &mapping,

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.