COIN-OR::LEMON - Graph Library

source: lemon-0.x/src/lemon/max_matching.h @ 1165:c5e56125959a

Last change on this file since 1165:c5e56125959a was 1165:c5e56125959a, checked in by jacint, 19 years ago

some minor changes, docs, etc.

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