- Gamma distributon random variable.
- Test file for random.h
3 * This file is a part of LEMON, a generic C++ optimization library
5 * Copyright (C) 2003-2007
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
7 * (Egervary Research Group on Combinatorial Optimization, EGRES).
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.
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
19 #ifndef LEMON_PR_BIPARTITE_MATCHING
20 #define LEMON_PR_BIPARTITE_MATCHING
22 #include <lemon/graph_utils.h>
23 #include <lemon/iterable_maps.h>
26 #include <lemon/elevator.h>
30 ///\brief Push-prelabel maximum matching algorithms in bipartite graphs.
34 ///Max cardinality matching algorithm based on push-relabel principle
37 ///Bipartite Max Cardinality Matching algorithm. This class uses the
38 ///push-relabel principle which in several cases has better runtime
39 ///performance than the augmenting path solutions.
41 ///\author Alpar Juttner
43 class PrBipartiteMatching {
44 typedef typename Graph::Node Node;
45 typedef typename Graph::ANodeIt ANodeIt;
46 typedef typename Graph::BNodeIt BNodeIt;
47 typedef typename Graph::UEdge UEdge;
48 typedef typename Graph::UEdgeIt UEdgeIt;
49 typedef typename Graph::IncEdgeIt IncEdgeIt;
56 typename Graph::template ANodeMap<typename Graph::UEdge> _matching;
57 Elevator<Graph,typename Graph::BNode> _levels;
58 typename Graph::template BNodeMap<int> _cov;
66 PrBipartiteMatching(const Graph &g) :
68 _node_num(countBNodes(g)),
75 /// \name Execution control
76 /// The simplest way to execute the algorithm is to use one of the
77 /// member functions called \c run(). \n
78 /// If you need more control on the execution, first
79 /// you must call \ref init() and then one variant of the start()
84 ///Initialize the data structures
86 ///This function constructs a prematching first, which is a
87 ///regular matching on the A-side of the graph, but on the B-side
88 ///each node could cover more matching edges. After that, the
89 ///B-nodes which multiple matched, will be pushed into the lowest
90 ///level of the Elevator. The remaning B-nodes will be pushed to
91 ///the consequent levels respect to a Bfs on following graph: the
92 ///nodes are the B-nodes of the original bipartite graph and two
93 ///nodes are adjacent if a node can pass over a matching edge to
94 ///an other node. The source of the Bfs are the lowest level
95 ///nodes. Last, the reached B-nodes without covered matching edge
99 _empty_level=_node_num;
100 for(ANodeIt n(_g);n!=INVALID;++n)
101 if((_matching[n]=IncEdgeIt(_g,n))!=INVALID)
102 ++_cov[_g.bNode(_matching[n])];
106 for(BNodeIt n(_g);n!=INVALID;++n)
108 _levels.initAddItem(n);
115 int nlev=_levels[n]+1;
116 for(IncEdgeIt e(_g,n);e!=INVALID;++e) {
117 Node m=_g.runningNode(e);
118 if(e==_matching[m]) {
119 for(IncEdgeIt f(_g,m);f!=INVALID;++f) {
120 Node r=_g.runningNode(f);
121 if(_levels[r]>nlev) {
122 for(;nlev>hlev;hlev++)
123 _levels.initNewLevel();
124 _levels.initAddItem(r);
131 _levels.initFinish();
132 for(BNodeIt n(_g);n!=INVALID;++n)
133 if(_cov[n]<1&&_levels[n]<_node_num)
137 ///Start the main phase of the algorithm
139 ///This algorithm calculates the maximum matching with the
140 ///push-relabel principle. This function should be called just
141 ///after the init() function which already set the initial
142 ///prematching, the level function on the B-nodes and the active,
143 ///ie. unmatched B-nodes.
145 ///The algorithm always takes highest active B-node, and it try to
146 ///find a B-node which is eligible to pass over one of it's
147 ///matching edge. This condition holds when the B-node is one
148 ///level lower, and the opposite node of it's matching edge is
149 ///adjacent to the highest active node. In this case the current
150 ///node steals the matching edge and becomes inactive. If there is
151 ///not eligible node then the highest active node should be lift
152 ///to the next proper level.
154 ///The nodes should not lift higher than the number of the
155 ///B-nodes, if a node reach this level it remains unmatched. If
156 ///during the execution one level becomes empty the nodes above it
157 ///can be deactivated and lift to the highest level.
161 Node last_activated=INVALID;
162 while((act=_levels.highestActive())!=INVALID) {
163 last_activated=INVALID;
164 int actlevel=_levels[act];
167 int nlevel=_node_num;
170 for(IncEdgeIt tbedge(_g,act);
171 tbedge!=INVALID && nlevel>=actlevel;
173 if((nnlevel=_levels[_g.bNode(_matching[_g.runningNode(tbedge)])])<
180 if(nlevel<_node_num) {
182 _levels.liftHighestActiveTo(nlevel+1);
183 bact=_g.bNode(_matching[_g.aNode(bedge)]);
185 _levels.activate(bact);
188 _matching[_g.aNode(bedge)]=bedge;
190 _levels.deactivate(act);
193 if(_node_num>actlevel)
194 _levels.liftHighestActiveTo(_node_num);
195 _levels.deactivate(act);
198 if(_levels.onLevel(actlevel)==0)
199 _levels.liftToTop(actlevel);
202 for(ANodeIt n(_g);n!=INVALID;++n) {
203 if (_matching[n]==INVALID)continue;
204 if (_cov[_g.bNode(_matching[n])]>1) {
205 _cov[_g.bNode(_matching[n])]--;
206 _matching[n]=INVALID;
213 ///Start the algorithm to find a perfect matching
215 ///This function is close to identical to the simple start()
216 ///member function but it calculates just perfect matching.
217 ///However, the perfect property is only checked on the B-side of
220 ///The main difference between the two function is the handling of
221 ///the empty levels. The simple start() function let the nodes
222 ///above the empty levels unmatched while this variant if it find
223 ///an empty level immediately terminates and gives back false
225 bool startPerfect() {
228 Node last_activated=INVALID;
229 while((act=_levels.highestActive())!=INVALID) {
230 last_activated=INVALID;
231 int actlevel=_levels[act];
234 int nlevel=_node_num;
237 for(IncEdgeIt tbedge(_g,act);
238 tbedge!=INVALID && nlevel>=actlevel;
240 if((nnlevel=_levels[_g.bNode(_matching[_g.runningNode(tbedge)])])<
247 if(nlevel<_node_num) {
249 _levels.liftHighestActiveTo(nlevel+1);
250 bact=_g.bNode(_matching[_g.aNode(bedge)]);
252 _levels.activate(bact);
255 _matching[_g.aNode(bedge)]=bedge;
257 _levels.deactivate(act);
260 if(_node_num>actlevel)
261 _levels.liftHighestActiveTo(_node_num);
262 _levels.deactivate(act);
265 if(_levels.onLevel(actlevel)==0)
266 _empty_level=actlevel;
269 _matching_size = _node_num;
273 ///Runs the algorithm
275 ///Just a shortcut for the next code:
285 ///Runs the algorithm to find a perfect matching
287 ///Just a shortcut for the next code:
293 ///\note If the two nodesets of the graph have different size then
294 ///this algorithm checks the perfect property on the B-side.
297 return startPerfect();
300 ///Runs the algorithm to find a perfect matching
302 ///Just a shortcut for the next code:
308 ///\note It checks that the size of the two nodesets are equal.
309 bool checkedRunPerfect() {
310 if (countANodes(_g) != _node_num) return false;
312 return startPerfect();
317 /// \name Query Functions
318 /// The result of the %Matching algorithm can be obtained using these
320 /// Before the use of these functions,
321 /// either run() or start() must be called.
324 ///Set true all matching uedge in the map.
326 ///Set true all matching uedge in the map. It does not change the
327 ///value mapped to the other uedges.
328 ///\return The number of the matching edges.
329 template <typename MatchingMap>
330 int quickMatching(MatchingMap& mm) const {
331 for (ANodeIt n(_g);n!=INVALID;++n) {
332 if (_matching[n]!=INVALID) mm.set(_matching[n],true);
334 return _matching_size;
337 ///Set true all matching uedge in the map and the others to false.
339 ///Set true all matching uedge in the map and the others to false.
340 ///\return The number of the matching edges.
341 template<class MatchingMap>
342 int matching(MatchingMap& mm) const {
343 for (UEdgeIt e(_g);e!=INVALID;++e) {
344 mm.set(e,e==_matching[_g.aNode(e)]);
346 return _matching_size;
349 ///Gives back the matching in an ANodeMap.
351 ///Gives back the matching in an ANodeMap. The parameter should
352 ///be a write ANodeMap of UEdge values.
353 ///\return The number of the matching edges.
354 template<class MatchingMap>
355 int aMatching(MatchingMap& mm) const {
356 for (ANodeIt n(_g);n!=INVALID;++n) {
357 mm.set(n,_matching[n]);
359 return _matching_size;
362 ///Gives back the matching in a BNodeMap.
364 ///Gives back the matching in a BNodeMap. The parameter should
365 ///be a write BNodeMap of UEdge values.
366 ///\return The number of the matching edges.
367 template<class MatchingMap>
368 int bMatching(MatchingMap& mm) const {
369 for (BNodeIt n(_g);n!=INVALID;++n) {
372 for (ANodeIt n(_g);n!=INVALID;++n) {
373 if (_matching[n]!=INVALID)
374 mm.set(_g.bNode(_matching[n]),_matching[n]);
376 return _matching_size;
380 ///Returns true if the given uedge is in the matching.
382 ///It returns true if the given uedge is in the matching.
384 bool matchingEdge(const UEdge& e) const {
385 return _matching[_g.aNode(e)]==e;
388 ///Returns the matching edge from the node.
390 ///Returns the matching edge from the node. If there is not such
391 ///edge it gives back \c INVALID.
392 ///\note If the parameter node is a B-node then the running time is
393 ///propotional to the degree of the node.
394 UEdge matchingEdge(const Node& n) const {
398 for (IncEdgeIt e(_g,n);e!=INVALID;++e)
399 if (e==_matching[_g.aNode(e)]) return e;
404 ///Gives back the number of the matching edges.
406 ///Gives back the number of the matching edges.
407 int matchingSize() const {
408 return _matching_size;
411 ///Gives back a barrier on the A-nodes
413 ///The barrier is s subset of the nodes on the same side of the
414 ///graph. If we tried to find a perfect matching and it failed
415 ///then the barrier size will be greater than its neighbours. If
416 ///the maximum matching searched then the barrier size minus its
417 ///neighbours will be exactly the unmatched nodes on the A-side.
418 ///\retval bar A WriteMap on the ANodes with bool value.
419 template<class BarrierMap>
420 void aBarrier(BarrierMap &bar) const
422 for(ANodeIt n(_g);n!=INVALID;++n)
423 bar.set(n,_matching[n]==INVALID ||
424 _levels[_g.bNode(_matching[n])]<_empty_level);
427 ///Gives back a barrier on the B-nodes
429 ///The barrier is s subset of the nodes on the same side of the
430 ///graph. If we tried to find a perfect matching and it failed
431 ///then the barrier size will be greater than its neighbours. If
432 ///the maximum matching searched then the barrier size minus its
433 ///neighbours will be exactly the unmatched nodes on the B-side.
434 ///\retval bar A WriteMap on the BNodes with bool value.
435 template<class BarrierMap>
436 void bBarrier(BarrierMap &bar) const
438 for(BNodeIt n(_g);n!=INVALID;++n) bar.set(n,_levels[n]>=_empty_level);
441 ///Returns a minimum covering of the nodes.
443 ///The minimum covering set problem is the dual solution of the
444 ///maximum bipartite matching. It provides a solution for this
445 ///problem what is proof of the optimality of the matching.
446 ///\param covering NodeMap of bool values, the nodes of the cover
447 ///set will set true while the others false.
448 ///\return The size of the cover set.
449 ///\note This function can be called just after the algorithm have
450 ///already found a matching.
451 template<class CoverMap>
452 int coverSet(CoverMap& covering) const {
454 for(BNodeIt n(_g);n!=INVALID;++n) {
455 if (_levels[n]<_empty_level) { covering.set(n,true); ++ret; }
456 else covering.set(n,false);
458 for(ANodeIt n(_g);n!=INVALID;++n) {
459 if (_matching[n]!=INVALID &&
460 _levels[_g.bNode(_matching[n])]>=_empty_level)
461 { covering.set(n,true); ++ret; }
462 else covering.set(n,false);
473 ///Maximum cardinality of the matchings in a bipartite graph
476 ///This function finds the maximum cardinality of the matchings
477 ///in a bipartite graph \c g.
478 ///\param g An undirected bipartite graph.
479 ///\return The cardinality of the maximum matching.
481 ///\note The the implementation is based
482 ///on the push-relabel principle.
483 template<class Graph>
484 int prBipartiteMatching(const Graph &g)
486 PrBipartiteMatching<Graph> bpm(g);
487 return bpm.matchingSize();
490 ///Maximum cardinality matching in a bipartite graph
493 ///This function finds a maximum cardinality matching
494 ///in a bipartite graph \c g.
495 ///\param g An undirected bipartite graph.
496 ///\retval matching A write ANodeMap of value type \c UEdge.
497 /// The found edges will be returned in this map,
498 /// i.e. for an \c ANode \c n the edge <tt>matching[n]</tt> is the one
499 /// that covers the node \c n.
500 ///\return The cardinality of the maximum matching.
502 ///\note The the implementation is based
503 ///on the push-relabel principle.
504 template<class Graph,class MT>
505 int prBipartiteMatching(const Graph &g,MT &matching)
507 PrBipartiteMatching<Graph> bpm(g);
509 bpm.aMatching(matching);
510 return bpm.matchingSize();
513 ///Maximum cardinality matching in a bipartite graph
516 ///This function finds a maximum cardinality matching
517 ///in a bipartite graph \c g.
518 ///\param g An undirected bipartite graph.
519 ///\retval matching A write ANodeMap of value type \c UEdge.
520 /// The found edges will be returned in this map,
521 /// i.e. for an \c ANode \c n the edge <tt>matching[n]</tt> is the one
522 /// that covers the node \c n.
523 ///\retval barrier A \c bool WriteMap on the BNodes. The map will be set
524 /// exactly once for each BNode. The nodes with \c true value represent
525 /// a barrier \e B, i.e. the cardinality of \e B minus the number of its
526 /// neighbor is equal to the number of the <tt>BNode</tt>s minus the
527 /// cardinality of the maximum matching.
528 ///\return The cardinality of the maximum matching.
530 ///\note The the implementation is based
531 ///on the push-relabel principle.
532 template<class Graph,class MT, class GT>
533 int prBipartiteMatching(const Graph &g,MT &matching,GT &barrier)
535 PrBipartiteMatching<Graph> bpm(g);
537 bpm.aMatching(matching);
538 bpm.bBarrier(barrier);
539 return bpm.matchingSize();
542 ///Perfect matching in a bipartite graph
545 ///This function checks whether the bipartite graph \c g
546 ///has a perfect matching.
547 ///\param g An undirected bipartite graph.
548 ///\return \c true iff \c g has a perfect matching.
550 ///\note The the implementation is based
551 ///on the push-relabel principle.
552 template<class Graph>
553 bool prPerfectBipartiteMatching(const Graph &g)
555 PrBipartiteMatching<Graph> bpm(g);
556 return bpm.runPerfect();
559 ///Perfect matching in a bipartite graph
562 ///This function finds a perfect matching in a bipartite graph \c g.
563 ///\param g An undirected bipartite graph.
564 ///\retval matching A write ANodeMap of value type \c UEdge.
565 /// The found edges will be returned in this map,
566 /// i.e. for an \c ANode \c n the edge <tt>matching[n]</tt> is the one
567 /// that covers the node \c n.
568 /// The values are unchanged if the graph
569 /// has no perfect matching.
570 ///\return \c true iff \c g has a perfect matching.
572 ///\note The the implementation is based
573 ///on the push-relabel principle.
574 template<class Graph,class MT>
575 bool prPerfectBipartiteMatching(const Graph &g,MT &matching)
577 PrBipartiteMatching<Graph> bpm(g);
578 bool ret = bpm.checkedRunPerfect();
579 if (ret) bpm.aMatching(matching);
583 ///Perfect matching in a bipartite graph
586 ///This function finds a perfect matching in a bipartite graph \c g.
587 ///\param g An undirected bipartite graph.
588 ///\retval matching A write ANodeMap of value type \c UEdge.
589 /// The found edges will be returned in this map,
590 /// i.e. for an \c ANode \c n the edge <tt>matching[n]</tt> is the one
591 /// that covers the node \c n.
592 /// The values are unchanged if the graph
593 /// has no perfect matching.
594 ///\retval barrier A \c bool WriteMap on the BNodes. The map will only
595 /// be set if \c g has no perfect matching. In this case it is set
596 /// exactly once for each BNode. The nodes with \c true value represent
597 /// a barrier, i.e. a subset \e B a of BNodes with the property that
598 /// the cardinality of \e B is greater than the number of its neighbors.
599 ///\return \c true iff \c g has a perfect matching.
601 ///\note The the implementation is based
602 ///on the push-relabel principle.
603 template<class Graph,class MT, class GT>
604 bool prPerfectBipartiteMatching(const Graph &g,MT &matching,GT &barrier)
606 PrBipartiteMatching<Graph> bpm(g);
607 bool ret=bpm.checkedRunPerfect();
609 bpm.aMatching(matching);
611 bpm.bBarrier(barrier);