[Boost-commit] svn:boost r50692 - in trunk: boost/graph boost/graph/detail libs/graph/test

Subject: [Boost-commit] svn:boost r50692 - in trunk: boost/graph boost/graph/detail libs/graph/test
From: jewillco_at_[hidden]
Date: 2009-01-20 16:37:32

Author: jewillco
Date: 2009-01-20 16:37:31 EST (Tue, 20 Jan 2009)
New Revision: 50692
URL: http://svn.boost.org/trac/boost/changeset/50692

Log:
Added d-ary heap implementation, changed Dijkstra performance test to use that rather than binary heap, changed m_decreased member in dijkstra_visitor to a local variable
Added:
   trunk/boost/graph/detail/d_ary_heap.hpp   (contents, props changed)
Text files modified: 
   trunk/boost/graph/dijkstra_shortest_paths.hpp       |    24 +++++++++++++++++-------                
   trunk/libs/graph/test/dijkstra_heap_performance.cpp |     7 ++++++-                                 
   2 files changed, 23 insertions(+), 8 deletions(-)

Added: trunk/boost/graph/detail/d_ary_heap.hpp
==============================================================================
--- (empty file)
+++ trunk/boost/graph/detail/d_ary_heap.hpp	2009-01-20 16:37:31 EST (Tue, 20 Jan 2009)
@@ -0,0 +1,260 @@
+//
+//=======================================================================
+// Copyright 2009 Trustees of Indiana University
+// Authors: Jeremiah J. Willcock, Andrew Lumsdaine
+//
+// Distributed under the Boost Software License, Version 1.0. (See
+// accompanying file LICENSE_1_0.txt or copy at
+// http://www.boost.org/LICENSE_1_0.txt)
+//=======================================================================
+//
+#ifndef BOOST_D_ARY_HEAP_HPP
+#define BOOST_D_ARY_HEAP_HPP
+
+#include <vector>
+#include <cstddef>
+#include <algorithm>
+#include <utility>
+#include <cassert>
+#include <boost/static_assert.hpp>
+#include <boost/property_map.hpp>
+
+namespace boost {
+
+  // Swap two elements in a property map without assuming they model
+  // LvaluePropertyMap -- currently not used
+  template <typename PropMap>
+  inline void property_map_swap(
+         PropMap prop_map,
+         const typename boost::property_traits<PropMap>::key_type& ka, 
+         const typename boost::property_traits<PropMap>::key_type& kb) {
+    typename boost::property_traits<PropMap>::value_type va = get(prop_map, ka);
+    put(prop_map, ka, get(prop_map, kb));
+    put(prop_map, kb, va);
+  }
+
+  // D-ary heap using an indirect compare operator (use identity_property_map
+  // as DistanceMap to get a direct compare operator).  This heap appears to be
+  // commonly used for Dijkstra's algorithm for its good practical performance
+  // on some platforms; asymptotically, it has an O(lg N) decrease-key
+  // operation while that can be done in constant time on a relaxed heap.  The
+  // implementation is mostly based on the binary heap page on Wikipedia and
+  // online sources that state that the operations are the same for d-ary
+  // heaps.  This code is not based on the old Boost d-ary heap code.
+  //
+  // - d_ary_heap_indirect is a model of UpdatableQueue as is needed for
+  //   dijkstra_shortest_paths.
+  //
+  // - Value must model Assignable.
+  // - Arity must be at least 2 (optimal value appears to be 4, both in my and
+  //   third-party experiments).
+  // - IndexInHeapMap must be a ReadWritePropertyMap from Value to
+  //   Container::size_type (to store the index of each stored value within the
+  //   heap for decrease-key aka update).
+  // - DistanceMap must be a ReadablePropertyMap from Value to something
+  //   (typedef'ed as distance_type).
+  // - Compare must be a BinaryPredicate used as a less-than operator on
+  //   distance_type.
+  // - Container must be a random-access, contiguous container (in practice,
+  //   the operations used probably require that it is std::vector<Value>).
+  //
+  template <typename Value,
+	    std::size_t Arity,
+            typename IndexInHeapPropertyMap,
+            typename DistanceMap,
+            typename Compare = std::less<Value>,
+            typename Container = std::vector<Value> >
+  class d_ary_heap_indirect {
+    BOOST_STATIC_ASSERT (Arity >= 2);
+
+    public:
+    typedef typename Container::size_type size_type;
+    typedef Value value_type;
+
+    d_ary_heap_indirect(DistanceMap distance,
+                        IndexInHeapPropertyMap index_in_heap,
+                        const Compare& compare = Compare())
+      : compare(compare), data(), distance(distance),
+        index_in_heap(index_in_heap) {}
+    /* Implicit copy constructor */
+    /* Implicit assignment operator */
+
+    size_type size() const {
+      return data.size();
+    }
+
+    bool empty() const {
+      return data.empty();
+    }
+
+    void push(const Value& v) {
+      size_type index = data.size();
+      data.push_back(v);
+      put(index_in_heap, v, index);
+      preserve_heap_property_up(index);
+      verify_heap();
+    }
+
+    Value& top() {
+      return data[0];
+    }
+
+    const Value& top() const {
+      return data[0];
+    }
+
+    void pop() {
+      data[0] = data.back();
+      put(index_in_heap, data[0], 0);
+      data.pop_back();
+      preserve_heap_property_down();
+      verify_heap();
+    }
+
+    // This function assumes the key has been updated (using an external write
+    // to the distance map or such)
+    // See http://coding.derkeiler.com/Archive/General/comp.theory/2007-05/msg00043.html
+    void update(const Value& v) { /* decrease-key */
+      size_type index = get(index_in_heap, v);
+      preserve_heap_property_up(index);
+      verify_heap();
+    }
+
+    private:
+    Compare compare;
+    Container data;
+    DistanceMap distance;
+    IndexInHeapPropertyMap index_in_heap;
+
+    // The distances being compared using compare and that are stored in the
+    // distance map
+    typedef typename boost::property_traits<DistanceMap>::value_type distance_type;
+
+    // Get the parent of a given node in the heap
+    static size_type parent(size_type index) {
+      return (index - 1) / Arity;
+    }
+
+    // Get the child_idx'th child of a given node; 0 <= child_idx < Arity
+    static size_type child(size_type index, std::size_t child_idx) {
+      return index * Arity + child_idx + 1;
+    }
+
+    // Swap two elements in the heap by index, updating index_in_heap
+    void swap_heap_elements(size_type index_a, size_type index_b) {
+      using std::swap;
+      Value value_a = data[index_a];
+      Value value_b = data[index_b];
+      data[index_a] = value_b;
+      data[index_b] = value_a;
+      put(index_in_heap, value_a, index_b);
+      put(index_in_heap, value_b, index_a);
+    }
+
+    // Emulate the indirect_cmp that is now folded into this heap class
+    bool compare_indirect(const Value& a, const Value& b) const {
+      return compare(get(distance, a), get(distance, b));
+    }
+
+    // Verify that the array forms a heap; commented out by default
+    void verify_heap() const {
+      // This is a very expensive test so it should be disabled even when
+      // NDEBUG is not defined
+#if 0
+      for (size_t i = 1; i < data.size(); ++i) {
+        if (compare_indirect(data[i], data[parent(i)])) {
+          assert (!"Element is smaller than its parent");
+        }
+      }
+#endif
+    }
+
+    // Starting at a node, move up the tree swapping elements to preserve the
+    // heap property
+    void preserve_heap_property_up(size_type index) {
+      size_type orig_index = index;
+      size_type num_levels_moved = 0;
+      // The first loop just saves swaps that need to be done in order to avoid
+      // aliasing issues in its search; there is a second loop that does the
+      // necessary swap operations
+      if (index == 0) return; // Do nothing on root
+      Value currently_being_moved = data[index];
+      distance_type currently_being_moved_dist =
+        get(distance, currently_being_moved);
+      for (;;) {
+        if (index == 0) break; // Stop at root
+        size_type parent_index = parent(index);
+        Value parent_value = data[parent_index];
+        if (compare(currently_being_moved_dist, get(distance, parent_value))) {
+          ++num_levels_moved;
+          index = parent_index;
+          continue;
+        } else {
+          break; // Heap property satisfied
+        }
+      }
+      // Actually do the moves -- move num_levels_moved elements down in the
+      // tree, then put currently_being_moved at the top
+      index = orig_index;
+      for (size_type i = 0; i < num_levels_moved; ++i) {
+        size_type parent_index = parent(index);
+        Value parent_value = data[parent_index];
+        put(index_in_heap, parent_value, index);
+        data[index] = parent_value;
+        index = parent_index;
+      }
+      data[index] = currently_being_moved;
+      put(index_in_heap, currently_being_moved, index);
+      verify_heap();
+    }
+
+    // From the root, swap elements (each one with its smallest child) if there
+    // are any parent-child pairs that violate the heap property
+    void preserve_heap_property_down() {
+      size_type index = 0;
+      Value currently_being_moved = data[0];
+      distance_type currently_being_moved_dist =
+        get(distance, currently_being_moved);
+      size_type heap_size = data.size();
+      Value* data_ptr = &data[0];
+      for (;;) {
+        size_type first_child_index = child(index, 0);
+        if (first_child_index >= heap_size) break; /* No children */
+        Value* child_base_ptr = data_ptr + first_child_index;
+        size_type smallest_child_index = 0;
+        distance_type smallest_child_dist = get(distance, child_base_ptr[smallest_child_index]);
+        if (first_child_index + Arity <= heap_size) {
+          // Special case for a statically known loop count (common case)
+          for (size_t i = 1; i < Arity; ++i) {
+            Value i_value = child_base_ptr[i];
+            distance_type i_dist = get(distance, i_value);
+            if (compare(i_dist, smallest_child_dist)) {
+              smallest_child_index = i;
+              smallest_child_dist = i_dist;
+            }
+          }
+        } else {
+          for (size_t i = 1; i < heap_size - first_child_index; ++i) {
+            distance_type i_dist = get(distance, child_base_ptr[i]);
+            if (compare(i_dist, smallest_child_dist)) {
+              smallest_child_index = i;
+              smallest_child_dist = i_dist;
+            }
+          }
+        }
+        if (compare(smallest_child_dist, currently_being_moved_dist)) {
+          swap_heap_elements(smallest_child_index + first_child_index, index);
+          index = smallest_child_index + first_child_index;
+          continue;
+        } else {
+          break; // Heap property satisfied
+        }
+      }
+      verify_heap();
+    }
+
+  };
+
+} // namespace boost
+
+#endif // BOOST_D_ARY_HEAP_HPP

Modified: trunk/boost/graph/dijkstra_shortest_paths.hpp
==============================================================================
--- trunk/boost/graph/dijkstra_shortest_paths.hpp	(original)
+++ trunk/boost/graph/dijkstra_shortest_paths.hpp	2009-01-20 16:37:31 EST (Tue, 20 Jan 2009)
@@ -17,6 +17,8 @@
 #include <boost/pending/indirect_cmp.hpp>
 #include <boost/graph/exception.hpp>
 #include <boost/pending/relaxed_heap.hpp>
+#include <boost/smart_ptr.hpp>
+#include <boost/graph/detail/d_ary_heap.hpp>
 
 #ifdef BOOST_GRAPH_DIJKSTRA_TESTING
 #  include <boost/pending/mutable_queue.hpp>
@@ -90,18 +92,18 @@
 
       template <class Edge, class Graph>
       void tree_edge(Edge e, Graph& g) {
-        m_decreased = relax(e, g, m_weight, m_predecessor, m_distance,
-                            m_combine, m_compare);
-        if (m_decreased)
+        bool decreased = relax(e, g, m_weight, m_predecessor, m_distance,
+                               m_combine, m_compare);
+        if (decreased)
           m_vis.edge_relaxed(e, g);
         else
           m_vis.edge_not_relaxed(e, g);
       }
       template <class Edge, class Graph>
       void gray_target(Edge e, Graph& g) {
-        m_decreased = relax(e, g, m_weight, m_predecessor, m_distance,
-                            m_combine, m_compare);
-        if (m_decreased) {
+        bool decreased = relax(e, g, m_weight, m_predecessor, m_distance,
+                               m_combine, m_compare);
+        if (decreased) {
           m_Q.update(target(e, g));
           m_vis.edge_relaxed(e, g);
         } else
@@ -134,7 +136,6 @@
       DistanceMap m_distance;
       BinaryFunction m_combine;
       BinaryPredicate m_compare;
-      bool m_decreased;
       D m_zero;
     };
 
@@ -182,10 +183,19 @@
 
 #ifdef BOOST_GRAPH_DIJKSTRA_TESTING
     if (!dijkstra_relaxed_heap) {
+#ifdef BOOST_GRAPH_TEST_D_ARY_HEAP
+      boost::scoped_array<std::size_t> index_in_heap_map(new std::size_t[num_vertices(g)]);
+      typedef boost::iterator_property_map<std::size_t*, IndexMap> IndexInHeapMap;
+      IndexInHeapMap index_in_heap(&index_in_heap_map[0], index_map);
+      typedef d_ary_heap_indirect<Vertex, 4, IndexInHeapMap, DistanceMap, Compare>
+        MutableQueue;
+      MutableQueue Q(distance, index_in_heap, compare);
+#else
       typedef mutable_queue<Vertex, std::vector<Vertex>, IndirectCmp, IndexMap>
         MutableQueue;
 
       MutableQueue Q(num_vertices(g), icmp, index_map);
+#endif
 
       detail::dijkstra_bfs_visitor<DijkstraVisitor, MutableQueue, WeightMap,
         PredecessorMap, DistanceMap, Combine, Compare>

Modified: trunk/libs/graph/test/dijkstra_heap_performance.cpp
==============================================================================
--- trunk/libs/graph/test/dijkstra_heap_performance.cpp	(original)
+++ trunk/libs/graph/test/dijkstra_heap_performance.cpp	2009-01-20 16:37:31 EST (Tue, 20 Jan 2009)
@@ -8,6 +8,7 @@
 //           Andrew Lumsdaine
 #ifndef BOOST_GRAPH_DIJKSTRA_TESTING_DIETMAR
 #  define BOOST_GRAPH_DIJKSTRA_TESTING
+#  define BOOST_GRAPH_TEST_D_ARY_HEAP
 #endif
 
 #include <boost/graph/dijkstra_shortest_paths.hpp>
@@ -105,8 +106,12 @@
   std::vector<double> binary_heap_distances(n);
   std::vector<double> relaxed_heap_distances(n);
 
-  // Run binary heap version
+  // Run binary or d-ary heap version
+#ifdef BOOST_GRAPH_TEST_D_ARY_HEAP
+  std::cout << "Running Dijkstra's with d-ary heap...";
+#else
   std::cout << "Running Dijkstra's with binary heap...";
+#endif
   std::cout.flush();
   timer t;
 #ifdef BOOST_GRAPH_DIJKSTRA_TESTING_DIETMAR