E061. 算法 Куна нахождения наибольшего паросочетания в двудольном 图е

Источник: e-maxx.ru/algo, страница PDF 184.

Дан двудольный 图

, содержащий

вершин и

рёбер. it is required find наибольшее паросочетание, т.е. выбрать как можно больше рёбер, чтобы ни одно выбранное edge не имело общей вершины ни с каким другим выбранным edgeм.

Описание 算法а

Необходимые определения

Паросочетанием

называется набор попарно несмежных рёбер 图а (иными словами, любой вершине

图а должно быть инцидентно не более одного ребра из множества

). Мощностью паросочетания назовём

количество рёбер в нём. Наибольшим (или максимальным) паросочетанием назовём паросочетание, мощность которого максимальна среди всех возможных паросочетаний в данном 图е. Все те вершины, у которых есть смежное edge из паросочетания (т.е. которые имеют степень ровно один в под图е, образованном

),

назовём насыщенными этим паросочетанием.

Цепью длины

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

содержащий

рёбер. Чередующейся цепью (в двудольном 图е, относительно некоторого паросочетания) назовём цепь, в которой рёбра поочередно принадлежат/не принадлежат паросочетанию. Увеличивающей цепью (в двудольном 图е, относительно некоторого паросочетания) назовём чередующуюся цепь, у которой начальная и конечная вершины не принадлежат паросочетанию.

Теорема Бержа

Формулировка. Паросочетание является максимальным тогда и только тогда, когда не существует увеличивающих относительно него цепей.

证明 необходимости. Покажем, что если паросочетание

максимально, то не

существует увеличивающей относительно него цепи. 证明 это будет конструктивным: мы покажем,

как увеличить с помощью этой увеличивающей цепи

мощность паросочетания

на единицу. Для этого выполним так называемое чередование паросочетания вдоль цепи

. Мы помним, что по определению

первое edge цепи

не принадлежит паросочетанию, второе — принадлежит, третье — снова не принадлежит, четвёртое — принадлежит, и т.д. Давайте поменяем состояние всех рёбер вдоль цепи

: те рёбра, которые не 输入или

в паросочетание (первое, третье и т.д. до последнего) включим в паросочетание, а рёбра, которые раньше 输入или в паросочетание (второе, четвёртое и т.д. до предпоследнего) — удалим из него. Понятно, что мощность паросочетания при этом увеличилась на единицу (потому что было добавлено на одно edge больше, чем удалено). Осталось проверить, что мы построили корректное паросочетание, т.е. что никакая vertex 图а не имеет сразу двух смежных рёбер из этого паросочетания. Для всех вершин чередующей цепи

,

кроме первой и последней, это следует из самого 算法а чередования: сначала мы у каждой такой вершины удалили смежное edge, потом добавили. Для первой и последней вершины цепи

также ничего не могло

нарушиться, поскольку до чередования они должны были быть ненасыщенными. Наконец, для всех остальных вершин,

— не 输入ящих в цепь

, — очевидно, ничего не поменялось. Таким образом, мы в самом деле построили паросочетание, и на единицу большей мощности, чем старое, что и завершает 证明 необходимости. 证明 достаточности. Докажем, что если относительно некоторого паросочетания

нет увеличивающих путей, то оно — максимально.

证明 проведём от противного. Пусть есть паросочетание

, имеющее бОльшую мощность, чем

. Рассмотрим симметрическую разность

этих двух паросочетаний, т.е. оставим все рёбра, 输入ящие в

или в

, но не в оба одновременно.

Понятно, что множество рёбер

— уже наверняка не паросочетание. Рассмотрим, какой вид это множество рёбер имеет; для удобства будем рассматривать его как 图. В этом 图е каждая vertex, очевидно, имеет степень не выше 2 (потому что каждая vertex может иметь максимум два смежных ребра — из одного паросочетания и из другого). Легко понять, что тогда этот 图 состоит только из циклов или путей, причём ни те, ни другие не пересекаются друг с другом.

Теперь заметим, что и пути в этом 图е

могут быть не любыми, а только чётной длины. В самом деле, в любом пути

в 图е

рёбра чередуются: после ребра из

идёт edge из

, и наоборот. Теперь, если мы рассмотрим какой-

то путь нечётной длины в 图е

, то получится, что в исходном 图е

это будет увеличивающей цепью либо

для паросочетания

, либо для

. Но этого быть не могло, потому что в случае паросочетания

это противоречит

с 题意м, а в случае

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

,

то

), но это 题意 не может выполниться в цикле нечётной длины — в нём обязательно найдутся два соседних ребра из одного паросочетания, что противоречит определению паросочетания.

Таким образом, все пути и циклы 图а

имеют чётную длину. Следовательно, 图

содержит равное количество рёбер из

и из

. Но, given, что в

содержатся все рёбра

и

,

за исключением их общих рёбер, то отсюда следует, что мощность

и

совпадают. Мы пришли к противоречию:

по предположению паросочетание

было не максимальным, значит, теорема доказана.

算法 Куна

算法 Куна — непосредственное применение теоремы Бержа. Его можно кратко описать так: сначала возьмём пустое паросочетание, а потом — пока в 图е удаётся find увеличивающую цепь, — будем выполнять чередование паросочетания вдоль этой цепи, и повторять процесс поиска увеличивающей цепи. Как только такую цепь find не удалось — процесс останавливаем, — текущее паросочетание и есть максимальное. Осталось детализировать способ нахождения увеличивающих цепей. 算法 Куна — просто ищет любую из таких цепей с помощью обхода в глубину или в ширину. 算法 Куна просматривает все вершины 图а по очереди, запуская из каждой обход, пытающийся find увеличивающую цепь, начинающуюся в этой вершине. Удобнее описывать этот 算法, считая, что 图 уже разбит на две доли (хотя на самом деле 算法 можно реализовать и так, чтобы ему не давался на 输入 图, явно разбитый на две доли).

算法 просматривает все вершины

первой доли 图а:

. Если текущая vertex

уже

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

текущей ненасыщенной вершине

первой доли. Просматриваем все рёбра из этой вершины, пусть текущее edge —

это edge

. Если vertex

ещё не насыщена паросочетанием, то, значит, мы смогли find

увеличивающую цепь: она состоит из единственного ребра

; в таком случае просто включаем это edge

в паросочетание и прекращаем поиск увеличивающей цепи из вершины

. Иначе, — если

уже насыщена каким-

то edgeм

, то попытаемся пройти вдоль этого ребра: тем самым мы попробуем find увеличивающую

цепь, проходящую через рёбра

,

. Для этого просто перейдём в нашем обходе в вершину

— теперь

мы уже пробуем find увеличивающую цепь из этой вершины. Можно понять, что в результате этот обход, запущенный из вершины

, либо найдёт увеличивающую цепь, и тем

самым насытит вершину

, либо же такой увеличивающей цепи не найдёт (и, следовательно, эта vertex

уже

не сможет стать насыщенной).

После того, как все вершины

будут просмотрены, текущее паросочетание будет максимальным.

运行时间

Итак, 算法 Куна можно представить как серию из

запусков обхода в глубину/ширину на всём 图е.

Следовательно, всего этот 算法 исполняется за время

, что в худшем случае есть

. Однако эту оценку можно немного улучшить. Оказывается, для 算法а Куна важно то, какая доля выбрана за первую, а какая — за вторую. В самом деле, в описанной выше реализации запуски обхода в глубину/ ширину происходят только из вершин первой доли, поэтому весь 算法 исполняется за время

, где

— number вершин первой доли. В худшем случае это составляет

(где

— number вершин второй

доли). Отсюда видно, что выгоднее, когда первая доля содержит меньшее number вершин, нежели вторая. На

очень несбалансированных 图ах (когда

и

сильно отличаются) это выливается в значительную разницу

времён работы.

实现

Приведём здесь реализацию вышеописанного 算法а, основанную на обходе в глубину, и принимающей двудольный 图 в виде явно разбитого на две доли 图а. Эта 实现 весьма лаконична, и, возможно, её стоит запомнить именно в таком виде.

Здесь

— number вершин в первой доле,

— во второй доле,

— список рёбер из вершины

первой доли (т.е.

список номеров вершин, в которые ведут эти рёбра из

). Вершины в обеих долях занумерованы независимо, т.е.

первая доля — с номерами

, вторая — с номерами

. Дальше идут два вспомогательных 数组а:

и

. Первый —

— содержит в себе информацию о

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

— это номер вершины первой доли, связанной edgeм с вершиной

второй доли (или

, если

никакого ребра паросочетания из

не 输出ит). Второй 数组 —

— обычный 数组 "посещённостей" вершин

в обходе в глубину (он нужен, просто чтобы обход в глубину не заходил в одну вершину дважды).

Функция

— и есть обход в глубину. Она returns

, если ей удалось find увеличивающую цепь

из вершины

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

Внутри функции просматриваются все рёбра, исходящие из вершины

первой доли, и затем проверяется: если это

edge ведёт в ненасыщенную вершину

, либо если эта vertex

насыщена, но удаётся find увеличивающую

цепь рекурсивным запуском из

, то мы говорим, что мы нашли увеличивающую цепь, и перед возвратом

из функции с результатом

производим чередование в текущем ребре: перенаправляем edge, смежное с

,

в вершину

. В основной программе сначала указывается, что текущее паросочетание — пустое (список

заполняется числами

). Затем перебирается vertex

первой доли, и из неё запускается обход в глубину

,

предварительно обнулив 数组

. Стоит заметить, что размер паросочетания легко получить как number вызовов

в основной

программе, вернувших результат

. Само искомое максимальное паросочетание содержится в 数组е

.

int n, k;

vector < vector<int> > g;

vector<int> mt;

vector<char> used;

bool try_kuhn (int v) {

if (used[v])  return false;

used[v] = true;

for (size_t i=0; i<g[v].size(); ++i) {

int to = g[v][i];

if (mt[to] == -1 || try_kuhn (mt[to])) {

mt[to] = v;

return true;

} }

return false;

}

int main() {

... чтение 图а ...

mt.assign (k, -1);

for (int v=0; v<n; ++v) {

used.assign (n, false);

try_kuhn (v);

}

for (int i=0; i<k; ++i)

if (mt[i] != -1)

printf ("%d %d\n", mt[i]+1, i+1);

} Ещё раз повторим, что 算法 Куна легко реализовать и так, чтобы он работал на 图ах, про которые известно, что они двудольные, но явное их разбиение на две доли не найдено. В этом случае придётся отказаться от удобного разбиения на две доли, и всю информацию хранить для всех вершин 图а. Для этого 数组 списков

теперь задаётся не только для вершин первой доли, а для всех вершин 图а (понятно, теперь вершины обеих

долей занумерованы в общей нумерации — от

до

). 数组ы

и

теперь также определены для

вершин обеих долей, и, соответственно, их нужно поддерживать в этом состоянии.

Улучшенная 实现

Модифицируем 算法 следующим образом. До основного цикла 算法а найдём каким-нибудь простым 算法ом произвольное паросочетание (простым эвристическим 算法ом), и лишь затем будем выполнять цикл с вызовами функции kuhn(), который будет улучшать это паросочетание. В результате 算法 будет работать заметно быстрее на случайных 图ах — потому что в большинстве 图ов можно легко набрать паросочетание достаточно большого веса с помощью эвристики, а потом улучшить найденное паросочетание до максимального уже обычным 算法ом Куна. Тем самым мы сэкономим на запусках обхода в глубину из тех вершин, которые мы уже включили с помощью эвристики в текущее паросочетание. На示例, можно просто перебрать все вершины первой доли, и для каждой из них find произвольное edge, которое можно добавить в паросочетание, и добавить его. Даже такая простая эвристика способна ускорить 算法 Куна в несколько раз. Следует обратить внимание на то, что основной цикл придётся немного модифицировать. Поскольку при вызове

функции

в основном цикле предполагается, что текущая vertex ещё не 输入ит в паросочетание, то нужно добавить соответствующую проверку. В реализации изменится только код в функции main():

int main() {

... чтение 图а ...

mt.assign (k, -1);

vector<char> used1 (n);

for (int i=0; i<n; ++i)

for (size_t j=0; j<g[i].size(); ++j)

if (mt[g[i][j]] == -1) {

mt[g[i][j]] = i;

used1[i] = true;

break;

}

for (int i=0; i<n; ++i) {

if (used1[i])  continue;

used.assign (n, false);

try_kuhn (i);

}

for (int i=0; i<k; ++i)

if (mt[i] != -1)

printf ("%d %d\n", mt[i]+1, i+1);

} Другой хорошей эвристикой является следующая. На каждом шаге будет искать вершину наименьшей степени (но не изолированную), из неё выбирать любое edge и добавлять его в паросочетание, затем удаляя обе эти вершины со всеми инцидентными им рёбрами из 图а. Такая жадность работает очень хорошо на случайных 图ах, даже в большинстве случаев строит максимальное паросочетание (хотя и против неё есть тест, на котором она найдёт паросочетание значительно меньшей величины, чем максимальное).

using System;
using System.Collections.Generic;
using System.Linq;

public static class AlgorithmDraft
{
    // Auto-generated C# draft from the original e-maxx C/C++ listing. Review before production use.
    int n, k;
    vector < List<int> > g;
    List<int> mt;
    List<char> used;
    bool try_kuhn (int v) {
            if (used[v])  return false;
            used[v] = true;
            for (size_t i=0; i<g[v].size(); ++i) {
                    int to = g[v][i];
                    if (mt[to] == -1 || try_kuhn (mt[to])) {
                            mt[to] = v;
                            return true;
                    }
            }
            return false;
    }
    int main() {
            ... чтение графа ...
            mt.assign (k, -1);
            for (int v=0; v<n; ++v) {
                    used.assign (n, false);
                    try_kuhn (v);
            }
            for (int i=0; i<k; ++i)
                    if (mt[i] != -1)
                            Console.Write ("%d %d\n", mt[i]+1, i+1);
    }
    int main() {
            ... чтение графа ...
            mt.assign (k, -1);
            List<char> used1 (n);
            for (int i=0; i<n; ++i)
                    for (size_t j=0; j<g[i].size(); ++j)
                            if (mt[g[i][j]] == -1) {
                                    mt[g[i][j]] = i;
                                    used1[i] = true;
                                    break;
                            }
            for (int i=0; i<n; ++i) {
                    if (used1[i])  continue;
                    used.assign (n, false);
                    try_kuhn (i);
            }
            for (int i=0; i<k; ++i)
                    if (mt[i] != -1)
                            Console.Write ("%d %d\n", mt[i]+1, i+1);
    }
}

int n, k;
vector < vector<int> > g;
vector<int> mt;
vector<char> used;
bool try_kuhn (int v) {
        if (used[v])  return false;
        used[v] = true;
        for (size_t i=0; i<g[v].size(); ++i) {
                int to = g[v][i];
                if (mt[to] == -1 || try_kuhn (mt[to])) {
                        mt[to] = v;
                        return true;
                }
        }
        return false;
}
int main() {
        ... чтение графа ...
        mt.assign (k, -1);
        for (int v=0; v<n; ++v) {
                used.assign (n, false);
                try_kuhn (v);
        }
        for (int i=0; i<k; ++i)
                if (mt[i] != -1)
                        printf ("%d %d\n", mt[i]+1, i+1);
}
int main() {
        ... чтение графа ...
        mt.assign (k, -1);
        vector<char> used1 (n);
        for (int i=0; i<n; ++i)
                for (size_t j=0; j<g[i].size(); ++j)
                        if (mt[g[i][j]] == -1) {
                                mt[g[i][j]] = i;
                                used1[i] = true;
                                break;
                        }
        for (int i=0; i<n; ++i) {
                if (used1[i])  continue;
                used.assign (n, false);
                try_kuhn (i);
        }
        for (int i=0; i<k; ++i)
                if (mt[i] != -1)
                        printf ("%d %d\n", mt[i]+1, i+1);
}

import java.util.*;
import java.math.*;

public class AlgorithmDraft {
    // Auto-generated Java draft from the original e-maxx C/C++ listing. Review before production use.
    int n, k;
    vector < ArrayList<Integer> > g;
    ArrayList<Integer> mt;
    ArrayList<Character> used;
    boolean try_kuhn (int v) {
            if (used[v])  return false;
            used[v] = true;
            for (size_t i=0; i<g[v].size(); ++i) {
                    int to = g[v][i];
                    if (mt[to] == -1 || try_kuhn (mt[to])) {
                            mt[to] = v;
                            return true;
                    }
            }
            return false;
    }
    int main() {
            ... чтение графа ...
            mt.assign (k, -1);
            for (int v=0; v<n; ++v) {
                    used.assign (n, false);
                    try_kuhn (v);
            }
            for (int i=0; i<k; ++i)
                    if (mt[i] != -1)
                            System.out.print ("%d %d\n", mt[i]+1, i+1);
    }
    int main() {
            ... чтение графа ...
            mt.assign (k, -1);
            ArrayList<Character> used1 (n);
            for (int i=0; i<n; ++i)
                    for (size_t j=0; j<g[i].size(); ++j)
                            if (mt[g[i][j]] == -1) {
                                    mt[g[i][j]] = i;
                                    used1[i] = true;
                                    break;
                            }
            for (int i=0; i<n; ++i) {
                    if (used1[i])  continue;
                    used.assign (n, false);
                    try_kuhn (i);
            }
            for (int i=0; i<k; ++i)
                    if (mt[i] != -1)
                            System.out.print ("%d %d\n", mt[i]+1, i+1);
    }
}