lemon-0.x: lemon/max_matching.h@9b56cca61e2e

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

    19 #ifndef LEMON_MAX_MATCHING_H

    20 #define LEMON_MAX_MATCHING_H

    22 #include <queue>

    23 #include <lemon/invalid.h>

    24 #include <lemon/unionfind.h>

    25 #include <lemon/graph_utils.h>

    27 ///\ingroup galgs

    28 ///\file

    29 ///\brief Maximum matching algorithm.

    31 namespace lemon {

    33   /// \addtogroup galgs

    34   /// @{

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

    38   ///This class provides Edmonds' alternating forest matching

    39   ///algorithm. The starting matching (if any) can be passed to the

    40   ///algorithm using read-in functions \ref readNMapNode, \ref

    41   ///readNMapEdge or \ref readEMapBool depending on the container. The

    42   ///resulting maximum matching can be attained by write-out functions

    43   ///\ref writeNMapNode, \ref writeNMapEdge or \ref writeEMapBool

    44   ///depending on the preferred container.

    45   ///

    46   ///The dual side of a matching is a map of the nodes to

    47   ///MaxMatching::pos_enum, having values D, A and C showing the

    48   ///Gallai-Edmonds decomposition of the graph. The nodes in D induce

    49   ///a graph with factor-critical components, the nodes in A form the

    50   ///barrier, and the nodes in C induce a graph having a perfect

    51   ///matching. This decomposition can be attained by calling \ref

    52   ///writePos after running the algorithm.

    53   ///

    54   ///\param Graph The undirected graph type the algorithm runs on.

    55   ///

    56   ///\author Jacint Szabo

    57   template <typename Graph>

    58   class MaxMatching {

    60   protected:

    62     typedef typename Graph::Node Node;

    63     typedef typename Graph::Edge Edge;

    64     typedef typename Graph::UEdge UEdge;

    65     typedef typename Graph::UEdgeIt UEdgeIt;

    66     typedef typename Graph::NodeIt NodeIt;

    67     typedef typename Graph::IncEdgeIt IncEdgeIt;

    69     typedef UnionFindEnum<Node, Graph::template NodeMap> 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       for(NodeIt v(g); v!=INVALID; ++v)

   118 	position.set(v,C);

   120       typename Graph::template NodeMap<Node> ear(g,INVALID);

   121       //undefined for the base nodes of the blossoms (i.e. for the

   122       //representative elements of UFE blossom) and for the nodes in C

   124       typename UFE::MapType blossom_base(g);

   125       UFE blossom(blossom_base);

   126       typename UFE::MapType tree_base(g);

   127       UFE tree(tree_base);

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

   129       //blossom_base and tree_base should be a member.

   131       for(NodeIt v(g); v!=INVALID; ++v) {

   132 	if ( position[v]==C && _mate[v]==INVALID ) {

   133 	  blossom.insert(v);

   134 	  tree.insert(v);

   135 	  position.set(v,D);

   136 	  if ( heur == 1 ) lateShrink( v, ear, blossom, tree );

   137 	  else normShrink( v, ear, blossom, tree );

   138 	}

   139       }

   140     }

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

   145     ///Starting form the actual matching stored, it finds a maximal

   146     ///greedy matching.

   147     void greedyMatching() {

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

   149 	if ( _mate[v]==INVALID ) {

   150 	  for( IncEdgeIt e(g,v); e!=INVALID ; ++e ) {

   151 	    Node y=g.runningNode(e);

   152 	    if ( _mate[y]==INVALID && y!=v ) {

   153 	      _mate.set(v,y);

   154 	      _mate.set(y,v);

   155 	      break;

   156 	    }

   157 	  }

   158 	}

   159     }

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

   163     ///Returns the size of the actual matching stored. After \ref

   164     ///run() it returns the size of a maximum matching in the graph.

   165     int size() const {

   166       int s=0;

   167       for(NodeIt v(g); v!=INVALID; ++v) {

   168 	if ( _mate[v]!=INVALID ) {

   169 	  ++s;

   170 	}

   171       }

   172       return s/2;

   173     }

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

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

   179     ///

   180     void resetMatching() {

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

   182 	_mate.set(v,INVALID);

   183     }

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

   187     ///Returns the mate of a \c node in the actual matching.

   188     ///Returns INVALID if the \c node is not covered by the actual matching.

   189     Node mate(Node& node) const {

   190       return _mate[node];

   191     }

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

   195     ///Reads a matching from a \c Node valued \c Node map. This map

   196     ///must be \e symmetric, i.e. if \c map[u]==v then \c map[v]==u

   197     ///must hold, and \c uv will be an edge of the matching.

   198     template<typename NMapN>

   199     void readNMapNode(NMapN& map) {

   200       for(NodeIt v(g); v!=INVALID; ++v) {

   201 	_mate.set(v,map[v]);

   202       }

   203     }

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

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

   208     ///resulting map will be \e symmetric, i.e. if \c map[u]==v then \c

   209     ///map[v]==u will hold, and now \c uv is an edge of the matching.

   210     template<typename NMapN>

   211     void writeNMapNode (NMapN& map) const {

   212       for(NodeIt v(g); v!=INVALID; ++v) {

   213 	map.set(v,_mate[v]);

   214       }

   215     }

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

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

   220     ///map[v] must be an \c UEdge incident to \c v. This map must

   221     ///have the property that if \c g.oppositeNode(u,map[u])==v then

   222     ///\c \c g.oppositeNode(v,map[v])==u holds, and now some edge

   223     ///joining \c u to \c v will be an edge of the matching.

   224     template<typename NMapE>

   225     void readNMapEdge(NMapE& map) {

   226      for(NodeIt v(g); v!=INVALID; ++v) {

   227        UEdge e=map[v];

   228 	if ( e!=INVALID )

   229 	  _mate.set(v,g.oppositeNode(v,e));

   230       }

   231     }

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

   235     ///Writes the stored matching to an \c UEdge valued \c Node

   236     ///map. \c map[v] will be an \c UEdge incident to \c v. This

   237     ///map will have the property that if \c g.oppositeNode(u,map[u])

   238     ///== v then \c map[u]==map[v] holds, and now this edge is an edge

   239     ///of the matching.

   240     template<typename NMapE>

   241     void writeNMapEdge (NMapE& map)  const {

   242       typename Graph::template NodeMap<bool> todo(g,true);

   243       for(NodeIt v(g); v!=INVALID; ++v) {

   244 	if ( todo[v] && _mate[v]!=INVALID ) {

   245 	  Node u=_mate[v];

   246 	  for(IncEdgeIt e(g,v); e!=INVALID; ++e) {

   247 	    if ( g.runningNode(e) == u ) {

   248 	      map.set(u,e);

   249 	      map.set(v,e);

   250 	      todo.set(u,false);

   251 	      todo.set(v,false);

   252 	      break;

   253 	    }

   254 	  }

   255 	}

   256       }

   257     }

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

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

   263     ///must have the property that there are no two incident edges \c

   264     ///e, \c f with \c map[e]==map[f]==true. The edges \c e with \c

   265     ///map[e]==true form the matching.

   266     template<typename EMapB>

   267     void readEMapBool(EMapB& map) {

   268       for(UEdgeIt e(g); e!=INVALID; ++e) {

   269 	if ( map[e] ) {

   270 	  Node u=g.source(e);

   271 	  Node v=g.target(e);

   272 	  _mate.set(u,v);

   273 	  _mate.set(v,u);

   274 	}

   275       }

   276     }

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

   281     ///Writes the matching stored to a \c bool valued \c Edge

   282     ///map. This map will have the property that there are no two

   283     ///incident edges \c e, \c f with \c map[e]==map[f]==true. The

   284     ///edges \c e with \c map[e]==true form the matching.

   285     template<typename EMapB>

   286     void writeEMapBool (EMapB& map) const {

   287       for(UEdgeIt e(g); e!=INVALID; ++e) map.set(e,false);

   289       typename Graph::template NodeMap<bool> todo(g,true);

   290       for(NodeIt v(g); v!=INVALID; ++v) {

   291 	if ( todo[v] && _mate[v]!=INVALID ) {

   292 	  Node u=_mate[v];

   293 	  for(IncEdgeIt e(g,v); e!=INVALID; ++e) {

   294 	    if ( g.runningNode(e) == u ) {

   295 	      map.set(e,true);

   296 	      todo.set(u,false);

   297 	      todo.set(v,false);

   298 	      break;

   299 	    }

   300 	  }

   301 	}

   302       }

   303     }

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

   307     ///the algorithm.

   309     ///After calling any run methods of the class, it writes the

   310     ///Gallai-Edmonds canonical decomposition of the graph. \c map

   311     ///must be a node map of \ref pos_enum 's.

   312     template<typename NMapEnum>

   313     void writePos (NMapEnum& map) const {

   314       for(NodeIt v(g); v!=INVALID; ++v)  map.set(v,position[v]);

   315     }

   317   private:

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

   321 		    UFE& blossom, UFE& tree);

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

   324 		    UFE& blossom, UFE& tree);

   326     bool noShrinkStep(Node x, typename Graph::template NodeMap<Node>& ear,

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

   329     void shrinkStep(Node& top, Node& middle, Node& bottom,

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

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

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

   334 		 UFE& blossom, UFE& tree);

   336   };

   339   // **********************************************************************

   340   //  IMPLEMENTATIONS

   341   // **********************************************************************

   344   template <typename Graph>

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

   346 				      UFE& blossom, UFE& tree) {

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

   349     Q.push(v);

   350     std::queue<Node> R;

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

   353     while ( !Q.empty() ) {

   354       Node x=Q.front();

   355       Q.pop();

   356       if ( noShrinkStep( x, ear, blossom, tree, Q ) ) return;

   357       else R.push(x);

   358     }

   360     while ( !R.empty() ) {

   361       Node x=R.front();

   362       R.pop();

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

   365 	Node y=g.runningNode(e);

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

   368 	  //x and y must be in the same tree

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

   372 	  Node b=blossom.find(x);

   373 	  path.set(b,true);

   374 	  b=_mate[b];

   375 	  while ( b!=INVALID ) {

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

   377 	    path.set(b,true);

   378 	    b=_mate[b];

   379 	  } //going till the root

   381 	  Node top=y;

   382 	  Node middle=blossom.find(top);

   383 	  Node bottom=x;

   384 	  while ( !path[middle] )

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

   387 	  Node base=middle;

   388 	  top=x;

   389 	  middle=blossom.find(top);

   390 	  bottom=y;

   391 	  Node blossom_base=blossom.find(base);

   392 	  while ( middle!=blossom_base )

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

   395 	  blossom.makeRep(base);

   396 	} // if shrink is needed

   398 	while ( !Q.empty() ) {

   399 	  Node x=Q.front();

   400 	  Q.pop();

   401 	  if ( noShrinkStep(x, ear, blossom, tree, Q) ) return;

   402 	  else R.push(x);

   403 	}

   404       } //for e

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

   406   }

   409   template <typename Graph>

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

   411 				      typename Graph::template

   412 				      NodeMap<Node>& ear,

   413 				      UFE& blossom, UFE& tree) {

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

   415     Q.push(v);

   416     while ( !Q.empty() ) {

   418       Node x=Q.front();

   419       Q.pop();

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

   422 	Node y=g.runningNode(e);

   424 	switch ( position[y] ) {

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

   427 	  if ( blossom.find(x) != blossom.find(y) ) { //shrink

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

   430 	    Node b=blossom.find(x);

   431 	    path.set(b,true);

   432 	    b=_mate[b];

   433 	    while ( b!=INVALID ) {

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

   435 	      path.set(b,true);

   436 	      b=_mate[b];

   437 	    } //going till the root

   439 	    Node top=y;

   440 	    Node middle=blossom.find(top);

   441 	    Node bottom=x;

   442 	    while ( !path[middle] )

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

   445 	    Node base=middle;

   446 	    top=x;

   447 	    middle=blossom.find(top);

   448 	    bottom=y;

   449 	    Node blossom_base=blossom.find(base);

   450 	    while ( middle!=blossom_base )

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

   453 	    blossom.makeRep(base);

   454 	  }

   455 	  break;

   456 	case C:

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

   459 	    ear.set(y,x);

   460 	    Node w=_mate[y];

   461 	    blossom.insert(w);

   462 	    position.set(y,A);

   463 	    position.set(w,D);

   464 	    tree.insert(y);

   465 	    tree.insert(w);

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

   467 	    tree.join(w,y);

   468 	    Q.push(w);

   469 	  } else {                 //augment

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

   471 	    _mate.set(x,y);

   472 	    _mate.set(y,x);

   473 	    return;

   474 	  } //if

   475 	  break;

   476 	default: break;

   477 	}

   478       }

   479     }

   480   }

   482   template <typename Graph>

   483   bool MaxMatching<Graph>::noShrinkStep(Node x,

   484 					typename Graph::template

   485 					NodeMap<Node>& ear,

   486 					UFE& blossom, UFE& tree,

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

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

   489       Node y=g.runningNode(e);

   491       if ( position[y]==C ) {

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

   493 	  ear.set(y,x);

   494 	  Node w=_mate[y];

   495 	  blossom.insert(w);

   496 	  position.set(y,A);

   497 	  position.set(w,D);

   498 	  tree.insert(y);

   499 	  tree.insert(w);

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

   501 	  tree.join(w,y);

   502 	  Q.push(w);

   503 	} else {                      //augment

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

   505 	  _mate.set(x,y);

   506 	  _mate.set(y,x);

   507 	  return true;

   508 	}

   509       }

   510     }

   511     return false;

   512   }

   514   template <typename Graph>

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

   516 				      typename Graph::template

   517 				      NodeMap<Node>& ear,

   518 				      UFE& blossom, UFE& tree,

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

   520     ear.set(top,bottom);

   521     Node t=top;

   522     while ( t!=middle ) {

   523       Node u=_mate[t];

   524       t=ear[u];

   525       ear.set(t,u);

   526     }

   527     bottom=_mate[middle];

   528     position.set(bottom,D);

   529     Q.push(bottom);

   530     top=ear[bottom];

   531     Node oldmiddle=middle;

   532     middle=blossom.find(top);

   533     tree.erase(bottom);

   534     tree.erase(oldmiddle);

   535     blossom.insert(bottom);

   536     blossom.join(bottom, oldmiddle);

   537     blossom.join(top, oldmiddle);

   538   }

   540   template <typename Graph>

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

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

   543 				   UFE& blossom, UFE& tree) {

   544     Node v=_mate[x];

   545     while ( v!=INVALID ) {

   547       Node u=ear[v];

   548       _mate.set(v,u);

   549       Node tmp=v;

   550       v=_mate[u];

   551       _mate.set(u,tmp);

   552     }

   553     typename UFE::ItemIt it;

   554     for (tree.first(it,blossom.find(x)); tree.valid(it); tree.next(it)) {

   555       if ( position[it] == D ) {

   556 	typename UFE::ItemIt b_it;

   557 	for (blossom.first(b_it,it); blossom.valid(b_it); blossom.next(b_it)) {

   558 	  position.set( b_it ,C);

   559 	}

   560 	blossom.eraseClass(it);

   561       } else position.set( it ,C);

   562     }

   563     tree.eraseClass(x);

   565   }

   567   /// @}

   569 } //END OF NAMESPACE LEMON

   571 #endif //LEMON_MAX_MATCHING_H

author	deba
	Mon, 27 Feb 2006 10:36:01 +0000
changeset 1986	9b56cca61e2e
parent 1909	2d806130e700
child 1993	2115143eceea
permissions	-rwxr-xr-x