1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/epoc32/include/stdapis/boost/graph/max_cardinality_matching.hpp Tue Mar 16 16:12:26 2010 +0000
1.3 @@ -0,0 +1,876 @@
1.4 +//=======================================================================
1.5 +// Copyright (c) 2005 Aaron Windsor
1.6 +//
1.7 +// Distributed under the Boost Software License, Version 1.0.
1.8 +// (See accompanying file LICENSE_1_0.txt or copy at
1.9 +// http://www.boost.org/LICENSE_1_0.txt)
1.10 +//
1.11 +//=======================================================================
1.12 +
1.13 +#ifndef BOOST_GRAPH_MAXIMUM_CARDINALITY_MATCHING_HPP
1.14 +#define BOOST_GRAPH_MAXIMUM_CARDINALITY_MATCHING_HPP
1.15 +
1.16 +#include <vector>
1.17 +#include <list>
1.18 +#include <deque>
1.19 +#include <algorithm> // for std::sort and std::stable_sort
1.20 +#include <utility> // for std::pair
1.21 +#include <boost/property_map.hpp>
1.22 +#include <boost/utility.hpp> // for boost::tie
1.23 +#include <boost/graph/graph_traits.hpp>
1.24 +#include <boost/graph/visitors.hpp>
1.25 +#include <boost/graph/depth_first_search.hpp>
1.26 +#include <boost/graph/filtered_graph.hpp>
1.27 +#include <boost/pending/disjoint_sets.hpp>
1.28 +#include <boost/assert.hpp>
1.29 +
1.30 +
1.31 +namespace boost
1.32 +{
1.33 + namespace graph { namespace detail {
1.34 + enum { V_EVEN, V_ODD, V_UNREACHED };
1.35 + } } // end namespace graph::detail
1.36 +
1.37 + template <typename Graph, typename MateMap, typename VertexIndexMap>
1.38 + typename graph_traits<Graph>::vertices_size_type
1.39 + matching_size(const Graph& g, MateMap mate, VertexIndexMap vm)
1.40 + {
1.41 + typedef typename graph_traits<Graph>::vertex_iterator vertex_iterator_t;
1.42 + typedef typename graph_traits<Graph>::vertex_descriptor
1.43 + vertex_descriptor_t;
1.44 + typedef typename graph_traits<Graph>::vertices_size_type v_size_t;
1.45 +
1.46 + v_size_t size_of_matching = 0;
1.47 + vertex_iterator_t vi, vi_end;
1.48 +
1.49 + for(tie(vi,vi_end) = vertices(g); vi != vi_end; ++vi)
1.50 + {
1.51 + vertex_descriptor_t v = *vi;
1.52 + if (get(mate,v) != graph_traits<Graph>::null_vertex()
1.53 + && get(vm,v) < get(vm,get(mate,v)))
1.54 + ++size_of_matching;
1.55 + }
1.56 + return size_of_matching;
1.57 + }
1.58 +
1.59 +
1.60 +
1.61 +
1.62 + template <typename Graph, typename MateMap>
1.63 + inline typename graph_traits<Graph>::vertices_size_type
1.64 + matching_size(const Graph& g, MateMap mate)
1.65 + {
1.66 + return matching_size(g, mate, get(vertex_index,g));
1.67 + }
1.68 +
1.69 +
1.70 +
1.71 +
1.72 + template <typename Graph, typename MateMap, typename VertexIndexMap>
1.73 + bool is_a_matching(const Graph& g, MateMap mate, VertexIndexMap vm)
1.74 + {
1.75 + typedef typename graph_traits<Graph>::vertex_descriptor
1.76 + vertex_descriptor_t;
1.77 + typedef typename graph_traits<Graph>::vertex_iterator vertex_iterator_t;
1.78 +
1.79 + vertex_iterator_t vi, vi_end;
1.80 + for( tie(vi,vi_end) = vertices(g); vi != vi_end; ++vi)
1.81 + {
1.82 + vertex_descriptor_t v = *vi;
1.83 + if (get(mate,v) != graph_traits<Graph>::null_vertex()
1.84 + && v != get(mate,get(mate,v)))
1.85 + return false;
1.86 + }
1.87 + return true;
1.88 + }
1.89 +
1.90 +
1.91 +
1.92 +
1.93 + template <typename Graph, typename MateMap>
1.94 + inline bool is_a_matching(const Graph& g, MateMap mate)
1.95 + {
1.96 + return is_a_matching(g, mate, get(vertex_index,g));
1.97 + }
1.98 +
1.99 +
1.100 +
1.101 +
1.102 + //***************************************************************************
1.103 + //***************************************************************************
1.104 + // Maximum Cardinality Matching Functors
1.105 + //***************************************************************************
1.106 + //***************************************************************************
1.107 +
1.108 + template <typename Graph, typename MateMap,
1.109 + typename VertexIndexMap = dummy_property_map>
1.110 + struct no_augmenting_path_finder
1.111 + {
1.112 + no_augmenting_path_finder(const Graph& g, MateMap mate, VertexIndexMap vm)
1.113 + { }
1.114 +
1.115 + inline bool augment_matching() { return false; }
1.116 +
1.117 + template <typename PropertyMap>
1.118 + void get_current_matching(PropertyMap p) {}
1.119 + };
1.120 +
1.121 +
1.122 +
1.123 +
1.124 + template <typename Graph, typename MateMap, typename VertexIndexMap>
1.125 + class edmonds_augmenting_path_finder
1.126 + {
1.127 + // This implementation of Edmonds' matching algorithm closely
1.128 + // follows Tarjan's description of the algorithm in "Data
1.129 + // Structures and Network Algorithms."
1.130 +
1.131 + public:
1.132 +
1.133 + //generates the type of an iterator property map from vertices to type X
1.134 + template <typename X>
1.135 + struct map_vertex_to_
1.136 + {
1.137 + typedef boost::iterator_property_map<typename std::vector<X>::iterator,
1.138 + VertexIndexMap> type;
1.139 + };
1.140 +
1.141 + typedef typename graph_traits<Graph>::vertex_descriptor
1.142 + vertex_descriptor_t;
1.143 + typedef typename std::pair< vertex_descriptor_t, vertex_descriptor_t >
1.144 + vertex_pair_t;
1.145 + typedef typename graph_traits<Graph>::edge_descriptor edge_descriptor_t;
1.146 + typedef typename graph_traits<Graph>::vertices_size_type v_size_t;
1.147 + typedef typename graph_traits<Graph>::edges_size_type e_size_t;
1.148 + typedef typename graph_traits<Graph>::vertex_iterator vertex_iterator_t;
1.149 + typedef typename graph_traits<Graph>::out_edge_iterator
1.150 + out_edge_iterator_t;
1.151 + typedef typename std::deque<vertex_descriptor_t> vertex_list_t;
1.152 + typedef typename std::vector<edge_descriptor_t> edge_list_t;
1.153 + typedef typename map_vertex_to_<vertex_descriptor_t>::type
1.154 + vertex_to_vertex_map_t;
1.155 + typedef typename map_vertex_to_<int>::type vertex_to_int_map_t;
1.156 + typedef typename map_vertex_to_<vertex_pair_t>::type
1.157 + vertex_to_vertex_pair_map_t;
1.158 + typedef typename map_vertex_to_<v_size_t>::type vertex_to_vsize_map_t;
1.159 + typedef typename map_vertex_to_<e_size_t>::type vertex_to_esize_map_t;
1.160 +
1.161 +
1.162 +
1.163 +
1.164 + edmonds_augmenting_path_finder(const Graph& arg_g, MateMap arg_mate,
1.165 + VertexIndexMap arg_vm) :
1.166 + g(arg_g),
1.167 + vm(arg_vm),
1.168 + n_vertices(num_vertices(arg_g)),
1.169 +
1.170 + mate_vector(n_vertices),
1.171 + ancestor_of_v_vector(n_vertices),
1.172 + ancestor_of_w_vector(n_vertices),
1.173 + vertex_state_vector(n_vertices),
1.174 + origin_vector(n_vertices),
1.175 + pred_vector(n_vertices),
1.176 + bridge_vector(n_vertices),
1.177 + ds_parent_vector(n_vertices),
1.178 + ds_rank_vector(n_vertices),
1.179 +
1.180 + mate(mate_vector.begin(), vm),
1.181 + ancestor_of_v(ancestor_of_v_vector.begin(), vm),
1.182 + ancestor_of_w(ancestor_of_w_vector.begin(), vm),
1.183 + vertex_state(vertex_state_vector.begin(), vm),
1.184 + origin(origin_vector.begin(), vm),
1.185 + pred(pred_vector.begin(), vm),
1.186 + bridge(bridge_vector.begin(), vm),
1.187 + ds_parent_map(ds_parent_vector.begin(), vm),
1.188 + ds_rank_map(ds_rank_vector.begin(), vm),
1.189 +
1.190 + ds(ds_rank_map, ds_parent_map)
1.191 + {
1.192 + vertex_iterator_t vi, vi_end;
1.193 + for(tie(vi,vi_end) = vertices(g); vi != vi_end; ++vi)
1.194 + mate[*vi] = get(arg_mate, *vi);
1.195 + }
1.196 +
1.197 +
1.198 +
1.199 +
1.200 + bool augment_matching()
1.201 + {
1.202 + //As an optimization, some of these values can be saved from one
1.203 + //iteration to the next instead of being re-initialized each
1.204 + //iteration, allowing for "lazy blossom expansion." This is not
1.205 + //currently implemented.
1.206 +
1.207 + e_size_t timestamp = 0;
1.208 + even_edges.clear();
1.209 +
1.210 + vertex_iterator_t vi, vi_end;
1.211 + for(tie(vi,vi_end) = vertices(g); vi != vi_end; ++vi)
1.212 + {
1.213 + vertex_descriptor_t u = *vi;
1.214 +
1.215 + origin[u] = u;
1.216 + pred[u] = u;
1.217 + ancestor_of_v[u] = 0;
1.218 + ancestor_of_w[u] = 0;
1.219 + ds.make_set(u);
1.220 +
1.221 + if (mate[u] == graph_traits<Graph>::null_vertex())
1.222 + {
1.223 + vertex_state[u] = graph::detail::V_EVEN;
1.224 + out_edge_iterator_t ei, ei_end;
1.225 + for(tie(ei,ei_end) = out_edges(u,g); ei != ei_end; ++ei)
1.226 + even_edges.push_back( *ei );
1.227 + }
1.228 + else
1.229 + vertex_state[u] = graph::detail::V_UNREACHED;
1.230 + }
1.231 +
1.232 + //end initializations
1.233 +
1.234 + vertex_descriptor_t v,w,w_free_ancestor,v_free_ancestor;
1.235 + w_free_ancestor = graph_traits<Graph>::null_vertex();
1.236 + v_free_ancestor = graph_traits<Graph>::null_vertex();
1.237 + bool found_alternating_path = false;
1.238 +
1.239 + while(!even_edges.empty() && !found_alternating_path)
1.240 + {
1.241 + // since we push even edges onto the back of the list as
1.242 + // they're discovered, taking them off the back will search
1.243 + // for augmenting paths depth-first.
1.244 + edge_descriptor_t current_edge = even_edges.back();
1.245 + even_edges.pop_back();
1.246 +
1.247 + v = source(current_edge,g);
1.248 + w = target(current_edge,g);
1.249 +
1.250 + vertex_descriptor_t v_prime = origin[ds.find_set(v)];
1.251 + vertex_descriptor_t w_prime = origin[ds.find_set(w)];
1.252 +
1.253 + // because of the way we put all of the edges on the queue,
1.254 + // v_prime should be labeled V_EVEN; the following is a
1.255 + // little paranoid but it could happen...
1.256 + if (vertex_state[v_prime] != graph::detail::V_EVEN)
1.257 + {
1.258 + std::swap(v_prime,w_prime);
1.259 + std::swap(v,w);
1.260 + }
1.261 +
1.262 + if (vertex_state[w_prime] == graph::detail::V_UNREACHED)
1.263 + {
1.264 + vertex_state[w_prime] = graph::detail::V_ODD;
1.265 + vertex_state[mate[w_prime]] = graph::detail::V_EVEN;
1.266 + out_edge_iterator_t ei, ei_end;
1.267 + for( tie(ei,ei_end) = out_edges(mate[w_prime], g); ei != ei_end; ++ei)
1.268 + even_edges.push_back(*ei);
1.269 + pred[w_prime] = v;
1.270 + }
1.271 + //w_prime == v_prime can happen below if we get an edge that has been
1.272 + //shrunk into a blossom
1.273 + else if (vertex_state[w_prime] == graph::detail::V_EVEN && w_prime != v_prime)
1.274 + {
1.275 + vertex_descriptor_t w_up = w_prime;
1.276 + vertex_descriptor_t v_up = v_prime;
1.277 + vertex_descriptor_t nearest_common_ancestor
1.278 + = graph_traits<Graph>::null_vertex();
1.279 + w_free_ancestor = graph_traits<Graph>::null_vertex();
1.280 + v_free_ancestor = graph_traits<Graph>::null_vertex();
1.281 +
1.282 + // We now need to distinguish between the case that
1.283 + // w_prime and v_prime share an ancestor under the
1.284 + // "parent" relation, in which case we've found a
1.285 + // blossom and should shrink it, or the case that
1.286 + // w_prime and v_prime both have distinct ancestors that
1.287 + // are free, in which case we've found an alternating
1.288 + // path between those two ancestors.
1.289 +
1.290 + ++timestamp;
1.291 +
1.292 + while (nearest_common_ancestor == graph_traits<Graph>::null_vertex() &&
1.293 + (v_free_ancestor == graph_traits<Graph>::null_vertex() ||
1.294 + w_free_ancestor == graph_traits<Graph>::null_vertex()
1.295 + )
1.296 + )
1.297 + {
1.298 + ancestor_of_w[w_up] = timestamp;
1.299 + ancestor_of_v[v_up] = timestamp;
1.300 +
1.301 + if (w_free_ancestor == graph_traits<Graph>::null_vertex())
1.302 + w_up = parent(w_up);
1.303 + if (v_free_ancestor == graph_traits<Graph>::null_vertex())
1.304 + v_up = parent(v_up);
1.305 +
1.306 + if (mate[v_up] == graph_traits<Graph>::null_vertex())
1.307 + v_free_ancestor = v_up;
1.308 + if (mate[w_up] == graph_traits<Graph>::null_vertex())
1.309 + w_free_ancestor = w_up;
1.310 +
1.311 + if (ancestor_of_w[v_up] == timestamp)
1.312 + nearest_common_ancestor = v_up;
1.313 + else if (ancestor_of_v[w_up] == timestamp)
1.314 + nearest_common_ancestor = w_up;
1.315 + else if (v_free_ancestor == w_free_ancestor &&
1.316 + v_free_ancestor != graph_traits<Graph>::null_vertex())
1.317 + nearest_common_ancestor = v_up;
1.318 + }
1.319 +
1.320 + if (nearest_common_ancestor == graph_traits<Graph>::null_vertex())
1.321 + found_alternating_path = true; //to break out of the loop
1.322 + else
1.323 + {
1.324 + //shrink the blossom
1.325 + link_and_set_bridges(w_prime, nearest_common_ancestor, std::make_pair(w,v));
1.326 + link_and_set_bridges(v_prime, nearest_common_ancestor, std::make_pair(v,w));
1.327 + }
1.328 + }
1.329 + }
1.330 +
1.331 + if (!found_alternating_path)
1.332 + return false;
1.333 +
1.334 + // retrieve the augmenting path and put it in aug_path
1.335 + reversed_retrieve_augmenting_path(v, v_free_ancestor);
1.336 + retrieve_augmenting_path(w, w_free_ancestor);
1.337 +
1.338 + // augment the matching along aug_path
1.339 + vertex_descriptor_t a,b;
1.340 + while (!aug_path.empty())
1.341 + {
1.342 + a = aug_path.front();
1.343 + aug_path.pop_front();
1.344 + b = aug_path.front();
1.345 + aug_path.pop_front();
1.346 + mate[a] = b;
1.347 + mate[b] = a;
1.348 + }
1.349 +
1.350 + return true;
1.351 +
1.352 + }
1.353 +
1.354 +
1.355 +
1.356 +
1.357 + template <typename PropertyMap>
1.358 + void get_current_matching(PropertyMap pm)
1.359 + {
1.360 + vertex_iterator_t vi,vi_end;
1.361 + for(tie(vi,vi_end) = vertices(g); vi != vi_end; ++vi)
1.362 + put(pm, *vi, mate[*vi]);
1.363 + }
1.364 +
1.365 +
1.366 +
1.367 +
1.368 + template <typename PropertyMap>
1.369 + void get_vertex_state_map(PropertyMap pm)
1.370 + {
1.371 + vertex_iterator_t vi,vi_end;
1.372 + for(tie(vi,vi_end) = vertices(g); vi != vi_end; ++vi)
1.373 + put(pm, *vi, vertex_state[origin[ds.find_set(*vi)]]);
1.374 + }
1.375 +
1.376 +
1.377 +
1.378 +
1.379 + private:
1.380 +
1.381 + vertex_descriptor_t parent(vertex_descriptor_t x)
1.382 + {
1.383 + if (vertex_state[x] == graph::detail::V_EVEN
1.384 + && mate[x] != graph_traits<Graph>::null_vertex())
1.385 + return mate[x];
1.386 + else if (vertex_state[x] == graph::detail::V_ODD)
1.387 + return origin[ds.find_set(pred[x])];
1.388 + else
1.389 + return x;
1.390 + }
1.391 +
1.392 +
1.393 +
1.394 +
1.395 + void link_and_set_bridges(vertex_descriptor_t x,
1.396 + vertex_descriptor_t stop_vertex,
1.397 + vertex_pair_t the_bridge)
1.398 + {
1.399 + for(vertex_descriptor_t v = x; v != stop_vertex; v = parent(v))
1.400 + {
1.401 + ds.union_set(v, stop_vertex);
1.402 + origin[ds.find_set(stop_vertex)] = stop_vertex;
1.403 +
1.404 + if (vertex_state[v] == graph::detail::V_ODD)
1.405 + {
1.406 + bridge[v] = the_bridge;
1.407 + out_edge_iterator_t oei, oei_end;
1.408 + for(tie(oei, oei_end) = out_edges(v,g); oei != oei_end; ++oei)
1.409 + even_edges.push_back(*oei);
1.410 + }
1.411 + }
1.412 + }
1.413 +
1.414 +
1.415 + // Since none of the STL containers support both constant-time
1.416 + // concatenation and reversal, the process of expanding an
1.417 + // augmenting path once we know one exists is a little more
1.418 + // complicated than it has to be. If we know the path is from v to
1.419 + // w, then the augmenting path is recursively defined as:
1.420 + //
1.421 + // path(v,w) = [v], if v = w
1.422 + // = concat([v, mate[v]], path(pred[mate[v]], w),
1.423 + // if v != w and vertex_state[v] == graph::detail::V_EVEN
1.424 + // = concat([v], reverse(path(x,mate[v])), path(y,w)),
1.425 + // if v != w, vertex_state[v] == graph::detail::V_ODD, and bridge[v] = (x,y)
1.426 + //
1.427 + // These next two mutually recursive functions implement this definition.
1.428 +
1.429 + void retrieve_augmenting_path(vertex_descriptor_t v, vertex_descriptor_t w)
1.430 + {
1.431 + if (v == w)
1.432 + aug_path.push_back(v);
1.433 + else if (vertex_state[v] == graph::detail::V_EVEN)
1.434 + {
1.435 + aug_path.push_back(v);
1.436 + aug_path.push_back(mate[v]);
1.437 + retrieve_augmenting_path(pred[mate[v]], w);
1.438 + }
1.439 + else //vertex_state[v] == graph::detail::V_ODD
1.440 + {
1.441 + aug_path.push_back(v);
1.442 + reversed_retrieve_augmenting_path(bridge[v].first, mate[v]);
1.443 + retrieve_augmenting_path(bridge[v].second, w);
1.444 + }
1.445 + }
1.446 +
1.447 +
1.448 + void reversed_retrieve_augmenting_path(vertex_descriptor_t v,
1.449 + vertex_descriptor_t w)
1.450 + {
1.451 +
1.452 + if (v == w)
1.453 + aug_path.push_back(v);
1.454 + else if (vertex_state[v] == graph::detail::V_EVEN)
1.455 + {
1.456 + reversed_retrieve_augmenting_path(pred[mate[v]], w);
1.457 + aug_path.push_back(mate[v]);
1.458 + aug_path.push_back(v);
1.459 + }
1.460 + else //vertex_state[v] == graph::detail::V_ODD
1.461 + {
1.462 + reversed_retrieve_augmenting_path(bridge[v].second, w);
1.463 + retrieve_augmenting_path(bridge[v].first, mate[v]);
1.464 + aug_path.push_back(v);
1.465 + }
1.466 + }
1.467 +
1.468 +
1.469 +
1.470 +
1.471 + //private data members
1.472 +
1.473 + const Graph& g;
1.474 + VertexIndexMap vm;
1.475 + v_size_t n_vertices;
1.476 +
1.477 + //storage for the property maps below
1.478 + std::vector<vertex_descriptor_t> mate_vector;
1.479 + std::vector<e_size_t> ancestor_of_v_vector;
1.480 + std::vector<e_size_t> ancestor_of_w_vector;
1.481 + std::vector<int> vertex_state_vector;
1.482 + std::vector<vertex_descriptor_t> origin_vector;
1.483 + std::vector<vertex_descriptor_t> pred_vector;
1.484 + std::vector<vertex_pair_t> bridge_vector;
1.485 + std::vector<vertex_descriptor_t> ds_parent_vector;
1.486 + std::vector<v_size_t> ds_rank_vector;
1.487 +
1.488 + //iterator property maps
1.489 + vertex_to_vertex_map_t mate;
1.490 + vertex_to_esize_map_t ancestor_of_v;
1.491 + vertex_to_esize_map_t ancestor_of_w;
1.492 + vertex_to_int_map_t vertex_state;
1.493 + vertex_to_vertex_map_t origin;
1.494 + vertex_to_vertex_map_t pred;
1.495 + vertex_to_vertex_pair_map_t bridge;
1.496 + vertex_to_vertex_map_t ds_parent_map;
1.497 + vertex_to_vsize_map_t ds_rank_map;
1.498 +
1.499 + vertex_list_t aug_path;
1.500 + edge_list_t even_edges;
1.501 + disjoint_sets< vertex_to_vsize_map_t, vertex_to_vertex_map_t > ds;
1.502 +
1.503 + };
1.504 +
1.505 +
1.506 +
1.507 +
1.508 + //***************************************************************************
1.509 + //***************************************************************************
1.510 + // Initial Matching Functors
1.511 + //***************************************************************************
1.512 + //***************************************************************************
1.513 +
1.514 + template <typename Graph, typename MateMap>
1.515 + struct greedy_matching
1.516 + {
1.517 + typedef typename graph_traits< Graph >::vertex_descriptor vertex_descriptor_t;
1.518 + typedef typename graph_traits< Graph >::vertex_iterator vertex_iterator_t;
1.519 + typedef typename graph_traits< Graph >::edge_descriptor edge_descriptor_t;
1.520 + typedef typename graph_traits< Graph >::edge_iterator edge_iterator_t;
1.521 +
1.522 + static void find_matching(const Graph& g, MateMap mate)
1.523 + {
1.524 + vertex_iterator_t vi, vi_end;
1.525 + for(tie(vi,vi_end) = vertices(g); vi != vi_end; ++vi)
1.526 + put(mate, *vi, graph_traits<Graph>::null_vertex());
1.527 +
1.528 + edge_iterator_t ei, ei_end;
1.529 + for( tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1.530 + {
1.531 + edge_descriptor_t e = *ei;
1.532 + vertex_descriptor_t u = source(e,g);
1.533 + vertex_descriptor_t v = target(e,g);
1.534 +
1.535 + if (get(mate,u) == get(mate,v))
1.536 + //only way equality can hold is if
1.537 + // mate[u] == mate[v] == null_vertex
1.538 + {
1.539 + put(mate,u,v);
1.540 + put(mate,v,u);
1.541 + }
1.542 + }
1.543 + }
1.544 + };
1.545 +
1.546 +
1.547 +
1.548 +
1.549 + template <typename Graph, typename MateMap>
1.550 + struct extra_greedy_matching
1.551 + {
1.552 + // The "extra greedy matching" is formed by repeating the
1.553 + // following procedure as many times as possible: Choose the
1.554 + // unmatched vertex v of minimum non-zero degree. Choose the
1.555 + // neighbor w of v which is unmatched and has minimum degree over
1.556 + // all of v's neighbors. Add (u,v) to the matching. Ties for
1.557 + // either choice are broken arbitrarily. This procedure takes time
1.558 + // O(m log n), where m is the number of edges in the graph and n
1.559 + // is the number of vertices.
1.560 +
1.561 + typedef typename graph_traits< Graph >::vertex_descriptor
1.562 + vertex_descriptor_t;
1.563 + typedef typename graph_traits< Graph >::vertex_iterator vertex_iterator_t;
1.564 + typedef typename graph_traits< Graph >::edge_descriptor edge_descriptor_t;
1.565 + typedef typename graph_traits< Graph >::edge_iterator edge_iterator_t;
1.566 + typedef std::pair<vertex_descriptor_t, vertex_descriptor_t> vertex_pair_t;
1.567 +
1.568 + struct select_first
1.569 + {
1.570 + inline static vertex_descriptor_t select_vertex(const vertex_pair_t p)
1.571 + {return p.first;}
1.572 + };
1.573 +
1.574 + struct select_second
1.575 + {
1.576 + inline static vertex_descriptor_t select_vertex(const vertex_pair_t p)
1.577 + {return p.second;}
1.578 + };
1.579 +
1.580 + template <class PairSelector>
1.581 + class less_than_by_degree
1.582 + {
1.583 + public:
1.584 + less_than_by_degree(const Graph& g): m_g(g) {}
1.585 + bool operator() (const vertex_pair_t x, const vertex_pair_t y)
1.586 + {
1.587 + return
1.588 + out_degree(PairSelector::select_vertex(x), m_g)
1.589 + < out_degree(PairSelector::select_vertex(y), m_g);
1.590 + }
1.591 + private:
1.592 + const Graph& m_g;
1.593 + };
1.594 +
1.595 +
1.596 + static void find_matching(const Graph& g, MateMap mate)
1.597 + {
1.598 + typedef std::vector<std::pair<vertex_descriptor_t, vertex_descriptor_t> >
1.599 + directed_edges_vector_t;
1.600 +
1.601 + directed_edges_vector_t edge_list;
1.602 + vertex_iterator_t vi, vi_end;
1.603 + for(tie(vi, vi_end) = vertices(g); vi != vi_end; ++vi)
1.604 + put(mate, *vi, graph_traits<Graph>::null_vertex());
1.605 +
1.606 + edge_iterator_t ei, ei_end;
1.607 + for(tie(ei, ei_end) = edges(g); ei != ei_end; ++ei)
1.608 + {
1.609 + edge_descriptor_t e = *ei;
1.610 + vertex_descriptor_t u = source(e,g);
1.611 + vertex_descriptor_t v = target(e,g);
1.612 + edge_list.push_back(std::make_pair(u,v));
1.613 + edge_list.push_back(std::make_pair(v,u));
1.614 + }
1.615 +
1.616 + //sort the edges by the degree of the target, then (using a
1.617 + //stable sort) by degree of the source
1.618 + std::sort(edge_list.begin(), edge_list.end(),
1.619 + less_than_by_degree<select_second>(g));
1.620 + std::stable_sort(edge_list.begin(), edge_list.end(),
1.621 + less_than_by_degree<select_first>(g));
1.622 +
1.623 + //construct the extra greedy matching
1.624 + for(typename directed_edges_vector_t::const_iterator itr = edge_list.begin(); itr != edge_list.end(); ++itr)
1.625 + {
1.626 + if (get(mate,itr->first) == get(mate,itr->second))
1.627 + //only way equality can hold is if mate[itr->first] == mate[itr->second] == null_vertex
1.628 + {
1.629 + put(mate, itr->first, itr->second);
1.630 + put(mate, itr->second, itr->first);
1.631 + }
1.632 + }
1.633 + }
1.634 + };
1.635 +
1.636 +
1.637 +
1.638 +
1.639 + template <typename Graph, typename MateMap>
1.640 + struct empty_matching
1.641 + {
1.642 + typedef typename graph_traits< Graph >::vertex_iterator vertex_iterator_t;
1.643 +
1.644 + static void find_matching(const Graph& g, MateMap mate)
1.645 + {
1.646 + vertex_iterator_t vi, vi_end;
1.647 + for(tie(vi,vi_end) = vertices(g); vi != vi_end; ++vi)
1.648 + put(mate, *vi, graph_traits<Graph>::null_vertex());
1.649 + }
1.650 + };
1.651 +
1.652 +
1.653 +
1.654 +
1.655 + //***************************************************************************
1.656 + //***************************************************************************
1.657 + // Matching Verifiers
1.658 + //***************************************************************************
1.659 + //***************************************************************************
1.660 +
1.661 + namespace detail
1.662 + {
1.663 +
1.664 + template <typename SizeType>
1.665 + class odd_components_counter : public dfs_visitor<>
1.666 + // This depth-first search visitor will count the number of connected
1.667 + // components with an odd number of vertices. It's used by
1.668 + // maximum_matching_verifier.
1.669 + {
1.670 + public:
1.671 + odd_components_counter(SizeType& c_count):
1.672 + m_count(c_count)
1.673 + {
1.674 + m_count = 0;
1.675 + }
1.676 +
1.677 + template <class Vertex, class Graph>
1.678 + void start_vertex(Vertex v, Graph&)
1.679 + {
1.680 + addend = -1;
1.681 + }
1.682 +
1.683 + template <class Vertex, class Graph>
1.684 + void discover_vertex(Vertex u, Graph&)
1.685 + {
1.686 + addend *= -1;
1.687 + m_count += addend;
1.688 + }
1.689 +
1.690 + protected:
1.691 + SizeType& m_count;
1.692 +
1.693 + private:
1.694 + SizeType addend;
1.695 +
1.696 + };
1.697 +
1.698 + }//namespace detail
1.699 +
1.700 +
1.701 +
1.702 +
1.703 + template <typename Graph, typename MateMap,
1.704 + typename VertexIndexMap = dummy_property_map>
1.705 + struct no_matching_verifier
1.706 + {
1.707 + inline static bool
1.708 + verify_matching(const Graph& g, MateMap mate, VertexIndexMap vm)
1.709 + { return true;}
1.710 + };
1.711 +
1.712 +
1.713 +
1.714 +
1.715 + template <typename Graph, typename MateMap, typename VertexIndexMap>
1.716 + struct maximum_cardinality_matching_verifier
1.717 + {
1.718 +
1.719 + template <typename X>
1.720 + struct map_vertex_to_
1.721 + {
1.722 + typedef boost::iterator_property_map<typename std::vector<X>::iterator,
1.723 + VertexIndexMap> type;
1.724 + };
1.725 +
1.726 + typedef typename graph_traits<Graph>::vertex_descriptor
1.727 + vertex_descriptor_t;
1.728 + typedef typename graph_traits<Graph>::vertices_size_type v_size_t;
1.729 + typedef typename graph_traits<Graph>::vertex_iterator vertex_iterator_t;
1.730 + typedef typename map_vertex_to_<int>::type vertex_to_int_map_t;
1.731 + typedef typename map_vertex_to_<vertex_descriptor_t>::type
1.732 + vertex_to_vertex_map_t;
1.733 +
1.734 + template <typename VertexStateMap>
1.735 + struct non_odd_vertex {
1.736 + //this predicate is used to create a filtered graph that
1.737 + //excludes vertices labeled "graph::detail::V_ODD"
1.738 + non_odd_vertex() : vertex_state(0) { }
1.739 + non_odd_vertex(VertexStateMap* arg_vertex_state)
1.740 + : vertex_state(arg_vertex_state) { }
1.741 + template <typename Vertex>
1.742 + bool operator()(const Vertex& v) const {
1.743 + BOOST_ASSERT(vertex_state);
1.744 + return get(*vertex_state, v) != graph::detail::V_ODD;
1.745 + }
1.746 + VertexStateMap* vertex_state;
1.747 + };
1.748 +
1.749 + static bool verify_matching(const Graph& g, MateMap mate, VertexIndexMap vm)
1.750 + {
1.751 + //For any graph G, let o(G) be the number of connected
1.752 + //components in G of odd size. For a subset S of G's vertex set
1.753 + //V(G), let (G - S) represent the subgraph of G induced by
1.754 + //removing all vertices in S from G. Let M(G) be the size of the
1.755 + //maximum cardinality matching in G. Then the Tutte-Berge
1.756 + //formula guarantees that
1.757 + //
1.758 + // 2 * M(G) = min ( |V(G)| + |U| + o(G - U) )
1.759 + //
1.760 + //where the minimum is taken over all subsets U of
1.761 + //V(G). Edmonds' algorithm finds a set U that achieves the
1.762 + //minimum in the above formula, namely the vertices labeled
1.763 + //"ODD." This function runs one iteration of Edmonds' algorithm
1.764 + //to find U, then verifies that the size of the matching given
1.765 + //by mate satisfies the Tutte-Berge formula.
1.766 +
1.767 + //first, make sure it's a valid matching
1.768 + if (!is_a_matching(g,mate,vm))
1.769 + return false;
1.770 +
1.771 + //We'll try to augment the matching once. This serves two
1.772 + //purposes: first, if we find some augmenting path, the matching
1.773 + //is obviously non-maximum. Second, running edmonds' algorithm
1.774 + //on a graph with no augmenting path will create the
1.775 + //Edmonds-Gallai decomposition that we need as a certificate of
1.776 + //maximality - we can get it by looking at the vertex_state map
1.777 + //that results.
1.778 + edmonds_augmenting_path_finder<Graph,MateMap,VertexIndexMap>
1.779 + augmentor(g,mate,vm);
1.780 + if (augmentor.augment_matching())
1.781 + return false;
1.782 +
1.783 + std::vector<int> vertex_state_vector(num_vertices(g));
1.784 + vertex_to_int_map_t vertex_state(vertex_state_vector.begin(), vm);
1.785 + augmentor.get_vertex_state_map(vertex_state);
1.786 +
1.787 + //count the number of graph::detail::V_ODD vertices
1.788 + v_size_t num_odd_vertices = 0;
1.789 + vertex_iterator_t vi, vi_end;
1.790 + for(tie(vi,vi_end) = vertices(g); vi != vi_end; ++vi)
1.791 + if (vertex_state[*vi] == graph::detail::V_ODD)
1.792 + ++num_odd_vertices;
1.793 +
1.794 + //count the number of connected components with odd cardinality
1.795 + //in the graph without graph::detail::V_ODD vertices
1.796 + non_odd_vertex<vertex_to_int_map_t> filter(&vertex_state);
1.797 + filtered_graph<Graph, keep_all, non_odd_vertex<vertex_to_int_map_t> > fg(g, keep_all(), filter);
1.798 +
1.799 + v_size_t num_odd_components;
1.800 + detail::odd_components_counter<v_size_t> occ(num_odd_components);
1.801 + depth_first_search(fg, visitor(occ).vertex_index_map(vm));
1.802 +
1.803 + if (2 * matching_size(g,mate,vm) == num_vertices(g) + num_odd_vertices - num_odd_components)
1.804 + return true;
1.805 + else
1.806 + return false;
1.807 + }
1.808 + };
1.809 +
1.810 +
1.811 +
1.812 +
1.813 + template <typename Graph,
1.814 + typename MateMap,
1.815 + typename VertexIndexMap,
1.816 + template <typename, typename, typename> class AugmentingPathFinder,
1.817 + template <typename, typename> class InitialMatchingFinder,
1.818 + template <typename, typename, typename> class MatchingVerifier>
1.819 + bool matching(const Graph& g, MateMap mate, VertexIndexMap vm)
1.820 + {
1.821 +
1.822 + InitialMatchingFinder<Graph,MateMap>::find_matching(g,mate);
1.823 +
1.824 + AugmentingPathFinder<Graph,MateMap,VertexIndexMap> augmentor(g,mate,vm);
1.825 + bool not_maximum_yet = true;
1.826 + while(not_maximum_yet)
1.827 + {
1.828 + not_maximum_yet = augmentor.augment_matching();
1.829 + }
1.830 + augmentor.get_current_matching(mate);
1.831 +
1.832 + return MatchingVerifier<Graph,MateMap,VertexIndexMap>::verify_matching(g,mate,vm);
1.833 +
1.834 + }
1.835 +
1.836 +
1.837 +
1.838 +
1.839 + template <typename Graph, typename MateMap, typename VertexIndexMap>
1.840 + inline bool checked_edmonds_maximum_cardinality_matching(const Graph& g, MateMap mate, VertexIndexMap vm)
1.841 + {
1.842 + return matching
1.843 + < Graph, MateMap, VertexIndexMap,
1.844 + edmonds_augmenting_path_finder, extra_greedy_matching, maximum_cardinality_matching_verifier>
1.845 + (g, mate, vm);
1.846 + }
1.847 +
1.848 +
1.849 +
1.850 +
1.851 + template <typename Graph, typename MateMap>
1.852 + inline bool checked_edmonds_maximum_cardinality_matching(const Graph& g, MateMap mate)
1.853 + {
1.854 + return checked_edmonds_maximum_cardinality_matching(g, mate, get(vertex_index,g));
1.855 + }
1.856 +
1.857 +
1.858 +
1.859 +
1.860 + template <typename Graph, typename MateMap, typename VertexIndexMap>
1.861 + inline void edmonds_maximum_cardinality_matching(const Graph& g, MateMap mate, VertexIndexMap vm)
1.862 + {
1.863 + matching < Graph, MateMap, VertexIndexMap,
1.864 + edmonds_augmenting_path_finder, extra_greedy_matching, no_matching_verifier>
1.865 + (g, mate, vm);
1.866 + }
1.867 +
1.868 +
1.869 +
1.870 +
1.871 + template <typename Graph, typename MateMap>
1.872 + inline void edmonds_maximum_cardinality_matching(const Graph& g, MateMap mate)
1.873 + {
1.874 + edmonds_maximum_cardinality_matching(g, mate, get(vertex_index,g));
1.875 + }
1.876 +
1.877 +}//namespace boost
1.878 +
1.879 +#endif //BOOST_GRAPH_MAXIMUM_CARDINALITY_MATCHING_HPP