Add support for maximum weight matching

This also adds more consistency with the heuristic version of the code.

src/graph/flow/Makefile.am

@@ -18,7 +18,6 @@ libgraph_tool_flow_la_SOURCES = \
     graph_edmonds_karp.cc \
     graph_push_relabel.cc \
     graph_kolmogorov.cc \
-    graph_maximum_cardinality_matching.cc \
     graph_minimum_cut.cc \
     graph_flow_bind.cc
 

src/graph/flow/graph_flow_bind.cc

@@ -26,7 +26,6 @@ void push_relabel_max_flow(GraphInterface& gi, size_t src, size_t sink,
                            boost::any capacity, boost::any res);
 void kolmogorov_max_flow(GraphInterface& gi, size_t src, size_t sink,
                          boost::any capacity, boost::any res);
-bool max_cardinality_matching(GraphInterface& gi, boost::any match);
 double min_cut(GraphInterface& gi, boost::any weight, boost::any part_map);
 void get_residual_graph(GraphInterface& gi, boost::any capacity, boost::any res,
                         boost::any oaugment);
@@ -40,7 +39,6 @@ BOOST_PYTHON_MODULE(libgraph_tool_flow)
     def("edmonds_karp_max_flow", &edmonds_karp_max_flow);
     def("push_relabel_max_flow", &push_relabel_max_flow);
     def("kolmogorov_max_flow", &kolmogorov_max_flow);
-    def("max_cardinality_matching", &max_cardinality_matching);
     def("min_cut", &min_cut);
     def("residual_graph", &get_residual_graph);
 }

src/graph/flow/graph_maximum_cardinality_matching.cc

-// graph-tool -- a general graph modification and manipulation thingy
-//
-// Copyright (C) 2006-2019 Tiago de Paula Peixoto <tiago@skewed.de>
-//
-// This program is free software; you can redistribute it and/or
-// modify it under the terms of the GNU General Public License
-// as published by the Free Software Foundation; either version 3
-// of the License, or (at your option) any later version.
-//
-// This program is distributed in the hope that it will be useful,
-// but WITHOUT ANY WARRANTY; without even the implied warranty of
-// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
-// GNU General Public License for more details.
-//
-// You should have received a copy of the GNU General Public License
-// along with this program. If not, see <http://www.gnu.org/licenses/>.
-
-#include "graph_filtering.hh"
-#include "graph_selectors.hh"
-#include "graph_properties.hh"
-#include "graph.hh"
-
-#include "graph_augment.hh"
-
-#include <boost/graph/max_cardinality_matching.hpp>
-
-#include <boost/bind.hpp>
-
-using namespace graph_tool;
-using namespace boost;
-
-struct get_max_cardinality_matching
-{
-    template <class Graph, class VertexIndex, class MatchMap>
-    void operator()(Graph& g, VertexIndex vertex_index, MatchMap match,
-                    bool &check) const
-    {
-        typedef typename graph_traits<Graph>::vertex_descriptor vertex_t;
-        unchecked_vector_property_map<vertex_t,VertexIndex>
-            mate(vertex_index, num_vertices(g));
-
-        typename graph_traits<Graph>::edge_iterator e, e_end;
-        for (tie(e, e_end) = edges(g); e != e_end; ++e)
-            match[*e] = false;
-
-        check =
-            checked_edmonds_maximum_cardinality_matching(g, mate, vertex_index);
-
-        for (tie(e, e_end) = edges(g); e != e_end; ++e)
-        {
-            if (mate[source(*e,g)] != graph_traits<Graph>::null_vertex() &&
-                mate[source(*e,g)] == target(*e,g))
-            {
-                // search for an already matched parallel edge
-                bool matched = false;
-                typename graph_traits<Graph>::out_edge_iterator oe, oe_end;
-                for (tie(oe, oe_end) = out_edges(source(*e, g),g); oe != oe_end;
-                     ++oe)
-                {
-                    if (match[*oe])
-                    {
-                        matched = true;
-                        break;
-                    }
-                }
-                if (!matched)
-                    match[*e] = true;
-            }
-        }
-    }
-};
-
-
-bool max_cardinality_matching(GraphInterface& gi, boost::any match)
-{
-    bool check;
-    run_action<graph_tool::detail::never_directed>()
-        (gi, std::bind(get_max_cardinality_matching(),
-                        std::placeholders::_1, gi.get_vertex_index(),
-                       std::placeholders::_2, std::ref(check)),
-         writable_edge_scalar_properties()) (match);
-    return check;
-}

src/graph/topology/Makefile.am

@@ -24,6 +24,7 @@ libgraph_tool_topology_la_SOURCES = \
     graph_dominator_tree.cc \
     graph_isomorphism.cc \
     graph_kcore.cc \
+    graph_matching.cc \
     graph_maximal_cliques.cc \
     graph_maximal_planar.cc \
     graph_maximal_vertex_set.cc \

src/graph/topology/graph_matching.cc

+// graph-tool -- a general graph modification and manipulation thingy
+//
+// Copyright (C) 2006-2019 Tiago de Paula Peixoto <tiago@skewed.de>
+//
+// This program is free software; you can redistribute it and/or
+// modify it under the terms of the GNU General Public License
+// as published by the Free Software Foundation; either version 3
+// of the License, or (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program. If not, see <http://www.gnu.org/licenses/>.
+
+#include "graph_filtering.hh"
+#include "graph.hh"
+#include "graph_properties.hh"
+
+#include <boost/graph/max_cardinality_matching.hpp>
+#include <boost/graph/maximum_weighted_matching.hpp>
+
+using namespace std;
+using namespace boost;
+using namespace graph_tool;
+
+
+void get_max_matching(GraphInterface& gi, std::string initial_matching,
+                      boost::any omatching)
+{
+    typedef typename vprop_map_t<int64_t>::type vprop_t;
+
+    vprop_t::unchecked_t matching = any_cast<vprop_t>(omatching).get_unchecked();
+
+    run_action<graph_tool::detail::never_directed>()
+        (gi,
+         [&](auto& g)
+         {
+             auto vindex = get(vertex_index, g);
+             typedef decltype(vindex) vindex_t;
+             typedef std::remove_reference_t<decltype(g)> g_t;
+
+             if (initial_matching == "empty")
+                 boost::matching<g_t, vprop_t::unchecked_t, vindex_t,
+                                 edmonds_augmenting_path_finder, empty_matching,
+                                 no_matching_verifier>
+                     (g, matching, vindex);
+             else if (initial_matching == "greedy")
+                 boost::matching<g_t, vprop_t::unchecked_t, vindex_t,
+                                 edmonds_augmenting_path_finder, greedy_matching,
+                                 no_matching_verifier>
+                     (g, matching, vindex);
+             else if (initial_matching == "extra_greedy")
+                 boost::matching<g_t, vprop_t::unchecked_t, vindex_t,
+                                 edmonds_augmenting_path_finder, extra_greedy_matching,
+                                 no_matching_verifier>
+                     (g, matching, vindex);
+             else
+                 throw ValueException("invalid initial matching: " +
+                                      initial_matching);
+
+             for (auto v : vertices_range(g))
+             {
+                 if (matching[v] == int64_t(graph_traits<g_t>::null_vertex()))
+                     matching[v] = std::numeric_limits<int64_t>::max();
+             }
+
+         })();
+}
+
+void get_max_weighted_matching(GraphInterface& gi, boost::any oweight,
+                               boost::any omatching, bool brute_force)
+{
+    typedef typename vprop_map_t<int64_t>::type vprop_t;
+
+    vprop_t::unchecked_t matching = any_cast<vprop_t>(omatching).get_unchecked();
+
+    run_action<graph_tool::detail::never_directed>()
+        (gi,
+         [&](auto& g, auto w)
+         {
+             typedef std::remove_reference_t<decltype(g)> g_t;
+
+             typedef typename graph_traits<g_t>::vertex_descriptor vertex_t;
+             typename vprop_map_t<vertex_t>::type match(get(vertex_index,g));
+
+             if (brute_force)
+                 brute_force_maximum_weighted_matching(g, w, match);
+             else
+                 maximum_weighted_matching(g, w, match);
+
+             for (auto v : vertices_range(g))
+             {
+                 if (match[v] == graph_traits<g_t>::null_vertex())
+                     matching[v] = std::numeric_limits<int64_t>::max();
+                 else
+                     matching[v] = match[v];
+             }
+         },
+         writable_edge_scalar_properties())(oweight);
+}
+
+void match_edges(GraphInterface& gi, boost::any omatching,
+                 boost::any oematching)
+{
+    typedef typename vprop_map_t<int64_t>::type vprop_t;
+    typedef typename eprop_map_t<uint8_t>::type eprop_t;
+
+    vprop_t::unchecked_t matching = any_cast<vprop_t>(omatching).get_unchecked();
+    eprop_t::unchecked_t ematching = any_cast<eprop_t>(oematching).get_unchecked();
+
+    run_action<graph_tool::detail::never_directed>()
+        (gi,
+         [&](auto& g)
+         {
+             for (auto v : vertices_range(g))
+             {
+                 auto u = matching[v];
+                 if (size_t(u) > num_vertices(g))
+                     continue;
+                 ematching[edge(v, u, g).first] = true;
+             }
+         })();
+}
+
+#include <boost/python.hpp>
+using namespace boost::python;
+
+void export_matching()
+{
+    def("get_max_matching", &get_max_matching);
+    def("get_max_weighted_matching", &get_max_weighted_matching);
+    def("match_edges", &match_edges);
+}

src/graph/topology/graph_random_matching.cc

@@ -43,8 +43,11 @@ struct do_random_matching
 
         vector<vertex_t> vlist;
         typename graph_traits<Graph>::vertex_iterator v, v_end;
-        for (tie(v, v_end) = vertices(g); v != v_end; ++v)
-            vlist.push_back(*v);
+        for (auto v : vertices_range(g))
+        {
+            vlist.push_back(v);
+            match[v] = numeric_limits<typename property_traits<MatchMap>::value_type>::max();
+        }
 
         unchecked_vector_property_map<uint8_t, VertexIndex>
             matched(vertex_index, num_vertices(g));
@@ -82,7 +85,8 @@ struct do_random_matching
             {
                 uniform_int_distribution<> sample(0, candidates.size() - 1);
                 size_t j = sample(rng);
-                match[candidates[j]] = true;
+                match[v] = target(candidates[j], g);
+                match[target(candidates[j], g)] = v;
                 matched[v] = true;
                 matched[target(candidates[j], g)] = true;
             }
@@ -103,7 +107,7 @@ void random_matching(GraphInterface& gi, boost::any weight, boost::any match,
     run_action<>()
         (gi, std::bind(do_random_matching(), std::placeholders::_1, gi.get_vertex_index(),
                        std::placeholders::_2, std::placeholders::_3, minimize, std::ref(rng)),
-         edge_props_t(), writable_edge_scalar_properties())(weight, match);
+         edge_props_t(), writable_vertex_scalar_properties())(weight, match);
 }
 
 void export_random_matching()

src/graph/topology/graph_topology.cc

@@ -65,6 +65,7 @@ void export_random_matching();
 void export_maximal_vertex_set();
 void export_vertex_similarity();
 void export_max_cliques();
+void export_matching();
 
 
 BOOST_PYTHON_MODULE(libgraph_tool_topology)
@@ -96,4 +97,5 @@ BOOST_PYTHON_MODULE(libgraph_tool_topology)
     export_maximal_vertex_set();
     export_vertex_similarity();
     export_max_cliques();
+    export_matching();
 }

src/graph_tool/draw/__init__.py

@@ -413,7 +413,7 @@ def _coarse_graph(g, vweight, eweight, mivs=False, groups=None):
             mivs = None
             m = max_cardinality_matching(GraphView(g, directed=False),
                                          heuristic=True, weight=eweight,
-                                         minimize=False)
+                                         minimize=False, edges=True)
             u = GraphView(g, efilt=m, directed=False)
             c = label_components(u)[0]
             c = g.own_property(c)

src/graph_tool/topology/__init__.py

@@ -2933,30 +2933,38 @@ def is_DAG(g):
     return is_DAG
 
 
-def max_cardinality_matching(g, heuristic=False, weight=None, minimize=True,
-                             match=None):
+def max_cardinality_matching(g, weight=None, init_match="extra_greedy",
+                             heuristic=False, minimize=False, edges=False,
+                             brute_force=False):
     r"""Find a maximum cardinality matching in the graph.
 
     Parameters
     ----------
     g : :class:`~graph_tool.Graph`
         Graph to be used.
-    heuristic : bool (optional, default: `False`)
-        If true, a random heuristic will be used, which runs in linear time.
     weight : :class:`~graph_tool.EdgePropertyMap` (optional, default: `None`)
-        If provided, the matching will minimize the edge weights (or maximize
-        if ``minimize == False``). This option has no effect if
-        ``heuristic == False``.
+        If provided, the matching will maximize the sum of edge weights.
+    init_match : string (optional, default: ``"extra_greedy"``)
+        Initial matching strategy. Can be one of: `"empty"`, `"greedy"`,
+        `"extra_greedy"`. Ignored if ``weight`` is given, or
+        ``heuristic == True``.
     minimize : bool (optional, default: `True`)
         If `True`, the matching will minimize the weights, otherwise they will
         be maximized. This option has no effect if ``heuristic == False``.
-    match : :class:`~graph_tool.EdgePropertyMap` (optional, default: `None`)
-        Edge property map where the matching will be specified.
+    heuristic : bool (optional, default: `False`)
+        If `True`, a random heuristic will be used, which runs in linear time.
+    edges : bool (optional, default: `False`)
+        If `True`, an edge property map will be returned, instead of a vertex
+        property map.
+    brute_force : bool (optional, default: `False`)
+        If `True`, and `weight` is not `None` and `heuristic` is `False`, a
+        slower, brute-force algorithm is used.
 
     Returns
     -------
-    match : :class:`~graph_tool.EdgePropertyMap`
-        Boolean edge property map where the matching is specified.
+    match : :class:`~graph_tool.VertexPropertyMap`
+        Vertex property map where the matching is specified. If ``edges ==
+        True`` a boolean-valued edge property map is returned instead.
 
     Notes
     -----
@@ -2964,19 +2972,24 @@ def max_cardinality_matching(g, heuristic=False, weight=None, minimize=True,
     share a common vertex. A *maximum cardinality matching* has maximum size
     over all matchings in the graph.
 
-    If the parameter ``weight`` is provided, as well as ``heuristic == True`` a
-    matching with maximum cardinality *and* maximum (or minimum) weight is
-    returned.
+    If the parameter ``weight`` is provided, a matching with maximum cardinality
+    *and* maximum weight is returned.
 
     If ``heuristic == True`` the algorithm does not necessarily return the
     maximum matching, instead the focus is to run on linear time.
 
     This algorithm runs in time :math:`O(EV\times\alpha(E,V))`, where
     :math:`\alpha(m,n)` is a slow growing function that is at most 4 for any
-    feasible input. If `heuristic == True`, the algorithm runs in time
-    :math:`O(V + E)`.
+    feasible input.
+
+    If weights are given, the algorithm runs in time :math:`O(V^3)`.
+
+    If `heuristic == True`, the algorithm runs in time :math:`O(V + E)`.
 
-    For a more detailed description, see [boost-max-matching]_.
+    If `brute_force == True`, the algorithm runs in time :math:`O(exp(E))`.
+
+    For a more detailed description, see [boost-max-matching]_ and
+    [boost-max-weighted-matching]_.
 
     Examples
     --------
@@ -2987,10 +3000,8 @@ def max_cardinality_matching(g, heuristic=False, weight=None, minimize=True,
        gt.seed_rng(43)
 
     >>> g = gt.GraphView(gt.price_network(300), directed=False)
-    >>> res = gt.max_cardinality_matching(g)
-    >>> print(res[1])
-    True
-    >>> w = res[0].copy("double")
+    >>> w = gt.max_cardinality_matching(g, edges=True)
+    >>> w = w.copy("double")
     >>> w.a = 2 * w.a + 2
     >>> gt.graph_draw(g, edge_color=res[0], edge_pen_width=w, vertex_fill_color="grey",
     ...               output="max_card_match.pdf")
@@ -3010,28 +3021,39 @@ def max_cardinality_matching(g, heuristic=False, weight=None, minimize=True,
     References
     ----------
     .. [boost-max-matching] http://www.boost.org/libs/graph/doc/maximum_matching.html
+    .. [boost-max-weighted-matching] http://www.boost.org/libs/graph/doc/maximum_weighted_matching.html
     .. [matching-heuristic] B. Hendrickson and R. Leland. "A Multilevel Algorithm
        for Partitioning Graphs." In S. Karin, editor, Proc. Supercomputing ’95,
        San Diego. ACM Press, New York, 1995, :doi:`10.1145/224170.224228`
 
     """
-    if match is None:
-        match = g.new_edge_property("bool")
-    _check_prop_scalar(match, "match")
-    _check_prop_writable(match, "match")
+    match = g.new_vp("int64_t")
     if weight is not None:
         _check_prop_scalar(weight, "weight")
 
     u = GraphView(g, directed=False)
     if not heuristic:
-        check = libgraph_tool_flow.\
-                max_cardinality_matching(u._Graph__graph, _prop("e", u, match))
-        return match, check
+        if weight is None:
+            libgraph_tool_topology.\
+                get_max_matching(u._Graph__graph, init_match,
+                                 _prop("v", u, match))
+        else:
+            libgraph_tool_topology.\
+                get_max_weighted_matching(u._Graph__graph,
+                                          _prop("e", u, weight),
+                                          _prop("v", u, match), brute_force)
     else:
-        libgraph_tool_topology.\
+         libgraph_tool_topology.\
                 random_matching(u._Graph__graph, _prop("e", u, weight),
-                                 _prop("e", u, match), minimize, _get_rng())
-        return match
+                                 _prop("v", u, match), minimize, _get_rng())
+
+    if edges:
+        ematch = g.new_ep("bool")
+        libgraph_tool_topology.match_edges(u._Graph__graph,
+                                           _prop("v", u, match),
+                                           _prop("e", u, ematch))
+        return ematch
+    return match
 
 
 def max_independent_vertex_set(g, high_deg=False, mivs=None):