Венгерский алгоритм решения задачи о назначениях

Венгерский алгоритм — алгоритм, решающий задачу о назначениях за полиномиальное время. Оригинальная версия была придумана и разработана Х. Куном в 1955 году и имела асимптотику [math] O(n^4) [/math], но позже Эдмонс и Карп (а также, независимо от них, Томидзава) показали, что можно улучшить ее до [math] O(n^3) [/math].

Постановка задачи

Пусть дан взвешенный полный двудольный граф c целыми весами ребер [math] K_{n, n} [/math], нужно найти в нем полное паросочетание минимального веса. Вес паросочетания определяется как сумма весов его ребер. Далее будем обозначать левую и правую доли графа за [math] X [/math] и [math] Y [/math] соответственно, вес ребра [math] xy [/math] — как [math] c(xy) [/math].

Некоторые полезные утверждения

Лемма:

Если веса всех ребер графа, инцидентных какой-либо вершине, изменить (увеличить или уменьшить) на одно и то же число, то в новом графе оптимальное паросочетание будет состоять из тех же ребер, что и в старом.

Доказательство:

[math]\triangleright[/math]

Полное паросочетание для каждой вершины содержит ровно одно ребро, инцидентное этой вершине. Указанная операция изменит на одно и то же число вес любого паросочетания. Значит, ребро, которое принадлежало оптимальному паросочетанию в старом графе, в новом графе тоже будет ему принадлежать.

[math]\triangleleft[/math]

Далее будем рассматривать только графы с неотрицательной весовой функцией, так как, согласно этой лемме, задачу о назначениях на остальных графах можно свести к задаче о назначениях на них.

Лемма:

Выделим в множествах [math]X[/math] и [math]Y[/math] подмножества [math]X', Y'[/math]. Пусть [math]d = \min \{c(xy)|\ x \in X \setminus X', y \in Y'\}[/math]. Прибавим [math] d [/math] ко всем весам ребер, инцидентных вершинам из [math]X'[/math]. Затем отнимем [math] d [/math] от всех весов ребер, инцидентных вершинам из [math]Y'[/math] (далее для краткости эта операция обозначается как [math] X' \uparrow\downarrow Y' [/math]). Тогда: Веса всех ребер графа останутся неотрицательными. Веса ребер вида [math]xy[/math], где [math]x \in X', y \in Y'[/math] или [math]x \in X \backslash X', y \in Y \backslash Y'[/math], не изменятся.

Доказательство:

[math]\triangleright[/math]

Рассмотрим матрицу весов графа. Не умаляя общности, можно сказать, что множества [math] X' [/math] и [math] Y' [/math] состоят из первых элементов множеств [math] X [/math] и [math] Y [/math] соответственно (мы упорядочиваем множества по номерам вершин). Тогда вся матрица делится на 4 блока: [math] Y' [/math] [math] Y \backslash Y' [/math] [math] X' [/math] [math] A + d - d [/math] [math] C + d [/math] [math] X \backslash X' [/math] [math] B - d [/math] [math] D [/math] Веса группы [math] A [/math] будут сначала увеличены, а потом уменьшены на [math] d [/math], поэтому они не изменятся, веса группы [math] D [/math] вообще изменяться не будут. Все веса группы [math] B [/math] будут уменьшены на [math] d [/math], но [math] d [/math] — минимум среди этих весов, поэтому они останутся неотрицательными.

[math]\triangleleft[/math]

Лемма:

Если веса всех ребер графа неотрицательны и некоторое полное паросочетание состоит из ребер нулевого веса, то оно является оптимальным.

Доказательство:

[math]\triangleright[/math]

Действительно, паросочетание с какими-то другими весами ребер имеет больший вес и оптимальным не является.

[math]\triangleleft[/math]

Общий метод

Доказанные ранее утверждения позволяют придумать схему алгоритма, решающего задачу о назначениях: нужно найти полное паросочетание из ребер нулевого веса в графе, полученном из исходного преобразованиями, описанными в первых двух леммах.

Алгоритм, решающий задачу, работает с графом, как с матрицей весов.

1. Вычитаем из каждой строки значение ее минимального элемента. Теперь в каждой строке есть хотя бы один нулевой элемент.
2. Вычитаем из каждого столбца значение его минимального элемента. Теперь в каждом столбце есть хотя бы один нулевой элемент.
3. Ищем в текущем графе полное паросочетание из ребер нулевого веса:
4. • Если оно найдено, то желаемый результат достигнут, алгоритм закончен.
   • В противном случае, покроем нули матрицы весов минимальным количеством строк и столбцов (это не что иное, как нахождение минимального вершинного покрытия в двудольном графе). Пусть [math] X_c [/math] и [math] Y_c [/math] — множества вершин минимального вершинного покрытия из левой и правой долей (то есть, строк и столбцов) соответственно, тогда применим преобразование [math] X_c \uparrow\downarrow (Y \setminus Y_c) [/math]. Для этого преобразования [math] d [/math] будет минимумом по всем ребрам между [math] X \setminus X_c [/math] и [math] Y \setminus Y_c [/math], то есть, ребер нулевого веса здесь нет, поэтому, после его выполнения в матрице весов появится новый нуль. После этого перейдем к шагу 1.

Анализ времени работы

Поиск максимального паросочетания или минимального вершинного покрытия в двудольном графе совершается за [math] O(n^3) [/math] операций. При каждом повторении шагов 1-4 в матрице весов появляется новый нуль. Этот нуль соответствует некоторому новому ребру между вершинами из множеств [math] X \setminus X_c [/math] и [math] Y \setminus Y_c [/math]. Всего в графе n^2 ребер, значит, всего будет совершено не более [math] O(n^2) [/math] итераций внешнего цикла. Поэтому, верхняя оценка времени работы данного метода — [math] O(n^5) [/math]. Более точная оценка довольно сложна и зависит от порядка чисел в матрице весов графа.

Алгоритм за [math]O(n^3)[/math]

Предложенный метод корректен, но не является оптимальным по времени работы, оно может быть улучшено до [math] O(n^3) [/math].

Заметим, что искать максимальное паросочетание каждый раз заново необязательно, можно просто постепенно изменять его, получая новые нулевые ребра и строя дополняющие цепи на них (примерно так же, как и в алгоритме Куна, но здесь мы фактически строим целое дерево возможных дополняющих цепей).

При их построении нам потребуется сохранять некоторую дополнительную информацию. Для каждого столбца, если в нем есть элемент из текущего паросочетания, будем хранить номер строки этого элемента. Также, при построении дополняющей цепи, будем для каждого элемента паросочетания хранить, на какой добавленный нуль из того же столбца его нужно заменить, если цепь будет проходить через него, а для каждого добавленного нуля — элемент из паросочетания, находящийся в той же строке (если он существует). Эти данные помогут нам восстановить дополняющую цепь.

Будем рассматривать строки матрицы последовательно. Пусть в данный момент мы рассматриваем [math] i [/math]-ую строку, а в первых [math] i - 1 [/math] строках максимальное паросочетание уже найдено.

Множество строк [math] C [/math], из которых мы будем вычитать минимумы, вначале состоит только из строки [math] i [/math], множество [math] R [/math] пусто. Применим к матрице операцию [math] R\uparrow\downarrow C [/math]. Если новый нуль с координатами [math] xy [/math] появился в столбце, не занятом паросочетанием, то дополняющая цепь найдена, [math] xy [/math] — ее последнее ребро. Если нет, то существует элемент [math] x'y [/math] из паросочетания. Тогда добавим строку [math] x' [/math] и столбец [math] y [/math] в множества [math] C [/math] и [math] R [/math] соответственно, отметим [math] xy [/math] как замену для [math] x'y [/math] и повторим операцию.

После того, как мы нашли дополняющую цепь, производим серию замен вдоль нее. Цепь состоит из элементов паросочетания и добавленных на текущей итерации нулей, причем для каждого элемента цепи (кроме, разумеется, последнего) известен следующий за ним. После того, как мы удалим из паросочетания все ребра цепи, принадлежащие ему и добавим не принадлежащие, суммарный размер паросочетания увеличится.

Краткий псевдокод данного алгоритма выглядит примерно так:

HungarianAlgorithm(A)
1 for i = 0..height-1
2   do
3     создать в матрице A новый нулевой элемент
4   while новый элемент приводит к коллизии
5   изменить метки элементов вдоль найденной цепочки
6 return найденное паросочетание

При добавлении строки с номером [math] i [/math] на поиск дополняющей цепи требуется [math] O(i) [/math] итераций, если каждая итерация требует [math] O(m) [/math] действий, то общее время работы алгоритма составляет [math] O(mn^2) [/math]. При наивной реализации [math] m = O(n^2) [/math]. Но можно не пересчитывать каждый раз все измененные элементы матрицы. Вместо них будем хранить и обновлять для каждой строки и каждого столбца минимум по данной строке/столбцу. Текущее значение элемента в каждый момент можно восстановить, используя эти вспомогательные массивы. При выполнении преобразования [math] R\uparrow\downarrow C [/math] изменим соответствующие элементы массивов на [math]d[/math]. Так можно выполнять все требуемые действия за [math] O(n) [/math] операций, тогда на выполнение всего алгоритма потребуется [math] O(n^3) [/math] действий.

Ссылки

• Венгерский алготитм в Википедии
• Визуализатор алгоритма
• Реализация венгерского алгоритма на C++

Литература

• Асанов М., Баранский В., Расин В. — Дискретная математика: Графы, матроиды, алгоритмы — 2010, 368 стр.