LEMON/LEMON-official Changeset - r962:07ec2b52e53d

Repositories » LEMON » LEMON-official

Summary
Changelog
Files
Switch To
- loading...
Options
- Lightweight changelog
- Search
Pull Requests

Changeset - r962:07ec2b52e53d

« 914:5b926c
« 961:d8ea85

963:7f6e2b »

r962:07ec2b52e53d

Wed, 17 Mar 2010 10:29:57

[Not Reviewed]

0 comments (0 inline)

alpar (Alpar Juttner)
alpar@cs.elte.hu

Merge #356

0 2 0

merge default

2 files changed with 71 insertions and 2 deletions:

lemon/matching.h

lemon/unionfind.h

↑ Collapse diff ↑

lemon/matching.h

Show inline comments

@@ -45,3378 +45,3438 @@
 		  /// It can be started from an arbitrary initial matching
 		  /// (the default is the empty one).
 		  ///
 		  /// The dual solution of the problem is a map of the nodes to
 		  /// \ref MaxMatching::Status "Status", having values \c EVEN (or \c D),
 		  /// \c ODD (or \c A) and \c MATCHED (or \c C) defining the Gallai-Edmonds
 		  /// decomposition of the graph. The nodes in \c EVEN/D induce a subgraph
 		  /// with factor-critical components, the nodes in \c ODD/A form the
 		  /// canonical barrier, and the nodes in \c MATCHED/C induce a graph having
 		  /// a perfect matching. The number of the factor-critical components
 		  /// minus the number of barrier nodes is a lower bound on the
 		  /// unmatched nodes, and the matching is optimal if and only if this bound is
 		  /// tight. This decomposition can be obtained using \ref status() or
 		  /// \ref statusMap() after running the algorithm.
 		  ///
 		  /// \tparam GR The undirected graph type the algorithm runs on.
 		  template <typename GR>
 		  class MaxMatching {
 		  public:
 		    /// The graph type of the algorithm
 		    typedef GR Graph;
 		    /// The type of the matching map
 		    typedef typename Graph::template NodeMap<typename Graph::Arc>
 		    MatchingMap;
 		    ///\brief Status constants for Gallai-Edmonds decomposition.
 		    ///
 		    ///These constants are used for indicating the Gallai-Edmonds
 		    ///decomposition of a graph. The nodes with status \c EVEN (or \c D)
 		    ///induce a subgraph with factor-critical components, the nodes with
 		    ///status \c ODD (or \c A) form the canonical barrier, and the nodes
 		    ///with status \c MATCHED (or \c C) induce a subgraph having a
 		    ///perfect matching.
 		    enum Status {
 		      EVEN = 1,       ///< = 1. (\c D is an alias for \c EVEN.)
 		      D = 1,
 		      MATCHED = 0,    ///< = 0. (\c C is an alias for \c MATCHED.)
 		      C = 0,
 		      ODD = -1,       ///< = -1. (\c A is an alias for \c ODD.)
 		      A = -1,
 		      UNMATCHED = -2  ///< = -2.
 		    };
 		    /// The type of the status map
 		    typedef typename Graph::template NodeMap<Status> StatusMap;
 		  private:
 		    TEMPLATE_GRAPH_TYPEDEFS(Graph);
 		    typedef UnionFindEnum<IntNodeMap> BlossomSet;
 		    typedef ExtendFindEnum<IntNodeMap> TreeSet;
 		    typedef RangeMap<Node> NodeIntMap;
 		    typedef MatchingMap EarMap;
 		    typedef std::vector<Node> NodeQueue;
 		    const Graph& _graph;
 		    MatchingMap* _matching;
 		    StatusMap* _status;
 		    EarMap* _ear;
 		    IntNodeMap* _blossom_set_index;
 		    BlossomSet* _blossom_set;
 		    NodeIntMap* _blossom_rep;
 		    IntNodeMap* _tree_set_index;
 		    TreeSet* _tree_set;
 		    NodeQueue _node_queue;
 		    int _process, _postpone, _last;
 		    int _node_num;
 		  private:
 		    void createStructures() {
 		      _node_num = countNodes(_graph);
 		      if (!_matching) {
 		        _matching = new MatchingMap(_graph);
+		      }
 		      if (!_status) {
 		        _status = new StatusMap(_graph);
+		      }
 		      if (!_ear) {
 		        _ear = new EarMap(_graph);
+		      }
 		      if (!_blossom_set) {
 		        _blossom_set_index = new IntNodeMap(_graph);
 		        _blossom_set = new BlossomSet(*_blossom_set_index);
+		      }
 		      if (!_blossom_rep) {
 		        _blossom_rep = new NodeIntMap(_node_num);
+		      }
 		      if (!_tree_set) {
 		        _tree_set_index = new IntNodeMap(_graph);
 		        _tree_set = new TreeSet(*_tree_set_index);
+		      }
 		      _node_queue.resize(_node_num);
+		    }
 		    void destroyStructures() {
 		      if (_matching) {
 		        delete _matching;
+		      }
 		      if (_status) {
 		        delete _status;
+		      }
 		      if (_ear) {
 		        delete _ear;
+		      }
 		      if (_blossom_set) {
 		        delete _blossom_set;
 		        delete _blossom_set_index;
+		      }
 		      if (_blossom_rep) {
 		        delete _blossom_rep;
+		      }
 		      if (_tree_set) {
 		        delete _tree_set_index;
 		        delete _tree_set;
+		      }
+		    }
 		    void processDense(const Node& n) {
 		      _process = _postpone = _last = 0;
 		      _node_queue[_last++] = n;
 		      while (_process != _last) {
 		        Node u = _node_queue[_process++];
 		        for (OutArcIt a(_graph, u); a != INVALID; ++a) {
 		          Node v = _graph.target(a);
 		          if ((*_status)[v] == MATCHED) {
 		            extendOnArc(a);
 		          } else if ((*_status)[v] == UNMATCHED) {
 		            augmentOnArc(a);
 		            return;
+		          }
+		        }
+		      }
 		      while (_postpone != _last) {
 		        Node u = _node_queue[_postpone++];
 		        for (OutArcIt a(_graph, u); a != INVALID ; ++a) {
 		          Node v = _graph.target(a);
 		          if ((*_status)[v] == EVEN) {
 		            if (_blossom_set->find(u) != _blossom_set->find(v)) {
 		              shrinkOnEdge(a);
+		            }
+		          }
 		          while (_process != _last) {
 		            Node w = _node_queue[_process++];
 		            for (OutArcIt b(_graph, w); b != INVALID; ++b) {
 		              Node x = _graph.target(b);
 		              if ((*_status)[x] == MATCHED) {
 		                extendOnArc(b);
 		              } else if ((*_status)[x] == UNMATCHED) {
 		                augmentOnArc(b);
 		                return;
+		              }
+		            }
+		          }
+		        }
+		      }
+		    }
 		    void processSparse(const Node& n) {
 		      _process = _last = 0;
 		      _node_queue[_last++] = n;
 		      while (_process != _last) {
 		        Node u = _node_queue[_process++];
 		        for (OutArcIt a(_graph, u); a != INVALID; ++a) {
 		          Node v = _graph.target(a);
 		          if ((*_status)[v] == EVEN) {
 		            if (_blossom_set->find(u) != _blossom_set->find(v)) {
 		              shrinkOnEdge(a);
+		            }
 		          } else if ((*_status)[v] == MATCHED) {
 		            extendOnArc(a);
 		          } else if ((*_status)[v] == UNMATCHED) {
 		            augmentOnArc(a);
 		            return;
+		          }
+		        }
+		      }
+		    }
 		    void shrinkOnEdge(const Edge& e) {
 		      Node nca = INVALID;
+		      {
 		        std::set<Node> left_set, right_set;
 		        Node left = (*_blossom_rep)[_blossom_set->find(_graph.u(e))];
 		        left_set.insert(left);
 		        Node right = (*_blossom_rep)[_blossom_set->find(_graph.v(e))];
 		        right_set.insert(right);
 		        while (true) {
 		          if ((*_matching)[left] == INVALID) break;
 		          left = _graph.target((*_matching)[left]);
 		          left = (*_blossom_rep)[_blossom_set->
 		                                 find(_graph.target((*_ear)[left]))];
 		          if (right_set.find(left) != right_set.end()) {
 		            nca = left;
 		            break;
+		          }
 		          left_set.insert(left);
 		          if ((*_matching)[right] == INVALID) break;
 		          right = _graph.target((*_matching)[right]);
 		          right = (*_blossom_rep)[_blossom_set->
 		                                  find(_graph.target((*_ear)[right]))];
 		          if (left_set.find(right) != left_set.end()) {
 		            nca = right;
 		            break;
+		          }
 		          right_set.insert(right);
+		        }
 		        if (nca == INVALID) {
 		          if ((*_matching)[left] == INVALID) {
 		            nca = right;
 		            while (left_set.find(nca) == left_set.end()) {
 		              nca = _graph.target((*_matching)[nca]);
 		              nca =(*_blossom_rep)[_blossom_set->
 		                                   find(_graph.target((*_ear)[nca]))];
+		            }
 		          } else {
 		            nca = left;
 		            while (right_set.find(nca) == right_set.end()) {
 		              nca = _graph.target((*_matching)[nca]);
 		              nca = (*_blossom_rep)[_blossom_set->
 		                                   find(_graph.target((*_ear)[nca]))];
+		            }
+		          }
+		        }
+		      }
+		      {
 		        Node node = _graph.u(e);
 		        Arc arc = _graph.direct(e, true);
 		        Node base = (*_blossom_rep)[_blossom_set->find(node)];
 		        while (base != nca) {
 		          (*_ear)[node] = arc;
 		          Node n = node;
 		          while (n != base) {
 		            n = _graph.target((*_matching)[n]);
 		            Arc a = (*_ear)[n];
 		            n = _graph.target(a);
 		            (*_ear)[n] = _graph.oppositeArc(a);
+		          }
 		          node = _graph.target((*_matching)[base]);
 		          _tree_set->erase(base);
 		          _tree_set->erase(node);
 		          _blossom_set->insert(node, _blossom_set->find(base));
 		          (*_status)[node] = EVEN;
 		          _node_queue[_last++] = node;
 		          arc = _graph.oppositeArc((*_ear)[node]);
 		          node = _graph.target((*_ear)[node]);
 		          base = (*_blossom_rep)[_blossom_set->find(node)];
 		          _blossom_set->join(_graph.target(arc), base);
+		        }
+		      }
 		      (*_blossom_rep)[_blossom_set->find(nca)] = nca;
+		      {
 		        Node node = _graph.v(e);
 		        Arc arc = _graph.direct(e, false);
 		        Node base = (*_blossom_rep)[_blossom_set->find(node)];
 		        while (base != nca) {
 		          (*_ear)[node] = arc;
 		          Node n = node;
 		          while (n != base) {
 		            n = _graph.target((*_matching)[n]);
 		            Arc a = (*_ear)[n];
 		            n = _graph.target(a);
 		            (*_ear)[n] = _graph.oppositeArc(a);
+		          }
 		          node = _graph.target((*_matching)[base]);
 		          _tree_set->erase(base);
 		          _tree_set->erase(node);
 		          _blossom_set->insert(node, _blossom_set->find(base));
 		          (*_status)[node] = EVEN;
 		          _node_queue[_last++] = node;
 		          arc = _graph.oppositeArc((*_ear)[node]);
 		          node = _graph.target((*_ear)[node]);
 		          base = (*_blossom_rep)[_blossom_set->find(node)];
 		          _blossom_set->join(_graph.target(arc), base);
+		        }
+		      }
 		      (*_blossom_rep)[_blossom_set->find(nca)] = nca;
+		    }
 		    void extendOnArc(const Arc& a) {
 		      Node base = _graph.source(a);
 		      Node odd = _graph.target(a);
 		      (*_ear)[odd] = _graph.oppositeArc(a);
 		      Node even = _graph.target((*_matching)[odd]);
 		      (*_blossom_rep)[_blossom_set->insert(even)] = even;
 		      (*_status)[odd] = ODD;
 		      (*_status)[even] = EVEN;
 		      int tree = _tree_set->find((*_blossom_rep)[_blossom_set->find(base)]);
 		      _tree_set->insert(odd, tree);
 		      _tree_set->insert(even, tree);
 		      _node_queue[_last++] = even;
+		    }
 		    void augmentOnArc(const Arc& a) {
 		      Node even = _graph.source(a);
 		      Node odd = _graph.target(a);
 		      int tree = _tree_set->find((*_blossom_rep)[_blossom_set->find(even)]);
 		      (*_matching)[odd] = _graph.oppositeArc(a);
 		      (*_status)[odd] = MATCHED;
 		      Arc arc = (*_matching)[even];
 		      (*_matching)[even] = a;
 		      while (arc != INVALID) {
 		        odd = _graph.target(arc);
 		        arc = (*_ear)[odd];
 		        even = _graph.target(arc);
 		        (*_matching)[odd] = arc;
 		        arc = (*_matching)[even];
 		        (*_matching)[even] = _graph.oppositeArc((*_matching)[odd]);
+		      }
 		      for (typename TreeSet::ItemIt it(*_tree_set, tree);
 		           it != INVALID; ++it) {
 		        if ((*_status)[it] == ODD) {
 		          (*_status)[it] = MATCHED;
 		        } else {
 		          int blossom = _blossom_set->find(it);
 		          for (typename BlossomSet::ItemIt jt(*_blossom_set, blossom);
 		               jt != INVALID; ++jt) {
 		            (*_status)[jt] = MATCHED;
+		          }
 		          _blossom_set->eraseClass(blossom);
+		        }
+		      }
 		      _tree_set->eraseClass(tree);
+		    }
 		  public:
 		    /// \brief Constructor
 		    ///
 		    /// Constructor.
 		    MaxMatching(const Graph& graph)
 		      : _graph(graph), _matching(0), _status(0), _ear(0),
 		        _blossom_set_index(0), _blossom_set(0), _blossom_rep(0),
 		        _tree_set_index(0), _tree_set(0) {}
 		    ~MaxMatching() {
 		      destroyStructures();
+		    }
 		    /// \name Execution Control
 		    /// The simplest way to execute the algorithm is to use the
 		    /// \c run() member function.\n
 		    /// If you need better control on the execution, you have to call
 		    /// one of the functions \ref init(), \ref greedyInit() or
 		    /// \ref matchingInit() first, then you can start the algorithm with
 		    /// \ref startSparse() or \ref startDense().
 		    ///@{
 		    /// \brief Set the initial matching to the empty matching.
 		    ///
 		    /// This function sets the initial matching to the empty matching.
 		    void init() {
 		      createStructures();
 		      for(NodeIt n(_graph); n != INVALID; ++n) {
 		        (*_matching)[n] = INVALID;
 		        (*_status)[n] = UNMATCHED;
+		      }
+		    }
 		    /// \brief Find an initial matching in a greedy way.
 		    ///
 		    /// This function finds an initial matching in a greedy way.
 		    void greedyInit() {
 		      createStructures();
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        (*_matching)[n] = INVALID;
 		        (*_status)[n] = UNMATCHED;
+		      }
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        if ((*_matching)[n] == INVALID) {
 		          for (OutArcIt a(_graph, n); a != INVALID ; ++a) {
 		            Node v = _graph.target(a);
 		            if ((*_matching)[v] == INVALID && v != n) {
 		              (*_matching)[n] = a;
 		              (*_status)[n] = MATCHED;
 		              (*_matching)[v] = _graph.oppositeArc(a);
 		              (*_status)[v] = MATCHED;
 		              break;
+		            }
+		          }
+		        }
+		      }
+		    }
 		    /// \brief Initialize the matching from a map.
 		    ///
 		    /// This function initializes the matching from a \c bool valued edge
 		    /// map. This map should have the property that there are no two incident
 		    /// edges with \c true value, i.e. it really contains a matching.
 		    /// \return \c true if the map contains a matching.
 		    template <typename MatchingMap>
 		    bool matchingInit(const MatchingMap& matching) {
 		      createStructures();
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        (*_matching)[n] = INVALID;
 		        (*_status)[n] = UNMATCHED;
+		      }
 		      for(EdgeIt e(_graph); e!=INVALID; ++e) {
 		        if (matching[e]) {
 		          Node u = _graph.u(e);
 		          if ((*_matching)[u] != INVALID) return false;
 		          (*_matching)[u] = _graph.direct(e, true);
 		          (*_status)[u] = MATCHED;
 		          Node v = _graph.v(e);
 		          if ((*_matching)[v] != INVALID) return false;
 		          (*_matching)[v] = _graph.direct(e, false);
 		          (*_status)[v] = MATCHED;
+		        }
+		      }
 		      return true;
+		    }
 		    /// \brief Start Edmonds' algorithm
 		    ///
 		    /// This function runs the original Edmonds' algorithm.
 		    ///
 		    /// \pre \ref init(), \ref greedyInit() or \ref matchingInit() must be
 		    /// called before using this function.
 		    void startSparse() {
 		      for(NodeIt n(_graph); n != INVALID; ++n) {
 		        if ((*_status)[n] == UNMATCHED) {
 		          (*_blossom_rep)[_blossom_set->insert(n)] = n;
 		          _tree_set->insert(n);
 		          (*_status)[n] = EVEN;
 		          processSparse(n);
+		        }
+		      }
+		    }
 		    /// \brief Start Edmonds' algorithm with a heuristic improvement
 		    /// for dense graphs
 		    ///
 		    /// This function runs Edmonds' algorithm with a heuristic of postponing
 		    /// shrinks, therefore resulting in a faster algorithm for dense graphs.
 		    ///
 		    /// \pre \ref init(), \ref greedyInit() or \ref matchingInit() must be
 		    /// called before using this function.
 		    void startDense() {
 		      for(NodeIt n(_graph); n != INVALID; ++n) {
 		        if ((*_status)[n] == UNMATCHED) {
 		          (*_blossom_rep)[_blossom_set->insert(n)] = n;
 		          _tree_set->insert(n);
 		          (*_status)[n] = EVEN;
 		          processDense(n);
+		        }
+		      }
+		    }
 		    /// \brief Run Edmonds' algorithm
 		    ///
 		    /// This function runs Edmonds' algorithm. An additional heuristic of
 		    /// postponing shrinks is used for relatively dense graphs
 		    /// (for which <tt>m>=2*n</tt> holds).
 		    void run() {
 		      if (countEdges(_graph) < 2 * countNodes(_graph)) {
 		        greedyInit();
 		        startSparse();
 		      } else {
 		        init();
 		        startDense();
+		      }
+		    }
 		    /// @}
 		    /// \name Primal Solution
 		    /// Functions to get the primal solution, i.e. the maximum matching.
 		    /// @{
 		    /// \brief Return the size (cardinality) of the matching.
 		    ///
 		    /// This function returns the size (cardinality) of the current matching.
 		    /// After run() it returns the size of the maximum matching in the graph.
 		    int matchingSize() const {
 		      int size = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        if ((*_matching)[n] != INVALID) {
 		          ++size;
+		        }
+		      }
 		      return size / 2;
+		    }
 		    /// \brief Return \c true if the given edge is in the matching.
 		    ///
 		    /// This function returns \c true if the given edge is in the current
 		    /// matching.
 		    bool matching(const Edge& edge) const {
 		      return edge == (*_matching)[_graph.u(edge)];
+		    }
 		    /// \brief Return the matching arc (or edge) incident to the given node.
 		    ///
 		    /// This function returns the matching arc (or edge) incident to the
 		    /// given node in the current matching or \c INVALID if the node is
 		    /// not covered by the matching.
 		    Arc matching(const Node& n) const {
 		      return (*_matching)[n];
+		    }
 		    /// \brief Return a const reference to the matching map.
 		    ///
 		    /// This function returns a const reference to a node map that stores
 		    /// the matching arc (or edge) incident to each node.
 		    const MatchingMap& matchingMap() const {
 		      return *_matching;
+		    }
 		    /// \brief Return the mate of the given node.
 		    ///
 		    /// This function returns the mate of the given node in the current
 		    /// matching or \c INVALID if the node is not covered by the matching.
 		    Node mate(const Node& n) const {
 		      return (*_matching)[n] != INVALID ?
 		        _graph.target((*_matching)[n]) : INVALID;
+		    }
 		    /// @}
 		    /// \name Dual Solution
 		    /// Functions to get the dual solution, i.e. the Gallai-Edmonds
 		    /// decomposition.
 		    /// @{
 		    /// \brief Return the status of the given node in the Edmonds-Gallai
 		    /// decomposition.
 		    ///
 		    /// This function returns the \ref Status "status" of the given node
 		    /// in the Edmonds-Gallai decomposition.
 		    Status status(const Node& n) const {
 		      return (*_status)[n];
+		    }
 		    /// \brief Return a const reference to the status map, which stores
 		    /// the Edmonds-Gallai decomposition.
 		    ///
 		    /// This function returns a const reference to a node map that stores the
 		    /// \ref Status "status" of each node in the Edmonds-Gallai decomposition.
 		    const StatusMap& statusMap() const {
 		      return *_status;
+		    }
 		    /// \brief Return \c true if the given node is in the barrier.
 		    ///
 		    /// This function returns \c true if the given node is in the barrier.
 		    bool barrier(const Node& n) const {
 		      return (*_status)[n] == ODD;
+		    }
 		    /// @}
 		  };
 		  /// \ingroup matching
 		  ///
 		  /// \brief Weighted matching in general graphs
 		  ///
 		  /// This class provides an efficient implementation of Edmond's
 		  /// maximum weighted matching algorithm. The implementation is based
 		  /// on extensive use of priority queues and provides
 		  /// \f$O(nm\log n)\f$ time complexity.
 		  ///
 		  /// The maximum weighted matching problem is to find a subset of the
 		  /// edges in an undirected graph with maximum overall weight for which
 		  /// each node has at most one incident edge.
 		  /// It can be formulated with the following linear program.
 		  /// \f[ \sum_{e \in \delta(u)}x_e \le 1 \quad \forall u\in V\f]
 		  /** \f[ \sum_{e \in \gamma(B)}x_e \le \frac{\vert B \vert - 1}{2}
 		      \quad \forall B\in\mathcal{O}\f] */
 		  /// \f[x_e \ge 0\quad \forall e\in E\f]
 		  /// \f[\max \sum_{e\in E}x_ew_e\f]
 		  /// where \f$\delta(X)\f$ is the set of edges incident to a node in
 		  /// \f$X\f$, \f$\gamma(X)\f$ is the set of edges with both ends in
 		  /// \f$X\f$ and \f$\mathcal{O}\f$ is the set of odd cardinality
 		  /// subsets of the nodes.
 		  ///
 		  /// The algorithm calculates an optimal matching and a proof of the
 		  /// optimality. The solution of the dual problem can be used to check
 		  /// the result of the algorithm. The dual linear problem is the
 		  /// following.
 		  /** \f[ y_u + y_v + \sum_{B \in \mathcal{O}, uv \in \gamma(B)}
 		      z_B \ge w_{uv} \quad \forall uv\in E\f] */
 		  /// \f[y_u \ge 0 \quad \forall u \in V\f]
 		  /// \f[z_B \ge 0 \quad \forall B \in \mathcal{O}\f]
 		  /** \f[\min \sum_{u \in V}y_u + \sum_{B \in \mathcal{O}}
 		      \frac{\vert B \vert - 1}{2}z_B\f] */
 		  ///
 		  /// The algorithm can be executed with the run() function.
 		  /// After it the matching (the primal solution) and the dual solution
 		  /// can be obtained using the query functions and the
 		  /// \ref MaxWeightedMatching::BlossomIt "BlossomIt" nested class,
 		  /// which is able to iterate on the nodes of a blossom.
 		  /// If the value type is integer, then the dual solution is multiplied
 		  /// by \ref MaxWeightedMatching::dualScale "4".
 		  ///
 		  /// \tparam GR The undirected graph type the algorithm runs on.
 		  /// \tparam WM The type edge weight map. The default type is
 		  /// \ref concepts::Graph::EdgeMap "GR::EdgeMap<int>".
 		#ifdef DOXYGEN
 		  template <typename GR, typename WM>
 		#else
 		  template <typename GR,
 		            typename WM = typename GR::template EdgeMap<int> >
 		#endif
 		  class MaxWeightedMatching {
 		  public:
 		    /// The graph type of the algorithm
 		    typedef GR Graph;
 		    /// The type of the edge weight map
 		    typedef WM WeightMap;
 		    /// The value type of the edge weights
 		    typedef typename WeightMap::Value Value;
 		    /// The type of the matching map
 		    typedef typename Graph::template NodeMap<typename Graph::Arc>
 		    MatchingMap;
 		    /// \brief Scaling factor for dual solution
 		    ///
 		    /// Scaling factor for dual solution. It is equal to 4 or 1
 		    /// according to the value type.
 		    static const int dualScale =
 		      std::numeric_limits<Value>::is_integer ? 4 : 1;
 		  private:
 		    TEMPLATE_GRAPH_TYPEDEFS(Graph);
 		    typedef typename Graph::template NodeMap<Value> NodePotential;
 		    typedef std::vector<Node> BlossomNodeList;
 		    struct BlossomVariable {
 		      int begin, end;
 		      Value value;
 		      BlossomVariable(int _begin, int _end, Value _value)
 		        : begin(_begin), end(_end), value(_value) {}
 		    };
 		    typedef std::vector<BlossomVariable> BlossomPotential;
 		    const Graph& _graph;
 		    const WeightMap& _weight;
 		    MatchingMap* _matching;
 		    NodePotential* _node_potential;
 		    BlossomPotential _blossom_potential;
 		    BlossomNodeList _blossom_node_list;
 		    int _node_num;
 		    int _blossom_num;
 		    typedef RangeMap<int> IntIntMap;
 		    enum Status {
 		      EVEN = -1, MATCHED = 0, ODD = 1
 		    };
 		    typedef HeapUnionFind<Value, IntNodeMap> BlossomSet;
 		    struct BlossomData {
 		      int tree;
 		      Status status;
 		      Arc pred, next;
 		      Value pot, offset;
 		      Node base;
 		    };
 		    IntNodeMap *_blossom_index;
 		    BlossomSet *_blossom_set;
 		    RangeMap<BlossomData>* _blossom_data;
 		    IntNodeMap *_node_index;
 		    IntArcMap *_node_heap_index;
 		    struct NodeData {
 		      NodeData(IntArcMap& node_heap_index)
 		        : heap(node_heap_index) {}
 		      int blossom;
 		      Value pot;
 		      BinHeap<Value, IntArcMap> heap;
 		      std::map<int, Arc> heap_index;
 		      int tree;
 		    };
 		    RangeMap<NodeData>* _node_data;
 		    typedef ExtendFindEnum<IntIntMap> TreeSet;
 		    IntIntMap *_tree_set_index;
 		    TreeSet *_tree_set;
 		    IntNodeMap *_delta1_index;
 		    BinHeap<Value, IntNodeMap> *_delta1;
 		    IntIntMap *_delta2_index;
 		    BinHeap<Value, IntIntMap> *_delta2;
 		    IntEdgeMap *_delta3_index;
 		    BinHeap<Value, IntEdgeMap> *_delta3;
 		    IntIntMap *_delta4_index;
 		    BinHeap<Value, IntIntMap> *_delta4;
 		    Value _delta_sum;
 		    int _unmatched;
 		    typedef MaxWeightedFractionalMatching<Graph, WeightMap> FractionalMatching;
 		    FractionalMatching *_fractional;
 		    void createStructures() {
 		      _node_num = countNodes(_graph);
 		      _blossom_num = _node_num * 3 / 2;
 		      if (!_matching) {
 		        _matching = new MatchingMap(_graph);
+		      }
 		      if (!_node_potential) {
 		        _node_potential = new NodePotential(_graph);
+		      }
 		      if (!_blossom_set) {
 		        _blossom_index = new IntNodeMap(_graph);
 		        _blossom_set = new BlossomSet(*_blossom_index);
 		        _blossom_data = new RangeMap<BlossomData>(_blossom_num);
 		      } else if (_blossom_data->size() != _blossom_num) {
 		        delete _blossom_data;
 		        _blossom_data = new RangeMap<BlossomData>(_blossom_num);
+		      }
 		      if (!_node_index) {
 		        _node_index = new IntNodeMap(_graph);
 		        _node_heap_index = new IntArcMap(_graph);
 		        _node_data = new RangeMap<NodeData>(_node_num,
-		                                              NodeData(*_node_heap_index));
 		                                            NodeData(*_node_heap_index));
 		      } else {
 		        delete _node_data;
 		        _node_data = new RangeMap<NodeData>(_node_num,
 		                                            NodeData(*_node_heap_index));
+		      }
 		      if (!_tree_set) {
 		        _tree_set_index = new IntIntMap(_blossom_num);
 		        _tree_set = new TreeSet(*_tree_set_index);
 		      } else {
 		        _tree_set_index->resize(_blossom_num);
+		      }
 		      if (!_delta1) {
 		        _delta1_index = new IntNodeMap(_graph);
 		        _delta1 = new BinHeap<Value, IntNodeMap>(*_delta1_index);
+		      }
 		      if (!_delta2) {
 		        _delta2_index = new IntIntMap(_blossom_num);
 		        _delta2 = new BinHeap<Value, IntIntMap>(*_delta2_index);
 		      } else {
 		        _delta2_index->resize(_blossom_num);
+		      }
 		      if (!_delta3) {
 		        _delta3_index = new IntEdgeMap(_graph);
 		        _delta3 = new BinHeap<Value, IntEdgeMap>(*_delta3_index);
+		      }
 		      if (!_delta4) {
 		        _delta4_index = new IntIntMap(_blossom_num);
 		        _delta4 = new BinHeap<Value, IntIntMap>(*_delta4_index);
 		      } else {
 		        _delta4_index->resize(_blossom_num);
+		      }
+		    }
 		    void destroyStructures() {
 		      if (_matching) {
 		        delete _matching;
+		      }
 		      if (_node_potential) {
 		        delete _node_potential;
+		      }
 		      if (_blossom_set) {
 		        delete _blossom_index;
 		        delete _blossom_set;
 		        delete _blossom_data;
+		      }
 		      if (_node_index) {
 		        delete _node_index;
 		        delete _node_heap_index;
 		        delete _node_data;
+		      }
 		      if (_tree_set) {
 		        delete _tree_set_index;
 		        delete _tree_set;
+		      }
 		      if (_delta1) {
 		        delete _delta1_index;
 		        delete _delta1;
+		      }
 		      if (_delta2) {
 		        delete _delta2_index;
 		        delete _delta2;
+		      }
 		      if (_delta3) {
 		        delete _delta3_index;
 		        delete _delta3;
+		      }
 		      if (_delta4) {
 		        delete _delta4_index;
 		        delete _delta4;
+		      }
+		    }
 		    void matchedToEven(int blossom, int tree) {
 		      if (_delta2->state(blossom) == _delta2->IN_HEAP) {
 		        _delta2->erase(blossom);
+		      }
 		      if (!_blossom_set->trivial(blossom)) {
 		        (*_blossom_data)[blossom].pot -=
 * (_delta_sum - (*_blossom_data)[blossom].offset);
+		      }
 		      for (typename BlossomSet::ItemIt n(*_blossom_set, blossom);
 		           n != INVALID; ++n) {
 		        _blossom_set->increase(n, std::numeric_limits<Value>::max());
 		        int ni = (*_node_index)[n];
 		        (*_node_data)[ni].heap.clear();
 		        (*_node_data)[ni].heap_index.clear();
 		        (*_node_data)[ni].pot += _delta_sum - (*_blossom_data)[blossom].offset;
 		        _delta1->push(n, (*_node_data)[ni].pot);
 		        for (InArcIt e(_graph, n); e != INVALID; ++e) {
 		          Node v = _graph.source(e);
 		          int vb = _blossom_set->find(v);
 		          int vi = (*_node_index)[v];
 		          Value rw = (*_node_data)[ni].pot + (*_node_data)[vi].pot -
 		            dualScale * _weight[e];
 		          if ((*_blossom_data)[vb].status == EVEN) {
 		            if (_delta3->state(e) != _delta3->IN_HEAP && blossom != vb) {
 		              _delta3->push(e, rw / 2);
+		            }
 		          } else {
 		            typename std::map<int, Arc>::iterator it =
 		              (*_node_data)[vi].heap_index.find(tree);
 		            if (it != (*_node_data)[vi].heap_index.end()) {
 		              if ((*_node_data)[vi].heap[it->second] > rw) {
 		                (*_node_data)[vi].heap.replace(it->second, e);
 		                (*_node_data)[vi].heap.decrease(e, rw);
 		                it->second = e;
+		              }
 		            } else {
 		              (*_node_data)[vi].heap.push(e, rw);
 		              (*_node_data)[vi].heap_index.insert(std::make_pair(tree, e));
+		            }
 		            if ((*_blossom_set)[v] > (*_node_data)[vi].heap.prio()) {
 		              _blossom_set->decrease(v, (*_node_data)[vi].heap.prio());
 		              if ((*_blossom_data)[vb].status == MATCHED) {
 		                if (_delta2->state(vb) != _delta2->IN_HEAP) {
 		                  _delta2->push(vb, _blossom_set->classPrio(vb) -
 		                               (*_blossom_data)[vb].offset);
 		                } else if ((*_delta2)[vb] > _blossom_set->classPrio(vb) -
 		                           (*_blossom_data)[vb].offset) {
 		                  _delta2->decrease(vb, _blossom_set->classPrio(vb) -
 		                                   (*_blossom_data)[vb].offset);
+		                }
+		              }
+		            }
+		          }
+		        }
+		      }
 		      (*_blossom_data)[blossom].offset = 0;
+		    }
 		    void matchedToOdd(int blossom) {
 		      if (_delta2->state(blossom) == _delta2->IN_HEAP) {
 		        _delta2->erase(blossom);
+		      }
 		      (*_blossom_data)[blossom].offset += _delta_sum;
 		      if (!_blossom_set->trivial(blossom)) {
 		        _delta4->push(blossom, (*_blossom_data)[blossom].pot / 2 +
 		                      (*_blossom_data)[blossom].offset);
+		      }
+		    }
 		    void evenToMatched(int blossom, int tree) {
 		      if (!_blossom_set->trivial(blossom)) {
 		        (*_blossom_data)[blossom].pot += 2 * _delta_sum;
+		      }
 		      for (typename BlossomSet::ItemIt n(*_blossom_set, blossom);
 		           n != INVALID; ++n) {
 		        int ni = (*_node_index)[n];
 		        (*_node_data)[ni].pot -= _delta_sum;
 		        _delta1->erase(n);
 		        for (InArcIt e(_graph, n); e != INVALID; ++e) {
 		          Node v = _graph.source(e);
 		          int vb = _blossom_set->find(v);
 		          int vi = (*_node_index)[v];
 		          Value rw = (*_node_data)[ni].pot + (*_node_data)[vi].pot -
 		            dualScale * _weight[e];
 		          if (vb == blossom) {
 		            if (_delta3->state(e) == _delta3->IN_HEAP) {
 		              _delta3->erase(e);
+		            }
 		          } else if ((*_blossom_data)[vb].status == EVEN) {
 		            if (_delta3->state(e) == _delta3->IN_HEAP) {
 		              _delta3->erase(e);
+		            }
 		            int vt = _tree_set->find(vb);
 		            if (vt != tree) {
 		              Arc r = _graph.oppositeArc(e);
 		              typename std::map<int, Arc>::iterator it =
 		                (*_node_data)[ni].heap_index.find(vt);
 		              if (it != (*_node_data)[ni].heap_index.end()) {
 		                if ((*_node_data)[ni].heap[it->second] > rw) {
 		                  (*_node_data)[ni].heap.replace(it->second, r);
 		                  (*_node_data)[ni].heap.decrease(r, rw);
 		                  it->second = r;
+		                }
 		              } else {
 		                (*_node_data)[ni].heap.push(r, rw);
 		                (*_node_data)[ni].heap_index.insert(std::make_pair(vt, r));
+		              }
 		              if ((*_blossom_set)[n] > (*_node_data)[ni].heap.prio()) {
 		                _blossom_set->decrease(n, (*_node_data)[ni].heap.prio());
 		                if (_delta2->state(blossom) != _delta2->IN_HEAP) {
 		                  _delta2->push(blossom, _blossom_set->classPrio(blossom) -
 		                               (*_blossom_data)[blossom].offset);
 		                } else if ((*_delta2)[blossom] >
 		                           _blossom_set->classPrio(blossom) -
 		                           (*_blossom_data)[blossom].offset){
 		                  _delta2->decrease(blossom, _blossom_set->classPrio(blossom) -
 		                                   (*_blossom_data)[blossom].offset);
+		                }
+		              }
+		            }
 		          } else {
 		            typename std::map<int, Arc>::iterator it =
 		              (*_node_data)[vi].heap_index.find(tree);
 		            if (it != (*_node_data)[vi].heap_index.end()) {
 		              (*_node_data)[vi].heap.erase(it->second);
 		              (*_node_data)[vi].heap_index.erase(it);
 		              if ((*_node_data)[vi].heap.empty()) {
 		                _blossom_set->increase(v, std::numeric_limits<Value>::max());
 		              } else if ((*_blossom_set)[v] < (*_node_data)[vi].heap.prio()) {
 		                _blossom_set->increase(v, (*_node_data)[vi].heap.prio());
+		              }
 		              if ((*_blossom_data)[vb].status == MATCHED) {
 		                if (_blossom_set->classPrio(vb) ==
 		                    std::numeric_limits<Value>::max()) {
 		                  _delta2->erase(vb);
 		                } else if ((*_delta2)[vb] < _blossom_set->classPrio(vb) -
 		                           (*_blossom_data)[vb].offset) {
 		                  _delta2->increase(vb, _blossom_set->classPrio(vb) -
 		                                   (*_blossom_data)[vb].offset);
+		                }
+		              }
+		            }
+		          }
+		        }
+		      }
+		    }
 		    void oddToMatched(int blossom) {
 		      (*_blossom_data)[blossom].offset -= _delta_sum;
 		      if (_blossom_set->classPrio(blossom) !=
 		          std::numeric_limits<Value>::max()) {
 		        _delta2->push(blossom, _blossom_set->classPrio(blossom) -
 		                      (*_blossom_data)[blossom].offset);
+		      }
 		      if (!_blossom_set->trivial(blossom)) {
 		        _delta4->erase(blossom);
+		      }
+		    }
 		    void oddToEven(int blossom, int tree) {
 		      if (!_blossom_set->trivial(blossom)) {
 		        _delta4->erase(blossom);
 		        (*_blossom_data)[blossom].pot -=
 * (2 * _delta_sum - (*_blossom_data)[blossom].offset);
+		      }
 		      for (typename BlossomSet::ItemIt n(*_blossom_set, blossom);
 		           n != INVALID; ++n) {
 		        int ni = (*_node_index)[n];
 		        _blossom_set->increase(n, std::numeric_limits<Value>::max());
 		        (*_node_data)[ni].heap.clear();
 		        (*_node_data)[ni].heap_index.clear();
 		        (*_node_data)[ni].pot +=
 * _delta_sum - (*_blossom_data)[blossom].offset;
 		        _delta1->push(n, (*_node_data)[ni].pot);
 		        for (InArcIt e(_graph, n); e != INVALID; ++e) {
 		          Node v = _graph.source(e);
 		          int vb = _blossom_set->find(v);
 		          int vi = (*_node_index)[v];
 		          Value rw = (*_node_data)[ni].pot + (*_node_data)[vi].pot -
 		            dualScale * _weight[e];
 		          if ((*_blossom_data)[vb].status == EVEN) {
 		            if (_delta3->state(e) != _delta3->IN_HEAP && blossom != vb) {
 		              _delta3->push(e, rw / 2);
+		            }
 		          } else {
 		            typename std::map<int, Arc>::iterator it =
 		              (*_node_data)[vi].heap_index.find(tree);
 		            if (it != (*_node_data)[vi].heap_index.end()) {
 		              if ((*_node_data)[vi].heap[it->second] > rw) {
 		                (*_node_data)[vi].heap.replace(it->second, e);
 		                (*_node_data)[vi].heap.decrease(e, rw);
 		                it->second = e;
+		              }
 		            } else {
 		              (*_node_data)[vi].heap.push(e, rw);
 		              (*_node_data)[vi].heap_index.insert(std::make_pair(tree, e));
+		            }
 		            if ((*_blossom_set)[v] > (*_node_data)[vi].heap.prio()) {
 		              _blossom_set->decrease(v, (*_node_data)[vi].heap.prio());
 		              if ((*_blossom_data)[vb].status == MATCHED) {
 		                if (_delta2->state(vb) != _delta2->IN_HEAP) {
 		                  _delta2->push(vb, _blossom_set->classPrio(vb) -
 		                               (*_blossom_data)[vb].offset);
 		                } else if ((*_delta2)[vb] > _blossom_set->classPrio(vb) -
 		                           (*_blossom_data)[vb].offset) {
 		                  _delta2->decrease(vb, _blossom_set->classPrio(vb) -
 		                                   (*_blossom_data)[vb].offset);
+		                }
+		              }
+		            }
+		          }
+		        }
+		      }
 		      (*_blossom_data)[blossom].offset = 0;
+		    }
 		    void alternatePath(int even, int tree) {
 		      int odd;
 		      evenToMatched(even, tree);
 		      (*_blossom_data)[even].status = MATCHED;
 		      while ((*_blossom_data)[even].pred != INVALID) {
 		        odd = _blossom_set->find(_graph.target((*_blossom_data)[even].pred));
 		        (*_blossom_data)[odd].status = MATCHED;
 		        oddToMatched(odd);
 		        (*_blossom_data)[odd].next = (*_blossom_data)[odd].pred;
 		        even = _blossom_set->find(_graph.target((*_blossom_data)[odd].pred));
 		        (*_blossom_data)[even].status = MATCHED;
 		        evenToMatched(even, tree);
 		        (*_blossom_data)[even].next =
 		          _graph.oppositeArc((*_blossom_data)[odd].pred);
+		      }
+		    }
 		    void destroyTree(int tree) {
 		      for (TreeSet::ItemIt b(*_tree_set, tree); b != INVALID; ++b) {
 		        if ((*_blossom_data)[b].status == EVEN) {
 		          (*_blossom_data)[b].status = MATCHED;
 		          evenToMatched(b, tree);
 		        } else if ((*_blossom_data)[b].status == ODD) {
 		          (*_blossom_data)[b].status = MATCHED;
 		          oddToMatched(b);
+		        }
+		      }
 		      _tree_set->eraseClass(tree);
+		    }
 		    void unmatchNode(const Node& node) {
 		      int blossom = _blossom_set->find(node);
 		      int tree = _tree_set->find(blossom);
 		      alternatePath(blossom, tree);
 		      destroyTree(tree);
 		      (*_blossom_data)[blossom].base = node;
 		      (*_blossom_data)[blossom].next = INVALID;
+		    }
 		    void augmentOnEdge(const Edge& edge) {
 		      int left = _blossom_set->find(_graph.u(edge));
 		      int right = _blossom_set->find(_graph.v(edge));
 		      int left_tree = _tree_set->find(left);
 		      alternatePath(left, left_tree);
 		      destroyTree(left_tree);
 		      int right_tree = _tree_set->find(right);
 		      alternatePath(right, right_tree);
 		      destroyTree(right_tree);
 		      (*_blossom_data)[left].next = _graph.direct(edge, true);
 		      (*_blossom_data)[right].next = _graph.direct(edge, false);
+		    }
 		    void augmentOnArc(const Arc& arc) {
 		      int left = _blossom_set->find(_graph.source(arc));
 		      int right = _blossom_set->find(_graph.target(arc));
 		      (*_blossom_data)[left].status = MATCHED;
 		      int right_tree = _tree_set->find(right);
 		      alternatePath(right, right_tree);
 		      destroyTree(right_tree);
 		      (*_blossom_data)[left].next = arc;
 		      (*_blossom_data)[right].next = _graph.oppositeArc(arc);
+		    }
 		    void extendOnArc(const Arc& arc) {
 		      int base = _blossom_set->find(_graph.target(arc));
 		      int tree = _tree_set->find(base);
 		      int odd = _blossom_set->find(_graph.source(arc));
 		      _tree_set->insert(odd, tree);
 		      (*_blossom_data)[odd].status = ODD;
 		      matchedToOdd(odd);
 		      (*_blossom_data)[odd].pred = arc;
 		      int even = _blossom_set->find(_graph.target((*_blossom_data)[odd].next));
 		      (*_blossom_data)[even].pred = (*_blossom_data)[even].next;
 		      _tree_set->insert(even, tree);
 		      (*_blossom_data)[even].status = EVEN;
 		      matchedToEven(even, tree);
+		    }
 		    void shrinkOnEdge(const Edge& edge, int tree) {
 		      int nca = -1;
 		      std::vector<int> left_path, right_path;
+		      {
 		        std::set<int> left_set, right_set;
 		        int left = _blossom_set->find(_graph.u(edge));
 		        left_path.push_back(left);
 		        left_set.insert(left);
 		        int right = _blossom_set->find(_graph.v(edge));
 		        right_path.push_back(right);
 		        right_set.insert(right);
 		        while (true) {
 		          if ((*_blossom_data)[left].pred == INVALID) break;
 		          left =
 		            _blossom_set->find(_graph.target((*_blossom_data)[left].pred));
 		          left_path.push_back(left);
 		          left =
 		            _blossom_set->find(_graph.target((*_blossom_data)[left].pred));
 		          left_path.push_back(left);
 		          left_set.insert(left);
 		          if (right_set.find(left) != right_set.end()) {
 		            nca = left;
 		            break;
+		          }
 		          if ((*_blossom_data)[right].pred == INVALID) break;
 		          right =
 		            _blossom_set->find(_graph.target((*_blossom_data)[right].pred));
 		          right_path.push_back(right);
 		          right =
 		            _blossom_set->find(_graph.target((*_blossom_data)[right].pred));
 		          right_path.push_back(right);
 		          right_set.insert(right);
 		          if (left_set.find(right) != left_set.end()) {
 		            nca = right;
 		            break;
+		          }
+		        }
 		        if (nca == -1) {
 		          if ((*_blossom_data)[left].pred == INVALID) {
 		            nca = right;
 		            while (left_set.find(nca) == left_set.end()) {
 		              nca =
 		                _blossom_set->find(_graph.target((*_blossom_data)[nca].pred));
 		              right_path.push_back(nca);
 		              nca =
 		                _blossom_set->find(_graph.target((*_blossom_data)[nca].pred));
 		              right_path.push_back(nca);
+		            }
 		          } else {
 		            nca = left;
 		            while (right_set.find(nca) == right_set.end()) {
 		              nca =
 		                _blossom_set->find(_graph.target((*_blossom_data)[nca].pred));
 		              left_path.push_back(nca);
 		              nca =
 		                _blossom_set->find(_graph.target((*_blossom_data)[nca].pred));
 		              left_path.push_back(nca);
+		            }
+		          }
+		        }
+		      }
 		      std::vector<int> subblossoms;
 		      Arc prev;
 		      prev = _graph.direct(edge, true);
 		      for (int i = 0; left_path[i] != nca; i += 2) {
 		        subblossoms.push_back(left_path[i]);
 		        (*_blossom_data)[left_path[i]].next = prev;
 		        _tree_set->erase(left_path[i]);
 		        subblossoms.push_back(left_path[i + 1]);
 		        (*_blossom_data)[left_path[i + 1]].status = EVEN;
 		        oddToEven(left_path[i + 1], tree);
 		        _tree_set->erase(left_path[i + 1]);
 		        prev = _graph.oppositeArc((*_blossom_data)[left_path[i + 1]].pred);
+		      }
 		      int k = 0;
 		      while (right_path[k] != nca) ++k;
 		      subblossoms.push_back(nca);
 		      (*_blossom_data)[nca].next = prev;
 		      for (int i = k - 2; i >= 0; i -= 2) {
 		        subblossoms.push_back(right_path[i + 1]);
 		        (*_blossom_data)[right_path[i + 1]].status = EVEN;
 		        oddToEven(right_path[i + 1], tree);
 		        _tree_set->erase(right_path[i + 1]);
 		        (*_blossom_data)[right_path[i + 1]].next =
 		          (*_blossom_data)[right_path[i + 1]].pred;
 		        subblossoms.push_back(right_path[i]);
 		        _tree_set->erase(right_path[i]);
+		      }
 		      int surface =
 		        _blossom_set->join(subblossoms.begin(), subblossoms.end());
 		      for (int i = 0; i < int(subblossoms.size()); ++i) {
 		        if (!_blossom_set->trivial(subblossoms[i])) {
 		          (*_blossom_data)[subblossoms[i]].pot += 2 * _delta_sum;
+		        }
 		        (*_blossom_data)[subblossoms[i]].status = MATCHED;
+		      }
 		      (*_blossom_data)[surface].pot = -2 * _delta_sum;
 		      (*_blossom_data)[surface].offset = 0;
 		      (*_blossom_data)[surface].status = EVEN;
 		      (*_blossom_data)[surface].pred = (*_blossom_data)[nca].pred;
 		      (*_blossom_data)[surface].next = (*_blossom_data)[nca].pred;
 		      _tree_set->insert(surface, tree);
 		      _tree_set->erase(nca);
+		    }
 		    void splitBlossom(int blossom) {
 		      Arc next = (*_blossom_data)[blossom].next;
 		      Arc pred = (*_blossom_data)[blossom].pred;
 		      int tree = _tree_set->find(blossom);
 		      (*_blossom_data)[blossom].status = MATCHED;
 		      oddToMatched(blossom);
 		      if (_delta2->state(blossom) == _delta2->IN_HEAP) {
 		        _delta2->erase(blossom);
+		      }
 		      std::vector<int> subblossoms;
 		      _blossom_set->split(blossom, std::back_inserter(subblossoms));
 		      Value offset = (*_blossom_data)[blossom].offset;
 		      int b = _blossom_set->find(_graph.source(pred));
 		      int d = _blossom_set->find(_graph.source(next));
 		      int ib = -1, id = -1;
 		      for (int i = 0; i < int(subblossoms.size()); ++i) {
 		        if (subblossoms[i] == b) ib = i;
 		        if (subblossoms[i] == d) id = i;
 		        (*_blossom_data)[subblossoms[i]].offset = offset;
 		        if (!_blossom_set->trivial(subblossoms[i])) {
 		          (*_blossom_data)[subblossoms[i]].pot -= 2 * offset;
+		        }
 		        if (_blossom_set->classPrio(subblossoms[i]) !=
 		            std::numeric_limits<Value>::max()) {
 		          _delta2->push(subblossoms[i],
 		                        _blossom_set->classPrio(subblossoms[i]) -
 		                        (*_blossom_data)[subblossoms[i]].offset);
+		        }
+		      }
 		      if (id > ib ? ((id - ib) % 2 == 0) : ((ib - id) % 2 == 1)) {
 		        for (int i = (id + 1) % subblossoms.size();
 		             i != ib; i = (i + 2) % subblossoms.size()) {
 		          int sb = subblossoms[i];
 		          int tb = subblossoms[(i + 1) % subblossoms.size()];
 		          (*_blossom_data)[sb].next =
 		            _graph.oppositeArc((*_blossom_data)[tb].next);
+		        }
 		        for (int i = ib; i != id; i = (i + 2) % subblossoms.size()) {
 		          int sb = subblossoms[i];
 		          int tb = subblossoms[(i + 1) % subblossoms.size()];
 		          int ub = subblossoms[(i + 2) % subblossoms.size()];
 		          (*_blossom_data)[sb].status = ODD;
 		          matchedToOdd(sb);
 		          _tree_set->insert(sb, tree);
 		          (*_blossom_data)[sb].pred = pred;
 		          (*_blossom_data)[sb].next =
 		            _graph.oppositeArc((*_blossom_data)[tb].next);
 		          pred = (*_blossom_data)[ub].next;
 		          (*_blossom_data)[tb].status = EVEN;
 		          matchedToEven(tb, tree);
 		          _tree_set->insert(tb, tree);
 		          (*_blossom_data)[tb].pred = (*_blossom_data)[tb].next;
+		        }
 		        (*_blossom_data)[subblossoms[id]].status = ODD;
 		        matchedToOdd(subblossoms[id]);
 		        _tree_set->insert(subblossoms[id], tree);
 		        (*_blossom_data)[subblossoms[id]].next = next;
 		        (*_blossom_data)[subblossoms[id]].pred = pred;
 		      } else {
 		        for (int i = (ib + 1) % subblossoms.size();
 		             i != id; i = (i + 2) % subblossoms.size()) {
 		          int sb = subblossoms[i];
 		          int tb = subblossoms[(i + 1) % subblossoms.size()];
 		          (*_blossom_data)[sb].next =
 		            _graph.oppositeArc((*_blossom_data)[tb].next);
+		        }
 		        for (int i = id; i != ib; i = (i + 2) % subblossoms.size()) {
 		          int sb = subblossoms[i];
 		          int tb = subblossoms[(i + 1) % subblossoms.size()];
 		          int ub = subblossoms[(i + 2) % subblossoms.size()];
 		          (*_blossom_data)[sb].status = ODD;
 		          matchedToOdd(sb);
 		          _tree_set->insert(sb, tree);
 		          (*_blossom_data)[sb].next = next;
 		          (*_blossom_data)[sb].pred =
 		            _graph.oppositeArc((*_blossom_data)[tb].next);
 		          (*_blossom_data)[tb].status = EVEN;
 		          matchedToEven(tb, tree);
 		          _tree_set->insert(tb, tree);
 		          (*_blossom_data)[tb].pred =
 		            (*_blossom_data)[tb].next =
 		            _graph.oppositeArc((*_blossom_data)[ub].next);
 		          next = (*_blossom_data)[ub].next;
+		        }
 		        (*_blossom_data)[subblossoms[ib]].status = ODD;
 		        matchedToOdd(subblossoms[ib]);
 		        _tree_set->insert(subblossoms[ib], tree);
 		        (*_blossom_data)[subblossoms[ib]].next = next;
 		        (*_blossom_data)[subblossoms[ib]].pred = pred;
+		      }
 		      _tree_set->erase(blossom);
+		    }
 		    void extractBlossom(int blossom, const Node& base, const Arc& matching) {
 		      if (_blossom_set->trivial(blossom)) {
 		        int bi = (*_node_index)[base];
 		        Value pot = (*_node_data)[bi].pot;
 		        (*_matching)[base] = matching;
 		        _blossom_node_list.push_back(base);
 		        (*_node_potential)[base] = pot;
 		      } else {
 		        Value pot = (*_blossom_data)[blossom].pot;
 		        int bn = _blossom_node_list.size();
 		        std::vector<int> subblossoms;
 		        _blossom_set->split(blossom, std::back_inserter(subblossoms));
 		        int b = _blossom_set->find(base);
 		        int ib = -1;
 		        for (int i = 0; i < int(subblossoms.size()); ++i) {
 		          if (subblossoms[i] == b) { ib = i; break; }
+		        }
 		        for (int i = 1; i < int(subblossoms.size()); i += 2) {
 		          int sb = subblossoms[(ib + i) % subblossoms.size()];
 		          int tb = subblossoms[(ib + i + 1) % subblossoms.size()];
 		          Arc m = (*_blossom_data)[tb].next;
 		          extractBlossom(sb, _graph.target(m), _graph.oppositeArc(m));
 		          extractBlossom(tb, _graph.source(m), m);
+		        }
 		        extractBlossom(subblossoms[ib], base, matching);
 		        int en = _blossom_node_list.size();
 		        _blossom_potential.push_back(BlossomVariable(bn, en, pot));
+		      }
+		    }
 		    void extractMatching() {
 		      std::vector<int> blossoms;
 		      for (typename BlossomSet::ClassIt c(*_blossom_set); c != INVALID; ++c) {
 		        blossoms.push_back(c);
+		      }
 		      for (int i = 0; i < int(blossoms.size()); ++i) {
 		        if ((*_blossom_data)[blossoms[i]].next != INVALID) {
 		          Value offset = (*_blossom_data)[blossoms[i]].offset;
 		          (*_blossom_data)[blossoms[i]].pot += 2 * offset;
 		          for (typename BlossomSet::ItemIt n(*_blossom_set, blossoms[i]);
 		               n != INVALID; ++n) {
 		            (*_node_data)[(*_node_index)[n]].pot -= offset;
+		          }
 		          Arc matching = (*_blossom_data)[blossoms[i]].next;
 		          Node base = _graph.source(matching);
 		          extractBlossom(blossoms[i], base, matching);
 		        } else {
 		          Node base = (*_blossom_data)[blossoms[i]].base;
 		          extractBlossom(blossoms[i], base, INVALID);
+		        }
+		      }
+		    }
 		  public:
 		    /// \brief Constructor
 		    ///
 		    /// Constructor.
 		    MaxWeightedMatching(const Graph& graph, const WeightMap& weight)
 		      : _graph(graph), _weight(weight), _matching(0),
 		        _node_potential(0), _blossom_potential(), _blossom_node_list(),
 		        _node_num(0), _blossom_num(0),
 		        _blossom_index(0), _blossom_set(0), _blossom_data(0),
 		        _node_index(0), _node_heap_index(0), _node_data(0),
 		        _tree_set_index(0), _tree_set(0),
 		        _delta1_index(0), _delta1(0),
 		        _delta2_index(0), _delta2(0),
 		        _delta3_index(0), _delta3(0),
 		        _delta4_index(0), _delta4(0),
 		        _delta_sum(), _unmatched(0),
 		        _fractional(0)
 		    {}
 		    ~MaxWeightedMatching() {
 		      destroyStructures();
 		      if (_fractional) {
 		        delete _fractional;
+		      }
+		    }
 		    /// \name Execution Control
 		    /// The simplest way to execute the algorithm is to use the
 		    /// \ref run() member function.
 		    ///@{
 		    /// \brief Initialize the algorithm
 		    ///
 		    /// This function initializes the algorithm.
 		    void init() {
 		      createStructures();
 		      _blossom_node_list.clear();
 		      _blossom_potential.clear();
 		      for (ArcIt e(_graph); e != INVALID; ++e) {
 		        (*_node_heap_index)[e] = BinHeap<Value, IntArcMap>::PRE_HEAP;
+		      }
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        (*_delta1_index)[n] = _delta1->PRE_HEAP;
+		      }
 		      for (EdgeIt e(_graph); e != INVALID; ++e) {
 		        (*_delta3_index)[e] = _delta3->PRE_HEAP;
+		      }
 		      for (int i = 0; i < _blossom_num; ++i) {
 		        (*_delta2_index)[i] = _delta2->PRE_HEAP;
 		        (*_delta4_index)[i] = _delta4->PRE_HEAP;
+		      }
 		      _unmatched = _node_num;
 		      _delta1->clear();
 		      _delta2->clear();
 		      _delta3->clear();
 		      _delta4->clear();
 		      _blossom_set->clear();
 		      _tree_set->clear();
 		      int index = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        Value max = 0;
 		        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
 		          if (_graph.target(e) == n) continue;
 		          if ((dualScale * _weight[e]) / 2 > max) {
 		            max = (dualScale * _weight[e]) / 2;
+		          }
+		        }
 		        (*_node_index)[n] = index;
 		        (*_node_data)[index].heap_index.clear();
 		        (*_node_data)[index].heap.clear();
 		        (*_node_data)[index].pot = max;
 		        _delta1->push(n, max);
 		        int blossom =
 		          _blossom_set->insert(n, std::numeric_limits<Value>::max());
 		        _tree_set->insert(blossom);
 		        (*_blossom_data)[blossom].status = EVEN;
 		        (*_blossom_data)[blossom].pred = INVALID;
 		        (*_blossom_data)[blossom].next = INVALID;
 		        (*_blossom_data)[blossom].pot = 0;
 		        (*_blossom_data)[blossom].offset = 0;
 		        ++index;
+		      }
 		      for (EdgeIt e(_graph); e != INVALID; ++e) {
 		        int si = (*_node_index)[_graph.u(e)];
 		        int ti = (*_node_index)[_graph.v(e)];
 		        if (_graph.u(e) != _graph.v(e)) {
 		          _delta3->push(e, ((*_node_data)[si].pot + (*_node_data)[ti].pot -
 		                            dualScale * _weight[e]) / 2);
+		        }
+		      }
+		    }
 		    /// \brief Initialize the algorithm with fractional matching
 		    ///
 		    /// This function initializes the algorithm with a fractional
 		    /// matching. This initialization is also called jumpstart heuristic.
 		    void fractionalInit() {
 		      createStructures();
 		      if (_fractional == 0) {
 		        _fractional = new FractionalMatching(_graph, _weight, false);
+		      }
 		      _fractional->run();
 		      for (ArcIt e(_graph); e != INVALID; ++e) {
 		        (*_node_heap_index)[e] = BinHeap<Value, IntArcMap>::PRE_HEAP;
+		      }
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        (*_delta1_index)[n] = _delta1->PRE_HEAP;
+		      }
 		      for (EdgeIt e(_graph); e != INVALID; ++e) {
 		        (*_delta3_index)[e] = _delta3->PRE_HEAP;
+		      }
 		      for (int i = 0; i < _blossom_num; ++i) {
 		        (*_delta2_index)[i] = _delta2->PRE_HEAP;
 		        (*_delta4_index)[i] = _delta4->PRE_HEAP;
+		      }
 		      _unmatched = 0;
 		      int index = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        Value pot = _fractional->nodeValue(n);
 		        (*_node_index)[n] = index;
 		        (*_node_data)[index].pot = pot;
 		        int blossom =
 		          _blossom_set->insert(n, std::numeric_limits<Value>::max());
 		        (*_blossom_data)[blossom].status = MATCHED;
 		        (*_blossom_data)[blossom].pred = INVALID;
 		        (*_blossom_data)[blossom].next = _fractional->matching(n);
 		        if (_fractional->matching(n) == INVALID) {
 		          (*_blossom_data)[blossom].base = n;
+		        }
 		        (*_blossom_data)[blossom].pot = 0;
 		        (*_blossom_data)[blossom].offset = 0;
 		        ++index;
+		      }
 		      typename Graph::template NodeMap<bool> processed(_graph, false);
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        if (processed[n]) continue;
 		        processed[n] = true;
 		        if (_fractional->matching(n) == INVALID) continue;
 		        int num = 1;
 		        Node v = _graph.target(_fractional->matching(n));
 		        while (n != v) {
 		          processed[v] = true;
 		          v = _graph.target(_fractional->matching(v));
 		          ++num;
+		        }
 		        if (num % 2 == 1) {
 		          std::vector<int> subblossoms(num);
 		          subblossoms[--num] = _blossom_set->find(n);
 		          _delta1->push(n, _fractional->nodeValue(n));
 		          v = _graph.target(_fractional->matching(n));
 		          while (n != v) {
 		            subblossoms[--num] = _blossom_set->find(v);
 		            _delta1->push(v, _fractional->nodeValue(v));
 		            v = _graph.target(_fractional->matching(v));
+		          }
 		          int surface =
 		            _blossom_set->join(subblossoms.begin(), subblossoms.end());
 		          (*_blossom_data)[surface].status = EVEN;
 		          (*_blossom_data)[surface].pred = INVALID;
 		          (*_blossom_data)[surface].next = INVALID;
 		          (*_blossom_data)[surface].pot = 0;
 		          (*_blossom_data)[surface].offset = 0;
 		          _tree_set->insert(surface);
 		          ++_unmatched;
+		        }
+		      }
 		      for (EdgeIt e(_graph); e != INVALID; ++e) {
 		        int si = (*_node_index)[_graph.u(e)];
 		        int sb = _blossom_set->find(_graph.u(e));
 		        int ti = (*_node_index)[_graph.v(e)];
 		        int tb = _blossom_set->find(_graph.v(e));
 		        if ((*_blossom_data)[sb].status == EVEN &&
 		            (*_blossom_data)[tb].status == EVEN && sb != tb) {
 		          _delta3->push(e, ((*_node_data)[si].pot + (*_node_data)[ti].pot -
 		                            dualScale * _weight[e]) / 2);
+		        }
+		      }
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        int nb = _blossom_set->find(n);
 		        if ((*_blossom_data)[nb].status != MATCHED) continue;
 		        int ni = (*_node_index)[n];
 		        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
 		          Node v = _graph.target(e);
 		          int vb = _blossom_set->find(v);
 		          int vi = (*_node_index)[v];
 		          Value rw = (*_node_data)[ni].pot + (*_node_data)[vi].pot -
 		            dualScale * _weight[e];
 		          if ((*_blossom_data)[vb].status == EVEN) {
 		            int vt = _tree_set->find(vb);
 		            typename std::map<int, Arc>::iterator it =
 		              (*_node_data)[ni].heap_index.find(vt);
 		            if (it != (*_node_data)[ni].heap_index.end()) {
 		              if ((*_node_data)[ni].heap[it->second] > rw) {
 		                (*_node_data)[ni].heap.replace(it->second, e);
 		                (*_node_data)[ni].heap.decrease(e, rw);
 		                it->second = e;
+		              }
 		            } else {
 		              (*_node_data)[ni].heap.push(e, rw);
 		              (*_node_data)[ni].heap_index.insert(std::make_pair(vt, e));
+		            }
+		          }
+		        }
 		        if (!(*_node_data)[ni].heap.empty()) {
 		          _blossom_set->decrease(n, (*_node_data)[ni].heap.prio());
 		          _delta2->push(nb, _blossom_set->classPrio(nb));
+		        }
+		      }
+		    }
 		    /// \brief Start the algorithm
 		    ///
 		    /// This function starts the algorithm.
 		    ///
 		    /// \pre \ref init() or \ref fractionalInit() must be called
 		    /// before using this function.
 		    void start() {
 		      enum OpType {
 		        D1, D2, D3, D4
 		      };
 		      while (_unmatched > 0) {
 		        Value d1 = !_delta1->empty() ?
 		          _delta1->prio() : std::numeric_limits<Value>::max();
 		        Value d2 = !_delta2->empty() ?
 		          _delta2->prio() : std::numeric_limits<Value>::max();
 		        Value d3 = !_delta3->empty() ?
 		          _delta3->prio() : std::numeric_limits<Value>::max();
 		        Value d4 = !_delta4->empty() ?
 		          _delta4->prio() : std::numeric_limits<Value>::max();
 		        _delta_sum = d3; OpType ot = D3;
 		        if (d1 < _delta_sum) { _delta_sum = d1; ot = D1; }
 		        if (d2 < _delta_sum) { _delta_sum = d2; ot = D2; }
 		        if (d4 < _delta_sum) { _delta_sum = d4; ot = D4; }
 		        switch (ot) {
 		        case D1:
+		          {
 		            Node n = _delta1->top();
 		            unmatchNode(n);
 		            --_unmatched;
+		          }
 		          break;
 		        case D2:
+		          {
 		            int blossom = _delta2->top();
 		            Node n = _blossom_set->classTop(blossom);
 		            Arc a = (*_node_data)[(*_node_index)[n]].heap.top();
 		            if ((*_blossom_data)[blossom].next == INVALID) {
 		              augmentOnArc(a);
 		              --_unmatched;
 		            } else {
 		              extendOnArc(a);
+		            }
+		          }
 		          break;
 		        case D3:
+		          {
 		            Edge e = _delta3->top();
 		            int left_blossom = _blossom_set->find(_graph.u(e));
 		            int right_blossom = _blossom_set->find(_graph.v(e));
 		            if (left_blossom == right_blossom) {
 		              _delta3->pop();
 		            } else {
 		              int left_tree = _tree_set->find(left_blossom);
 		              int right_tree = _tree_set->find(right_blossom);
 		              if (left_tree == right_tree) {
 		                shrinkOnEdge(e, left_tree);
 		              } else {
 		                augmentOnEdge(e);
 		                _unmatched -= 2;
+		              }
+		            }
 		          } break;
 		        case D4:
 		          splitBlossom(_delta4->top());
 		          break;
+		        }
+		      }
 		      extractMatching();
+		    }
 		    /// \brief Run the algorithm.
 		    ///
 		    /// This method runs the \c %MaxWeightedMatching algorithm.
 		    ///
 		    /// \note mwm.run() is just a shortcut of the following code.
 		    /// \code
 		    ///   mwm.fractionalInit();
 		    ///   mwm.start();
 		    /// \endcode
 		    void run() {
 		      fractionalInit();
 		      start();
+		    }
 		    /// @}
 		    /// \name Primal Solution
 		    /// Functions to get the primal solution, i.e. the maximum weighted
 		    /// matching.\n
 		    /// Either \ref run() or \ref start() function should be called before
 		    /// using them.
 		    /// @{
 		    /// \brief Return the weight of the matching.
 		    ///
 		    /// This function returns the weight of the found matching.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Value matchingWeight() const {
 		      Value sum = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        if ((*_matching)[n] != INVALID) {
 		          sum += _weight[(*_matching)[n]];
+		        }
+		      }
 		      return sum / 2;
+		    }
 		    /// \brief Return the size (cardinality) of the matching.
 		    ///
 		    /// This function returns the size (cardinality) of the found matching.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    int matchingSize() const {
 		      int num = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        if ((*_matching)[n] != INVALID) {
 		          ++num;
+		        }
+		      }
 		      return num /= 2;
+		    }
 		    /// \brief Return \c true if the given edge is in the matching.
 		    ///
 		    /// This function returns \c true if the given edge is in the found
 		    /// matching.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    bool matching(const Edge& edge) const {
 		      return edge == (*_matching)[_graph.u(edge)];
+		    }
 		    /// \brief Return the matching arc (or edge) incident to the given node.
 		    ///
 		    /// This function returns the matching arc (or edge) incident to the
 		    /// given node in the found matching or \c INVALID if the node is
 		    /// not covered by the matching.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Arc matching(const Node& node) const {
 		      return (*_matching)[node];
+		    }
 		    /// \brief Return a const reference to the matching map.
 		    ///
 		    /// This function returns a const reference to a node map that stores
 		    /// the matching arc (or edge) incident to each node.
 		    const MatchingMap& matchingMap() const {
 		      return *_matching;
+		    }
 		    /// \brief Return the mate of the given node.
 		    ///
 		    /// This function returns the mate of the given node in the found
 		    /// matching or \c INVALID if the node is not covered by the matching.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Node mate(const Node& node) const {
 		      return (*_matching)[node] != INVALID ?
 		        _graph.target((*_matching)[node]) : INVALID;
+		    }
 		    /// @}
 		    /// \name Dual Solution
 		    /// Functions to get the dual solution.\n
 		    /// Either \ref run() or \ref start() function should be called before
 		    /// using them.
 		    /// @{
 		    /// \brief Return the value of the dual solution.
 		    ///
 		    /// This function returns the value of the dual solution.
 		    /// It should be equal to the primal value scaled by \ref dualScale
 		    /// "dual scale".
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Value dualValue() const {
 		      Value sum = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        sum += nodeValue(n);
+		      }
 		      for (int i = 0; i < blossomNum(); ++i) {
 		        sum += blossomValue(i) * (blossomSize(i) / 2);
+		      }
 		      return sum;
+		    }
 		    /// \brief Return the dual value (potential) of the given node.
 		    ///
 		    /// This function returns the dual value (potential) of the given node.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Value nodeValue(const Node& n) const {
 		      return (*_node_potential)[n];
+		    }
 		    /// \brief Return the number of the blossoms in the basis.
 		    ///
 		    /// This function returns the number of the blossoms in the basis.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    /// \see BlossomIt
 		    int blossomNum() const {
 		      return _blossom_potential.size();
+		    }
 		    /// \brief Return the number of the nodes in the given blossom.
 		    ///
 		    /// This function returns the number of the nodes in the given blossom.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    /// \see BlossomIt
 		    int blossomSize(int k) const {
 		      return _blossom_potential[k].end - _blossom_potential[k].begin;
+		    }
 		    /// \brief Return the dual value (ptential) of the given blossom.
 		    ///
 		    /// This function returns the dual value (ptential) of the given blossom.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Value blossomValue(int k) const {
 		      return _blossom_potential[k].value;
+		    }
 		    /// \brief Iterator for obtaining the nodes of a blossom.
 		    ///
 		    /// This class provides an iterator for obtaining the nodes of the
 		    /// given blossom. It lists a subset of the nodes.
 		    /// Before using this iterator, you must allocate a
 		    /// MaxWeightedMatching class and execute it.
 		    class BlossomIt {
 		    public:
 		      /// \brief Constructor.
 		      ///
 		      /// Constructor to get the nodes of the given variable.
 		      ///
 		      /// \pre Either \ref MaxWeightedMatching::run() "algorithm.run()" or
 		      /// \ref MaxWeightedMatching::start() "algorithm.start()" must be
 		      /// called before initializing this iterator.
 		      BlossomIt(const MaxWeightedMatching& algorithm, int variable)
 		        : _algorithm(&algorithm)
+		      {
 		        _index = _algorithm->_blossom_potential[variable].begin;
 		        _last = _algorithm->_blossom_potential[variable].end;
+		      }
 		      /// \brief Conversion to \c Node.
 		      ///
 		      /// Conversion to \c Node.
 		      operator Node() const {
 		        return _algorithm->_blossom_node_list[_index];
+		      }
 		      /// \brief Increment operator.
 		      ///
 		      /// Increment operator.
 		      BlossomIt& operator++() {
 		        ++_index;
 		        return *this;
+		      }
 		      /// \brief Validity checking
 		      ///
 		      /// Checks whether the iterator is invalid.
 		      bool operator==(Invalid) const { return _index == _last; }
 		      /// \brief Validity checking
 		      ///
 		      /// Checks whether the iterator is valid.
 		      bool operator!=(Invalid) const { return _index != _last; }
 		    private:
 		      const MaxWeightedMatching* _algorithm;
 		      int _last;
 		      int _index;
 		    };
 		    /// @}
 		  };
 		  /// \ingroup matching
 		  ///
 		  /// \brief Weighted perfect matching in general graphs
 		  ///
 		  /// This class provides an efficient implementation of Edmond's
 		  /// maximum weighted perfect matching algorithm. The implementation
 		  /// is based on extensive use of priority queues and provides
 		  /// \f$O(nm\log n)\f$ time complexity.
 		  ///
 		  /// The maximum weighted perfect matching problem is to find a subset of
 		  /// the edges in an undirected graph with maximum overall weight for which
 		  /// each node has exactly one incident edge.
 		  /// It can be formulated with the following linear program.
 		  /// \f[ \sum_{e \in \delta(u)}x_e = 1 \quad \forall u\in V\f]
 		  /** \f[ \sum_{e \in \gamma(B)}x_e \le \frac{\vert B \vert - 1}{2}
 		      \quad \forall B\in\mathcal{O}\f] */
 		  /// \f[x_e \ge 0\quad \forall e\in E\f]
 		  /// \f[\max \sum_{e\in E}x_ew_e\f]
 		  /// where \f$\delta(X)\f$ is the set of edges incident to a node in
 		  /// \f$X\f$, \f$\gamma(X)\f$ is the set of edges with both ends in
 		  /// \f$X\f$ and \f$\mathcal{O}\f$ is the set of odd cardinality
 		  /// subsets of the nodes.
 		  ///
 		  /// The algorithm calculates an optimal matching and a proof of the
 		  /// optimality. The solution of the dual problem can be used to check
 		  /// the result of the algorithm. The dual linear problem is the
 		  /// following.
 		  /** \f[ y_u + y_v + \sum_{B \in \mathcal{O}, uv \in \gamma(B)}z_B \ge
 		      w_{uv} \quad \forall uv\in E\f] */
 		  /// \f[z_B \ge 0 \quad \forall B \in \mathcal{O}\f]
 		  /** \f[\min \sum_{u \in V}y_u + \sum_{B \in \mathcal{O}}
 		      \frac{\vert B \vert - 1}{2}z_B\f] */
 		  ///
 		  /// The algorithm can be executed with the run() function.
 		  /// After it the matching (the primal solution) and the dual solution
 		  /// can be obtained using the query functions and the
 		  /// \ref MaxWeightedPerfectMatching::BlossomIt "BlossomIt" nested class,
 		  /// which is able to iterate on the nodes of a blossom.
 		  /// If the value type is integer, then the dual solution is multiplied
 		  /// by \ref MaxWeightedMatching::dualScale "4".
 		  ///
 		  /// \tparam GR The undirected graph type the algorithm runs on.
 		  /// \tparam WM The type edge weight map. The default type is
 		  /// \ref concepts::Graph::EdgeMap "GR::EdgeMap<int>".
 		#ifdef DOXYGEN
 		  template <typename GR, typename WM>
 		#else
 		  template <typename GR,
 		            typename WM = typename GR::template EdgeMap<int> >
 		#endif
 		  class MaxWeightedPerfectMatching {
 		  public:
 		    /// The graph type of the algorithm
 		    typedef GR Graph;
 		    /// The type of the edge weight map
 		    typedef WM WeightMap;
 		    /// The value type of the edge weights
 		    typedef typename WeightMap::Value Value;
 		    /// \brief Scaling factor for dual solution
 		    ///
 		    /// Scaling factor for dual solution, it is equal to 4 or 1
 		    /// according to the value type.
 		    static const int dualScale =
 		      std::numeric_limits<Value>::is_integer ? 4 : 1;
 		    /// The type of the matching map
 		    typedef typename Graph::template NodeMap<typename Graph::Arc>
 		    MatchingMap;
 		  private:
 		    TEMPLATE_GRAPH_TYPEDEFS(Graph);
 		    typedef typename Graph::template NodeMap<Value> NodePotential;
 		    typedef std::vector<Node> BlossomNodeList;
 		    struct BlossomVariable {
 		      int begin, end;
 		      Value value;
 		      BlossomVariable(int _begin, int _end, Value _value)
 		        : begin(_begin), end(_end), value(_value) {}
 		    };
 		    typedef std::vector<BlossomVariable> BlossomPotential;
 		    const Graph& _graph;
 		    const WeightMap& _weight;
 		    MatchingMap* _matching;
 		    NodePotential* _node_potential;
 		    BlossomPotential _blossom_potential;
 		    BlossomNodeList _blossom_node_list;
 		    int _node_num;
 		    int _blossom_num;
 		    typedef RangeMap<int> IntIntMap;
 		    enum Status {
 		      EVEN = -1, MATCHED = 0, ODD = 1
 		    };
 		    typedef HeapUnionFind<Value, IntNodeMap> BlossomSet;
 		    struct BlossomData {
 		      int tree;
 		      Status status;
 		      Arc pred, next;
 		      Value pot, offset;
 		    };
 		    IntNodeMap *_blossom_index;
 		    BlossomSet *_blossom_set;
 		    RangeMap<BlossomData>* _blossom_data;
 		    IntNodeMap *_node_index;
 		    IntArcMap *_node_heap_index;
 		    struct NodeData {
 		      NodeData(IntArcMap& node_heap_index)
 		        : heap(node_heap_index) {}
 		      int blossom;
 		      Value pot;
 		      BinHeap<Value, IntArcMap> heap;
 		      std::map<int, Arc> heap_index;
 		      int tree;
 		    };
 		    RangeMap<NodeData>* _node_data;
 		    typedef ExtendFindEnum<IntIntMap> TreeSet;
 		    IntIntMap *_tree_set_index;
 		    TreeSet *_tree_set;
 		    IntIntMap *_delta2_index;
 		    BinHeap<Value, IntIntMap> *_delta2;
 		    IntEdgeMap *_delta3_index;
 		    BinHeap<Value, IntEdgeMap> *_delta3;
 		    IntIntMap *_delta4_index;
 		    BinHeap<Value, IntIntMap> *_delta4;
 		    Value _delta_sum;
 		    int _unmatched;
 		    typedef MaxWeightedPerfectFractionalMatching<Graph, WeightMap>
 		    FractionalMatching;
 		    FractionalMatching *_fractional;
 		    void createStructures() {
 		      _node_num = countNodes(_graph);
 		      _blossom_num = _node_num * 3 / 2;
 		      if (!_matching) {
 		        _matching = new MatchingMap(_graph);
+		      }
 		      if (!_node_potential) {
 		        _node_potential = new NodePotential(_graph);
+		      }
 		      if (!_blossom_set) {
 		        _blossom_index = new IntNodeMap(_graph);
 		        _blossom_set = new BlossomSet(*_blossom_index);
 		        _blossom_data = new RangeMap<BlossomData>(_blossom_num);
 		      } else if (_blossom_data->size() != _blossom_num) {
 		        delete _blossom_data;
 		        _blossom_data = new RangeMap<BlossomData>(_blossom_num);
+		      }
 		      if (!_node_index) {
 		        _node_index = new IntNodeMap(_graph);
 		        _node_heap_index = new IntArcMap(_graph);
 		        _node_data = new RangeMap<NodeData>(_node_num,
 		                                            NodeData(*_node_heap_index));
 		      } else if (_node_data->size() != _node_num) {
 		        delete _node_data;
 		        _node_data = new RangeMap<NodeData>(_node_num,
 		                                            NodeData(*_node_heap_index));
+		      }
 		      if (!_tree_set) {
 		        _tree_set_index = new IntIntMap(_blossom_num);
 		        _tree_set = new TreeSet(*_tree_set_index);
 		      } else {
 		        _tree_set_index->resize(_blossom_num);
+		      }
 		      if (!_delta2) {
 		        _delta2_index = new IntIntMap(_blossom_num);
 		        _delta2 = new BinHeap<Value, IntIntMap>(*_delta2_index);
 		      } else {
 		        _delta2_index->resize(_blossom_num);
+		      }
 		      if (!_delta3) {
 		        _delta3_index = new IntEdgeMap(_graph);
 		        _delta3 = new BinHeap<Value, IntEdgeMap>(*_delta3_index);
+		      }
 		      if (!_delta4) {
 		        _delta4_index = new IntIntMap(_blossom_num);
 		        _delta4 = new BinHeap<Value, IntIntMap>(*_delta4_index);
 		      } else {
 		        _delta4_index->resize(_blossom_num);
+		      }
+		    }
 		    void destroyStructures() {
 		      if (_matching) {
 		        delete _matching;
+		      }
 		      if (_node_potential) {
 		        delete _node_potential;
+		      }
 		      if (_blossom_set) {
 		        delete _blossom_index;
 		        delete _blossom_set;
 		        delete _blossom_data;
+		      }
 		      if (_node_index) {
 		        delete _node_index;
 		        delete _node_heap_index;
 		        delete _node_data;
+		      }
 		      if (_tree_set) {
 		        delete _tree_set_index;
 		        delete _tree_set;
+		      }
 		      if (_delta2) {
 		        delete _delta2_index;
 		        delete _delta2;
+		      }
 		      if (_delta3) {
 		        delete _delta3_index;
 		        delete _delta3;
+		      }
 		      if (_delta4) {
 		        delete _delta4_index;
 		        delete _delta4;
+		      }
+		    }
 		    void matchedToEven(int blossom, int tree) {
 		      if (_delta2->state(blossom) == _delta2->IN_HEAP) {
 		        _delta2->erase(blossom);
+		      }
 		      if (!_blossom_set->trivial(blossom)) {
 		        (*_blossom_data)[blossom].pot -=
 * (_delta_sum - (*_blossom_data)[blossom].offset);
+		      }
 		      for (typename BlossomSet::ItemIt n(*_blossom_set, blossom);
 		           n != INVALID; ++n) {
 		        _blossom_set->increase(n, std::numeric_limits<Value>::max());
 		        int ni = (*_node_index)[n];
 		        (*_node_data)[ni].heap.clear();
 		        (*_node_data)[ni].heap_index.clear();
 		        (*_node_data)[ni].pot += _delta_sum - (*_blossom_data)[blossom].offset;
 		        for (InArcIt e(_graph, n); e != INVALID; ++e) {
 		          Node v = _graph.source(e);
 		          int vb = _blossom_set->find(v);
 		          int vi = (*_node_index)[v];
 		          Value rw = (*_node_data)[ni].pot + (*_node_data)[vi].pot -
 		            dualScale * _weight[e];
 		          if ((*_blossom_data)[vb].status == EVEN) {
 		            if (_delta3->state(e) != _delta3->IN_HEAP && blossom != vb) {
 		              _delta3->push(e, rw / 2);
+		            }
 		          } else {
 		            typename std::map<int, Arc>::iterator it =
 		              (*_node_data)[vi].heap_index.find(tree);
 		            if (it != (*_node_data)[vi].heap_index.end()) {
 		              if ((*_node_data)[vi].heap[it->second] > rw) {
 		                (*_node_data)[vi].heap.replace(it->second, e);
 		                (*_node_data)[vi].heap.decrease(e, rw);
 		                it->second = e;
+		              }
 		            } else {
 		              (*_node_data)[vi].heap.push(e, rw);
 		              (*_node_data)[vi].heap_index.insert(std::make_pair(tree, e));
+		            }
 		            if ((*_blossom_set)[v] > (*_node_data)[vi].heap.prio()) {
 		              _blossom_set->decrease(v, (*_node_data)[vi].heap.prio());
 		              if ((*_blossom_data)[vb].status == MATCHED) {
 		                if (_delta2->state(vb) != _delta2->IN_HEAP) {
 		                  _delta2->push(vb, _blossom_set->classPrio(vb) -
 		                               (*_blossom_data)[vb].offset);
 		                } else if ((*_delta2)[vb] > _blossom_set->classPrio(vb) -
 		                           (*_blossom_data)[vb].offset){
 		                  _delta2->decrease(vb, _blossom_set->classPrio(vb) -
 		                                   (*_blossom_data)[vb].offset);
+		                }
+		              }
+		            }
+		          }
+		        }
+		      }
 		      (*_blossom_data)[blossom].offset = 0;
+		    }
 		    void matchedToOdd(int blossom) {
 		      if (_delta2->state(blossom) == _delta2->IN_HEAP) {
 		        _delta2->erase(blossom);
+		      }
 		      (*_blossom_data)[blossom].offset += _delta_sum;
 		      if (!_blossom_set->trivial(blossom)) {
 		        _delta4->push(blossom, (*_blossom_data)[blossom].pot / 2 +
 		                     (*_blossom_data)[blossom].offset);
+		      }
+		    }
 		    void evenToMatched(int blossom, int tree) {
 		      if (!_blossom_set->trivial(blossom)) {
 		        (*_blossom_data)[blossom].pot += 2 * _delta_sum;
+		      }
 		      for (typename BlossomSet::ItemIt n(*_blossom_set, blossom);
 		           n != INVALID; ++n) {
 		        int ni = (*_node_index)[n];
 		        (*_node_data)[ni].pot -= _delta_sum;
 		        for (InArcIt e(_graph, n); e != INVALID; ++e) {
 		          Node v = _graph.source(e);
 		          int vb = _blossom_set->find(v);
 		          int vi = (*_node_index)[v];
 		          Value rw = (*_node_data)[ni].pot + (*_node_data)[vi].pot -
 		            dualScale * _weight[e];
 		          if (vb == blossom) {
 		            if (_delta3->state(e) == _delta3->IN_HEAP) {
 		              _delta3->erase(e);
+		            }
 		          } else if ((*_blossom_data)[vb].status == EVEN) {
 		            if (_delta3->state(e) == _delta3->IN_HEAP) {
 		              _delta3->erase(e);
+		            }
 		            int vt = _tree_set->find(vb);
 		            if (vt != tree) {
 		              Arc r = _graph.oppositeArc(e);
 		              typename std::map<int, Arc>::iterator it =
 		                (*_node_data)[ni].heap_index.find(vt);
 		              if (it != (*_node_data)[ni].heap_index.end()) {
 		                if ((*_node_data)[ni].heap[it->second] > rw) {
 		                  (*_node_data)[ni].heap.replace(it->second, r);
 		                  (*_node_data)[ni].heap.decrease(r, rw);
 		                  it->second = r;
+		                }
 		              } else {
 		                (*_node_data)[ni].heap.push(r, rw);
 		                (*_node_data)[ni].heap_index.insert(std::make_pair(vt, r));
+		              }
 		              if ((*_blossom_set)[n] > (*_node_data)[ni].heap.prio()) {
 		                _blossom_set->decrease(n, (*_node_data)[ni].heap.prio());
 		                if (_delta2->state(blossom) != _delta2->IN_HEAP) {
 		                  _delta2->push(blossom, _blossom_set->classPrio(blossom) -
 		                               (*_blossom_data)[blossom].offset);
 		                } else if ((*_delta2)[blossom] >
 		                           _blossom_set->classPrio(blossom) -
 		                           (*_blossom_data)[blossom].offset){
 		                  _delta2->decrease(blossom, _blossom_set->classPrio(blossom) -
 		                                   (*_blossom_data)[blossom].offset);
+		                }
+		              }
+		            }
 		          } else {
 		            typename std::map<int, Arc>::iterator it =
 		              (*_node_data)[vi].heap_index.find(tree);
 		            if (it != (*_node_data)[vi].heap_index.end()) {
 		              (*_node_data)[vi].heap.erase(it->second);
 		              (*_node_data)[vi].heap_index.erase(it);
 		              if ((*_node_data)[vi].heap.empty()) {
 		                _blossom_set->increase(v, std::numeric_limits<Value>::max());
 		              } else if ((*_blossom_set)[v] < (*_node_data)[vi].heap.prio()) {
 		                _blossom_set->increase(v, (*_node_data)[vi].heap.prio());
+		              }
 		              if ((*_blossom_data)[vb].status == MATCHED) {
 		                if (_blossom_set->classPrio(vb) ==
 		                    std::numeric_limits<Value>::max()) {
 		                  _delta2->erase(vb);
 		                } else if ((*_delta2)[vb] < _blossom_set->classPrio(vb) -
 		                           (*_blossom_data)[vb].offset) {
 		                  _delta2->increase(vb, _blossom_set->classPrio(vb) -
 		                                   (*_blossom_data)[vb].offset);
+		                }
+		              }
+		            }
+		          }
+		        }
+		      }
+		    }
 		    void oddToMatched(int blossom) {
 		      (*_blossom_data)[blossom].offset -= _delta_sum;
 		      if (_blossom_set->classPrio(blossom) !=
 		          std::numeric_limits<Value>::max()) {
 		        _delta2->push(blossom, _blossom_set->classPrio(blossom) -
 		                       (*_blossom_data)[blossom].offset);
+		      }
 		      if (!_blossom_set->trivial(blossom)) {
 		        _delta4->erase(blossom);
+		      }
+		    }
 		    void oddToEven(int blossom, int tree) {
 		      if (!_blossom_set->trivial(blossom)) {
 		        _delta4->erase(blossom);
 		        (*_blossom_data)[blossom].pot -=
 * (2 * _delta_sum - (*_blossom_data)[blossom].offset);
+		      }
 		      for (typename BlossomSet::ItemIt n(*_blossom_set, blossom);
 		           n != INVALID; ++n) {
 		        int ni = (*_node_index)[n];
 		        _blossom_set->increase(n, std::numeric_limits<Value>::max());
 		        (*_node_data)[ni].heap.clear();
 		        (*_node_data)[ni].heap_index.clear();
 		        (*_node_data)[ni].pot +=
 * _delta_sum - (*_blossom_data)[blossom].offset;
 		        for (InArcIt e(_graph, n); e != INVALID; ++e) {
 		          Node v = _graph.source(e);
 		          int vb = _blossom_set->find(v);
 		          int vi = (*_node_index)[v];
 		          Value rw = (*_node_data)[ni].pot + (*_node_data)[vi].pot -
 		            dualScale * _weight[e];
 		          if ((*_blossom_data)[vb].status == EVEN) {
 		            if (_delta3->state(e) != _delta3->IN_HEAP && blossom != vb) {
 		              _delta3->push(e, rw / 2);
+		            }
 		          } else {
 		            typename std::map<int, Arc>::iterator it =
 		              (*_node_data)[vi].heap_index.find(tree);
 		            if (it != (*_node_data)[vi].heap_index.end()) {
 		              if ((*_node_data)[vi].heap[it->second] > rw) {
 		                (*_node_data)[vi].heap.replace(it->second, e);
 		                (*_node_data)[vi].heap.decrease(e, rw);
 		                it->second = e;
+		              }
 		            } else {
 		              (*_node_data)[vi].heap.push(e, rw);
 		              (*_node_data)[vi].heap_index.insert(std::make_pair(tree, e));
+		            }
 		            if ((*_blossom_set)[v] > (*_node_data)[vi].heap.prio()) {
 		              _blossom_set->decrease(v, (*_node_data)[vi].heap.prio());
 		              if ((*_blossom_data)[vb].status == MATCHED) {
 		                if (_delta2->state(vb) != _delta2->IN_HEAP) {
 		                  _delta2->push(vb, _blossom_set->classPrio(vb) -
 		                               (*_blossom_data)[vb].offset);
 		                } else if ((*_delta2)[vb] > _blossom_set->classPrio(vb) -
 		                           (*_blossom_data)[vb].offset) {
 		                  _delta2->decrease(vb, _blossom_set->classPrio(vb) -
 		                                   (*_blossom_data)[vb].offset);
+		                }
+		              }
+		            }
+		          }
+		        }
+		      }
 		      (*_blossom_data)[blossom].offset = 0;
+		    }
 		    void alternatePath(int even, int tree) {
 		      int odd;
 		      evenToMatched(even, tree);
 		      (*_blossom_data)[even].status = MATCHED;
 		      while ((*_blossom_data)[even].pred != INVALID) {
 		        odd = _blossom_set->find(_graph.target((*_blossom_data)[even].pred));
 		        (*_blossom_data)[odd].status = MATCHED;
 		        oddToMatched(odd);
 		        (*_blossom_data)[odd].next = (*_blossom_data)[odd].pred;
 		        even = _blossom_set->find(_graph.target((*_blossom_data)[odd].pred));
 		        (*_blossom_data)[even].status = MATCHED;
 		        evenToMatched(even, tree);
 		        (*_blossom_data)[even].next =
 		          _graph.oppositeArc((*_blossom_data)[odd].pred);
+		      }
+		    }
 		    void destroyTree(int tree) {
 		      for (TreeSet::ItemIt b(*_tree_set, tree); b != INVALID; ++b) {
 		        if ((*_blossom_data)[b].status == EVEN) {
 		          (*_blossom_data)[b].status = MATCHED;
 		          evenToMatched(b, tree);
 		        } else if ((*_blossom_data)[b].status == ODD) {
 		          (*_blossom_data)[b].status = MATCHED;
 		          oddToMatched(b);
+		        }
+		      }
 		      _tree_set->eraseClass(tree);
+		    }
 		    void augmentOnEdge(const Edge& edge) {
 		      int left = _blossom_set->find(_graph.u(edge));
 		      int right = _blossom_set->find(_graph.v(edge));
 		      int left_tree = _tree_set->find(left);
 		      alternatePath(left, left_tree);
 		      destroyTree(left_tree);
 		      int right_tree = _tree_set->find(right);
 		      alternatePath(right, right_tree);
 		      destroyTree(right_tree);
 		      (*_blossom_data)[left].next = _graph.direct(edge, true);
 		      (*_blossom_data)[right].next = _graph.direct(edge, false);
+		    }
 		    void extendOnArc(const Arc& arc) {
 		      int base = _blossom_set->find(_graph.target(arc));
 		      int tree = _tree_set->find(base);
 		      int odd = _blossom_set->find(_graph.source(arc));
 		      _tree_set->insert(odd, tree);
 		      (*_blossom_data)[odd].status = ODD;
 		      matchedToOdd(odd);
 		      (*_blossom_data)[odd].pred = arc;
 		      int even = _blossom_set->find(_graph.target((*_blossom_data)[odd].next));
 		      (*_blossom_data)[even].pred = (*_blossom_data)[even].next;
 		      _tree_set->insert(even, tree);
 		      (*_blossom_data)[even].status = EVEN;
 		      matchedToEven(even, tree);
+		    }
 		    void shrinkOnEdge(const Edge& edge, int tree) {
 		      int nca = -1;
 		      std::vector<int> left_path, right_path;
+		      {
 		        std::set<int> left_set, right_set;
 		        int left = _blossom_set->find(_graph.u(edge));
 		        left_path.push_back(left);
 		        left_set.insert(left);
 		        int right = _blossom_set->find(_graph.v(edge));
 		        right_path.push_back(right);
 		        right_set.insert(right);
 		        while (true) {
 		          if ((*_blossom_data)[left].pred == INVALID) break;
 		          left =
 		            _blossom_set->find(_graph.target((*_blossom_data)[left].pred));
 		          left_path.push_back(left);
 		          left =
 		            _blossom_set->find(_graph.target((*_blossom_data)[left].pred));
 		          left_path.push_back(left);
 		          left_set.insert(left);
 		          if (right_set.find(left) != right_set.end()) {
 		            nca = left;
 		            break;
+		          }
 		          if ((*_blossom_data)[right].pred == INVALID) break;
 		          right =
 		            _blossom_set->find(_graph.target((*_blossom_data)[right].pred));
 		          right_path.push_back(right);
 		          right =
 		            _blossom_set->find(_graph.target((*_blossom_data)[right].pred));
 		          right_path.push_back(right);
 		          right_set.insert(right);
 		          if (left_set.find(right) != left_set.end()) {
 		            nca = right;
 		            break;
+		          }
+		        }
 		        if (nca == -1) {
 		          if ((*_blossom_data)[left].pred == INVALID) {
 		            nca = right;
 		            while (left_set.find(nca) == left_set.end()) {
 		              nca =
 		                _blossom_set->find(_graph.target((*_blossom_data)[nca].pred));
 		              right_path.push_back(nca);
 		              nca =
 		                _blossom_set->find(_graph.target((*_blossom_data)[nca].pred));
 		              right_path.push_back(nca);
+		            }
 		          } else {
 		            nca = left;
 		            while (right_set.find(nca) == right_set.end()) {
 		              nca =
 		                _blossom_set->find(_graph.target((*_blossom_data)[nca].pred));
 		              left_path.push_back(nca);
 		              nca =
 		                _blossom_set->find(_graph.target((*_blossom_data)[nca].pred));
 		              left_path.push_back(nca);
+		            }
+		          }
+		        }
+		      }
 		      std::vector<int> subblossoms;
 		      Arc prev;
 		      prev = _graph.direct(edge, true);
 		      for (int i = 0; left_path[i] != nca; i += 2) {
 		        subblossoms.push_back(left_path[i]);
 		        (*_blossom_data)[left_path[i]].next = prev;
 		        _tree_set->erase(left_path[i]);
 		        subblossoms.push_back(left_path[i + 1]);
 		        (*_blossom_data)[left_path[i + 1]].status = EVEN;
 		        oddToEven(left_path[i + 1], tree);
 		        _tree_set->erase(left_path[i + 1]);
 		        prev = _graph.oppositeArc((*_blossom_data)[left_path[i + 1]].pred);
+		      }
 		      int k = 0;
 		      while (right_path[k] != nca) ++k;
 		      subblossoms.push_back(nca);
 		      (*_blossom_data)[nca].next = prev;
 		      for (int i = k - 2; i >= 0; i -= 2) {
 		        subblossoms.push_back(right_path[i + 1]);
 		        (*_blossom_data)[right_path[i + 1]].status = EVEN;
 		        oddToEven(right_path[i + 1], tree);
 		        _tree_set->erase(right_path[i + 1]);
 		        (*_blossom_data)[right_path[i + 1]].next =
 		          (*_blossom_data)[right_path[i + 1]].pred;
 		        subblossoms.push_back(right_path[i]);
 		        _tree_set->erase(right_path[i]);
+		      }
 		      int surface =
 		        _blossom_set->join(subblossoms.begin(), subblossoms.end());
 		      for (int i = 0; i < int(subblossoms.size()); ++i) {
 		        if (!_blossom_set->trivial(subblossoms[i])) {
 		          (*_blossom_data)[subblossoms[i]].pot += 2 * _delta_sum;
+		        }
 		        (*_blossom_data)[subblossoms[i]].status = MATCHED;
+		      }
 		      (*_blossom_data)[surface].pot = -2 * _delta_sum;
 		      (*_blossom_data)[surface].offset = 0;
 		      (*_blossom_data)[surface].status = EVEN;
 		      (*_blossom_data)[surface].pred = (*_blossom_data)[nca].pred;
 		      (*_blossom_data)[surface].next = (*_blossom_data)[nca].pred;
 		      _tree_set->insert(surface, tree);
 		      _tree_set->erase(nca);
+		    }
 		    void splitBlossom(int blossom) {
 		      Arc next = (*_blossom_data)[blossom].next;
 		      Arc pred = (*_blossom_data)[blossom].pred;
 		      int tree = _tree_set->find(blossom);
 		      (*_blossom_data)[blossom].status = MATCHED;
 		      oddToMatched(blossom);
 		      if (_delta2->state(blossom) == _delta2->IN_HEAP) {
 		        _delta2->erase(blossom);
+		      }
 		      std::vector<int> subblossoms;
 		      _blossom_set->split(blossom, std::back_inserter(subblossoms));
 		      Value offset = (*_blossom_data)[blossom].offset;
 		      int b = _blossom_set->find(_graph.source(pred));
 		      int d = _blossom_set->find(_graph.source(next));
 		      int ib = -1, id = -1;
 		      for (int i = 0; i < int(subblossoms.size()); ++i) {
 		        if (subblossoms[i] == b) ib = i;
 		        if (subblossoms[i] == d) id = i;
 		        (*_blossom_data)[subblossoms[i]].offset = offset;
 		        if (!_blossom_set->trivial(subblossoms[i])) {
 		          (*_blossom_data)[subblossoms[i]].pot -= 2 * offset;
+		        }
 		        if (_blossom_set->classPrio(subblossoms[i]) !=
 		            std::numeric_limits<Value>::max()) {
 		          _delta2->push(subblossoms[i],
 		                        _blossom_set->classPrio(subblossoms[i]) -
 		                        (*_blossom_data)[subblossoms[i]].offset);
+		        }
+		      }
 		      if (id > ib ? ((id - ib) % 2 == 0) : ((ib - id) % 2 == 1)) {
 		        for (int i = (id + 1) % subblossoms.size();
 		             i != ib; i = (i + 2) % subblossoms.size()) {
 		          int sb = subblossoms[i];
 		          int tb = subblossoms[(i + 1) % subblossoms.size()];
 		          (*_blossom_data)[sb].next =
 		            _graph.oppositeArc((*_blossom_data)[tb].next);
+		        }
 		        for (int i = ib; i != id; i = (i + 2) % subblossoms.size()) {
 		          int sb = subblossoms[i];
 		          int tb = subblossoms[(i + 1) % subblossoms.size()];
 		          int ub = subblossoms[(i + 2) % subblossoms.size()];
 		          (*_blossom_data)[sb].status = ODD;
 		          matchedToOdd(sb);
 		          _tree_set->insert(sb, tree);
 		          (*_blossom_data)[sb].pred = pred;
 		          (*_blossom_data)[sb].next =
 		                           _graph.oppositeArc((*_blossom_data)[tb].next);
 		          pred = (*_blossom_data)[ub].next;
 		          (*_blossom_data)[tb].status = EVEN;
 		          matchedToEven(tb, tree);
 		          _tree_set->insert(tb, tree);
 		          (*_blossom_data)[tb].pred = (*_blossom_data)[tb].next;
+		        }
 		        (*_blossom_data)[subblossoms[id]].status = ODD;
 		        matchedToOdd(subblossoms[id]);
 		        _tree_set->insert(subblossoms[id], tree);
 		        (*_blossom_data)[subblossoms[id]].next = next;
 		        (*_blossom_data)[subblossoms[id]].pred = pred;
 		      } else {
 		        for (int i = (ib + 1) % subblossoms.size();
 		             i != id; i = (i + 2) % subblossoms.size()) {
 		          int sb = subblossoms[i];
 		          int tb = subblossoms[(i + 1) % subblossoms.size()];
 		          (*_blossom_data)[sb].next =
 		            _graph.oppositeArc((*_blossom_data)[tb].next);
+		        }
 		        for (int i = id; i != ib; i = (i + 2) % subblossoms.size()) {
 		          int sb = subblossoms[i];
 		          int tb = subblossoms[(i + 1) % subblossoms.size()];
 		          int ub = subblossoms[(i + 2) % subblossoms.size()];
 		          (*_blossom_data)[sb].status = ODD;
 		          matchedToOdd(sb);
 		          _tree_set->insert(sb, tree);
 		          (*_blossom_data)[sb].next = next;
 		          (*_blossom_data)[sb].pred =
 		            _graph.oppositeArc((*_blossom_data)[tb].next);
 		          (*_blossom_data)[tb].status = EVEN;
 		          matchedToEven(tb, tree);
 		          _tree_set->insert(tb, tree);
 		          (*_blossom_data)[tb].pred =
 		            (*_blossom_data)[tb].next =
 		            _graph.oppositeArc((*_blossom_data)[ub].next);
 		          next = (*_blossom_data)[ub].next;
+		        }
 		        (*_blossom_data)[subblossoms[ib]].status = ODD;
 		        matchedToOdd(subblossoms[ib]);
 		        _tree_set->insert(subblossoms[ib], tree);
 		        (*_blossom_data)[subblossoms[ib]].next = next;
 		        (*_blossom_data)[subblossoms[ib]].pred = pred;
+		      }
 		      _tree_set->erase(blossom);
+		    }
 		    void extractBlossom(int blossom, const Node& base, const Arc& matching) {
 		      if (_blossom_set->trivial(blossom)) {
 		        int bi = (*_node_index)[base];
 		        Value pot = (*_node_data)[bi].pot;
 		        (*_matching)[base] = matching;
 		        _blossom_node_list.push_back(base);
 		        (*_node_potential)[base] = pot;
 		      } else {
 		        Value pot = (*_blossom_data)[blossom].pot;
 		        int bn = _blossom_node_list.size();
 		        std::vector<int> subblossoms;
 		        _blossom_set->split(blossom, std::back_inserter(subblossoms));
 		        int b = _blossom_set->find(base);
 		        int ib = -1;
 		        for (int i = 0; i < int(subblossoms.size()); ++i) {
 		          if (subblossoms[i] == b) { ib = i; break; }
+		        }
 		        for (int i = 1; i < int(subblossoms.size()); i += 2) {
 		          int sb = subblossoms[(ib + i) % subblossoms.size()];
 		          int tb = subblossoms[(ib + i + 1) % subblossoms.size()];
 		          Arc m = (*_blossom_data)[tb].next;
 		          extractBlossom(sb, _graph.target(m), _graph.oppositeArc(m));
 		          extractBlossom(tb, _graph.source(m), m);
+		        }
 		        extractBlossom(subblossoms[ib], base, matching);
 		        int en = _blossom_node_list.size();
 		        _blossom_potential.push_back(BlossomVariable(bn, en, pot));
+		      }
+		    }
 		    void extractMatching() {
 		      std::vector<int> blossoms;
 		      for (typename BlossomSet::ClassIt c(*_blossom_set); c != INVALID; ++c) {
 		        blossoms.push_back(c);
+		      }
 		      for (int i = 0; i < int(blossoms.size()); ++i) {
 		        Value offset = (*_blossom_data)[blossoms[i]].offset;
 		        (*_blossom_data)[blossoms[i]].pot += 2 * offset;
 		        for (typename BlossomSet::ItemIt n(*_blossom_set, blossoms[i]);
 		             n != INVALID; ++n) {
 		          (*_node_data)[(*_node_index)[n]].pot -= offset;
+		        }
 		        Arc matching = (*_blossom_data)[blossoms[i]].next;
 		        Node base = _graph.source(matching);
 		        extractBlossom(blossoms[i], base, matching);
+		      }
+		    }
 		  public:
 		    /// \brief Constructor
 		    ///
 		    /// Constructor.
 		    MaxWeightedPerfectMatching(const Graph& graph, const WeightMap& weight)
 		      : _graph(graph), _weight(weight), _matching(0),
 		        _node_potential(0), _blossom_potential(), _blossom_node_list(),
 		        _node_num(0), _blossom_num(0),
 		        _blossom_index(0), _blossom_set(0), _blossom_data(0),
 		        _node_index(0), _node_heap_index(0), _node_data(0),
 		        _tree_set_index(0), _tree_set(0),
 		        _delta2_index(0), _delta2(0),
 		        _delta3_index(0), _delta3(0),
 		        _delta4_index(0), _delta4(0),
 		        _delta_sum(), _unmatched(0),
 		        _fractional(0)
 		    {}
 		    ~MaxWeightedPerfectMatching() {
 		      destroyStructures();
 		      if (_fractional) {
 		        delete _fractional;
+		      }
+		    }
 		    /// \name Execution Control
 		    /// The simplest way to execute the algorithm is to use the
 		    /// \ref run() member function.
 		    ///@{
 		    /// \brief Initialize the algorithm
 		    ///
 		    /// This function initializes the algorithm.
 		    void init() {
 		      createStructures();
 		      _blossom_node_list.clear();
 		      _blossom_potential.clear();
 		      for (ArcIt e(_graph); e != INVALID; ++e) {
 		        (*_node_heap_index)[e] = BinHeap<Value, IntArcMap>::PRE_HEAP;
+		      }
 		      for (EdgeIt e(_graph); e != INVALID; ++e) {
 		        (*_delta3_index)[e] = _delta3->PRE_HEAP;
+		      }
 		      for (int i = 0; i < _blossom_num; ++i) {
 		        (*_delta2_index)[i] = _delta2->PRE_HEAP;
 		        (*_delta4_index)[i] = _delta4->PRE_HEAP;
+		      }
 		      _unmatched = _node_num;
 		      _delta2->clear();
 		      _delta3->clear();
 		      _delta4->clear();
 		      _blossom_set->clear();
 		      _tree_set->clear();
 		      int index = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        Value max = - std::numeric_limits<Value>::max();
 		        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
 		          if (_graph.target(e) == n) continue;
 		          if ((dualScale * _weight[e]) / 2 > max) {
 		            max = (dualScale * _weight[e]) / 2;
+		          }
+		        }
 		        (*_node_index)[n] = index;
 		        (*_node_data)[index].heap_index.clear();
 		        (*_node_data)[index].heap.clear();
 		        (*_node_data)[index].pot = max;
 		        int blossom =
 		          _blossom_set->insert(n, std::numeric_limits<Value>::max());
 		        _tree_set->insert(blossom);
 		        (*_blossom_data)[blossom].status = EVEN;
 		        (*_blossom_data)[blossom].pred = INVALID;
 		        (*_blossom_data)[blossom].next = INVALID;
 		        (*_blossom_data)[blossom].pot = 0;
 		        (*_blossom_data)[blossom].offset = 0;
 		        ++index;
+		      }
 		      for (EdgeIt e(_graph); e != INVALID; ++e) {
 		        int si = (*_node_index)[_graph.u(e)];
 		        int ti = (*_node_index)[_graph.v(e)];
 		        if (_graph.u(e) != _graph.v(e)) {
 		          _delta3->push(e, ((*_node_data)[si].pot + (*_node_data)[ti].pot -
 		                            dualScale * _weight[e]) / 2);
+		        }
+		      }
+		    }
 		    /// \brief Initialize the algorithm with fractional matching
 		    ///
 		    /// This function initializes the algorithm with a fractional
 		    /// matching. This initialization is also called jumpstart heuristic.
 		    void fractionalInit() {
 		      createStructures();
 		      if (_fractional == 0) {
 		        _fractional = new FractionalMatching(_graph, _weight, false);
+		      }
 		      if (!_fractional->run()) {
 		        _unmatched = -1;
 		        return;
+		      }
 		      for (ArcIt e(_graph); e != INVALID; ++e) {
 		        (*_node_heap_index)[e] = BinHeap<Value, IntArcMap>::PRE_HEAP;
+		      }
 		      for (EdgeIt e(_graph); e != INVALID; ++e) {
 		        (*_delta3_index)[e] = _delta3->PRE_HEAP;
+		      }
 		      for (int i = 0; i < _blossom_num; ++i) {
 		        (*_delta2_index)[i] = _delta2->PRE_HEAP;
 		        (*_delta4_index)[i] = _delta4->PRE_HEAP;
+		      }
 		      _unmatched = 0;
 		      int index = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        Value pot = _fractional->nodeValue(n);
 		        (*_node_index)[n] = index;
 		        (*_node_data)[index].pot = pot;
 		        int blossom =
 		          _blossom_set->insert(n, std::numeric_limits<Value>::max());
 		        (*_blossom_data)[blossom].status = MATCHED;
 		        (*_blossom_data)[blossom].pred = INVALID;
 		        (*_blossom_data)[blossom].next = _fractional->matching(n);
 		        (*_blossom_data)[blossom].pot = 0;
 		        (*_blossom_data)[blossom].offset = 0;
 		        ++index;
+		      }
 		      typename Graph::template NodeMap<bool> processed(_graph, false);
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        if (processed[n]) continue;
 		        processed[n] = true;
 		        if (_fractional->matching(n) == INVALID) continue;
 		        int num = 1;
 		        Node v = _graph.target(_fractional->matching(n));
 		        while (n != v) {
 		          processed[v] = true;
 		          v = _graph.target(_fractional->matching(v));
 		          ++num;
+		        }
 		        if (num % 2 == 1) {
 		          std::vector<int> subblossoms(num);
 		          subblossoms[--num] = _blossom_set->find(n);
 		          v = _graph.target(_fractional->matching(n));
 		          while (n != v) {
 		            subblossoms[--num] = _blossom_set->find(v);
 		            v = _graph.target(_fractional->matching(v));
+		          }
 		          int surface =
 		            _blossom_set->join(subblossoms.begin(), subblossoms.end());
 		          (*_blossom_data)[surface].status = EVEN;
 		          (*_blossom_data)[surface].pred = INVALID;
 		          (*_blossom_data)[surface].next = INVALID;
 		          (*_blossom_data)[surface].pot = 0;
 		          (*_blossom_data)[surface].offset = 0;
 		          _tree_set->insert(surface);
 		          ++_unmatched;
+		        }
+		      }
 		      for (EdgeIt e(_graph); e != INVALID; ++e) {
 		        int si = (*_node_index)[_graph.u(e)];
 		        int sb = _blossom_set->find(_graph.u(e));
 		        int ti = (*_node_index)[_graph.v(e)];
 		        int tb = _blossom_set->find(_graph.v(e));
 		        if ((*_blossom_data)[sb].status == EVEN &&
 		            (*_blossom_data)[tb].status == EVEN && sb != tb) {
 		          _delta3->push(e, ((*_node_data)[si].pot + (*_node_data)[ti].pot -
 		                            dualScale * _weight[e]) / 2);
+		        }
+		      }
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        int nb = _blossom_set->find(n);
 		        if ((*_blossom_data)[nb].status != MATCHED) continue;
 		        int ni = (*_node_index)[n];
 		        for (OutArcIt e(_graph, n); e != INVALID; ++e) {
 		          Node v = _graph.target(e);
 		          int vb = _blossom_set->find(v);
 		          int vi = (*_node_index)[v];
 		          Value rw = (*_node_data)[ni].pot + (*_node_data)[vi].pot -
 		            dualScale * _weight[e];
 		          if ((*_blossom_data)[vb].status == EVEN) {
 		            int vt = _tree_set->find(vb);
 		            typename std::map<int, Arc>::iterator it =
 		              (*_node_data)[ni].heap_index.find(vt);
 		            if (it != (*_node_data)[ni].heap_index.end()) {
 		              if ((*_node_data)[ni].heap[it->second] > rw) {
 		                (*_node_data)[ni].heap.replace(it->second, e);
 		                (*_node_data)[ni].heap.decrease(e, rw);
 		                it->second = e;
+		              }
 		            } else {
 		              (*_node_data)[ni].heap.push(e, rw);
 		              (*_node_data)[ni].heap_index.insert(std::make_pair(vt, e));
+		            }
+		          }
+		        }
 		        if (!(*_node_data)[ni].heap.empty()) {
 		          _blossom_set->decrease(n, (*_node_data)[ni].heap.prio());
 		          _delta2->push(nb, _blossom_set->classPrio(nb));
+		        }
+		      }
+		    }
 		    /// \brief Start the algorithm
 		    ///
 		    /// This function starts the algorithm.
 		    ///
 		    /// \pre \ref init() or \ref fractionalInit() must be called before
 		    /// using this function.
 		    bool start() {
 		      enum OpType {
 		        D2, D3, D4
 		      };
 		      if (_unmatched == -1) return false;
 		      while (_unmatched > 0) {
 		        Value d2 = !_delta2->empty() ?
 		          _delta2->prio() : std::numeric_limits<Value>::max();
 		        Value d3 = !_delta3->empty() ?
 		          _delta3->prio() : std::numeric_limits<Value>::max();
 		        Value d4 = !_delta4->empty() ?
 		          _delta4->prio() : std::numeric_limits<Value>::max();
 		        _delta_sum = d3; OpType ot = D3;
 		        if (d2 < _delta_sum) { _delta_sum = d2; ot = D2; }
 		        if (d4 < _delta_sum) { _delta_sum = d4; ot = D4; }
 		        if (_delta_sum == std::numeric_limits<Value>::max()) {
 		          return false;
+		        }
 		        switch (ot) {
 		        case D2:
+		          {
 		            int blossom = _delta2->top();
 		            Node n = _blossom_set->classTop(blossom);
 		            Arc e = (*_node_data)[(*_node_index)[n]].heap.top();
 		            extendOnArc(e);
+		          }
 		          break;
 		        case D3:
+		          {
 		            Edge e = _delta3->top();
 		            int left_blossom = _blossom_set->find(_graph.u(e));
 		            int right_blossom = _blossom_set->find(_graph.v(e));
 		            if (left_blossom == right_blossom) {
 		              _delta3->pop();
 		            } else {
 		              int left_tree = _tree_set->find(left_blossom);
 		              int right_tree = _tree_set->find(right_blossom);
 		              if (left_tree == right_tree) {
 		                shrinkOnEdge(e, left_tree);
 		              } else {
 		                augmentOnEdge(e);
 		                _unmatched -= 2;
+		              }
+		            }
 		          } break;
 		        case D4:
 		          splitBlossom(_delta4->top());
 		          break;
+		        }
+		      }
 		      extractMatching();
 		      return true;
+		    }
 		    /// \brief Run the algorithm.
 		    ///
 		    /// This method runs the \c %MaxWeightedPerfectMatching algorithm.
 		    ///
 		    /// \note mwpm.run() is just a shortcut of the following code.
 		    /// \code
 		    ///   mwpm.fractionalInit();
 		    ///   mwpm.start();
 		    /// \endcode
 		    bool run() {
 		      fractionalInit();
 		      return start();
+		    }
 		    /// @}
 		    /// \name Primal Solution
 		    /// Functions to get the primal solution, i.e. the maximum weighted
 		    /// perfect matching.\n
 		    /// Either \ref run() or \ref start() function should be called before
 		    /// using them.
 		    /// @{
 		    /// \brief Return the weight of the matching.
 		    ///
 		    /// This function returns the weight of the found matching.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Value matchingWeight() const {
 		      Value sum = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        if ((*_matching)[n] != INVALID) {
 		          sum += _weight[(*_matching)[n]];
+		        }
+		      }
 		      return sum / 2;
+		    }
 		    /// \brief Return \c true if the given edge is in the matching.
 		    ///
 		    /// This function returns \c true if the given edge is in the found
 		    /// matching.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    bool matching(const Edge& edge) const {
 		      return static_cast<const Edge&>((*_matching)[_graph.u(edge)]) == edge;
+		    }
 		    /// \brief Return the matching arc (or edge) incident to the given node.
 		    ///
 		    /// This function returns the matching arc (or edge) incident to the
 		    /// given node in the found matching or \c INVALID if the node is
 		    /// not covered by the matching.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Arc matching(const Node& node) const {
 		      return (*_matching)[node];
+		    }
 		    /// \brief Return a const reference to the matching map.
 		    ///
 		    /// This function returns a const reference to a node map that stores
 		    /// the matching arc (or edge) incident to each node.
 		    const MatchingMap& matchingMap() const {
 		      return *_matching;
+		    }
 		    /// \brief Return the mate of the given node.
 		    ///
 		    /// This function returns the mate of the given node in the found
 		    /// matching or \c INVALID if the node is not covered by the matching.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Node mate(const Node& node) const {
 		      return _graph.target((*_matching)[node]);
+		    }
 		    /// @}
 		    /// \name Dual Solution
 		    /// Functions to get the dual solution.\n
 		    /// Either \ref run() or \ref start() function should be called before
 		    /// using them.
 		    /// @{
 		    /// \brief Return the value of the dual solution.
 		    ///
 		    /// This function returns the value of the dual solution.
 		    /// It should be equal to the primal value scaled by \ref dualScale
 		    /// "dual scale".
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Value dualValue() const {
 		      Value sum = 0;
 		      for (NodeIt n(_graph); n != INVALID; ++n) {
 		        sum += nodeValue(n);
+		      }
 		      for (int i = 0; i < blossomNum(); ++i) {
 		        sum += blossomValue(i) * (blossomSize(i) / 2);
+		      }
 		      return sum;
+		    }
 		    /// \brief Return the dual value (potential) of the given node.
 		    ///
 		    /// This function returns the dual value (potential) of the given node.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Value nodeValue(const Node& n) const {
 		      return (*_node_potential)[n];
+		    }
 		    /// \brief Return the number of the blossoms in the basis.
 		    ///
 		    /// This function returns the number of the blossoms in the basis.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    /// \see BlossomIt
 		    int blossomNum() const {
 		      return _blossom_potential.size();
+		    }
 		    /// \brief Return the number of the nodes in the given blossom.
 		    ///
 		    /// This function returns the number of the nodes in the given blossom.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    /// \see BlossomIt
 		    int blossomSize(int k) const {
 		      return _blossom_potential[k].end - _blossom_potential[k].begin;
+		    }
 		    /// \brief Return the dual value (ptential) of the given blossom.
 		    ///
 		    /// This function returns the dual value (ptential) of the given blossom.
 		    ///
 		    /// \pre Either run() or start() must be called before using this function.
 		    Value blossomValue(int k) const {
 		      return _blossom_potential[k].value;
+		    }
 		    /// \brief Iterator for obtaining the nodes of a blossom.
 		    ///
 		    /// This class provides an iterator for obtaining the nodes of the
 		    /// given blossom. It lists a subset of the nodes.
 		    /// Before using this iterator, you must allocate a
 		    /// MaxWeightedPerfectMatching class and execute it.
 		    class BlossomIt {
 		    public:
 		      /// \brief Constructor.
 		      ///
 		      /// Constructor to get the nodes of the given variable.
 		      ///
 		      /// \pre Either \ref MaxWeightedPerfectMatching::run() "algorithm.run()"
 		      /// or \ref MaxWeightedPerfectMatching::start() "algorithm.start()"
 		      /// must be called before initializing this iterator.
 		      BlossomIt(const MaxWeightedPerfectMatching& algorithm, int variable)
 		        : _algorithm(&algorithm)
+		      {
 		        _index = _algorithm->_blossom_potential[variable].begin;
 		        _last = _algorithm->_blossom_potential[variable].end;
+		      }
 		      /// \brief Conversion to \c Node.
 		      ///
 		      /// Conversion to \c Node.
 		      operator Node() const {
 		        return _algorithm->_blossom_node_list[_index];
+		      }
 		      /// \brief Increment operator.
 		      ///
 		      /// Increment operator.
 		      BlossomIt& operator++() {
 		        ++_index;
 		        return *this;
+		      }
 		      /// \brief Validity checking
 		      ///
 		      /// This function checks whether the iterator is invalid.
 		      bool operator==(Invalid) const { return _index == _last; }
 		      /// \brief Validity checking
 		      ///
 		      /// This function checks whether the iterator is valid.
 		      bool operator!=(Invalid) const { return _index != _last; }
 		    private:
 		      const MaxWeightedPerfectMatching* _algorithm;
 		      int _last;
 		      int _index;
 		    };
 		    /// @}
 		  };
 		} //END OF NAMESPACE LEMON
 		#endif //LEMON_MATCHING_H

lemon/unionfind.h

Show inline comments

@@ -523,1293 +523,1302 @@
 		      ///
 		      /// It steps to the next representant item.
 		      ClassIt& operator++() {
 		        cdx = unionFind->classes[cdx].next;
 		        return *this;
+		      }
 		      /// \brief Conversion operator
 		      ///
 		      /// It converts the iterator to the current representant item.
 		      operator int() const {
 		        return cdx;
+		      }
 		      /// \brief Equality operator
 		      ///
 		      /// Equality operator
 		      bool operator==(const ClassIt& i) {
 		        return i.cdx == cdx;
+		      }
 		      /// \brief Inequality operator
 		      ///
 		      /// Inequality operator
 		      bool operator!=(const ClassIt& i) {
 		        return i.cdx != cdx;
+		      }
 		    private:
 		      const UnionFindEnum* unionFind;
 		      int cdx;
 		    };
 		    /// \brief LEMON style iterator for the items of a component.
 		    ///
 		    /// ClassIt is a lemon style iterator for the components. It iterates
 		    /// on the items of a class. By example if you want to iterate on
 		    /// each items of each classes then you may write the next code.
 		    ///\code
 		    /// for (ClassIt cit(ufe); cit != INVALID; ++cit) {
 		    ///   std::cout << "Class: ";
 		    ///   for (ItemIt iit(ufe, cit); iit != INVALID; ++iit) {
 		    ///     std::cout << toString(iit) << ' ' << std::endl;
 		    ///   }
 		    ///   std::cout << std::endl;
 		    /// }
 		    ///\endcode
 		    class ItemIt {
 		    public:
 		      /// \brief Constructor of the iterator
 		      ///
 		      /// Constructor of the iterator. The iterator iterates
 		      /// on the class of the \c item.
 		      ItemIt(const UnionFindEnum& ufe, int cls) : unionFind(&ufe) {
 		        fdx = idx = unionFind->classes[cls].firstItem;
+		      }
 		      /// \brief Constructor to get invalid iterator
 		      ///
 		      /// Constructor to get invalid iterator
 		      ItemIt(Invalid) : unionFind(0), idx(-1) {}
 		      /// \brief Increment operator
 		      ///
 		      /// It steps to the next item in the class.
 		      ItemIt& operator++() {
 		        idx = unionFind->items[idx].next;
 		        if (idx == fdx) idx = -1;
 		        return *this;
+		      }
 		      /// \brief Conversion operator
 		      ///
 		      /// It converts the iterator to the current item.
 		      operator const Item&() const {
 		        return unionFind->items[idx].item;
+		      }
 		      /// \brief Equality operator
 		      ///
 		      /// Equality operator
 		      bool operator==(const ItemIt& i) {
 		        return i.idx == idx;
+		      }
 		      /// \brief Inequality operator
 		      ///
 		      /// Inequality operator
 		      bool operator!=(const ItemIt& i) {
 		        return i.idx != idx;
+		      }
 		    private:
 		      const UnionFindEnum* unionFind;
 		      int idx, fdx;
 		    };
 		  };
 		  /// \ingroup auxdat
 		  ///
 		  /// \brief A \e Extend-Find data structure implementation which
 		  /// is able to enumerate the components.
 		  ///
 		  /// The class implements an \e Extend-Find data structure which is
 		  /// able to enumerate the components and the items in a
 		  /// component. The data structure is a simplification of the
 		  /// Union-Find structure, and it does not allow to merge two components.
 		  ///
 		  /// \pre You need to add all the elements by the \ref insert()
 		  /// method.
 		  template <typename IM>
 		  class ExtendFindEnum {
 		  public:
 		    ///\e
 		    typedef IM ItemIntMap;
 		    ///\e
 		    typedef typename ItemIntMap::Key Item;
 		  private:
 		    ItemIntMap& index;
 		    struct ItemT {
 		      int cls;
 		      Item item;
 		      int next, prev;
 		    };
 		    std::vector<ItemT> items;
 		    int firstFreeItem;
 		    struct ClassT {
 		      int firstItem;
 		      int next, prev;
 		    };
 		    std::vector<ClassT> classes;
 		    int firstClass, firstFreeClass;
 		    int newClass() {
 		      if (firstFreeClass != -1) {
 		        int cdx = firstFreeClass;
 		        firstFreeClass = classes[cdx].next;
 		        return cdx;
 		      } else {
 		        classes.push_back(ClassT());
 		        return classes.size() - 1;
+		      }
+		    }
 		    int newItem() {
 		      if (firstFreeItem != -1) {
 		        int idx = firstFreeItem;
 		        firstFreeItem = items[idx].next;
 		        return idx;
 		      } else {
 		        items.push_back(ItemT());
 		        return items.size() - 1;
+		      }
+		    }
 		  public:
 		    /// \brief Constructor
 		    ExtendFindEnum(ItemIntMap& _index)
 		      : index(_index), items(), firstFreeItem(-1),
 		        classes(), firstClass(-1), firstFreeClass(-1) {}
 		    /// \brief Inserts the given element into a new component.
 		    ///
 		    /// This method creates a new component consisting only of the
 		    /// given element.
 		    int insert(const Item& item) {
 		      int cdx = newClass();
 		      classes[cdx].prev = -1;
 		      classes[cdx].next = firstClass;
 		      if (firstClass != -1) {
 		        classes[firstClass].prev = cdx;
+		      }
 		      firstClass = cdx;
 		      int idx = newItem();
 		      items[idx].item = item;
 		      items[idx].cls = cdx;
 		      items[idx].prev = idx;
 		      items[idx].next = idx;
 		      classes[cdx].firstItem = idx;
 		      index.set(item, idx);
 		      return cdx;
+		    }
 		    /// \brief Inserts the given element into the given component.
 		    ///
 		    /// This methods inserts the element \e item a into the \e cls class.
 		    void insert(const Item& item, int cls) {
 		      int idx = newItem();
 		      int rdx = classes[cls].firstItem;
 		      items[idx].item = item;
 		      items[idx].cls = cls;
 		      items[idx].prev = rdx;
 		      items[idx].next = items[rdx].next;
 		      items[items[rdx].next].prev = idx;
 		      items[rdx].next = idx;
 		      index.set(item, idx);
+		    }
 		    /// \brief Clears the union-find data structure
 		    ///
 		    /// Erase each item from the data structure.
 		    void clear() {
 		      items.clear();
 		      classes.clear();
 		      firstClass = firstFreeClass = firstFreeItem = -1;
+		    }
 		    /// \brief Gives back the class of the \e item.
 		    ///
 		    /// Gives back the class of the \e item.
 		    int find(const Item &item) const {
 		      return items[index[item]].cls;
+		    }
 		    /// \brief Gives back a representant item of the component.
 		    ///
 		    /// Gives back a representant item of the component.
 		    Item item(int cls) const {
 		      return items[classes[cls].firstItem].item;
+		    }
 		    /// \brief Removes the given element from the structure.
 		    ///
 		    /// Removes the element from its component and if the component becomes
 		    /// empty then removes that component from the component list.
 		    ///
 		    /// \warning It is an error to remove an element which is not in
 		    /// the structure.
 		    void erase(const Item &item) {
 		      int idx = index[item];
 		      int cdx = items[idx].cls;
 		      if (idx == items[idx].next) {
 		        if (classes[cdx].prev != -1) {
 		          classes[classes[cdx].prev].next = classes[cdx].next;
 		        } else {
 		          firstClass = classes[cdx].next;
+		        }
 		        if (classes[cdx].next != -1) {
 		          classes[classes[cdx].next].prev = classes[cdx].prev;
+		        }
 		        classes[cdx].next = firstFreeClass;
 		        firstFreeClass = cdx;
 		      } else {
 		        classes[cdx].firstItem = items[idx].next;
 		        items[items[idx].next].prev = items[idx].prev;
 		        items[items[idx].prev].next = items[idx].next;
+		      }
 		      items[idx].next = firstFreeItem;
 		      firstFreeItem = idx;
+		    }
 		    /// \brief Removes the component of the given element from the structure.
 		    ///
 		    /// Removes the component of the given element from the structure.
 		    ///
 		    /// \warning It is an error to give an element which is not in the
 		    /// structure.
 		    void eraseClass(int cdx) {
 		      int idx = classes[cdx].firstItem;
 		      items[items[idx].prev].next = firstFreeItem;
 		      firstFreeItem = idx;
 		      if (classes[cdx].prev != -1) {
 		        classes[classes[cdx].prev].next = classes[cdx].next;
 		      } else {
 		        firstClass = classes[cdx].next;
+		      }
 		      if (classes[cdx].next != -1) {
 		        classes[classes[cdx].next].prev = classes[cdx].prev;
+		      }
 		      classes[cdx].next = firstFreeClass;
 		      firstFreeClass = cdx;
+		    }
 		    /// \brief LEMON style iterator for the classes.
 		    ///
 		    /// ClassIt is a lemon style iterator for the components. It iterates
 		    /// on the ids of classes.
 		    class ClassIt {
 		    public:
 		      /// \brief Constructor of the iterator
 		      ///
 		      /// Constructor of the iterator
 		      ClassIt(const ExtendFindEnum& ufe) : extendFind(&ufe) {
 		        cdx = extendFind->firstClass;
+		      }
 		      /// \brief Constructor to get invalid iterator
 		      ///
 		      /// Constructor to get invalid iterator
 		      ClassIt(Invalid) : extendFind(0), cdx(-1) {}
 		      /// \brief Increment operator
 		      ///
 		      /// It steps to the next representant item.
 		      ClassIt& operator++() {
 		        cdx = extendFind->classes[cdx].next;
 		        return *this;
+		      }
 		      /// \brief Conversion operator
 		      ///
 		      /// It converts the iterator to the current class id.
 		      operator int() const {
 		        return cdx;
+		      }
 		      /// \brief Equality operator
 		      ///
 		      /// Equality operator
 		      bool operator==(const ClassIt& i) {
 		        return i.cdx == cdx;
+		      }
 		      /// \brief Inequality operator
 		      ///
 		      /// Inequality operator
 		      bool operator!=(const ClassIt& i) {
 		        return i.cdx != cdx;
+		      }
 		    private:
 		      const ExtendFindEnum* extendFind;
 		      int cdx;
 		    };
 		    /// \brief LEMON style iterator for the items of a component.
 		    ///
 		    /// ClassIt is a lemon style iterator for the components. It iterates
 		    /// on the items of a class. By example if you want to iterate on
 		    /// each items of each classes then you may write the next code.
 		    ///\code
 		    /// for (ClassIt cit(ufe); cit != INVALID; ++cit) {
 		    ///   std::cout << "Class: ";
 		    ///   for (ItemIt iit(ufe, cit); iit != INVALID; ++iit) {
 		    ///     std::cout << toString(iit) << ' ' << std::endl;
 		    ///   }
 		    ///   std::cout << std::endl;
 		    /// }
 		    ///\endcode
 		    class ItemIt {
 		    public:
 		      /// \brief Constructor of the iterator
 		      ///
 		      /// Constructor of the iterator. The iterator iterates
 		      /// on the class of the \c item.
 		      ItemIt(const ExtendFindEnum& ufe, int cls) : extendFind(&ufe) {
 		        fdx = idx = extendFind->classes[cls].firstItem;
+		      }
 		      /// \brief Constructor to get invalid iterator
 		      ///
 		      /// Constructor to get invalid iterator
 		      ItemIt(Invalid) : extendFind(0), idx(-1) {}
 		      /// \brief Increment operator
 		      ///
 		      /// It steps to the next item in the class.
 		      ItemIt& operator++() {
 		        idx = extendFind->items[idx].next;
 		        if (fdx == idx) idx = -1;
 		        return *this;
+		      }
 		      /// \brief Conversion operator
 		      ///
 		      /// It converts the iterator to the current item.
 		      operator const Item&() const {
 		        return extendFind->items[idx].item;
+		      }
 		      /// \brief Equality operator
 		      ///
 		      /// Equality operator
 		      bool operator==(const ItemIt& i) {
 		        return i.idx == idx;
+		      }
 		      /// \brief Inequality operator
 		      ///
 		      /// Inequality operator
 		      bool operator!=(const ItemIt& i) {
 		        return i.idx != idx;
+		      }
 		    private:
 		      const ExtendFindEnum* extendFind;
 		      int idx, fdx;
 		    };
 		  };
 		  /// \ingroup auxdat
 		  ///
 		  /// \brief A \e Union-Find data structure implementation which
 		  /// is able to store a priority for each item and retrieve the minimum of
 		  /// each class.
 		  ///
 		  /// A \e Union-Find data structure implementation which is able to
 		  /// store a priority for each item and retrieve the minimum of each
 		  /// class. In addition, it supports the joining and splitting the
 		  /// components. If you don't need this feature then you makes
 		  /// better to use the \ref UnionFind class which is more efficient.
 		  ///
 		  /// The union-find data strcuture based on a (2, 16)-tree with a
 		  /// tournament minimum selection on the internal nodes. The insert
 		  /// operation takes O(1), the find, set, decrease and increase takes
 		  /// O(log(n)), where n is the number of nodes in the current
 		  /// component.  The complexity of join and split is O(log(n)*k),
 		  /// where n is the sum of the number of the nodes and k is the
 		  /// number of joined components or the number of the components
 		  /// after the split.
 		  ///
 		  /// \pre You need to add all the elements by the \ref insert()
 		  /// method.
 		  template <typename V, typename IM, typename Comp = std::less<V> >
 		  class HeapUnionFind {
 		  public:
 		    ///\e
 		    typedef V Value;
 		    ///\e
 		    typedef typename IM::Key Item;
 		    ///\e
 		    typedef IM ItemIntMap;
 		    ///\e
 		    typedef Comp Compare;
 		  private:
 		    static const int cmax = 16;
 		    ItemIntMap& index;
 		    struct ClassNode {
 		      int parent;
 		      int depth;
 		      int left, right;
 		      int next, prev;
 		    };
 		    int first_class;
 		    int first_free_class;
 		    std::vector<ClassNode> classes;
 		    int newClass() {
 		      if (first_free_class < 0) {
 		        int id = classes.size();
 		        classes.push_back(ClassNode());
 		        return id;
 		      } else {
 		        int id = first_free_class;
 		        first_free_class = classes[id].next;
 		        return id;
+		      }
+		    }
 		    void deleteClass(int id) {
 		      classes[id].next = first_free_class;
 		      first_free_class = id;
+		    }
 		    struct ItemNode {
 		      int parent;
 		      Item item;
 		      Value prio;
 		      int next, prev;
 		      int left, right;
 		      int size;
 		    };
 		    int first_free_node;
 		    std::vector<ItemNode> nodes;
 		    int newNode() {
 		      if (first_free_node < 0) {
 		        int id = nodes.size();
 		        nodes.push_back(ItemNode());
 		        return id;
 		      } else {
 		        int id = first_free_node;
 		        first_free_node = nodes[id].next;
 		        return id;
+		      }
+		    }
 		    void deleteNode(int id) {
 		      nodes[id].next = first_free_node;
 		      first_free_node = id;
+		    }
 		    Comp comp;
 		    int findClass(int id) const {
 		      int kd = id;
 		      while (kd >= 0) {
 		        kd = nodes[kd].parent;
+		      }
 		      return ~kd;
+		    }
 		    int leftNode(int id) const {
 		      int kd = ~(classes[id].parent);
 		      for (int i = 0; i < classes[id].depth; ++i) {
 		        kd = nodes[kd].left;
+		      }
 		      return kd;
+		    }
 		    int nextNode(int id) const {
 		      int depth = 0;
 		      while (id >= 0 && nodes[id].next == -1) {
 		        id = nodes[id].parent;
 		        ++depth;
+		      }
 		      if (id < 0) {
 		        return -1;
+		      }
 		      id = nodes[id].next;
 		      while (depth--) {
 		        id = nodes[id].left;
+		      }
 		      return id;
+		    }
 		    void setPrio(int id) {
 		      int jd = nodes[id].left;
 		      nodes[id].prio = nodes[jd].prio;
 		      nodes[id].item = nodes[jd].item;
 		      jd = nodes[jd].next;
 		      while (jd != -1) {
 		        if (comp(nodes[jd].prio, nodes[id].prio)) {
 		          nodes[id].prio = nodes[jd].prio;
 		          nodes[id].item = nodes[jd].item;
+		        }
 		        jd = nodes[jd].next;
+		      }
+		    }
 		    void push(int id, int jd) {
 		      nodes[id].size = 1;
 		      nodes[id].left = nodes[id].right = jd;
 		      nodes[jd].next = nodes[jd].prev = -1;
 		      nodes[jd].parent = id;
+		    }
 		    void pushAfter(int id, int jd) {
 		      int kd = nodes[id].parent;
 		      if (nodes[id].next != -1) {
 		        nodes[nodes[id].next].prev = jd;
 		        if (kd >= 0) {
 		          nodes[kd].size += 1;
+		        }
 		      } else {
 		        if (kd >= 0) {
 		          nodes[kd].right = jd;
 		          nodes[kd].size += 1;
+		        }
+		      }
 		      nodes[jd].next = nodes[id].next;
 		      nodes[jd].prev = id;
 		      nodes[id].next = jd;
 		      nodes[jd].parent = kd;
+		    }
 		    void pushRight(int id, int jd) {
 		      nodes[id].size += 1;
 		      nodes[jd].prev = nodes[id].right;
 		      nodes[jd].next = -1;
 		      nodes[nodes[id].right].next = jd;
 		      nodes[id].right = jd;
 		      nodes[jd].parent = id;
+		    }
 		    void popRight(int id) {
 		      nodes[id].size -= 1;
 		      int jd = nodes[id].right;
 		      nodes[nodes[jd].prev].next = -1;
 		      nodes[id].right = nodes[jd].prev;
+		    }
 		    void splice(int id, int jd) {
 		      nodes[id].size += nodes[jd].size;
 		      nodes[nodes[id].right].next = nodes[jd].left;
 		      nodes[nodes[jd].left].prev = nodes[id].right;
 		      int kd = nodes[jd].left;
 		      while (kd != -1) {
 		        nodes[kd].parent = id;
 		        kd = nodes[kd].next;
+		      }
 		      nodes[id].right = nodes[jd].right;
+		    }
 		    void split(int id, int jd) {
 		      int kd = nodes[id].parent;
 		      nodes[kd].right = nodes[id].prev;
 		      nodes[nodes[id].prev].next = -1;
 		      nodes[jd].left = id;
 		      nodes[id].prev = -1;
 		      int num = 0;
 		      while (id != -1) {
 		        nodes[id].parent = jd;
 		        nodes[jd].right = id;
 		        id = nodes[id].next;
 		        ++num;
+		      }
 		      nodes[kd].size -= num;
 		      nodes[jd].size = num;
+		    }
 		    void pushLeft(int id, int jd) {
 		      nodes[id].size += 1;
 		      nodes[jd].next = nodes[id].left;
 		      nodes[jd].prev = -1;
 		      nodes[nodes[id].left].prev = jd;
 		      nodes[id].left = jd;
 		      nodes[jd].parent = id;
+		    }
 		    void popLeft(int id) {
 		      nodes[id].size -= 1;
 		      int jd = nodes[id].left;
 		      nodes[nodes[jd].next].prev = -1;
 		      nodes[id].left = nodes[jd].next;
+		    }
 		    void repairLeft(int id) {
 		      int jd = ~(classes[id].parent);
 		      while (nodes[jd].left != -1) {
 		        int kd = nodes[jd].left;
 		        if (nodes[jd].size == 1) {
 		          if (nodes[jd].parent < 0) {
 		            classes[id].parent = ~kd;
 		            classes[id].depth -= 1;
 		            nodes[kd].parent = ~id;
 		            deleteNode(jd);
 		            jd = kd;
 		          } else {
 		            int pd = nodes[jd].parent;
 		            if (nodes[nodes[jd].next].size < cmax) {
 		              pushLeft(nodes[jd].next, nodes[jd].left);
 		              if (less(jd, nodes[jd].next) ||
 		                  nodes[jd].item == nodes[pd].item) {
 		                nodes[nodes[jd].next].prio = nodes[jd].prio;
 		                nodes[nodes[jd].next].item = nodes[jd].item;
+		              }
 		              popLeft(pd);
 		              deleteNode(jd);
 		              jd = pd;
 		            } else {
 		              int ld = nodes[nodes[jd].next].left;
 		              popLeft(nodes[jd].next);
 		              pushRight(jd, ld);
 		              if (less(ld, nodes[jd].left) ||
 		                  nodes[ld].item == nodes[pd].item) {
 		                nodes[jd].item = nodes[ld].item;
 		                nodes[jd].prio = nodes[ld].prio;
+		              }
 		              if (nodes[nodes[jd].next].item == nodes[ld].item) {
 		                setPrio(nodes[jd].next);
+		              }
 		              jd = nodes[jd].left;
+		            }
+		          }
 		        } else {
 		          jd = nodes[jd].left;
+		        }
+		      }
+		    }
 		    void repairRight(int id) {
 		      int jd = ~(classes[id].parent);
 		      while (nodes[jd].right != -1) {
 		        int kd = nodes[jd].right;
 		        if (nodes[jd].size == 1) {
 		          if (nodes[jd].parent < 0) {
 		            classes[id].parent = ~kd;
 		            classes[id].depth -= 1;
 		            nodes[kd].parent = ~id;
 		            deleteNode(jd);
 		            jd = kd;
 		          } else {
 		            int pd = nodes[jd].parent;
 		            if (nodes[nodes[jd].prev].size < cmax) {
 		              pushRight(nodes[jd].prev, nodes[jd].right);
 		              if (less(jd, nodes[jd].prev) ||
 		                  nodes[jd].item == nodes[pd].item) {
 		                nodes[nodes[jd].prev].prio = nodes[jd].prio;
 		                nodes[nodes[jd].prev].item = nodes[jd].item;
+		              }
 		              popRight(pd);
 		              deleteNode(jd);
 		              jd = pd;
 		            } else {
 		              int ld = nodes[nodes[jd].prev].right;
 		              popRight(nodes[jd].prev);
 		              pushLeft(jd, ld);
 		              if (less(ld, nodes[jd].right) ||
 		                  nodes[ld].item == nodes[pd].item) {
 		                nodes[jd].item = nodes[ld].item;
 		                nodes[jd].prio = nodes[ld].prio;
+		              }
 		              if (nodes[nodes[jd].prev].item == nodes[ld].item) {
 		                setPrio(nodes[jd].prev);
+		              }
 		              jd = nodes[jd].right;
+		            }
+		          }
 		        } else {
 		          jd = nodes[jd].right;
+		        }
+		      }
+		    }
 		    bool less(int id, int jd) const {
 		      return comp(nodes[id].prio, nodes[jd].prio);
+		    }
 		  public:
 		    /// \brief Returns true when the given class is alive.
 		    ///
 		    /// Returns true when the given class is alive, ie. the class is
 		    /// not nested into other class.
 		    bool alive(int cls) const {
 		      return classes[cls].parent < 0;
+		    }
 		    /// \brief Returns true when the given class is trivial.
 		    ///
 		    /// Returns true when the given class is trivial, ie. the class
 		    /// contains just one item directly.
 		    bool trivial(int cls) const {
 		      return classes[cls].left == -1;
+		    }
 		    /// \brief Constructs the union-find.
 		    ///
 		    /// Constructs the union-find.
 		    /// \brief _index The index map of the union-find. The data
 		    /// structure uses internally for store references.
 		    HeapUnionFind(ItemIntMap& _index)
 		      : index(_index), first_class(-1),
 		        first_free_class(-1), first_free_node(-1) {}
 		    /// \brief Clears the union-find data structure
 		    ///
 		    /// Erase each item from the data structure.
 		    void clear() {
 		      nodes.clear();
 		      classes.clear();
 		      first_free_node = first_free_class = first_class = -1;
+		    }
 		    /// \brief Insert a new node into a new component.
 		    ///
 		    /// Insert a new node into a new component.
 		    /// \param item The item of the new node.
 		    /// \param prio The priority of the new node.
 		    /// \return The class id of the one-item-heap.
 		    int insert(const Item& item, const Value& prio) {
 		      int id = newNode();
 		      nodes[id].item = item;
 		      nodes[id].prio = prio;
 		      nodes[id].size = 0;
 		      nodes[id].prev = -1;
 		      nodes[id].next = -1;
 		      nodes[id].left = -1;
 		      nodes[id].right = -1;
 		      nodes[id].item = item;
 		      index[item] = id;
 		      int class_id = newClass();
 		      classes[class_id].parent = ~id;
 		      classes[class_id].depth = 0;
 		      classes[class_id].left = -1;
 		      classes[class_id].right = -1;
 		      if (first_class != -1) {
 		        classes[first_class].prev = class_id;
+		      }
 		      classes[class_id].next = first_class;
 		      classes[class_id].prev = -1;
 		      first_class = class_id;
 		      nodes[id].parent = ~class_id;
 		      return class_id;
+		    }
 		    /// \brief The class of the item.
 		    ///
 		    /// \return The alive class id of the item, which is not nested into
 		    /// other classes.
 		    ///
 		    /// The time complexity is O(log(n)).
 		    int find(const Item& item) const {
 		      return findClass(index[item]);
+		    }
 		    /// \brief Joins the classes.
 		    ///
 		    /// The current function joins the given classes. The parameter is
 		    /// an STL range which should be contains valid class ids. The
 		    /// time complexity is O(log(n)*k) where n is the overall number
 		    /// of the joined nodes and k is the number of classes.
 		    /// \return The class of the joined classes.
 		    /// \pre The range should contain at least two class ids.
 		    template <typename Iterator>
 		    int join(Iterator begin, Iterator end) {
 		      std::vector<int> cs;
 		      for (Iterator it = begin; it != end; ++it) {
 		        cs.push_back(*it);
+		      }
 		      int class_id = newClass();
 		      { // creation union-find
 		        if (first_class != -1) {
 		          classes[first_class].prev = class_id;
+		        }
 		        classes[class_id].next = first_class;
 		        classes[class_id].prev = -1;
 		        first_class = class_id;
 		        classes[class_id].depth = classes[cs[0]].depth;
 		        classes[class_id].parent = classes[cs[0]].parent;
 		        nodes[~(classes[class_id].parent)].parent = ~class_id;
 		        int l = cs[0];
 		        classes[class_id].left = l;
 		        classes[class_id].right = l;
 		        if (classes[l].next != -1) {
 		          classes[classes[l].next].prev = classes[l].prev;
+		        }
 		        classes[classes[l].prev].next = classes[l].next;
 		        classes[l].prev = -1;
 		        classes[l].next = -1;
 		        classes[l].depth = leftNode(l);
 		        classes[l].parent = class_id;
+		      }
 		      { // merging of heap
 		        int l = class_id;
 		        for (int ci = 1; ci < int(cs.size()); ++ci) {
 		          int r = cs[ci];
 		          int rln = leftNode(r);
 		          if (classes[l].depth > classes[r].depth) {
 		            int id = ~(classes[l].parent);
 		            for (int i = classes[r].depth + 1; i < classes[l].depth; ++i) {
 		              id = nodes[id].right;
+		            }
 		            while (id >= 0 && nodes[id].size == cmax) {
 		              int new_id = newNode();
 		              int right_id = nodes[id].right;
 		              popRight(id);
 		              if (nodes[id].item == nodes[right_id].item) {
 		                setPrio(id);
+		              }
 		              push(new_id, right_id);
 		              pushRight(new_id, ~(classes[r].parent));
 		              if (less(~classes[r].parent, right_id)) {
 		                nodes[new_id].item = nodes[~classes[r].parent].item;
 		                nodes[new_id].prio = nodes[~classes[r].parent].prio;
 		              } else {
 		                nodes[new_id].item = nodes[right_id].item;
 		                nodes[new_id].prio = nodes[right_id].prio;
+		              }
 		              id = nodes[id].parent;
 		              classes[r].parent = ~new_id;
+		            }
 		            if (id < 0) {
 		              int new_parent = newNode();
 		              nodes[new_parent].next = -1;
 		              nodes[new_parent].prev = -1;
 		              nodes[new_parent].parent = ~l;
 		              push(new_parent, ~(classes[l].parent));
 		              pushRight(new_parent, ~(classes[r].parent));
 		              setPrio(new_parent);
 		              classes[l].parent = ~new_parent;
 		              classes[l].depth += 1;
 		            } else {
 		              pushRight(id, ~(classes[r].parent));
 		              while (id >= 0 && less(~(classes[r].parent), id)) {
 		                nodes[id].prio = nodes[~(classes[r].parent)].prio;
 		                nodes[id].item = nodes[~(classes[r].parent)].item;
 		                id = nodes[id].parent;
+		              }
+		            }
 		          } else if (classes[r].depth > classes[l].depth) {
 		            int id = ~(classes[r].parent);
 		            for (int i = classes[l].depth + 1; i < classes[r].depth; ++i) {
 		              id = nodes[id].left;
+		            }
 		            while (id >= 0 && nodes[id].size == cmax) {
 		              int new_id = newNode();
 		              int left_id = nodes[id].left;
 		              popLeft(id);
 		              if (nodes[id].prio == nodes[left_id].prio) {
 		                setPrio(id);
+		              }
 		              push(new_id, left_id);
 		              pushLeft(new_id, ~(classes[l].parent));
 		              if (less(~classes[l].parent, left_id)) {
 		                nodes[new_id].item = nodes[~classes[l].parent].item;
 		                nodes[new_id].prio = nodes[~classes[l].parent].prio;
 		              } else {
 		                nodes[new_id].item = nodes[left_id].item;
 		                nodes[new_id].prio = nodes[left_id].prio;
+		              }
 		              id = nodes[id].parent;
 		              classes[l].parent = ~new_id;
+		            }
 		            if (id < 0) {
 		              int new_parent = newNode();
 		              nodes[new_parent].next = -1;
 		              nodes[new_parent].prev = -1;
 		              nodes[new_parent].parent = ~l;
 		              push(new_parent, ~(classes[r].parent));
 		              pushLeft(new_parent, ~(classes[l].parent));
 		              setPrio(new_parent);
 		              classes[r].parent = ~new_parent;
 		              classes[r].depth += 1;
 		            } else {
 		              pushLeft(id, ~(classes[l].parent));
 		              while (id >= 0 && less(~(classes[l].parent), id)) {
 		                nodes[id].prio = nodes[~(classes[l].parent)].prio;
 		                nodes[id].item = nodes[~(classes[l].parent)].item;
 		                id = nodes[id].parent;
+		              }
+		            }
 		            nodes[~(classes[r].parent)].parent = ~l;
 		            classes[l].parent = classes[r].parent;
 		            classes[l].depth = classes[r].depth;
 		          } else {
 		            if (classes[l].depth != 0 &&
 		                nodes[~(classes[l].parent)].size +
 		                nodes[~(classes[r].parent)].size <= cmax) {
 		              splice(~(classes[l].parent), ~(classes[r].parent));
 		              deleteNode(~(classes[r].parent));
 		              if (less(~(classes[r].parent), ~(classes[l].parent))) {
 		                nodes[~(classes[l].parent)].prio =
 		                  nodes[~(classes[r].parent)].prio;
 		                nodes[~(classes[l].parent)].item =
 		                  nodes[~(classes[r].parent)].item;
+		              }
 		            } else {
 		              int new_parent = newNode();
 		              nodes[new_parent].next = nodes[new_parent].prev = -1;
 		              push(new_parent, ~(classes[l].parent));
 		              pushRight(new_parent, ~(classes[r].parent));
 		              setPrio(new_parent);
 		              classes[l].parent = ~new_parent;
 		              classes[l].depth += 1;
 		              nodes[new_parent].parent = ~l;
+		            }
+		          }
 		          if (classes[r].next != -1) {
 		            classes[classes[r].next].prev = classes[r].prev;
+		          }
 		          classes[classes[r].prev].next = classes[r].next;
 		          classes[r].prev = classes[l].right;
 		          classes[classes[l].right].next = r;
 		          classes[l].right = r;
 		          classes[r].parent = l;
 		          classes[r].next = -1;
 		          classes[r].depth = rln;
+		        }
+		      }
 		      return class_id;
+		    }
 		    /// \brief Split the class to subclasses.
 		    ///
 		    /// The current function splits the given class. The join, which
 		    /// made the current class, stored a reference to the
 		    /// subclasses. The \c splitClass() member restores the classes
 		    /// and creates the heaps. The parameter is an STL output iterator
 		    /// which will be filled with the subclass ids. The time
 		    /// complexity is O(log(n)*k) where n is the overall number of
 		    /// nodes in the splitted classes and k is the number of the
 		    /// classes.
 		    template <typename Iterator>
 		    void split(int cls, Iterator out) {
 		      std::vector<int> cs;
 		      { // splitting union-find
 		        int id = cls;
 		        int l = classes[id].left;
 		        classes[l].parent = classes[id].parent;
 		        classes[l].depth = classes[id].depth;
 		        nodes[~(classes[l].parent)].parent = ~l;
 		        *out++ = l;
 		        while (l != -1) {
 		          cs.push_back(l);
 		          l = classes[l].next;
+		        }
 		        classes[classes[id].right].next = first_class;
 		        classes[first_class].prev = classes[id].right;
 		        first_class = classes[id].left;
 		        if (classes[id].next != -1) {
 		          classes[classes[id].next].prev = classes[id].prev;
+		        }
 		        classes[classes[id].prev].next = classes[id].next;
 		        deleteClass(id);
+		      }
+		      {
 		        for (int i = 1; i < int(cs.size()); ++i) {
 		          int l = classes[cs[i]].depth;
 		          while (nodes[nodes[l].parent].left == l) {
 		            l = nodes[l].parent;
+		          }
 		          int r = l;
 		          while (nodes[l].parent >= 0) {
 		            l = nodes[l].parent;
 		            int new_node = newNode();
 		            nodes[new_node].prev = -1;
 		            nodes[new_node].next = -1;
 		            split(r, new_node);
 		            pushAfter(l, new_node);
 		            setPrio(l);
 		            setPrio(new_node);
 		            r = new_node;
+		          }
 		          classes[cs[i]].parent = ~r;
 		          classes[cs[i]].depth = classes[~(nodes[l].parent)].depth;
 		          nodes[r].parent = ~cs[i];
 		          nodes[l].next = -1;
 		          nodes[r].prev = -1;
 		          repairRight(~(nodes[l].parent));
 		          repairLeft(cs[i]);
 		          *out++ = cs[i];
+		        }
+		      }
+		    }
 		    /// \brief Gives back the priority of the current item.
 		    ///
 		    /// Gives back the priority of the current item.
 		    const Value& operator[](const Item& item) const {
 		      return nodes[index[item]].prio;
+		    }
 		    /// \brief Sets the priority of the current item.
 		    ///
 		    /// Sets the priority of the current item.
 		    void set(const Item& item, const Value& prio) {
 		      if (comp(prio, nodes[index[item]].prio)) {
 		        decrease(item, prio);
 		      } else if (!comp(prio, nodes[index[item]].prio)) {
 		        increase(item, prio);
+		      }
+		    }
 		    /// \brief Increase the priority of the current item.
 		    ///
 		    /// Increase the priority of the current item.
 		    void increase(const Item& item, const Value& prio) {
 		      int id = index[item];
 		      int kd = nodes[id].parent;
 		      nodes[id].prio = prio;
 		      while (kd >= 0 && nodes[kd].item == item) {
 		        setPrio(kd);
 		        kd = nodes[kd].parent;
+		      }
+		    }
 		    /// \brief Increase the priority of the current item.
 		    ///
 		    /// Increase the priority of the current item.
 		    void decrease(const Item& item, const Value& prio) {
 		      int id = index[item];
 		      int kd = nodes[id].parent;
 		      nodes[id].prio = prio;
 		      while (kd >= 0 && less(id, kd)) {
 		        nodes[kd].prio = prio;
 		        nodes[kd].item = item;
 		        kd = nodes[kd].parent;
+		      }
+		    }
 		    /// \brief Gives back the minimum priority of the class.
 		    ///
 		    /// Gives back the minimum priority of the class.
 		    const Value& classPrio(int cls) const {
 		      return nodes[~(classes[cls].parent)].prio;
+		    }
 		    /// \brief Gives back the minimum priority item of the class.
 		    ///
 		    /// \return Gives back the minimum priority item of the class.
 		    const Item& classTop(int cls) const {
 		      return nodes[~(classes[cls].parent)].item;
+		    }
 		    /// \brief Gives back a representant item of the class.
 		    ///
 		    /// Gives back a representant item of the class.
 		    /// The representant is indpendent from the priorities of the
 		    /// items.
 		    const Item& classRep(int id) const {
 		      int parent = classes[id].parent;
 		      return nodes[parent >= 0 ? classes[id].depth : leftNode(id)].item;
+		    }
 		    /// \brief LEMON style iterator for the items of a class.
 		    ///
 		    /// ClassIt is a lemon style iterator for the components. It iterates
 		    /// on the items of a class. By example if you want to iterate on
 		    /// each items of each classes then you may write the next code.
 		    ///\code
 		    /// for (ClassIt cit(huf); cit != INVALID; ++cit) {
 		    ///   std::cout << "Class: ";
 		    ///   for (ItemIt iit(huf, cit); iit != INVALID; ++iit) {
 		    ///     std::cout << toString(iit) << ' ' << std::endl;
 		    ///   }
 		    ///   std::cout << std::endl;
 		    /// }
 		    ///\endcode
 		    class ItemIt {
 		    private:
 		      const HeapUnionFind* _huf;
 		      int _id, _lid;
 		    public:
 		      /// \brief Default constructor
 		      ///
 		      /// Default constructor
 		      ItemIt() {}
 		      ItemIt(const HeapUnionFind& huf, int cls) : _huf(&huf) {
 		        int id = cls;
 		        int parent = _huf->classes[id].parent;
 		        if (parent >= 0) {
 		          _id = _huf->classes[id].depth;
 		          if (_huf->classes[id].next != -1) {
 		            _lid = _huf->classes[_huf->classes[id].next].depth;
 		          } else {
 		            _lid = -1;
+		          }
 		        } else {
 		          _id = _huf->leftNode(id);
 		          _lid = -1;
+		        }
+		      }
 		      /// \brief Increment operator
 		      ///
 		      /// It steps to the next item in the class.
 		      ItemIt& operator++() {
 		        _id = _huf->nextNode(_id);
 		        return *this;
+		      }
 		      /// \brief Conversion operator
 		      ///
 		      /// It converts the iterator to the current item.
 		      operator const Item&() const {
 		        return _huf->nodes[_id].item;
+		      }
 		      /// \brief Equality operator
 		      ///
 		      /// Equality operator
 		      bool operator==(const ItemIt& i) {
 		        return i._id == _id;
+		      }
 		      /// \brief Inequality operator
 		      ///
 		      /// Inequality operator
 		      bool operator!=(const ItemIt& i) {
 		        return i._id != _id;
+		      }
 		      /// \brief Equality operator
 		      ///
 		      /// Equality operator
 		      bool operator==(Invalid) {
 		        return _id == _lid;
+		      }
 		      /// \brief Inequality operator
 		      ///
 		      /// Inequality operator
 		      bool operator!=(Invalid) {
 		        return _id != _lid;
+		      }
 		    };
 		    /// \brief Class iterator
 		    ///
 		    /// The iterator stores
 		    class ClassIt {
 		    private:
 		      const HeapUnionFind* _huf;
 		      int _id;
 		    public:
 		      ClassIt(const HeapUnionFind& huf)
 		        : _huf(&huf), _id(huf.first_class) {}
 		      ClassIt(const HeapUnionFind& huf, int cls)
 		        : _huf(&huf), _id(huf.classes[cls].left) {}
 		      ClassIt(Invalid) : _huf(0), _id(-1) {}
 		      const ClassIt& operator++() {
 		        _id = _huf->classes[_id].next;
 		        return *this;
+		      }
 		      /// \brief Equality operator
 		      ///
 		      /// Equality operator
 		      bool operator==(const ClassIt& i) {
 		        return i._id == _id;
+		      }
 		      /// \brief Inequality operator
 		      ///
 		      /// Inequality operator
 		      bool operator!=(const ClassIt& i) {
 		        return i._id != _id;
+		      }
 		      operator int() const {
 		        return _id;
+		      }
 		    };
 		  };
 		  //! @}
 		} //namespace lemon
 		#endif //LEMON_UNION_FIND_H

0 comments (0 inline)

RhodeCode

Login to your account