Source code for tryalgo.kuhn_munkres

#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""\
Maximum profit bipartite matching by Kuhn-Munkres
jill-jenn vie, christoph durr and samuel tardieu - 2014-2019

primal LP

    max sum_{u,v} w[u,v] * x[u,v]

    such that
    for all u in U: sum_v x[u,v] == 1   (l[u])

    for all v in V: sum_u x[u,v] <= 1   (l[v])

    for all u,v: x[u,v] >= 0


dual LP

    min sum_u l[u] + sum_v l[v]

    such that
    for all u,v:  l[u] + l[v] >= w[u,v]   (*)

    for all u in U: l[u] is arbitrary (free variable)
    for all v in V: l[v] >= 0


primal-dual algorithm:

    Start with trivial solution l for dual and with trivial
    non-solution x for primal.

    Iteratively improve primal or dual solution, maintaining complementary
    slackness conditions.

"""


# snip{
# pylint: disable=too-many-locals, too-many-branches
[docs]def kuhn_munkres(G, TOLERANCE=1e-6): """Maximum profit bipartite matching by Kuhn-Munkres :param G: weight matrix where G[u][v] is the weight of edge (u,v), :param TOLERANCE: a value with absolute value below tolerance is considered as being zero. If G consists of integer or fractional values then TOLERANCE can be chosen 0. :requires: graph (U,V,E) is complete bi-partite graph with len(U) <= len(V) float('-inf') or float('inf') entries in G are allowed but not None. :returns: matching table from U to V, value of matching :complexity: :math:`O(|U|^2 |V|)` """ nU = len(G) U = range(nU) nV = len(G[0]) V = range(nV) assert nU <= nV mu = [None] * nU # empty matching mv = [None] * nV lu = [max(row) for row in G] # trivial labels lv = [0] * nV for root in U: # build an alternate tree au = [False] * nU # au, av mark nodes... au[root] = True # ... covered by the tree Av = [None] * nV # Av[v] successor of v in the tree # for every vertex u, slack[u] := (val, v) such that # val is the smallest slack on the constraints (*) # with fixed u and v being the corresponding vertex slack = [(lu[root] + lv[v] - G[root][v], root) for v in V] while True: (delta, u), v = min((slack[v], v) for v in V if Av[v] is None) assert au[u] if delta > TOLERANCE: # tree is full for u0 in U: # improve labels if au[u0]: lu[u0] -= delta for v0 in V: if Av[v0] is not None: lv[v0] += delta else: (val, arg) = slack[v0] slack[v0] = (val - delta, arg) assert abs(lu[u] + lv[v] - G[u][v]) <= TOLERANCE # equality Av[v] = u # add (u, v) to A if mv[v] is None: break # alternating path found u1 = mv[v] assert not au[u1] au[u1] = True # add (u1, v) to A for v1 in V: if Av[v1] is None: # update margins alt = (lu[u1] + lv[v1] - G[u1][v1], u1) if slack[v1] > alt: slack[v1] = alt while v is not None: # ... alternating path found u = Av[v] # along path to root prec = mu[u] mv[v] = u # augment matching mu[u] = v v = prec return (mu, sum(lu) + sum(lv))
# snip}