Fix Stoer Wagner

include/boost/graph/stoer_wagner_min_cut.hpp

@@ -29,94 +29,130 @@ namespace boost
 
 namespace detail
 {
-    template < typename ParityMap, typename WeightMap, typename IndexMap >
-    class mas_min_cut_visitor : public boost::default_mas_visitor
+    /**
+     * \brief Performs a phase of the Stoer-Wagner min-cut algorithm
+     *
+     * Performs a phase of the Stoer-Wagner min-cut algorithm.
+     *
+     * As described by Stoer & Wagner (1997), a phase is simply a maximum
+     * adjacency search (also called a maximum cardinality search), which
+     * results in the selection of two vertices \em s and \em t, and, as a side
+     * product, a minimum <em>s</em>-<em>t</em> cut of the input graph. Here,
+     * the input graph is basically \p g, but some vertices are virtually
+     * assigned to others as a way of viewing \p g as a graph with some sets of
+     * vertices merged together.
+     *
+     * This implementation is a translation of pseudocode by Professor Uri
+     * Zwick, School of Computer Science, Tel Aviv University.
+     *
+     * \pre \p g is a connected, undirected graph
+     * \param[in] g the input graph
+     * \param[in] assignments a read/write property map from each vertex to the
+     *                        vertex that it is assigned to
+     * \param[in] assignedVertices a list of vertices that are assigned to
+     *                             others
+     * \param[in] weights a readable property map from each edge to its
+     *                    weight (a non-negative value)
+     * \param[out] pq a keyed, updatable max-priority queue
+     * \returns a tuple (\em s, \em t, \em w) of the "<em>s</em>" and
+     *          "<em>t</em>" of the minimum <em>s</em>-<em>t</em> cut and the
+     *          cut weight \em w of the minimum <em>s</em>-<em>t</em> cut.
+     * \see http://www.cs.tau.ac.il/~zwick/grad-algo-08/gmc.pdf
+     *
+     * \author Daniel Trebbien
+     * \date 2010-09-11
+     */
+    template < class UndirectedGraph, class VertexAssignmentMap,
+        class WeightMap, class KeyedUpdatablePriorityQueue >
+    boost::tuple<
+        typename boost::graph_traits< UndirectedGraph >::vertex_descriptor,
+        typename boost::graph_traits< UndirectedGraph >::vertex_descriptor,
+        typename boost::property_traits< WeightMap >::value_type >
+    stoer_wagner_phase(const UndirectedGraph& g,
+        VertexAssignmentMap assignments,
+        const std::set< typename boost::graph_traits<
+            UndirectedGraph >::vertex_descriptor >& assignedVertices,
+        WeightMap weights, KeyedUpdatablePriorityQueue& pq)
     {
-        typedef one_bit_color_map< IndexMap > InternalParityMap;
+        typedef
+            typename boost::graph_traits< UndirectedGraph >::vertex_descriptor
+                vertex_descriptor;
         typedef typename boost::property_traits< WeightMap >::value_type
             weight_type;
 
-    public:
-        template < typename Graph >
-        mas_min_cut_visitor(const Graph& g, ParityMap parity,
-            weight_type& cutweight, const WeightMap& weight_map,
-            IndexMap index_map)
-        : m_bestParity(parity)
-        , m_parity(make_one_bit_color_map(num_vertices(g), index_map))
-        , m_bestWeight(cutweight)
-        , m_cutweight(0)
-        , m_visited(0)
-        , m_weightMap(weight_map)
-        {
-            // set here since the init list sets the reference
-            m_bestWeight = (std::numeric_limits< weight_type >::max)();
-        }
-
-        template < typename Vertex, typename Graph >
-        void initialize_vertex(Vertex u, const Graph& g)
-        {
-            typedef typename boost::property_traits< ParityMap >::value_type
-                parity_type;
-            typedef
-                typename boost::property_traits< InternalParityMap >::value_type
-                    internal_parity_type;
-
-            put(m_parity, u, internal_parity_type(0));
-            put(m_bestParity, u, parity_type(0));
-        }
+        BOOST_ASSERT(pq.empty());
+        typename KeyedUpdatablePriorityQueue::key_map keys = pq.keys();
 
-        template < typename Edge, typename Graph >
-        void examine_edge(Edge e, const Graph& g)
+        BGL_FORALL_VERTICES_T(v, g, UndirectedGraph)
         {
-            weight_type w = get(m_weightMap, e);
+            if (v == get(assignments, v))
+            { // foreach u \in V do
+                put(keys, v, weight_type(0));
 
-            // if the target of e is already marked then decrease cutweight
-            // otherwise, increase it
-            if (get(m_parity, boost::target(e, g)))
-            {
-                m_cutweight -= w;
-            }
-            else
-            {
-                m_cutweight += w;
+                pq.push(v);
             }
         }
 
-        template < typename Vertex, typename Graph >
-        void finish_vertex(Vertex u, const Graph& g)
-        {
-            typedef
-                typename boost::property_traits< InternalParityMap >::value_type
-                    internal_parity_type;
+        BOOST_ASSERT(pq.size() >= 2);
 
-            ++m_visited;
-            put(m_parity, u, internal_parity_type(1));
+        vertex_descriptor s
+            = boost::graph_traits< UndirectedGraph >::null_vertex();
+        vertex_descriptor t
+            = boost::graph_traits< UndirectedGraph >::null_vertex();
+        weight_type w;
+        while (!pq.empty())
+        { // while PQ \neq {} do
+            const vertex_descriptor u = pq.top(); // u = extractmax(PQ)
+            w = get(keys, u);
+            pq.pop();
+
+            s = t;
+            t = u;
+
+            BGL_FORALL_OUTEDGES_T(u, e, g, UndirectedGraph)
+            { // foreach (u, v) \in E do
+                const vertex_descriptor v = get(assignments, target(e, g));
+
+                if (pq.contains(v))
+                { // if v \in PQ then
+                    put(keys, v,
+                        get(keys, v)
+                            + get(weights,
+                                e)); // increasekey(PQ, v, wA(v) + w(u, v))
+                    pq.update(v);
+                }
+            }
 
-            if (m_cutweight < m_bestWeight && m_visited < num_vertices(g))
+            typename std::set< vertex_descriptor >::const_iterator
+                assignedVertexIt,
+                assignedVertexEnd = assignedVertices.end();
+            for (assignedVertexIt = assignedVertices.begin();
+                 assignedVertexIt != assignedVertexEnd; ++assignedVertexIt)
             {
-                m_bestWeight = m_cutweight;
-                BGL_FORALL_VERTICES_T(i, g, Graph)
+                const vertex_descriptor uPrime = *assignedVertexIt;
+
+                if (get(assignments, uPrime) == u)
                 {
-                    put(m_bestParity, i, get(m_parity, i));
+                    BGL_FORALL_OUTEDGES_T(uPrime, e, g, UndirectedGraph)
+                    { // foreach (u, v) \in E do
+                        const vertex_descriptor v
+                            = get(assignments, target(e, g));
+
+                        if (pq.contains(v))
+                        { // if v \in PQ then
+                            put(keys, v,
+                                get(keys, v)
+                                    + get(weights, e)); // increasekey(PQ, v,
+                                                        // wA(v) + w(u, v))
+                            pq.update(v);
+                        }
+                    }
                 }
             }
         }
 
-        inline void clear()
-        {
-            m_bestWeight = (std::numeric_limits< weight_type >::max)();
-            m_visited = 0;
-            m_cutweight = 0;
-        }
-
-    private:
-        ParityMap m_bestParity;
-        InternalParityMap m_parity;
-        weight_type& m_bestWeight;
-        weight_type m_cutweight;
-        unsigned m_visited;
-        const WeightMap& m_weightMap;
-    };
+        return boost::make_tuple(s, t, w);
+    }
 
     /**
      * \brief Computes a min-cut of the input graph
@@ -156,44 +192,71 @@ namespace detail
                 vertex_descriptor;
         typedef typename boost::property_traits< WeightMap >::value_type
             weight_type;
+        typedef
+            typename boost::graph_traits< UndirectedGraph >::vertices_size_type
+                vertices_size_type;
+        typedef typename boost::property_traits< ParityMap >::value_type
+            parity_type;
 
-        typename graph_traits< UndirectedGraph >::vertex_iterator u_iter, u_end;
+        vertices_size_type n = num_vertices(g);
 
-        weight_type bestW = (std::numeric_limits< weight_type >::max)();
-        weight_type bestThisTime = (std::numeric_limits< weight_type >::max)();
-        vertex_descriptor bestStart
-            = boost::graph_traits< UndirectedGraph >::null_vertex();
+        std::set< vertex_descriptor > assignedVertices;
+
+        // initialize `assignments` (all vertices are initially assigned to
+        // themselves)
+        BGL_FORALL_VERTICES_T(v, g, UndirectedGraph) { put(assignments, v, v); }
 
-        detail::mas_min_cut_visitor< ParityMap, WeightMap, IndexMap > vis(
-            g, parities, bestThisTime, weights, index_map);
+        vertex_descriptor s, t;
+        weight_type bestW;
 
-        // for each node in the graph,
-        for (boost::tie(u_iter, u_end) = vertices(g); u_iter != u_end; ++u_iter)
+        boost::tie(s, t, bestW) = boost::detail::stoer_wagner_phase(
+            g, assignments, assignedVertices, weights, pq);
+        BOOST_ASSERT(s != t);
+        BGL_FORALL_VERTICES_T(v, g, UndirectedGraph)
         {
-            // run the MAS and find the min cut
-            vis.clear();
-            boost::maximum_adjacency_search(g,
-                boost::weight_map(weights)
-                    .visitor(vis)
-                    .root_vertex(*u_iter)
-                    .vertex_assignment_map(assignments)
-                    .max_priority_queue(pq));
-            if (bestThisTime < bestW)
+            put(parities, v, parity_type(v == t ? 1 : 0));
+        }
+        put(assignments, t, s);
+        assignedVertices.insert(t);
+        --n;
+
+        for (; n >= 2; --n)
+        {
+            weight_type w;
+            boost::tie(s, t, w) = boost::detail::stoer_wagner_phase(
+                g, assignments, assignedVertices, weights, pq);
+            BOOST_ASSERT(s != t);
+
+            if (w < bestW)
+            {
+                BGL_FORALL_VERTICES_T(v, g, UndirectedGraph)
+                {
+                    put(parities, v,
+                        parity_type(get(assignments, v) == t ? 1 : 0));
+
+                    if (get(assignments, v)
+                        == t) // all vertices that were assigned to t are now
+                              // assigned to s
+                        put(assignments, v, s);
+                }
+
+                bestW = w;
+            }
+            else
             {
-                bestW = bestThisTime;
-                bestStart = *u_iter;
+                BGL_FORALL_VERTICES_T(v, g, UndirectedGraph)
+                {
+                    if (get(assignments, v)
+                        == t) // all vertices that were assigned to t are now
+                              // assigned to s
+                        put(assignments, v, s);
+                }
             }
+            put(assignments, t, s);
+            assignedVertices.insert(t);
         }
 
-        // Run one more time, starting from the best start location, to
-        // ensure the visitor has the best values.
-        vis.clear();
-        boost::maximum_adjacency_search(g,
-            boost::vertex_assignment_map(assignments)
-                .weight_map(weights)
-                .visitor(vis)
-                .root_vertex(bestStart)
-                .max_priority_queue(pq));
+        BOOST_ASSERT(pq.empty());
 
         return bestW;
     }

test/stoer_wagner_test.cpp

@@ -184,6 +184,32 @@ void test4()
     BOOST_TEST_EQ(parity2, get(parities, 7));
 }
 
+// Non regression test for github.com/boostorg/graph/issues/286
+void test5()
+{
+    edge_t edges[] = { { 0, 1 }, { 0, 2 }, { 0, 3 }, { 1, 2 }, { 1, 3 },
+        { 2, 3 }, { 4, 5 }, { 4, 6 }, { 4, 7 }, { 5, 6 }, { 5, 7 }, { 6, 7 },
+        { 0, 4 } };
+    weight_type ws[] = { 3, 3, 3, 2, 2, 2, 3, 3, 3, 2, 2, 2, 6 };
+    undirected_graph g(edges, edges + 13, ws, 8, 13);
+
+    weight_map_type weights = get(boost::edge_weight, g);
+    std::map< int, bool > parity;
+    boost::associative_property_map< std::map< int, bool > > parities(parity);
+    int w
+        = boost::stoer_wagner_min_cut(g, weights, boost::parity_map(parities));
+    BOOST_TEST_EQ(w, 6);
+    const bool parity0 = get(parities, 0);
+    BOOST_TEST_EQ(parity0, get(parities, 1));
+    BOOST_TEST_EQ(parity0, get(parities, 2));
+    BOOST_TEST_EQ(parity0, get(parities, 3));
+    const bool parity4 = get(parities, 4);
+    BOOST_TEST_NE(parity0, parity4);
+    BOOST_TEST_EQ(parity4, get(parities, 5));
+    BOOST_TEST_EQ(parity4, get(parities, 6));
+    BOOST_TEST_EQ(parity4, get(parities, 7));
+}
+
 // The input for the `test_prgen` family of tests comes from a program, named
 // `prgen`, that comes with a package of min-cut solvers by Chandra Chekuri,
 // Andrew Goldberg, David Karger, Matthew Levine, and Cliff Stein. `prgen` was
@@ -285,13 +311,15 @@ void test_prgen_50_70_2()
 
 int main(int argc, char* argv[])
 {
-    if (BOOST_TEST(argc == 2)) {
+    if (BOOST_TEST(argc == 2))
+    {
         test_dir = argv[1];
         test0();
         test1();
         test2();
         test3();
         test4();
+        test5();
         test_prgen_20_70_2();
         test_prgen_50_70_2();
     }