Changeset 2505:1bb471764ab8 in lemon-0.x
- Timestamp:
- 10/30/07 21:21:10 (17 years ago)
- Branch:
- default
- Phase:
- public
- Convert:
- svn:c9d7d8f5-90d6-0310-b91f-818b3a526b0e/lemon/trunk@3351
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
lemon/max_matching.h
r2391 r2505 31 31 namespace lemon { 32 32 33 /// 34 35 /// Edmonds' alternating forest maximum matching algorithm.36 33 ///\ingroup matching 34 35 ///\brief Edmonds' alternating forest maximum matching algorithm. 36 /// 37 37 ///This class provides Edmonds' alternating forest matching 38 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. 39 ///algorithm using some of init functions. 44 40 /// 45 41 ///The dual side of a matching is a map of the nodes to 46 ///MaxMatching:: pos_enum, having values D, A and C showing the47 /// Gallai-Edmonds decomposition of the graph. The nodes in D induce48 /// a graph with factor-critical components, the nodes in A form the49 /// barrier, and the nodes in C induce a graph having a perfect50 /// matching. This decomposition can be attained by calling \ref51 /// writePos after running the algorithm.42 ///MaxMatching::DecompType, having values \c D, \c A and \c C 43 ///showing the Gallai-Edmonds decomposition of the graph. The nodes 44 ///in \c D induce a graph with factor-critical components, the nodes 45 ///in \c A form the barrier, and the nodes in \c C induce a graph 46 ///having a perfect matching. This decomposition can be attained by 47 ///calling \c decomposition() after running the algorithm. 52 48 /// 53 49 ///\param Graph The undirected graph type the algorithm runs on. … … 68 64 typedef typename Graph::template NodeMap<int> UFECrossRef; 69 65 typedef UnionFindEnum<UFECrossRef> UFE; 66 typedef std::vector<Node> NV; 67 68 typedef typename Graph::template NodeMap<int> EFECrossRef; 69 typedef ExtendFindEnum<EFECrossRef> EFE; 70 70 71 71 public: 72 72 73 /// Indicates the Gallai-Edmonds decomposition of the graph.74 73 ///\brief Indicates the Gallai-Edmonds decomposition of the graph. 74 /// 75 75 ///Indicates the Gallai-Edmonds decomposition of the graph, which 76 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-critical77 ///nodes with DecompType \c D induce a graph with factor-critical 78 78 ///components, the nodes in \c A form the canonical barrier, and the 79 79 ///nodes in \c C induce a graph having a perfect matching. 80 enum pos_enum{80 enum DecompType { 81 81 D=0, 82 82 A=1, … … 89 89 const Graph& g; 90 90 typename Graph::template NodeMap<Node> _mate; 91 typename Graph::template NodeMap< pos_enum> position;91 typename Graph::template NodeMap<DecompType> position; 92 92 93 93 public: 94 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() { 95 MaxMatching(const Graph& _g) 96 : g(_g), _mate(_g), position(_g) {} 97 98 ///\brief Sets the actual matching to the empty matching. 99 /// 100 ///Sets the actual matching to the empty matching. 101 /// 102 void init() { 103 for(NodeIt v(g); v!=INVALID; ++v) { 104 _mate.set(v,INVALID); 105 position.set(v,C); 106 } 107 } 108 109 ///\brief Finds a greedy matching for initial matching. 110 /// 111 ///For initial matchig it finds a maximal greedy matching. 112 void greedyInit() { 113 for(NodeIt v(g); v!=INVALID; ++v) { 114 _mate.set(v,INVALID); 115 position.set(v,C); 116 } 156 117 for(NodeIt v(g); v!=INVALID; ++v) 157 118 if ( _mate[v]==INVALID ) { … … 167 128 } 168 129 169 ///Returns the size of the actual matching stored. 170 130 ///\brief Initialize the matching from each nodes' mate. 131 /// 132 ///Initialize the matching from a \c Node valued \c Node map. This 133 ///map must be \e symmetric, i.e. if \c map[u]==v then \c 134 ///map[v]==u must hold, and \c uv will be an edge of the initial 135 ///matching. 136 template <typename MateMap> 137 void mateMapInit(MateMap& map) { 138 for(NodeIt v(g); v!=INVALID; ++v) { 139 _mate.set(v,map[v]); 140 position.set(v,C); 141 } 142 } 143 144 ///\brief Initialize the matching from a node map with the 145 ///incident matching edges. 146 /// 147 ///Initialize the matching from an \c UEdge valued \c Node map. \c 148 ///map[v] must be an \c UEdge incident to \c v. This map must have 149 ///the property that if \c g.oppositeNode(u,map[u])==v then \c \c 150 ///g.oppositeNode(v,map[v])==u holds, and now some edge joining \c 151 ///u to \c v will be an edge of the matching. 152 template<typename MatchingMap> 153 void matchingMapInit(MatchingMap& map) { 154 for(NodeIt v(g); v!=INVALID; ++v) { 155 position.set(v,C); 156 UEdge e=map[v]; 157 if ( e!=INVALID ) 158 _mate.set(v,g.oppositeNode(v,e)); 159 else 160 _mate.set(v,INVALID); 161 } 162 } 163 164 ///\brief Initialize the matching from the map containing the 165 ///undirected matching edges. 166 /// 167 ///Initialize the matching from a \c bool valued \c UEdge map. This 168 ///map must have the property that there are no two incident edges 169 ///\c e, \c f with \c map[e]==map[f]==true. The edges \c e with \c 170 ///map[e]==true form the matching. 171 template <typename MatchingMap> 172 void matchingInit(MatchingMap& map) { 173 for(NodeIt v(g); v!=INVALID; ++v) { 174 _mate.set(v,INVALID); 175 position.set(v,C); 176 } 177 for(UEdgeIt e(g); e!=INVALID; ++e) { 178 if ( map[e] ) { 179 Node u=g.source(e); 180 Node v=g.target(e); 181 _mate.set(u,v); 182 _mate.set(v,u); 183 } 184 } 185 } 186 187 188 ///\brief Runs Edmonds' algorithm. 189 /// 190 ///Runs Edmonds' algorithm for sparse graphs (number of edges < 191 ///2*number of nodes), and a heuristical Edmonds' algorithm with a 192 ///heuristic of postponing shrinks for dense graphs. 193 void run() { 194 if (countUEdges(g) < HEUR_density * countNodes(g)) { 195 greedyInit(); 196 startSparse(); 197 } else { 198 init(); 199 startDense(); 200 } 201 } 202 203 204 ///\brief Starts Edmonds' algorithm. 205 /// 206 ///If runs the original Edmonds' algorithm. 207 void startSparse() { 208 209 typename Graph::template NodeMap<Node> ear(g,INVALID); 210 //undefined for the base nodes of the blossoms (i.e. for the 211 //representative elements of UFE blossom) and for the nodes in C 212 213 UFECrossRef blossom_base(g); 214 UFE blossom(blossom_base); 215 NV rep(countNodes(g)); 216 217 EFECrossRef tree_base(g); 218 EFE tree(tree_base); 219 220 //If these UFE's would be members of the class then also 221 //blossom_base and tree_base should be a member. 222 223 //We build only one tree and the other vertices uncovered by the 224 //matching belong to C. (They can be considered as singleton 225 //trees.) If this tree can be augmented or no more 226 //grow/augmentation/shrink is possible then we return to this 227 //"for" cycle. 228 for(NodeIt v(g); v!=INVALID; ++v) { 229 if (position[v]==C && _mate[v]==INVALID) { 230 rep[blossom.insert(v)] = v; 231 tree.insert(v); 232 position.set(v,D); 233 normShrink(v, ear, blossom, rep, tree); 234 } 235 } 236 } 237 238 ///\brief Starts Edmonds' algorithm. 239 /// 240 ///It runs Edmonds' algorithm with a heuristic of postponing 241 ///shrinks, giving a faster algorithm for dense graphs. 242 void startDense() { 243 244 typename Graph::template NodeMap<Node> ear(g,INVALID); 245 //undefined for the base nodes of the blossoms (i.e. for the 246 //representative elements of UFE blossom) and for the nodes in C 247 248 UFECrossRef blossom_base(g); 249 UFE blossom(blossom_base); 250 NV rep(countNodes(g)); 251 252 EFECrossRef tree_base(g); 253 EFE tree(tree_base); 254 255 //If these UFE's would be members of the class then also 256 //blossom_base and tree_base should be a member. 257 258 //We build only one tree and the other vertices uncovered by the 259 //matching belong to C. (They can be considered as singleton 260 //trees.) If this tree can be augmented or no more 261 //grow/augmentation/shrink is possible then we return to this 262 //"for" cycle. 263 for(NodeIt v(g); v!=INVALID; ++v) { 264 if ( position[v]==C && _mate[v]==INVALID ) { 265 rep[blossom.insert(v)] = v; 266 tree.insert(v); 267 position.set(v,D); 268 lateShrink(v, ear, blossom, rep, tree); 269 } 270 } 271 } 272 273 274 275 ///\brief Returns the size of the actual matching stored. 276 /// 171 277 ///Returns the size of the actual matching stored. After \ref 172 278 ///run() it returns the size of a maximum matching in the graph. … … 182 288 183 289 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 290 ///\brief Returns the mate of a node in the actual matching. 291 /// 195 292 ///Returns the mate of a \c node in the actual matching. 196 293 ///Returns INVALID if the \c node is not covered by the actual matching. 197 Node mate( Node& node) const {294 Node mate(const Node& node) const { 198 295 return _mate[node]; 199 296 } 200 297 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 } 298 ///\brief Returns the matching edge incident to the given node. 299 /// 300 ///Returns the matching edge of a \c node in the actual matching. 301 ///Returns INVALID if the \c node is not covered by the actual matching. 302 UEdge matchingEdge(const Node& node) const { 303 if (_mate[node] == INVALID) return INVALID; 304 Node n = node < _mate[node] ? node : _mate[node]; 305 for (IncEdgeIt e(g, n); e != INVALID; ++e) { 306 if (g.oppositeNode(n, e) == _mate[n]) { 307 return e; 308 } 309 } 310 return INVALID; 211 311 } 212 213 ///Writes the stored matching to a \c Node valued \c Node map. 214 312 313 /// \brief Returns the class of the node in the Edmonds-Gallai 314 /// decomposition. 315 /// 316 /// Returns the class of the node in the Edmonds-Gallai 317 /// decomposition. 318 DecompType decomposition(const Node& n) { 319 return position[n] == A; 320 } 321 322 /// \brief Returns true when the node is in the barrier. 323 /// 324 /// Returns true when the node is in the barrier. 325 bool barrier(const Node& n) { 326 return position[n] == A; 327 } 328 329 ///\brief Gives back the matching in a \c Node of mates. 330 /// 215 331 ///Writes the stored matching to a \c Node valued \c Node map. The 216 332 ///resulting map will be \e symmetric, i.e. if \c map[u]==v then \c 217 333 ///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 {334 template <typename MateMap> 335 void mateMap(MateMap& map) const { 220 336 for(NodeIt v(g); v!=INVALID; ++v) { 221 337 map.set(v,_mate[v]); 222 338 } 223 339 } 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 340 341 ///\brief Gives back the matching in an \c UEdge valued \c Node 342 ///map. 343 /// 243 344 ///Writes the stored matching to an \c UEdge valued \c Node 244 345 ///map. \c map[v] will be an \c UEdge incident to \c v. This … … 246 347 ///== v then \c map[u]==map[v] holds, and now this edge is an edge 247 348 ///of the matching. 248 template <typename NMapE>249 void writeNMapEdge (NMapE& map) const {349 template <typename MatchingMap> 350 void matchingMap(MatchingMap& map) const { 250 351 typename Graph::template NodeMap<bool> todo(g,true); 251 352 for(NodeIt v(g); v!=INVALID; ++v) { 252 if ( todo[v] && _mate[v]!=INVALID) {353 if (_mate[v]!=INVALID && v < _mate[v]) { 253 354 Node u=_mate[v]; 254 355 for(IncEdgeIt e(g,v); e!=INVALID; ++e) { … … 266 367 267 368 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 369 ///\brief Gives back the matching in a \c bool valued \c UEdge 370 ///map. 371 /// 289 372 ///Writes the matching stored to a \c bool valued \c Edge 290 373 ///map. This map will have the property that there are no two 291 374 ///incident edges \c e, \c f with \c map[e]==map[f]==true. The 292 375 ///edges \c e with \c map[e]==true form the matching. 293 template<typename EMapB>294 void writeEMapBool (EMapB& map) const {376 template<typename MatchingMap> 377 void matching(MatchingMap& map) const { 295 378 for(UEdgeIt e(g); e!=INVALID; ++e) map.set(e,false); 296 379 … … 312 395 313 396 314 /// Writes the canonical decomposition of the graph after running397 ///\brief Returns the canonical decomposition of the graph after running 315 398 ///the algorithm. 316 399 /// 317 400 ///After calling any run methods of the class, it writes the 318 401 ///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]); 402 ///must be a node map of \ref DecompType 's. 403 template <typename DecompositionMap> 404 void decomposition(DecompositionMap& map) const { 405 for(NodeIt v(g); v!=INVALID; ++v) map.set(v,position[v]); 406 } 407 408 ///\brief Returns a barrier on the nodes. 409 /// 410 ///After calling any run methods of the class, it writes a 411 ///canonical barrier on the nodes. The odd component number of the 412 ///remaining graph minus the barrier size is a lower bound for the 413 ///uncovered nodes in the graph. The \c map must be a node map of 414 ///bools. 415 template <typename BarrierMap> 416 void barrier(BarrierMap& barrier) { 417 for(NodeIt v(g); v!=INVALID; ++v) barrier.set(v,position[v] == A); 323 418 } 324 419 … … 327 422 328 423 void lateShrink(Node v, typename Graph::template NodeMap<Node>& ear, 329 UFE& blossom, UFE& tree); 424 UFE& blossom, NV& rep, EFE& tree) { 425 //We have one tree which we grow, and also shrink but only if it 426 //cannot be postponed. If we augment then we return to the "for" 427 //cycle of runEdmonds(). 428 429 std::queue<Node> Q; //queue of the totally unscanned nodes 430 Q.push(v); 431 std::queue<Node> R; 432 //queue of the nodes which must be scanned for a possible shrink 433 434 while ( !Q.empty() ) { 435 Node x=Q.front(); 436 Q.pop(); 437 for( IncEdgeIt e(g,x); e!= INVALID; ++e ) { 438 Node y=g.runningNode(e); 439 //growOrAugment grows if y is covered by the matching and 440 //augments if not. In this latter case it returns 1. 441 if (position[y]==C && 442 growOrAugment(y, x, ear, blossom, rep, tree, Q)) return; 443 } 444 R.push(x); 445 } 446 447 while ( !R.empty() ) { 448 Node x=R.front(); 449 R.pop(); 450 451 for( IncEdgeIt e(g,x); e!=INVALID ; ++e ) { 452 Node y=g.runningNode(e); 453 454 if ( position[y] == D && blossom.find(x) != blossom.find(y) ) 455 //Recall that we have only one tree. 456 shrink( x, y, ear, blossom, rep, tree, Q); 457 458 while ( !Q.empty() ) { 459 Node z=Q.front(); 460 Q.pop(); 461 for( IncEdgeIt f(g,z); f!= INVALID; ++f ) { 462 Node w=g.runningNode(f); 463 //growOrAugment grows if y is covered by the matching and 464 //augments if not. In this latter case it returns 1. 465 if (position[w]==C && 466 growOrAugment(w, z, ear, blossom, rep, tree, Q)) return; 467 } 468 R.push(z); 469 } 470 } //for e 471 } // while ( !R.empty() ) 472 } 330 473 331 474 void normShrink(Node v, typename Graph::template NodeMap<Node>& ear, 332 UFE& blossom, UFE& tree); 475 UFE& blossom, NV& rep, EFE& tree) { 476 //We have one tree, which we grow and shrink. If we augment then we 477 //return to the "for" cycle of runEdmonds(). 478 479 std::queue<Node> Q; //queue of the unscanned nodes 480 Q.push(v); 481 while ( !Q.empty() ) { 482 483 Node x=Q.front(); 484 Q.pop(); 485 486 for( IncEdgeIt e(g,x); e!=INVALID; ++e ) { 487 Node y=g.runningNode(e); 488 489 switch ( position[y] ) { 490 case D: //x and y must be in the same tree 491 if ( blossom.find(x) != blossom.find(y)) 492 //x and y are in the same tree 493 shrink(x, y, ear, blossom, rep, tree, Q); 494 break; 495 case C: 496 //growOrAugment grows if y is covered by the matching and 497 //augments if not. In this latter case it returns 1. 498 if (growOrAugment(y, x, ear, blossom, rep, tree, Q)) return; 499 break; 500 default: break; 501 } 502 } 503 } 504 } 333 505 334 506 void shrink(Node x,Node y, typename Graph::template NodeMap<Node>& ear, 335 UFE& blossom, UFE& tree,std::queue<Node>& Q); 507 UFE& blossom, NV& rep, EFE& tree,std::queue<Node>& Q) { 508 //x and y are the two adjacent vertices in two blossoms. 509 510 typename Graph::template NodeMap<bool> path(g,false); 511 512 Node b=rep[blossom.find(x)]; 513 path.set(b,true); 514 b=_mate[b]; 515 while ( b!=INVALID ) { 516 b=rep[blossom.find(ear[b])]; 517 path.set(b,true); 518 b=_mate[b]; 519 } //we go until the root through bases of blossoms and odd vertices 520 521 Node top=y; 522 Node middle=rep[blossom.find(top)]; 523 Node bottom=x; 524 while ( !path[middle] ) 525 shrinkStep(top, middle, bottom, ear, blossom, rep, tree, Q); 526 //Until we arrive to a node on the path, we update blossom, tree 527 //and the positions of the odd nodes. 528 529 Node base=middle; 530 top=x; 531 middle=rep[blossom.find(top)]; 532 bottom=y; 533 Node blossom_base=rep[blossom.find(base)]; 534 while ( middle!=blossom_base ) 535 shrinkStep(top, middle, bottom, ear, blossom, rep, tree, Q); 536 //Until we arrive to a node on the path, we update blossom, tree 537 //and the positions of the odd nodes. 538 539 rep[blossom.find(base)] = base; 540 } 336 541 337 542 void shrinkStep(Node& top, Node& middle, Node& bottom, 338 543 typename Graph::template NodeMap<Node>& ear, 339 UFE& blossom, UFE& tree, std::queue<Node>& Q); 544 UFE& blossom, NV& rep, EFE& tree, std::queue<Node>& Q) { 545 //We traverse a blossom and update everything. 546 547 ear.set(top,bottom); 548 Node t=top; 549 while ( t!=middle ) { 550 Node u=_mate[t]; 551 t=ear[u]; 552 ear.set(t,u); 553 } 554 bottom=_mate[middle]; 555 position.set(bottom,D); 556 Q.push(bottom); 557 top=ear[bottom]; 558 Node oldmiddle=middle; 559 middle=rep[blossom.find(top)]; 560 tree.erase(bottom); 561 tree.erase(oldmiddle); 562 blossom.insert(bottom); 563 blossom.join(bottom, oldmiddle); 564 blossom.join(top, oldmiddle); 565 } 566 567 340 568 341 569 bool growOrAugment(Node& y, Node& x, typename Graph::template 342 NodeMap<Node>& ear, UFE& blossom, UFE& tree, 343 std::queue<Node>& Q); 570 NodeMap<Node>& ear, UFE& blossom, NV& rep, EFE& tree, 571 std::queue<Node>& Q) { 572 //x is in a blossom in the tree, y is outside. If y is covered by 573 //the matching we grow, otherwise we augment. In this case we 574 //return 1. 575 576 if ( _mate[y]!=INVALID ) { //grow 577 ear.set(y,x); 578 Node w=_mate[y]; 579 rep[blossom.insert(w)] = w; 580 position.set(y,A); 581 position.set(w,D); 582 int t = tree.find(rep[blossom.find(x)]); 583 tree.insert(y,t); 584 tree.insert(w,t); 585 Q.push(w); 586 } else { //augment 587 augment(x, ear, blossom, rep, tree); 588 _mate.set(x,y); 589 _mate.set(y,x); 590 return true; 591 } 592 return false; 593 } 344 594 345 595 void augment(Node x, typename Graph::template NodeMap<Node>& ear, 346 UFE& blossom, UFE& tree); 596 UFE& blossom, NV& rep, EFE& tree) { 597 Node v=_mate[x]; 598 while ( v!=INVALID ) { 599 600 Node u=ear[v]; 601 _mate.set(v,u); 602 Node tmp=v; 603 v=_mate[u]; 604 _mate.set(u,tmp); 605 } 606 int y = tree.find(rep[blossom.find(x)]); 607 for (typename EFE::ItemIt tit(tree, y); tit != INVALID; ++tit) { 608 if ( position[tit] == D ) { 609 int b = blossom.find(tit); 610 for (typename UFE::ItemIt bit(blossom, b); bit != INVALID; ++bit) { 611 position.set(bit, C); 612 } 613 blossom.eraseClass(b); 614 } else position.set(tit, C); 615 } 616 tree.eraseClass(y); 617 618 } 347 619 348 620 }; 349 350 351 // **********************************************************************352 // IMPLEMENTATIONS353 // **********************************************************************354 355 356 template <typename Graph>357 void MaxMatching<Graph>::lateShrink(Node v, typename Graph::template358 NodeMap<Node>& ear, UFE& blossom,359 UFE& tree) {360 //We have one tree which we grow, and also shrink but only if it cannot be361 //postponed. If we augment then we return to the "for" cycle of362 //runEdmonds().363 364 std::queue<Node> Q; //queue of the totally unscanned nodes365 Q.push(v);366 std::queue<Node> R;367 //queue of the nodes which must be scanned for a possible shrink368 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 and375 //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 }381 382 while ( !R.empty() ) {383 Node x=R.front();384 R.pop();385 386 for( IncEdgeIt e(g,x); e!=INVALID ; ++e ) {387 Node y=g.runningNode(e);388 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);392 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 and399 //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 e406 } // while ( !R.empty() )407 }408 409 410 template <typename Graph>411 void MaxMatching<Graph>::normShrink(Node v,412 typename Graph::template413 NodeMap<Node>& ear,414 UFE& blossom, UFE& tree) {415 //We have one tree, which we grow and shrink. If we augment then we416 //return to the "for" cycle of runEdmonds().417 418 std::queue<Node> Q; //queue of the unscanned nodes419 Q.push(v);420 while ( !Q.empty() ) {421 422 Node x=Q.front();423 Q.pop();424 425 for( IncEdgeIt e(g,x); e!=INVALID; ++e ) {426 Node y=g.runningNode(e);427 428 switch ( position[y] ) {429 case D: //x and y must be in the same tree430 if ( blossom.find(x) != blossom.find(y) )431 //x and y are in the same tree432 shrink( x, y, ear, blossom, tree, Q);433 break;434 case C:435 //growOrAugment grows if y is covered by the matching and436 //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 }444 445 446 template <typename Graph>447 void MaxMatching<Graph>::shrink(Node x,Node y, typename448 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.451 452 typename Graph::template NodeMap<bool> path(g,false);453 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 vertices462 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, tree469 //and the positions of the odd nodes.470 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, tree479 //and the positions of the odd nodes.480 481 blossom.makeRep(base);482 }483 484 485 486 template <typename Graph>487 void MaxMatching<Graph>::shrinkStep(Node& top, Node& middle, Node& bottom,488 typename Graph::template489 NodeMap<Node>& ear,490 UFE& blossom, UFE& tree,491 std::queue<Node>& Q) {492 //We traverse a blossom and update everything.493 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 }513 514 515 template <typename Graph>516 bool MaxMatching<Graph>::growOrAugment(Node& y, Node& x, typename Graph::template517 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 by520 //the matching we grow, otherwise we augment. In this case we521 //return 1.522 523 if ( _mate[y]!=INVALID ) { //grow524 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 { //augment535 augment(x, ear, blossom, tree);536 _mate.set(x,y);537 _mate.set(y,x);538 return true;539 }540 return false;541 }542 543 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 ) {550 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);567 568 }569 570 621 571 622 } //END OF NAMESPACE LEMON -
lemon/unionfind.h
r2427 r2505 179 179 private: 180 180 181 ItemIntMap& index; 182 181 183 // If the parent stores negative value for an item then that item 182 // is root item and it has -items[it].parent component size. Else184 // is root item and it has ~(items[it].parent) component id. Else 183 185 // the items[it].parent contains the index of the parent. 184 186 // 185 // The \c nextItem and \c prevItem provides the double-linked 186 // cyclic list of one component's items. The \c prevClass and 187 // \c nextClass gives the double linked list of the representant 188 // items. 187 // The \c next and \c prev provides the double-linked 188 // cyclic list of one component's items. 189 189 struct ItemT { 190 190 int parent; 191 191 Item item; 192 192 193 int nextItem, prevItem; 194 int nextClass, prevClass; 193 int next, prev; 195 194 }; 196 195 197 196 std::vector<ItemT> items; 198 ItemIntMap& index; 199 200 int firstClass; 197 int firstFreeItem; 198 199 struct ClassT { 200 int size; 201 int firstItem; 202 int next, prev; 203 }; 204 205 std::vector<ClassT> classes; 206 int firstClass, firstFreeClass; 207 208 int newClass() { 209 if (firstFreeClass == -1) { 210 int cdx = classes.size(); 211 classes.push_back(ClassT()); 212 return cdx; 213 } else { 214 int cdx = firstFreeClass; 215 firstFreeClass = classes[firstFreeClass].next; 216 return cdx; 217 } 218 } 219 220 int newItem() { 221 if (firstFreeItem == -1) { 222 int idx = items.size(); 223 items.push_back(ItemT()); 224 return idx; 225 } else { 226 int idx = firstFreeItem; 227 firstFreeItem = items[firstFreeItem].next; 228 return idx; 229 } 230 } 201 231 202 232 … … 218 248 } 219 249 220 void unlaceClass(int k) { 221 if (items[k].prevClass != -1) { 222 items[items[k].prevClass].nextClass = items[k].nextClass; 250 int classIndex(int idx) const { 251 return ~(items[repIndex(idx)].parent); 252 } 253 254 void singletonItem(int idx) { 255 items[idx].next = idx; 256 items[idx].prev = idx; 257 } 258 259 void laceItem(int idx, int rdx) { 260 items[idx].prev = rdx; 261 items[idx].next = items[rdx].next; 262 items[items[rdx].next].prev = idx; 263 items[rdx].next = idx; 264 } 265 266 void unlaceItem(int idx) { 267 items[items[idx].prev].next = items[idx].next; 268 items[items[idx].next].prev = items[idx].prev; 269 270 items[idx].next = firstFreeItem; 271 firstFreeItem = idx; 272 } 273 274 void spliceItems(int ak, int bk) { 275 items[items[ak].prev].next = bk; 276 items[items[bk].prev].next = ak; 277 int tmp = items[ak].prev; 278 items[ak].prev = items[bk].prev; 279 items[bk].prev = tmp; 280 281 } 282 283 void laceClass(int cls) { 284 if (firstClass != -1) { 285 classes[firstClass].prev = cls; 286 } 287 classes[cls].next = firstClass; 288 classes[cls].prev = -1; 289 firstClass = cls; 290 } 291 292 void unlaceClass(int cls) { 293 if (classes[cls].prev != -1) { 294 classes[classes[cls].prev].next = classes[cls].next; 223 295 } else { 224 firstClass = items[k].nextClass; 225 } 226 if (items[k].nextClass != -1) { 227 items[items[k].nextClass].prevClass = items[k].prevClass; 228 } 296 firstClass = classes[cls].next; 297 } 298 if (classes[cls].next != -1) { 299 classes[classes[cls].next].prev = classes[cls].prev; 300 } 301 302 classes[cls].next = firstFreeClass; 303 firstFreeClass = cls; 229 304 } 230 305 231 void spliceItems(int ak, int bk) {232 items[items[ak].prevItem].nextItem = bk;233 items[items[bk].prevItem].nextItem = ak;234 int tmp = items[ak].prevItem;235 items[ak].prevItem = items[bk].prevItem;236 items[bk].prevItem = tmp;237 238 }239 240 306 public: 241 307 242 308 UnionFindEnum(ItemIntMap& _index) 243 : items(), index(_index), firstClass(-1) {} 309 : index(_index), items(), firstFreeItem(-1), 310 firstClass(-1), firstFreeClass(-1) {} 244 311 245 312 /// \brief Inserts the given element into a new component. … … 248 315 /// given element. 249 316 /// 250 void insert(const Item& item) { 251 ItemT t; 252 253 int idx = items.size(); 317 int insert(const Item& item) { 318 int idx = newItem(); 319 254 320 index.set(item, idx); 255 321 256 t.nextItem = idx; 257 t.prevItem = idx; 258 t.item = item; 259 t.parent = -1; 260 261 t.nextClass = firstClass; 262 if (firstClass != -1) { 263 items[firstClass].prevClass = idx; 264 } 265 t.prevClass = -1; 266 firstClass = idx; 267 268 items.push_back(t); 322 singletonItem(idx); 323 items[idx].item = item; 324 325 int cdx = newClass(); 326 327 items[idx].parent = ~cdx; 328 329 laceClass(cdx); 330 classes[cdx].size = 1; 331 classes[cdx].firstItem = idx; 332 333 firstClass = cdx; 334 335 return cdx; 269 336 } 270 337 … … 273 340 /// This methods inserts the element \e a into the component of the 274 341 /// element \e comp. 275 void insert(const Item& item, const Item& comp) { 276 int k = repIndex(index[comp]); 277 ItemT t; 278 279 int idx = items.size(); 342 void insert(const Item& item, int cls) { 343 int rdx = classes[cls].firstItem; 344 int idx = newItem(); 345 280 346 index.set(item, idx); 281 347 282 t.prevItem = k; 283 t.nextItem = items[k].nextItem; 284 items[items[k].nextItem].prevItem = idx; 285 items[k].nextItem = idx; 286 287 t.item = item; 288 t.parent = k; 289 290 --items[k].parent; 291 292 items.push_back(t); 348 laceItem(idx, rdx); 349 350 items[idx].item = item; 351 items[idx].parent = rdx; 352 353 ++classes[~(items[rdx].parent)].size; 293 354 } 294 355 … … 299 360 items.clear(); 300 361 firstClass = -1; 301 } 302 303 /// \brief Finds the leader of the component of the given element. 304 /// 305 /// The method returns the leader of the component of the given element. 306 const Item& find(const Item &item) const { 307 return items[repIndex(index[item])].item; 362 firstFreeItem = -1; 363 } 364 365 /// \brief Finds the component of the given element. 366 /// 367 /// The method returns the component id of the given element. 368 int find(const Item &item) const { 369 return ~(items[repIndex(index[item])].parent); 308 370 } 309 371 … … 313 375 /// Joins the component of element \e a and component of 314 376 /// element \e b. If \e a and \e b are in the same component then 315 /// returns false else returns true.316 booljoin(const Item& a, const Item& b) {377 /// returns -1 else returns the remaining class. 378 int join(const Item& a, const Item& b) { 317 379 318 380 int ak = repIndex(index[a]); … … 320 382 321 383 if (ak == bk) { 322 return false; 323 } 324 325 if ( items[ak].parent < items[bk].parent ) { 326 unlaceClass(bk); 327 items[ak].parent += items[bk].parent; 384 return -1; 385 } 386 387 int acx = ~(items[ak].parent); 388 int bcx = ~(items[bk].parent); 389 390 int rcx; 391 392 if (classes[acx].size > classes[bcx].size) { 393 classes[acx].size += classes[bcx].size; 328 394 items[bk].parent = ak; 395 unlaceClass(bcx); 396 rcx = acx; 329 397 } else { 330 unlaceClass(ak); 331 items[bk].parent += items[ak].parent; 398 classes[bcx].size += classes[acx].size; 332 399 items[ak].parent = bk; 400 unlaceClass(acx); 401 rcx = bcx; 333 402 } 334 403 spliceItems(ak, bk); 335 404 336 return true; 337 } 338 339 /// \brief Returns the size of the component of element \e a. 340 /// 341 /// Returns the size of the component of element \e a. 342 int size(const Item &item) const { 343 return - items[repIndex(index[item])].parent; 344 } 345 346 /// \brief Splits up the component of the element. 347 /// 348 /// Splitting the component of the element into sigleton 349 /// components (component of size one). 350 void split(const Item &item) { 351 int k = repIndex(index[item]); 352 int idx = items[k].nextItem; 353 while (idx != k) { 354 int next = items[idx].nextItem; 405 return rcx; 406 } 407 408 /// \brief Returns the size of the class. 409 /// 410 /// Returns the size of the class. 411 int size(int cls) const { 412 return classes[cls].size; 413 } 414 415 /// \brief Splits up the component. 416 /// 417 /// Splitting the component into singleton components (component 418 /// of size one). 419 void split(int cls) { 420 int fdx = classes[cls].firstItem; 421 int idx = items[fdx].next; 422 while (idx != fdx) { 423 int next = items[idx].next; 424 425 singletonItem(idx); 426 427 int cdx = newClass(); 428 items[idx].parent = ~cdx; 429 430 laceClass(cdx); 431 classes[cdx].size = 1; 432 classes[cdx].firstItem = idx; 355 433 356 items[idx].parent = -1;357 items[idx].prevItem = idx;358 items[idx].nextItem = idx;359 360 items[idx].nextClass = firstClass;361 items[firstClass].prevClass = idx;362 firstClass = idx;363 364 434 idx = next; 365 435 } 366 436 367 items[idx].parent = -1; 368 items[idx].prevItem = idx; 369 items[idx].nextItem = idx; 370 371 items[firstClass].prevClass = -1; 372 } 373 374 /// \brief Sets the given element to the leader element of its component. 375 /// 376 /// Sets the given element to the leader element of its component. 377 void makeRep(const Item &item) { 378 int nk = index[item]; 379 int k = repIndex(nk); 380 if (nk == k) return; 381 382 if (items[k].prevClass != -1) { 383 items[items[k].prevClass].nextClass = nk; 384 } else { 385 firstClass = nk; 386 } 387 if (items[k].nextClass != -1) { 388 items[items[k].nextClass].prevClass = nk; 389 } 390 391 int idx = items[k].nextItem; 392 while (idx != k) { 393 items[idx].parent = nk; 394 idx = items[idx].nextItem; 395 } 396 397 items[nk].parent = items[k].parent; 398 items[k].parent = nk; 437 items[idx].prev = idx; 438 items[idx].next = idx; 439 440 classes[~(items[idx].parent)].size = 1; 441 399 442 } 400 443 … … 406 449 /// \warning It is an error to remove an element which is not in 407 450 /// the structure. 408 void erase(const Item &item) { 451 /// \warning This running time of this operation is proportional to the 452 /// number of the items in this class. 453 void erase(const Item& item) { 409 454 int idx = index[item]; 410 if (rep(idx)) { 411 int k = idx; 412 if (items[k].parent == -1) { 413 unlaceClass(idx); 414 return; 415 } else { 416 int nk = items[k].nextItem; 417 if (items[k].prevClass != -1) { 418 items[items[k].prevClass].nextClass = nk; 419 } else { 420 firstClass = nk; 421 } 422 if (items[k].nextClass != -1) { 423 items[items[k].nextClass].prevClass = nk; 424 } 425 426 int l = items[k].nextItem; 427 while (l != k) { 428 items[l].parent = nk; 429 l = items[l].nextItem; 430 } 455 int fdx = items[idx].next; 456 457 int cdx = classIndex(idx); 458 if (idx == fdx) { 459 unlaceClass(cdx); 460 items[idx].next = firstFreeItem; 461 firstFreeItem = idx; 462 return; 463 } else { 464 classes[cdx].firstItem = fdx; 465 --classes[cdx].size; 466 items[fdx].parent = ~cdx; 467 468 unlaceItem(idx); 469 idx = items[fdx].next; 470 while (idx != fdx) { 471 items[idx].parent = fdx; 472 idx = items[idx].next; 473 } 431 474 432 items[nk].parent = items[k].parent + 1; 433 } 434 } else { 435 436 int k = repIndex(idx); 437 idx = items[k].nextItem; 438 while (idx != k) { 439 items[idx].parent = k; 440 idx = items[idx].nextItem; 441 } 442 443 ++items[k].parent; 444 } 445 446 idx = index[item]; 447 items[items[idx].prevItem].nextItem = items[idx].nextItem; 448 items[items[idx].nextItem].prevItem = items[idx].prevItem; 475 } 449 476 450 477 } 451 452 /// \brief Moves the given element to another component.453 ///454 /// This method moves the element \e a from its component455 /// to the component of \e comp.456 /// If \e a and \e comp are in the same component then457 /// it returns false otherwise it returns true.458 bool move(const Item &item, const Item &comp) {459 if (repIndex(index[item]) == repIndex(index[comp])) return false;460 erase(item);461 insert(item, comp);462 return true;463 }464 465 478 466 479 /// \brief Removes the component of the given element from the structure. … … 470 483 /// \warning It is an error to give an element which is not in the 471 484 /// structure. 472 void eraseClass(const Item &item) { 473 unlaceClass(repIndex(index[item])); 485 void eraseClass(int cls) { 486 int fdx = classes[cls].firstItem; 487 unlaceClass(cls); 488 items[items[fdx].prev].next = firstFreeItem; 489 firstFreeItem = fdx; 474 490 } 475 491 … … 477 493 /// 478 494 /// ClassIt is a lemon style iterator for the components. It iterates 479 /// on the representant items of the classes.495 /// on the ids of the classes. 480 496 class ClassIt { 481 497 public: … … 484 500 /// Constructor of the iterator 485 501 ClassIt(const UnionFindEnum& ufe) : unionFind(&ufe) { 486 idx = unionFind->firstClass;502 cdx = unionFind->firstClass; 487 503 } 488 504 … … 490 506 /// 491 507 /// Constructor to get invalid iterator 492 ClassIt(Invalid) : unionFind(0), idx(-1) {}508 ClassIt(Invalid) : unionFind(0), cdx(-1) {} 493 509 494 510 /// \brief Increment operator … … 496 512 /// It steps to the next representant item. 497 513 ClassIt& operator++() { 498 idx = unionFind->items[idx].nextClass;514 cdx = unionFind->classes[cdx].next; 499 515 return *this; 500 516 } … … 503 519 /// 504 520 /// It converts the iterator to the current representant item. 505 operator const Item&() const {506 return unionFind->items[idx].item;521 operator int() const { 522 return cdx; 507 523 } 508 524 … … 511 527 /// Equality operator 512 528 bool operator==(const ClassIt& i) { 513 return i. idx == idx;529 return i.cdx == cdx; 514 530 } 515 531 … … 518 534 /// Inequality operator 519 535 bool operator!=(const ClassIt& i) { 520 return i. idx != idx;536 return i.cdx != cdx; 521 537 } 522 538 523 539 private: 524 540 const UnionFindEnum* unionFind; 525 int idx;541 int cdx; 526 542 }; 527 543 … … 546 562 /// Constructor of the iterator. The iterator iterates 547 563 /// on the class of the \c item. 548 ItemIt(const UnionFindEnum& ufe, const Item& item) : unionFind(&ufe) {549 idx = unionFind->repIndex(unionFind->index[item]);564 ItemIt(const UnionFindEnum& ufe, int cls) : unionFind(&ufe) { 565 fdx = idx = unionFind->classes[cls].firstItem; 550 566 } 551 567 … … 559 575 /// It steps to the next item in the class. 560 576 ItemIt& operator++() { 561 idx = unionFind->items[idx].next Item;562 if ( unionFind->rep(idx)) idx = -1;577 idx = unionFind->items[idx].next; 578 if (idx == fdx) idx = -1; 563 579 return *this; 564 580 } … … 587 603 private: 588 604 const UnionFindEnum* unionFind; 589 int idx ;605 int idx, fdx; 590 606 }; 591 607 592 608 }; 593 609 610 /// \ingroup auxdat 611 /// 612 /// \brief A \e Extend-Find data structure implementation which 613 /// is able to enumerate the components. 614 /// 615 /// The class implements an \e Extend-Find data structure which is 616 /// able to enumerate the components and the items in a 617 /// component. The data structure is a simplification of the 618 /// Union-Find structure, and it does not allow to merge two components. 619 /// 620 /// \pre You need to add all the elements by the \ref insert() 621 /// method. 622 template <typename _ItemIntMap> 623 class ExtendFindEnum { 624 public: 625 626 typedef _ItemIntMap ItemIntMap; 627 typedef typename ItemIntMap::Key Item; 628 629 private: 630 631 ItemIntMap& index; 632 633 struct ItemT { 634 int cls; 635 Item item; 636 int next, prev; 637 }; 638 639 std::vector<ItemT> items; 640 int firstFreeItem; 641 642 struct ClassT { 643 int firstItem; 644 int next, prev; 645 }; 646 647 std::vector<ClassT> classes; 648 649 int firstClass, firstFreeClass; 650 651 int newClass() { 652 if (firstFreeClass != -1) { 653 int cdx = firstFreeClass; 654 firstFreeClass = classes[cdx].next; 655 return cdx; 656 } else { 657 classes.push_back(ClassT()); 658 return classes.size() - 1; 659 } 660 } 661 662 int newItem() { 663 if (firstFreeItem != -1) { 664 int idx = firstFreeItem; 665 firstFreeItem = items[idx].next; 666 return idx; 667 } else { 668 items.push_back(ItemT()); 669 return items.size() - 1; 670 } 671 } 672 673 public: 674 675 /// \brief Constructor 676 ExtendFindEnum(ItemIntMap& _index) 677 : index(_index), items(), firstFreeItem(-1), 678 classes(), firstClass(-1), firstFreeClass(-1) {} 679 680 /// \brief Inserts the given element into a new component. 681 /// 682 /// This method creates a new component consisting only of the 683 /// given element. 684 int insert(const Item& item) { 685 int cdx = newClass(); 686 classes[cdx].prev = -1; 687 classes[cdx].next = firstClass; 688 firstClass = cdx; 689 690 int idx = newItem(); 691 items[idx].item = item; 692 items[idx].cls = cdx; 693 items[idx].prev = idx; 694 items[idx].next = idx; 695 696 classes[cdx].firstItem = idx; 697 698 index.set(item, idx); 699 700 return cdx; 701 } 702 703 /// \brief Inserts the given element into the given component. 704 /// 705 /// This methods inserts the element \e item a into the \e cls class. 706 void insert(const Item& item, int cls) { 707 int idx = newItem(); 708 int rdx = classes[cls].firstItem; 709 items[idx].item = item; 710 items[idx].cls = cls; 711 712 items[idx].prev = rdx; 713 items[idx].next = items[rdx].next; 714 items[items[rdx].next].prev = idx; 715 items[rdx].next = idx; 716 717 index.set(item, idx); 718 } 719 720 /// \brief Clears the union-find data structure 721 /// 722 /// Erase each item from the data structure. 723 void clear() { 724 items.clear(); 725 classes.clear; 726 firstClass = firstFreeClass = firstFreeItem = -1; 727 } 728 729 /// \brief Gives back the class of the \e item. 730 /// 731 /// Gives back the class of the \e item. 732 int find(const Item &item) const { 733 return items[index[item]].cls; 734 } 735 736 /// \brief Removes the given element from the structure. 737 /// 738 /// Removes the element from its component and if the component becomes 739 /// empty then removes that component from the component list. 740 /// 741 /// \warning It is an error to remove an element which is not in 742 /// the structure. 743 void erase(const Item &item) { 744 int idx = index[item]; 745 int cdx = items[idx].cls; 746 747 if (idx == items[idx].next) { 748 if (classes[cdx].prev != -1) { 749 classes[classes[cdx].prev].next = classes[cdx].next; 750 } else { 751 firstClass = classes[cdx].next; 752 } 753 if (classes[cdx].next != -1) { 754 classes[classes[cdx].next].prev = classes[cdx].prev; 755 } 756 classes[cdx].next = firstFreeClass; 757 firstFreeClass = cdx; 758 } else { 759 classes[cdx].firstItem = items[idx].next; 760 items[items[idx].next].prev = items[idx].prev; 761 items[items[idx].prev].next = items[idx].next; 762 } 763 items[idx].next = firstFreeItem; 764 firstFreeItem = idx; 765 766 } 767 768 769 /// \brief Removes the component of the given element from the structure. 770 /// 771 /// Removes the component of the given element from the structure. 772 /// 773 /// \warning It is an error to give an element which is not in the 774 /// structure. 775 void eraseClass(int cdx) { 776 int idx = classes[cdx].firstItem; 777 items[items[idx].prev].next = firstFreeItem; 778 firstFreeItem = idx; 779 780 if (classes[cdx].prev != -1) { 781 classes[classes[cdx].prev].next = classes[cdx].next; 782 } else { 783 firstClass = classes[cdx].next; 784 } 785 if (classes[cdx].next != -1) { 786 classes[classes[cdx].next].prev = classes[cdx].prev; 787 } 788 classes[cdx].next = firstFreeClass; 789 firstFreeClass = cdx; 790 } 791 792 /// \brief Lemon style iterator for the classes. 793 /// 794 /// ClassIt is a lemon style iterator for the components. It iterates 795 /// on the ids of classes. 796 class ClassIt { 797 public: 798 /// \brief Constructor of the iterator 799 /// 800 /// Constructor of the iterator 801 ClassIt(const ExtendFindEnum& ufe) : extendFind(&ufe) { 802 cdx = extendFind->firstClass; 803 } 804 805 /// \brief Constructor to get invalid iterator 806 /// 807 /// Constructor to get invalid iterator 808 ClassIt(Invalid) : extendFind(0), cdx(-1) {} 809 810 /// \brief Increment operator 811 /// 812 /// It steps to the next representant item. 813 ClassIt& operator++() { 814 cdx = extendFind->classes[cdx].next; 815 return *this; 816 } 817 818 /// \brief Conversion operator 819 /// 820 /// It converts the iterator to the current class id. 821 operator int() const { 822 return cdx; 823 } 824 825 /// \brief Equality operator 826 /// 827 /// Equality operator 828 bool operator==(const ClassIt& i) { 829 return i.cdx == cdx; 830 } 831 832 /// \brief Inequality operator 833 /// 834 /// Inequality operator 835 bool operator!=(const ClassIt& i) { 836 return i.cdx != cdx; 837 } 838 839 private: 840 const ExtendFindEnum* extendFind; 841 int cdx; 842 }; 843 844 /// \brief Lemon style iterator for the items of a component. 845 /// 846 /// ClassIt is a lemon style iterator for the components. It iterates 847 /// on the items of a class. By example if you want to iterate on 848 /// each items of each classes then you may write the next code. 849 ///\code 850 /// for (ClassIt cit(ufe); cit != INVALID; ++cit) { 851 /// std::cout << "Class: "; 852 /// for (ItemIt iit(ufe, cit); iit != INVALID; ++iit) { 853 /// std::cout << toString(iit) << ' ' << std::endl; 854 /// } 855 /// std::cout << std::endl; 856 /// } 857 ///\endcode 858 class ItemIt { 859 public: 860 /// \brief Constructor of the iterator 861 /// 862 /// Constructor of the iterator. The iterator iterates 863 /// on the class of the \c item. 864 ItemIt(const ExtendFindEnum& ufe, int cls) : extendFind(&ufe) { 865 fdx = idx = extendFind->classes[cls].firstItem; 866 } 867 868 /// \brief Constructor to get invalid iterator 869 /// 870 /// Constructor to get invalid iterator 871 ItemIt(Invalid) : extendFind(0), idx(-1) {} 872 873 /// \brief Increment operator 874 /// 875 /// It steps to the next item in the class. 876 ItemIt& operator++() { 877 idx = extendFind->items[idx].next; 878 if (fdx == idx) idx = -1; 879 return *this; 880 } 881 882 /// \brief Conversion operator 883 /// 884 /// It converts the iterator to the current item. 885 operator const Item&() const { 886 return extendFind->items[idx].item; 887 } 888 889 /// \brief Equality operator 890 /// 891 /// Equality operator 892 bool operator==(const ItemIt& i) { 893 return i.idx == idx; 894 } 895 896 /// \brief Inequality operator 897 /// 898 /// Inequality operator 899 bool operator!=(const ItemIt& i) { 900 return i.idx != idx; 901 } 902 903 private: 904 const ExtendFindEnum* extendFind; 905 int idx, fdx; 906 }; 907 908 }; 594 909 595 910 //! @} -
test/max_matching_test.cc
r2391 r2505 67 67 68 68 MaxMatching<Graph> max_matching(g); 69 max_matching.runEdmonds(0); 69 max_matching.init(); 70 max_matching.startDense(); 70 71 71 72 int s=0; 72 73 Graph::NodeMap<Node> mate(g,INVALID); 73 max_matching. writeNMapNode(mate);74 max_matching.mateMap(mate); 74 75 for(NodeIt v(g); v!=INVALID; ++v) { 75 76 if ( mate[v]!=INVALID ) ++s; … … 83 84 check ( size == max_matching.size(), "mate() returns a different size matching than max_matching.size()" ); 84 85 85 Graph::NodeMap<MaxMatching<Graph>:: pos_enum> pos0(g);86 max_matching. writePos(pos0);86 Graph::NodeMap<MaxMatching<Graph>::DecompType> pos0(g); 87 max_matching.decomposition(pos0); 87 88 88 max_matching. resetMatching();89 max_matching. runEdmonds(1);89 max_matching.init(); 90 max_matching.startSparse(); 90 91 s=0; 91 max_matching. writeNMapNode(mate);92 max_matching.mateMap(mate); 92 93 for(NodeIt v(g); v!=INVALID; ++v) { 93 94 if ( mate[v]!=INVALID ) ++s; … … 95 96 check ( int(s/2) == size, "The size does not equal!" ); 96 97 97 Graph::NodeMap<MaxMatching<Graph>:: pos_enum> pos1(g);98 max_matching. writePos(pos1);98 Graph::NodeMap<MaxMatching<Graph>::DecompType> pos1(g); 99 max_matching.decomposition(pos1); 99 100 100 101 max_matching.run(); 101 102 s=0; 102 max_matching. writeNMapNode(mate);103 max_matching.mateMap(mate); 103 104 for(NodeIt v(g); v!=INVALID; ++v) { 104 105 if ( mate[v]!=INVALID ) ++s; … … 106 107 check ( int(s/2) == size, "The size does not equal!" ); 107 108 108 Graph::NodeMap<MaxMatching<Graph>:: pos_enum> pos2(g);109 max_matching. writePos(pos2);109 Graph::NodeMap<MaxMatching<Graph>::DecompType> pos2(g); 110 max_matching.decomposition(pos2); 110 111 111 max_matching.resetMatching();112 112 max_matching.run(); 113 113 s=0; 114 max_matching. writeNMapNode(mate);114 max_matching.mateMap(mate); 115 115 for(NodeIt v(g); v!=INVALID; ++v) { 116 116 if ( mate[v]!=INVALID ) ++s; … … 118 118 check ( int(s/2) == size, "The size does not equal!" ); 119 119 120 Graph::NodeMap<MaxMatching<Graph>:: pos_enum> pos(g);121 max_matching. writePos(pos);120 Graph::NodeMap<MaxMatching<Graph>::DecompType> pos(g); 121 max_matching.decomposition(pos); 122 122 123 123 bool ismatching=true; … … 140 140 bool noedge=true; 141 141 for(UEdgeIt e(g); e!=INVALID; ++e) { 142 if ( (pos[g.target(e)]==max_matching.C && pos[g.source(e)]==max_matching.D) || 143 (pos[g.target(e)]==max_matching.D && pos[g.source(e)]==max_matching.C) ) 142 if ( (pos[g.target(e)]==max_matching.C && 143 pos[g.source(e)]==max_matching.D) || 144 (pos[g.target(e)]==max_matching.D && 145 pos[g.source(e)]==max_matching.C) ) 144 146 noedge=false; 145 147 } -
test/unionfind_test.cc
r2391 r2505 50 50 U.insert(2); 51 51 52 check(U.join(1,2) ,"Test failed.");52 check(U.join(1,2) != -1,"Test failed."); 53 53 54 54 U.insert(3); … … 58 58 U.insert(7); 59 59 60 check(U.join(1,4),"Test failed.");61 check(!U.join(2,4),"Test failed.");62 check(U.join(3,5),"Test failed.");63 60 64 U.insert(8,5); 61 check(U.join(1,4) != -1,"Test failed."); 62 check(U.join(2,4) == -1,"Test failed."); 63 check(U.join(3,5) != -1,"Test failed."); 65 64 66 check(U.size(4) == 3,"Test failed."); 67 check(U.size(5) == 3,"Test failed."); 68 check(U.size(6) == 1,"Test failed."); 69 check(U.size(2) == 3,"Test failed."); 65 66 U.insert(8,U.find(5)); 67 68 69 check(U.size(U.find(4)) == 3,"Test failed."); 70 check(U.size(U.find(5)) == 3,"Test failed."); 71 check(U.size(U.find(6)) == 1,"Test failed."); 72 check(U.size(U.find(2)) == 3,"Test failed."); 73 70 74 71 75 U.insert(9); 72 U.insert(10, 9);76 U.insert(10,U.find(9)); 73 77 74 check(U.join(8,10),"Test failed.");75 78 76 check(U.move(9,4),"Test failed."); 77 check(!U.move(9,2),"Test failed."); 79 check(U.join(8,10) != -1,"Test failed."); 78 80 79 check(U.size(4) == 4,"Test failed.");80 check(U.size(9) == 4,"Test failed.");81 81 82 check(U.move(5,6),"Test failed."); 82 check(U.size(U.find(4)) == 3,"Test failed."); 83 check(U.size(U.find(9)) == 5,"Test failed."); 83 84 84 check(U.size(5) == 2,"Test failed."); 85 check(U.size(8) == 3,"Test failed."); 86 87 check(U.move(7,10),"Test failed."); 88 check(U.size(7) == 4,"Test failed."); 85 check(U.size(U.find(8)) == 5,"Test failed."); 89 86 90 87 U.erase(9); 91 88 U.erase(1); 92 89 93 check(U.size( 4) == 2,"Test failed.");94 check(U.size( 2) == 2,"Test failed.");90 check(U.size(U.find(10)) == 4,"Test failed."); 91 check(U.size(U.find(2)) == 2,"Test failed."); 95 92 96 93 U.erase(6); 97 U.split(8); 98 99 check(U.size(4) == 2,"Test failed."); 100 check(U.size(3) == 1,"Test failed."); 101 check(U.size(2) == 2,"Test failed."); 102 103 check(U.join(3,4),"Test failed."); 104 check(!U.join(2,4),"Test failed."); 105 106 check(U.size(4) == 3,"Test failed."); 107 check(U.size(3) == 3,"Test failed."); 108 check(U.size(2) == 3,"Test failed."); 109 110 U.makeRep(4); 111 U.makeRep(3); 112 U.makeRep(2); 113 114 check(U.size(4) == 3,"Test failed."); 115 check(U.size(3) == 3,"Test failed."); 116 check(U.size(2) == 3,"Test failed."); 94 U.split(U.find(8)); 117 95 118 96 119 U.eraseClass(4); 120 U.eraseClass(7); 97 check(U.size(U.find(4)) == 2,"Test failed."); 98 check(U.size(U.find(3)) == 1,"Test failed."); 99 check(U.size(U.find(2)) == 2,"Test failed."); 121 100 101 102 check(U.join(3,4) != -1,"Test failed."); 103 check(U.join(2,4) == -1,"Test failed."); 104 105 106 check(U.size(U.find(4)) == 3,"Test failed."); 107 check(U.size(U.find(3)) == 3,"Test failed."); 108 check(U.size(U.find(2)) == 3,"Test failed."); 109 110 U.eraseClass(U.find(4)); 111 U.eraseClass(U.find(7)); 112 113 return 0; 122 114 }
Note: See TracChangeset
for help on using the changeset viewer.