Slide 1

Slide 1 text

Optimal control problems in density spaces Wuchen Li Level Set Collective, 2017 Joint work with Yat-Tin Chow, Stanley Osher and Wotao Yin.

Slide 2

Slide 2 text

Motivation Optimal control problem of densities (histograms) play critical roles in image processing and mean field games, which are widely used in social Network, Biology species, Virus, Trading, Cancer and Congestion etc. 2

Slide 3

Slide 3 text

In this talk, we will design fast numerics towards the Mean field game system, focus on the following examples: Mean field optimal control problem; Earth Mover’s distance; Schr¨ odinger bridge problem. 3

Slide 4

Slide 4 text

Mean field: One to all, all for one. Strategy set: S = {C, D}; Players: Infinity, i.e. players form (ρC , ρD ) with ρC + ρD = 1; Payoffs: F(ρ) = (FC (ρ), FD (ρ))T = Wρ, where W = 3 0 2 2 , meaning a Deer worthing 6, a rabbit worthing 2. 4

Slide 5

Slide 5 text

Static game Population games by extending finite player games, model the strategic interactions in large populations of small, anonymous agents. E.g. Discrete static Strategy set S; Players (Simplex) P(S) = {(ρ(x))x∈S ∈ R|S| : x∈S ρ(x) = 1 , ρ(x) ≥ 0} ; Payoff function to strategy x, F(x, ·) : P(S) → R. Nash Equilibrium (NE): Players have no unilateral incentive to deviate from their current strategies. ρ∗ = (ρ∗(x))x∈S is a NE if ρ∗(x) > 0 implies that F(x, ρ∗) ≥ F(y, ρ∗) for all y ∈ S. 5

Slide 6

Slide 6 text

Variational approach A particular type of game, named potential games, are widely considered: There exists a potential F : P(S) → R, such that ∇ρ F(ρ) = F(ρ) . In potential games, from KKT condition, NE is the critical point of max ρ F(ρ) : ρ ∈ P(S) . Similar games can be formulated into differential games. 6

Slide 7

Slide 7 text

Differential games 7

Slide 8

Slide 8 text

Finite player potential games All players minimize the potential: inf X,u 1 N t 0 N i=1 L(Xi (s), ui (s)) − F(X1 (s), · · · , XN (s))ds + G(X1 (0), · · · , XN (0)) , where F, G are given potential, terminal functions, and the infimum is taken among all player i’s controls (strategy) vectors ui (s) and position Xi (s) d ds Xi = ui (s) , 0 ≤ s ≤ t , Xi (t) = xi . 8

Slide 9

Slide 9 text

Mean field potential games If the number of players goes to infinity, and F, G satisfy certain symmetric properties, then one approximates the game by the following minimization problem: inf ρ,u t 0 { Td L(x, u(s, x))ρ(s, x)dx − F(ρ(s, ·)}ds + G(ρ(0, ·)) , where the infimum is taken among all vector fields u(s, x) and density ρ(s, x): ∂ρ ∂s + ∇ · (ρu) = 0 , 0 ≤ s ≤ t , ρ(t, ·) = ρ(·) . 9

Slide 10

Slide 10 text

Analogs E.g. t = 0 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 0.01 0.02 0.03 0.04 0.05 t = 1 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0.5 1 1.5 2 2.5 3 x 10−5 In above two systems, many similar structures have been discovered: Primal dual PDEs [Larsy, Lions]; Hamilton-Jacobi equation in probability set [Gangbo]. 10

Slide 11

Slide 11 text

Goal We plan to numerically solve the mean field optimal control problems. Difficulties Curse of dimensionality (Infinite dimension); Structure keeping spatial discretization (Time reversible). Main tools: Hopf-Lax formula overcome the curse of dimensionality1; Optimal transport on finite graphs2. 1Y.T. Chow, J. Darbon, S. Osher and W. Yin, Algorithm for Overcoming the Curse of Dimensionality For Time-Dependent Non-convex Hamilton-Jacobi Equations Arising From Optimal Control and Differential Games Problems, 2016. 2W. Li, E. Ryu, S. Osher, W. Yin and W. Gangbo, a parallel method for earth mover’s distance, 2017. 11

Slide 12

Slide 12 text

Discrete strategy set Strategy graph G = (S, E), S is the finite strategy set, E is the edge set; Probability set P(G) = {(ρi )i∈S | i∈S ρi = 1, ρi ≥ 0}; Discrete potential energy and Terminal condition: F, G : P(G) → R . 12

Slide 13

Slide 13 text

Minimal flux problem Denote m(s, x) = ρ(s, x)u(s, x). The variational problem forms inf ρ,u t 0 { Td L(x, m(s, x) ρ(s, x) )ρ(s, x)dx − F(ρ(s, ·)}ds + G(ρ(0, ·)) , where the infimum is taken among all flux function m(s, x) and density ρ(s, x): ∂ρ ∂s + ∇ · m = 0 , 0 ≤ s ≤ t , ρ(t, ·) = ρ(·) . 13

Slide 14

Slide 14 text

Transport on finite graphs To mimic the minimal flux problem, we consider the discrete flux function div(m)|i = 1 ∆x d v=1 (mi+ 1 2 ev − mi− 1 2 ev ) , and the cost functional L(m, ρ) =        i+ ev 2 ∈E L m i+ 1 2 ev g i+ 1 2 ev gi+ 1 2 ev if gi+ ev 2 > 0 ; 0 if gi+ ev 2 = 0 and mi+ ev 2 = 0 ; +∞ Otherwise . where gi+ 1 2 ev := 1 2 (ρi + ρi+ev ) is the discrete probability on the edge i + ev 2 ∈ E. The time interval [0, 1] is divided into N interval, ∆t = 1 N . 14

Slide 15

Slide 15 text

Discrete strategy Mean field games Consider the discrete optimal control system: ˜ U(t, ρ) := inf m,ρ N n=1 L(mn, ρn) − N n=1 F(ρn) + G(ρ0) where the minimizer is taken among {ρ}n i , {m}n i+ ev 2 , such that ρn+1 i − ρn i + ∆t · div(m)|i = 0 , ρN i = ρi . 15

Slide 16

Slide 16 text

Primal-Dual structure sup Φ inf m,ρ n L(mn, ρn)∆t − n ∆tF({ρ}n i ) + G({ρ}0 i ) + n i Φn i ρn+1 i − ρn i + ∆t · div(m)|i = sup Φ inf ρ inf m n  L(mn, ρn) + i+ ev 2 ∈E 1 ∆x (Φn i − Φn i+ev )mi+ 1 2 ev   ∆t − n ∆tF({ρ}n i ) + G({ρ}0 i ) + n i Φn i ρn+1 i − ρn i = sup Φ inf ρ − n i+ ev 2 ∈E H 1 ∆x (Φn i − Φn i+ev ) gi+ 1 2 ev ∆t − n ∆tF({ρ}n i ) + G({ρ}0 i ) + n i Φn i ρn+1 i − ρn i where H is the Legendre transform of L. 16

Slide 17

Slide 17 text

Example 1: Kinetic energy A typical Lagrangian is the kinetic energy L(x, u) = u 2 . Consider inf m,ρ t 0 Td m2(s, x) ρ(s, x) dx − F(ρ(s, ·))ds + G(ρ(0, ·)) such that ∂ρ(s, x) ∂s + ∇ · (m(s, x)) = 0 , ρ(t, ·) = ρ . 17

Slide 18

Slide 18 text

Hopf formula Following the primal-dual structure, we arrive at the Hopf formula (Application of state-dependent Hopf formula3): sup {Φi} i ΦN−1 i ρi − n ∆t F(ρn) − i ∂ ∂ρi F(ρn)ρn i − G∗({Φ}0 i ) s.t. ρn+1 i − ρn i + ∆t j∼i [∂H]( Φn i −Φn j ∆x )gn ij = 0 Φn i − Φn−1 i + ∆t 4 j∼i H( Φn i −Φn j ∆x ) + ∂ ∂ρi F(ρn) = 0 ρN+1 i = ˜ ρi ΦN i = Φi We apply the gradient descent method towards it. 3Y.T. Chow, J. Darbon, S. Osher and W. Yin, Algorithm for Overcoming the Curse of Dimensionality for State-dependent Hamilton-Jacobi equations, 2017. 18

Slide 19

Slide 19 text

Case 1 ρ, Φoptimal, ∇x Φoptimal −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0.5 1 1.5 2 2.5 3 x 10−5 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 −0.5 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 19

Slide 20

Slide 20 text

Case 1: Evolution of Density [t = 0] t = 0 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 0.01 0.02 0.03 0.04 0.05 t = 0.2 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 2 4 6 8 10 12 14 x 10−3 t = 0.4 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 0.5 1 1.5 2 2.5 3 x 10−3 t = 0.6 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 1 2 3 4 5 6 7 x 10−4 t = 0.8 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 2 4 6 8 10 12 14 x 10−5 t = 1 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0.5 1 1.5 2 2.5 3 x 10−5 20

Slide 21

Slide 21 text

Case 2 ρ, Φoptimal, ∇x Φoptimal −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 0.5 1 1.5 2 2.5 3 3.5 4 x 10−5 −3 −2 −1 0 1 2 −3 −2 −1 0 1 2 −0.5 0 0.5 1 1.5 2 2.5 3 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 21

Slide 22

Slide 22 text

Case 2: Evolution of Density t = 0 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 1 2 3 4 5 6 7 8 x 10−3 t = 0.2 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 0.5 1 1.5 2 2.5 3 3.5 x 10−3 t = 0.4 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 2 4 6 8 10 12 14 x 10−4 t = 0.6 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 x 10−4 t = 0.8 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 5 10 15 x 10−5 t = 1 −3 −2 −1 0 1 2 3 −3 −2 −1 0 1 2 3 0 0.5 1 1.5 2 2.5 3 3.5 4 x 10−5 22

Slide 23

Slide 23 text

Example 2: Earth Mover’s distance A special attention is paid into the homogenous degree one Lagrangian L(x, u) = u . Consider inf m,ρ 1 0 Td m(t, x) dxdt such that ∂ρ(t, x) ∂t + ∇ · (m(t, x)) = 0 , ρ(0, ·) = ρ0 , ρ(1, ·) = ρ1 . By Jensen’s inequality in time. Let ˜ m(x) = 1 0 m(t, x)dt, one minimizer is attached at a time independent optimization: inf ˜ m { Td ˜ m(x) dx: ∇ · ˜ m(x) + ρ1(x) − ρ0(x) = 0} This is an L1 minimization problem, which shares many similarities to the one in compressed sensing. 23

Slide 24

Slide 24 text

L1 Primal Dual system In this setting, the discretized minimization problem forms minimize m m subject to div(m) + p1 − p0 = 0 , We solve it by looking at its saddle point structure. Denote Φ = (Φi )N i=1 as a Lagrange multiplier: min m max Φ m + ΦT (div(m) + p1 − p0) . The iteration steps are as follows (using Chambolle and Pock): mk+1 = arg minm m + (Φk)T div(m) + m−mk 2 2 2µ ; Φk+1 = arg maxΦ ΦT div(2mk+1 − mk + p1 − p0) − Φ−Φk 2 2 2τ . 24

Slide 25

Slide 25 text

Algorithm: 2 line codes Primal-dual method for EMD 1. For k = 1, 2, · · · Iterates until convergence 2. mk+1 i+ 1 2 = shrink2 (mk i+ 1 2 + µ∇Φk i+ 1 2 , µ) ; 3. Φk+1 i = Φk i + τ{div(2mk+1 i − mk i ) + p1 i − p0 i } ; 4. End Here the shrink2 operator for the Euclidean metric is shrink2 (y, α) := y y 2 max{ y 2 − α, 0} , where y ∈ R2 . 25

Slide 26

Slide 26 text

Optimal flux I (c) (d) 26

Slide 27

Slide 27 text

Comparison Grids size EMD CUDA EMD CPU Ling Pele 32 × 32 0.012s 0.08s 0.007s 2.74s 64 × 64 0.063s 0.9s 0.009s N/A 128 × 128 0.336s 12.9s 2.3s N/A 256 × 256 6.8s 245.5s 80.8s N/A Table: Runtime of algorithms. 27

Slide 28

Slide 28 text

Example 3: Scrh¨ odinger bridge problem What is the optimal way to transport under white noise perturbations? 28

Slide 29

Slide 29 text

History remark 29

Slide 30

Slide 30 text

Problem formulation Schr¨ odinger in 1931 proposed one type of Mean field games: inf b 1 0 Rd 1 2 b2ρdx dt , where the infimum is among all drift vector fields b(t, x), such that ∂ρ ∂t + ∇ · (ρb) = β∆ρ , ρ(0) = ρ0 , ρ(1) = ρ1 . 30

Slide 31

Slide 31 text

Fisher Regularization The key idea (inherit from Nelson) is from the change of variables v = b − β∇ log ρ . Substituting the new v into the problem, inf v { 1 0 Rd 1 2 v2ρdx + β2 2 I(ρ) + β · D(ρ1|ρ0) dt : ∂ρ ∂t + ∇ · (ρv) = 0} , where D(ρ1|ρ0) = ρ1 log ρ1 − ρ0 log ρ0dx and the functional I(ρ) = (∇ log ρ)2ρdx , is called Fisher information. 31

Slide 32

Slide 32 text

Minimization The discrete minimization problem forms min m,p N n=1 i+ ev 2 ∈E { (mn i+ ev 2 )2 gn i+ ev 2 + β2 ∆x2 (log ρn i ρn i+ev )2gn i+ ev 2 } subject to ρn+1 i − ρn i ∆t + 1 ∆x d v=1 (mn i+ 1 2 ev − mn i− 1 2 ev ) = 0 ; ρi,0 = ρ0 i , ρi,N+1 = ρ1 i . Importance of Fisher information regularization: Boundary repeller; Enforces the strictly convexity. These two properties allow us to apply a simple Newton’s method. 32

Slide 33

Slide 33 text

Case 1 33

Slide 34

Slide 34 text

Case 2 34

Slide 35

Slide 35 text

References Yat-Tin Chow, Wuchen Li, Stanley Osher and Wotao Yin Numerics towards Hamilton-Jacobi equation in probability spaces, 2017. Y.T. Chow, J. Darbon, S. Osher and W. Yin Algorithm for Overcoming the Curse of Dimensionality for State-dependent Hamilton-Jacobi equations, 2017. Wuchen Li, Ryu Ernest, Stanley Osher, Wotao Yin and Wilfrid Gangbo A parallel algorithm to Earth Mover’s distance, 2017. Wuchen Li, Penghang Yin, and Stanley Osher Computations of optimal transport distance with Fisher information regularization, 2017. 35

Slide 36

Slide 36 text

Thanks! 36