]> https://gitweb.dealii.org/ - dealii.git/commitdiff
DoFTools::map_dofs_to_support_points
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 23 Apr 2001 17:11:30 +0000 (17:11 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 23 Apr 2001 17:11:30 +0000 (17:11 +0000)
git-svn-id: https://svn.dealii.org/trunk@4471 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/source/dofs/dof_tools.cc

index 4a9d0c7e95d45e6613d1ca636621729fa61380df..c4106f24f75ff0cbd7a9f03fcb062dae9c76e807 100644 (file)
@@ -27,7 +27,7 @@ template <int dim> class DoFHandler;
 template <int dim> class MGDoFHandler;
 class ConstraintMatrix;
 template <template <int> class GridClass, int dim> class InterGridMap;
-
+template <int dim> class Mapping;
 
 
 /**
@@ -137,7 +137,7 @@ template <template <int> class GridClass, int dim> class InterGridMap;
  * don't do it here) and not relevant in this context.
  *
  *
- * @author Wolfgang Bangerth and others, 1998, 1999, 2000
+ * @author Wolfgang Bangerth and others, 1998, 1999, 2000, 2001
  */
 class DoFTools
 {
@@ -870,8 +870,35 @@ class DoFTools
     map_dof_to_boundary_indices (const DoFHandler<dim>         &dof_handler,
                                 const std::set<unsigned char> &boundary_indicators,
                                 std::vector<unsigned int>     &mapping);
+
+                                    /**
+                                     * Return a list of support
+                                     * points for all the degrees of
+                                     * freedom handled by this DoF
+                                     * handler object. This function,
+                                     * of course, only works if the
+                                     * finite element object used by
+                                     * the DoF handler object
+                                     * actually provides support
+                                     * points, i.e. no edge elements
+                                     * or the like. Otherwise, an
+                                     * exception is thrown.
+                                     *
+                                     * The given array must have a
+                                     * length of as many elements as
+                                     * there are degrees of freedom.
+                                     */
+    template <int dim>
+    static void
+    map_dofs_to_support_points (const Mapping<dim>       &mapping,
+                               const DoFHandler<dim>    &dof_handler,
+                               std::vector<Point<dim> > &support_points);
+
+                                    /**
+                                     * Exception
+                                     */
+    DeclException0 (ExcFEHasNoSupportPoints);
     
-                                  
                                     /**
                                      * Exception
                                      */
index 74210838d3bb7b7cab30356c9b89e704c0ae01ca..cdb2b675408328ecde9fe1e467ef8eeafbb4f514 100644 (file)
@@ -14,6 +14,7 @@
 
 #include <base/multithread_info.h>
 #include <base/thread_management.h>
+#include <base/quadrature_lib.h>
 #include <grid/tria.h>
 #include <grid/tria_iterator.h>
 #include <grid/intergrid_map.h>
@@ -23,6 +24,7 @@
 #include <multigrid/mg_dof_handler.h>
 #include <multigrid/mg_dof_accessor.h>
 #include <fe/fe.h>
+#include <fe/fe_values.h>
 #include <dofs/dof_tools.h>
 #include <lac/sparsity_pattern.h>
 #include <lac/block_sparsity_pattern.h>
@@ -1861,6 +1863,47 @@ void DoFTools::map_dof_to_boundary_indices (const DoFHandler<dim>         &dof_h
 
 
 
+template <int dim>
+void
+DoFTools::map_dofs_to_support_points (const Mapping<dim>       &mapping,
+                                     const DoFHandler<dim>    &dof_handler,
+                                     std::vector<Point<dim> > &support_points)
+{
+  const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell;
+  
+                                  // check whether fe has support
+                                  // points
+  Assert (dof_handler.get_fe().has_support_points(),
+         ExcFEHasNoSupportPoints());
+  Assert (support_points.size() == dof_handler.n_dofs(),
+         ExcWrongSize (support_points.size(), dof_handler.n_dofs()));
+
+                                  // now loop over all cells and
+                                  // enquire the support points on
+                                  // each of these
+  QMidpoint<dim> q_dummy;
+  FEValues<dim> fe_values (mapping, dof_handler.get_fe(),
+                          q_dummy, update_support_points);
+  typename DoFHandler<dim>::active_cell_iterator
+    cell = dof_handler.begin_active(),
+    endc = dof_handler.end();
+
+  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);
+//TODO:[WB] Uncomment code once FEValues::get_support_points is available. Add a ChangeLog entry then.
+      Assert (false, ExcNotImplemented());
+//        const std::vector<Point<dim> > & support_points
+//     = fe_values.get_support_points ();
+//        for (unsigned int i=0; i<dofs_per_cell; ++i)
+//     support_points[local_dof_indices[i]] = support_points[i];
+    };
+};
+
+
+
 // explicit instantiations
 #if deal_II_dimension > 1
 template void
@@ -1962,3 +2005,11 @@ DoFTools::map_dof_to_boundary_indices (const DoFHandler<deal_II_dimension> &,
                                       std::vector<unsigned int> &);
 
 #endif 
+
+
+
+template
+void
+DoFTools::map_dofs_to_support_points (const Mapping<deal_II_dimension>       &mapping,
+                                     const DoFHandler<deal_II_dimension>    &dof_handler,
+                                     std::vector<Point<deal_II_dimension> > &support_points);

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.