lemon-0.x: lemon/max_matching.h@e4f8367beb41

     1 /* -*- C++ -*-

     2  *

     3  * This file is a part of LEMON, a generic C++ optimization library

     4  *

     5  * Copyright (C) 2003-2007

     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  */

    19 #ifndef LEMON_MAX_MATCHING_H

    20 #define LEMON_MAX_MATCHING_H

    22 #include <queue>

    23 #include <lemon/bits/invalid.h>

    24 #include <lemon/unionfind.h>

    25 #include <lemon/graph_utils.h>

    27 ///\ingroup matching

    28 ///\file

    29 ///\brief Maximum matching algorithm in undirected graph.

    31 namespace lemon {

    33   /// \ingroup matching

    35   ///Edmonds' alternating forest maximum matching algorithm.

    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 {

    59   protected:

    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;

    68     typedef typename Graph::template NodeMap<int> UFECrossRef;

    69     typedef UnionFindEnum<UFECrossRef> UFE;

    71   public:

    73     ///Indicates the Gallai-Edmonds decomposition of the graph.

    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     };

    86   protected:

    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;

    93   public:

    95     MaxMatching(const Graph& _g) : g(_g), _mate(_g,INVALID), position(_g) {}

    97     ///Runs Edmonds' algorithm.

    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     }

   110     ///Runs Edmonds' algorithm.

   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 ) {

   117       //each vertex is put to C

   118       for(NodeIt v(g); v!=INVALID; ++v)

   119 	position.set(v,C);

   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

   125       UFECrossRef blossom_base(g);

   126       UFE blossom(blossom_base);

   128       UFECrossRef tree_base(g);

   129       UFE tree(tree_base);

   131       //If these UFE's would be members of the class then also

   132       //blossom_base and tree_base should be a member.

   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     }

   151     ///Finds a greedy matching starting from the actual matching.

   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     }

   169     ///Returns the size of the actual matching stored.

   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     }

   184     ///Resets the actual matching to the empty matching.

   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     }

   193     ///Returns the mate of a node in the actual matching.

   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     }

   201     ///Reads a matching from a \c Node valued \c Node map.

   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     }

   213     ///Writes the stored matching to a \c Node valued \c Node map.

   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     }

   225     ///Reads a matching from an \c UEdge valued \c Node map.

   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     }

   241     ///Writes the matching stored to an \c UEdge valued \c Node map.

   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     }

   268     ///Reads a matching from a \c bool valued \c Edge map.

   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     }

   287     ///Writes the matching stored to a \c bool valued \c Edge map.

   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);

   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     }

   314     ///Writes the canonical decomposition of the graph after running

   315     ///the algorithm.

   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     }

   325   private:

   328     void lateShrink(Node v, typename Graph::template NodeMap<Node>& ear,

   329 		    UFE& blossom, UFE& tree);

   331     void normShrink(Node v, typename Graph::template NodeMap<Node>& ear,

   332 		    UFE& blossom, UFE& tree);

   334     void shrink(Node x,Node y, typename Graph::template NodeMap<Node>& ear,

   335 		UFE& blossom, UFE& tree,std::queue<Node>& Q);

   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);

   341     bool growOrAugment(Node& y, Node& x, typename Graph::template

   342 		       NodeMap<Node>& ear, UFE& blossom, UFE& tree,

   343 		       std::queue<Node>& Q);

   345     void augment(Node x, typename Graph::template NodeMap<Node>& ear,

   346 		 UFE& blossom, UFE& tree);

   348   };

   351   // **********************************************************************

   352   //  IMPLEMENTATIONS

   353   // **********************************************************************

   356   template <typename Graph>

   357   void MaxMatching<Graph>::lateShrink(Node v, typename Graph::template

   358 				      NodeMap<Node>& ear, UFE& blossom,

   359                                       UFE& tree) {

   360     //We have one tree which we grow, and also shrink but only if it cannot be

   361     //postponed. If we augment then we return to the "for" cycle of

   362     //runEdmonds().

   364     std::queue<Node> Q;   //queue of the totally unscanned nodes

   365     Q.push(v);

   366     std::queue<Node> R;

   367     //queue of the nodes which must be scanned for a possible shrink

   369     while ( !Q.empty() ) {

   370       Node x=Q.front();

   371       Q.pop();

   372       for( IncEdgeIt e(g,x); e!= INVALID; ++e ) {

   373 	Node y=g.runningNode(e);

   374 	//growOrAugment grows if y is covered by the matching and

   375 	//augments if not. In this latter case it returns 1.

   376 	if ( position[y]==C && growOrAugment(y, x, ear, blossom, tree, Q) )

   377           return;

   378       }

   379       R.push(x);

   380     }

   382     while ( !R.empty() ) {

   383       Node x=R.front();

   384       R.pop();

   386       for( IncEdgeIt e(g,x); e!=INVALID ; ++e ) {

   387 	Node y=g.runningNode(e);

   389 	if ( position[y] == D && blossom.find(x) != blossom.find(y) )

   390 	  //Recall that we have only one tree.

   391 	  shrink( x, y, ear, blossom, tree, Q);

   393 	while ( !Q.empty() ) {

   394 	  Node z=Q.front();

   395 	  Q.pop();

   396 	  for( IncEdgeIt f(g,z); f!= INVALID; ++f ) {

   397 	    Node w=g.runningNode(f);

   398 	    //growOrAugment grows if y is covered by the matching and

   399 	    //augments if not. In this latter case it returns 1.

   400 	    if ( position[w]==C && growOrAugment(w, z, ear, blossom, tree, Q) )

   401               return;

   402 	  }

   403 	  R.push(z);

   404 	}

   405       } //for e

   406     } // while ( !R.empty() )

   407   }

   410   template <typename Graph>

   411   void MaxMatching<Graph>::normShrink(Node v,

   412 				      typename Graph::template

   413 				      NodeMap<Node>& ear,

   414 				      UFE& blossom, UFE& tree) {

   415     //We have one tree, which we grow and shrink. If we augment then we

   416     //return to the "for" cycle of runEdmonds().

   418     std::queue<Node> Q;   //queue of the unscanned nodes

   419     Q.push(v);

   420     while ( !Q.empty() ) {

   422       Node x=Q.front();

   423       Q.pop();

   425       for( IncEdgeIt e(g,x); e!=INVALID; ++e ) {

   426 	Node y=g.runningNode(e);

   428 	switch ( position[y] ) {

   429 	case D:          //x and y must be in the same tree

   430 	  if ( blossom.find(x) != blossom.find(y) )

   431 	    //x and y are in the same tree

   432 	    shrink( x, y, ear, blossom, tree, Q);

   433 	  break;

   434 	case C:

   435 	  //growOrAugment grows if y is covered by the matching and

   436 	  //augments if not. In this latter case it returns 1.

   437 	  if ( growOrAugment(y, x, ear, blossom, tree, Q) ) return;

   438 	  break;

   439 	default: break;

   440        	}

   441       }

   442     }

   443   }

   446   template <typename Graph>

   447     void MaxMatching<Graph>::shrink(Node x,Node y, typename

   448 				    Graph::template NodeMap<Node>& ear,

   449 				    UFE& blossom, UFE& tree, std::queue<Node>& Q) {

   450     //x and y are the two adjacent vertices in two blossoms.

   452     typename Graph::template NodeMap<bool> path(g,false);

   454     Node b=blossom.find(x);

   455     path.set(b,true);

   456     b=_mate[b];

   457     while ( b!=INVALID ) {

   458       b=blossom.find(ear[b]);

   459       path.set(b,true);

   460       b=_mate[b];

   461     } //we go until the root through bases of blossoms and odd vertices

   463     Node top=y;

   464     Node middle=blossom.find(top);

   465     Node bottom=x;

   466     while ( !path[middle] )

   467       shrinkStep(top, middle, bottom, ear, blossom, tree, Q);

   468     //Until we arrive to a node on the path, we update blossom, tree

   469     //and the positions of the odd nodes.

   471     Node base=middle;

   472     top=x;

   473     middle=blossom.find(top);

   474     bottom=y;

   475     Node blossom_base=blossom.find(base);

   476     while ( middle!=blossom_base )

   477       shrinkStep(top, middle, bottom, ear, blossom, tree, Q);

   478     //Until we arrive to a node on the path, we update blossom, tree

   479     //and the positions of the odd nodes.

   481     blossom.makeRep(base);

   482   }

   486   template <typename Graph>

   487   void MaxMatching<Graph>::shrinkStep(Node& top, Node& middle, Node& bottom,

   488 				      typename Graph::template

   489 				      NodeMap<Node>& ear,

   490 				      UFE& blossom, UFE& tree,

   491 				      std::queue<Node>& Q) {

   492     //We traverse a blossom and update everything.

   494     ear.set(top,bottom);

   495     Node t=top;

   496     while ( t!=middle ) {

   497       Node u=_mate[t];

   498       t=ear[u];

   499       ear.set(t,u);

   500     }

   501     bottom=_mate[middle];

   502     position.set(bottom,D);

   503     Q.push(bottom);

   504     top=ear[bottom];

   505     Node oldmiddle=middle;

   506     middle=blossom.find(top);

   507     tree.erase(bottom);

   508     tree.erase(oldmiddle);

   509     blossom.insert(bottom);

   510     blossom.join(bottom, oldmiddle);

   511     blossom.join(top, oldmiddle);

   512   }

   515   template <typename Graph>

   516   bool MaxMatching<Graph>::growOrAugment(Node& y, Node& x, typename Graph::template

   517 					 NodeMap<Node>& ear, UFE& blossom, UFE& tree,

   518 					 std::queue<Node>& Q) {

   519     //x is in a blossom in the tree, y is outside. If y is covered by

   520     //the matching we grow, otherwise we augment. In this case we

   521     //return 1.

   523     if ( _mate[y]!=INVALID ) {       //grow

   524       ear.set(y,x);

   525       Node w=_mate[y];

   526       blossom.insert(w);

   527       position.set(y,A);

   528       position.set(w,D);

   529       tree.insert(y);

   530       tree.insert(w);

   531       tree.join(y,blossom.find(x));

   532       tree.join(w,y);

   533       Q.push(w);

   534     } else {                      //augment

   535       augment(x, ear, blossom, tree);

   536       _mate.set(x,y);

   537       _mate.set(y,x);

   538       return true;

   539     }

   540     return false;

   541   }

   544   template <typename Graph>

   545   void MaxMatching<Graph>::augment(Node x,

   546 				   typename Graph::template NodeMap<Node>& ear,

   547 				   UFE& blossom, UFE& tree) {

   548     Node v=_mate[x];

   549     while ( v!=INVALID ) {

   551       Node u=ear[v];

   552       _mate.set(v,u);

   553       Node tmp=v;

   554       v=_mate[u];

   555       _mate.set(u,tmp);

   556     }

   557     Node y=blossom.find(x);

   558     for (typename UFE::ItemIt tit(tree, y); tit != INVALID; ++tit) {

   559       if ( position[tit] == D ) {

   560 	for (typename UFE::ItemIt bit(blossom, tit); bit != INVALID; ++bit) {

   561 	  position.set( bit ,C);

   562 	}

   563 	blossom.eraseClass(tit);

   564       } else position.set( tit ,C);

   565     }

   566     tree.eraseClass(y);

   568   }

   571 } //END OF NAMESPACE LEMON

   573 #endif //LEMON_MAX_MATCHING_H

author	ladanyi
	Sat, 13 Oct 2007 08:48:07 +0000
changeset 2495	e4f8367beb41
parent 2386	81b47fc5c444
child 2505	1bb471764ab8
permissions	-rwxr-xr-x