]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
DoFTools::map_dofs_to_support_points now works within the hp framekwork.
authorgoll <goll@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 2 Feb 2012 09:25:53 +0000 (09:25 +0000)
committergoll <goll@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 2 Feb 2012 09:25:53 +0000 (09:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@24975 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/examples/step-27/step-27.cc
deal.II/include/deal.II/dofs/dof_tools.h
deal.II/source/dofs/dof_tools.cc
deal.II/source/dofs/dof_tools.inst.in

index d5d7e084b87e15be92def7690fd4d7fd88c21ae4..b63d73834a31cfef61f22d5ed8b433e1d4cdf879 100644 (file)
@@ -81,6 +81,10 @@ enabled due to a missing include file in file
 <h3>Specific improvements</h3>
 
 <ol>
+<li> Improved: DoFTools::map_dofs_to_support_points() now also works within the hp framework.
+<br>
+(Christian Goll, 2012/02/02)
+
 <li> Improved: There is now a constructor for FESystem that allows to
 create collections of finite elements of arbitrary length.
 <br>
index a3c091558d07b14ee1233ce6f469deae310d4b0d..8c8d5ec7b3c709abc144bb9318077dd2df04c2f0 100644 (file)
@@ -454,6 +454,15 @@ namespace Step27
                                        solution,
                                        estimated_error_per_cell);
 
+    std::vector<dealii::Point<2> > support_points(dof_handler.n_dofs());
+
+    DoFTools::map_dofs_to_support_points(hp::StaticMappingQ1<dim>::mapping_collection,
+        dof_handler, support_points);
+
+
+    dealii::Point<dim> p(0.5,0.75);
+    std::cout << VectorTools::point_value(dof_handler, solution, p) << std::endl;
+
     Vector<float> smoothness_indicators (triangulation.n_active_cells());
     estimate_smoothness (smoothness_indicators);
 
index e5202f3d3627f91f5fb6591c5a3993ab4f41842a..d8d3d4e6a1be708f8fc80cfe38bed4b3ef320d98 100644 (file)
@@ -23,6 +23,7 @@
 #include <deal.II/dofs/function_map.h>
 #include <deal.II/dofs/dof_handler.h>
 #include <deal.II/fe/fe.h>
+#include <deal.II/hp/mapping_collection.h>
 
 #include <vector>
 #include <set>
@@ -39,6 +40,7 @@ template <int dim, int spacedim> class DoFHandler;
 namespace hp
 {
   template <int dim, int spacedim> class DoFHandler;
+  template <int dim, int spacedim> class MappingCollection;
 }
 template <int dim, int spacedim> class MGDoFHandler;
 class ConstraintMatrix;
@@ -2192,6 +2194,16 @@ namespace DoFTools
                              const DoFHandler<dim,spacedim>    &dof_handler,
                              std::vector<Point<spacedim> >     &support_points);
 
+            /**
+             * Same as above for the hp case.
+             */
+
+  template <int dim, int spacedim>
+  void
+  map_dofs_to_support_points (const dealii::hp::MappingCollection<dim,spacedim>   &mapping,
+        const hp::DoFHandler<dim,spacedim>    &dof_handler,
+        std::vector<Point<spacedim> > &support_points);
+
                                   /**
                                    * This is the opposite function
                                    * to the one above. It generates
index a071a393d99797d83f7b3fe7e40e21db2a675d6d..5adc4567168fe81cc03f2f8ad8caec4e4c442fc1 100644 (file)
@@ -36,6 +36,8 @@
 #include <deal.II/fe/fe_values.h>
 #include <deal.II/fe/fe_tools.h>
 #include <deal.II/hp/fe_collection.h>
+#include <deal.II/hp/q_collection.h>
+#include <deal.II/hp/fe_values.h>
 #include <deal.II/dofs/dof_tools.h>
 #include <deal.II/numerics/vectors.h>
 
@@ -5688,51 +5690,103 @@ namespace DoFTools
                     dof_handler.n_boundary_dofs (boundary_indicators));
   }
 
+  namespace internal
+  {
+    namespace
+    {
+      template<class DH>
+      void  map_dofs_to_support_points_internal(
+          const hp::MappingCollection<DH::dimension, DH::space_dimension> & mapping,
+          const DH  &dof_handler,
+          std::vector<Point<DH::space_dimension> >  &support_points)
+      {
+          AssertDimension(support_points.size(), dof_handler.n_dofs());
+
+          const unsigned int dim = DH::dimension;
+          const unsigned int spacedim = DH::space_dimension;
+
+          hp::FECollection<dim, spacedim> fe_collection(dof_handler.get_fe());
+          hp::QCollection<dim> q_coll_dummy;
+
+          for (unsigned int fe_index = 0; fe_index < fe_collection.size(); ++fe_index)
+            {
+              // check whether every fe in the collection
+              // has support points
+              Assert(fe_collection[fe_index].has_support_points(),
+                  typename FiniteElement<dim>::ExcFEHasNoSupportPoints());
+              q_coll_dummy.push_back(
+                  Quadrature<dim> (
+                      fe_collection[fe_index].get_unit_support_points()));
+            }
 
+          // now loop over all cells and
+          // enquire the support points on
+          // each of these. we use dummy
+          // quadrature formulas where the
+          // quadrature points are located at
+          // the unit support points to
+          // enquire the location of the
+          // support points in real space
+          //
+          // the weights of the quadrature
+          // rule have been set to invalid values
+          // by the used constructor.
+          hp::FEValues<dim, spacedim> hp_fe_values(mapping, fe_collection,
+              q_coll_dummy, update_quadrature_points);
+          typename DH::active_cell_iterator cell =
+              dof_handler.begin_active(), endc = dof_handler.end();
+
+          std::vector<unsigned int> local_dof_indices;
+          for (; cell != endc; ++cell)
+            {
+              hp_fe_values.reinit(cell);
+              const FEValues<dim, spacedim> &fe_values = hp_fe_values.get_present_fe_values();
+
+              local_dof_indices.resize(cell->get_fe().dofs_per_cell);
+              cell->get_dof_indices(local_dof_indices);
+
+              const std::vector<Point<spacedim> > & points =
+                  fe_values.get_quadrature_points();
+              for (unsigned int i = 0; i < cell->get_fe().dofs_per_cell; ++i)
+                support_points[local_dof_indices[i]] = points[i];
+            };
+      }
+    }
+  }
 
   template <int dim, int spacedim>
   void
   map_dofs_to_support_points (const Mapping<dim,spacedim>       &mapping,
-                             const DoFHandler<dim,spacedim>    &dof_handler,
-                             std::vector<Point<spacedim> > &support_points)
+      const DoFHandler<dim,spacedim>    &dof_handler,
+      std::vector<Point<spacedim> > &support_points)
   {
-    const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
+    AssertDimension(support_points.size(), dof_handler.n_dofs());
 
-                                    // check whether fe has support
-                                    // points
-    Assert (dof_handler.get_fe().has_support_points(),
-           typename FiniteElement<dim>::ExcFEHasNoSupportPoints());
-    AssertDimension (support_points.size(), dof_handler.n_dofs());
+    //Let the internal function do all the work,
+    //just make sure that it gets a MappingCollection
+    const hp::MappingCollection<dim, spacedim> mapping_collection(mapping);
 
-                                    // now loop over all cells and
-                                    // enquire the support points on
-                                    // each of these. use a dummy
-                                    // quadrature formula where the
-                                    // quadrature points are located at
-                                    // the unit support points to
-                                    // enquire the location of the
-                                    // support points in real space
-                                    //
-                                    // the weights of the quadrature
-                                    // rule are set to invalid values
-                                    // by the used constructor.
-    Quadrature<dim> q_dummy(dof_handler.get_fe().get_unit_support_points());
-    FEValues<dim,spacedim> fe_values (mapping, dof_handler.get_fe(),
-                                     q_dummy, update_quadrature_points);
-    typename DoFHandler<dim,spacedim>::active_cell_iterator
-      cell = dof_handler.begin_active(),
-      endc = dof_handler.end();
+    internal::map_dofs_to_support_points_internal(
+        mapping_collection,
+        dof_handler, support_points);
+  }
 
-    std::vector<unsigned int> local_dof_indices (dofs_per_cell);
-    for (; cell!=endc; ++cell)
-      {
-       fe_values.reinit (cell);
-       cell->get_dof_indices (local_dof_indices);
-       const std::vector<Point<spacedim> > & points
-         = fe_values.get_quadrature_points ();
-       for (unsigned int i=0; i<dofs_per_cell; ++i)
-         support_points[local_dof_indices[i]] = points[i];
-      };
+
+  template<int dim, int spacedim>
+  void
+  map_dofs_to_support_points(
+      const hp::MappingCollection<dim, spacedim> &mapping,
+      const hp::DoFHandler<dim, spacedim> &dof_handler,
+      std::vector<Point<spacedim> > &support_points)
+  {
+    AssertDimension(support_points.size(), dof_handler.n_dofs());
+
+    //Let the internal function do all the work,
+    //just make sure that it gets a MappingCollection
+    internal::map_dofs_to_support_points_internal(
+        mapping,
+        dof_handler,
+        support_points);
   }
 
 
index 996ba002a9c558c0574327558755801d926da015..af7cdb324f3421c99d34f447b107f812736f045b 100644 (file)
@@ -618,6 +618,14 @@ DoFTools::map_dofs_to_support_points<deal_II_dimension>
 (const Mapping<deal_II_dimension,deal_II_dimension>&,
  const DoFHandler<deal_II_dimension>&,
  std::vector<Point<deal_II_dimension> >&);
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension>
+(const hp::MappingCollection<deal_II_dimension,deal_II_dimension>&,
+ const hp::DoFHandler<deal_II_dimension>&,
+ std::vector<Point<deal_II_dimension> >&);
 
 #if deal_II_dimension < 3
 
@@ -627,6 +635,13 @@ DoFTools::map_dofs_to_support_points<deal_II_dimension,deal_II_dimension+1>
 (const Mapping<deal_II_dimension,deal_II_dimension+1>&,
  const DoFHandler<deal_II_dimension, deal_II_dimension+1>&,
  std::vector<Point<deal_II_dimension+1> >&);
+template
+void
+DoFTools::map_dofs_to_support_points<deal_II_dimension,deal_II_dimension+1>
+(const hp::MappingCollection<deal_II_dimension,deal_II_dimension+1>&,
+ const hp::DoFHandler<deal_II_dimension, deal_II_dimension+1>&,
+ std::vector<Point<deal_II_dimension+1> >&);
 
 template
 void

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.