]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement DoF-wise downstream numbering.
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 19 Mar 2010 13:57:56 +0000 (13:57 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 19 Mar 2010 13:57:56 +0000 (13:57 +0000)
git-svn-id: https://svn.dealii.org/trunk@20850 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/dofs/dof_renumbering.h
deal.II/deal.II/include/dofs/dof_tools.h
deal.II/deal.II/source/dofs/dof_renumbering.cc
deal.II/doc/news/changes.h

index c88d012540b26b399c5bc096c4a76a651e3c9621..7291f1fcda8012ccc77b417f2b399fe3a40ef2c6 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -72,7 +72,7 @@ DEAL_II_NAMESPACE_OPEN
  * published in 1984), but certainly not with those used in this
  * library, featuring several 10,000 to a few 100,000 elements.
  *
- * 
+ *
  * <h4>Implementation of renumbering schemes</h4>
  *
  * The renumbering algorithms need quite a lot of memory, since they
@@ -92,7 +92,7 @@ DEAL_II_NAMESPACE_OPEN
  * function. This will then give inferior results, since knots in the
  * graph (representing dofs) are not found to be neighbors even if
  * they would be after condensation.
- * 
+ *
  * The renumbering algorithms work on a purely algebraic basis, due to
  * the isomorphism between the graph theoretical groundwork underlying
  * the algorithms and binary matrices (matrices of which the entries
@@ -109,7 +109,7 @@ DEAL_II_NAMESPACE_OPEN
  * deleted. If no index is allowed, the algorithm will search for its
  * own starting point.
  *
- * 
+ *
  * <h4>Results of renumbering</h4>
  *
  * The renumbering schemes mentioned above do not lead to optimal
@@ -364,7 +364,7 @@ DEAL_II_NAMESPACE_OPEN
  * @ingroup dofs
  * @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2004, 2007, 2008
  */
-namespace DoFRenumbering 
+namespace DoFRenumbering
 {
                                   /**
                                    * Direction based comparator for
@@ -383,7 +383,7 @@ namespace DoFRenumbering
                                        */
       CompareDownstream (const Point<dim> &dir)
                      :
-                     dir(dir) 
+                     dir(dir)
        {}
                                       /**
                                        * Return true if c1 less c2.
@@ -393,14 +393,52 @@ namespace DoFRenumbering
          const Point<dim> diff = c2->center() - c1->center();
          return (diff*dir > 0);
        }
-      
+
+    private:
+                                      /**
+                                       * Flow direction.
+                                       */
+      const Point<dim> dir;
+  };
+
+
+                                  /**
+                                   * Point based comparator for downstream
+                                   * directions: it returns @p true if the
+                                   * second point is downstream of the first
+                                   * one with respect to the direction given
+                                   * to the constructor. If the points are
+                                   * the same with respect to the downstream
+                                   * direction, the point with the lower DoF
+                                   * number is considered smaller.
+                                   */
+  template <int dim>
+  struct ComparePointwiseDownstream
+  {
+                                      /**
+                                       * Constructor.
+                                       */
+      ComparePointwiseDownstream (const Point<dim> &dir)
+                     :
+                     dir(dir)
+       {}
+                                      /**
+                                       * Return true if c1 less c2.
+                                       */
+      bool operator () (const std::pair<Point<dim>,unsigned int> &c1,
+                       const std::pair<Point<dim>,unsigned int> &c2) const
+        {
+         const Point<dim> diff = c2.first-c1.first;
+         return (diff*dir > 0 || (diff*dir==0 && c1.second<c2.second));
+       }
+
     private:
                                       /**
                                        * Flow direction.
                                        */
       const Point<dim> dir;
   };
-  
+
                                   /**
                                    * A namespace for the
                                    * implementation of some
@@ -447,7 +485,7 @@ namespace DoFRenumbering
                                      * DoFRenumbering namespace.
                                      */
     template <class DH>
-    void 
+    void
     Cuthill_McKee (DH&                              dof_handler,
                   const bool                       reversed_numbering = false,
                   const bool                       use_constraints    = false);
@@ -460,7 +498,7 @@ namespace DoFRenumbering
                                      * the DoFHandler dofs but
                                      * returns the renumbering
                                      * vector.
-                                     */    
+                                     */
     template <class DH>
     void
     compute_Cuthill_McKee (std::vector<unsigned int>& new_dof_indices,
@@ -494,7 +532,7 @@ namespace DoFRenumbering
                                      * step-22.
                                      */
     template <class DH>
-    void 
+    void
     king_ordering (DH&                              dof_handler,
                   const bool                       reversed_numbering = false,
                   const bool                       use_constraints    = false);
@@ -539,7 +577,7 @@ namespace DoFRenumbering
                                      * DoFRenumbering namespace.
                                      */
     template <class DH>
-    void 
+    void
     minimum_degree (DH&                              dof_handler,
                    const bool                       reversed_numbering = false,
                    const bool                       use_constraints    = false);
@@ -559,7 +597,7 @@ namespace DoFRenumbering
                            const bool reversed_numbering = false,
                            const bool use_constraints    = false);
   }
-    
+
                                   /**
                                    * Renumber the degrees of
                                    * freedom according to the
@@ -579,7 +617,7 @@ namespace DoFRenumbering
                                    * DoFRenumbering namespace.
                                    */
   template <class DH>
-  void 
+  void
   Cuthill_McKee (DH&                              dof_handler,
                 const bool                       reversed_numbering = false,
                 const bool                       use_constraints    = false,
@@ -593,7 +631,7 @@ namespace DoFRenumbering
                                    * the DoFHandler dofs but
                                    * returns the renumbering
                                    * vector.
-                                   */    
+                                   */
   template <class DH>
   void
   compute_Cuthill_McKee (std::vector<unsigned int>& new_dof_indices,
@@ -626,7 +664,7 @@ namespace DoFRenumbering
                                    * the different methods.
                                    */
   template <int dim>
-  void 
+  void
   Cuthill_McKee (MGDoFHandler<dim>          &dof_handler,
                 const unsigned int          level,
                 const bool                  reversed_numbering = false,
@@ -681,7 +719,7 @@ namespace DoFRenumbering
                                    * are associated with the first
                                    * vector component to which they
                                    * belong.
-                                   * 
+                                   *
                                    * For finite elements with only
                                    * one component, or a single
                                    * non-primitive base element,
@@ -737,7 +775,7 @@ namespace DoFRenumbering
   void
   component_wise (MGDoFHandler<dim>&               dof_handler,
                  const std::vector<unsigned int>& target_component = std::vector<unsigned int>());
-    
+
                                   /**
                                    * Computes the renumbering
                                    * vector needed by the
@@ -746,7 +784,7 @@ namespace DoFRenumbering
                                    * the renumbering on the
                                    * DoFHandler dofs but returns
                                    * the renumbering vector.
-                                   */    
+                                   */
   template <int dim, class ITERATOR, class ENDITERATOR>
   static unsigned int
   compute_component_wise (std::vector<unsigned int>& new_dof_indices,
@@ -830,7 +868,7 @@ namespace DoFRenumbering
   cell_wise (MGDoFHandler<dim>   &dof_handler,
             const unsigned int   level,
             const std::vector<typename MGDoFHandler<dim>::cell_iterator> &cell_order);
-    
+
                                   /**
                                    * @deprecated Use cell_wise() instead.
                                    *
@@ -844,7 +882,7 @@ namespace DoFRenumbering
   cell_wise_dg (MGDoFHandler<dim>   &dof_handler,
                const unsigned int   level,
                const std::vector<typename MGDoFHandler<dim>::cell_iterator> &cell_order);
-    
+
                                   /**
                                    * Computes the renumbering
                                    * vector needed by the
@@ -882,9 +920,18 @@ namespace DoFRenumbering
 
 
                                   /**
-                                   * Cell-wise downstream numbering
-                                   * with respect to a constant
-                                   * flow direction.
+                                   * Downstream numbering with respect to a
+                                   * constant flow direction. If the
+                                   * additional argument @p
+                                   * dof_wise_renumbering is set to @p false,
+                                   * the numbering is performed cell-wise. If
+                                   * it is set to @p true, the dofs are
+                                   * renumbered for finite elements that are
+                                   * based on support points. Note that in
+                                   * this case, the numbering of points will
+                                   * be unaffected in case two points are
+                                   * located the same with respect to the
+                                   * downstream direction.
                                    *
                                    * The cells are sorted such that
                                    * the centers of higher numbers
@@ -911,7 +958,8 @@ namespace DoFRenumbering
   template <class DH, int dim>
   void
   downstream (DH&               dof_handler,
-             const Point<dim>& direction);
+             const Point<dim>& direction,
+             const bool        dof_wise_renumbering = false);
 
                                   /**
                                    * @deprecated Use downstream() instead.
@@ -930,18 +978,20 @@ namespace DoFRenumbering
                                    */
   template <int dim>
   void
-  downstream (MGDoFHandler<dim>  &dof_handler,
+  downstream (MGDoFHandler<dim> &dof_handler,
              const unsigned int level,
-             const Point<dim> &direction);
+             const Point<dim>  &direction,
+             const bool         dof_wise_renumbering = false);
+
                                   /**
                                    * @deprecated Use downstream()
                                    * instead.
                                    */
   template <int dim>
   void
-  downstream_dg (MGDoFHandler<dim>  &dof_handler,
+  downstream_dg (MGDoFHandler<dim> &dof_handler,
                 const unsigned int level,
-                const Point<dim> &direction);
+                const Point<dim>  &direction);
 
                                   /**
                                    * @deprecated The new function
@@ -979,8 +1029,9 @@ namespace DoFRenumbering
   compute_downstream (std::vector<unsigned int>& new_dof_indices,
                      std::vector<unsigned int>& reverse,
                      const DH&                  dof_handler,
-                     const Point<dim>&          direction);
-    
+                     const Point<dim>&          direction,
+                     const bool                 dof_wise_renumbering);
+
                                   /**
                                    * @deprecated Use
                                    * compute_downstream() instead
@@ -1007,8 +1058,9 @@ namespace DoFRenumbering
                      std::vector<unsigned int>& reverse,
                      const MGDoFHandler<dim>&   dof_handler,
                      const unsigned int         level,
-                     const Point<dim>&          direction);
-    
+                     const Point<dim>&          direction,
+                     const bool                 dof_wise_renumbering);
+
                                   /**
                                    * @deprecated Use
                                    * compute_downstream() instead
@@ -1117,7 +1169,7 @@ namespace DoFRenumbering
                                    * renumbering on the DoFHandler
                                    * dofs but returns the
                                    * renumbering vector.
-                                   */   
+                                   */
   template <class DH>
   void
   compute_random (std::vector<unsigned int> &new_dof_indices,
@@ -1170,12 +1222,12 @@ namespace DoFRenumbering
                                    * renumbering on the @p
                                    * DoFHandler dofs but returns
                                    * the renumbering vector.
-                                   */   
+                                   */
   template <class DH>
   void
   compute_subdomain_wise (std::vector<unsigned int> &new_dof_indices,
                          const DH                  &dof_handler);
-    
+
                                   /**
                                    * Exception
                                    *
index 297c42b1694b2e5eb8f3e14d23ba5903aece6d18..4823207c444289bedb3c2c319c5d97d1b8d9e6f5 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -1735,7 +1735,7 @@ class DoFTools
     static void
     map_dofs_to_support_points (const Mapping<dim,spacedim>       &mapping,
                                const DoFHandler<dim,spacedim>    &dof_handler,
-                               std::vector<Point<spacedim> > &support_points);
+                               std::vector<Point<spacedim> >     &support_points);
 
                                     /**
                                      * This is the opposite function
@@ -1767,7 +1767,7 @@ class DoFTools
     template <class DH, class Comp>
     static void
     map_support_points_to_dofs (const Mapping<DH::dimension, DH::space_dimension> &mapping,
-                               const DH                     &dof_handler,
+                               const DH                                          &dof_handler,
                                std::map<Point<DH::space_dimension>, unsigned int, Comp> &point_to_index_map);
 
                                     /**
@@ -2047,8 +2047,8 @@ DoFTools::make_sparsity_pattern (const DH                              &dof,
 template <class DH, class Comp>
 void
 DoFTools::map_support_points_to_dofs (
-  const Mapping<DH::dimension,DH::space_dimension> &mapping,
-  const DH                     &dof_handler,
+  const Mapping<DH::dimension,DH::space_dimension>         &mapping,
+  const DH                                                 &dof_handler,
   std::map<Point<DH::space_dimension>, unsigned int, Comp> &point_to_index_map)
 {
                                   // let the checking of arguments be
index 5316c36a22eb2b49b649abac096f98d882104216..81d3fbad6e659a8c4e3c025109cea61f8cf57a52 100644 (file)
@@ -2,7 +2,7 @@
 //    $id$
 //    Version: $name$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
 
 #include <base/thread_management.h>
 #include <base/utilities.h>
+#include <base/quadrature_lib.h>
+
 #include <lac/sparsity_pattern.h>
 #include <lac/sparsity_tools.h>
 #include <lac/compressed_simple_sparsity_pattern.h>
 #include <lac/constraint_matrix.h>
+
 #include <dofs/dof_accessor.h>
+#include <dofs/dof_handler.h>
+#include <dofs/dof_tools.h>
+#include <dofs/dof_renumbering.h>
+
 #include <grid/tria_iterator.h>
 #include <grid/tria.h>
-#include <dofs/dof_handler.h>
+
+#include <fe/fe.h>
 #include <hp/dof_handler.h>
 #include <hp/fe_collection.h>
-#include <dofs/dof_tools.h>
-#include <fe/fe.h>
-#include <dofs/dof_renumbering.h>
+#include <hp/fe_values.h>
 
 #include <multigrid/mg_dof_handler.h>
 #include <multigrid/mg_dof_accessor.h>
@@ -1154,11 +1160,13 @@ namespace DoFRenumbering
 
   template <class DH, int dim>
   void
-  downstream (DH& dof, const Point<dim>& direction)
+  downstream (DH& dof, const Point<dim>& direction,
+             const bool dof_wise_renumbering)
   {
     std::vector<unsigned int> renumbering(dof.n_dofs());
     std::vector<unsigned int> reverse(dof.n_dofs());
-    compute_downstream(renumbering, reverse, dof, direction);
+    compute_downstream(renumbering, reverse, dof, direction,
+                      dof_wise_renumbering);
 
     dof.renumber_dofs(renumbering);
   }
@@ -1216,19 +1224,82 @@ namespace DoFRenumbering
     std::vector<unsigned int>& new_indices,
     std::vector<unsigned int>& reverse,
     const DH& dof,
-    const Point<dim>& direction)
+    const Point<dim>& direction,
+    const bool dof_wise_renumbering)
   {
-    std::vector<typename DH::cell_iterator>
-      ordered_cells(dof.get_tria().n_active_cells());
-    const CompareDownstream<typename DH::cell_iterator, dim> comparator(direction);
+    if (dof_wise_renumbering == false)
+      {
+       std::vector<typename DH::cell_iterator>
+         ordered_cells(dof.get_tria().n_active_cells());
+       const CompareDownstream<typename DH::cell_iterator, dim> comparator(direction);
 
-    typename DH::active_cell_iterator begin = dof.begin_active();
-    typename DH::active_cell_iterator end = dof.end();
+       typename DH::active_cell_iterator begin = dof.begin_active();
+       typename DH::active_cell_iterator end = dof.end();
 
-    copy (begin, end, ordered_cells.begin());
-    std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
+       copy (begin, end, ordered_cells.begin());
+       std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
+
+       compute_cell_wise(new_indices, reverse, dof, ordered_cells);
+      }
+    else
+      {
+                               // similar code as for
+                               // DoFTools::map_dofs_to_support_points, but
+                               // need to do this for general DH classes and
+                               // want to be able to sort the result
+                               // (otherwise, could use something like
+                               // DoFTools::map_support_points_to_dofs)
+       const unsigned int n_dofs = dof.n_dofs();
+       std::vector<std::pair<Point<dim>,unsigned int> > support_point_list
+         (n_dofs);
+
+       const hp::FECollection<dim> fe_collection (dof.get_fe ());
+       Assert (fe_collection[0].has_support_points(),
+               DoFTools::ExcFEHasNoSupportPoints());
+       hp::QCollection<dim> quadrature_collection;
+       for (unsigned int comp=0; comp<fe_collection.size(); ++comp)
+         {
+           Assert (fe_collection[comp].has_support_points(),
+                   DoFTools::ExcFEHasNoSupportPoints());
+           quadrature_collection.push_back
+             (Quadrature<DH::dimension> (fe_collection[comp].
+                                         get_unit_support_points()));
+         }
+       hp::FEValues<DH::dimension,DH::space_dimension>
+         hp_fe_values (fe_collection, quadrature_collection,
+                       update_quadrature_points);
+
+       std::vector<bool> already_touched (n_dofs, false);
+
+       std::vector<unsigned int> local_dof_indices;
+       typename DH::active_cell_iterator begin = dof.begin_active();
+       typename DH::active_cell_iterator end = dof.end();
+       for ( ; begin != end; ++begin)
+         {
+           const unsigned int dofs_per_cell = begin->get_fe().dofs_per_cell;
+           local_dof_indices.resize (dofs_per_cell);
+           hp_fe_values.reinit (begin);
+           const FEValues<dim> &fe_values =
+             hp_fe_values.get_present_fe_values ();
+           begin->get_dof_indices(local_dof_indices);
+           const std::vector<Point<DH::space_dimension> > & points
+             = fe_values.get_quadrature_points ();
+           for (unsigned int i=0; i<dofs_per_cell; ++i)
+             if (!already_touched[local_dof_indices[i]])
+               {
+                 support_point_list[local_dof_indices[i]].first = points[i];
+                 support_point_list[local_dof_indices[i]].second =
+                   local_dof_indices[i];
+                 already_touched[local_dof_indices[i]] = true;
+               }
+         }
 
-    compute_cell_wise(new_indices, reverse, dof, ordered_cells);
+       ComparePointwiseDownstream<dim> comparator (direction);
+       std::sort (support_point_list.begin(), support_point_list.end(),
+                  comparator);
+       for (unsigned int i=0; i<n_dofs; ++i)
+         new_indices[support_point_list[i].second] = i;
+      }
   }
 
 
@@ -1236,7 +1307,7 @@ namespace DoFRenumbering
   template <int dim>
   void downstream_dg (MGDoFHandler<dim>& dof,
                      const unsigned int level,
-                     const Point<dim>& direction)
+                     const Point<dim>&  direction)
   {
     std::vector<unsigned int> renumbering(dof.n_dofs(level));
     std::vector<unsigned int> reverse(dof.n_dofs(level));
@@ -1250,11 +1321,13 @@ namespace DoFRenumbering
   template <int dim>
   void downstream (MGDoFHandler<dim>& dof,
                   const unsigned int level,
-                  const Point<dim>& direction)
+                  const Point<dim>&  direction,
+                  const bool         dof_wise_renumbering)
   {
     std::vector<unsigned int> renumbering(dof.n_dofs(level));
     std::vector<unsigned int> reverse(dof.n_dofs(level));
-    compute_downstream(renumbering, reverse, dof, level, direction);
+    compute_downstream(renumbering, reverse, dof, level, direction,
+                      dof_wise_renumbering);
 
     dof.renumber_dofs(level, renumbering);
   }
@@ -1272,7 +1345,8 @@ namespace DoFRenumbering
   {
     std::vector<typename MGDoFHandler<dim>::cell_iterator>
       ordered_cells(dof.get_tria().n_cells(level));
-    const CompareDownstream<typename MGDoFHandler<dim>::cell_iterator, dim> comparator(direction);
+    const CompareDownstream<typename MGDoFHandler<dim>::cell_iterator, dim>
+      comparator(direction);
 
     typename MGDoFHandler<dim>::cell_iterator begin = dof.begin(level);
     typename MGDoFHandler<dim>::cell_iterator end = dof.end(level);
@@ -1292,19 +1366,63 @@ namespace DoFRenumbering
     std::vector<unsigned int>& reverse,
     const MGDoFHandler<dim>& dof,
     const unsigned int level,
-    const Point<dim>& direction)
+    const Point<dim>& direction,
+    const bool dof_wise_renumbering)
   {
-    std::vector<typename MGDoFHandler<dim>::cell_iterator>
-      ordered_cells(dof.get_tria().n_cells(level));
-    const CompareDownstream<typename MGDoFHandler<dim>::cell_iterator, dim> comparator(direction);
+    if (dof_wise_renumbering == false)
+      {
+       std::vector<typename MGDoFHandler<dim>::cell_iterator>
+         ordered_cells(dof.get_tria().n_cells(level));
+       const CompareDownstream<typename MGDoFHandler<dim>::cell_iterator, dim> comparator(direction);
 
-    typename MGDoFHandler<dim>::cell_iterator begin = dof.begin(level);
-    typename MGDoFHandler<dim>::cell_iterator end = dof.end(level);
+       typename MGDoFHandler<dim>::cell_iterator begin = dof.begin(level);
+       typename MGDoFHandler<dim>::cell_iterator end = dof.end(level);
 
-    std::copy (begin, end, ordered_cells.begin());
-    std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
+       std::copy (begin, end, ordered_cells.begin());
+       std::sort (ordered_cells.begin(), ordered_cells.end(), comparator);
 
-    compute_cell_wise(new_indices, reverse, dof, level, ordered_cells);
+       compute_cell_wise(new_indices, reverse, dof, level, ordered_cells);
+      }
+    else
+      {
+       Assert (dof.get_fe().has_support_points(),
+               DoFTools::ExcFEHasNoSupportPoints());
+       const unsigned int n_dofs = dof.n_dofs(level);
+       std::vector<std::pair<Point<dim>,unsigned int> > support_point_list
+         (n_dofs);
+
+       Quadrature<dim>   q_dummy(dof.get_fe().get_unit_support_points());
+       FEValues<dim,dim> fe_values (dof.get_fe(), q_dummy,
+                                    update_quadrature_points);
+
+       std::vector<bool> already_touched (dof.n_dofs(), false);
+
+       const unsigned int dofs_per_cell = dof.get_fe().dofs_per_cell;
+       std::vector<unsigned int> local_dof_indices (dofs_per_cell);
+       typename MGDoFHandler<dim>::cell_iterator begin = dof.begin(level);
+       typename MGDoFHandler<dim>::cell_iterator end = dof.end(level);
+       for ( ; begin != end; ++begin)
+         {
+           begin->get_mg_dof_indices(local_dof_indices);
+           fe_values.reinit (begin);
+           const std::vector<Point<dim> > & points
+             = fe_values.get_quadrature_points ();
+           for (unsigned int i=0; i<dofs_per_cell; ++i)
+             if (!already_touched[local_dof_indices[i]])
+               {
+                 support_point_list[local_dof_indices[i]].first = points[i];
+                 support_point_list[local_dof_indices[i]].second =
+                   local_dof_indices[i];
+                 already_touched[local_dof_indices[i]] = true;
+               }
+         }
+
+       ComparePointwiseDownstream<dim> comparator (direction);
+       std::sort (support_point_list.begin(), support_point_list.end(),
+                  comparator);
+       for (unsigned int i=0; i<n_dofs; ++i)
+         new_indices[support_point_list[i].second] = i;
+      }
   }
 
 
@@ -1693,12 +1811,14 @@ namespace DoFRenumbering
 
   template void
   downstream<DoFHandler<deal_II_dimension> >
-  (DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&);
+  (DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&,
+   const bool);
 
   template void
   compute_downstream<DoFHandler<deal_II_dimension> >
   (std::vector<unsigned int>&,std::vector<unsigned int>&,
-   const DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&);
+   const DoFHandler<deal_II_dimension>&, const Point<deal_II_dimension>&,
+   const bool);
 
   template
   void
@@ -1755,13 +1875,15 @@ namespace DoFRenumbering
   template void
   downstream<hp::DoFHandler<deal_II_dimension> >
   (hp::DoFHandler<deal_II_dimension>&,
-   const Point<deal_II_dimension>&);
+   const Point<deal_II_dimension>&,
+   const bool);
 
   template void
   compute_downstream<hp::DoFHandler<deal_II_dimension> >
   (std::vector<unsigned int>&,std::vector<unsigned int>&,
    const hp::DoFHandler<deal_II_dimension>&,
-   const Point<deal_II_dimension>&);
+   const Point<deal_II_dimension>&,
+   const bool);
 
   template
   void
@@ -1787,7 +1909,8 @@ namespace DoFRenumbering
   template
   void downstream<deal_II_dimension>
   (MGDoFHandler<deal_II_dimension>&, const unsigned int,
-   const Point<deal_II_dimension>&);
+   const Point<deal_II_dimension>&,
+   const bool);
 
   template
   void clockwise_dg<deal_II_dimension>
index 77558246eca9fbd46700b0f20bcdf87fe7153a19..cb36bc9a8747fd9b57c4e46c76374a4e856ee1f7 100644 (file)
@@ -589,12 +589,26 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
+  <li>
+  <p>
+  New: The function DoFRenumbering::downstream has now an additional bool
+  argument. If enabled, the downstream comparison is performed on a DoF
+  basis, as opposed to the cell-based comparison that is used with a false
+  argument.
+  <br>
+  (Martin Kronbichler 2010/03/19)
+  </p>
+  </li>
 
-<li><p> Improved: DoFHandler now has a BlockInfo object, automatically updated after DoFHandler::distribute_dofs() and accessible by DoFHandler::block_info(). This object can be used to initialize block vectors and obliterates the necessity to count dofs per block (or dofs per component in legacy code) in application programs.
-<br>
-(GK 2010/03/18)
-</p>
-</li>
+  <li><p> Improved: DoFHandler now has a BlockInfo object, automatically
+  updated after DoFHandler::distribute_dofs() and accessible by
+  DoFHandler::block_info(). This object can be used to initialize block
+  vectors and obliterates the necessity to count dofs per block (or dofs
+  per component in legacy code) in application programs.
+  <br>
+  (GK 2010/03/18)
+  </p>
+  </li>
 
   <li>
   <p>

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.