Рандомизированный распределенный адаптивный алгоритм решения задачи о максимальном потоке Текст научной статьи по специальности «Математика»

Компьютерные инструменты в образовании, 2016 № 5: 46-62 УДК: 519.178 http://ipo.spb.ru/journal

РАНДОМИЗИРОВАННЫЙ РАСПРЕДЕЛЕННЫЙ АДАПТИВНЫЙ АЛГОРИТМ РЕШЕНИЯ ЗАДАЧИ О МАКСИМАЛЬНОМ ПОТОКЕ*

Мальковский Н.В.1 1СПбГУ, Санкт-Петербург, Россия.

Аннотация

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

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

Цитирование: Мальковский Н.В. Рандомизированный распределенный адаптивный алгоритм решения задачи о максимальном потоке // Компьютерные инструменты в образовании, 2016. № 5. С. 46-62 .

1. ВВЕДЕНИЕ

Задача о максимальном потоке является одной из классических задач оптимизации. Она имеет многие приложения в сетевых системах, а также часто возникает как техническая подзадача в других задачах оптимизации. Задача изучается уже более полувека, наиболее свежий обзор алгоритмов её решения можно найти в [1]. В этой работе представлен метод решения задачи о максимальном потоке, не являющийся самым эффективным по времени работы, однако обладающий некоторыми прочими полезными свойствами: возможностью распределенной / мультиагентной реализации (то есть если задача возникает в распределенной сети, при этом вычислительные устройства этой сети образуют вершины графа, а каналы связи образуют дуги графа, то алгоритм имеет

*Работа выполнена при поддержке гранта СПбГУ 6.37.181.2014.

«естественную» реализацию на узлах этой сети) и адаптивности (если пропускная способность дуг меняется по мере вычислений, то алгоритм сохраняет близость к максимальному потоку при условии, что эти изменения ограничены и достаточно малы). За основу взят метод балансирования дуг, описанный в [2]: оказывается, этот метод можно представить как проективный градиентный спуск [3]. В рамках этой работы будут исследованы две вариации метода: синхронный (обычный проективный градиентный спуск) и рандомизированный асинхронный (специальный вид рандомизированного алгоритма стохастической аппроксимации [4]). Второй вариант больше походит на метод, описанный в [2], с той лишь разницей, что используется случайный выбор дуг.

Основной целью этой статьи является разработка специфических алгоритмов решения задачи о максимальном потоке. Стоит отметить, что результаты работы также обобщают идеи проективного градиентного спуска [3], рандомизированного алгоритма стохастической аппроксимации [4, 5], алгоритмов «сплетен» [6] и протоколов консенсуса [7]. По сравнению с классическими результатами для проективного градиентного спуска, в этой работе получена экспоненциальная скорость сходимости для квадратичной функции, не являющейся сильно выпуклой. По сравнению с рандомизированными алгоритмами стохастической аппроксимации, рассматриваемыми в [4], в этой работе решается задача с ограничениями в виде неравенств. Алгоритмы консенсуса по общему стилю схожи с синхронной версией разработанного алгоритма, а алгоритмы сплетен — с рандомизированной версией. Связь между разработанными методами и сплетнями / консенсусом состоит в следующем: консенсус и сплетни минимизируют функционал утВВту, в то время как разработанный алгоритм минимизирует функционал хтВтВх (в обоих случаях В — матрица инцидентности некоторого графа). Ключевым различием является наличие дополнительных ограничений в задаче, заданных неравенствами. Лемма 3 является ключевой для получения оценок сходимости, аналогичных [6, 7].

2. ЗАДАЧА О МАКСИМА1ЬНОМ ПОТОКЕ

Пусть задан некоторый ориентированный граф О = {V,Е) с выделенными вершинами истока з и стока Ь. Традиционно задача о максимальном потоке ставится в следующем виде:

максимизировать val(x) = £eein(l) xe eEout(i) xe,

T.eein(i) xe -Y.eeout(i) xe = 0, i £ V, i Ф S, t, (1)

0 ^ Xe ^ Ce,

при условии

где переменными являются величины потока по дугам xi,...,xm (здесь и далее x = (x1,...,xm)T), in(i) — множество входящих в i дуг, out(i) — множество дуг, выходящих из i. Первое условие означает, что поток не задерживается в промежуточных вершинах. Единственные две вершины с неотрицательной разностью входящих / исходящих потоков — s и t. Более того,

xe - xe = - xe - x

eein(t) eeout(t) \eein(s) eeout(s)

Задачу (1) можно записать в более компактной форме, используя матрицу инцидентности B:

{1, e выходит из i, -1, e входит в i, 0 в остальных случаях

Обозначив вектор %st: %sst = -1, %stt = 1, xf = i Ф t, t задачу (1) перепишем в виде:

максимизировать v,

i Bx = vxst, при условии < o

[ 0 ^ xe ^ ce,

3. ОБЩАЯ СХЕМА МЕТОДА БАЛАНСИРОВАНИЯ ДУГ

Метод балансирования дуг состоит из инициализации и последовательном применении некоторых локальных изменений потока:

Function Initialized, c) for (i, j) e E do

if i = t или j = t then

Lxijcij;

_ Xij ^ 0; return x;

Function Move(G,c,x, e)

Xe - 1[BTBx]e

xe

if xe > ce then

xe ce

L

if xe < 0 then

|_ xe ^ 0;

Function ArcBalancing(G, c, stopping criterion) x — Initialize(G, c);

while stopping criterion is not satisfied do pick e = (i, j) e E: i, j + s, t; Move(G, c, x, e);

return x;

Величину [Bx]; = £eein{i) xe ee0ut(i) xe обычно принято называть избытком вершины i. Величина [BTBx];,j) = [Bx]j - [Bx]; есть разница избытков вершин i, j. Отсюда название алгоритма: функция Move пытается изменить значение потока на дуге e так, чтобы разница между избытками на концах дуги была минимальной при условии, что 0 ^ xe ^ ce. При этом от способа выбора дуги зависит сходимость функции ArcBalancing. Обычно выбираются только те дуги, для которых функция Move изменит значение потока, однако если выбирать дуги таким образом, а останавливаться в случае, когда таких дуг нет, то нет никакой гарантии, что ArcBalancing закончит работу.

Одна из основных идей, предложенных в этой работе, заключается в том, что при равномерном выборе дуг для операции Move возможно доказать асимптотическую сходимость, при этом можно достаточно легко извлечь максимальный поток / минимальный разрез из предельного состояния. Далее в работе будут проанализированы следующие варианты метода ArcBalancing:RandomizedArcBalancing и ConcurrentArcBalancing.

Algorithm 1: RandomizedArcBalancing(G, c, stopping criterion) x — Initialize(G, c);

while stopping criterion is not satisfied do

pick random e = (i, j) e E': i, j + s, t w.p. pi; Move(G, c, x, e);

return x;

Algorithm 2: ConcurrentArcBalancing(G, c, stopping criterion)

x — 0m;

while stopping criterion is not satisfied do

x — x - 2d(B'TB'x + B'Tb);

for e e E' do

if xe > ce then

xe — ce;

if xe < 0 then

|_ xe — 0;

return x;

В описании этих методов используются обозначения E' = {(г, j) е E | г Ф х, j Ф Ь}, d- = |{(г, j) е Е'}|, d+ = ^^, г) е Е'}|, d > 0 — некоторый скалярный параметр, В' — мат-

, е' соответствует дугам

рица инцидентности графа О' = {V, Е'). Положим, что с = из Е', с" соответствует дугам из Е \ Е', Ь = В "с".

4. СХОДИМОСТЬ МЕТОДОВ

Для начала отметим, что оба метода очень схожи с проективным градиентным спуском, более того, если обозначить /(х) = ^(2хтВ'тВ'х+хтВ'тЬ), V/(х) = ^(В'тВ'х+В'тЬ), то алгоритм СопсиггеШАгсБаЫпащ в точности является проективным градиентным спуском для функции /. Общий метод проективного градиентного спуска был впервые сформулирован в работе [3]. Метод RandomizedArcBalancing обновляет поток похожим образом: на каждом шаге случайным образом выбирается координата и делается шаг градиентного спуска только по этой координате. Прежде чем переходить к анализу сходимости, выпишем некоторые ключевые свойства задачи.

Лемма 1. Если Б — диагональная матрица с неотрицательными элементами, а В е ипхт — произвольная вещественная матрица, то матрица БВтВ имеет только вещественные неотрицательные собственные числа. Более того, если все диагональные элементы Б строго положительны, то

БВтВх = 0т о Вх = 0П.

Доказательство. Так как Б — диагональная матрица с неотрицательными элементами, то БВтВ имеет только вещественные собственные числа, для матрицы Б допустимо разложение

Б = \Тбт \[Б,

где у/Б = diag {%/Д[ъ • ••>%/ Бтт}. С другой стороны, у матриц ББТБ и БББТ одинаковый набор ненулевых собственных чисел:

Ах = ББТБ х ^ БББТ (Б х) = Б (ББТБ х) = Б (Ах) = ЦБ х),

а значит, при Бх Ф 0п А является собственным числом БББТ. Если же Бх = 0п, то ББТБх=0 т, а значит, А = 0. Аналогично в обратную сторону:

АхТ = хТБББТ ^ (хТБ)ББТБ = (хБББТ)Б = (АхТ)Б = А(хТБ),

при хТБ Ф 0п А является собственным числом ББТБ, при хТБ = 0п хТБББТ = 0т, а значит, А = 0. В силу

хТБББТх = \\^ББТх\\2 & 0,

получаем, что ББТБ имеет только неотрицательные собственные числа. Наконец, если все диагональные элементы Б строго положительны, то Б обратима, а значит,

ББТБх = 0т ^ 0 = хБ-1 ББТБх = \\Бх\\2 о Бх = 0п•

Лемма 2. Если С — граф без петель и кратных дуг, Бее * тах{1/2^- + d+ ),1/(ndj )} при е = а, ]), то а(ББТБ') * 2.

Доказательство. Для е = Ц, ]) обозначим de = dij = Бее, dij = 0 при Ц, ]) £ ЕТак как у атри ББ 'ТБ' и Б'ББ'Т одинаковый набор ненулевых собственных чисел, то у них одинаковый спектральный радиус. Непосредственным вычислением получаем, что

[Б 'ББ 'Т

].. = \ ~(dij + dji), г Ф j, lj \ T.kФidik + ^> i = j•

гк + '

Если dij * 1/2(dj + d++), то

X+ ^) * I -1 +) + I -1 +) *

jФ1 а,j)еЕ' 2(а1 + ) (Л)еЕ + а1 )

V — V — <

а,j)еЕ 2ai (jЛ)еЕ 2ai

Последнее неравенство обращается в равенство, если обе суммы содержат хотя бы одно слагаемое. Если же de * 1/(ndj), то

^ + £ 1 ^ jФ1 а,j)еЕ"(Л)еЕ ndj

I Щ + dji) * I — + I — * - + * 1

па па. п п

Как и в предыдущем случае, предпоследнее неравенство переходит в равенство, если соответствующие суммы содержат хотя бы одно слагаемое. В последнем неравенстве мы воспользовались тем, что в графе без петель и кратных дуг У? : d+ * п - 1. Наконец, в силу критерия Гершгорина, мы получаем, что все собственные числа матрицы Б 'ББ 'Т удовлетворяют соотношению \А - тах.Е.ф(dij + dji)\ * max^jФi(dij + dji), из чего следует, что \А -1\ * 1, а значит, а(Б,Т Б 'Б) < 2. □

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

Лемма 3. Пусть д = {хе Кп | Ах < ф, Рд(х) = аг^туед ||х-у||, Ь1;...,Ьп — некоторый ортогональный базис, х, у е Кп, а 1,..., ап и в1,..., вп — разложения х, у по базису Ь, а^,..., аП и в1,..., вп — разложения Рд (х), Рд (у) по базису Ь, тогда

а - в;-1<|а,- - в11,1 < I < п.

Доказательство. Во-первых, отметим, что Рд представляется в виде композиции проекций Рд1 о... о Рдт, где Qi = {х е Шп | х < }, поэтому достаточно доказать указанное свойство для полупространства. Пусть д = {х е Кп | аТх < й}. Рассмотрим некоторый ортогональный базис, содержащий вектор а: а = а1; а2,..., ап, g1,..., gn и Н1,..., кп — разложения х,у по этому базису. Очевидно, что разложение по этому базису для Рд(х), Рд(у) задается формулами:

^ = т*пI к1 = т1пК пар)

g'i = gi, i > 1 1 Ц = кг, i > 1

з чего следует

^ 1- 1 - к1| gг• - к = gi - к, г > 1.

Так как Ь — ортогональный базис, то

а

||Ьг||= хгЬ = £ gjа] Ь.

1=1

Остальные коэффициенты получаются аналогичным образом. Наконец,

К - в\ | =

|| Ь г || - к1||аТ Ьг | 1

пп

£ g1 а]Ь - а]Ьг

1=1 1=1

|| Ь г

|| Ь г ||

£ gj а] Ьг- £ а] Ьг

lg1 - к1||а]Ьг| < ||Ь г || ^

= |аг -в11.

1=1

1=1

Далее мы воспользуемся следующим очевидным следствием этой леммы: если М

у + уу

х, иу, и^, иу е М, ух, Уу, уХ, уу е М1, то

линейное подпространство Кп, х = их + ух, у = иу + у у, Рд (х) = иХ + уХ , Рд (у) = и'у + Уу,

||иХ - иу||<||их - иу||,

||уХ - уу ||<||Ух - Уу ||.

4.1. Сходимость синхронного алгоритма

Во всей этой секции мы будем считать, что матрица || 1В'ТВ || ^ 2 (параметр й можно выбирать как угодно, поэтому достаточно легко гарантировать это свойство, выбрав й, используя, например, лемму (2)). Как и раньше, обозначим f (х) = хТВ'ТВ'х + хТВ'ТЬ),

1

V / (х) = 1(Б 'Т Б 'х + Б'ТЬ). Метод СопсиггеШАгсБа1апа^ генерирует последовательность {хк0 по правилу

хк+! = PQ |хк - 2 VГ(хк)), (2)

при этом начальное приближение выбирается как х0 = 0т. Если бы матрица Б,ТБ' была строго положительно определенной, то можно было бы вывести экспоненциальную сходимость из результатов [3]. К сожалению, это происходит только в случае, если С' не содержит циклов (без учета ориентации). Тем не менее, для квадратичных функций также можно доказать экспоненциальную сходимость. Обозначим А = 1Б'ТБ', ш = 1Б,ТЬ.

Теорема 1. Пусть хк — последовательность, генерируемая по правилу (2) с любым начальным приближением х0 е Q, \\ А\\ * 2, у = тахх^ \\ Ах + ш|\, А2(А) — минимальное отличное от нуля собственное число матрицы А, ц = [1 - ЩА), тогда

1. Последовательность хк сходится к некоторой точке минимума х* / на Q, причем

цк/2

- у/Щ

цк/2 ,-

\\хк - х*\\<^— л/г\\х0 - х* !!• 1-х Гц у

2. /(хк) убывает, причем

/(хк) - /(х*) < уцк\\х0 -х*!\,

3. Если П — ядро А, то есть П = {у е Кт \ Ау = 0т}, то

б(хк -х*,П) < цк\\х0 -х*\\, где 6(-, ■) — евклидова метрика. Доказательство.

/(хк+:) - Г(хк) = 1 хк+1ТАхк+! + шТхк+1 - 1 хкТАхк - шТхк = 22

!(хк+1 -хк)ТА(хк+: -хк) + (Ахк + ш)Т(хк+: -хк) *

(Ахк + ш)Т (хк+1 - хк) + \\хк+1 - хк \\2 =

, |хк +1 -хк \ \ 2 - ^ - 2(Ахк + -хк+.,Т^ -хк, - 2,хк+. -хк,Т^ -хк, *

- \ \хк+1 - хк \ \ г.

В последнем равенстве использовалось свойство проекции Ух е Кп, у е Q : (х - PQ(х)) (у - PQ(х)) * 0. Таким образом, /(хк) убывает. Так как / ограничена снизу, то существует предел /(хк), а значит, \ \ хк+1 -хк \ \ 2 ^ 0. Пусть х* - любая точка минимума / на Q, тогда, повторив аналогичное рассуждение для х* и х' = PQ(х - 2 V/(х*)), получаем

0 * /(х') - /(х*) * - \\ х*-х' \ \ 2 ^ х*=

Далее, рассмотрим последовательности ик е П, ук е П1: хк = ик + ук (эти последовательности задаются однозначно по хк). Отметим, что /(хк) = /(ук):

/(хк) - /(ук) = 1 (хк - ук) ТА(хк + ук) + утт(хк-ук\ = 0^ =икА=0т =0 по лемме 1

г

Пусть и* е П,у* е П1, и* + у* = х*, тогда

/(хк) - f (х*) = /(ук) - f (у*) = 2(Ук - у*)тА(ук + у*) + шт(ук - у*) ^

1 А(ук + у*) + ш

ук у*

2

В последнем выражении левый сомножитель ограничен на ((, в силу выбора у, выполняется неравенство

/(хк) - /(х*) ^ Г||ук -у*||.

Далее, применяя лемму 3 к точкам хк - 1V/(хк) и х* - 2V/(х*) и учитывая ик,и* е П,

(1 _ )|

ук, у* е П1 получаем

||ук+1 -у*||^||(7- 1 А)(ук -у*)||^ [1 - ^^ ||ук -у*||.

Таким образом,

и

6(хк -х*,П) ^ дк||х0 -х*||

/(хк) - /(х*) ^ Тдк||х0 -х*||.

Теперь заметим следующее: раз f (хк) ^ /(х*), то £то0 1|хг+1 -хг ||2 ^ ^^0 f (хг) - /(хг+1) = f (х0) -/(х*). Для сходимости хк достаточно сходимости ряда £ТО0 ||хг+1 -хг ||, что в общем

случае не следует из сходимости £ТО0 ||хг+1 -хг ||2. Однако в нашем случае сходимость

г к

экспоненциальная, что позволяет вывести сходимость хк:

||хк -хк+1||2 ^ /(хк) -/(хк+1) < /(хк) -/(х*),

следовательно,

||хк -хк+1|| ^ qkl2 ^у||х0 -х*||,

а значит, хк сходится. Не умаляя общности, можно считать, что хк ^ х*. Наконец, выве-

к

дем скорость сходимости для хк:

то

1|хк -х*||^£||хг' -хг'+1||^

г=к

то /- ^к/2 I-

^й112 Йк/Чг||х0 -х* || = ^ Vг||х0 -х* ||.

\г=0 1 Vй

Замечания. Для процесса хк+1 = хк - 1 (Ахк + ш) легко показать, что и* = и0 = и1 = ..., что облегчает доказательство, однако при наличии проекции это свойство не выполняется, для чего и понадобилась лемма 3. Величину Л2(А) принято называть числом Фидле-ра или алгебраической связностью, оно всегда фигурирует в качестве параметра скорости сходимости лапласиановых систем. Пункт 3 дает оценку расстояния до ближайшей точки минимума, хотя нет гарантии, что последовательность сходится именно к ней.

Так как удалось установить экспоненциальную сходимость для ||хк - х* || и при этом размер шага постоянен, то без особого труда можно доказать, что в случае если на каждой итерации пропускные способности немного изменяются, то ||хк -х* || будет экспоненциально сходиться к некоторому шару с центром в х*. Пусть ск — вектор пропускных

способностей на итерации к, Qk = {х е Кт \ 0 * хе * ск}, /к(х) = 2хТАх + хТшк, шк = Б'ТЬк, Ьк = Б"с"к. В этих обозначениях мы полагаем, что в случае изменения пропускных способностей последовательность оценок генерируется по правилу:

хк+1 = PQk (х - 2 V/к (хк)), (3)

Обозначим также ^ = и^ 1Qk.

Теорема 2. Пусть хк — последовательность, генерируемая по правилу (3) с любым на-чальнымприближениемх0 е Q0, приэтом \\ ск \\ * а, \\ ск-ск-1 \\ * в, Т = 8иРхе^,к \\Ах+шк \\, ц = 1 - А2(А)/2, V = в(1 + £§»), тогда

1. Ближайшая к хк точка минимума х*к функции /к находится на расстоянии не бо-

лее

цк1 \ \х0 - х*0 \\- — + V

1 - ц 1 - ц

2. Близость значения функции /к к оптимальному значению оценивается следующим неравенством

к(хк) - к(х*к) * г [цк ( \ \х0 -х*0 \ \ - ^) + ^

Доказательство. Доказательство практически повторяет теорему 1, если обозначить х*к,и*к,у*к — точка минимума /к и её разложение по П,П1 соответственно, то все выкладки вплоть до

/к(хк) - к(х*к) * г\ \ Ук -у*к \ \

идентичны. Далее, из-за изменения пропускных способностей точка минимума сдвигается не более, чем на в из-за изменения множества Qk, и не более чем на \ \ Б 'ТБ " \ \ в/А2(Б 'ТБ') из-за изменения функции /к. Таким образом,

\ \ ук +1 - у*к+1 \ \ *\ \ (I - 1 А)(ук - У*к ) \ \+ V *

ц \ \ук - ук*\ \ + V•

Вычтем из обоих частей

(1-ц):

\ \ ук+1 -у*к+1 \\-—Ц * ц\\ук -у*к\\+ V-—^ = ц| \ \ук -у*к\\ - —-—-(1 - ц) (1 - ц) (1 - ц)

Таким образом,

\ \ук - у*к \ \* цк1 \ \у0 - у*0 \ \--Ц + ,

1 - ц 1 - ц

' (у4 Иг-ц)) * ц* (их" - х-и-

Так как множество Пк = {х*к + у \ у е П} состоит только из точек минимума /к, то расстояние до ближайшей точки минимума и есть \ \ ук - у*к \ \ .

□

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

4.2. Сходимость рандомизированного алгоритма

Как уже было сказано, рандомизированная версия алгоритма представляет специальный вид рандомизированной стохастической аппроксимации. Обозначим /к(х) = 2хтВ,тВ'х + хтВ,тЬк, Дк — рандомизированное пробное возмущение на итерации к, определяемое по правилу

Дк _ ±4 с. в. Р, (4)

где е1т — т-мерный г-ый орт. Далее, Ак _ ДкД1^В'тВ', А _ БА1, wгk _ е1тв1гТВ'тЪк, wk = Бwгk, у+ _ /(хк + аДк) + (2к и У- _ /(хк - аДк) + (2к+1 — зашумленные измерения на итерации к, при этом (к — погрешность, возникающая при измерении Ьк (который, в свою очередь, зависит от пропускных способностей), (к — случайная величина, распределение которой неизвестно, ((к _ {х е Ит | 0 ^ хе ^ ск}. В этих обозначениях RandomizedArcBalancing генерирует последовательность оценок по правилу

хк+1 = ^к

к ук+ - ук-хк - Дк к к

4а

(5)

Теорема 3. Пусть хк — последовательность, генерируемая по правилу (5) с любым начальным приближением х0 е ((0, Б (к ^ а2> множество и к (к ограничено и у _ 8ирхе^ к II Ах+wk ||, ||ск - ск+1|| ^ в, ^2(А) — минимальное отличное от нуля собственное число матрицы А, д = - ^^], а > 0, V _ в(1 + я2ВВ'ВВ'))' Я _ у + ' тогда

1. Если х*к — любая точка минимума /к на ((к, то

Б8(хк -х*к, П) ^ дк/2(||х0 -х*01| - Я) + Я,

где 6 — евклидова метрика.

2. Если х*к — любая точка минимума /к на ((к, то

Б{/(хк) - /(х*к)} ^ Г(дк/2(||х0 -х*01| - Я) + Я).

3. Если помехи в измерениях отсутствуют, а пропускные способности не изменяются (то есть Я _ 0, ((к _ ((к+1, /к = /к+1) и хк сходится к х*, то

дк/4 ,-

Б||хк -х*|| ^ у у||х0 -х*0||.

1 - у

Доказательство. Во-первых, отметим, что для любого возможного Дк имеет место ||Ак|| ^ 2. Непосредственным вычислением получаем, что при выборе дуги I _ (г,]) выполняется

.I „1тв

(к, ] )еБ' (],Н)еБ' (к,г )еБ' (г ,Н)еБ'

ететВ Вх _ ^ хк ^ хк ^ хк + ^ хк.

Каждая компонента вектора х может входить в любую из этих сумм, а значит любая компонента входит с коэффициентом -2,-1,0,1 или 2, из чего и следует ||Ак|| ^ 2. Из выбора у следует

/(хк) - /(х*) ^ Г||ук - V* ||,

где V и и — разложения х по П1, П. Из леммы 1 следует, что ядро А и В'т В' совпадают. Так как 0 Ц Ак Ц 21, то

Ак т Ак Ак А

Е{(I- —)т(I- —)} ЦЕ{I- —} = I-

Из определения Дк следует, что ЕДк = 0т, ДкДТДк = Дк, а значит, ЕДкДТДк = 0т и, следовательно, ЕДкАк = 0. Таким образом, применяя лемму 3 для точек хк - Дк(у+ - у-)/4а и х*к, получаем

Е{|^к+1 - v*k||2 | хк} ^ Е{||(I- — )(vk -V*) + Дк<2к+1 <2к ||2} =

2 4а

Ак Ак (vk - v*)TE{(I- —)т(I- — )}^к - v*) +

е {(^2к+1 - ^2к )дТ (I - Aг)}(vk - v*k)+4а Е {ДТ Дк (<^+1 - <2к)} ^

А2( А) V,, „*, ,2. V2

а значит,

1--^^ llvk - V*||2 + —,

2 4а

Е{||vk+1 - v*k+1|| | хк} ^ Е{|^ - v*k|| | хк} + V ^ (Е^1 -v*k||2 | хк})1/2 + V ^ [1 - ^ 12 ||vk -v*k|| + V + а

2 ) 2 у/а

Отсюда выводим

По выбору у

Е{|^к -v*k||} - Я ^ qk/2(||х0 - х*0||- Я).

Е/к(хк) - /к(х*к) ^ ТЕ||vk - v*k|| < у^Шх0 -х*0|| - Я) + Я).

Если же а = 0, /к = /к+1, (¿к = Qk+1, то Я = 0, множество точек минимума не изменяется, и расстояние до ближайшей точки минимума в среднем сходится с экспоненциальной скоростью.

Так как ¿к — прямоугольный параллелепипед, а на одном шаге изменения выполняются только по одной координате, то существует такой набор скаляров 0 ^ фк ^ 1, что

хк+1 = хк + фк(хк - 1(Акхк + wk)) -хк) = хк - фк 1(Акхк + wk). (6)

Из выпуклости /к следует

/к (хк+1) - /к (хк) ^

/к (хк) + фк\М хк - 1 (Акхк + wk ) | - /к (хк ) | | - /к (хк) =

Фк /к (хк - 2(Акхк + шк)) - /к(хк)). (7)

Далее, заметим,что /к(х) = ^хТВ,тВ'х+В,тЬк = 2(11В'х+Ьк||2-||Ьк||2). В случае отсутствия помех (5) делает оптимальный шаг в выбранном направлении: выбирается дуга е = (г, ])

и уравниваются значения [В'хк+Ьк]; и [В'хк + Ьк] •, при этом прочие координаты вектора В'хк + Ьк не изменяются. Используя этот факт, можно вычислить Д(хк+1) - Д(хк). Для удобства обозначим ъ _ В'хк Ьк, тогда

Д( хк - 1( Ак хк + wk )1 - Д (хк) _ + ^

22

• Г-( ^) _- 2 с - - >2.

С другой стороны,

а значит,

||1 (Акхк + wk) ||2 _ 11 [В'хк + Ьк] г- - [В'хк + Ь]) |2 _ 1С - Zj• )2, 2 2 2

/к (хк - 1 (Агхк + wk) 1 - /к(хк) _ -|| 2 (Агхк + wk) ||2.

Следовательно, в силу (6) и (7), получаем

|хк+1 - хк ||2 ^ Д (хк) - Д (хк+1)

Если Я _ 0 и хк сходится к х*, то

то то

Б||хк - х* || ^ Б||хг - хг+11| ^ Б(/к (хг) - /к (х*к)1/2 ^

г_к г_к

то /- дк/4 ,-

X у||х0 -х*0|| _ у/Г||х0 -х*0||

г_к

Замечания. А _ ЮВ'тВ', где Б _ diag..., рт}. Для того чтобы выполнялось || А|| ^ 2, годится выбор ре _ 1/(^-), е _ (г,•) или же просто ре _ 1/|Б'|. По сравнению с синхронной версией алгоритма, сходимость немного более чем в два раза медленней, но при этом на одну итерацию гораздо меньше действий.

В пункте 1 оценивается расстояние до ближайшей точки минимума. В случае отсутствия помех и изменений функций /к (пропускных способностей) Я _ 0, имеет место экспоненциальная сходимость к множеству точек минимума /к. В пункте 3 можно вывести сходимость ряда £ ТО 0 Б ||хг -хг+1||, а также сходимость ряда£ ТО 0 ||хг -хг+1||2 для любой возможной последовательности хк, что, тем не менее, ничего не говорит о сходимости хк.

5. ИЗВЛЕЧЕНИЕ РЕШЕНИЯ ЗОДАЧИ О МАКСИМАЛЬНОМ ПОТОКЕ

Из теорем 1 и 3 следует, что оба рассматриваемых алгоритма решают следующую задачу оптимизации:

минимизировать 1 хт В 'тВ 'х + хт В,т Ь, (8)

при условии 0 ^ хе ^ се, е е Б'.

При этом на исходном графе О мы получаем псевдопоток, удовлетворяющий пропускным способностям, в котором насыщены все дуги, инцидентные истоку и стоку. Обозначим Щх) _ {г е V \{х, г} | [В'тх + В'тЬ]г _ 0}, У-(х) _ {г е V \{х, г} | [В'тх + В'тЬ]г < 0}, V+ (х) _ {г е V \ {х, г} | [В'тх + В'тЬ] г > 0}. Ключевым моментом в сведении задачи о максимальном потоке к задаче (8) является следующая теорема:

Теорема 4. Если х* — решение задачи (8), то

1. Множества А = {х} и У+(х*) и Щх*) и V \ А образуют минимальный вЬ-разрез графа в.

2. Для получения максимального потока достаточно провести частичную декомпо-

*

зицию потока х .

Доказательство. Во-первых, применяя условия Каруша-Куна-Такера к задаче (8), получаем

[В'х* + Ь]] - [В'х* + Ь] г + Ле - Хе = 0,

¡Л, Х & 0т,

0 = Ле(се - хе) = хехе>

где е = (г,]), це — множитель Лагранжа, соответствующий ограничению хе ^ се, а ве — множитель Лагранжа, соответствующий 0 ^ хе. В данном случае как ве, так и це играют вспомогательную роль, и от них можно избавиться следующим образом: пусть Ех = {е = (г,]) е Е' | се - хе > 0} и {(], г) | е = (г,]) е Е', хе > 0}, тогда условия ККТ для (8) можно переписать в виде:

[В'х*+ Ь]г ^ [В'х*+ Ь]], (г,]) е Ех*. (9)

Проще говоря, если избыток вершины г (величина [В'х+Ь] г) больше избытка вершины ], то (г, ]) £ Ех*. Теперь вернемся к исходному графу в. Избыток потока для вершин в графе

G задается вектором B

, по построению для вершин, отличных от s, t, га совпадает

с B'x+b, для истока s — это величина (-£(s,;)eEce) ^ 0, для стока — величина £(г-,t)eEce ^ 0. По теореме о декомпозиции потока существует конечный набор путей, каждый из которых ведет из вершины с отрицательным избытком в вершину с положительным избытком и проходит только по дугам с положительным потоком. В силу (9), для вершин с положительным избытком кроме t такой путь полностью проходит только по вершинам из A, а значит, начинается в s, для вершин с отрицательным избытком кроме s такой путь полностью проходит только по вершинам из V \ A, а значит, заканчивается в t. Пусть, y — поток, получившийся в результате удаления потока по соответствующим путям из s в V+(x*) и путям из V-(x*) в t. По построению y является допустимой точкой задачи о максимальном потоке (1), при этом не существует дуг (i, j) e Ey, i e A, j e V \ A, а значит, из теоремы Форда-Фалкерсона y - максимальный поток, а A, V \ A - минимальный st-разрез. □

Замечание. В этой теореме также можно было взять A = {s} и V+(x*). В работе [2] также отмечено, что, после того как в общей процедуре балансирования дуг ни одна из операций Move не дает изменения потока больше, чем на 1/n2, возможно точное извлечение максимального потока / минимального разреза, что позволяет использовать приближенное решение (8), а не точное.

x*

6. РАСПРЕДЕЛЕННАЯ РЕМИЗАЦИЯ М1ГОРИТМОВ

Оба алгоритма СопстгеШАгсБаЫпащ и RandomizedArcBalancing допускают распределенную реализацию в сети, если топология этой сети задана тем же графом, что и задача о максимальном потоке. Также можно считать, что вычислительные узлы соответствуют не всем вершинам, а только V \ {х, Ь}. Для ускорения вычислений имеет смысл непосредственным образом вычислять избыток вершины. Для обоих алгоритмов мы будем предполагать следующее:

• Алгоритмы оперируют с величинами xe, e е E', yi, i е V \ {s, t}.

• Узел i хранит в локальной памяти величины xe и ce, где e = (i, j) е E'.

• В синхронной версии используются дополнительные промежуточные переменные yi,...,yn, в рандомизировнной версии каждый узел хранит дополнительно

четыре переменных для промежуточных вычислений z+, z +, z-, z-.

i j i j

Далее представлены распределенные версии алгоритмов ConcurrentArcBalancing и RandomizedArcBalancing.

Фактически, они ничем не отличаются от представленных ранее, но имеют более явную форму. В алгоритме ConcurrentArcBalancing предполагаем, что в сети есть синхронизация времени, все операции двух внутренних циклов for происходят одновременно по тактам системы. В этом случае, на каждой итерации достаточно хранить только значения xk,xk+1, yk, yk+1. В обоих алгоритмах предполагается, что узел i должен иметь возможность запросить информацию соседних узлов, то есть таких j, что (i, j) е E, и при этом передать информацию на этот узел. Другими словами, условие двустороннего взаимодействия необходимо для реализации этого алгоритма (однако это еще не значит, что для соответствующей задачи о максимальном потоке пропускные способности должны быть двусторонними).

Algorithm 3: Distributed ConcurrentArcBalancingG, c, stopping criterion)

x0 - 0m, к — 0;

for e = (s, i) е E do

L

x° <- c,

e;

for e = (i, t) е E do

x° <- c,

e;

for i e V \ {s, t} do

L y°i ^^(j,i)eEx°e -E(i,j

/* End of initialization */

while stopping criterion is not satisfied do for e = (i, j) e E' do

xk+1 xk _ yj -yt . x~ * x „

2

if xk+1 > ce then I xk+1 c

L xe — be. if xe < 0 then

l_ xk+1 - 0.

for i е V \ {s, t} do

I yk +1 v xk v- vk.

L yi — (j,i)еExe (i,j)еExe;

к — к + 1;

return x;

В рандомизированном алгоритме предполагается, что в системе нет синхронизации, каждый из узлов совершает Move относительно случайной исходящей дуги. Таким образом, можно считать, что сначала случайным образом выбирается узел, который совершает действие, после чего случайно выбирается исходящая дуга, что дает вероятность выбора дуги e = (i, j): pe = —!р. Из леммы 2 следует, что такой выбор вероятности га-

рантирует выполнение неравенства ||А|| ^ 2 в соответствующих теоремах сходимости. Результатом обоих алгоритмов является решение (8), извлечение решения задачи о максимальном потоке было обсуждено в предыдущем разделе, оно заключается в частичной декомпозиции потока и представляет собой гораздо менее сложную задачу, нежели сама задача о максимальном потоке. Извлечение минимального разреза возможно в распределенном виде: если > 0, то узел г находится в одной части разреза, если же у I ^ 0, то в другой.

Algorithm 4: Distributed RandomizedArcBalancing(G, c, stopping criterion) x — 0m; while stopping criterion is not satisfied do pick random node i uniformly;

pick random outcoming arc e = (i, j) e E' of i uniformly; let A — ±1 w. p. 1/2;

z+ — (T.(k,i)eE' xe - X(i,k)eE' xe + s,i)eEce - X(i,t)eE' xe + aA); z+ — (E(k,j)eE' xe -Y.(j,k)eE' xe + T.(s,j)eEce -Y.(j,t)eE' xe + aA); Z- — (T.(k,i)eE' xe - T.(i,k)eE' xe + T.(s,i)eEce - T.(i,t)eE' xe - aA); z- — (£(k,j)eE' xe (j,k)eE' xe + Z(s,j)eEce -Y.(j,t)eE' xe - aA);

(z+ )2+(z+)2-(z-)2-(z-)2 xe — xe - A 80 ; if xe > ce then

L xe — ce •

if xe < 0 then

L xe — 0.

return x;

7. ЗАОЮЧЕНИЕ

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

Список литературы

1. Goldberg A., TarjanR. Efficient maximum flow algorithms // Communications of the ACM, 2014. Vol. 57, №. 8. P. 82-89.

2. Tarjan R., Ward J., Zhang B., Zhou Y., Mao J. Balancing applied to maximum network flow problems // European Symposium on Algorithms, 2006. Springer, P. 612-623.

3. Левитин E.C., ПолякБ.Т. Методы минимизации при наличии ограничений // Журнал вычислительной математики и математической физики, 1966. Т. 6, №. 5. С. 787-823.

4. Граничин О.Н. Поисковые алгоритмы стохастической аппроксимации с рандомизацией на входе // Автоматика и телемеханика, 2015. №. 5. С. 43-59.

5. Granichin O., Volkovich V., Toledano-Kitai D. Randomized algorithms in automatic control and data mining. Springer, 2014.

6. Boyd S., Ghosh A., Prabhakar B., Shah D. Randomized gossip algorithms //Information Theory, IEEE Transactions on, 2006. Vol. 52, №. 6. P. 2508-2530.

7. Olfati-Saber R., Murray R.M. Consensus problems in networks of agents with switching topology and time-delays // Automatic Control, IEEE Transactions on, 2004. Vol. 49, №. 9. P. 1520-1533.

Поступила в редакцию 10.08.2016, окончательный вариант — 11.09.2016.

Randomized Distributed Adaptive Algorithm for a Maximum Flow Problem

Computer tools in education, 2016 № 5:46-62

http://ipo.spb.ru/journal

RANDOMIZED DISTRIBUTED ADAPTIVE ALGORITHM FOR A MAXIMUM FLOW PROBLEM

Мальковский Н.В.1 1SPbSU, Saint-Petersburg, Russia

Abstract

This paper presents an algorithm for solving one of the fundamental optimization problems — maximum flow problem. The algorithm is based on the ideas of arc balancing procedure and projected gradient method. Although the algorithm does not improve asymptotic complexity over existing methods it still possess some usefull properties: natural implementation in distributed systems and convergence to an optimal solution even if network parameters are changing slightly over time. Two versions of the algorithm are considered: randomized and concurrent. Both versions mainain adaptability but randomized version is preferable due to better convergence rate and simpler implementation in distributed systems.

Keywords: constraint optimization, gradient descent, distributed algorithms, randomized algorithms, network flow.

Citation: Malkovskiy, N., 2016. "Randomizirovannyi raspredelennyi adaptivnyi algoritm resheniya zadachi o maksimal'nom potoke" ["Randomized Distributed Adaptive Algorithm for a Maximum Flow Problem"]. Computer tools in education, no. 5, pp. 46-61.

Received 10.08.2016, the Anal version — 11.09.2016.

Мальковский Николай Владимирович, магистрант математико-механического факультета СПбГУ; 198504 Санкт-Петербург, Старый Петергоф, Университетский пр. 28, математико-механический факультет, [email protected].

© Наши авторы, 2016. Our authors, 2016.

Мальковский Николай Владимирович, магистрант математико-механического факультета СПбГУ; 198504 Санкт-Петербург, Старый Петергоф, Университетский пр. 28, математико-механический факультет, [email protected].