In this paper, we consider a large class of hierarchical congestion population games. One can show that the equilibrium in a game of such type can be described as a minimum point in a properly constructed multi-level convex optimization problem. We propose a fast primal-dual composite gradient method and apply it to the problem, which is dual to the problem describing the equilibrium in the considered class of games. We prove that this method allows to nd an approximate solution of the initial problem without increasing the complexity.
In this subsection, we brie y describe a variational principle for equilibrium description in hierarchical congestion population games. In particular, we consider a multistage model of traffic ows. Further details can be found in [1].
We consider the traffic network described by the directed graph 1 = ⟨V 1; E1⟩. Some of its vertices O1 V 1 are sources (origins), and some are sinks (destinations) D1 V 1. We denote a set of source-sink pairs by OD1 O1 D1. Let us assume that for each pair w1 2 OD1 there is a ow of network users of the amount of d1w1 := d1w1 M , where M ≫ 1, per unit time who moves from the origin of w1 to its destination. We call the pair w1; d1w1 as correspondence.
Let edges 1 be partitioned into two types E1 = E~1 ⨿E1. The edges of type E~1 are characterized by non-decreasing functions of expenses e11 (fe11 ) := e11 (fe11 =M ). Expenses e11 (fe11 ) are incurred by those users who use in their path an edge e1 2 E~1, the ow of users on this edge being equal to fe11 . The pairs of vertices setting the edges of type E1 are in turn a source-sink pairs OD2 (with correspondences d2w2 = fe11 , w2 = e1 2 E1) in a traffic network of the second level 2 = ⟨V 2; E2⟩whose edges are partitioned in turn into two types E2 = E~2 ⨿E2. The edges having type E~2 are characterized by non-decreasing functions of expenses e22 (fe22 ) := e22 (fe22 =M ). Expenses e22 (fe22 ) are incurred by those users who use in their path an edge e2 2 E~2, the ow of users on this edge being equal to fe22 .
The pairs of vertices setting the edges having type E2 are in turn source-sink pairs OD3 (with correspondences d3w3 = fe22 , w3 = e2 2 E2) in a traffic network of a higher level 3 = ⟨V 3; E3⟩, etc. We assume that in total there are m levels: E~m = Em. Usually, in applications, the number m is small and varies from 2 to 10.
Let Pw11 be the set of all paths in 1 which correspond to a correspondence w1. Each user in the graph 1 chooses a path p1w1 2 Pw11 (a consecutive set of the edges passed by the user) corresponding to his correspondence w1 2 OD1. Having de ned a path p1w1 , it is possible to restore unambiguously the edges having type E1 which belong to this path. On each of these edges w2 2 E1, user can choose a path pw2 2 Pw22 (Pw22 is a set of all paths corresponding in 2 the graph 2 to the correspondence w2), etc. Let us assume that each user have made the choice.
We denote by x1p1 the size of the ow of users on a path p1 2 P 1 = w12⨿OD1 Pw11 , x2p2 the size of the ow of users on a path p2 2 P 2 = w22⨿OD2 Pw22 , etc. Let us notice that
k xpkwk and that 0; pwk 2 Pwkk ; k
∑ pkwk 2Pwkk xpkwk = dkwk ; wk 2 ODk; k = 1; :::; m
k wk+1 (= ek)2 ODk+1 (= Ek); dkw+k+11 = fekk ; k = 1; :::; m 1:
For all k = 1; :::; m, we introduce for the graph a matrix k = ekpk ek2Ek;pk2P k ; ekpk =
k and the set of paths P k { 1; ek 2 pk 0; ek 2= pk : Then, for all k = 1; :::; m, the vector f k of ows on the edges of thegraph de ned in a unique way by the vector of ows on the paths xk = {xpkk }pk2P k : k is f k = kxk: We introduce the following notation x = {xk}m k=1 ; k=1 ; f = {f k}m = diag { k}m k=1 :
Let us now describe the probabilistic model for the choice of the path by a network user. We assume that each user l of a traffic network who uses a correspondence wk 2 ODk at a level k (and simultaniously the edge ek 1(= wk) 2 Ek 1 at the level k 1) chooses to use a path pk 2 Pwkk if pk = arg qkm2aPxwkk f gqkk (t) + qkk;lg; where qkk;l are iid random variables with double exponential distribution (also known as Gumbel's distribution) with cumulative distribution function where E 0:5772 is Euler{Mascheroni constant. In this case
P ( qkk;l < ) = expf e = k E ;
g M [ qkk;l] = 0;
D[ qkk;l] = ( k)2 2=6: Also, it turns out that, when the number of agents on each correspondence wk 2 ODk, k = 1; :::; m tends to in nity, i. e. M ! 1, the limiting distribution of users among paths is the Gibbs's distribution (also known as logit distribution) xpkk = dwk k
∑
exp( gpk~k (t)= k) ; pk 2 Pwkk ; wk 2 ODk; k = 1; :::; m:
(
[ For the sake of convenience we introduce the graph
m = ⨿ k=1
⟨ k =
V; E = m ⨿ E~k k=1 ⟩ and denote te = e(fe), e 2 E.
Assume that, for a given vector of expenses t on edges E, which is identical
to all users, each user chooses the shortest path at each level based on noisy
information and averaging of the information from the higher levels. Then, in
the limit number of users tending to in nity, such behavior of users leads to
the description of distribution of users on paths/edges given in (
Theorem 1 (Variational principle). The described above xed point equilibrium t can be found as a solution of the following problem (here and below we denote by dom the effective domain of the function conjugated to a function ) mf;ixnf (x; f ) : f = x; x 2 Xg =
min
t2dom
{ 1 1(t= 1) + ∑ e (te)}; (
} e(z)dz ; e(z)dz
= fe : te = e(fe); e 2 E; empm tem =
empm tem ; ∑ em2Em d e (te) = d dte
{ max fete dte fe ∑ ek2E~k gpmm (t) = ekpk tek fe ∫ 0 ∑ em2E~m ∑ ek2Ek pk2Pwkk 1(t) = gpkk (t) = ekpk k+1 ekk+1(t= k+1); k = 1; :::; m 1; wkk (t) = ln ( ∑ exp( gpkk (t))); k = 1; :::; m;
∑ In this subsection, we describe one of our contributions made by this paper, namely a general accelerated primal-dual gradient method for composite minimization problems.
We consider the following convex composite optimization problem [3]:
min [ϕ(x) := f (x) + (x)]:
x2Q
(
Here Q E is a closed convex set, the function f is differentiable and convex on Q, and function is closed and convex on Q (not necessarily differentiable).
In what follows we assume that f is Lf -smooth on Q: ∥∇f (x) ∇f (y)∥
Lf ∥x y∥ ; 8x; y 2 Q: We stress that the constant Lf > 0 arises only in theoretical analysis and not in the actual implementation of the proposed method. Moreover, we assume that the set Q is unbounded and that Lf can be unbounded on the set Q.
The space E is endowed with a norm ∥ ∥ (which can be arbitrary). The corresponding dual norm is ∥g∥ := maxx2Ef⟨g; x⟩ : ∥x∥ 1g; g 2 E . For mirror descent, we need to introduce the Bregman divergence. Let ! : Q ! R be a distance generating function, i.e. a 1-strongly convex function on Q in the ∥ ∥-norm: !(y) !(x) + ⟨!′(w); y
1 x⟩ + 2 ∥y x∥2 ; 8x; y 2 Q: Then, the corresponding Bregman divergence is de ned as
Vx(y) := !(y) !(x) ⟨!′(x); y x⟩; x; y 2 Q:
Finally, we generalize the Grad and Mirr operators from [2] to composite functions:
{ GradL(x) := argmin ⟨∇f (x); y y2Q
{ Mirrz (g) := argmin ⟨g; y y2Q z⟩ +
L x⟩ + 2 ∥y 1
} Vz(y) + (y) ;
}
x∥2 + (y) ; x 2 Q;
g 2 E ; z 2 Q:
(
Algorithm description Below is the proposed scheme of the new method. The main differences between this algorithm and the algorithm of [2] are as follows: 1) now the Grad and Mirr operators contain the (y) term inside; 2) now the algorithm does not require the actual Lipschitz constant Lf , instead it requires an arbitrary number L01 and automatically adapts the Lipschitz constant in iterations; 3) now we need to use a different formula for k+1 to guarantee convergence (see next section).
Note that Algorihtm 1 if well-de ned in the sense that it is always guaranteed that k 2 [0; 1] and, hence, xk+1 2 Q as a convex combination of points from Q. Indeed, from the formula for k+1 we have therefore k = k+1Lk+1 1 k+1Lk+1 1.
(√
1 4L2k+1 +
1 ) 2Lk+1
Lk+1 = 1;
(
Lk+1 2Lk+1 end while zk+1 Mirrzkk+1 (∇f (xk+1)) end for
return yT Therefore k+1⟨∇f (xk+1); zk
u⟩ = k+1⟨∇f (xk+1); zk
k+1⟨∇f (xk+1); zk + k+1⟨ ′(zk+1); u
zk+1⟩ ( k+1⟨∇f (xk+1); zk zk+1⟩
k+1 (zk+1))
+ ⟨Vz′k (zk+1); u
zk+1⟩ + k+1 (u);
zk+1⟩ + k+1⟨∇f (xk+1); zk+1
zk+1⟩ + ⟨Vz′k (zk+1); u
zk+1⟩
u⟩
(
Using the triangle equality of the Bregman divergence, ⟨Vx′(y); u y⟩ = Vx(u)
Vy(u)
we get
⟨Vz′k (zk+1); u zk+1⟩ = Vzk (u) Vzk+1(u) Vzk (zk+1)
1
Vzk (u) Vzk+1(u) 2 ∥zk+1 zk∥2 ;
(
So we have
12 ∥zk+1 zk∥2 in the last inequality. k+1⟨∇f (xk+1); zk (
u⟩
+ ( k2+1Lk+1
ϕ(xk+1) ϕ(yk+1)
(xk+1):
(
⊔⊓ (19) (20) (21) ⊔⊓ Lemma 2. For any u 2 Q and k = Now we apply Lemma 1 to bound the last term, group the terms and get k+1( (xk+1) (u)) + k+1⟨∇f (xk+1); xk+1 u ⟩ k+1ϕ(xk+1) ( k2+1Lk+1)ϕ(yk+1)
Now we are ready to prove the convergence theorem for Algorithm 1.
Theorem 2. For the sequence fykgk 0 in Algorithm 1 we have ( T2 LT )ϕ(yT )
k=1 { T min ∑ k (f (xk) + ⟨∇f (xk); u x2Q and, hence, the following rate of convergence: xk⟩ + (u)) + Vz0 (u) Proof. Note that the special choice of f kgk 0 in Algorithm 1 gives us ϕ(yT ) ϕ(x ) 4Lf R2
T 2
: 2 k+1Lk+1 k+1 = k2Lk; k 0: } (22) (23) (24) } (26) (27) (28) Therefore, taking the sum over k = 0; : : : ; T VzT (u) 0 we get, for any u 2 Q, ( T2 LT )ϕ(yT )
T ∑ k (f (xk) + ⟨∇f (xk); u k=1 1 in (19) and using that 0 = 0, xk⟩ + (u)) + Vz0 (u) (25) and (22) is straightforward. At the same time, using the convexity of f (x), the de nition of ϕ(x), and u = x = argminx2Q ϕ(x), we obtain ( T2 LT )ϕ(yT )
k=1 { T min ∑ k (f (xk) + ⟨∇f (xk); u x2Q ( T ) ∑ k ϕ(x ) + Vz0 (x ): k=1 From (24) it follows that ∑kT=1 k = T2 LT , so xk⟩ + (u)) + Vz0 (u) ϕ(yT ) ϕ(x ) +
Vz0 (x ): Now it remains to estimate the rate of growth of coefficients Ak := this we use the technique from [3]. Note that from (24) we have k2Lk. For Ak+1
Ak = 1 T2 LT √ Ak+1
Lk+1 Rearranging and using (a + b)2 2a2 + 2b2 and Ak
Ak+1, we get Ak+1 = Lk+1(Ak+1
Ak)2 = Lk+1 (√Ak+1 + p
Ak)2 (√Ak+1 p
Ak )2 4Lk+1Ak+1 (√Ak+1 p
Ak
)2 From this it follows that
k
1 ∑
2
i=0
p
1
Li
:
Note that according to (
2Lf . Hence, √Ak+1 k + 1 2√2Lf ()
Ak+1 (k + 1)2 8Lf : Thus, combining (31) and (27) with Vz0 (x ) =: R22 , we have proved (23).
Using the same arguments to [3], it is also possible to prove that the average number of evaluations of the function f per iteration in Algorithm 1 equals 4. Theorem 3. Let Nk be the total number of evaluations of the function f in Algorithm 1 after the rst k iterations. Then for any k 0 we have Nk 4(k + 1) + 2 log2 L0
: Lf 3
Application to the Equilibrium Problem
In this section, we apply Algorithm 1 to solve the dual problem in (
The inequality (22) leads to the fact that Algorithm 1 is primal-dual [6{9], which means that the sequences ftig (which is in the role of fxkg) and ft~ig (which is in the role of fykg) generated by this method have the following property: (30) (31)
⊔⊓ (32) min t2dom { 1
T 1 1(t~T = 1) + ∑ e2E
e (t~eT ) ∑ [ i( 1 1(ti= 1) + ⟨∇ 1(ti= 1); t ti⟩)] + ∑ e2E e (te) } (33) where 4L2R22
T 2
; k )
∑ with lw1 being the total number of edges (among all of the levels) in the longest path for correspondence w1, R22 = maxfR~22; R^22g;
R~22 = (1=2) ∥t 4L2R22
T 2
; where xpkk;i = dwk k ∑ exp( gpk~k (ti)= k) ; xi = {xpkk;i}k=1;:::;m
pk2Pwkk ;wk2ODk ; pk 2 Pwkk ; wk 2 ODk; k = 1; :::; m; f T =
Theorem 2 provides the bound for the number of iterations in order to solve
the problem (
In general, it is worth noting that the approach of adding some arti cial vertices, edges, sources, sinks is very useful in different applications [12{14].
Acknowledgements. The research was supported by RFBR, research project No. 15-31-70001 mol a mos and No. 15-31-20571 mol a ved.