COIN-OR::LEMON - Graph Library

source: lemon-0.x/lemon/max_matching.h @ 2229:4dbb6dd2dd4b

Last change on this file since 2229:4dbb6dd2dd4b was 2205:c20b0eb92a33, checked in by Balazs Dezso, 13 years ago

UnionFind?
Changing the representation of the union-find
it has the same running time but it takes just 2/3 space
! does not auto insert items /performance/

UnionFindEnum?
Changing the interface - more convenient to UnionFind?
Does not based on the stl data structures /it could be disadvantage/

=> does not use singular iterator assignment /not stl conform, but always work/

Just new iterator interface

MaxMatching? + UnionFindTest?
Using new iterator interface instead of the old

  • Property exe set to *
File size: 16.8 KB
Line 
1/* -*- C++ -*-
2 *
3 * This file is a part of LEMON, a generic C++ optimization library
4 *
5 * Copyright (C) 2003-2006
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
8 *
9 * Permission to use, modify and distribute this software is granted
10 * provided that this copyright notice appears in all copies. For
11 * precise terms see the accompanying LICENSE file.
12 *
13 * This software is provided "AS IS" with no warranty of any kind,
14 * express or implied, and with no claim as to its suitability for any
15 * purpose.
16 *
17 */
18
19#ifndef LEMON_MAX_MATCHING_H
20#define LEMON_MAX_MATCHING_H
21
22#include <queue>
23#include <lemon/bits/invalid.h>
24#include <lemon/unionfind.h>
25#include <lemon/graph_utils.h>
26
27///\ingroup matching
28///\file
29///\brief Maximum matching algorithm in undirected graph.
30
31namespace lemon {
32
33  /// \ingroup matching
34
35  ///Edmonds' alternating forest maximum matching algorithm.
36
37  ///This class provides Edmonds' alternating forest matching
38  ///algorithm. The starting matching (if any) can be passed to the
39  ///algorithm using read-in functions \ref readNMapNode, \ref
40  ///readNMapEdge or \ref readEMapBool depending on the container. The
41  ///resulting maximum matching can be attained by write-out functions
42  ///\ref writeNMapNode, \ref writeNMapEdge or \ref writeEMapBool
43  ///depending on the preferred container.
44  ///
45  ///The dual side of a matching is a map of the nodes to
46  ///MaxMatching::pos_enum, having values D, A and C showing the
47  ///Gallai-Edmonds decomposition of the graph. The nodes in D induce
48  ///a graph with factor-critical components, the nodes in A form the
49  ///barrier, and the nodes in C induce a graph having a perfect
50  ///matching. This decomposition can be attained by calling \ref
51  ///writePos after running the algorithm.
52  ///
53  ///\param Graph The undirected graph type the algorithm runs on.
54  ///
55  ///\author Jacint Szabo 
56  template <typename Graph>
57  class MaxMatching {
58
59  protected:
60
61    typedef typename Graph::Node Node;
62    typedef typename Graph::Edge Edge;
63    typedef typename Graph::UEdge UEdge;
64    typedef typename Graph::UEdgeIt UEdgeIt;
65    typedef typename Graph::NodeIt NodeIt;
66    typedef typename Graph::IncEdgeIt IncEdgeIt;
67
68    typedef typename Graph::template NodeMap<int> UFECrossRef;
69    typedef UnionFindEnum<Node, UFECrossRef> UFE;
70
71  public:
72   
73    ///Indicates the Gallai-Edmonds decomposition of the graph.
74
75    ///Indicates the Gallai-Edmonds decomposition of the graph, which
76    ///shows an upper bound on the size of a maximum matching. The
77    ///nodes with pos_enum \c D induce a graph with factor-critical
78    ///components, the nodes in \c A form the canonical barrier, and the
79    ///nodes in \c C induce a graph having a perfect matching.
80    enum pos_enum {
81      D=0,
82      A=1,
83      C=2
84    };
85
86  protected:
87
88    static const int HEUR_density=2;
89    const Graph& g;
90    typename Graph::template NodeMap<Node> _mate;
91    typename Graph::template NodeMap<pos_enum> position;
92     
93  public:
94   
95    MaxMatching(const Graph& _g) : g(_g), _mate(_g,INVALID), position(_g) {}
96
97    ///Runs Edmonds' algorithm.
98
99    ///Runs Edmonds' algorithm for sparse graphs (number of edges <
100    ///2*number of nodes), and a heuristical Edmonds' algorithm with a
101    ///heuristic of postponing shrinks for dense graphs.
102    void run() {
103      if ( countUEdges(g) < HEUR_density*countNodes(g) ) {
104        greedyMatching();
105        runEdmonds(0);
106      } else runEdmonds(1);
107    }
108
109
110    ///Runs Edmonds' algorithm.
111   
112    ///If heur=0 it runs Edmonds' algorithm. If heur=1 it runs
113    ///Edmonds' algorithm with a heuristic of postponing shrinks,
114    ///giving a faster algorithm for dense graphs. 
115    void runEdmonds( int heur = 1 ) {
116     
117      //each vertex is put to C
118      for(NodeIt v(g); v!=INVALID; ++v)
119        position.set(v,C);     
120     
121      typename Graph::template NodeMap<Node> ear(g,INVALID);
122      //undefined for the base nodes of the blossoms (i.e. for the
123      //representative elements of UFE blossom) and for the nodes in C
124     
125      UFECrossRef blossom_base(g);
126      UFE blossom(blossom_base);
127
128      UFECrossRef tree_base(g);
129      UFE tree(tree_base);
130
131      //If these UFE's would be members of the class then also
132      //blossom_base and tree_base should be a member.
133     
134      //We build only one tree and the other vertices uncovered by the
135      //matching belong to C. (They can be considered as singleton
136      //trees.) If this tree can be augmented or no more
137      //grow/augmentation/shrink is possible then we return to this
138      //"for" cycle.
139      for(NodeIt v(g); v!=INVALID; ++v) {
140        if ( position[v]==C && _mate[v]==INVALID ) {
141          blossom.insert(v);
142          tree.insert(v);
143          position.set(v,D);
144          if ( heur == 1 ) lateShrink( v, ear, blossom, tree );
145          else normShrink( v, ear, blossom, tree );
146        }
147      }
148    }
149
150
151    ///Finds a greedy matching starting from the actual matching.
152   
153    ///Starting form the actual matching stored, it finds a maximal
154    ///greedy matching.
155    void greedyMatching() {
156      for(NodeIt v(g); v!=INVALID; ++v)
157        if ( _mate[v]==INVALID ) {
158          for( IncEdgeIt e(g,v); e!=INVALID ; ++e ) {
159            Node y=g.runningNode(e);
160            if ( _mate[y]==INVALID && y!=v ) {
161              _mate.set(v,y);
162              _mate.set(y,v);
163              break;
164            }
165          }
166        }
167    }
168
169    ///Returns the size of the actual matching stored.
170
171    ///Returns the size of the actual matching stored. After \ref
172    ///run() it returns the size of a maximum matching in the graph.
173    int size() const {
174      int s=0;
175      for(NodeIt v(g); v!=INVALID; ++v) {
176        if ( _mate[v]!=INVALID ) {
177          ++s;
178        }
179      }
180      return s/2;
181    }
182
183
184    ///Resets the actual matching to the empty matching.
185
186    ///Resets the actual matching to the empty matching. 
187    ///
188    void resetMatching() {
189      for(NodeIt v(g); v!=INVALID; ++v)
190        _mate.set(v,INVALID);     
191    }
192
193    ///Returns the mate of a node in the actual matching.
194
195    ///Returns the mate of a \c node in the actual matching.
196    ///Returns INVALID if the \c node is not covered by the actual matching.
197    Node mate(Node& node) const {
198      return _mate[node];
199    }
200
201    ///Reads a matching from a \c Node valued \c Node map.
202
203    ///Reads a matching from a \c Node valued \c Node map. This map
204    ///must be \e symmetric, i.e. if \c map[u]==v then \c map[v]==u
205    ///must hold, and \c uv will be an edge of the matching.
206    template<typename NMapN>
207    void readNMapNode(NMapN& map) {
208      for(NodeIt v(g); v!=INVALID; ++v) {
209        _mate.set(v,map[v]);   
210      }
211    }
212   
213    ///Writes the stored matching to a \c Node valued \c Node map.
214
215    ///Writes the stored matching to a \c Node valued \c Node map. The
216    ///resulting map will be \e symmetric, i.e. if \c map[u]==v then \c
217    ///map[v]==u will hold, and now \c uv is an edge of the matching.
218    template<typename NMapN>
219    void writeNMapNode (NMapN& map) const {
220      for(NodeIt v(g); v!=INVALID; ++v) {
221        map.set(v,_mate[v]);   
222      }
223    }
224
225    ///Reads a matching from an \c UEdge valued \c Node map.
226
227    ///Reads a matching from an \c UEdge valued \c Node map. \c
228    ///map[v] must be an \c UEdge incident to \c v. This map must
229    ///have the property that if \c g.oppositeNode(u,map[u])==v then
230    ///\c \c g.oppositeNode(v,map[v])==u holds, and now some edge
231    ///joining \c u to \c v will be an edge of the matching.
232    template<typename NMapE>
233    void readNMapEdge(NMapE& map) {
234      for(NodeIt v(g); v!=INVALID; ++v) {
235        UEdge e=map[v];
236        if ( e!=INVALID )
237          _mate.set(v,g.oppositeNode(v,e));
238      }
239    }
240   
241    ///Writes the matching stored to an \c UEdge valued \c Node map.
242
243    ///Writes the stored matching to an \c UEdge valued \c Node
244    ///map. \c map[v] will be an \c UEdge incident to \c v. This
245    ///map will have the property that if \c g.oppositeNode(u,map[u])
246    ///== v then \c map[u]==map[v] holds, and now this edge is an edge
247    ///of the matching.
248    template<typename NMapE>
249    void writeNMapEdge (NMapE& map)  const {
250      typename Graph::template NodeMap<bool> todo(g,true);
251      for(NodeIt v(g); v!=INVALID; ++v) {
252        if ( todo[v] && _mate[v]!=INVALID ) {
253          Node u=_mate[v];
254          for(IncEdgeIt e(g,v); e!=INVALID; ++e) {
255            if ( g.runningNode(e) == u ) {
256              map.set(u,e);
257              map.set(v,e);
258              todo.set(u,false);
259              todo.set(v,false);
260              break;
261            }
262          }
263        }
264      }
265    }
266
267
268    ///Reads a matching from a \c bool valued \c Edge map.
269   
270    ///Reads a matching from a \c bool valued \c Edge map. This map
271    ///must have the property that there are no two incident edges \c
272    ///e, \c f with \c map[e]==map[f]==true. The edges \c e with \c
273    ///map[e]==true form the matching.
274    template<typename EMapB>
275    void readEMapBool(EMapB& map) {
276      for(UEdgeIt e(g); e!=INVALID; ++e) {
277        if ( map[e] ) {
278          Node u=g.source(e);     
279          Node v=g.target(e);
280          _mate.set(u,v);
281          _mate.set(v,u);
282        }
283      }
284    }
285
286
287    ///Writes the matching stored to a \c bool valued \c Edge map.
288
289    ///Writes the matching stored to a \c bool valued \c Edge
290    ///map. This map will have the property that there are no two
291    ///incident edges \c e, \c f with \c map[e]==map[f]==true. The
292    ///edges \c e with \c map[e]==true form the matching.
293    template<typename EMapB>
294    void writeEMapBool (EMapB& map) const {
295      for(UEdgeIt e(g); e!=INVALID; ++e) map.set(e,false);
296
297      typename Graph::template NodeMap<bool> todo(g,true);
298      for(NodeIt v(g); v!=INVALID; ++v) {
299        if ( todo[v] && _mate[v]!=INVALID ) {
300          Node u=_mate[v];
301          for(IncEdgeIt e(g,v); e!=INVALID; ++e) {
302            if ( g.runningNode(e) == u ) {
303              map.set(e,true);
304              todo.set(u,false);
305              todo.set(v,false);
306              break;
307            }
308          }
309        }
310      }
311    }
312
313
314    ///Writes the canonical decomposition of the graph after running
315    ///the algorithm.
316
317    ///After calling any run methods of the class, it writes the
318    ///Gallai-Edmonds canonical decomposition of the graph. \c map
319    ///must be a node map of \ref pos_enum 's.
320    template<typename NMapEnum>
321    void writePos (NMapEnum& map) const {
322      for(NodeIt v(g); v!=INVALID; ++v)  map.set(v,position[v]);
323    }
324
325  private:
326
327 
328    void lateShrink(Node v, typename Graph::template NodeMap<Node>& ear, 
329                    UFE& blossom, UFE& tree);
330
331    void normShrink(Node v, typename Graph::template NodeMap<Node>& ear, 
332                    UFE& blossom, UFE& tree);
333
334    void shrink(Node x,Node y, typename Graph::template NodeMap<Node>& ear, 
335                UFE& blossom, UFE& tree,std::queue<Node>& Q);
336
337    void shrinkStep(Node& top, Node& middle, Node& bottom,
338                    typename Graph::template NodeMap<Node>& ear, 
339                    UFE& blossom, UFE& tree, std::queue<Node>& Q);
340
341    bool growOrAugment(Node& y, Node& x, typename Graph::template
342                       NodeMap<Node>& ear, UFE& blossom, UFE& tree,
343                       std::queue<Node>& Q);
344
345    void augment(Node x, typename Graph::template NodeMap<Node>& ear, 
346                 UFE& blossom, UFE& tree);
347
348  };
349
350
351  // **********************************************************************
352  //  IMPLEMENTATIONS
353  // **********************************************************************
354
355
356  template <typename Graph>
357  void MaxMatching<Graph>::lateShrink(Node v, typename Graph::template
358                                      NodeMap<Node>& ear, UFE& blossom, UFE& tree) {
359    //We have one tree which we grow, and also shrink but only if it cannot be
360    //postponed. If we augment then we return to the "for" cycle of
361    //runEdmonds().
362
363    std::queue<Node> Q;   //queue of the totally unscanned nodes
364    Q.push(v); 
365    std::queue<Node> R;   
366    //queue of the nodes which must be scanned for a possible shrink
367     
368    while ( !Q.empty() ) {
369      Node x=Q.front();
370      Q.pop();
371      for( IncEdgeIt e(g,x); e!= INVALID; ++e ) {
372        Node y=g.runningNode(e);
373        //growOrAugment grows if y is covered by the matching and
374        //augments if not. In this latter case it returns 1.
375        if ( position[y]==C && growOrAugment(y, x, ear, blossom, tree, Q) ) return;
376      }
377      R.push(x);
378    }
379     
380    while ( !R.empty() ) {
381      Node x=R.front();
382      R.pop();
383       
384      for( IncEdgeIt e(g,x); e!=INVALID ; ++e ) {
385        Node y=g.runningNode(e);
386
387        if ( position[y] == D && blossom.find(x) != blossom.find(y) ) 
388          //Recall that we have only one tree.
389          shrink( x, y, ear, blossom, tree, Q);
390       
391        while ( !Q.empty() ) {
392          Node x=Q.front();
393          Q.pop();
394          for( IncEdgeIt e(g,x); e!= INVALID; ++e ) {
395            Node y=g.runningNode(e);
396            //growOrAugment grows if y is covered by the matching and
397            //augments if not. In this latter case it returns 1.
398            if ( position[y]==C && growOrAugment(y, x, ear, blossom, tree, Q) ) return;
399          }
400          R.push(x);
401        }
402      } //for e
403    } // while ( !R.empty() )
404  }
405
406
407  template <typename Graph>
408  void MaxMatching<Graph>::normShrink(Node v,
409                                      typename Graph::template
410                                      NodeMap<Node>& ear, 
411                                      UFE& blossom, UFE& tree) {
412    //We have one tree, which we grow and shrink. If we augment then we
413    //return to the "for" cycle of runEdmonds().
414   
415    std::queue<Node> Q;   //queue of the unscanned nodes
416    Q.push(v); 
417    while ( !Q.empty() ) {
418
419      Node x=Q.front();
420      Q.pop();
421       
422      for( IncEdgeIt e(g,x); e!=INVALID; ++e ) {
423        Node y=g.runningNode(e);
424             
425        switch ( position[y] ) {
426        case D:          //x and y must be in the same tree
427          if ( blossom.find(x) != blossom.find(y) )
428            //x and y are in the same tree
429            shrink( x, y, ear, blossom, tree, Q);
430          break;
431        case C:
432          //growOrAugment grows if y is covered by the matching and
433          //augments if not. In this latter case it returns 1.
434          if ( growOrAugment(y, x, ear, blossom, tree, Q) ) return;
435          break;
436        default: break;
437        }
438      }
439    }
440  }
441 
442
443  template <typename Graph>
444    void MaxMatching<Graph>::shrink(Node x,Node y, typename
445                                    Graph::template NodeMap<Node>& ear, 
446                                    UFE& blossom, UFE& tree, std::queue<Node>& Q) {
447    //x and y are the two adjacent vertices in two blossoms.
448   
449    typename Graph::template NodeMap<bool> path(g,false);
450   
451    Node b=blossom.find(x);
452    path.set(b,true);
453    b=_mate[b];
454    while ( b!=INVALID ) {
455      b=blossom.find(ear[b]);
456      path.set(b,true);
457      b=_mate[b];
458    } //we go until the root through bases of blossoms and odd vertices
459   
460    Node top=y;
461    Node middle=blossom.find(top);
462    Node bottom=x;
463    while ( !path[middle] )
464      shrinkStep(top, middle, bottom, ear, blossom, tree, Q);
465    //Until we arrive to a node on the path, we update blossom, tree
466    //and the positions of the odd nodes.
467   
468    Node base=middle;
469    top=x;
470    middle=blossom.find(top);
471    bottom=y;
472    Node blossom_base=blossom.find(base);
473    while ( middle!=blossom_base )
474      shrinkStep(top, middle, bottom, ear, blossom, tree, Q);
475    //Until we arrive to a node on the path, we update blossom, tree
476    //and the positions of the odd nodes.
477   
478    blossom.makeRep(base);
479  }
480
481
482
483  template <typename Graph>
484  void MaxMatching<Graph>::shrinkStep(Node& top, Node& middle, Node& bottom,
485                                      typename Graph::template
486                                      NodeMap<Node>& ear, 
487                                      UFE& blossom, UFE& tree,
488                                      std::queue<Node>& Q) {
489    //We traverse a blossom and update everything.
490   
491    ear.set(top,bottom);
492    Node t=top;
493    while ( t!=middle ) {
494      Node u=_mate[t];
495      t=ear[u];
496      ear.set(t,u);
497    }
498    bottom=_mate[middle];
499    position.set(bottom,D);
500    Q.push(bottom);
501    top=ear[bottom];           
502    Node oldmiddle=middle;
503    middle=blossom.find(top);
504    tree.erase(bottom);
505    tree.erase(oldmiddle);
506    blossom.insert(bottom);
507    blossom.join(bottom, oldmiddle);
508    blossom.join(top, oldmiddle);
509  }
510
511
512  template <typename Graph>
513  bool MaxMatching<Graph>::growOrAugment(Node& y, Node& x, typename Graph::template
514                                         NodeMap<Node>& ear, UFE& blossom, UFE& tree,
515                                         std::queue<Node>& Q) {
516    //x is in a blossom in the tree, y is outside. If y is covered by
517    //the matching we grow, otherwise we augment. In this case we
518    //return 1.
519   
520    if ( _mate[y]!=INVALID ) {       //grow
521      ear.set(y,x);
522      Node w=_mate[y];
523      blossom.insert(w);
524      position.set(y,A);
525      position.set(w,D);
526      tree.insert(y);
527      tree.insert(w);
528      tree.join(y,blossom.find(x)); 
529      tree.join(w,y); 
530      Q.push(w);
531    } else {                      //augment
532      augment(x, ear, blossom, tree);
533      _mate.set(x,y);
534      _mate.set(y,x);
535      return true;
536    }
537    return false;
538  }
539 
540
541  template <typename Graph>
542  void MaxMatching<Graph>::augment(Node x,
543                                   typename Graph::template NodeMap<Node>& ear, 
544                                   UFE& blossom, UFE& tree) {
545    Node v=_mate[x];
546    while ( v!=INVALID ) {
547       
548      Node u=ear[v];
549      _mate.set(v,u);
550      Node tmp=v;
551      v=_mate[u];
552      _mate.set(u,tmp);
553    }
554    Node y=blossom.find(x);
555    for (typename UFE::ItemIt tit(tree, y); tit != INVALID; ++tit) {   
556      if ( position[tit] == D ) {
557        for (typename UFE::ItemIt bit(blossom, tit); bit != INVALID; ++bit) { 
558          position.set( bit ,C);
559        }
560        blossom.eraseClass(tit);
561      } else position.set( tit ,C);
562    }
563    tree.eraseClass(y);
564
565  }
566
567 
568} //END OF NAMESPACE LEMON
569
570#endif //LEMON_MAX_MATCHING_H
Note: See TracBrowser for help on using the repository browser.