Bipartite Graph Max Cardinality Matching (Hopcroft-Karp)
authordeba
Fri, 07 Apr 2006 09:51:23 +0000
changeset 2040c7bd55c0d820
parent 2039 dacc4ce9474d
child 2041 28df5272df99
Bipartite Graph Max Cardinality Matching (Hopcroft-Karp)
Test for it

Some BpUgraph improvments
lemon/Makefile.am
lemon/bipartite_matching.h
lemon/bpugraph_adaptor.h
test/Makefile.am
test/bipartite_matching_test.cc
     1.1 --- a/lemon/Makefile.am	Thu Apr 06 09:33:29 2006 +0000
     1.2 +++ b/lemon/Makefile.am	Fri Apr 07 09:51:23 2006 +0000
     1.3 @@ -31,6 +31,8 @@
     1.4  	dfs.h \
     1.5  	bin_heap.h \
     1.6  	bpugraph_adaptor.h \
     1.7 +	bipartite_matching.h \
     1.8 +	bucket_heap.h \
     1.9  	color.h \
    1.10  	config.h \
    1.11  	counter.h \
    1.12 @@ -53,7 +55,6 @@
    1.13  	iterable_maps.h \
    1.14  	johnson.h \
    1.15  	kruskal.h \
    1.16 -	bucket_heap.h \
    1.17  	list_graph.h \
    1.18  	lp.h \
    1.19  	lp_base.h \
     2.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     2.2 +++ b/lemon/bipartite_matching.h	Fri Apr 07 09:51:23 2006 +0000
     2.3 @@ -0,0 +1,466 @@
     2.4 +/* -*- C++ -*-
     2.5 + *
     2.6 + * This file is a part of LEMON, a generic C++ optimization library
     2.7 + *
     2.8 + * Copyright (C) 2003-2006
     2.9 + * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
    2.10 + * (Egervary Research Group on Combinatorial Optimization, EGRES).
    2.11 + *
    2.12 + * Permission to use, modify and distribute this software is granted
    2.13 + * provided that this copyright notice appears in all copies. For
    2.14 + * precise terms see the accompanying LICENSE file.
    2.15 + *
    2.16 + * This software is provided "AS IS" with no warranty of any kind,
    2.17 + * express or implied, and with no claim as to its suitability for any
    2.18 + * purpose.
    2.19 + *
    2.20 + */
    2.21 +
    2.22 +#ifndef LEMON_BIPARTITE_MATCHING
    2.23 +#define LEMON_BIPARTITE_MATCHING
    2.24 +
    2.25 +#include <lemon/bpugraph_adaptor.h>
    2.26 +#include <lemon/bfs.h>
    2.27 +
    2.28 +#include <iostream>
    2.29 +
    2.30 +///\ingroup matching
    2.31 +///\file
    2.32 +///\brief Maximum matching algorithms in bipartite graphs.
    2.33 +
    2.34 +namespace lemon {
    2.35 +
    2.36 +  /// \ingroup matching
    2.37 +  ///
    2.38 +  /// \brief Bipartite Max Cardinality Matching algorithm
    2.39 +  ///
    2.40 +  /// Bipartite Max Cardinality Matching algorithm. This class implements
    2.41 +  /// the Hopcroft-Karp algorithm wich has \f$ O(e\sqrt{n}) \f$ time
    2.42 +  /// complexity.
    2.43 +  template <typename BpUGraph>
    2.44 +  class MaxBipartiteMatching {
    2.45 +  protected:
    2.46 +
    2.47 +    typedef BpUGraph Graph;
    2.48 +
    2.49 +    typedef typename Graph::Node Node;
    2.50 +    typedef typename Graph::ANodeIt ANodeIt;
    2.51 +    typedef typename Graph::BNodeIt BNodeIt;
    2.52 +    typedef typename Graph::UEdge UEdge;
    2.53 +    typedef typename Graph::UEdgeIt UEdgeIt;
    2.54 +    typedef typename Graph::IncEdgeIt IncEdgeIt;
    2.55 +
    2.56 +    typedef typename BpUGraph::template ANodeMap<UEdge> ANodeMatchingMap;
    2.57 +    typedef typename BpUGraph::template BNodeMap<UEdge> BNodeMatchingMap;
    2.58 +
    2.59 +
    2.60 +  public:
    2.61 +
    2.62 +    /// \brief Constructor.
    2.63 +    ///
    2.64 +    /// Constructor of the algorithm. 
    2.65 +    MaxBipartiteMatching(const BpUGraph& _graph) 
    2.66 +      : anode_matching(_graph), bnode_matching(_graph), graph(&_graph) {}
    2.67 +
    2.68 +    /// \name Execution control
    2.69 +    /// The simplest way to execute the algorithm is to use
    2.70 +    /// one of the member functions called \c run().
    2.71 +    /// \n
    2.72 +    /// If you need more control on the execution,
    2.73 +    /// first you must call \ref init() or one alternative for it.
    2.74 +    /// Finally \ref start() will perform the matching computation or
    2.75 +    /// with step-by-step execution you can augment the solution.
    2.76 +
    2.77 +    /// @{
    2.78 +
    2.79 +    /// \brief Initalize the data structures.
    2.80 +    ///
    2.81 +    /// It initalizes the data structures and creates an empty matching.
    2.82 +    void init() {
    2.83 +      for (ANodeIt it(*graph); it != INVALID; ++it) {
    2.84 +        anode_matching[it] = INVALID;
    2.85 +      }
    2.86 +      for (BNodeIt it(*graph); it != INVALID; ++it) {
    2.87 +        bnode_matching[it] = INVALID;
    2.88 +      }
    2.89 +      matching_value = 0;
    2.90 +    }
    2.91 +
    2.92 +    /// \brief Initalize the data structures.
    2.93 +    ///
    2.94 +    /// It initalizes the data structures and creates a greedy
    2.95 +    /// matching.  From this matching sometimes it is faster to get
    2.96 +    /// the matching than from the initial empty matching.
    2.97 +    void greedyInit() {
    2.98 +      matching_value = 0;
    2.99 +      for (BNodeIt it(*graph); it != INVALID; ++it) {
   2.100 +        bnode_matching[it] = INVALID;
   2.101 +      }
   2.102 +      for (ANodeIt it(*graph); it != INVALID; ++it) {
   2.103 +        anode_matching[it] = INVALID;
   2.104 +        for (IncEdgeIt jt(*graph, it); jt != INVALID; ++jt) {
   2.105 +          if (bnode_matching[graph->bNode(jt)] == INVALID) {
   2.106 +            anode_matching[it] = jt;
   2.107 +            bnode_matching[graph->bNode(jt)] = jt;
   2.108 +            ++matching_value;
   2.109 +            break;
   2.110 +          }
   2.111 +        }
   2.112 +      }
   2.113 +    }
   2.114 +
   2.115 +    /// \brief Initalize the data structures with an initial matching.
   2.116 +    ///
   2.117 +    /// It initalizes the data structures with an initial matching.
   2.118 +    template <typename MatchingMap>
   2.119 +    void matchingInit(const MatchingMap& matching) {
   2.120 +      for (ANodeIt it(*graph); it != INVALID; ++it) {
   2.121 +        anode_matching[it] = INVALID;
   2.122 +      }
   2.123 +      for (BNodeIt it(*graph); it != INVALID; ++it) {
   2.124 +        bnode_matching[it] = INVALID;
   2.125 +      }
   2.126 +      matching_value = 0;
   2.127 +      for (UEdgeIt it(*graph); it != INVALID; ++it) {
   2.128 +        if (matching[it]) {
   2.129 +          ++matching_value;
   2.130 +          anode_matching[graph->aNode(it)] = it;
   2.131 +          bnode_matching[graph->bNode(it)] = it;
   2.132 +        }
   2.133 +      }
   2.134 +    }
   2.135 +
   2.136 +    /// \brief Initalize the data structures with an initial matching.
   2.137 +    ///
   2.138 +    /// It initalizes the data structures with an initial matching.
   2.139 +    /// \return %True when the given map contains really a matching.
   2.140 +    template <typename MatchingMap>
   2.141 +    void checkedMatchingInit(const MatchingMap& matching) {
   2.142 +      for (ANodeIt it(*graph); it != INVALID; ++it) {
   2.143 +        anode_matching[it] = INVALID;
   2.144 +      }
   2.145 +      for (BNodeIt it(*graph); it != INVALID; ++it) {
   2.146 +        bnode_matching[it] = INVALID;
   2.147 +      }
   2.148 +      matching_value = 0;
   2.149 +      for (UEdgeIt it(*graph); it != INVALID; ++it) {
   2.150 +        if (matching[it]) {
   2.151 +          ++matching_value;
   2.152 +          if (anode_matching[graph->aNode(it)] != INVALID) {
   2.153 +            return false;
   2.154 +          }
   2.155 +          anode_matching[graph->aNode(it)] = it;
   2.156 +          if (bnode_matching[graph->aNode(it)] != INVALID) {
   2.157 +            return false;
   2.158 +          }
   2.159 +          bnode_matching[graph->bNode(it)] = it;
   2.160 +        }
   2.161 +      }
   2.162 +      return false;
   2.163 +    }
   2.164 +
   2.165 +    /// \brief An augmenting phase of the Hopcroft-Karp algorithm
   2.166 +    ///
   2.167 +    /// It runs an augmenting phase of the Hopcroft-Karp
   2.168 +    /// algorithm. The phase finds maximum count of edge disjoint
   2.169 +    /// augmenting paths and augments on these paths. The algorithm
   2.170 +    /// consists at most of \f$ O(\sqrt{n}) \f$ phase and one phase is
   2.171 +    /// \f$ O(e) \f$ long.
   2.172 +    bool augment() {
   2.173 +
   2.174 +      typename Graph::template ANodeMap<bool> areached(*graph, false);
   2.175 +      typename Graph::template BNodeMap<bool> breached(*graph, false);
   2.176 +
   2.177 +      typename Graph::template BNodeMap<UEdge> bpred(*graph, INVALID);
   2.178 +
   2.179 +      std::vector<Node> queue, bqueue;
   2.180 +      for (ANodeIt it(*graph); it != INVALID; ++it) {
   2.181 +        if (anode_matching[it] == INVALID) {
   2.182 +          queue.push_back(it);
   2.183 +          areached[it] = true;
   2.184 +        }
   2.185 +      }
   2.186 +
   2.187 +      bool success = false;
   2.188 +
   2.189 +      while (!success && !queue.empty()) {
   2.190 +        std::vector<Node> newqueue;
   2.191 +        for (int i = 0; i < (int)queue.size(); ++i) {
   2.192 +          Node anode = queue[i];
   2.193 +          for (IncEdgeIt jt(*graph, anode); jt != INVALID; ++jt) {
   2.194 +            Node bnode = graph->bNode(jt);
   2.195 +            if (breached[bnode]) continue;
   2.196 +            breached[bnode] = true;
   2.197 +            bpred[bnode] = jt;
   2.198 +            if (bnode_matching[bnode] == INVALID) {
   2.199 +              bqueue.push_back(bnode);
   2.200 +              success = true;
   2.201 +            } else {           
   2.202 +              Node newanode = graph->aNode(bnode_matching[bnode]);
   2.203 +              if (!areached[newanode]) {
   2.204 +                areached[newanode] = true;
   2.205 +                newqueue.push_back(newanode);
   2.206 +              }
   2.207 +            }
   2.208 +          }
   2.209 +        }
   2.210 +        queue.swap(newqueue);
   2.211 +      }
   2.212 +
   2.213 +      if (success) {
   2.214 +
   2.215 +        typename Graph::template ANodeMap<bool> aused(*graph, false);
   2.216 +        
   2.217 +        for (int i = 0; i < (int)bqueue.size(); ++i) {
   2.218 +          Node bnode = bqueue[i];
   2.219 +
   2.220 +          bool used = false;
   2.221 +
   2.222 +          while (bnode != INVALID) {
   2.223 +            UEdge uedge = bpred[bnode];
   2.224 +            Node anode = graph->aNode(uedge);
   2.225 +            
   2.226 +            if (aused[anode]) {
   2.227 +              used = true;
   2.228 +              break;
   2.229 +            }
   2.230 +            
   2.231 +            bnode = anode_matching[anode] != INVALID ? 
   2.232 +              graph->bNode(anode_matching[anode]) : INVALID;
   2.233 +            
   2.234 +          }
   2.235 +          
   2.236 +          if (used) continue;
   2.237 +
   2.238 +          bnode = bqueue[i];
   2.239 +          while (bnode != INVALID) {
   2.240 +            UEdge uedge = bpred[bnode];
   2.241 +            Node anode = graph->aNode(uedge);
   2.242 +            
   2.243 +            bnode_matching[bnode] = uedge;
   2.244 +
   2.245 +            bnode = anode_matching[anode] != INVALID ? 
   2.246 +              graph->bNode(anode_matching[anode]) : INVALID;
   2.247 +            
   2.248 +            anode_matching[anode] = uedge;
   2.249 +
   2.250 +            aused[anode] = true;
   2.251 +          }
   2.252 +          ++matching_value;
   2.253 +
   2.254 +        }
   2.255 +      }
   2.256 +      return success;
   2.257 +    }
   2.258 +
   2.259 +    /// \brief An augmenting phase of the Ford-Fulkerson algorithm
   2.260 +    ///
   2.261 +    /// It runs an augmenting phase of the Ford-Fulkerson
   2.262 +    /// algorithm. The phase finds only one augmenting path and 
   2.263 +    /// augments only on this paths. The algorithm consists at most 
   2.264 +    /// of \f$ O(n) \f$ simple phase and one phase is at most 
   2.265 +    /// \f$ O(e) \f$ long.
   2.266 +    bool simpleAugment() {
   2.267 +
   2.268 +      typename Graph::template ANodeMap<bool> areached(*graph, false);
   2.269 +      typename Graph::template BNodeMap<bool> breached(*graph, false);
   2.270 +
   2.271 +      typename Graph::template BNodeMap<UEdge> bpred(*graph, INVALID);
   2.272 +
   2.273 +      std::vector<Node> queue;
   2.274 +      for (ANodeIt it(*graph); it != INVALID; ++it) {
   2.275 +        if (anode_matching[it] == INVALID) {
   2.276 +          queue.push_back(it);
   2.277 +          areached[it] = true;
   2.278 +        }
   2.279 +      }
   2.280 +
   2.281 +      while (!queue.empty()) {
   2.282 +        std::vector<Node> newqueue;
   2.283 +        for (int i = 0; i < (int)queue.size(); ++i) {
   2.284 +          Node anode = queue[i];
   2.285 +          for (IncEdgeIt jt(*graph, anode); jt != INVALID; ++jt) {
   2.286 +            Node bnode = graph->bNode(jt);
   2.287 +            if (breached[bnode]) continue;
   2.288 +            breached[bnode] = true;
   2.289 +            bpred[bnode] = jt;
   2.290 +            if (bnode_matching[bnode] == INVALID) {
   2.291 +              while (bnode != INVALID) {
   2.292 +                UEdge uedge = bpred[bnode];
   2.293 +                anode = graph->aNode(uedge);
   2.294 +                
   2.295 +                bnode_matching[bnode] = uedge;
   2.296 +                
   2.297 +                bnode = anode_matching[anode] != INVALID ? 
   2.298 +                  graph->bNode(anode_matching[anode]) : INVALID;
   2.299 +                
   2.300 +                anode_matching[anode] = uedge;
   2.301 +                
   2.302 +              }
   2.303 +              ++matching_value;
   2.304 +              return true;
   2.305 +            } else {           
   2.306 +              Node newanode = graph->aNode(bnode_matching[bnode]);
   2.307 +              if (!areached[newanode]) {
   2.308 +                areached[newanode] = true;
   2.309 +                newqueue.push_back(newanode);
   2.310 +              }
   2.311 +            }
   2.312 +          }
   2.313 +        }
   2.314 +        queue.swap(newqueue);
   2.315 +      }
   2.316 +      
   2.317 +      return false;
   2.318 +    }
   2.319 +
   2.320 +    /// \brief Starts the algorithm.
   2.321 +    ///
   2.322 +    /// Starts the algorithm. It runs augmenting phases until the optimal
   2.323 +    /// solution reached.
   2.324 +    void start() {
   2.325 +      while (augment()) {}
   2.326 +    }
   2.327 +
   2.328 +    /// \brief Runs the algorithm.
   2.329 +    ///
   2.330 +    /// It just initalize the algorithm and then start it.
   2.331 +    void run() {
   2.332 +      init();
   2.333 +      start();
   2.334 +    }
   2.335 +
   2.336 +    /// @}
   2.337 +
   2.338 +    /// \name Query Functions
   2.339 +    /// The result of the %Matching algorithm can be obtained using these
   2.340 +    /// functions.\n
   2.341 +    /// Before the use of these functions,
   2.342 +    /// either run() or start() must be called.
   2.343 +    
   2.344 +    ///@{
   2.345 +
   2.346 +    /// \brief Returns an minimum covering of the nodes.
   2.347 +    ///
   2.348 +    /// The minimum covering set problem is the dual solution of the
   2.349 +    /// maximum bipartite matching. It provides an solution for this
   2.350 +    /// problem what is proof of the optimality of the matching.
   2.351 +    /// \return The size of the cover set.
   2.352 +    template <typename CoverMap>
   2.353 +    int coverSet(CoverMap& covering) {
   2.354 +
   2.355 +      typename Graph::template ANodeMap<bool> areached(*graph, false);
   2.356 +      typename Graph::template BNodeMap<bool> breached(*graph, false);
   2.357 +      
   2.358 +      std::vector<Node> queue;
   2.359 +      for (ANodeIt it(*graph); it != INVALID; ++it) {
   2.360 +        if (anode_matching[it] == INVALID) {
   2.361 +          queue.push_back(it);
   2.362 +        }
   2.363 +      }
   2.364 +
   2.365 +      while (!queue.empty()) {
   2.366 +        std::vector<Node> newqueue;
   2.367 +        for (int i = 0; i < (int)queue.size(); ++i) {
   2.368 +          Node anode = queue[i];
   2.369 +          for (IncEdgeIt jt(*graph, anode); jt != INVALID; ++jt) {
   2.370 +            Node bnode = graph->bNode(jt);
   2.371 +            if (breached[bnode]) continue;
   2.372 +            breached[bnode] = true;
   2.373 +            if (bnode_matching[bnode] != INVALID) {
   2.374 +              Node newanode = graph->aNode(bnode_matching[bnode]);
   2.375 +              if (!areached[newanode]) {
   2.376 +                areached[newanode] = true;
   2.377 +                newqueue.push_back(newanode);
   2.378 +              }
   2.379 +            }
   2.380 +          }
   2.381 +        }
   2.382 +        queue.swap(newqueue);
   2.383 +      }
   2.384 +
   2.385 +      int size = 0;
   2.386 +      for (ANodeIt it(*graph); it != INVALID; ++it) {
   2.387 +        covering[it] = !areached[it] && anode_matching[it] != INVALID;
   2.388 +        if (!areached[it] && anode_matching[it] != INVALID) {
   2.389 +          ++size;
   2.390 +        }
   2.391 +      }
   2.392 +      for (BNodeIt it(*graph); it != INVALID; ++it) {
   2.393 +        covering[it] = breached[it];
   2.394 +        if (breached[it]) {
   2.395 +          ++size;
   2.396 +        }
   2.397 +      }
   2.398 +      return size;
   2.399 +    }
   2.400 +
   2.401 +    /// \brief Set true all matching uedge in the map.
   2.402 +    /// 
   2.403 +    /// Set true all matching uedge in the map. It does not change the
   2.404 +    /// value mapped to the other uedges.
   2.405 +    /// \return The number of the matching edges.
   2.406 +    template <typename MatchingMap>
   2.407 +    int quickMatching(MatchingMap& matching) {
   2.408 +      for (ANodeIt it(*graph); it != INVALID; ++it) {
   2.409 +        if (anode_matching[it] != INVALID) {
   2.410 +          matching[anode_matching[it]] = true;
   2.411 +        }
   2.412 +      }
   2.413 +      return matching_value;
   2.414 +    }
   2.415 +
   2.416 +    /// \brief Set true all matching uedge in the map and the others to false.
   2.417 +    /// 
   2.418 +    /// Set true all matching uedge in the map and the others to false.
   2.419 +    /// \return The number of the matching edges.
   2.420 +    template <typename MatchingMap>
   2.421 +    int matching(MatchingMap& matching) {
   2.422 +      for (UEdgeIt it(*graph); it != INVALID; ++it) {
   2.423 +        matching[it] = it == anode_matching[graph->aNode(it)];
   2.424 +      }
   2.425 +      return matching_value;
   2.426 +    }
   2.427 +
   2.428 +
   2.429 +    /// \brief Return true if the given uedge is in the matching.
   2.430 +    /// 
   2.431 +    /// It returns true if the given uedge is in the matching.
   2.432 +    bool matchingEdge(const UEdge& edge) {
   2.433 +      return anode_matching[graph->aNode(edge)] == edge;
   2.434 +    }
   2.435 +
   2.436 +    /// \brief Returns the matching edge from the node.
   2.437 +    /// 
   2.438 +    /// Returns the matching edge from the node. If there is not such
   2.439 +    /// edge it gives back \c INVALID.
   2.440 +    UEdge matchingEdge(const Node& node) {
   2.441 +      if (graph->aNode(node)) {
   2.442 +        return anode_matching[node];
   2.443 +      } else {
   2.444 +        return bnode_matching[node];
   2.445 +      }
   2.446 +    }
   2.447 +
   2.448 +    /// \brief Gives back the number of the matching edges.
   2.449 +    ///
   2.450 +    /// Gives back the number of the matching edges.
   2.451 +    int matchingValue() const {
   2.452 +      return matching_value;
   2.453 +    }
   2.454 +
   2.455 +    /// @}
   2.456 +
   2.457 +  private:
   2.458 +
   2.459 +    ANodeMatchingMap anode_matching;
   2.460 +    BNodeMatchingMap bnode_matching;
   2.461 +    const Graph *graph;
   2.462 +
   2.463 +    int matching_value;
   2.464 +  
   2.465 +  };
   2.466 +
   2.467 +}
   2.468 +
   2.469 +#endif
     3.1 --- a/lemon/bpugraph_adaptor.h	Thu Apr 06 09:33:29 2006 +0000
     3.2 +++ b/lemon/bpugraph_adaptor.h	Fri Apr 07 09:51:23 2006 +0000
     3.3 @@ -103,6 +103,12 @@
     3.4      Node source(const Edge& e) const { return graph->source(e); }
     3.5      Node target(const Edge& e) const { return graph->target(e); }
     3.6  
     3.7 +    Node aNode(const UEdge& e) const { return graph->aNode(e); }
     3.8 +    Node bNode(const UEdge& e) const { return graph->bNode(e); }
     3.9 +
    3.10 +    bool aNode(const Node& n) const { return graph->aNode(n); }
    3.11 +    bool bNode(const Node& n) const { return graph->bNode(n); }
    3.12 +
    3.13      typedef NodeNumTagIndicator<Graph> NodeNumTag;
    3.14      int nodeNum() const { return graph->nodeNum(); }
    3.15      int aNodeNum() const { return graph->aNodeNum(); }
    3.16 @@ -122,7 +128,8 @@
    3.17        return graph->findUEdge(source, target, prev);
    3.18      }
    3.19    
    3.20 -    Node addNode() const { return graph->addNode(); }
    3.21 +    Node addANode() const { return graph->addANode(); }
    3.22 +    Node addBNode() const { return graph->addBNode(); }
    3.23      UEdge addEdge(const Node& source, const Node& target) const { 
    3.24        return graph->addEdge(source, target); 
    3.25      }
    3.26 @@ -281,6 +288,11 @@
    3.27    };
    3.28  
    3.29    /// \ingroup graph_adaptors
    3.30 +  ///
    3.31 +  /// \brief Trivial Bipartite Undirected Graph Adaptor
    3.32 +  ///
    3.33 +  /// Trivial Bipartite Undirected Graph Adaptor. It does not change
    3.34 +  /// the functionality of the adapted graph.
    3.35    template <typename _BpUGraph>
    3.36    class BpUGraphAdaptor 
    3.37      : public BpUGraphAdaptorExtender< BpUGraphAdaptorBase<_BpUGraph> > { 
    3.38 @@ -310,6 +322,17 @@
    3.39      typedef typename Parent::Node Node;
    3.40      typedef typename Parent::BNode ANode;
    3.41      typedef typename Parent::ANode BNode;
    3.42 +    typedef typename Parent::Edge Edge;
    3.43 +    typedef typename Parent::UEdge UEdge;
    3.44 +
    3.45 +    bool direction(const Edge& e) const { return !Parent::direction(e); }
    3.46 +    Edge direct(const UEdge& e, bool b) const { return Parent::direct(e, !b); }
    3.47 +
    3.48 +    Node aNode(const UEdge& e) const { return Parent::bNode(e); }
    3.49 +    Node bNode(const UEdge& e) const { return Parent::aNode(e); }
    3.50 +
    3.51 +    bool aNode(const Node& n) const { return Parent::bNode(n); }
    3.52 +    bool bNode(const Node& n) const { return Parent::aNode(n); }
    3.53  
    3.54      void firstANode(Node& i) const { Parent::firstBNode(i); }
    3.55      void firstBNode(Node& i) const { Parent::firstANode(i); }
    3.56 @@ -407,6 +430,126 @@
    3.57    };
    3.58  
    3.59  
    3.60 +  template <typename _BpUGraph, 
    3.61 +            typename _ANMatchingMap, typename _BNMatchingMap>
    3.62 +  class MatchingBpUGraphAdaptorBase 
    3.63 +    : public BpUGraphAdaptorBase<const _BpUGraph>
    3.64 +  {
    3.65 +  public:
    3.66 +    
    3.67 +    typedef _BpUGraph Graph;
    3.68 +    typedef _ANMatchingMap ANodeMatchingMap;
    3.69 +    typedef _BNMatchingMap BNodeMatchingMap;
    3.70 +
    3.71 +    typedef BpUGraphAdaptorBase<const _BpUGraph> Parent;
    3.72 +
    3.73 +  protected:
    3.74 +    
    3.75 +    MatchingBpUGraphAdaptorBase() {}
    3.76 +
    3.77 +    void setANodeMatchingMap(ANodeMatchingMap& _anode_matching) {
    3.78 +      anode_matching = &_anode_matching;
    3.79 +    }
    3.80 +
    3.81 +    void setBNodeMatchingMap(BNodeMatchingMap& _bnode_matching) {
    3.82 +      bnode_matching = &_bnode_matching;
    3.83 +    }
    3.84 +
    3.85 +  public:
    3.86 +
    3.87 +    typedef typename Parent::Node Node;
    3.88 +    typedef typename Parent::Edge Edge;
    3.89 +
    3.90 +    void firstOut(Edge& edge, const Node& node) const {
    3.91 +      if (Parent::aNode(node)) {
    3.92 +        Parent::firstOut(edge, node);
    3.93 +        if (edge != INVALID && (*anode_matching)[node] == edge) {
    3.94 +          Parent::nextOut(edge);
    3.95 +        }
    3.96 +      } else {
    3.97 +        edge = (*bnode_matching)[node] != INVALID ?
    3.98 +          Parent::direct((*bnode_matching)[node], false) : INVALID;
    3.99 +      }
   3.100 +    }
   3.101 +
   3.102 +    void nextOut(Edge& edge) const {
   3.103 +      if (Parent::aNode(Parent::source(edge))) {
   3.104 +        Parent::nextOut(edge);
   3.105 +        if (edge != INVALID && (*anode_matching)[Parent::aNode(edge)] == edge) {
   3.106 +          Parent::nextOut(edge);
   3.107 +        }
   3.108 +      } else {
   3.109 +        edge = INVALID;
   3.110 +      }
   3.111 +    }
   3.112 + 
   3.113 +    void firstIn(Edge& edge, const Node& node) const {
   3.114 +      if (Parent::aNode(node)) {
   3.115 +        edge = (*bnode_matching)[node] != INVALID ?
   3.116 +          Parent::direct((*anode_matching)[node], false) : INVALID;
   3.117 +      } else {
   3.118 +        Parent::firstIn(edge, node);
   3.119 +        if (edge != INVALID && (*bnode_matching)[node] == edge) {
   3.120 +          Parent::nextIn(edge);
   3.121 +        }
   3.122 +      }
   3.123 +    }
   3.124 +
   3.125 +    void nextIn(Edge& edge) const {
   3.126 +      if (Parent::aNode(Parent::target(edge))) {
   3.127 +        edge = INVALID;
   3.128 +      } else {
   3.129 +        Parent::nextIn(edge);
   3.130 +        if (edge != INVALID && (*bnode_matching)[Parent::bNode(edge)] == edge) {
   3.131 +          Parent::nextIn(edge);
   3.132 +        }
   3.133 +      }
   3.134 +    } 
   3.135 +
   3.136 +  private:
   3.137 +    ANodeMatchingMap* anode_matching;
   3.138 +    BNodeMatchingMap* bnode_matching;
   3.139 +  };
   3.140 +
   3.141 +
   3.142 +  /// \ingroup graph_adaptors
   3.143 +  ///
   3.144 +  /// \brief Bipartite graph adaptor to implement matching algorithms.
   3.145 +  ///
   3.146 +  /// Bipartite graph adaptor to implement matchings. It implements
   3.147 +  /// the residual graph of the matching.
   3.148 +  template <typename _BpUGraph, 
   3.149 +            typename _ANMatchingMap, typename _BNMatchingMap>
   3.150 +  class MatchingBpUGraphAdaptor 
   3.151 +    : public BpUGraphAdaptorExtender<
   3.152 +    MatchingBpUGraphAdaptorBase<_BpUGraph, _ANMatchingMap, _BNMatchingMap> > 
   3.153 +  { 
   3.154 +  public:
   3.155 +
   3.156 +    typedef _BpUGraph BpUGraph;
   3.157 +    typedef _BpUGraph Graph;
   3.158 +    typedef _ANMatchingMap ANodeMatchingMap;
   3.159 +    typedef _BNMatchingMap BNodeMatchingMap;
   3.160 +    typedef BpUGraphAdaptorExtender<
   3.161 +      MatchingBpUGraphAdaptorBase<BpUGraph, 
   3.162 +                                  ANodeMatchingMap, BNodeMatchingMap> > 
   3.163 +    Parent;
   3.164 +
   3.165 +  protected:
   3.166 +    MatchingBpUGraphAdaptor() : Parent() {}
   3.167 +
   3.168 +  public:
   3.169 +
   3.170 +    MatchingBpUGraphAdaptor(const Graph& _graph, 
   3.171 +                            ANodeMatchingMap& _anode_matching,
   3.172 +                            BNodeMatchingMap& _bnode_matching) {
   3.173 +      setGraph(_graph);
   3.174 +      setANodeMatchingMap(_anode_matching);
   3.175 +      setBNodeMatchingMap(_bnode_matching);
   3.176 +    }
   3.177 +
   3.178 +  };
   3.179 +
   3.180  }
   3.181  
   3.182  #endif
     4.1 --- a/test/Makefile.am	Thu Apr 06 09:33:29 2006 +0000
     4.2 +++ b/test/Makefile.am	Fri Apr 07 09:51:23 2006 +0000
     4.3 @@ -15,6 +15,7 @@
     4.4  	all_pairs_shortest_path_test \
     4.5  	arborescence_test \
     4.6  	bfs_test \
     4.7 +	bipartite_matching_test \
     4.8  	counter_test \
     4.9  	dfs_test \
    4.10  	dijkstra_test \
    4.11 @@ -55,6 +56,7 @@
    4.12  all_pairs_shortest_path_test_SOURCES = all_pairs_shortest_path_test.cc
    4.13  arborescence_test_SOURCES = arborescence_test.cc
    4.14  bfs_test_SOURCES = bfs_test.cc
    4.15 +bipartite_matching_test_SOURCES = bipartite_matching_test.cc
    4.16  counter_test_SOURCES = counter_test.cc
    4.17  dfs_test_SOURCES = dfs_test.cc
    4.18  dijkstra_test_SOURCES = dijkstra_test.cc
     5.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     5.2 +++ b/test/bipartite_matching_test.cc	Fri Apr 07 09:51:23 2006 +0000
     5.3 @@ -0,0 +1,138 @@
     5.4 +#include <cstdlib>
     5.5 +#include <iostream>
     5.6 +#include <sstream>
     5.7 +
     5.8 +#include <lemon/list_graph.h>
     5.9 +
    5.10 +#include <lemon/bpugraph_adaptor.h>
    5.11 +#include <lemon/bipartite_matching.h>
    5.12 +
    5.13 +#include <lemon/graph_utils.h>
    5.14 +
    5.15 +#include "test_tools.h"
    5.16 +
    5.17 +using namespace std;
    5.18 +using namespace lemon;
    5.19 +
    5.20 +typedef ListBpUGraph Graph;
    5.21 +BPUGRAPH_TYPEDEFS(Graph);
    5.22 +
    5.23 +int main(int argc, const char *argv[]) {
    5.24 +  Graph graph;
    5.25 +  vector<Node> aNodes;
    5.26 +  vector<Node> bNodes;
    5.27 +  int n = argc > 1 ? atoi(argv[1]) : 100;
    5.28 +  int m = argc > 2 ? atoi(argv[2]) : 100;
    5.29 +  int e = argc > 3 ? atoi(argv[3]) : (int)((n+m) * log(n+m));
    5.30 +
    5.31 +  for (int i = 0; i < n; ++i) {
    5.32 +    Node node = graph.addANode();
    5.33 +    aNodes.push_back(node);
    5.34 +  }
    5.35 +  for (int i = 0; i < m; ++i) {
    5.36 +    Node node = graph.addBNode();
    5.37 +    bNodes.push_back(node);
    5.38 +  }
    5.39 +  for (int i = 0; i < e; ++i) {
    5.40 +    Node aNode = aNodes[urandom(n)];
    5.41 +    Node bNode = bNodes[urandom(m)];
    5.42 +    graph.addEdge(aNode, bNode);
    5.43 +  }
    5.44 +
    5.45 +  {
    5.46 +    MaxBipartiteMatching<Graph> bpmatch(graph);
    5.47 +
    5.48 +    bpmatch.run();
    5.49 +
    5.50 +    Graph::UEdgeMap<bool> mm(graph);
    5.51 +    Graph::NodeMap<bool> cs(graph);
    5.52 +    
    5.53 +    check(bpmatch.coverSet(cs) == bpmatch.matching(mm), "INVALID PRIMAL-DUAL");
    5.54 +    
    5.55 +    for (UEdgeIt it(graph); it != INVALID; ++it) {
    5.56 +      check(cs[graph.aNode(it)] || cs[graph.bNode(it)], "INVALID DUAL");
    5.57 +    }
    5.58 +    
    5.59 +    for (ANodeIt it(graph); it != INVALID; ++it) {
    5.60 +      int num = 0;
    5.61 +      for (IncEdgeIt jt(graph, it); jt != INVALID; ++jt) {
    5.62 +        if (mm[jt]) ++num;
    5.63 +    }
    5.64 +      check(num <= 1, "INVALID PRIMAL");
    5.65 +    }
    5.66 +  }
    5.67 +
    5.68 +  {
    5.69 +    MaxBipartiteMatching<Graph> bpmatch(graph);
    5.70 +
    5.71 +    bpmatch.greedyInit();
    5.72 +    bpmatch.start();
    5.73 +
    5.74 +    Graph::UEdgeMap<bool> mm(graph);
    5.75 +    Graph::NodeMap<bool> cs(graph);
    5.76 +
    5.77 +    check(bpmatch.coverSet(cs) == bpmatch.matching(mm), "INVALID PRIMAL-DUAL");
    5.78 +    
    5.79 +    for (UEdgeIt it(graph); it != INVALID; ++it) {
    5.80 +      check(cs[graph.aNode(it)] || cs[graph.bNode(it)], "INVALID DUAL");
    5.81 +    }
    5.82 +    
    5.83 +    for (ANodeIt it(graph); it != INVALID; ++it) {
    5.84 +      int num = 0;
    5.85 +      for (IncEdgeIt jt(graph, it); jt != INVALID; ++jt) {
    5.86 +        if (mm[jt]) ++num;
    5.87 +    }
    5.88 +      check(num <= 1, "INVALID PRIMAL");
    5.89 +    }
    5.90 +  }
    5.91 +
    5.92 +  {
    5.93 +    MaxBipartiteMatching<Graph> bpmatch(graph);
    5.94 +
    5.95 +    bpmatch.greedyInit();
    5.96 +    while (bpmatch.simpleAugment());
    5.97 +    
    5.98 +    Graph::UEdgeMap<bool> mm(graph);
    5.99 +    Graph::NodeMap<bool> cs(graph);
   5.100 +    
   5.101 +    check(bpmatch.coverSet(cs) == bpmatch.matching(mm), "INVALID PRIMAL-DUAL");
   5.102 +    
   5.103 +    for (UEdgeIt it(graph); it != INVALID; ++it) {
   5.104 +      check(cs[graph.aNode(it)] || cs[graph.bNode(it)], "INVALID DUAL");
   5.105 +    }
   5.106 +    
   5.107 +    for (ANodeIt it(graph); it != INVALID; ++it) {
   5.108 +      int num = 0;
   5.109 +      for (IncEdgeIt jt(graph, it); jt != INVALID; ++jt) {
   5.110 +        if (mm[jt]) ++num;
   5.111 +    }
   5.112 +      check(num <= 1, "INVALID PRIMAL");
   5.113 +    }
   5.114 +  }
   5.115 +
   5.116 +  {
   5.117 +    SwapBpUGraphAdaptor<Graph> swapgraph(graph);
   5.118 +    MaxBipartiteMatching<SwapBpUGraphAdaptor<Graph> > bpmatch(swapgraph);
   5.119 +
   5.120 +    bpmatch.run();
   5.121 +
   5.122 +    Graph::UEdgeMap<bool> mm(graph);
   5.123 +    Graph::NodeMap<bool> cs(graph);
   5.124 +    
   5.125 +    check(bpmatch.coverSet(cs) == bpmatch.matching(mm), "INVALID PRIMAL-DUAL");
   5.126 +    
   5.127 +    for (UEdgeIt it(graph); it != INVALID; ++it) {
   5.128 +      check(cs[graph.aNode(it)] || cs[graph.bNode(it)], "INVALID DUAL");
   5.129 +    }
   5.130 +    
   5.131 +    for (ANodeIt it(graph); it != INVALID; ++it) {
   5.132 +      int num = 0;
   5.133 +      for (IncEdgeIt jt(graph, it); jt != INVALID; ++jt) {
   5.134 +        if (mm[jt]) ++num;
   5.135 +    }
   5.136 +      check(num <= 1, "INVALID PRIMAL");
   5.137 +    }
   5.138 +  }
   5.139 +
   5.140 +  return 0;
   5.141 +}