В этом посте я расскажу о четырех основных и простых алгоритмах поиска максимального потока, к каждому алгоритму прилагается описание работы, псевдокод, нестрогое объяснение корректности, реализация на C++11, и ссылки на задачи.
I. Задача
Пусть дан граф $$$G=<V,E>$$$ и две вершины $$$s,t \in V$$$, которые в дальнейшем будем называть исток и сток, и каждому ребру из вершины $$$a$$$ в вершину $$$b$$$ сопоставлена пропускная способность $$$c_{a,b}$$$. Определим сопоставление каждому ребру $$$<a,b>$$$ числа $$$g_{a,b}$$$. Назовем это сопоставление потоком, если для каждого ребра $$$<a,b>$$$ верны два утверждения:
1) $$$0 \le g_{a,b} \le c_{a,b}$$$
2) $$$ \sum_{i \in V} g_{i,a} = \sum_{i \in V} g_{a,i} \forall a \in V \setminus \{ s,t \} $$$
Величина потока равна $$$ \sum_{i \in V} g_{s,i} = \sum_{i \in V} g_{i,t} = |F|$$$. Задача — вычислить $$$|F|$$$ и все $$$g_{a,b}$$$.
II. Метод проталкивания потока
a) Алгоритм Форда-Фалкерсона
Пусть дан граф с пропускными способностями. Будем хранить его пропускные способности в матрице смежности $$$c[][]$$$, если ребра $$$<a,b>$$$ нет в графе, $$$c[a][b]=0$$$. Определим операцию проталкивания потока: пусть мы выбрали путь $$$p_{1},p_{2}, \dots ,p_{k}$$$ из истока $$$s$$$ в сток $$$t$$$ такой, что пропускные способности между соседними вершинами положительны: $$$c_{p_{i-1},p_{i}} > 0 \forall i \in [2,k]$$$. Пусть $$$AugFlow$$$ — минимальная пропускная способность на ребре на пути: $$$AugFlow = min_{i \in [2,k]} c_{p_{i-1},p_{i}}$$$. Пройдем по этому пути и вычтем из пропускных способностей $$$AugFlow$$$ по направлению к стоку и добавим против стока. Итак, если у нас есть путь ненулевой пропускной способности, мы можем увеличить поток. Можно ли, вычисляя такие пути, всегда найти максимальный поток? Работает следующая теорема:
Теорема 1.1 Для того, чтобы поток был максимальным, необходимо и достаточно, чтобы отсутствовал путь из истока в сток ненулевой пропускной способности.
Алгоритмпоток:=0;
пока истина:
найти любой путь из истока в сток p[];
если путь не нашелся:
выйти из цикла;
иначе:
дополняющий_поток:=min(c[p[i-1]][p[i]]);
поток+=дополняющий_поток;
для всех i:=1 до p.size-1:
c[p[i-1]][p[i]]-=дополняющий_поток;
c[p[i]][p[i-1]]+=дополняющий_поток;
Асимптотика Если пропускные способности целочисленные, то дополняющий поток всегда не меньше единицы. Тогда, так как асимптотика обхода графа $$$O(E)$$$, асимптотика алгоритма составляет $$$O(E|F|)$$$. К сожалению, это не так для графов с вещественными пропускными способностями. Но если использовать обход графа в ширину, действует следующая теорема:
Теорема 1.2 В течении алгоритма с обходом в ширину расстояние от истока до любой вершины по ребрам ненулевой пропускной способности не уменьшается.
Это дает еще одну асимптотику: $$$O(VE^{2})$$$. Данный алгоритм называется алгоритмом Эдмондса-Карпа. Его реализация здесь дана.
Реализация/**
* Дано количество вершин и ребер, исток и сток, затем ребра с их пропускными способностями.
* Везде 0-индексация. Необходимо вывести максимальный поток и матрицу остаточных пропускных способностей.
*/
#include <stdio.h>
#include <vector>
#include <queue>
using namespace std;
const int inf=1e9;
int n,m,s,t;//Количество вершин и ребер, исток и сток
int c[1000][1000];//Матрица пропускных способностей
int main(){
scanf("%d %d %d %d",&n,&m,&s,&t);
for(int i=0;i<m;i++){
int a,b,cap;
scanf("%d %d %d",&a,&b,&cap);
c[a][b]=cap;
}
int MaxFlow=0;//Искомый максимальный поток
while(true){
/**< Поиск в ширину */
vector <int> parent(n,-1);
vector <bool> used(n,false);
queue <int> q;
used[s]=true;
q.push(s);
while(!q.empty()){
int v=q.front();
q.pop();
for(int i=0;i<n;i++){
if(!used[i] && c[v][i]>0){
parent[i]=v;
used[i]=true;
q.push(i);
}
}
}
if(!used[t])//Не дошли до стока - поток уже максимальный
break;
int AugFlow=inf;//Дополнительный поток
/**< Бежим по пути и ищем ребро с минимальной пропускной способностью */
int ptr=t;
while(ptr!=s){
AugFlow=min(AugFlow,c[parent[ptr]][ptr]);
ptr=parent[ptr];
}
/**< Изменяем пропускные способности */
ptr=t;
while(ptr!=s){
c[parent[ptr]][ptr]-=AugFlow;//По пути вычитаем
c[ptr][parent[ptr]]+=AugFlow;//Против пути добавляем
ptr=parent[ptr];
}
MaxFlow+=AugFlow;
}
/**< Выводим ответ */
printf("%d\n",MaxFlow);
for(int i=0;i<n;i++){
for(int j=0;j<n;j++){
printf("%3d ",c[i][j]);
}
printf("\n");
}
return 0;
}
Задачи 1 2 3
b) Алгоритм Диница
Алгоритм Форда-Фалкерсона можно ускорить, использую следующую идею. Выполним поиск в ширину из истока и, игнорируя ребра с нулевой пропускной способностью, вычислим расстояния $$$d[]$$$ до каждой вершины. Будем говорить, что мы построили слоистую сеть. Затем будем искать такие дополняющие пути $$$p[]$$$, что $$$d[p[1]]-d[p[0]]=d[p[2]]-d[p[1]]= \dots =1$$$. Все эти пути можно найти довольно быстро — за $$$O(VE)$$$: очевидно, что в силу ограничения, наложенного на путь, если по ребру $$$<a,b>$$$ не удалось протолкнуть путь, в течении фазы по нему нельзя будет протолкнуть путь. Поэтому будем поддерживать у каждой вершины первую нерассмотренную вершину для каждой фазы, и искать следующие пути в этой фазе, начиная с первой нерассмотренной вершины и увеличивать это значение, если путь с этой вершины не нашелся. После того, как мы нашли все пути, мы приступаем к следующей фазе. При этом работает следующая теорема:
Теорема 2.1 Расстояние до стока после каждой фазы строго возрастает.
Тогда, так как путь из истока в сток не больше $$$V$$$, количество фаз $$$O(V)$$$.
Алгоритмc[][]
left[];
d[];
функция решить():
поток:=0;
пока истина:
вычислить расстояния до каждой вершины d[];
если пути до стока нет:
выйти из цикла;
заполнить(left[],V,0);
пока истина:
дополняющий_поток=протолкнуть(s,inf);
если дополняющий_поток==0:
выйти из цикла;
иначе:
поток+=дополняющий_поток;
функция протолкнуть(вершина,поток):
если вершина==сток:
вернуть поток;
пока left[вершина]<V:
если d[left[вершина]]+1==d[вершина] и c[вершина][left[вершина]]>0:
дополняющий_поток=протолкнуть(left[вершина],min(поток,c[вершина][left[вершина]]));
если дополняющий_поток>0:
c[вершина][left[вершина]]-=дополняющий_поток;
c[left[вершина]][вершина]+=дополняющий_поток;
вернуть дополняющий_поток;
вернуть 0;
Асимптотика $$$O(V)$$$ фаз и $$$O(VE)$$$ на каждую фазу — итого $$$O(V^{2}E)$$$.
Частный случай использования этого алгоритма для поиска максимального паросочетания в двудольном графе называется алгоритмом Хопкрофта-Карпа. Он имеет более сильную асимптотику $$$O(\sqrt{V}E)$$$.
Реализация/**
* Дано количество вершин и ребер, исток и сток, затем ребра с их пропускными способностями.
* Везде 0-индексация. Необходимо вывести максимальный поток и матрицу остаточных пропускных способностей.
*/
#include <stdio.h>
#include <vector>
#include <queue>
using namespace std;
const int inf=1e9;
int n,m,s,t;//Количество вершин и ребер, исток и сток
int c[1000][1000];//Матрица пропускных способностей
int d[1000];//Расстояния
int left[1000];//Самый левый нерассмотренный сосед
/**< Построить слоистую сеть обходом в ширину */
void buildNet(){
for(int i=0;i<n;i++)
d[i]=-1;//Пути до вершины нет
queue <int> q;
d[s]=0;
q.push(s);
while(!q.empty()){
int v=q.front();
q.pop();
for(int i=0;i<n;i++){
if(d[i]==-1 && c[v][i]>0){
d[i]=d[v]+1;
q.push(i);
}
}
}
}
/**< Протолкнуть поток обходом в глубину */
int push(int v,int flow){
/**< Уже прошли какой-то путь из s в v, минимальная пропускная способность - flow */
if(v==t){//Если зашли в сток, нашли дополняющий путь, возвращаем минимальную пропускную способность на пути, завершаем обход
return flow;
}
for(;left[v]<n;left[v]++){//Перебираем вершины не от начала (без этого алгоритм деградирует!)
int to=left[v];
if(d[to]==d[v]+1 && c[v][to]>0){//Если вершина на следующем слое
int AugFlow=push(to,min(flow,c[v][to]));//Проталкиваем от нее, путь увеличился на ребро c[v][to], поэтому flow возможно уменьшился
if(AugFlow>0){//Если протолкнули - завершаем обход
c[v][to]-=AugFlow;//По пути вычитаем
c[to][v]+=AugFlow;//Против пути добавляем
return AugFlow;
}
}
}
return 0;
}
int main(){
scanf("%d %d %d %d",&n,&m,&s,&t);
for(int i=0;i<m;i++){
int a,b,cap;
scanf("%d %d %d",&a,&b,&cap);
c[a][b]=cap;
}
int MaxFlow=0;//Искомый максимальный поток
while(true){
buildNet();//Строим новую слоистую сеть
if(d[t]==-1)//Если пути до стока нет, поток максимальный
break;
for(int i=0;i<n;i++)//Начинаем новые проталкивания
left[i]=0;
while(true){
int AugFlow=push(0,inf);//Проталкиваем из истока, пока минимальная пропускная способность на пути отсутствует
if(AugFlow==0)//Если не протолкнули, строим новую слоистую сеть
break;
else
MaxFlow+=AugFlow;//Добавляем дополняющий поток
}
}
/**< Выводим ответ */
printf("%d\n",MaxFlow);
for(int i=0;i<n;i++){
for(int j=0;j<n;j++){
printf("%3d ",c[i][j]);
}
printf("\n");
}
return 0;
}
Задачи 4 5
III. Метод проталкивания предпотока
a) Разрядить, пока есть что разрядить
Назовем предпотоком сопоставление каждому ребру в графе числа, при котором выполняется:
1) $$$0 \le g_{a,b} \le c_{a,b}$$$
2) $$$ \sum_{i \in V} g_{i,a} \ge \sum_{i \in V} g_{a,i} \forall a \in V \setminus \{ s,t \} $$$
Определим переполнение $$$e_{i}$$$ вершины как значение $$$ \sum_{i \in V} g_{i,a} - \sum_{i \in V} g_{a,i} $$$, и высоту $$$h_{i}$$$ вершины. Будем выполнять две операции:
1) Проталкивание из вершины $$$a$$$ в вершину $$$b$$$. Допустимо тогда, когда $$$h_{a}=h_{b}+1$$$. Проталкиваем столько потока, на сколько переполнена $$$a$$$, но не более, чем может вместить ребро $$$<a,b>$$$: $$$pushed=min(e_{a},c_{a,b})$$$.
2) Поднятие вершины $$$a$$$: допустимо, если $$$a$$$ переполнена, но нет такой вершины $$$b$$$, что $$$h_{a}=h_{b}+1$$$.
Алгоритм заключается в использовании этих операций в любом порядке, пока хотя бы одна из них допустима. Можно для удобства определить операцию разрядить вершину, которая выполняет проталкивания из нее и поднятия, пока она переполнена (это не повлияет на асимптотику).
Перед основным циклом необходимо выполнить следующую инициализацию: проталкиваем тривиально из истока во все соседние вершины.
Алгоритмc[][];
e[];
h[];
функция решить():
e[s]=inf;
h[s]=n;
для всех i:=1 до n:
h[i]:=0;
если c[s][i]>0:
протолкнуть(s,i);
пока есть переполненная вершина i, которая не исток или сток:
разрядить(i);
функция протолкнуть(откуда,куда):
величина:=min(e[откуда],c[откуда][куда]);
e[откуда]-=величина;
e[куда]+=величина;
c[откуда][куда]-=величина;
c[откуда][куда]+=величина;
функция поднять(вершина):
h[вершина]:=1+min(h[i]) для такой i что c[вершина][i]>0, если такой i существует;
функция разрядить(вершина):
пока вершина переполнена:
для всех i:=1 до n:
если c[вершина][i]>0 и h[вершина]==h[i]+1:
протолкнуть(вершина,i);
поднять(вершина);
Теорема 3.1 Подъем вершины строго увеличивает ее высоту.
Теорема 3.2 Если ни одно действие недопустимо — предпоток есть максимальный поток.
Асимптотика Высота каждой вершины во время работы алгоритма не превышает $$$2V$$$, откуда количество вызовов операции поднять есть $$$O(V^{2})$$$, а количество вызовов операции протолкнуть есть $$$O(V^{2}E)$$$. Итоговая асимптотика — $$$O(V^{2}E)$$$.
Реализация/**
* Дано количество вершин и ребер, исток и сток, затем ребра с их пропускными способностями.
* Везде 0-индексация. Необходимо вывести максимальный поток и матрицу остаточных пропускных способностей.
*/
#include <stdio.h>
#include <vector>
#include <queue>
using namespace std;
const int inf=1e9;
int n,m,s,t;//Количество вершин и ребер, исток и сток
int c[1000][1000];//Матрица пропускных способностей
int e[1000],h[1000];//Переполнение и высота
int MaxFlow=0;//Максимальный поток
/**< Протолкнуть по ребру */
int push(int from,int to){
int pushed=min(e[from],c[from][to]);//Проталкиваем на столько, на сколько переполнен from, но не более, чем позволит протолкнуть ребро
if(to==t)//Если протолкнули в сток, добавляем в ответ
MaxFlow+=pushed;
e[from]-=pushed;
e[to]+=pushed;
c[from][to]-=pushed;
c[to][from]+=pushed;
return pushed;
}
/**< Поднять вершину */
void lift(int v){
/**< Ищем вершину с минимальной высотой и устанавливаем высоту на единицу больше, если таковая есть */
int New=inf;
for(int i=0;i<n;i++){
if(c[v][i]>0){
New=min(New,h[i]);
}
}
if(New!=inf)
h[v]=New+1;
}
/**< Разрядить */
void discharge(int v){
/**< Пока вершина переполнена, проталкиваем во все вершины, затем поднимаем */
while(e[v]>0){
for(int i=0;i<n;i++){
if(c[v][i]>0 && h[v]==h[i]+1){
push(v,i);
}
}
lift(v);
}
}
int main(){
scanf("%d %d %d %d",&n,&m,&s,&t);
for(int i=0;i<m;i++){
int a,b,cap;
scanf("%d %d %d",&a,&b,&cap);
c[a][b]=cap;
}
/**< Начальная инициализация */
e[s]=inf;
h[s]=n;
for(int i=0;i<n;i++){
if(c[s][i]>0){
push(s,i);
}
}
while(true){
bool good=false;//Есть вершины, которые нужно разрядить
for(int i=0;i<n;i++){
if(i!=s && i!=t && e[i]>0){
discharge(i);
good=true;
}
}
if(!good)//Все вершины разряжены - предпоток есть максимальный поток
break;
}
/**< Выводим ответ */
printf("%d\n",MaxFlow);
for(int i=0;i<n;i++){
for(int j=0;j<n;j++){
printf("%3d ",c[i][j]);
}
printf("\n");
}
return 0;
}
Задачи 4 5
Те же, что и для Диница
б) Разрядить в очереди
Краткая идея: будем хранить все переполненные вершины в очереди, и разряжать первую вершину.
Алгоритмc[][];
e[];
h[];
очередь;
функция решить():
e[s]=inf;
h[s]=n;
для всех i:=1 до n:
h[i]:=0;
если c[s][i]>0:
протолкнуть(s,i);
пока очередь не пуста:
разрядить(очередь.начало);
очередь.удалить_начало;
функция протолкнуть(откуда,куда):
величина:=min(e[откуда],c[откуда][куда]);
если e[куда]==0 и куда не исток или сток:
положить куда в очередь;
e[откуда]-=величина;
e[куда]+=величина;
c[откуда][куда]-=величина;
c[откуда][куда]+=величина;
функция поднять(вершина):
h[вершина]:=1+min(h[i]) для такой i что c[вершина][i]>0, если такой i существует;
функция разрядить(вершина):
пока вершина переполнена:
для всех i:=1 до n:
если c[вершина][i]>0 и h[вершина]==h[i]+1:
протолкнуть(вершина,i);
поднять(вершина);
Асимптотика $$$O(V^{3})$$$.
Реализация/**
* Дано количество вершин и ребер, исток и сток, затем ребра с их пропускными способностями.
* Везде 0-индексация. Необходимо вывести максимальный поток и матрицу остаточных пропускных способностей.
*/
#include <stdio.h>
#include <vector>
#include <queue>
using namespace std;
const int inf=1e9;
int n,m,s,t;//Количество вершин и ребер, исток и сток
int c[1000][1000];//Матрица пропускных способностей
int e[1000],h[1000];//Переполнение и высота
int MaxFlow=0;//Максимальный поток
queue <int> q;//Очередь
/**< Протолкнуть по ребру */
int push(int from,int to){
int pushed=min(e[from],c[from][to]);//Проталкиваем на столько, на сколько переполнен from, но не более, чем позволит протолкнуть ребро
if(to==t)//Если протолкнули в сток, добавляем в ответ
MaxFlow+=pushed;
if(e[to]==0 && to!=t)//Вершина to не была переполненной, а теперь переполненна
q.push(to);
e[from]-=pushed;
e[to]+=pushed;
c[from][to]-=pushed;
c[to][from]+=pushed;
return pushed;
}
/**< Поднять вершину */
void lift(int v){
/**< Ищем вершину с минимальной высотой и устанавливаем высоту на единицу больше, если таковая есть */
int New=inf;
for(int i=0;i<n;i++){
if(c[v][i]>0){
New=min(New,h[i]);
}
}
if(New!=inf)
h[v]=New+1;
}
/**< Разрядить */
void discharge(int v){
/**< Пока вершина переполнена, проталкиваем во все вершины, затем поднимаем */
while(e[v]>0){
for(int i=0;i<n;i++){
if(c[v][i]>0 && h[v]==h[i]+1){
push(v,i);
}
}
lift(v);
}
}
int main(){
scanf("%d %d %d %d",&n,&m,&s,&t);
for(int i=0;i<m;i++){
int a,b,cap;
scanf("%d %d %d",&a,&b,&cap);
c[a][b]=cap;
}
/**< Начальная инициализация */
e[s]=inf;
h[s]=n;
for(int i=0;i<n;i++){
if(c[s][i]>0){
push(s,i);
}
}
/**< Пока есть переполненные вершины - разрядить */
while(!q.empty()){
int v=q.front();
q.pop();
discharge(v);
}
/**< Выводим ответ */
printf("%d\n",MaxFlow);
for(int i=0;i<n;i++){
for(int j=0;j<n;j++){
printf("%3d ",c[i][j]);
}
printf("\n");
}
return 0;
}
Задачи 6
К сожалению, к данной задаче принимается Диниц, а не должен.