diff -r d9a221218ea4 -r 08652c1763f6 lemon/bipartite_matching.h --- a/lemon/bipartite_matching.h Fri Apr 14 18:05:02 2006 +0000 +++ b/lemon/bipartite_matching.h Fri Apr 14 18:07:33 2006 +0000 @@ -19,8 +19,10 @@ #ifndef LEMON_BIPARTITE_MATCHING #define LEMON_BIPARTITE_MATCHING -#include -#include +#include + +#include +#include #include @@ -35,7 +37,7 @@ /// \brief Bipartite Max Cardinality Matching algorithm /// /// Bipartite Max Cardinality Matching algorithm. This class implements - /// the Hopcroft-Karp algorithm wich has \f$ O(e\sqrt{n}) \f$ time + /// the Hopcroft-Karp algorithm which has \f$ O(e\sqrt{n}) \f$ time /// complexity. template class MaxBipartiteMatching { @@ -83,7 +85,7 @@ for (BNodeIt it(*graph); it != INVALID; ++it) { bnode_matching[it] = INVALID; } - matching_value = 0; + matching_size = 0; } /// \brief Initalize the data structures. @@ -92,7 +94,7 @@ /// matching. From this matching sometimes it is faster to get /// the matching than from the initial empty matching. void greedyInit() { - matching_value = 0; + matching_size = 0; for (BNodeIt it(*graph); it != INVALID; ++it) { bnode_matching[it] = INVALID; } @@ -102,7 +104,7 @@ if (bnode_matching[graph->bNode(jt)] == INVALID) { anode_matching[it] = jt; bnode_matching[graph->bNode(jt)] = jt; - ++matching_value; + ++matching_size; break; } } @@ -120,10 +122,10 @@ for (BNodeIt it(*graph); it != INVALID; ++it) { bnode_matching[it] = INVALID; } - matching_value = 0; + matching_size = 0; for (UEdgeIt it(*graph); it != INVALID; ++it) { if (matching[it]) { - ++matching_value; + ++matching_size; anode_matching[graph->aNode(it)] = it; bnode_matching[graph->bNode(it)] = it; } @@ -142,10 +144,10 @@ for (BNodeIt it(*graph); it != INVALID; ++it) { bnode_matching[it] = INVALID; } - matching_value = 0; + matching_size = 0; for (UEdgeIt it(*graph); it != INVALID; ++it) { if (matching[it]) { - ++matching_value; + ++matching_size; if (anode_matching[graph->aNode(it)] != INVALID) { return false; } @@ -246,7 +248,7 @@ aused[anode] = true; } - ++matching_value; + ++matching_size; } } @@ -297,7 +299,7 @@ anode_matching[anode] = uedge; } - ++matching_value; + ++matching_size; return true; } else { Node newanode = graph->aNode(bnode_matching[bnode]); @@ -407,7 +409,7 @@ matching[anode_matching[it]] = true; } } - return matching_value; + return matching_size; } /// \brief Set true all matching uedge in the map and the others to false. @@ -419,7 +421,7 @@ for (UEdgeIt it(*graph); it != INVALID; ++it) { matching[it] = it == anode_matching[graph->aNode(it)]; } - return matching_value; + return matching_size; } @@ -445,8 +447,8 @@ /// \brief Gives back the number of the matching edges. /// /// Gives back the number of the matching edges. - int matchingValue() const { - return matching_value; + int matchingSize() const { + return matching_size; } /// @} @@ -457,7 +459,1022 @@ BNodeMatchingMap bnode_matching; const Graph *graph; - int matching_value; + int matching_size; + + }; + + /// \brief Default traits class for weighted bipartite matching algoritms. + /// + /// Default traits class for weighted bipartite matching algoritms. + /// \param _BpUGraph The bipartite undirected graph type. + /// \param _WeightMap Type of weight map. + template + struct WeightedBipartiteMatchingDefaultTraits { + /// \brief The type of the weight of the undirected edges. + typedef typename _WeightMap::Value Value; + + /// The undirected bipartite graph type the algorithm runs on. + typedef _BpUGraph BpUGraph; + + /// The map of the edges weights + typedef _WeightMap WeightMap; + + /// \brief The cross reference type used by heap. + /// + /// The cross reference type used by heap. + /// Usually it is \c Graph::NodeMap. + typedef typename BpUGraph::template NodeMap HeapCrossRef; + + /// \brief Instantiates a HeapCrossRef. + /// + /// This function instantiates a \ref HeapCrossRef. + /// \param graph is the graph, to which we would like to define the + /// HeapCrossRef. + static HeapCrossRef *createHeapCrossRef(const BpUGraph &graph) { + return new HeapCrossRef(graph); + } + + /// \brief The heap type used by weighted matching algorithms. + /// + /// The heap type used by weighted matching algorithms. It should + /// minimize the priorities and the heap's key type is the graph's + /// anode graph's node. + /// + /// \sa BinHeap + typedef BinHeap Heap; + + /// \brief Instantiates a Heap. + /// + /// This function instantiates a \ref Heap. + /// \param crossref The cross reference of the heap. + static Heap *createHeap(HeapCrossRef& crossref) { + return new Heap(crossref); + } + + }; + + + /// \ingroup matching + /// + /// \brief Bipartite Max Weighted Matching algorithm + /// + /// This class implements the bipartite Max Weighted Matching + /// algorithm. It uses the successive shortest path algorithm to + /// calculate the maximum weighted matching in the bipartite + /// graph. The algorithm can be used also to calculate the maximum + /// cardinality maximum weighted matching. The time complexity + /// of the algorithm is \f$ O(ne\log(n)) \f$ with the default binary + /// heap implementation but this can be improved to + /// \f$ O(n^2\log(n)+ne) \f$ if we use fibonacci heaps. + /// + /// The algorithm also provides a potential function on the nodes + /// which a dual solution of the matching algorithm and it can be + /// used to proof the optimality of the given pimal solution. +#ifdef DOXYGEN + template +#else + template , + typename _Traits = WeightedBipartiteMatchingDefaultTraits<_BpUGraph, _WeightMap> > +#endif + class MaxWeightedBipartiteMatching { + public: + + typedef _Traits Traits; + typedef typename Traits::BpUGraph BpUGraph; + typedef typename Traits::WeightMap WeightMap; + typedef typename Traits::Value Value; + + protected: + + typedef typename Traits::HeapCrossRef HeapCrossRef; + typedef typename Traits::Heap Heap; + + + typedef typename BpUGraph::Node Node; + typedef typename BpUGraph::ANodeIt ANodeIt; + typedef typename BpUGraph::BNodeIt BNodeIt; + typedef typename BpUGraph::UEdge UEdge; + typedef typename BpUGraph::UEdgeIt UEdgeIt; + typedef typename BpUGraph::IncEdgeIt IncEdgeIt; + + typedef typename BpUGraph::template ANodeMap ANodeMatchingMap; + typedef typename BpUGraph::template BNodeMap BNodeMatchingMap; + + typedef typename BpUGraph::template ANodeMap ANodePotentialMap; + typedef typename BpUGraph::template BNodeMap BNodePotentialMap; + + + public: + + /// \brief \ref Exception for uninitialized parameters. + /// + /// This error represents problems in the initialization + /// of the parameters of the algorithms. + class UninitializedParameter : public lemon::UninitializedParameter { + public: + virtual const char* exceptionName() const { + return "lemon::MaxWeightedBipartiteMatching::UninitializedParameter"; + } + }; + + ///\name Named template parameters + + ///@{ + + template + struct DefHeapTraits : public Traits { + typedef CR HeapCrossRef; + typedef H Heap; + static HeapCrossRef *createHeapCrossRef(const BpUGraph &) { + throw UninitializedParameter(); + } + static Heap *createHeap(HeapCrossRef &) { + throw UninitializedParameter(); + } + }; + + /// \brief \ref named-templ-param "Named parameter" for setting heap + /// and cross reference type + /// + /// \ref named-templ-param "Named parameter" for setting heap and cross + /// reference type + template > + struct DefHeap + : public MaxWeightedBipartiteMatching > { + typedef MaxWeightedBipartiteMatching > Create; + }; + + template + struct DefStandardHeapTraits : public Traits { + typedef CR HeapCrossRef; + typedef H Heap; + static HeapCrossRef *createHeapCrossRef(const BpUGraph &graph) { + return new HeapCrossRef(graph); + } + static Heap *createHeap(HeapCrossRef &crossref) { + return new Heap(crossref); + } + }; + + /// \brief \ref named-templ-param "Named parameter" for setting heap and + /// cross reference type with automatic allocation + /// + /// \ref named-templ-param "Named parameter" for setting heap and cross + /// reference type. It can allocate the heap and the cross reference + /// object if the cross reference's constructor waits for the graph as + /// parameter and the heap's constructor waits for the cross reference. + template > + struct DefStandardHeap + : public MaxWeightedBipartiteMatching > { + typedef MaxWeightedBipartiteMatching > + Create; + }; + + ///@} + + + /// \brief Constructor. + /// + /// Constructor of the algorithm. + MaxWeightedBipartiteMatching(const BpUGraph& _graph, + const WeightMap& _weight) + : graph(&_graph), weight(&_weight), + anode_matching(_graph), bnode_matching(_graph), + anode_potential(_graph), bnode_potential(_graph), + _heap_cross_ref(0), local_heap_cross_ref(false), + _heap(0), local_heap(0) {} + + /// \brief Destructor. + /// + /// Destructor of the algorithm. + ~MaxWeightedBipartiteMatching() { + destroyStructures(); + } + + /// \brief Sets the heap and the cross reference used by algorithm. + /// + /// Sets the heap and the cross reference used by algorithm. + /// If you don't use this function before calling \ref run(), + /// it will allocate one. The destuctor deallocates this + /// automatically allocated map, of course. + /// \return \c (*this) + MaxWeightedBipartiteMatching& heap(Heap& heap, HeapCrossRef &crossRef) { + if(local_heap_cross_ref) { + delete _heap_cross_ref; + local_heap_cross_ref = false; + } + _heap_cross_ref = &crossRef; + if(local_heap) { + delete _heap; + local_heap = false; + } + _heap = &heap; + return *this; + } + + /// \name Execution control + /// The simplest way to execute the algorithm is to use + /// one of the member functions called \c run(). + /// \n + /// If you need more control on the execution, + /// first you must call \ref init() or one alternative for it. + /// Finally \ref start() will perform the matching computation or + /// with step-by-step execution you can augment the solution. + + /// @{ + + /// \brief Initalize the data structures. + /// + /// It initalizes the data structures and creates an empty matching. + void init() { + initStructures(); + for (ANodeIt it(*graph); it != INVALID; ++it) { + anode_matching[it] = INVALID; + anode_potential[it] = 0; + } + for (BNodeIt it(*graph); it != INVALID; ++it) { + bnode_matching[it] = INVALID; + bnode_potential[it] = 0; + for (IncEdgeIt jt(*graph, it); jt != INVALID; ++jt) { + if (-(*weight)[jt] <= bnode_potential[it]) { + bnode_potential[it] = -(*weight)[jt]; + } + } + } + matching_value = 0; + matching_size = 0; + } + + + /// \brief An augmenting phase of the weighted matching algorithm + /// + /// It runs an augmenting phase of the weighted matching + /// algorithm. The phase finds the best augmenting path and + /// augments only on this paths. + /// + /// The algorithm consists at most + /// of \f$ O(n) \f$ phase and one phase is \f$ O(n\log(n)+e) \f$ + /// long with Fibonacci heap or \f$ O((n+e)\log(n)) \f$ long + /// with binary heap. + /// \param decrease If the given parameter true the matching value + /// can be decreased in the augmenting phase. If we would like + /// to calculate the maximum cardinality maximum weighted matching + /// then we should let the algorithm to decrease the matching + /// value in order to increase the number of the matching edges. + bool augment(bool decrease = false) { + + typename BpUGraph::template BNodeMap bdist(*graph); + typename BpUGraph::template BNodeMap bpred(*graph, INVALID); + + Node bestNode = INVALID; + Value bestValue = 0; + + _heap->clear(); + for (ANodeIt it(*graph); it != INVALID; ++it) { + (*_heap_cross_ref)[it] = Heap::PRE_HEAP; + } + + for (ANodeIt it(*graph); it != INVALID; ++it) { + if (anode_matching[it] == INVALID) { + _heap->push(it, 0); + } + } + + Value bdistMax = 0; + while (!_heap->empty()) { + Node anode = _heap->top(); + Value avalue = _heap->prio(); + _heap->pop(); + for (IncEdgeIt jt(*graph, anode); jt != INVALID; ++jt) { + if (jt == anode_matching[anode]) continue; + Node bnode = graph->bNode(jt); + Value bvalue = avalue + anode_potential[anode] - + bnode_potential[bnode] - (*weight)[jt]; + if (bpred[bnode] == INVALID || bvalue < bdist[bnode]) { + bdist[bnode] = bvalue; + bpred[bnode] = jt; + } + if (bvalue > bdistMax) { + bdistMax = bvalue; + } + if (bnode_matching[bnode] != INVALID) { + Node newanode = graph->aNode(bnode_matching[bnode]); + switch (_heap->state(newanode)) { + case Heap::PRE_HEAP: + _heap->push(newanode, bvalue); + break; + case Heap::IN_HEAP: + if (bvalue < (*_heap)[newanode]) { + _heap->decrease(newanode, bvalue); + } + break; + case Heap::POST_HEAP: + break; + } + } else { + if (bestNode == INVALID || + - bvalue - bnode_potential[bnode] > bestValue) { + bestValue = - bvalue - bnode_potential[bnode]; + bestNode = bnode; + } + } + } + } + + if (bestNode == INVALID || (!decrease && bestValue < 0)) { + return false; + } + + matching_value += bestValue; + ++matching_size; + + for (BNodeIt it(*graph); it != INVALID; ++it) { + if (bpred[it] != INVALID) { + bnode_potential[it] += bdist[it]; + } else { + bnode_potential[it] += bdistMax; + } + } + for (ANodeIt it(*graph); it != INVALID; ++it) { + if (anode_matching[it] != INVALID) { + Node bnode = graph->bNode(anode_matching[it]); + if (bpred[bnode] != INVALID) { + anode_potential[it] += bdist[bnode]; + } else { + anode_potential[it] += bdistMax; + } + } + } + + while (bestNode != INVALID) { + UEdge uedge = bpred[bestNode]; + Node anode = graph->aNode(uedge); + + bnode_matching[bestNode] = uedge; + if (anode_matching[anode] != INVALID) { + bestNode = graph->bNode(anode_matching[anode]); + } else { + bestNode = INVALID; + } + anode_matching[anode] = uedge; + } + + + return true; + } + + /// \brief Starts the algorithm. + /// + /// Starts the algorithm. It runs augmenting phases until the + /// optimal solution reached. + /// + /// \param maxCardinality If the given value is true it will + /// calculate the maximum cardinality maximum matching instead of + /// the maximum matching. + void start(bool maxCardinality = false) { + while (augment(maxCardinality)) {} + } + + /// \brief Runs the algorithm. + /// + /// It just initalize the algorithm and then start it. + /// + /// \param maxCardinality If the given value is true it will + /// calculate the maximum cardinality maximum matching instead of + /// the maximum matching. + void run(bool maxCardinality = false) { + init(); + start(maxCardinality); + } + + /// @} + + /// \name Query Functions + /// The result of the %Matching algorithm can be obtained using these + /// functions.\n + /// Before the use of these functions, + /// either run() or start() must be called. + + ///@{ + + /// \brief Gives back the potential in the NodeMap + /// + /// Gives back the potential in the NodeMap. The potential + /// is feasible if \f$ \pi(a) - \pi(b) - w(ab) = 0 \f$ for + /// each matching edges and \f$ \pi(a) - \pi(b) - w(ab) \ge 0 \f$ + /// for each edges. + template + void potential(PotentialMap& potential) { + for (ANodeIt it(*graph); it != INVALID; ++it) { + potential[it] = anode_potential[it]; + } + for (BNodeIt it(*graph); it != INVALID; ++it) { + potential[it] = bnode_potential[it]; + } + } + + /// \brief Set true all matching uedge in the map. + /// + /// Set true all matching uedge in the map. It does not change the + /// value mapped to the other uedges. + /// \return The number of the matching edges. + template + int quickMatching(MatchingMap& matching) { + for (ANodeIt it(*graph); it != INVALID; ++it) { + if (anode_matching[it] != INVALID) { + matching[anode_matching[it]] = true; + } + } + return matching_size; + } + + /// \brief Set true all matching uedge in the map and the others to false. + /// + /// Set true all matching uedge in the map and the others to false. + /// \return The number of the matching edges. + template + int matching(MatchingMap& matching) { + for (UEdgeIt it(*graph); it != INVALID; ++it) { + matching[it] = it == anode_matching[graph->aNode(it)]; + } + return matching_size; + } + + + /// \brief Return true if the given uedge is in the matching. + /// + /// It returns true if the given uedge is in the matching. + bool matchingEdge(const UEdge& edge) { + return anode_matching[graph->aNode(edge)] == edge; + } + + /// \brief Returns the matching edge from the node. + /// + /// Returns the matching edge from the node. If there is not such + /// edge it gives back \c INVALID. + UEdge matchingEdge(const Node& node) { + if (graph->aNode(node)) { + return anode_matching[node]; + } else { + return bnode_matching[node]; + } + } + + /// \brief Gives back the sum of weights of the matching edges. + /// + /// Gives back the sum of weights of the matching edges. + Value matchingValue() const { + return matching_value; + } + + /// \brief Gives back the number of the matching edges. + /// + /// Gives back the number of the matching edges. + int matchingSize() const { + return matching_size; + } + + /// @} + + private: + + void initStructures() { + if (!_heap_cross_ref) { + local_heap_cross_ref = true; + _heap_cross_ref = Traits::createHeapCrossRef(*graph); + } + if (!_heap) { + local_heap = true; + _heap = Traits::createHeap(*_heap_cross_ref); + } + } + + void destroyStructures() { + if (local_heap_cross_ref) delete _heap_cross_ref; + if (local_heap) delete _heap; + } + + + private: + + const BpUGraph *graph; + const WeightMap* weight; + + ANodeMatchingMap anode_matching; + BNodeMatchingMap bnode_matching; + + ANodePotentialMap anode_potential; + BNodePotentialMap bnode_potential; + + Value matching_value; + int matching_size; + + HeapCrossRef *_heap_cross_ref; + bool local_heap_cross_ref; + + Heap *_heap; + bool local_heap; + + }; + + /// \brief Default traits class for minimum cost bipartite matching + /// algoritms. + /// + /// Default traits class for minimum cost bipartite matching + /// algoritms. + /// + /// \param _BpUGraph The bipartite undirected graph + /// type. + /// + /// \param _CostMap Type of cost map. + template + struct MinCostMaxBipartiteMatchingDefaultTraits { + /// \brief The type of the cost of the undirected edges. + typedef typename _CostMap::Value Value; + + /// The undirected bipartite graph type the algorithm runs on. + typedef _BpUGraph BpUGraph; + + /// The map of the edges costs + typedef _CostMap CostMap; + + /// \brief The cross reference type used by heap. + /// + /// The cross reference type used by heap. + /// Usually it is \c Graph::NodeMap. + typedef typename BpUGraph::template NodeMap HeapCrossRef; + + /// \brief Instantiates a HeapCrossRef. + /// + /// This function instantiates a \ref HeapCrossRef. + /// \param graph is the graph, to which we would like to define the + /// HeapCrossRef. + static HeapCrossRef *createHeapCrossRef(const BpUGraph &graph) { + return new HeapCrossRef(graph); + } + + /// \brief The heap type used by costed matching algorithms. + /// + /// The heap type used by costed matching algorithms. It should + /// minimize the priorities and the heap's key type is the graph's + /// anode graph's node. + /// + /// \sa BinHeap + typedef BinHeap Heap; + + /// \brief Instantiates a Heap. + /// + /// This function instantiates a \ref Heap. + /// \param crossref The cross reference of the heap. + static Heap *createHeap(HeapCrossRef& crossref) { + return new Heap(crossref); + } + + }; + + + /// \ingroup matching + /// + /// \brief Bipartite Min Cost Matching algorithm + /// + /// This class implements the bipartite Min Cost Matching algorithm. + /// It uses the successive shortest path algorithm to calculate the + /// minimum cost maximum matching in the bipartite graph. The time + /// complexity of the algorithm is \f$ O(ne\log(n)) \f$ with the + /// default binary heap implementation but this can be improved to + /// \f$ O(n^2\log(n)+ne) \f$ if we use fibonacci heaps. + /// + /// The algorithm also provides a potential function on the nodes + /// which a dual solution of the matching algorithm and it can be + /// used to proof the optimality of the given pimal solution. +#ifdef DOXYGEN + template +#else + template , + typename _Traits = MinCostMaxBipartiteMatchingDefaultTraits<_BpUGraph, _CostMap> > +#endif + class MinCostMaxBipartiteMatching { + public: + + typedef _Traits Traits; + typedef typename Traits::BpUGraph BpUGraph; + typedef typename Traits::CostMap CostMap; + typedef typename Traits::Value Value; + + protected: + + typedef typename Traits::HeapCrossRef HeapCrossRef; + typedef typename Traits::Heap Heap; + + + typedef typename BpUGraph::Node Node; + typedef typename BpUGraph::ANodeIt ANodeIt; + typedef typename BpUGraph::BNodeIt BNodeIt; + typedef typename BpUGraph::UEdge UEdge; + typedef typename BpUGraph::UEdgeIt UEdgeIt; + typedef typename BpUGraph::IncEdgeIt IncEdgeIt; + + typedef typename BpUGraph::template ANodeMap ANodeMatchingMap; + typedef typename BpUGraph::template BNodeMap BNodeMatchingMap; + + typedef typename BpUGraph::template ANodeMap ANodePotentialMap; + typedef typename BpUGraph::template BNodeMap BNodePotentialMap; + + + public: + + /// \brief \ref Exception for uninitialized parameters. + /// + /// This error represents problems in the initialization + /// of the parameters of the algorithms. + class UninitializedParameter : public lemon::UninitializedParameter { + public: + virtual const char* exceptionName() const { + return "lemon::MinCostMaxBipartiteMatching::UninitializedParameter"; + } + }; + + ///\name Named template parameters + + ///@{ + + template + struct DefHeapTraits : public Traits { + typedef CR HeapCrossRef; + typedef H Heap; + static HeapCrossRef *createHeapCrossRef(const BpUGraph &) { + throw UninitializedParameter(); + } + static Heap *createHeap(HeapCrossRef &) { + throw UninitializedParameter(); + } + }; + + /// \brief \ref named-templ-param "Named parameter" for setting heap + /// and cross reference type + /// + /// \ref named-templ-param "Named parameter" for setting heap and cross + /// reference type + template > + struct DefHeap + : public MinCostMaxBipartiteMatching > { + typedef MinCostMaxBipartiteMatching > Create; + }; + + template + struct DefStandardHeapTraits : public Traits { + typedef CR HeapCrossRef; + typedef H Heap; + static HeapCrossRef *createHeapCrossRef(const BpUGraph &graph) { + return new HeapCrossRef(graph); + } + static Heap *createHeap(HeapCrossRef &crossref) { + return new Heap(crossref); + } + }; + + /// \brief \ref named-templ-param "Named parameter" for setting heap and + /// cross reference type with automatic allocation + /// + /// \ref named-templ-param "Named parameter" for setting heap and cross + /// reference type. It can allocate the heap and the cross reference + /// object if the cross reference's constructor waits for the graph as + /// parameter and the heap's constructor waits for the cross reference. + template > + struct DefStandardHeap + : public MinCostMaxBipartiteMatching > { + typedef MinCostMaxBipartiteMatching > + Create; + }; + + ///@} + + + /// \brief Constructor. + /// + /// Constructor of the algorithm. + MinCostMaxBipartiteMatching(const BpUGraph& _graph, + const CostMap& _cost) + : graph(&_graph), cost(&_cost), + anode_matching(_graph), bnode_matching(_graph), + anode_potential(_graph), bnode_potential(_graph), + _heap_cross_ref(0), local_heap_cross_ref(false), + _heap(0), local_heap(0) {} + + /// \brief Destructor. + /// + /// Destructor of the algorithm. + ~MinCostMaxBipartiteMatching() { + destroyStructures(); + } + + /// \brief Sets the heap and the cross reference used by algorithm. + /// + /// Sets the heap and the cross reference used by algorithm. + /// If you don't use this function before calling \ref run(), + /// it will allocate one. The destuctor deallocates this + /// automatically allocated map, of course. + /// \return \c (*this) + MinCostMaxBipartiteMatching& heap(Heap& heap, HeapCrossRef &crossRef) { + if(local_heap_cross_ref) { + delete _heap_cross_ref; + local_heap_cross_ref = false; + } + _heap_cross_ref = &crossRef; + if(local_heap) { + delete _heap; + local_heap = false; + } + _heap = &heap; + return *this; + } + + /// \name Execution control + /// The simplest way to execute the algorithm is to use + /// one of the member functions called \c run(). + /// \n + /// If you need more control on the execution, + /// first you must call \ref init() or one alternative for it. + /// Finally \ref start() will perform the matching computation or + /// with step-by-step execution you can augment the solution. + + /// @{ + + /// \brief Initalize the data structures. + /// + /// It initalizes the data structures and creates an empty matching. + void init() { + initStructures(); + for (ANodeIt it(*graph); it != INVALID; ++it) { + anode_matching[it] = INVALID; + anode_potential[it] = 0; + } + for (BNodeIt it(*graph); it != INVALID; ++it) { + bnode_matching[it] = INVALID; + bnode_potential[it] = 0; + } + matching_cost = 0; + matching_size = 0; + } + + + /// \brief An augmenting phase of the costed matching algorithm + /// + /// It runs an augmenting phase of the matching algorithm. The + /// phase finds the best augmenting path and augments only on this + /// paths. + /// + /// The algorithm consists at most + /// of \f$ O(n) \f$ phase and one phase is \f$ O(n\log(n)+e) \f$ + /// long with Fibonacci heap or \f$ O((n+e)\log(n)) \f$ long + /// with binary heap. + bool augment() { + + typename BpUGraph::template BNodeMap bdist(*graph); + typename BpUGraph::template BNodeMap bpred(*graph, INVALID); + + Node bestNode = INVALID; + Value bestValue = 0; + + _heap->clear(); + for (ANodeIt it(*graph); it != INVALID; ++it) { + (*_heap_cross_ref)[it] = Heap::PRE_HEAP; + } + + for (ANodeIt it(*graph); it != INVALID; ++it) { + if (anode_matching[it] == INVALID) { + _heap->push(it, 0); + } + } + + while (!_heap->empty()) { + Node anode = _heap->top(); + Value avalue = _heap->prio(); + _heap->pop(); + for (IncEdgeIt jt(*graph, anode); jt != INVALID; ++jt) { + if (jt == anode_matching[anode]) continue; + Node bnode = graph->bNode(jt); + Value bvalue = avalue + (*cost)[jt] + + anode_potential[anode] - bnode_potential[bnode]; + if (bpred[bnode] == INVALID || bvalue < bdist[bnode]) { + bdist[bnode] = bvalue; + bpred[bnode] = jt; + } + if (bnode_matching[bnode] != INVALID) { + Node newanode = graph->aNode(bnode_matching[bnode]); + switch (_heap->state(newanode)) { + case Heap::PRE_HEAP: + _heap->push(newanode, bvalue); + break; + case Heap::IN_HEAP: + if (bvalue < (*_heap)[newanode]) { + _heap->decrease(newanode, bvalue); + } + break; + case Heap::POST_HEAP: + break; + } + } else { + if (bestNode == INVALID || + bvalue + bnode_potential[bnode] < bestValue) { + bestValue = bvalue + bnode_potential[bnode]; + bestNode = bnode; + } + } + } + } + + if (bestNode == INVALID) { + return false; + } + + matching_cost += bestValue; + ++matching_size; + + for (BNodeIt it(*graph); it != INVALID; ++it) { + if (bpred[it] != INVALID) { + bnode_potential[it] += bdist[it]; + } + } + for (ANodeIt it(*graph); it != INVALID; ++it) { + if (anode_matching[it] != INVALID) { + Node bnode = graph->bNode(anode_matching[it]); + if (bpred[bnode] != INVALID) { + anode_potential[it] += bdist[bnode]; + } + } + } + + while (bestNode != INVALID) { + UEdge uedge = bpred[bestNode]; + Node anode = graph->aNode(uedge); + + bnode_matching[bestNode] = uedge; + if (anode_matching[anode] != INVALID) { + bestNode = graph->bNode(anode_matching[anode]); + } else { + bestNode = INVALID; + } + anode_matching[anode] = uedge; + } + + + return true; + } + + /// \brief Starts the algorithm. + /// + /// Starts the algorithm. It runs augmenting phases until the + /// optimal solution reached. + void start() { + while (augment()) {} + } + + /// \brief Runs the algorithm. + /// + /// It just initalize the algorithm and then start it. + void run() { + init(); + start(); + } + + /// @} + + /// \name Query Functions + /// The result of the %Matching algorithm can be obtained using these + /// functions.\n + /// Before the use of these functions, + /// either run() or start() must be called. + + ///@{ + + /// \brief Gives back the potential in the NodeMap + /// + /// Gives back the potential in the NodeMap. The potential + /// is feasible if \f$ \pi(a) - \pi(b) + w(ab) = 0 \f$ for + /// each matching edges and \f$ \pi(a) - \pi(b) + w(ab) \ge 0 \f$ + /// for each edges. + template + void potential(PotentialMap& potential) { + for (ANodeIt it(*graph); it != INVALID; ++it) { + potential[it] = anode_potential[it]; + } + for (BNodeIt it(*graph); it != INVALID; ++it) { + potential[it] = bnode_potential[it]; + } + } + + /// \brief Set true all matching uedge in the map. + /// + /// Set true all matching uedge in the map. It does not change the + /// value mapped to the other uedges. + /// \return The number of the matching edges. + template + int quickMatching(MatchingMap& matching) { + for (ANodeIt it(*graph); it != INVALID; ++it) { + if (anode_matching[it] != INVALID) { + matching[anode_matching[it]] = true; + } + } + return matching_size; + } + + /// \brief Set true all matching uedge in the map and the others to false. + /// + /// Set true all matching uedge in the map and the others to false. + /// \return The number of the matching edges. + template + int matching(MatchingMap& matching) { + for (UEdgeIt it(*graph); it != INVALID; ++it) { + matching[it] = it == anode_matching[graph->aNode(it)]; + } + return matching_size; + } + + + /// \brief Return true if the given uedge is in the matching. + /// + /// It returns true if the given uedge is in the matching. + bool matchingEdge(const UEdge& edge) { + return anode_matching[graph->aNode(edge)] == edge; + } + + /// \brief Returns the matching edge from the node. + /// + /// Returns the matching edge from the node. If there is not such + /// edge it gives back \c INVALID. + UEdge matchingEdge(const Node& node) { + if (graph->aNode(node)) { + return anode_matching[node]; + } else { + return bnode_matching[node]; + } + } + + /// \brief Gives back the sum of costs of the matching edges. + /// + /// Gives back the sum of costs of the matching edges. + Value matchingCost() const { + return matching_cost; + } + + /// \brief Gives back the number of the matching edges. + /// + /// Gives back the number of the matching edges. + int matchingSize() const { + return matching_size; + } + + /// @} + + private: + + void initStructures() { + if (!_heap_cross_ref) { + local_heap_cross_ref = true; + _heap_cross_ref = Traits::createHeapCrossRef(*graph); + } + if (!_heap) { + local_heap = true; + _heap = Traits::createHeap(*_heap_cross_ref); + } + } + + void destroyStructures() { + if (local_heap_cross_ref) delete _heap_cross_ref; + if (local_heap) delete _heap; + } + + + private: + + const BpUGraph *graph; + const CostMap* cost; + + ANodeMatchingMap anode_matching; + BNodeMatchingMap bnode_matching; + + ANodePotentialMap anode_potential; + BNodePotentialMap bnode_potential; + + Value matching_cost; + int matching_size; + + HeapCrossRef *_heap_cross_ref; + bool local_heap_cross_ref; + + Heap *_heap; + bool local_heap; };