MAXimal

добавлено: 13 Aug 2010 19:48
редактировано: 14 Jun 2012 5:01

Алгоритм Диница

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

Пусть дана сеть, т.е. ориентированный граф G, в котором каждому ребру (u,v) приписана пропускная способность c_{uv}, а также выделены две вершины — исток s и сток t.

Требуется найти в этой сети поток f_{uv} из истока s в сток t максимальной величины.

Немного истории

Этот алгоритм был опубликован советским (израильским) учёным Ефимом Диницем (Yefim Dinic, иногда пишется как "Dinitz") в 1970 г., т.е. даже на два года раньше опубликования алгоритма Эдмондса-Карпа (впрочем, оба алгоритма были независимо открыты в 1968 г.).

Кроме того, следует отметить, что некоторые упрощения алгоритма были произведены Шимоном Ивеном (Shimon Even) и его учеником Алоном Итаи (Alon Itai) в 1979 г. Именно благодаря им алгоритм получил свой современный облик: они применили к идее Диница концепцию блокирующих потоков Александра Карзанова (Alexander Karzanov, 1974 г.), а также переформулировали алгоритм к той комбинации обхода в ширину и в глубину, в которой сейчас этот алгоритм и излагается везде.

Развитие идей по отношению к потоковым алгоритмам крайне интересно рассматривать, учитывая "железный занавес" тех лет, разделявший СССР и Запад. Видно, как иногда похожие идеи появлялись почти одновременно (как в случае алгоритма Диница и алгоритма Эдмондса-Карпа), правда, имея при этом разную эффективность (алгоритм Диница на один порядок быстрее); иногда же, наоборот, появление идеи по одну сторону "занавеса" опережало аналогичный ход по другую сторону более чем на десятилетие (как алгоритм Карзанова проталкивания в 1974 г. и алгоритм Гольдберга (Goldberg) проталкивания в 1985 г.).

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

Введём три необходимых определения (каждое из них является независимым от остальных), которые затем будут использоваться в алгоритме Диница.

Остаточной сетью G^R по отношению к сети G и некоторому потоку f в ней называется сеть, в которой каждому ребру (u,v) \in G с пропускной способностью c_{uv} и потоком f_{uv} соответствуют два ребра:

  • (u,v) с пропускной способностью c_{uv}^R = c_{uv} - f_{uv}

  • (v,u) с пропускной способностью c_{vu}^R = f_{uv}

Стоит отметить, что при таком определении в остаточной сети могут появляться кратные рёбра: если в исходной сети было как ребро (u,v), так и (v,u).

Остаточное ребро можно интуитивно понимать как меру того, насколько ещё можно увеличить поток вдоль какого-то ребра. В самом деле, если по ребру (u,v) с пропускной способностью c_{uv} протекает поток f_{uv}, то потенциально по нему можно пропустить ещё c_{uv}-f_{uv} единиц потока, а в обратную сторону можно пропустить до f_{uv} единиц потока, что будет означать отмену потока в первоначальном направлении.

Блокирующим потоком в данной сети называется такой поток, что любой путь из истока s в сток t содержит насыщенное этим потоком ребро. Иными словами, в данной сети не найдётся такого пути из истока в сток, вдоль которого можно беспрепятственно увеличить поток.

Блокирующий поток не обязательно максимален. Теорема Форда-Фалкерсона говорит о том, что поток будет максимальным тогда и только тогда, когда в остаточной сети не найдётся s-t пути; в блокирующем же потоке ничего не утверждается о существовании пути по рёбрам, появляющимся в остаточной сети.

Слоистая сеть для данной сети строится следующим образом. Сначала определяются длины кратчайших путей из истока s до всех остальных вершин; назовём уровнем {\rm level}[v] вершины её расстояние от истока. Тогда в слоистую сеть включают все те рёбра (u,v) исходной сети, которые ведут с одного уровня на какой-либо другой, более поздний, уровень, т.е. {\rm level}[u] + 1 = {\rm level}[v] (почему в этом случае разница расстояний не может превосходить единицы, следует из свойства кратчайших расстояний). Таким образом, удаляются все рёбра, расположенные целиком внутри уровней, а также рёбра, ведущие назад, к предыдущим уровням.

Очевидно, слоистая сеть ациклична. Кроме того, любой s-t путь в слоистой сети является кратчайшим путём в исходной сети.

Построить слоистую сеть по данной сети очень легко: для этого надо запустить обход в ширину по рёбрам этой сети, посчитав тем самым для каждой вершины величину {\rm level}[], и затем внести в слоистую сеть все подходящие рёбра.

Примечание. Термин "слоистая сеть" в русскоязычной литературе не употребляется; обычно эта конструкция называется просто "вспомогательным графом". Впрочем, на английском языке обычно используется термин "layered network".

Алгоритм

Схема алгоритма

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

Этот алгоритм схож с алгоритмом Эдмондса-Карпа, но основное отличие можно понимать так: на каждой итерации поток увеличивается не вдоль одного кратчайшего s-t пути, а вдоль целого набора таких путей (ведь именно такими путями и являются пути в блокирующем потоке слоистой сети).

Корректность алгоритма

Покажем, что если алгоритм завершается, то на выходе у него получается поток именно максимальной величины.

В самом деле, предположим, что в какой-то момент в слоистой сети, построенной для остаточной сети, не удалось найти блокирующий поток. Это означает, что сток t вообще не достижим в слоистой сети из истока s. Но поскольку слоистая сеть содержит в себе все кратчайшие пути из истока в остаточной сети, это в свою очередь означает, что в остаточной сети нет пути из истока в сток. Следовательно, применяя теорему Форда-Фалкерсона, получаем, что текущий поток в самом деле максимален.

Оценка числа фаз

Покажем, что алгоритм Диница всегда выполняет менее n фаз. Для этого докажем две леммы:

Лемма 1. Кратчайшее расстояние от истока до каждой вершины не уменьшается с выполнением каждой итерации, т.е.

 {\rm level}_{i+1}[v] \ge {\rm level}_i[v]

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

Доказательство. Зафиксируем произвольную фазу i и произвольную вершину v и рассмотрим любой кратчайший s-v путь P в сети G^R_{i+1} (напомним, так мы обозначаем остаточную сеть, взятую перед выполнением i+1-ой фазы). Очевидно, длина пути P равна {\rm level}_{i+1}[v].

Заметим, что в остаточную сеть G^R_{i+1} могут входить только рёбра из G^R, а также рёбра, обратные рёбрам из G^R (это следует из определения остаточной сети). Рассмотрим два случая:

  • Путь P содержит только рёбра из G^R. Тогда, понятно, длина пути P больше либо равна {\rm level}_i[v] (потому что {\rm level}_i[v] по определению — длина кратчайшего пути), что и означает выполнение неравенства.
  • Путь P содержит как минимум одно ребро, не содержащееся в G^R (но обратное какому-то ребру из G^R). Рассмотрим первое такое ребро; пусть это будет ребро (u,w).

     s \Longrightarrow u \rightarrow w \Longrightarrow[...]

    Мы можем применить нашу лемму к вершине u, потому что она подпадает под первый случай; итак, мы получаем неравенство {\rm level}_{i+1}[u] \ge {\rm level}_i[u].

    Теперь заметим, что поскольку ребро (u,w) появилось в остаточной сети только после выполнения i-ой фазы, то отсюда следует, что вдоль ребра (w,u) был дополнительно пропущен какой-то поток; следовательно, ребро (w,u) принадлежало слоистой сети перед i-ой фазой, а потому {\rm level}_i[u] = {\rm level}_i[w] + 1. Учтём, что по свойству кратчайших путей {\rm level}_{i+1}[w] = {\rm level}_{i+1}[u] + 1, и объединяя это равенство с двумя предыдущими неравенствами, получаем:

     {\rm level}_{i+1}[w] \ge {\rm level}_i[w] + 2.

    Теперь мы можем применять те же самые рассуждения ко всему оставшемуся пути до v (т.е. что каждое инвертированное ребро добавляет к \rm level как минимум два), и в итоге получим требуемое неравенство.

Лемма 2. Расстояние между истоком и стоком строго увеличивается после каждой фазы алгоритма, т.е.:

 {\rm level}^\prime[t] > {\rm level}[t],

где штрихом помечено значение, полученное на следующей фазе алгоритма.

Доказательство: от противного. Предположим, что после выполнения текущей фазы оказалось, что  {\rm level}^\prime[t] = {\rm level}[t] . Рассмотрим кратчайший путь из истока в сток; по предположению, его длина должна сохраниться неизменной. Однако остаточная сеть на следующей фазе содержит только рёбра остаточной сети перед выполнением текущей фазы, либо обратные к ним. Таким образом, пришли к противоречию: нашёлся s-t путь, который не содержит насыщенных рёбер, и имеет ту же длину, что и кратчайший путь. Этот путь должен был быть "заблокирован" блокирующим потоком, чего не произошло, в чём и заключается противоречие, что и требовалось доказать.

Эту лемму интуитивно можно понимать следующим образом: на i-ой фазе алгоритм Диница выявляет и насыщает все s-t пути длины i.

Поскольку длина кратчайшего пути из s в t не может превосходить n-1, то, следовательно, алгоритм Диница совершает не более n-1 фазы.

Поиск блокирующего потока

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

Мы рассмотрим три возможных варианта реализации поиска блокирующего потока:

  • Искать s-t пути по одному, пока такие пути находятся. Путь можно найти за O(m) обходом в глубину, а всего таких путей будет O(m) (поскольку каждый путь насыщает как минимум одно ребро). Итоговая асимптотика поиска одного блокирующего потока составит O(m^2).

  • Аналогично предыдущей идее, однако удалять в процессе обхода в глубину из графа все "лишние" рёбра, т.е. рёбра, вдоль которых не получится дойти до стока.

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

    Оценим асимптотику этого решения. Каждый обход в глубину завершается либо насыщением как минимум одного ребра (если этот обход достиг стока), либо продвижением вперёд как минимум одного указателя (в противном случае). Можно понять, что один запуск обхода в глубину из основной программы работает за O (k + n), где k — число продвижений указателей. Учитывая, что всего запусков обхода в глубину в рамках поиска одного блокирующего потока будет O (p), где p — число рёбер, насыщенных этим блокирующим потоком, то весь алгоритм поиска блокирующего потока отработает за O (p k + p n), что, учитывая, что все указатели в сумме прошли расстояние O (m), даёт асимптотику O (m + pn). В худшем случае, когда блокирующий поток насыщает все рёбра, асимптотика получается O (n m); эта асимптотика и будет использоваться далее.

    Можно сказать, что этот способ нахождения блокирующего потока чрезвычайно эффективен в том смысле, что на поиск одного увеличивающего пути он тратит O (n) операций в среднем. Именно в этом и кроется разность на целый порядок эффективностей алгоритма Диница и Эдмондса-Карпа (который ищет один увеличивающий путь за O (m)).

    Этот способ решения является по-прежнему простым для реализации, но достаточно эффективным, и потому наиболее часто применяется на практике.

  • Можно применить специальные структуры данных — динамические деревья Слетора (Sleator) и Тарьяна (Tarjan)). Тогда каждый блокирующий поток можно найти за время O (m \log n).

Асимптотика

Таким образом, весь алгоритм Диница выполняется за O (n^2 m), если блокирующий поток искать описанным выше способом за O (n m). Реализация с использованием динамических деревьев Слетора и Тарьяна будет работать за время O (n m \log n).

Единичные сети

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

Этот случай является достаточно важным, поскольку в задаче поиска максимального паросочетания построенная сеть является именно единичной.

Докажем, что на единичных сетях алгоритм Диница даже в простой реализации (которая на произвольных графах отрабатывает за O (n^2 m)) работает за время O (m \sqrt{n}), достигая на задаче поиска наибольшего паросочетания один из лучших известных алгоритмов — алгоритм Хопкрофта-Карпа.

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

Во-вторых, оценим общее количество фаз, которое могло произойти в случае единичных сетей.

Пусть уже было произведено \sqrt{n} фаз алгоритма Диница; тогда все увеличивающие пути длины не более \sqrt{n} уже обнаружены. Пусть f — текущий найденный поток, а f^* — искомый максимальный поток; рассмотрим их разность: f^* - f. Она представляет из себя поток в остаточной сети G^R. Этот поток имеет величину |f^*| - |f|, и вдоль каждого ребра равен нулю или единице. Его можно декомпозировать на набор из |f^*| - |f| путей из s в t и, возможно, циклов. Поскольку сеть единична, то все эти пути не могут иметь общих вершин, поэтому, учитывая вышесказанное, суммарное количество вершин в них cnt можно оценить как:

 cnt \ge (|f^*| - |f|) \cdot \sqrt{n}.

С другой стороны, учитывая, что cnt \le n, мы получаем отсюда:

 |f^*| - |f| \le \sqrt{n},

что означает, что ещё через \sqrt{n} фаз алгоритма Диница гарантированно найдётся максимальный поток.

Следовательно, общее число фаз алгоритма Диница, выполняемое на единичных сетях, можно оценить как 2 \sqrt{n}, что и требовалось доказать.

Реализация

Приведём две реализации алгоритма за O (n^2 m), работающие на сетях, заданных матрицами смежности и списками смежности соответственно.

Реализация над графами в виде матриц смежности

const int MAXN = ...; // число вершин
const int INF = 1000000000; // константа-бесконечность
 
int n, c[MAXN][MAXN], f[MAXN][MAXN], s, t, d[MAXN], ptr[MAXN], q[MAXN];
 
bool bfs() {
	int qh=0, qt=0;
	q[qt++] = s;
	memset (d, -1, n * sizeof d[0]);
	d[s] = 0;
	while (qh < qt) {
		int v = q[qh++];
		for (int to=0; to<n; ++to)
			if (d[to] == -1 && f[v][to] < c[v][to]) {
				q[qt++] = to;
				d[to] = d[v] + 1;
			}
	}
	return d[t] != -1;
}
 
int dfs (int v, int flow) {
	if (!flow)  return 0;
	if (v == t)  return flow;
	for (int & to=ptr[v]; to<n; ++to) {
		if (d[to] != d[v] + 1)  continue;
		int pushed = dfs (to, min (flow, c[v][to] - f[v][to]));
		if (pushed) {
			f[v][to] += pushed;
			f[to][v] -= pushed;
			return pushed;
		}
	}
	return 0;
}
 
int dinic() {
	int flow = 0;
	for (;;) {
		if (!bfs())  break;
		memset (ptr, 0, n * sizeof ptr[0]);
		while (int pushed = dfs (s, INF))
			flow += pushed;
	}
	return flow;
}

Сеть должна быть предварительно считана: должны быть заданы переменные n, s, t, а также считана матрица пропускных способностей c[][]. Основная функция решения — \rm dinic(), которая возвращает величину найденного максимального потока.

Реализация над графами в виде списков смежности

const int MAXN = ...; // число вершин
const int INF = 1000000000; // константа-бесконечность
 
struct edge {
	int a, b, cap, flow;
};
 
int n, s, t, d[MAXN], ptr[MAXN], q[MAXN];
vector<edge> e;
vector<int> g[MAXN];
 
void add_edge (int a, int b, int cap) {
	edge e1 = { a, b, cap, 0 };
	edge e2 = { b, a, 0, 0 };
	g[a].push_back ((int) e.size());
	e.push_back (e1);
	g[b].push_back ((int) e.size());
	e.push_back (e2);
}
 
bool bfs() {
	int qh=0, qt=0;
	q[qt++] = s;
	memset (d, -1, n * sizeof d[0]);
	d[s] = 0;
	while (qh < qt && d[t] == -1) {
		int v = q[qh++];
		for (size_t i=0; i<g[v].size(); ++i) {
			int id = g[v][i],
				to = e[id].b;
			if (d[to] == -1 && e[id].flow < e[id].cap) {
				q[qt++] = to;
				d[to] = d[v] + 1;
			}
		}
	}
	return d[t] != -1;
}
 
int dfs (int v, int flow) {
	if (!flow)  return 0;
	if (v == t)  return flow;
	for (; ptr[v]<(int)g[v].size(); ++ptr[v]) {
		int id = g[v][ptr[v]],
			to = e[id].b;
		if (d[to] != d[v] + 1)  continue;
		int pushed = dfs (to, min (flow, e[id].cap - e[id].flow));
		if (pushed) {
			e[id].flow += pushed;
			e[id^1].flow -= pushed;
			return pushed;
		}
	}
	return 0;
}
 
int dinic() {
	int flow = 0;
	for (;;) {
		if (!bfs())  break;
		memset (ptr, 0, n * sizeof ptr[0]);
		while (int pushed = dfs (s, INF))
			flow += pushed;
	}
	return flow;
}

Сеть должна быть предварительно считана: должны быть заданы переменные n, s, t, а также добавлены все рёбра (ориентированные) с помощью вызовов функции \rm add\_edge. Основная функция решения — \rm dinic(), которая возвращает величину найденного максимального потока.