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

УДК 51.77

А. В. Гасников1'2'3, Е.В. Гасникова3, П.Е. Двуреченский2'3'4, Е. И. Ершов2,

A.A. Лагуновская3'5

1 Центр исследований транспортной политики, Институт экономики транспорта и транспортной

политики НИУ ВШЭ 2Институт проблем передачи информации им. А.А. Харкевича РАН 3Лаборатория структурных методов анализа данных в предсказательном моделировании при

МФТИ(ГУ) (ПреМоЛаб) 4Weierstrass Institute for Applied Analysis and Stochastics, Berlin 5Институт прикладной математики им. М.В. Келдыша РАН

Поиск стохастических равновесий в транспортных моделях равновесного распределения потоков

В работе предложены эффективные способы поиска стохастических равновесий в популяционных играх загрузок. Поиск равновесия Нэша в таких играх всегда сводится к задаче оптимизации. Мы рассматриваем модели равновесного распределения потоков по путям Бэкмана и Нестерова-де Пальмы. Поиск стохастических равновесий Нэша(-Вардропа) приводит к энтропийной регуляризации выпуклого функционала, отвечающего этим моделям. Данная работа посвящена тому, как эффективно решать такого рода задачи. В основе подхода лежит идея композитной оптимизации и особенность постановки, что функционал имеет вид суммы (сепарабельный функционал). Это обстоятельство вместе с неограниченностью константы Липшица градиента функционала мотивирует переформулировку исходной задачи оптимизации таким образом, чтобы этот сепарабельный функционал стал композитным членом. Рассматриваются и развиваются также и классические способы решения отмеченной задачи с помощью аппарата характеристических функций на графе.

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

A. V. Gasnikov1'2'3, E. V. Gasnikova3, P. E. Dvurechensky2'3'4, E. I. Ershov2,

A. A. Lagunovskaya3'5

Research Centre for Transport Policy Studies, Institute for Transport Economics and Transport Policy Studies, National Research University Higher School of Economics

2Institute for Information Transmission Problems of the Russian Academy of Sciences

(Kharkevich Institute)

3Research laboratory in Predictive Modeling and Optimization at PhysTech (PreMoLab) 4Weierstrass Institute for Applied Analysis and Stochastics, Berlin 5Keldysh Institute of Applied Mathematics

Search for the stochastic equilibria in the transport models of equilibrium flow distribution

In this paper, we propose efficient ways to search stochastic equilibria in congestion population games. A Nash equilibrium search in these games can always be reduced to an optimization problem. We consider the model of equilibrium flow distribution over paths, which is proposed by Beckman and Nesterov and de Palma. The stochastic Nash equilibria (-Vardrop) search leads to the entropic regularization of a convex functional corresponding to these models. This paper is devoted to effective approaches to such problems. The approach is based on the idea of composite optimization and a well- posed problem that the functional has the form of sum (i.e. a separable functional). This fact together with the unlimited Lipschitz constant of the functional gradient motivates reformulation of the original optimization problem so that the separable functional becomes a composite term. The classical methods for solving the problem in question using the characteristic functions on the graph are also considered and developed.

Key words: composite fast gradient method, sparsity, equilibrium flow distribution.

1. Введение

В работе [1] было анонсировано, что в цикле статей планируется описать эффективные численные методы поиска равновесий (и стохастических равновесий) в популярных моделях распределения транспортных потоков по путям. В частности, это планировалось сделать для модели Бэкмана и модели стабильной динамики (Нестерова-де Пальмы). В работе [2] мы сконцентрировались на поиске обычных равновесий (Нэша). Данная работа посвящена эффективным численным методам поиска стохастических равновесий в этих моделях. Стохастические равновесия играют важную роль в практических приложениях этих моделей [3, 4]. В популяционной теории игр исследуемые нами стохастические равновесия иногда называют логит-равновесиями [5]. Такие равновесия получаются, если вместо динамики наилучших ответов (приводящей к равновесию Нэша при условиях потенциальности игры с выпуклым потенциалом [5]) и положительно коррелирующих с ней динамик (также приводящих к равновесию Нэша) рассматривается логит-динамика (гиббсовская динамика), отражающая поведение ограниченно рациональных агентов, т.е. агентов, действующих в условиях не полной информации (немного подробнее об этом будет написано в п. 2, посвященном постановке задачи и в п. 5). Последнее обстоятельство позволяет связать такие равновесия с достаточно популярным разделом современной экономической теории - с теорией дискретного выбора [6].

В пп. 2-4 приводятся новые подходы к поиску стохастических равновесий, в основе которых лежат современные результаты из композитной оптимизации [7,8]. Ожидается, что эти методы будут неплохо работать в случае, когда число возможных путей в графе транспортной сети не очень большое по сравнению с числом ребер. В общем же случае число возможных путей может быть «экспоненциально» велико по сравнению с числом ребер. Тогда методы из пп. 2-4 доминируются методом, базирующимся на аппарате характеристических функций на графе [4, 9]. Дополнительно к описанию подхода работ [4, 9] мы приводим в п. 5 способ ускорения описанного там метода для графов специальной (упорядоченной) структуры, а также приводим новый рандомизированный метод, тесно связанный с алгоритмом из п. 3 [2].

В данной работе мы сосредоточились на двух моделях равновесного распределения транспортных потоков по путям: модели Бэкмана и модели стабильной динамики. Однако легко перенести все, что далее будет написано, например, на смешанные модели [2]. Также заложенные в этой статье и в [2] подходы впоследствии должны позволить разработать эффективные численные методы для поиска равновесий в многостадийных транспортных моделях, одним из блоков которых являются изучаемые в данной статье модели распределения транспортных потоков по путям [1,10-12].

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

Рассмотрим транспортную сеть, которую будем представлять ориентированным графом {V, Е), где V - множество вершин (как правило, можно считать, что | /4 ^ | V| ^ |^|), а Е

- множество ребер (|^| = т). Обозначим множество пар и> = (I,]) источник-сток через ОИ,

(1Ш - корреспонденция, отвечающая паре ■ш, хр - поток по пути р; Р,ш - множество путей,

отвечающих корреспонденции и> (начинающихся в г и заканчивающихся в ]), Р = У Рт

-шеоо

- множество всех путей (|Р| = п). Затраты на прохождение ребра е € Е описываются функцией те (/е), где /е - поток по ребру е.

Опишем марковскую логит-динамику (также говорят гиббсовскую динамику) в повторяющейся игре загрузки графа транспортной сети [5]. Пусть каждой корреспонденции отвечает М агентов (М » 1),те (/е) := те (/е/М).

Пусть имеется ТЖ шагов (Ж ^ 1). Каждый агент независимо от остальных на шаге £ + 1 выбирает с вероятностью

Л ехр

" Д ехр ИЛ)

путь р е Р,ш, где Ср - затраты на пути р на шаге £ (С® = 0), ас вероятность 1 — ЛЖ путь, который использовал на шаге ¿. Такая динамика отражает ограниченную рациональность агентов (см. п. 5), и часто используется в популяционной теории игр [5] и теории дискретного выбора [6]. Оказывается, эта марковская динамика в пределе N ^ ж превращается в марковскую динамику в непрерывном времени (вырождающуюся при 7 ^ 0+ в динамику наилучших ответов [5]), которая в свою очередь допускает два предельных перехода (обоснование перестановочности этих пределов см. в [13]): Т ^ ж, М ^ ж или М ^ ж, Т ^ ж. При первом порядке переходов мы сначала (Т ^ ж) согласно эрго-дической теореме для марковских процессов (в нашем случае марковский процесс - модель стохастической химической кинетики с унарными реакциями в условиях детального баланса [14-16]) приходим к финальной (стационарной) вероятностной мере, имеющей в основе мультиномиальное распределение. С ростом числа агентов (М ^ ж) эта мера концентрируется около наиболее вероятного состояния, поиск которого сводится к решению задачи (1) ниже. Функционал в этой задаче оптимизации с точностью до потенцирования и мультипликативных и аддитивных констант соответствует исследуемой стационарной мере, то есть это функционал Санова [16,17]. При обратном порядке переходов мы сначала осуществляем так называемый канонический скейлинг [13,14], приводящий к детерминированной кинетической динамике, описываемой СОДУ на х, а затем (Т ^ ж) ищем аттрактор получившейся СОДУ. Глобальным аттрактором оказывается неподвижная точка, которая определяется решением задачи (1) ниже. Более того, функционал, стоящий в (1), является функцией Ляпунова полученной кинетической динамики (то есть функционалом Больцмана). Последнее утверждение - достаточно общий факт (функционал Санова, является функционалом Больцмана), верный при намного более общих условиях [16]:

У^е (/е)+7 У2 У2 ХР 1п (хР/(1™) ™П ^ (1)

евЕ т€ОПр€Рт

и

где 7 > 0,ое (/е) = / те (г) йг - выпуклые функции; о

© = II^ер\\е£Е,р£.Р = Ц©^ = > ^ = ^ ^

важно заметить, что матрица © сильно разрежена по столбцам (в - среднее число ребер в пути, т.е. среднее число ненулевых элементов в столбце матрицы ©, как правило, можно считать, что в ^ л/т );

X = <х ^ 0: = &ш ^ е ОБ

{ реР„

- прямое произведение симплексов; в качестве те (/е) обычно выбирают БРИ-функции [1,4,18,19]:

Те (/е)=*е • (1+р • Це/¡е)^ . В пределе модели стабильной динамики [1]

(/е)=й • (1+Р • ,

Нт т? (/е) = { ^

^0+ [ [t е,

задача перепишется как [1]

0 < fe < fe

e ( ' e)^[ie, «>) , fe = U

У2 fete + 7 У У xp In (xp/dw) ^ min (2)

eeE weODpePw f /j

Данную задачу можно переписать как задачу энтропийно-линейного программирования и использовать соответствующие численные методы [20,21]. Это приводит к следующей оценке числа арифметических операций: О [вПл/Ъ^Щ/е 1, где L2 = max

Цв^Ц;; = Н (6<р>) -

р-й столбец матрицы 6 , Н - максимальное число ребер в пути, как правило, можно считать, что Н = О (y/m), R2 - евклидов размер решения двойственной задачи (итерационный процесс в двойственном пространстве стартует с нуля), е - точность по функции (в прямой задаче). При этом ограничения выполняются с точностью е/R2. Интересно было бы попробовать использовать для решения двойственной задачи к (2) не быстрый градиентный метод FGM, как в [21], а его покомпонентный вариант ACRCD, также обладающий прямо-двойственной структурой [22]. Тем не менее на данный момент нам не удалось теоретически обосновать, что это может принести какие-то дивиденды по сравнению с использованием обычного FGM (даже при дополнительных предположениях типа разреженности матрицы 6).

3. Композитный подход для задачи (1)

Прежде всего, заметим, что f (x) = 6x и функции oe (fe) гладкие (с ограниченной константой липшица градиента) на компакте {/ = 6x : x £ X}. Следовательно, к задаче (1), записанной в следующей форме:

^ Oe (6^x) + 7 ^ Y,xp In (xp/dw) ^ min (3)

eG.E weODpePw

можно применять композитный быстрый градиентный метод [7,8], считая композитом энтропию 7 ^ ^ хр \п(хр/йт), выбирая норму теоореРы

llxll = \\{x}vepw 111

|x|| =./ > \lix}pGPw ||1

weOD

в прямом пространстве и 1-сильно выпуклую в этой норме прокс-функцию

(х) = ^ йт •(Т ({хр}р€р^ , теоо

({Хр}реРш) = (т \п |РШ | + ^ Хр \п(хр/йт),

реРы

если 7 > 0 мало (не учитываем сильную выпуклость функционала), иначе (учитываем сильную выпуклость функционала) выбираем 1-сильно выпуклую прокс-функцию

(2 (х) = ^ • —{хр}реры

теоо а'ш 1

ttw -

2 In |Pw| 21n|Pw| - 1'

В случае, когда 7 > 0 мало, имеем (согласно [7]) следующую оценку трудоемкость метода в терминах общего числа арифметических операций (здесь использовано, что О (вN) -

с

стоимость одной итерации, поскольку шаг композитного метода осуществляется по явным формулам, и основная сложность заключается в вычислении / = ©ж):

Т)

где

' (4)

L =

f=&x, х£

max max (h, 0 diag \т e (fe)\ 0h) < max т„ (fe) ■ max 0^

&x,xex im, <1 \ 6 1 eyjeJ>eeE / ^ f=ex,xex eKJe' PeP

евЕ

R2 = max^ (x) = V^ d2w ln IPW |.

x x

xex

weOD

В случае, когда ^ > 0 немало, имеем (согласно [7]) следующую оценку трудоемкости метода (здесь использовано, что О (вп) - стоимость одной итерации [21]):

(' ^ ^1П ("

0\snJLX ■ 1nl ЦХ l|x° -x^)) , (5)

где

1

ц = 7 ■ I max dw ) \weOD J

- константа сильной выпуклости композита в введенной норме на X, а множитель X < 2 max ln |PW | возникает из-за использования техники рестартов в прямом простран-

wGOD

стве и прокс-функции d2 (x) (см., например, с. 21-22 [23]1 и п. 5 работы [24]). Из формул (3), (4) можно сделать вывод, что учитывать 1-сильную выпуклость композита на X стоит, если R2/е » х/ц, т.е. при

\ -1

7 » е ■ I у ] dw

(Edw) .

\weOD J

В частном случае |О = 1 описанный здесь подход ранее предлагался в работе [21]. Можно использовать и другие способы решения задачи (3) из [21], например, базирующиеся на рандомизированном двойственном покомпонентном методе.

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

Помочь в решении отмеченной проблемы может следующее наблюдение. В действительности, функции е ( е) имеют штрафную природу: при приближении е к пропускной способности /е наблюдается резкий рост. На самом решении (в равновесии) и в его окрестностях не стоит ожидать, что пропускные способности сильно превышены. Это позволяет надеяться, что по ходу итерационного процесса это свойство также будет выполнено. Для этого нужно еще стартовать с разумной точки, в которой нет сильных нарушений пропускных способностей. По условию мы считаем, что это возможно, иначе искомое равновесие будет «одной большой пробкой». Если предположить, что есть метод, который настраивается на реальную константу Липшица градиента, а не на ее оценку сверху Ь, то во многом отмеченная проблема будет решена. Оказывается, можно предложить универсальный вариант [8,25] описанного выше композитного метода (в том числе с учетом сильной

1 Заметим, что в сильно выпуклом случае мы используем сильную выпуклость композита (этим объясняется некоторое противоречие с тем, что написано на с. 21-22 [23] - там речь идет об одной и той же функции, когда говорится, что не имеет смысл рассматривать в сильно выпуклом случае прокс-структуры, отличные от евклидовой: это верно для не композитных постановок).

2

2

выпуклости и возникающей неточности при вычислении шага градиентного отображения в сильно выпуклом случае [26]), который оптимально самонастраивается на текущую гладкость оптимизируемой функции (за вычетом композита). Другой способ борьбы с большим Ь будет описан далее.

4. Штрафное раздутие задач

Задачи (1) и (2) можно переписать следующим образом (см. [21] в связи с тем, как следует выбирать Л):

1 ||вх - /У22 + Ai ■ I (Д) + 7 ln (xp/dw) ) ^ XXX (6)

\ееЕ weODpeP, х

1 ||вх-/||2 + A2 ■ (У fete +7 У У хр ln (xp/d-) ) ^ min _. (7)

2 \е€Е weODpePw J xeX'f

Мы выбираем норму в прямом пространстве (х, /) следующим образом:

II (х, / )|| = ||х|| + ||/112.

Прокс-функция выбирается: di(x) + 21|/1|2. Далее мы имеем композитную постановку с композитом вида

Yhae (fe )+7 Y1 YhXP ln(xP/dw )

ееЕ -шеоореРы

С точностью до симплексных ограничений на х это сепарабельный композит. Далее к нему следует применить (должным образом модернезировав) методы из работы [21].

Приведем, например, оценку трудоемкости композитного варианта быстрого градиентного метода для задач (6), (7):

O

I /33Г .LR

V J

L = max (h, втeh) + max (f, f) = max 7 \\f \\2<i p^p

e{p)

+ 1,

2

2 1 -R = dWw ln iPw| + -II/

2

теоо

обеспечивающую -решение (по функции) задачам (1), (2) и выполнение неравенства ||©х — /Ц2 ^ е. Отметим, что каждый шаг композитного метода считается по явным формулам, поскольку переменные распадаются на две группы: по х - это стандартное экспоненциальное взвешивание, а по потребуется решать уравнение четвертой степени для БРИ-функций (см. п. 2), что также можно сделать по явным формулам.

2

5. Использование аппарата характеристических функций на графе

Запишем двойственную задачу к (1) [1,2,4] (далее мы используем обозначение ёоша* - область определения сопряженной к а функции):

mm < (fe) + 7 y y xp 1п(жр/dw): f = ©x, x gx

f'x \eee weodpepw

mm < У max [ fete - a*e (ie)] + 7 У У xp In (xp/dw) : f = ©x, x G X

f,x\ — ie€dom a* L—' L—'

edE

weODpePw

max < mm

iGdom a* I f,x

fete + 7 ^ Xp In (xp/dw) : f = ©x, x G X

E

weODpePw

min j71 (i/7) + (ie)l

iedom a* z—'

I eEE )

- £ ^e

E

где

1 (t) = £ dwiw (t) , iw (t) = In I £ exM - £ 5epte I , weod \pePw \ eeE J J

^-1 E Me)

f = -V1 (t/i), Xp = dw-, eeE-Y, p G pw,

E exp ( -1 E Septe)

ре P,

e E

для

Te (fe) =te • + P ^ J имеем [1] (BPR-функция, см. п. 2, получается при v = |)

© 0,

(Г*е (ie) = SUp I (te - te) • fe - ^

U>0\ 1+P

V fe

P • —

1+ 7

= fe

fe"

(te - te \

Vie •P/

t e te\ (t e te)

P J 1 + V '

Собственно, формула (5) есть не что иное, как отражение формулы f = -Vip(t/7) и связи te = Te( fe),e G ^. Действительно, по формуле Демьянова-Данскина-Рубинова [1]

Va*e(te) = Vie max I te fe - f Te (z) dz^ = fe : te = Te( fe).

В свою очередь, формула f = -У^ (£/7) может интерпретироваться как следствие соотношений f = ©ж и формулы распределения Гиббса (логит-распределения):

ехР ( -1 Е Ме )

= , _V е / р

Жр --/ \ , Р ^ ^ ад.

Е ехЫ -1 £ 5еptе ) ре V еее /

Обратим внимание, что вектор распределения потоков по путям ж при поиске стохастического равновесия (7 > 0) получается не разреженным в отличие от поиска обычного равновесия Нэша(-Вардропа) [2] (7 = 0). Как следствие, чтобы вычислить этот вектор требуются затраты, существенно зависящие от потенциально огромной размерности те. К счастью, в приложениях, как правило, не требуется знание этого вектора, достаточно определить вектор потоков на ребрах /, который, как мы увидим ниже, может быть вычислен намного эффективнее (в частности, с затратами независящими от те).

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

ребра которого взвешены вектором ¿, равны др (£) = ^ 5ерЬе, последняя формула для хр

ееЕ

есть не что иное, как отражение следующего принципа поведения (ограниченной рациональности агентов [6]): каждый агент к (пользователь транспортной сети), отвечающий корреспонденции — € ОИ, выбирает маршрут следования р € Р,ш, если

Р = arg max { -gq (i) + ,

где независимые случайные величины ^, имеют одинаковое двойное экспоненциальное распределение, также называемое распределением Гумбеля3 [5,6]:

Р < с) =ехр {-е-с/"'-Е } ,7> 0.

Отметим также, что если взять Е & 0.5772 - константа Эйлера, то М £к = 0, И £к = 72к2/6. Распределение Гиббса получается в пределе, когда число агентов на каждой корреспонденции стремится к бесконечности (случайность исчезает и описание переходит на средние величины). Заметим, что при 7 ^ 0+ распределение водителей по путям вырождается, и все водители (агенты) будут использовать только кратчайшие пути:

- Ит 7фт {Ъ/7) = тт др ф.

реР-ш

Полезно также иметь в виду, что [7]

7^ (t/7) = МШреРт

В пределе стабильной динамики задача (5) вырождается в задачу

7Ф (t/7) + {f,t - t) ^ min . (9)

max {-9p № + (p}

pZPw

t>t

В пределе 7 ^ 0+ все приведенные выше формулы переходят в соответствующие формулы для модели Бэкмана и модели стабильной динамики [2].

Задачи (5), (8) можно не различать по сложности при композитном подходе к ним [7]. В задаче (8) сепарабельный композит проще сепарабельного композита задачи (5), но зато дополнительно добавляется сепарабельное ограничение простой структуры (для BPR-функций все также сводится к решению уравнения четвертой степени, см. п. 2). В любом случае основные затраты на каждой итерации при использовании быстрого градиентного метода в композитном варианте [7] связаны с необходимостью расчета градиента гладкой функции 7Ф (t/7), при этом константа Липшица градиента этой функции в 2-норме может быть оценена сверху числом [20,21]

1 2

L2 = - У dw max 0<р> .

7 ТТ^г, peP™ 2

weOD

Заметим, что возникающие при расчете V^ (t/7) отношения экспонент стоит сразу же приводить (для большей вычислительной устойчивости) к дробям с числителем, равным 1

ea 1

еа + еЬ + ... 1 + еЬ-а + .

Распределение Гумбеля можно объяснить исходя из идемпотентного аналога центральной предельной теоремы (вместо суммы случайных величин - максимум) для независимых случайных величин с экспоненциальным и более быстро убывающим правым хвостом [27]. Распределение Гумбеля возникает в данном контексте, например, если при принятии решения водитель собирает информацию с большого числа разных (независимых) зашумленных источников, ориентируясь на худшие прогнозы по каждому из путей.

Используя технику быстрого автоматического дифференцирования [28], можно показать, что сложность вычисления градиента функции ф(t/7) не более чем в четыре раза дороже сложности вычисления значения функции, однако при этом возрастают затраты памяти, поскольку приходится хранить все «дерево вычислений» ф (¿/7).

Приведем, следуя [4,9], сглаженный идемпотентный аналог метода Форда-Беллмана [29,30], позволяющий эффективно рассчитывать значение характеристической функции ф (£/7). Для этого предположим, что любые движения по ребрам графа с учетом их ориентации являются допустимыми, т.е. множество путей, соединяющих заданные две вершины (источник и сток), - это множество всевозможных способов добраться из источника в сток по ребрам имеющегося графа с учетом их ориентации. Этого всегда можно добиться раздутием исходного графа в несколько раз за счет введения дополнительных вершин и ребер. Такое раздутие заведомо можно сделать за О (т).

Будем считать, что число ребер в любом пути не больше Н = О (л/т) (диаметр графа с манхетенской структурой, то есть квадратной решетки). Введем классы путей: Р^ -множество всех путей из г в j, состоящих ровно из I ребер, Р^- - множество всех путей из I в ], состоящих из не более чем I ребер. Зафиксируем источник (вершину) г £ V и введем следующие функции для ] £ V,1 = 1,..., Н:

alj (í) = 1ФРi. (í/7) = 7In I Е ехр -Е 5epte 7

ipep¿

Е Septe!7)

Ь1г1 (t) = 1ФР1] (t/7) =7 ln | Е ехр Е Me/7) j .

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

(t) = Н, (t) =

- te, е = (г --ж, е = (г

Э) /Е,

>j) /Е,

4+1 (í)=7ln

Е ехр ((<4 (t) - te) /7)

,k: e=(k^-j)eE

Í+1 (t) =7 ln (ехр (ъ1г1 (í)/7) + ехр (ai+1 (t)/7) )

j/V, l = l,...,H - 1.

На каждом шаге l необходимо сделать О (т) арифметических операций. Следовательно, для вычисления ф (t/j) необходимо сделать О (SHm) арифметических операций, где S = |О| - число источников, как правило, можно считать S ^ т. Причем вычисление функции ф (t/7) (и ее градиента) может быть распараллелено на S процессорах. При 7 — 0+ процедура вырождается в известный метод Форда-Беллмана (динамическое программирование). В процедуре Форда-Беллмана требуется посчитать H-степень матрицы А = \\aij\\i,jev,

aij = te,e = (i — j) / Е; aij = ж,е = (i — j) / Е,

в идемпотентной математике (вместо обычного поля используется тропическое полуполе [31] со следующими операциями: сложение a фЬ = min{a, Ъ} , произведение a <&Ъ = a + b). Учитывая, что

a фЬ = min{a, b} = - lim 7 ln (ехр (-a/7 ) + ехр(-6/7)),

7^0+

a = a + b = - lim 7 Ы(ехр(-a/7) ■ ехр(-Ь/7))

7^0+

можно посчитать обычную (над обычным полем) H-степень матрицы А1 = ||a11|íjev,

al = е-= (i ^ j) £ Е;

al = 0, е = (i ^ i) £ Е,

и применить поэлементно к полученной матрице —7 ln(-). В пределе 7 ^ 0+ получим метод Форда-Беллмана. Однако если не делать предельный переход, то получается нужный нам сглаженный вариант этого алгоритма с такой же временной сложностью. Некоторый аналог этого сглаженного варианта, по сути, и был описан выше.

Оказывается, при определенных условиях (возможности упорядочивания) оценку О (SHm) можно редуцировать до О (Sm). Опишем эти условия. Для простоты обозначений будем далее считать, что S = |D| ^ m - число стоков. Переход от источников к стокам также является подготовительным шагом для последующего описания рандомизированного подхода. Предположим, что при любом фиксированном стоке j £ D можно так пронумеровать вершины, что сам сток имеет наибольший номер v = |V| (далее j = v), и для любых двух вершин k, г £ V из того, что к > г следует, что не существует ребра из к в . Выполнение этого свойства можно попытаться добиться раздутием исходного графа в несколько раз за счет введения дополнительных вершин и ребер. Такое раздутие и нумерацию можно осуществить за О (m). Далее, следуя п. 5.4 [32], опишем процедуру расчета Фiv (t/7) ,i £ V:

фт (t/7) = 0, ф™ (t/7) =7ln I £ exP ((tpkv (t/7) — te)/7)

\k: e=(í^k)eE

Отсюда видно, что для вычисления ф^ (t/7) ,i £ V, необходимо сделать О(т) арифметических операций, а для вычисления ф (t/7) соответственно О (Sm).

Таким образом, используя композитный быстрый градиентный метод [7,8], можно рассчитывать на следующее число арифметических операций (в действительности, лучше использовать композитный вариант универсального метода [8, 25, 26], поскольку явное использование в методе полученной выше оценки константы Липшица градиента L2 приводит к «лишним итерациям», ввиду того, что полученная оценка L2 - завышена; при универсальном подходе метод сам настраивается на подходящую гладкость - адаптивно и не требует никаких констант на входе):

О ISm\'L2Щ

L

Или в общем случае

О I SHm ',L2ie>

Здесь Щ - евклидов размер решения задачи (5) (или (8) в зависимости от контекста).

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

7Ф^к/7) + (¿k)

E

к

. 1

— min < —

íGdom a* Ак

-Ы {tlh) + (f/7),t — f

a .¿=0

+ £ ^e (te)

<

С L2R2

I Ак

E

где константа С < 10, ак ~ к, Ак = Е¿=оа», Ак ~ к2, Щ - евклидов размер решения задачи (5). Для метода зеркального спуска правая часть этого неравенства записывается по-другому (см., например, [2]). Положим

ехр(-1 Е^ ерУе) . к 1 к

Г = -Уф ^/7) , х; = ^-^^-/ , р е Вш, ¡к = £ хк = Ар Е агхг

Е ехр ( -1 Е «У ^ Ак Ак -о

V ееЕ /

1

реРы

Учитывая, что

- к

^аг (Уф (¿г/7) ,1

. 1

Ш1П < —

4€ёош ст* Ак

г=о

1

4€ёош <г* | \ Ак

£*е* (*е )}

ееЕ )

= - та^ — Е аг/г, <\ - Е ■* (*е)} = - Цк

=о Е Е

7ф (¿г/7) - (уф (¿г/7) , ¿г> = 7ф (¿г/7) + </г, ?) = - Е Е хр 1П {х;/й,ш),

р р

™еовреРт

к

Е ■ (- 7 А Е аг Е Е хр 1П (хр/йад) <

е€Е к г=0 теОИр^Ры

< - Е ■^ (- 7 Е Е хк 1П (хк/^),

ееЕ теовреРт

получаем

0 < [уф /7) + Е ■* ) } + | Е ■е ) +7 Е Е хкр 1П (хк/Лт) [ < ^,

I ееЕ ) \ееЕ теОИр^Ры ' к

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

Обратим внимание, что эта оценка не зависит от потенциально экспоненциально большого числа путей |Р| = п (п ^ 2л/™, например, для манхетенских сетей), в отличие от оценок, приведенных в пп. 2-4. Также обратим внимание, что эта оценка чувствительна к предельному переходу 7 ^ 0+, поскольку ¿2 ~ 7-1. В случае малых 7 можно воспользоваться подходом п. 3 [2], который можно представить для данной задачи следующим образом. Используем рандомизированный метод зеркального спуска, считая вместо

Уф(Ь/7) стохастический градиент: согласно распределению вероятностей / Е вы-

теои

бираем ■ е ОБ. Первый раз это делается за О (|ОБ|), а все последующие разы за О (1пт), затем независимо считаем стохастический градиент фт(¿/7). Это можно сделать за О(т) при условиях возможности упорядочивания (см. выше). Опишем, следуя п. 5.4 [32], соответствующую процедуру. Пусть ■ = (г, V). Считаем за О (т) (см. выше) ф^ (¿/7) ,к е У. Далее используем следующий подход, базирующийся на формуле полной вероятности. Стартуем из вершины , выбирая ребро, по которому пойдем согласно распределению вероятностей

Р (е = Ц ^ к)) = ехр«у/^)/7> ,к : (, ^ к) е Е.

ехр (фгУ ^/7)/7)

По аналогичному принципу выбираем ребро, по которому пойдем из вершины к и т.д.

Таким образом, мы сгенерируем случайный путь р е Ргь. Поток ^ по этому пути

теои

определяет потоки на ребрах, входящих в этот путь (на остальных ребрах потоки нулевые). Так задаваемый вектор потоков по ребрам и будет несмещенной оценкой -Уф (¿/7).

Оценка сложности (в среднем) метода получается аналогично п. 3 [2]: О [тМ^Щ/е2^, где = Н • I Е dw I , а О(-) = О(-) с точностью до логарифмического множителя. Более

\weOD )

того, эта оценка имеет и такой же вид, как оценка для аналога этого метода, приведенного в п. 3 [2] при 7 = 0. Описанный подход конкурентно способен с приведенным выше в этом пункте подходом лишь при малых значениях параметра 7 > 0. Можно явно выписывать здесь оценку на 7 > 0 подобно п. 3, но далее мы это сделаем в более интересном случае.

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

7 = 2 £ dw 1п ^ |

weOD

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

выше) О(Бту^ К/е2). Мы немного завысили настоящую оценку, чтобы можно было привести получающуюся оценку к максимально удобному для сравнения виду. В зависимости от того, что больше: О БтМ^Щ/е2^ или О (^Бт^М^К2/, следует принимать решение о том вводить или нет регуляризацию.

Исследование в п. 5 выполнено ИППИ РАН за счет гранта Российского научного фонда (проект № 14-50-00150). Исследование в других пунктах выполнено при поддержке грантов РФФИ 15-31-70001-мол_а_мос, 15-31-20571-мол_а_вед.

Литература

1. Гасников А.В., Дорн Ю.В., Нестеров Ю.Е, Шпирко С.В. О трехстадийной версии модели стационарной динамики транспортных потоков // Математическое моделирование. 2014. Т. 26:6. C. 34-70. arXiv:1405.7630

2. Гасников А.В., Двуреченский П.Е., Дорн Ю.В., Максимов Ю.В. Численные методы поиска равновесного распределения потоков в модели Бэкмана и модели стабильной динамики // Математическое моделирование. 2016. Т. 28. (принята к печати). arXiv:1506.00293

3. Sheffi Y. Urban transportation networks: Equilibrium analysis with mathematical programming methods. N.J.: Prentice-Hall Inc., Englewood Cliffs, 1985.

4. Введение в математическое моделирование транспортных потоков. 2-е изд. / под ред. А.В. Гасникова. М.: МЦНМО, 2013. 427 с. http://www.mou.mipt.ru/gasnikov1129.pdf

5. Sandholm W. Population games and Evolutionary dynamics. Economic Learning and Social Evolution. MIT Press; Cambridge, 2010.

6. Andersen S.P., de Palma A., Thisse J.-F. Discrete choice theory of product differentiation. MIT Press; Cambridge, 1992.

7. Nesterov Yu. Gradient methods for minimizing composite functions // Math. Prog. 2013. V. 140. N 1. P. 125-161.

http://www.uclouvain.be/cps/ucl/doc/core/documents/Composit.pdf

8. Nesterov Yu., Nemirovski A. On first order algorithms for l\/ nuclear norm minimization // Acta Numerica. 2013. V. 22. P. 509-575.

9. Nesterov Y. Characteristic functions of directed graphs and applications to stochastic equilibrium problems // Optim. Engineering. 2007. V. 8. P. 193-214.

10. Гасников А.В. Об эффективной вычислимости конкурентных равновесий в транспортно-экономических моделях // Математическое моделирование. 2015. Т. 27, № 12. С. 121-136. (принята к печати). arXiv:1410.3123

11. Бабичева Т.С., Гасников А.В., Лагуновская А.А., Мендель М.А. Двухстадийная модель равновесного распределения транспортных потоков // Труды МФТИ. 2015. Т. 7, № 3. С. 31-41.

12. Гасников А.В., Двуреченский П.Е., Камзолов Д.И., Нестеров Ю.Е., Спокойный В.Г., Стецюк П.И., Суворикова А.Л., Чернов А.В. Поиск равновесий в многостадийныйх транспортных моделях // Труды МФТИ. 2015. Т. 7, № 4. С. 143-155. (принята к печати). arXiv:1506.00292

13. Ethier N.S., Kurtz T.G. Markov processes. Wiley Series in Probability and Mathematical Statistics: Probability and Mathematical Statistics. John Wiley & Sons Inc., New York, 1986.

14. Малышев В.А., Пирогов С.А. Обратимость и необратимость в стохастической химической кинетике // Успехи мат. наук. 2008. Т. 63, вып. 1(379). С. 4-36.

15. Гасников А.В., Гасникова Е.В. Об энтропийно-подобных функционалах, возникающих в стохастической химической кинетике при концентрации инвариантной меры и в качестве функций Ляпунова динамики квазисредних // Математические заметки. 2013. Т. 94, № 6. С. 816-824.

16. Баймурзина Д.Р.,Гасников А.В., Гасникова Е.В. Теория макросистем с точки зрения стохастической химической кинетики // Труды МФТИ. 2015. Т. 7, № 4. C. 95-103. (принята к печати).

17. Санов И.Н. О вероятности больших отклонений случайных величин // Матем. сб. 1957. Т. 42(84):1. C. 11-44.

18. Patriksson M. The traffic assignment problem. Models and methods. Utrecht, Netherlands: VSP, 1994.

19. Ortuzar J.D., Willumsen L.G. Modelling transport. JohnWilley & Sons, 2011.

20. Гасников А.В., Гасникова Е.В., Нестеров Ю.Е., Чернов А.В. Об эффективных численных методах решения задач энтропийно-линейного программирования // ЖВМ и МФ. 2016. Т. 56, № 4. arXiv:1410.7719 (принята к печати).

21. Anikin A., Dvurechensky P., Gasnikov A., Golov A., Gornov A., Maximov Yu., Mendel M., Spokoiny V. Modern efficient numerical approaches to regularized regression problems in application to traffic demands matrix calculation from link loads // Proceedings of International conference ITAS-2015. Russia, Sochi, September, 2015. arXiv:1508.00858

22. Гасников А.В., Двуреченский П.Е., Усманова И.Н. О нетривиальности быстрых (ускоренных) рандомизированных методов. e-print, 2015. arXiv:1508.02182

23. Гасников А.В., Двуреченский П.Е., Нестеров Ю.Е. Стохастические градиентные методы с неточным оракулом // e-print, 2015. arXiv:1411.4218

24. Allen-Zhu Z., Orecchia L. Linear coupling: An ultimate unification of gradient and mirror descent // e-print, 2014. arXiv:1407.1537

25. Nesterov Yu. Universal gradient methods for convex optimization problems // Mathematical Programming. 2015. V. 152. С. 381-404.

26. Гасников А.В, Двуреченский П.Е., Камзолов Д.И. Градиентные и прямые методы с неточным оракулом для задач стохастической оптимизации // Труды SDCP-2014. Екатеринбург, 2015. arXiv:1502.06259

27. Лидбеттер М., Линдгрен Г., Ротсен Х. Экстремумы случайных последовательностей и процессов. М.: Мир, 1989.

28. Ким К., Нестеров Ю., Скоков В., Черкасский Б. Эффективные алгоритмы для дифференцирования и задачи экстремали // Экономика и математические методы. 1984. Т. 20. С. 309-318.

29. Форд Л., Фалкерсон Д. Потоки в сетях. М.: Мир, 1966.

30. Ahuja R.K., Magnati T.L., Orlin J.B. Network flows: Theory, algorithms and applications. Prentice Hall, 1993.

31. Kolokoltsov V.N., Maslov V.P. Idempotent analysis and applications. Dordrecht: Kluwer Acad. Publ., 1997.

32. Lugosi G., Cesa-Bianchi N. Prediction, learning and games. New York: Cambridge University Press, 2006.

References

1. Gasnikov A., Dorn Yu, Nesterov Yu., Shpirko S. On the three-stage version of stable dynamic model. Matem. Mod., V. 26:6 (2014). C. 34-70. arXiv:1405.7630

2. Gasnikov A.V., Dvurechensky P.E., Dorn Y.V., Maksimov Y.V. Searching equillibriums in Beckmann's and Nesterov-de Palma's models. Mathematical modelling. 2016. V. 28. (In Russian). arXiv:1506.00293 (accepted for publication).

3. Sheffi Y. Urban transportation networks: Equilibrium analysis with mathematical programming methods. N.J.: Prentice-Hall Inc., Englewood Cliffs, 1985.

4. Introduction to the mathematical modeling of traffic flows. Ed. Gasnikov A.V., with a foreword by hand of the Department of Transport in Moscow Liksutov M.S. M.: MCCME, 2013. (In Russian). http://www.mou.mipt.ru/gasnikov1129.pdf

5. Sandholm W. Population games and Evolutionary dynamics. Economic Learning and Social Evolution. MIT Press; Cambridge, 2010.

6. Andersen S.P., de Palma A., Thisse J.-F. Discrete choice theory of product differentiation. MIT Press; Cambridge, 1992.

7. Nesterov Yu. Gradient methods for minimizing composite functions. Math. Prog. 2013. V. 140, N 1. P. 125-161.

http://www.uclouvain.be/cps/ucl/doc/core/documents/Composit.pdf

8. Nesterov Yu., Nemirovski A. On first order algorithms for l\/ nuclear norm minimization // Acta Numerica. 2013. V. 22. P. 509-575.

9. Nesterov Y. Characteristic functions of directed graphs and applications to stochastic equilibrium problems. Optim. Engineering. 2007. V. 8. P. 193-214.

10. Gasnikov A.V. About reduction of searching competetive equillibrium to the minimax problem in application to different network problems. Mathematical modelling. 2015. V. 27, N 12. P. 121-136. (In Russian). arXiv:1410.3123 (accepted for publication).

11. Babicheva T.S., Gasnikov A.V, Lagunovskaya A.A., Mendel M.A. Two-stage model of equilibrium distributions of traffic flows. Proceedings of MIPT. 2015. V. 4. P. 31-41.

12. Gasnikov A.V., Dvurechensky P.E., Kamzolov D.I., Nesterov Y.E., Spokoyny V.G., Stetsyuk P.I., Suvorikova A.L., Chernov A.V. Searching equillibriums in multistage transport problems. Proceedings of MIPT. 2015. V. 4. P. 143-155. (accepted for publication). arXiv:1506.00292

13. Ethier N.S., Kurtz T.G. Markov processes. Wiley Series in Probability and Mathematical Statistics: Probability and Mathematical Statistics. John Wiley & Sons Inc., New York, 1986.

14. Malyshev V.A., Pirogo S.A. Reversibility and irreversibility in stochastic chemical kinetics. Uspekhi Mat. Nauk, 63:1(379) (2008), 3-36.

15. Gasnikov A.V., Gasnikova E.V. On Entropy-Type Functionals Arising in Stochastic Chemical Kinetics Related to the Con-centration of the Invariant Measure and Playing the Role of Lyapunov Functions in the Dy-namics of Quasiaverages. Mat. Zametki, 94:6 (2013), 819-827.

16. Baimurzina D.R., Gasnikov A.V., Gasnikova E.V. Theory of macrosystems from the point of view of stochastic chemical kinetics. Proceedings of MIPT. 2015. V. 7, N 4. 95-103. (accepted for publication).

17. Sanov I.N. On the probability of large deviations of random magnitudes. Mat. Sb. (N.S.), 42(84):1 (1957), 11-44.

18. Patriksson M. The traffic assignment problem. Models and methods. Utrecht, Netherlands: VSP, 1994.

19. Ortuzar J.D., Willumsen L.G. Modelling transport. JohnWilley & Sons, 2011.

20. Gasnikov A.V., Gasnikova E.V., Nesterov Y.E., Chernov A.V. Entropy linear programming. Computational Mathematics and Mathematical Physics. 2016. V. 56, N 4. (In Russian). arXiv:1410.7719 (accepted for publication)

21. Anikin A., Dvurechensky P., Gasnikov A, Golov A., Gornov A., Maximov Yu., Mendel M., Spokoiny V. Modern efficient numerical approaches to regularized regression problems in application to traffic demands matrix calculation from link loads. Proceedings of International conference ITAS-2015. Russia, Sochi, September, 2015. arXiv:1508.00858

22. Gasnikov A.V., Dvurechensky P.E., Usmanova I.N. About accelerated randomized methods. e-print, 2015. arXiv:1508.02182

23. Gasnikov A., Dvurechensky P., Nesterov Yu. Stochastic gradient methods with inexact oracle. e-print, 2015. arXiv:1411.4218

24. Allen-Zhu Z, Orecchia L. Linear coupling: An ultimate unification of gradient and mirror descent. e-print, 2014. arXiv:1407.1537

25. Nesterov Yu. Universal gradient methods for convex optimization problems. Mathematical Programming. 2015. V. 152. P. 381-404.

26. Gasnikov A.V., Dvurechensky P.E., Kamzolov D.I. Gradient and direct methods with inexact oracle for stochastic optimization problems. Proceedings of SDCP-2014. Ekaterinburg, 2015. arXiv:1502.06259

27. Leadbetter M.R., Lindgren G., Rootzen H. Extremes and Related Properties of Random Sequences and Processes. Springer, 1983.

28. Kim K., Nesterov Yu., Skokov V., Cherkassky B. Efficient algorithms for the derivation and extremaly problem. Economics and mathematical methods. 1984. V. 20, 309-318.

29. Ford L.R., Fulkerson D.R. Flows in Networks. Princeton University Press, 1962.

30. Ahuja R.K., Magnati T.L., Orlin J.B. Network flows: Theory, algorithms and applications. Prentice Hall, 1993.

31. Kolokoltsov V.N., Maslov V.P. Idempotent analysis and applications. Dordrecht: Kluwer Acad. Publ., 1997.

32. Lugosi G., Cesa-Bianchi N. Prediction, learning and games. New York: Cambridge University Press, 2006.

Поступила в редакцию 16.11.2015.