Update README wrt C++ run time

Split into separate files concatenable_queue.hpp and priority_queue.hpp.
Move classes into namespace.
Separate unit testse.

.gitignore

@@ -2,6 +2,8 @@ __pycache__/
 .coverage
 
 cpp/test_mwmatching
+cpp/test_concatenable_queue
+cpp/test_priority_queue
 cpp/run_matching
 cpp/run_matching_dbg
 tests/lemon/lemon_matching

README.md

@@ -42,10 +42,6 @@ print(matching)  # prints [(1, 4), (2, 3)]
 
 The folder [cpp/](cpp/) contains a header-only C++ implementation of maximum weighted matching.
 
-**NOTE:**
-The C++ code currently implements a slower algorithm that runs in _O(n<sup>3</sup>)_ steps.
-I plan to eventually update the C++ code to implement the faster _O(n m log n)_ algorithm.
-
 The C++ code is self-contained and can easily be linked into an application.
 It is also reasonably efficient.
 

cpp/Makefile

@@ -6,22 +6,34 @@ CXXFLAGS = -std=c++11 -Wall -Wextra $(OPTFLAGS) $(DBGFLAGS)
 LIB_BOOST_TEST = -l:libboost_unit_test_framework.a
 
 .PHONY: all
-all: run_matching run_matching_dbg test_mwmatching
+all: run_matching run_matching_dbg test_mwmatching test_concatenable_queue test_priority_queue
 
-run_matching: run_matching.cpp mwmatching.h
+run_matching: run_matching.cpp mwmatching.hpp concatenable_queue.hpp priority_queue.hpp
 	$(CXX) $(CPPFLAGS) $(CXXFLAGS) $(LDFLAGS) -o $@ $<
 
 run_matching_dbg: OPTFLAGS = -Og
 run_matching_dbg: DBGFLAGS = -g -fsanitize=address -fsanitize=undefined
-run_matching_dbg: run_matching.cpp mwmatching.h
+run_matching_dbg: run_matching.cpp mwmatching.hpp concatenable_queue.hpp priority_queue.hpp
 	$(CXX) $(CPPFLAGS) $(CXXFLAGS) $(LDFLAGS) -o $@ $<
 
 test_mwmatching: OPTFLAGS = -O1
 test_mwmatching: DBGFLAGS = -fsanitize=address -fsanitize=undefined
-test_mwmatching: test_mwmatching.cpp mwmatching.h
+test_mwmatching: test_mwmatching.cpp mwmatching.hpp concatenable_queue.hpp priority_queue.hpp
+	$(CXX) $(CPPFLAGS) $(CXXFLAGS) $(LDFLAGS) -o $@ $< $(LIB_BOOST_TEST)
+
+test_concatenable_queue: OPTFLAGS = -O1
+test_concatenable_queue: DBGFLAGS = -fsanitize=address -fsanitize=undefined
+test_concatenable_queue: test_concatenable_queue.cpp concatenable_queue.hpp
+	$(CXX) $(CPPFLAGS) $(CXXFLAGS) $(LDFLAGS) -o $@ $< $(LIB_BOOST_TEST)
+
+test_priority_queue: OPTFLAGS = -O1
+test_priority_queue: DBGFLAGS = -fsanitize=address -fsanitize=undefined
+test_priority_queue: test_priority_queue.cpp priority_queue.hpp
 	$(CXX) $(CPPFLAGS) $(CXXFLAGS) $(LDFLAGS) -o $@ $< $(LIB_BOOST_TEST)
 
 .PHONY: clean
 clean:
-	$(RM) run_matching run_matching_dbg test_mwmatching
+	$(RM) run_matching run_matching_dbg
+	$(RM) test_mwmatching
+	$(RM) test_concatenable_queue test_priority_queue
 

cpp/concatenable_queue.hpp

@@ -0,0 +1,841 @@
+/**
+ * Concatenable queue data structure.
+ */
+
+#ifndef MWMATCHING_CONCATENABLE_QUEUE_H_
+#define MWMATCHING_CONCATENABLE_QUEUE_H_
+
+#include <algorithm>
+#include <cassert>
+#include <tuple>
+
+
+namespace mwmatching {
+
+
+/**
+ * Priority queue supporting efficient merge and split operations.
+ *
+ * Conceptually, this is a combination of a disjoint set and a priority queue.
+ * Each queue has a "name".
+ * Each element has associated "data".
+ * Each element has a priority.
+ * Each element is contained in at most one queue at any time.
+ *
+ * The following operations can be done efficiently:
+ *  - Find the name of the queue that contains a given element.
+ *  - Change the priority of a given element.
+ *  - Find the element with lowest priority in a given queue.
+ *  - Merge two or more queues.
+ *  - Undo a previous merge step.
+ *
+ * A ConcatenableQueue instance may be destructed if it is empty and not
+ * currently merged into a larger queue. Alternatively, a group of related
+ * ConcatenableQueue instances and their Node instances may be destructed
+ * together, even if non-empty. In this case, the objects may be destructed
+ * in any order. No other interactions with the objects are allowed once
+ * destruction has started.
+ *
+ * This data structure is implemented as an AVL tree, with minimum-priority
+ * tracking added to it.
+ * See also
+ *   https://en.wikipedia.org/wiki/Avl_tree
+ * and
+ *   G. Blelloch, D. Ferizovic, Y. Sun, Parallel Ordered Sets Using Join,
+ *   https://arxiv.org/abs/1602.02120
+ */
+template <typename PrioType,
+          typename NameType,
+          typename DataType>
+class ConcatenableQueue
+{
+public:
+
+    /**
+     * A Node instance represents an element in a concatenable queue.
+     *
+     * A Node instance must remain valid while it is contained in a queue.
+     * The containing queue holds a pointer to the Node.
+     *
+     * A Node instance may be destructed if it is not contained in any queue.
+     * Alternatively, a Node instance may be destructed just before or after
+     * destructing the containing queue. In this case, no intermediate
+     * interactions with the node or queue are allowed.
+     */
+    class Node
+    {
+      public:
+        /** Construct an unused node, not yet contained in any queue. */
+        Node()
+          : owner_(nullptr)
+          , parent_(nullptr)
+          , left_child_(nullptr)
+          , right_child_(nullptr)
+          , height_(0)
+        { }
+
+        // Prevent copying.
+        Node(const Node&) = delete;
+        Node& operator=(const Node&) = delete;
+
+        /**
+         * Return the name of the queue that contains this element.
+         *
+         * The node must be contained in a queue.
+         * This function takes time O(log(n)).
+         */
+        NameType find() const
+        {
+            const Node* node = this;
+            while (node->parent_) {
+                node = node->parent_;
+            }
+            assert(node->owner_);
+            return node->owner_->name();
+        }
+
+        /**
+         * Return the priority of this element.
+         *
+         * The node must be contained in a queue.
+         * This function takes time O(1).
+         */
+        PrioType prio() const
+        {
+            assert(height_ != 0);
+            return prio_;
+        }
+
+        /**
+         * Change the priority of this element.
+         *
+         * The node must be contained in a queue.
+         * This function takes time O(log(n)).
+         */
+        void set_prio(PrioType new_prio)
+        {
+            assert(height_ != 0);
+
+            prio_ = new_prio;
+
+            Node* node = this;
+            while (node) {
+                PrioType min_prio = node->prio_;
+                DataType min_data = node->data_;
+                for (Node* child : { node->left_child_,
+                                     node->right_child_ }) {
+                    if (child && child->min_prio_ < min_prio) {
+                        min_prio = child->min_prio_;
+                        min_data = child->min_data_;
+                    }
+                }
+                node->min_prio_ = min_prio;
+                node->min_data_ = min_data;
+                node = node->parent_;
+            }
+        }
+
+      private:
+        PrioType            prio_;
+        DataType            data_;
+        PrioType            min_prio_;
+        DataType            min_data_;
+        ConcatenableQueue*  owner_;
+        Node*               parent_;
+        Node*               left_child_;
+        Node*               right_child_;
+        unsigned int        height_;
+
+        friend class ConcatenableQueue;
+    };
+
+    /** Construct an empty queue. */
+    explicit ConcatenableQueue(const NameType& name)
+      : name_(name)
+      , tree_(nullptr)
+      , first_node_(nullptr)
+      , first_subqueue_(nullptr)
+      , next_subqueue_(nullptr)
+    { }
+
+    // Prevent copying.
+    ConcatenableQueue(const ConcatenableQueue&) = delete;
+    ConcatenableQueue& operator=(const ConcatenableQueue&) = delete;
+
+    /** Return the name of this queue. */
+    NameType name() const
+    {
+        return name_;
+    }
+
+    /**
+     * Insert the specified Node into the queue.
+     *
+     * The Node instance must be unused (not yet contained in any queue).
+     *
+     * The queue must be empty. Only one element can be inserted into
+     * a queue in this way. Larger queues can only result from merging.
+     *
+     * The queue stores a pointer to the Node instance. The Node instance must
+     * remain valid for as long as it is contained in any queue.
+     *
+     * This function takes time O(1).
+     */
+    void insert(Node* node, PrioType prio, const DataType& data)
+    {
+        assert(! tree_);
+        assert(! first_node_);
+        assert(node->height_ == 0);
+        assert(! node->parent_);
+        assert(! node->left_child_);
+        assert(! node->right_child_);
+
+        node->prio_ = prio;
+        node->data_ = data;
+        node->min_prio_ = prio;
+        node->min_data_ = data;
+        node->owner_ = this;
+        node->height_ = 1;
+
+        tree_ = node;
+        first_node_ = node;
+    }
+
+    /**
+     * Return the minimum priority of any element in the queue.
+     *
+     * The queue must be non-empty.
+     * This function takes time O(1).
+     */
+    PrioType min_prio() const
+    {
+        assert(tree_);
+        return tree_->min_prio_;
+    }
+
+    /**
+     * Return the element with minimum priority.
+     *
+     * The queue must be non-empty.
+     * This function takes time O(1).
+     */
+    DataType min_elem() const
+    {
+        assert(tree_);
+        return tree_->min_data_;
+    }
+
+    /**
+     * Merge the specified queues.
+     *
+     * This queue instance must inititially be empty.
+     * All specified sub-queues must initially be non-empty.
+     *
+     * This function removes all elements from the specified sub-queues
+     * and adds them to this queue.
+     *
+     * After merging, this queue retains references to the sub-queues.
+     * These may be used later to split (undo the merge step).
+     *
+     * This function takes time O(n_sub_queues * log(n)).
+     */
+    template <typename QueueIterator>
+    void merge(QueueIterator sub_queues_begin, QueueIterator sub_queues_end)
+    {
+        assert(! tree_);
+        assert(! first_node_);
+        assert(sub_queues_begin != sub_queues_end);
+
+        // Pull the root node from the first sub-queue.
+        ConcatenableQueue* sub = *sub_queues_begin;
+        assert(sub->tree_);
+        Node* merged_tree = sub->tree_;
+        sub->tree_ = nullptr;
+
+        // Clear owner pointer from tree.
+        assert(merged_tree->owner_ == sub);
+        merged_tree->owner_ = nullptr;
+
+        // Copy first node to this queue.
+        assert(sub->first_node_);
+        first_node_ = sub->first_node_;
+
+        // Build linked list of sub-queues, starting with the first sub-queue.
+        assert(! sub->next_subqueue_);
+        first_subqueue_ = sub;
+
+        // Merge remaining sub-queues.
+        QueueIterator it = sub_queues_begin;
+        ++it;
+        while (it != sub_queues_end) {
+            ConcatenableQueue* prev_sub = sub;
+            sub = *it;
+
+            // Clear owner pointer and root node from the sub-queue.
+            assert(sub->tree_);
+            assert(sub->tree_->owner_ == sub);
+            sub->tree_->owner_ = nullptr;
+
+            // Merge sub-queue tree into our current tree.
+            merged_tree = merge_tree(merged_tree, sub->tree_);
+
+            // Clear root pointer from sub-queue.
+            sub->tree_ = nullptr;
+
+            // Add sub-queue to linked list of sub-queues.
+            assert(! sub->next_subqueue_);
+            prev_sub->next_subqueue_ = sub;
+
+            ++it;
+        }
+
+        // Put owner pointer in the root node of the tree.
+        merged_tree->owner_ = this;
+        tree_ = merged_tree;
+    }
+
+    /**
+     * Undo the merge step that created this queue.
+     *
+     * Remove all elements from this queue and put them back in
+     * the sub-queues from which they came.
+     *
+     * After splitting, this queue will be empty.
+     *
+     * This function takes time O(n_sub_queues * log(n)).
+     */
+    void split()
+    {
+        assert(tree_);
+        assert(first_subqueue_);
+
+        // Clear the owner pointer from the root node.
+        assert(tree_->owner_ == this);
+        tree_->owner_ = nullptr;
+
+        Node* rtree = tree_;
+        ConcatenableQueue* sub = first_subqueue_;
+
+        // Repeatedly split the tree to reconstruct each sub-queue.
+        while (sub->next_subqueue_) {
+            ConcatenableQueue* next_sub = sub->next_subqueue_;
+            sub->next_subqueue_ = nullptr;
+
+            // Split the tree on the first node of the next sub-queue.
+            assert(next_sub->first_node_);
+            Node* ltree;
+            std::tie(ltree, rtree) = split_tree(next_sub->first_node_);
+
+            // Put the left half of the tree in the current sub-queue.
+            sub->tree_ = ltree;
+            ltree->owner_ = sub;
+
+            // Continue to next sub-queue.
+            sub = next_sub;
+        }
+
+        // Put the remaining part of the tree in the last sub-queue.
+        rtree->owner_ = sub;
+        sub->tree_ = rtree;
+
+        // Clean up this instance.
+        tree_ = nullptr;
+        first_node_ = nullptr;
+        first_subqueue_ = nullptr;
+    }
+
+private:
+    /**
+     * Merge two trees and rebalance.
+     *
+     * This function takes time O(log(n)).
+     */
+    static Node* merge_tree(Node* ltree, Node* rtree)
+    {
+        // Remove the last node from the left tree.
+        Node* node;
+        std::tie(ltree, node) = split_last_node(ltree);
+
+        // Join the trees.
+        return join(ltree, node, rtree);
+    }
+
+    /**
+     * Remove the last node from the tree and rebalance.
+     *
+     * This function takes time O(log(n)).
+     */
+    static std::tuple<Node*, Node*> split_last_node(Node* tree)
+    {
+        assert(tree);
+
+        /*
+         * Descend down the right spine of the tree to find the last node.
+         *
+         *      tree
+         *      /  \
+         *          X
+         *         / \
+         *            X  <-- parent
+         *           / \
+         *             last_node
+         *              /
+         *             X  <-- ltree
+         */
+
+        // Find last node.
+        Node* last_node = tree;
+        while (last_node->right_child_) {
+            last_node = last_node->right_child_;
+        }
+
+        // Detach its left subtree.
+        Node* ltree = last_node->left_child_;
+        last_node->left_child_ = nullptr;
+        if (ltree) {
+            ltree->parent_ = nullptr;
+        }
+
+        // Detach from the parent.
+        Node* parent = last_node->parent_;
+        last_node->parent_ = nullptr;
+        if (parent) {
+            parent->right_child_ = nullptr;
+        }
+
+        // Reconstruct along the right spine of the original tree.
+        while (parent) {
+
+            // Step to parent.
+            Node* cur = parent;
+            parent = cur->parent_;
+
+            // Detach from its own parent.
+            cur->parent_ = nullptr;
+            if (parent) {
+                parent->right_child_ = nullptr;
+            }
+
+            // Join the current node to the reconstructed tree.
+            ltree = join(cur->left_child_, cur, ltree);
+        }
+
+        return std::make_tuple(ltree, last_node);
+    }
+
+    /**
+     * Split a tree on a specified node.
+     *
+     * Two new trees will be constructed.
+     * Leaf nodes to the left of "split_node" will go to the left tree.
+     * Leaf nodes to the right of "split_node", and "split_node" itself,
+     * will go to the right tree.
+     *
+     * This function takes time O(log(n)).
+     */
+    static std::tuple<Node*, Node*> split_tree(Node* split_node)
+    {
+        assert(split_node);
+
+        // All left-descendants of "split_node" belong in the left tree.
+        // Detach it from "split_node".
+        Node* ltree = split_node->left_child_;
+        split_node->left_child_ = nullptr;
+        if (ltree) {
+            ltree->parent_ = nullptr;
+        }
+
+        // Start with an empty right tree.
+        Node* rtree = nullptr;
+
+        // Start at "split_node".
+        // Take note of the fact that the "cut" between the two halves of
+        // the tree runs down its left branch, since "split_node" itself
+        // belongs to the right tree.
+        Node* node = split_node;
+        bool left_branch = true;
+
+        // Work upwards through the tree, assigning each node either to
+        // the new left tree or the new right tree.
+        //
+        // This loop runs for O(log(n)) iterations.
+        // Each iteration calls join() once, taking time proportional
+        // to the difference in height between the intermediate trees.
+        // The total run time of all join() calls together is O(log(n)).
+        while (node) {
+
+            // Detach the current node from its parent.
+            // Remember to which branch of the parent it was attached.
+            Node* parent = node->parent_;
+            node->parent_ = nullptr;
+            bool parent_left_branch = true;
+            if (parent) {
+                parent_left_branch = (parent->left_child_ == node);
+                if (parent_left_branch) {
+                    parent->left_child_ = nullptr;
+                } else {
+                    assert(parent->right_child_ == node);
+                    parent->right_child_ = nullptr;
+                }
+            }
+
+            // Join the current node and its remaining descendents either
+            // to the left tree or to the right tree.
+            if (left_branch) {
+                /*
+                 * "node" belongs to the right tree.
+                 * Join like this:
+                 *
+                 *             (node)  <--- new rtree
+                 *             /    \
+                 *        (rtree)  (node->right_child)
+                 */
+                assert(! node->left_child_);
+                rtree = join(rtree, node, node->right_child_);
+            } else {
+                /*
+                 * "node" belongs to the left tree.
+                 * Join like this:
+                 *
+                 *                   (node)  <--- new ltree
+                 *                   /    \
+                 *   (node->left_child)   (ltree)
+                 */
+                assert(! node->right_child_);
+                ltree = join(node->left_child_, node, ltree);
+            }
+
+            // Continue with the parent node.
+            node = parent;
+            left_branch = parent_left_branch;
+        }
+
+        // Done. All that remains of the old tree is the two new halves.
+        return std::make_tuple(ltree, rtree);
+    }
+
+    /** Return node height, or 0 if node == nullptr. */
+    static unsigned int get_node_height(const Node* node)
+    {
+        return node ? node->height_ : 0;
+    }
+
+    /**
+     * Repair the height and min-priority information of a node
+     * after modifying its children.
+     *
+     * After repairing a node, it is typically necessary to also repair
+     * its ancestors.
+     */
+    static void repair_node(Node* node)
+    {
+        Node* lchild = node->left_child_;
+        Node* rchild = node->right_child_;
+
+        // Repair node height.
+        node->height_ = 1 + std::max(get_node_height(lchild),
+                                     get_node_height(rchild));
+
+        // Repair min-priority.
+        PrioType min_prio = node->prio_;
+        DataType min_data = node->data_;
+        for (Node* child : { lchild, rchild }) {
+            if (child && child->min_prio_ < min_prio) {
+                min_prio = child->min_prio_;
+                min_data = child->min_data_;
+            }
+        }
+        node->min_prio_ = min_prio;
+        node->min_data_ = min_data;
+    }
+
+    /** Rotate the subtree to the left and return the new root of the subtree. */
+    static Node* rotate_left(Node* node)
+    {
+        /*
+         *     N                C
+         *    / \              / \
+         *   A   C    --->    N   D
+         *      / \          / \
+         *     B   D        A   B
+         */
+        Node* parent = node->parent_;
+        Node* new_top = node->right_child_;
+        assert(new_top);
+
+        Node* nb = new_top->left_child_;
+        node->right_child_ = nb;
+        if (nb) {
+            nb->parent_ = node;
+        }
+
+        new_top->left_child_ = node;
+        node->parent_ = new_top;
+
+        new_top->parent_ = parent;
+
+        if (parent) {
+            if (parent->left_child_ == node) {
+                parent->left_child_ = new_top;
+            } else {
+                assert(parent->right_child_ == node);
+                parent->right_child_ = new_top;
+            }
+        }
+
+        repair_node(node);
+        repair_node(new_top);
+
+        return new_top;
+    }
+
+    /** Rotate the subtree to the right and return the new root of the subtree. */
+    static Node* rotate_right(Node* node)
+    {
+        /*
+         *       N                B
+         *      / \              / \
+         *     B   D    --->    A   N
+         *    / \                  / \
+         *   A   C                C   D
+         */
+        Node* parent = node->parent_;
+        Node* new_top = node->left_child_;
+        assert(new_top);
+
+        Node* nc = new_top->right_child_;
+        node->left_child_ = nc;
+        if (nc) {
+            nc->parent_ = node;
+        }
+
+        new_top->right_child_ = node;
+        node->parent_ = new_top;
+
+        new_top->parent_ = parent;
+
+        if (parent) {
+            if (parent->left_child_ == node) {
+                parent->left_child_ = new_top;
+            } else {
+                assert(parent->right_child_ == node);
+                parent->right_child_ = new_top;
+            }
+        }
+
+        repair_node(node);
+        repair_node(new_top);
+
+        return new_top;
+    }
+
+    /**
+     * Join a left subtree, middle node and right subtree together.
+     *
+     * The left subtree is higher than the right subtree.
+     */
+    static Node* join_right(Node* ltree, Node* node, Node* rtree)
+    {
+        assert(ltree);
+        unsigned int lh = ltree->height_;
+        unsigned int rh = get_node_height(rtree);
+        assert(lh > rh + 1);
+
+        /*
+         * Descend down the right spine of "ltree".
+         * Stop at a node with compatible height, then insert "node"
+         * and attach "rtree".
+         *
+         *     ltree
+         *      /  \
+         *          X
+         *         / \
+         *            X  <-- cur
+         *           / \
+         *             node
+         *              / \
+         *             X   rtree
+         */
+
+        // Descend to a point with compatible height.
+        Node* cur = ltree;
+        while (cur->right_child_ && (cur->right_child_->height_ > rh + 1)) {
+            cur = cur->right_child_;
+        }
+
+        // Insert "node" and "rtree".
+        node->left_child_ = cur->right_child_;
+        node->right_child_ = rtree;
+        if (node->left_child_) {
+            node->left_child_->parent_ = node;
+        }
+        if (rtree) {
+            rtree->parent_ = node;
+        }
+        cur->right_child_ = node;
+        node->parent_ = cur;
+
+        // A double rotation may be necessary.
+        if ((! cur->left_child_) || (cur->left_child_->height_ <= rh)) {
+            node = rotate_right(node);
+            cur = rotate_left(cur);
+        } else {
+            repair_node(node);
+            repair_node(cur);
+        }
+
+        // Ascend from "cur" to the root of the tree; repair and rebalance.
+        while (cur->parent_) {
+            cur = cur->parent_;
+            assert(cur->left_child_);
+            assert(cur->right_child_);
+
+            if (cur->left_child_->height_ + 1 < cur->right_child_->height_) {
+                cur = rotate_left(cur);
+            } else {
+                repair_node(cur);
+            }
+        }
+
+        return cur;
+    }
+
+    /**
+     * Join a left subtree, middle node and right subtree together.
+     *
+     * The right subtree is higher than the left subtree.
+     */
+    static Node* join_left(Node* ltree, Node* node, Node* rtree)
+    {
+        assert(rtree);
+        unsigned int lh = get_node_height(ltree);
+        unsigned int rh = rtree->height_;
+        assert(lh + 1 < rh);
+
+        /*
+         * Descend down the left spine of "rtree".
+         * Stop at a node with compatible height, then insert "node"
+         * and attach "ltree".
+         *
+         *             rtree
+         *             /  \
+         *            X
+         *           / \
+         * cur -->  X
+         *         / \
+         *      node
+         *       / \
+         *  ltree   X
+         */
+
+        // Descend to a point with compatible height.
+        Node* cur = rtree;
+        while (cur->left_child_ && (cur->left_child_->height_ > lh + 1)) {
+            cur = cur->left_child_;
+        }
+
+        // Insert "node" and "ltree".
+        node->left_child_ = ltree;
+        node->right_child_ = cur->left_child_;
+        if (ltree) {
+            ltree->parent_ = node;
+        }
+        if (node->right_child_) {
+            node->right_child_->parent_ = node;
+        }
+        cur->left_child_ = node;
+        node->parent_ = cur;
+
+        // A double rotation may be necessary.
+        if ((! cur->right_child_) || (cur->right_child_->height_ <= lh)) {
+            node = rotate_left(node);
+            cur = rotate_right(cur);
+        } else {
+            repair_node(node);
+            repair_node(cur);
+        }
+
+        // Ascend from "cur" to the root of the tree; repair and rebalance.
+        while (cur->parent_) {
+            cur = cur->parent_;
+            assert(cur->left_child_);
+            assert(cur->right_child_);
+
+            if (cur->left_child_->height_ > cur->right_child_->height_ + 1) {
+                cur = rotate_right(cur);
+            } else {
+                repair_node(cur);
+            }
+        }
+
+        return cur;
+    }
+
+    /**
+     * Join a left subtree, middle node and right subtree together.
+     *
+     * The left or right subtree may initially be a child of the middle
+     * node; such links will be broken as needed.
+     *
+     * The left and right subtrees must be consistent, semi-balanced trees.
+     * Height and priority of the middle node may initially be inconsistent;
+     * this function will repair it.
+     *
+     * @return Root node of the joined tree.
+     */
+    static Node* join(Node* ltree, Node* node, Node* rtree)
+    {
+        unsigned int lh = get_node_height(ltree);
+        unsigned int rh = get_node_height(rtree);
+
+        if (lh > rh + 1) {
+            assert(ltree);
+            ltree->parent_ = nullptr;
+            return join_right(ltree, node, rtree);
+        } else if (lh + 1 < rh) {
+            assert(rtree);
+            rtree->parent_ = nullptr;
+            return join_left(ltree, node, rtree);
+        } else {
+            /*
+             * Subtree heights are compatible. Just join them.
+             *
+             *       node
+             *       /  \
+             *   ltree  rtree
+             *    / \    / \
+             */
+            node->parent_ = nullptr;
+            node->left_child_ = ltree;
+            if (ltree) {
+                ltree->parent_ = node;
+            }
+            node->right_child_ = rtree;
+            if (rtree) {
+                rtree->parent_ = node;
+            }
+            repair_node(node);
+            return node;
+        }
+    }
+
+    /** Name of this queue. */
+    const NameType name_;
+
+    /** Root node of the tree, or "nullptr" if the queue is empty. */
+    Node* tree_;
+
+    /** Left-most node that belongs to this queue. */
+    Node* first_node_;
+
+    /** Pointer to first sub-queue that got merged into this instance. */
+    ConcatenableQueue* first_subqueue_;
+
+    /** Linked list of sub-queues that were merged to build this queue. */
+    ConcatenableQueue* next_subqueue_;
+};
+
+
+} // namespace mwmatching
+
+#endif  // MWMATCHING_CONCATENABLE_QUEUE_H_

cpp/priority_queue.hpp

@@ -0,0 +1,266 @@
+/**
+ * Priority queue data structure.
+ */
+
+#ifndef MWMATCHING_PRIORITY_QUEUE_H_
+#define MWMATCHING_PRIORITY_QUEUE_H_
+
+#include <cassert>
+#include <limits>
+#include <vector>
+
+
+namespace mwmatching {
+
+
+/**
+ * Min-priority queue implemented as a binary heap.
+ *
+ * Elements in a heap have a priority and associated "data".
+ *
+ * The following operations can be done efficiently:
+ *  - Insert an element into the queue.
+ *  - Remove an element from the queue.
+ *  - Change the priority of a given element.
+ *  - Find the element with lowest priority in the queue.
+ */
+template <typename PrioType,
+          typename DataType>
+class PriorityQueue
+{
+public:
+    typedef unsigned int IndexType;
+    static constexpr IndexType INVALID_INDEX =
+        std::numeric_limits<IndexType>::max();
+
+    /**
+     * A Node instance represents an element in a PriorityQueue.
+     *
+     * A Node instance must remain valid while it is contained in a queue.
+     * The containing queue holds a pointer to the Node.
+     *
+     * A Node instance may be destructed if it is not contained in any queue.
+     * Alternatively, a Node instance may be destructed after its containing
+     * queue instance has been destructed.
+     */
+    class Node
+    {
+      public:
+        /** Construct an invalid node, not contained in any queue. */
+        Node()
+          : index_(INVALID_INDEX)
+        { }
+
+        // Prevent copying.
+        Node(const Node&) = delete;
+        Node& operator=(const Node&) = delete;
+
+        /** Return true if this node is contained in a queue. */
+        bool valid() const
+        {
+            return (index_ != INVALID_INDEX);
+        }
+
+        /**
+         * Return the priority of this node in the queue.
+         *
+         * The node must be contained in a queue.
+         * This function takes time O(1).
+         */
+        PrioType prio() const
+        {
+            assert(index_ != INVALID_INDEX);
+            return prio_;
+        }
+
+      private:
+        IndexType       index_;
+        PrioType        prio_;
+        DataType        data_;
+
+        friend class PriorityQueue;
+    };
+
+    /** Construct an empty queue. */
+    PriorityQueue()
+    { }
+
+    // Prevent copying.
+    PriorityQueue(const PriorityQueue&) = delete;
+    PriorityQueue& operator=(const PriorityQueue&) = delete;
+
+    /**
+     * Remove all elements from the queue.
+     *
+     * This function takes time O(n).
+     */
+    void clear()
+    {
+        for (Node* node : heap_) {
+            node->index_ = INVALID_INDEX;
+        }
+        heap_.clear();
+    }
+
+    /** Return true if the queue is empty. */
+    bool empty() const
+    {
+        return heap_.empty();
+    }
+
+    /**
+     * Return the minimum priority of any element in the queue.
+     *
+     * The queue must be non-empty.
+     * This function takes time O(1).
+     */
+    PrioType min_prio() const
+    {
+        assert(! heap_.empty());
+        Node* top = heap_.front();
+        return top->prio_;
+    }
+
+    /**
+     * Return the element with minimum priority.
+     *
+     * The queue must be non-empty.
+     * This function takes time O(1).
+     */
+    DataType min_elem() const
+    {
+        assert(! heap_.empty());
+        Node* top = heap_.front();
+        return top->data_;
+    }
+
+    /**
+     * Insert the given node into the queue with associated data.
+     *
+     * The node must not currently be contained in any queue.
+     * This function takes amortized time O(log(n)).
+     */
+    void insert(Node* node, PrioType prio, const DataType& data)
+    {
+        assert(node->index_ == INVALID_INDEX);
+
+        node->index_ = heap_.size();
+        node->prio_ = prio;
+        node->data_ = data;
+
+        heap_.push_back(node);
+        sift_up(node->index_);
+    }
+
+    /**
+     * Update priority of an existing node.
+     *
+     * This function takes time O(log(n)).
+     */
+    void set_prio(Node* node, PrioType prio)
+    {
+        IndexType index = node->index_;
+        assert(index != INVALID_INDEX);
+        assert(heap_[index] == node);
+
+        PrioType prev_prio = node->prio_;
+        node->prio_ = prio;
+
+        if (prio < prev_prio) {
+            sift_up(index);
+        } else if (prio > prev_prio) {
+            sift_down(index);
+        }
+    }
+
+    /**
+     * Remove the specified element from the queue.
+     *
+     * This function takes time O(log(n)).
+     */
+    void remove(Node* node)
+    {
+        IndexType index = node->index_;
+        assert(index != INVALID_INDEX);
+        assert(heap_[index] == node);
+
+        node->index_ = INVALID_INDEX;
+
+        Node* move_node = heap_.back();
+        heap_.pop_back();
+
+        if (index < heap_.size()) {
+            heap_[index] = move_node;
+            move_node->index_ = index;
+            if (move_node->prio_ < node->prio_) {
+                sift_up(index);
+            } else if (move_node->prio_ > node->prio_) {
+                sift_down(index);
+            }
+        }
+    }
+
+private:
+    /** Repair the heap along an ascending path to the root. */
+    void sift_up(IndexType index)
+    {
+        Node* node = heap_[index];
+        PrioType prio = node->prio_;
+
+        while (index > 0) {
+            IndexType next_index = (index - 1) / 2;
+            Node* next_node = heap_[next_index];
+            if (next_node->prio_ <= prio) {
+                break;
+            }
+            heap_[index] = next_node;
+            next_node->index_ = index;
+            index = next_index;
+        }
+
+        node->index_ = index;
+        heap_[index] = node;
+    }
+
+    /** Repair the heap along a descending path. */
+    void sift_down(IndexType index)
+    {
+        Node* node = heap_[index];
+        PrioType prio = node->prio_;
+
+        IndexType num_elem = heap_.size();
+
+        while (index < num_elem / 2) {
+            IndexType next_index = 2 * index + 1;
+            Node* next_node = heap_[next_index];
+
+            if (next_index + 1 < num_elem) {
+                Node* tmp_node = heap_[next_index + 1];
+                if (tmp_node->prio_ <= next_node->prio_) {
+                    ++next_index;
+                    next_node = tmp_node;
+                }
+            }
+
+            if (next_node->prio_ >= prio) {
+                break;
+            }
+
+            heap_[index] = next_node;
+            next_node->index_ = index;
+
+            index = next_index;
+        }
+
+        heap_[index] = node;
+        node->index_ = index;
+    }
+
+    /** Heap data structure. */
+    std::vector<Node*> heap_;
+};
+
+
+} // namespace mwmatching
+
+#endif  // MWMATCHING_PRIORITY_QUEUE_H_

cpp/run_matching.cpp

@@ -13,7 +13,7 @@
 #include <utility>
 #include <vector>
 
-#include "mwmatching.h"
+#include "mwmatching.hpp"
 
 
 namespace {  // anonymous namespace

cpp/test_concatenable_queue.cpp

@@ -0,0 +1,403 @@
+/*
+ * Unit tests for ConcatenableQueue data structure.
+ *
+ * Depends on the Boost.Test unit test framework.
+ * Tested with Boost v1.74, available from https://www.boost.org/
+ */
+
+#include <algorithm>
+#include <memory>
+#include <random>
+#include <set>
+#include <string>
+#include <unordered_map>
+#include <vector>
+
+#define BOOST_TEST_MODULE datastruct
+#include <boost/test/unit_test.hpp>
+
+#include "concatenable_queue.hpp"
+
+
+BOOST_AUTO_TEST_SUITE(test_concatenable_queue)
+
+BOOST_AUTO_TEST_CASE(test_single)
+{
+    using Queue = mwmatching::ConcatenableQueue<int, std::string, std::string>;
+    Queue q("Q");
+    BOOST_TEST(q.name() == std::string("Q"));
+
+    Queue::Node n;
+    q.insert(&n, 4, "a");
+
+    BOOST_TEST(n.find() == std::string("Q"));
+    BOOST_TEST(n.prio() == 4);
+    BOOST_TEST(q.min_prio() == 4);
+    BOOST_TEST(q.min_elem() == std::string("a"));
+
+    n.set_prio(8);
+
+    BOOST_TEST(n.prio() == 8);
+    BOOST_TEST(n.find() == std::string("Q"));
+    BOOST_TEST(q.min_prio() == 8);
+    BOOST_TEST(q.min_elem() == std::string("a"));
+}
+
+BOOST_AUTO_TEST_CASE(test_simple)
+{
+    using Queue = mwmatching::ConcatenableQueue<int, std::string, char>;
+    Queue::Node n1, n2, n3, n4, n5;
+
+    Queue q1("A");
+    q1.insert(&n1, 5, 'a');
+
+    Queue q2("B");
+    q2.insert(&n2, 6, 'b');
+
+    Queue q3("C");
+    q3.insert(&n3, 7, 'c');
+
+    Queue q4("D");
+    q4.insert(&n4, 4, 'd');
+
+    Queue q5("E");
+    q5.insert(&n5, 3, 'e');
+
+    auto m345 = {&q3, &q4, &q5};
+    Queue q345("P");
+    q345.merge(m345.begin(), m345.end());
+
+    BOOST_TEST(n1.find() == std::string("A"));
+    BOOST_TEST(n2.find() == std::string("B"));
+    BOOST_TEST(n3.find() == std::string("P"));
+    BOOST_TEST(n4.find() == std::string("P"));
+    BOOST_TEST(n5.find() == std::string("P"));
+    BOOST_TEST(q345.min_prio() == 3);
+    BOOST_TEST(q345.min_elem() == 'e');
+
+    n5.set_prio(6);
+    BOOST_TEST(q345.min_prio() == 4);
+    BOOST_TEST(q345.min_elem() == 'd');
+
+    auto m12 = {&q1, &q2};
+    Queue q12("Q");
+    q12.merge(m12.begin(), m12.end());
+
+    BOOST_TEST(n1.find() == std::string("Q"));
+    BOOST_TEST(n2.find() == std::string("Q"));
+    BOOST_TEST(q12.min_prio() == 5);
+    BOOST_TEST(q12.min_elem() == 'a');
+
+    auto m12345 = {&q12, &q345};
+    Queue q12345("R");
+    q12345.merge(m12345.begin(), m12345.end());
+
+    BOOST_TEST(n1.find() == std::string("R"));
+    BOOST_TEST(n2.find() == std::string("R"));
+    BOOST_TEST(n3.find() == std::string("R"));
+    BOOST_TEST(n4.find() == std::string("R"));
+    BOOST_TEST(n5.find() == std::string("R"));
+    BOOST_TEST(q12345.min_prio() == 4);
+    BOOST_TEST(q12345.min_elem() == 'd');
+
+    n4.set_prio(8);
+
+    BOOST_TEST(q12345.min_prio() == 5);
+    BOOST_TEST(q12345.min_elem() == 'a');
+
+    n3.set_prio(2);
+
+    BOOST_TEST(q12345.min_prio() == 2);
+    BOOST_TEST(q12345.min_elem() == 'c');
+
+    q12345.split();
+
+    BOOST_TEST(n1.find() == std::string("Q"));
+    BOOST_TEST(n2.find() == std::string("Q"));
+    BOOST_TEST(n3.find() == std::string("P"));
+    BOOST_TEST(n4.find() == std::string("P"));
+    BOOST_TEST(n5.find() == std::string("P"));
+    BOOST_TEST(q12.min_prio() == 5);
+    BOOST_TEST(q12.min_elem() == 'a');
+    BOOST_TEST(q345.min_prio() == 2);
+    BOOST_TEST(q345.min_elem() == 'c');
+
+    q12.split();
+    BOOST_TEST(n1.find() == std::string("A"));
+    BOOST_TEST(n2.find() == std::string("B"));
+
+    q345.split();
+    BOOST_TEST(n3.find() == std::string("C"));
+    BOOST_TEST(n4.find() == std::string("D"));
+    BOOST_TEST(n5.find() == std::string("E"));
+    BOOST_TEST(q3.min_prio() == 2);
+    BOOST_TEST(q3.min_elem() == 'c');
+}
+
+BOOST_AUTO_TEST_CASE(test_medium)
+{
+    using Queue = mwmatching::ConcatenableQueue<int, char, char>;
+
+    std::vector<char> queue_names(19);
+    for (int i = 0; i < 14; i++) {
+        queue_names[i] = 'A' + i;
+    }
+    queue_names[14] = 'P';
+    queue_names[15] = 'Q';
+    queue_names[16] = 'R';
+    queue_names[17] = 'S';
+    queue_names[18] = 'Z';
+
+    std::vector<Queue> queues(queue_names.begin(), queue_names.end());
+    Queue::Node nodes[14];
+
+    int prios[14] = {3, 8, 6, 2, 9, 4, 6, 8, 1, 5, 9, 4, 7, 8};
+
+    auto min_prio = [&prios](int begin, int end) -> int {
+        int p = begin;
+        int m = prios[p];
+        ++p;
+        while (p != end) {
+            m = std::min(m, prios[p]);
+            ++p;
+        }
+        return m;
+    };
+
+    for (int i = 0; i < 14; i++) {
+        BOOST_TEST(queues[i].name() == 'A' + i);
+        queues[i].insert(&nodes[i], prios[i], 'a' + i);
+    }
+
+    auto m14 = {&queues[0], &queues[1]};
+    queues[14].merge(m14.begin(), m14.end());
+    BOOST_TEST(queues[14].name() == 'P');
+    BOOST_TEST(queues[14].min_prio() == min_prio(0, 2));
+
+    auto m15 = {&queues[2], &queues[3], &queues[4]};
+    queues[15].merge(m15.begin(), m15.end());
+    BOOST_TEST(queues[15].name() == 'Q');
+    BOOST_TEST(queues[15].min_prio() == min_prio(2, 5));
+
+    auto m16 = {&queues[5], &queues[6], &queues[7], &queues[8]};
+    queues[16].merge(m16.begin(), m16.end());
+    BOOST_TEST(queues[16].name() == 'R');
+    BOOST_TEST(queues[16].min_prio() == min_prio(5, 9));
+
+    auto m17 = {&queues[9], &queues[10], &queues[11], &queues[12],
+                &queues[13]};
+    queues[17].merge(m17.begin(), m17.end());
+    BOOST_TEST(queues[17].name() == 'S');
+    BOOST_TEST(queues[17].min_prio() == min_prio(9, 14));
+
+    for (int i = 0; i < 2; i++) {
+        BOOST_TEST(nodes[i].find() == 'P');
+    }
+    for (int i = 2; i < 5; i++) {
+        BOOST_TEST(nodes[i].find() == 'Q');
+    }
+    for (int i = 5; i < 9; i++) {
+        BOOST_TEST(nodes[i].find() == 'R');
+    }
+    for (int i = 9; i < 14; i++) {
+        BOOST_TEST(nodes[i].find() == 'S');
+    }
+
+    auto m18 = {&queues[14], &queues[15], &queues[16], &queues[17]};
+    queues[18].merge(m18.begin(), m18.end());
+    BOOST_TEST(queues[18].name() == 'Z');
+    BOOST_TEST(queues[18].min_prio() == 1);
+    BOOST_TEST(queues[18].min_elem() == 'i');
+
+    for (int i = 0; i < 14; i++) {
+        BOOST_TEST(nodes[i].find() == 'Z');
+    }
+
+    prios[8] = 5;
+    nodes[8].set_prio(prios[8]);
+    BOOST_TEST(queues[18].min_prio() == 2);
+    BOOST_TEST(queues[18].min_elem() == 'd');
+
+    queues[18].split();
+
+    for (int i = 0; i < 2; i++) {
+        BOOST_TEST(nodes[i].find() == 'P');
+    }
+    for (int i = 2; i < 5; i++) {
+        BOOST_TEST(nodes[i].find() == 'Q');
+    }
+    for (int i = 5; i < 9; i++) {
+        BOOST_TEST(nodes[i].find() == 'R');
+    }
+    for (int i = 9; i < 14; i++) {
+        BOOST_TEST(nodes[i].find() == 'S');
+    }
+
+    BOOST_TEST(queues[14].min_prio() == min_prio(0, 2));
+    BOOST_TEST(queues[15].min_prio() == min_prio(2, 5));
+    BOOST_TEST(queues[16].min_prio() == min_prio(5, 9));
+    BOOST_TEST(queues[17].min_prio() == min_prio(9, 14));
+
+    for (int i = 14; i < 18; i++) {
+        queues[i].split();
+    }
+
+    for (int i = 0; i < 14; i++) {
+        BOOST_TEST(nodes[i].find() == 'A' + i);
+        BOOST_TEST(queues[i].min_prio() == prios[i]);
+        BOOST_TEST(queues[i].min_elem() == 'a' + i);
+    }
+}
+
+BOOST_AUTO_TEST_CASE(test_random)
+{
+    using Queue = mwmatching::ConcatenableQueue<double, int, int>;
+
+    constexpr int NUM_NODES = 4000;
+
+    std::mt19937 rng(23456);
+    std::uniform_real_distribution<> prio_distribution;
+    std::uniform_int_distribution<> node_distribution(0, NUM_NODES - 1);
+
+    double prios[NUM_NODES];
+    Queue::Node nodes[NUM_NODES];
+    std::unordered_map<int, std::unique_ptr<Queue>> queues;
+    std::unordered_map<int, std::set<int>> queue_nodes;
+    std::unordered_map<int, std::set<int>> queue_subs;
+    std::set<int> live_queues;
+    std::set<int> live_merged_queues;
+
+    // Make trivial queues.
+    for (int i = 0; i < NUM_NODES; i++) {
+        int name = 10000 + i;
+        prios[i] = prio_distribution(rng);
+        queues[name] = std::unique_ptr<Queue>(new Queue(name));
+        queues[name]->insert(&nodes[i], prios[i], i);
+        queue_nodes[name].insert(i);
+        live_queues.insert(name);
+    }
+
+    // Run modifications.
+    for (int i = 0; i < 4000; i++) {
+
+        // Find top-level queue of few nodes.
+        for (int k = 0; k < 200; k++) {
+            int t = node_distribution(rng);
+            int name = nodes[t].find();
+            BOOST_TEST(queue_nodes[name].count(t) == 1);
+        }
+
+        // Change priority of a few nodes.
+        for (int k = 0; k < 10; k++) {
+            int t = node_distribution(rng);
+            int name = nodes[t].find();
+            BOOST_TEST(live_queues.count(name) == 1);
+            BOOST_TEST(queue_nodes[name].count(t) == 1);
+            BOOST_TEST(nodes[t].prio() == prios[t]);
+            double p = prio_distribution(rng);
+            prios[t] = p;
+            nodes[t].set_prio(p);
+            for (int tt : queue_nodes[name]) {
+                if (prios[tt] < p) {
+                    t = tt;
+                    p = prios[tt];
+                }
+            }
+            BOOST_TEST(queues[name]->min_prio() == p);
+            BOOST_TEST(queues[name]->min_elem() == t);
+        }
+
+        if (live_queues.size() > 100) {
+
+            // Choose number of queues to merge between 2 and 100.
+            int k = std::uniform_int_distribution<>(2, 100)(rng);
+            k = std::uniform_int_distribution<>(2, k)(rng);
+
+            // Choose queues to merge.
+            std::vector<int> live_queue_vec(live_queues.begin(),
+                                            live_queues.end());
+            std::vector<int> sub_names;
+            std::vector<Queue*> sub_queues;
+            for (int ki = 0; ki < k; ki++) {
+                int t = std::uniform_int_distribution<>(
+                    0, live_queue_vec.size() - 1)(rng);
+                int name = live_queue_vec[t];
+                sub_names.push_back(name);
+                sub_queues.push_back(queues[name].get());
+                live_queue_vec[t] = live_queue_vec.back();
+                live_queue_vec.pop_back();
+                live_queues.erase(name);
+                live_merged_queues.erase(name);
+            }
+
+            // Create new queue by merging selected queues.
+            int name = 20000 + i;
+            queues[name] = std::unique_ptr<Queue>(new Queue(name));
+            queues[name]->merge(sub_queues.begin(), sub_queues.end());
+            for (int nn : sub_names) {
+                queue_nodes[name].insert(queue_nodes[nn].begin(),
+                                         queue_nodes[nn].end());
+            }
+            queue_subs[name].insert(sub_names.begin(), sub_names.end());
+            live_queues.insert(name);
+            live_merged_queues.insert(name);
+
+            // Check new queue.
+            {
+                double p = 2;
+                int t = 0;
+                for (int tt : queue_nodes[name]) {
+                    if (prios[tt] < p) {
+                        t = tt;
+                        p = prios[tt];
+                    }
+                }
+                BOOST_TEST(queues[name]->min_prio() == p);
+                BOOST_TEST(queues[name]->min_elem() == t);
+            }
+        }
+
+        if ((live_queues.size() <= 100)
+                || (live_merged_queues.size() >= 100)) {
+
+            // Choose a random queue to split.
+            std::vector<int> live_queue_vec(live_merged_queues.begin(),
+                                            live_merged_queues.end());
+            int k = std::uniform_int_distribution<>(
+                0, live_queue_vec.size() - 1)(rng);
+            int name = live_queue_vec[k];
+
+            queues[name]->split();
+
+            for (int nn : queue_subs[name]) {
+                // Check reconstructed sub-queue.
+                double p = 2;
+                int t = 0;
+                for (int tt : queue_nodes[nn]) {
+                    if (prios[tt] < p) {
+                        t = tt;
+                        p = prios[tt];
+                    }
+                }
+                BOOST_TEST(queues[nn]->min_prio() == p);
+                BOOST_TEST(queues[nn]->min_elem() == t);
+
+                // Mark sub-queue as live.
+                live_queues.insert(nn);
+                if (queue_subs.count(nn) == 1) {
+                    live_merged_queues.insert(nn);
+                }
+            }
+
+            live_merged_queues.erase(name);
+            live_queues.erase(name);
+
+            queues.erase(name);
+            queue_nodes.erase(name);
+            queue_subs.erase(name);
+        }
+
+    }
+}
+
+BOOST_AUTO_TEST_SUITE_END()

cpp/test_mwmatching.cpp

@@ -13,7 +13,7 @@
 #define BOOST_TEST_MODULE mwmatching
 #include <boost/test/unit_test.hpp>
 
-#include "mwmatching.h"
+#include "mwmatching.hpp"
 
 
 using EdgeVectorLong = std::vector<mwmatching::Edge<long>>;
@@ -327,19 +327,40 @@ BOOST_AUTO_TEST_CASE(test51_augment_blossom_nested)
     BOOST_TEST(mwmatching::maximum_weight_matching(edges) == expect);
 }
 
-BOOST_AUTO_TEST_CASE(test61_unsigned)
+BOOST_AUTO_TEST_CASE(test52_augment_blossom_nested2)
 {
-    std::vector<mwmatching::Edge<unsigned int>> edges = {{1, 2, 5}, {2, 3, 11}, {3, 4, 5}};
-    Matching expect = {{2, 3}};
+    /*
+     *
+     *   [4]--15     19--[2]
+     *    |     \   /     |
+     *   16      [1]     21
+     *    |     /   \     |
+     *   [5]--17     19--[3]
+     *    |
+     *   10
+     *    |
+     *   [6]
+     */
+    EdgeVectorLong edges = {
+        {1, 2, 19}, {1, 3, 19}, {1, 4, 15}, {1, 5, 17}, {2, 3, 21}, {4, 5, 16}, {5, 6, 10}};
+    Matching expect = {{1, 4}, {2, 3}, {5, 6}};
     BOOST_TEST(mwmatching::maximum_weight_matching(edges) == expect);
 }
 
-BOOST_AUTO_TEST_CASE(test62_unsigned2)
+BOOST_AUTO_TEST_CASE(test61_triangles_n9)
 {
-    std::vector<mwmatching::Edge<unsigned int>> edges = {
-        {1, 2, 8}, {1, 3, 9}, {2, 3, 10}, {3, 4, 7}, {1, 6, 5}, {4, 5, 6}};
-    Matching expect = {{2, 3}, {1, 6}, {4, 5}};
-    BOOST_TEST(mwmatching::maximum_weight_matching(edges) == expect);
+    /*
+     * t.f 9 nodes
+     *
+     *  [1]------[4]      [7]
+     *   | \      | \      | \
+     *   |  [3]   |  [6]   |  [9]
+     *   | /      | /      | /
+     *  [2]      [5]------[8]
+     */
+    EdgeVectorLong edges = {
+        {1, 2, 1}, {1, 3, 1}, {2, 3, 1}, {4, 5, 1}, {4, 6, 1}, {5, 6, 1}, {7, 8, 1}, {7, 9, 1}, {8, 9, 1}, {1, 4, 1}, {5, 8, 1}};
+    BOOST_TEST(mwmatching::maximum_weight_matching(edges).size() == 4);
 }
 
 BOOST_AUTO_TEST_CASE(test_fail_bad_input)
@@ -749,3 +770,72 @@ BOOST_AUTO_TEST_CASE(test_blossom_not_full)
 }
 
 BOOST_AUTO_TEST_SUITE_END()
+
+
+/* **********  Test graphs that force big values for dual / slack  ********** */
+
+BOOST_AUTO_TEST_SUITE(test_value_range)
+
+BOOST_AUTO_TEST_CASE(test_big_blossom_dual)
+{
+    /*
+     * Force modified blossom dual close to 2*maxweight.
+     *
+     *              [3]
+     *              / \
+     *          W-4/   \W
+     *      7     /     \
+     *  [1]-----[2]-----[4]
+     *          |   W-4
+     *         5|
+     *          |   1
+     *         [5]-----[6]
+     */
+    long w = std::numeric_limits<long>::max() / 6;
+    EdgeVectorLong edges = {{1, 2, 7}, {2, 3, w - 4}, {2, 4, w - 4}, {2, 5, 5}, {3, 4, w}, {5, 6, 1}};
+    Matching expect = {{1, 2}, {3, 4}, {5, 6}};
+    BOOST_TEST(mwmatching::maximum_weight_matching(edges) == expect);
+}
+
+BOOST_AUTO_TEST_CASE(test_negative_blossom_dual)
+{
+    /*
+     * Force modified blossom dual close to -maxweight.
+     *
+     *              [3]
+     *              / \
+     *            5/   \7
+     *       1    /     \
+     *  [1]-----[2]-----[4]
+     *           |   5
+     *          1|
+     *           |   W
+     *          [5]-----[6]
+     */
+    long w = std::numeric_limits<long>::max() / 6;
+    EdgeVectorLong edges = {{1, 2, 1}, {2, 3, 5}, {2, 4, 5}, {2, 5, 1}, {3, 4, 7}, {5, 6, w}};
+    Matching expect = {{1, 2}, {3, 4}, {5, 6}};
+    BOOST_TEST(mwmatching::maximum_weight_matching(edges) == expect);
+}
+
+BOOST_AUTO_TEST_CASE(test_big_edge_slack)
+{
+    /*
+     * Force modified edge slack close to 3*maxweight.
+     *
+     *      6       W      W-2      5
+     * [1]-----[2]-----[3]-----[4]-----[5]
+     *  |               |
+     *  |1              |1
+     *  |               |
+     * [6]-----[7]-----[8]-----[9]-----[10]
+     *      6       W      W-2      5
+     */
+    long w = std::numeric_limits<long>::max() / 6;
+    EdgeVectorLong edges = {{1, 2, 6}, {1, 6, 1}, {2, 3, w}, {3, 4, w - 2}, {3, 8, 1}, {4, 5, 5},
+                            {6, 7, 6}, {7, 8, w}, {8, 9, w - 2}, {9, 10, 5}};
+    Matching expect = {{1, 6}, {2, 3}, {4, 5}, {7, 8}, {9, 10}};
+    BOOST_TEST(mwmatching::maximum_weight_matching(edges) == expect);
+}
+
+BOOST_AUTO_TEST_SUITE_END()

cpp/test_priority_queue.cpp

@@ -0,0 +1,252 @@
+/*
+ * Unit tests for PriorityQueue data structure.
+ *
+ * Depends on the Boost.Test unit test framework.
+ * Tested with Boost v1.74, available from https://www.boost.org/
+ */
+
+#include <algorithm>
+#include <climits>
+#include <memory>
+#include <random>
+#include <string>
+#include <tuple>
+#include <vector>
+
+#define BOOST_TEST_MODULE datastruct
+#include <boost/test/unit_test.hpp>
+
+#include "priority_queue.hpp"
+
+
+BOOST_AUTO_TEST_SUITE(test_priority_queue)
+
+BOOST_AUTO_TEST_CASE(test_empty)
+{
+    using Queue = mwmatching::PriorityQueue<int, std::string>;
+    Queue q;
+    BOOST_TEST(q.empty() == true);
+}
+
+BOOST_AUTO_TEST_CASE(test_single)
+{
+    using Queue = mwmatching::PriorityQueue<int, std::string>;
+    Queue q;
+
+    Queue::Node n1;
+    BOOST_TEST(n1.valid() == false);
+
+    q.insert(&n1, 5, "a");
+
+    BOOST_TEST(n1.valid() == true);
+    BOOST_TEST(n1.prio() == 5);
+    BOOST_TEST(q.empty() == false);
+    BOOST_TEST(q.min_prio() == 5);
+    BOOST_TEST(q.min_elem() == std::string("a"));
+
+    q.set_prio(&n1, 3);
+    BOOST_TEST(n1.prio() == 3);
+    BOOST_TEST(q.min_prio() == 3);
+    BOOST_TEST(q.min_elem() == std::string("a"));
+
+    q.remove(&n1);
+    BOOST_TEST(n1.valid() == false);
+    BOOST_TEST(q.empty() == true);
+}
+
+BOOST_AUTO_TEST_CASE(test_simple)
+{
+    using Queue = mwmatching::PriorityQueue<int, char>;
+    Queue q;
+    Queue::Node nodes[10];
+
+    q.insert(&nodes[0], 9, 'a');
+    BOOST_TEST(q.min_prio() == 9);
+    BOOST_TEST(q.min_elem() == 'a');
+
+    q.insert(&nodes[1], 4, 'b');
+    BOOST_TEST(q.min_prio() == 4);
+    BOOST_TEST(q.min_elem() == 'b');
+
+    q.insert(&nodes[2], 7, 'c');
+    BOOST_TEST(q.min_prio() == 4);
+    BOOST_TEST(q.min_elem() == 'b');
+
+    q.insert(&nodes[3], 5, 'd');
+    BOOST_TEST(q.min_prio() == 4);
+    BOOST_TEST(q.min_elem() == 'b');
+
+    q.insert(&nodes[4], 8, 'e');
+    BOOST_TEST(q.min_prio() == 4);
+    BOOST_TEST(q.min_elem() == 'b');
+
+    q.insert(&nodes[5], 6, 'f');
+    BOOST_TEST(q.min_prio() == 4);
+    BOOST_TEST(q.min_elem() == 'b');
+
+    q.insert(&nodes[6], 4, 'g');
+    q.insert(&nodes[7], 5, 'h');
+    q.insert(&nodes[8], 2, 'i');
+    BOOST_TEST(q.min_prio() == 2);
+    BOOST_TEST(q.min_elem() == 'i');
+
+    q.insert(&nodes[9], 6, 'j');
+    BOOST_TEST(q.min_prio() == 2);
+    BOOST_TEST(q.min_elem() == 'i');
+
+    q.set_prio(&nodes[2], 1);
+    BOOST_TEST(q.min_prio() == 1);
+    BOOST_TEST(q.min_elem() == 'c');
+
+    q.set_prio(&nodes[4], 3);
+    BOOST_TEST(q.min_prio() == 1);
+    BOOST_TEST(q.min_elem() == 'c');
+
+    q.remove(&nodes[2]);
+    BOOST_TEST(q.min_prio() == 2);
+    BOOST_TEST(q.min_elem() == 'i');
+
+    q.remove(&nodes[8]);
+    BOOST_TEST(q.min_prio() == 3);
+    BOOST_TEST(q.min_elem() == 'e');
+
+    q.remove(&nodes[4]);
+    q.remove(&nodes[1]);
+    BOOST_TEST(q.min_prio() == 4);
+    BOOST_TEST(q.min_elem() == 'g');
+
+    q.remove(&nodes[3]);
+    q.remove(&nodes[9]);
+    BOOST_TEST(q.min_prio() == 4);
+    BOOST_TEST(q.min_elem() == 'g');
+
+    q.remove(&nodes[6]);
+    BOOST_TEST(q.min_prio() == 5);
+    BOOST_TEST(q.min_elem() == 'h');
+
+    BOOST_TEST(q.empty() == false);
+    q.clear();
+    BOOST_TEST(q.empty() == true);
+}
+
+BOOST_AUTO_TEST_CASE(test_increase_prio)
+{
+    using Queue = mwmatching::PriorityQueue<int, char>;
+
+    Queue q;
+    Queue::Node n1;
+
+    q.insert(&n1, 5, 'a');
+    q.set_prio(&n1, 8);
+    BOOST_TEST(n1.prio() == 8);
+    BOOST_TEST(q.min_prio() == 8);
+
+    q.clear();
+
+    Queue::Node n2, n3, n4;
+
+    q.insert(&n1, 9, 'A');
+    q.insert(&n2, 4, 'b');
+    q.insert(&n3, 7, 'c');
+    q.insert(&n4, 5, 'd');
+    BOOST_TEST(q.min_prio() == 4);
+    BOOST_TEST(q.min_elem() == 'b');
+
+    q.set_prio(&n2, 8);
+    BOOST_TEST(n2.prio() == 8);
+    BOOST_TEST(q.min_elem() == 'd');
+    BOOST_TEST(q.min_prio() == 5);
+
+    q.set_prio(&n3, 10);
+    BOOST_TEST(n3.prio() == 10);
+    BOOST_TEST(q.min_elem() == 'd');
+
+    q.remove(&n4);
+    BOOST_TEST(q.min_elem() == 'b');
+
+    q.remove(&n2);
+    BOOST_TEST(q.min_prio() == 9);
+    BOOST_TEST(q.min_elem() == 'A');
+
+    q.remove(&n1);
+    BOOST_TEST(q.min_elem() == 'c');
+    BOOST_TEST(q.min_prio() == 10);
+
+    q.remove(&n3);
+    BOOST_TEST(q.empty() == true);
+}
+
+BOOST_AUTO_TEST_CASE(test_random)
+{
+    using Queue = mwmatching::PriorityQueue<int, int>;
+    Queue q;
+
+    const int num_elem = 1000;
+    std::vector<std::tuple<std::unique_ptr<Queue::Node>, int, int>> elems;
+    int next_data = 0;
+
+    std::mt19937 rng(34567);
+
+    auto check = [&q,&elems]() {
+        int min_prio = q.min_prio();
+        int min_data = q.min_elem();
+        int best_prio = INT_MAX;
+        bool found = false;
+        for (const auto& v : elems) {
+            int this_prio = std::get<1>(v);
+            int this_data = std::get<2>(v);
+            best_prio = std::min(best_prio, this_prio);
+            if ((this_prio == min_prio) && (this_data == min_data)) {
+                found = true;
+            }
+        }
+        BOOST_TEST(found == true);
+        BOOST_TEST(min_prio == best_prio);
+    };
+
+    for (int i = 0; i < num_elem; ++i) {
+        ++next_data;
+        int prio = std::uniform_int_distribution<>(0, 1000000)(rng);
+        std::unique_ptr<Queue::Node> nptr(new Queue::Node);
+        q.insert(nptr.get(), prio, next_data);
+        elems.push_back(std::make_tuple(std::move(nptr), prio, next_data));
+        check();
+    }
+
+    for (int i = 0; i < 10000; ++i) {
+        int p = std::uniform_int_distribution<>(0, num_elem - 1)(rng);
+        Queue::Node* node = std::get<0>(elems[p]).get();
+        int prio = std::get<1>(elems[p]);
+        prio = std::uniform_int_distribution<>(0, 1000000)(rng);
+        q.set_prio(node, prio);
+        std::get<1>(elems[p]) = prio;
+        check();
+
+        p = std::uniform_int_distribution<>(0, num_elem - 1)(rng);
+        node = std::get<0>(elems[p]).get();
+        q.remove(node);
+        elems.erase(elems.begin() + p);
+        check();
+
+        ++next_data;
+        prio = std::uniform_int_distribution<>(0, 1000000)(rng);
+        std::unique_ptr<Queue::Node> nptr(new Queue::Node);
+        q.insert(nptr.get(), prio, next_data);
+        elems.push_back(std::make_tuple(std::move(nptr), prio, next_data));
+        check();
+    }
+
+    for (int i = 0; i < num_elem; ++i) {
+        int p = std::uniform_int_distribution<>(0, num_elem - 1 - i)(rng);
+        Queue::Node* node = std::get<0>(elems[p]).get();
+        q.remove(node);
+        elems.erase(elems.begin() + p);
+        if (! elems.empty()) {
+            check();
+        }
+    }
+
+    BOOST_TEST(q.empty() == true);
+}
+
+BOOST_AUTO_TEST_SUITE_END()

pylintrc

@@ -13,6 +13,7 @@ disable=consider-using-in,
         too-many-lines,
         too-many-locals,
         too-many-nested-blocks,
+        too-many-positional-arguments,
         too-many-public-methods,
         too-many-return-statements,
         too-many-statements,

python/mwmatching/algorithm.py

@@ -1688,24 +1688,16 @@ class MatchingContext:
                 (p, q, _w) = edges[e]
                 y = p if p != x else q
 
-                # Consider the edge between vertices "x" and "y".
-                # Update delta2 or delta3 tracking accordingly.
-                #
-                # We don't actually use the edge right now to extend
-                # the alternating tree or create a blossom or alternating path.
-                # If appropriate, insert this edge into delta2 or delta3
-                # tracking.
-                # Insert this edge into delta2 or delta3 tracking
-                # Try to pull this edge into an alternating tree.
-
                 # Ignore edges that are internal to a blossom.
                 by = self.top_level_blossom(y)
                 if bx is by:
                     continue
 
                 if by.label == LABEL_S:
+                    # Edge between S-vertices.
                     self.delta3_add_edge(e)
                 else:
+                    # Edge to T-vertex or unlabeled vertex.
                     self.delta2_add_edge(e, y, by)
 
         self.scan_queue.clear()

test_cpp.sh

@@ -9,8 +9,11 @@ echo
 g++ --version
 echo
 
-make -C cpp run_matching test_mwmatching
+make -C cpp clean
+make -C cpp run_matching test_mwmatching test_concatenable_queue test_priority_queue
 cpp/test_mwmatching
+cpp/test_concatenable_queue
+cpp/test_priority_queue
 
 echo
 echo ">> Running test graphs"

tox.ini

@@ -10,6 +10,7 @@ env_list =
     py310
     py311
     py312
+    py313
     pypy3
 
 [testenv]