]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Replace a quadratic algorithm by a linear one. This greatly accelerates a good number...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 8 Jul 2011 01:36:31 +0000 (01:36 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 8 Jul 2011 01:36:31 +0000 (01:36 +0000)
git-svn-id: https://svn.dealii.org/trunk@23935 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/doc/news/changes.h
deal.II/include/deal.II/grid/grid_reordering_internal.h
deal.II/source/grid/grid_reordering.cc

index 1dbab530f86f31c1c6d0bf37fecf6f6ed427e78d..049f03b4f4a30424cdece02e340f8c547abc0d23 100644 (file)
@@ -216,6 +216,13 @@ should be fixed now.
 <h3>Specific improvements</h3>
 
 <ol>
+<li> Fixed: The 2d grid reordering algorithm that is used by all grid readers had
+a component that was quadratic in its complexity, sometimes leading to cases
+where reading in a mesh in debug mode could take minutes for just a few tens
+of thousands of cells. This has now been fixed.
+<br>
+(Wolfgang Bangerth 2011/07/07)
+
 <li> New: The function DoFTools::count_dofs_per_component now also works
 for objects of type hp::DoFHandler.
 <br>
index c7e5261d194c7d66248433ae28e20871c1047262..4f926379a17025dbce9d161f5673b661d7baf381 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2003, 2004, 2005, 2006 by the deal.II authors
+//    Copyright (C) 2003, 2004, 2005, 2006, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -142,22 +142,6 @@ namespace internal
                                          * data of this object.
                                          */
        CellData<2>  original_cell_data;
-       
-                                        /**
-                                         * Makes an MQuad from the
-                                         * given CellData and MSide
-                                         * list.  Is derived from
-                                         * binary_function to be
-                                         * usable with STL
-                                         * containers.
-                                         *
-                                         * Also assumes that the
-                                         * edges listed present in
-                                         * the CellData are already
-                                         * present in the elist
-                                         * vector.
-                                         */ 
-       struct MakeQuad;
     };
 
 /**
index 779213bb2bfa5800b9bf5238ddf2053efcfea4e1..dbaa9f719001bb95143885d327511658579c045c 100644 (file)
@@ -266,48 +266,74 @@ namespace internal
     }
 
 
-    struct MQuad::MakeQuad : public std::binary_function<CellData<2>,
-                 std::vector<MSide>,
-                 MQuad>
+    namespace
     {
-       MQuad operator()(const CellData<2> &q,
-                        const std::vector<MSide> &elist) const
-         {
-                                            // compute the indices of the
-                                            // four sides that bound this
-                                            // quad
-           unsigned int edges[4] = { numbers::invalid_unsigned_int,
-                                     numbers::invalid_unsigned_int,
-                                     numbers::invalid_unsigned_int,
-                                     numbers::invalid_unsigned_int };
-
-           unsigned int index = 0;
-           for (std::vector<MSide>::const_iterator
-                  edge = elist.begin(); edge != elist.end(); ++edge, ++index)
-             for (unsigned int i=0; i<4; ++i)
-                                                // for each of the four
-                                                // sides, find the place in
-                                                // the list with the
-                                                // corresponding line and
-                                                // store it. shortcut the
-                                                // comparison if it has
-                                                // already been found
-               if ((edges[i] == numbers::invalid_unsigned_int)
-                   &&
-                   (MSide::SideSortLess()(*edge, quadside(q,i)) == false))
-                 edges[i] = index;
-
-           for (unsigned int i=0; i<4; ++i)
-             Assert (edges[i] != numbers::invalid_unsigned_int,
-                     ExcInternalError());
+                                      /**
+                                       * A replacement for std::lower_bound
+                                       * if we know that the input sequence
+                                       * is in fact sorted. In that case, we
+                                       * can do a binary search.
+                                       */
+      template<typename Iterator, typename T, typename Comp>
+      Iterator
+      lower_bound_in_sorted_array (Iterator first, Iterator last,
+                                  const T& val, Comp comp)
+      {
+       unsigned int length = std::distance(first, last);
+       unsigned int half;
+       Iterator middle;
 
-           return MQuad(q.vertices[0],q.vertices[1], q.vertices[2], q.vertices[3],
-                        edges[0], edges[1], edges[2], edges[3],
-                        q);
+       while (length > 0)
+         {
+           half = length >> 1;
+           middle = first;
+           std::advance(middle, half);
+           if (comp(*middle, val))
+             {
+               first = middle;
+               ++first;
+               length = length - half - 1;
+             }
+           else
+             length = half;
          }
+       return first;
+      }
+    
 
-    };
-
+                                      /**
+                                       * Create an MQuad object from the
+                                       * indices of the four vertices by
+                                       * looking up the indices of the four
+                                       * sides.
+                                       */
+      MQuad build_quad_from_vertices(const CellData<2> &q,
+                                    const std::vector<MSide> &elist)
+      {
+                                        // compute the indices of the four
+                                        // sides that bound this quad. note
+                                        // that the incoming list elist is
+                                        // sorted with regard to the
+                                        // MSide::SideSortLess criterion
+       unsigned int edges[4] = { numbers::invalid_unsigned_int,
+                                 numbers::invalid_unsigned_int,
+                                 numbers::invalid_unsigned_int,
+                                 numbers::invalid_unsigned_int };
+
+       for (unsigned int i=0; i<4; ++i)
+         edges[i] = (lower_bound_in_sorted_array (elist.begin(),
+                                                  elist.end(),
+                                                  quadside(q,i),
+                                                  MSide::SideSortLess())
+                     -
+                     elist.begin());
+
+       return MQuad(q.vertices[0],q.vertices[1], q.vertices[2], q.vertices[3],
+                    edges[0], edges[1], edges[2], edges[3],
+                    q);
+      }
+    }
+    
 
 
     void
@@ -324,7 +350,6 @@ namespace internal
     {
                                       //Reserve some space
       sides.reserve(4*inquads.size());
-      mquads.reserve(inquads.size());
 
                                       //Insert all the sides into the side vector
       for (int i = 0;i<4;++i)
@@ -348,11 +373,16 @@ namespace internal
                                       // Swap trick to shrink the
                                       // side vector
       std::vector<MSide>(sides).swap(sides);
-
-                                      //Assigns the correct sides to
-                                      //each quads
-      std::transform(inquads.begin(),inquads.end(), std::back_inserter(mquads),
-                    std::bind2nd(MQuad::MakeQuad(),sides) );
+       
+                                      // Now assign the correct sides to
+                                      // each quads
+      mquads.reserve(inquads.size());
+      std::transform(inquads.begin(),
+                    inquads.end(),
+                    std::back_inserter(mquads),
+                    std_cxx1x::bind(build_quad_from_vertices,
+                                    std_cxx1x::_1,
+                                    std_cxx1x::cref(sides)) );
 
                                       // Assign the quads to their sides also.
       int qctr = 0;

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.