Toward Feedback Control of Densities in Nonlinear Systems Abhishek - PowerPoint PPT Presentation

Toward Feedback Control of Densities in Nonlinear Systems Abhishek Halder Department of Applied Mathematics University of California, Santa Cruz Santa Cruz, CA 95064 Joint work with Kenneth F. Caluya (UC Santa Cruz) CDC 2019, Workshop: Uncertainty Synthesis, Nice

Overarching Theme Systems-control theory for densities

What is density?

Probability Density Fn. n 
 l ! # x $ ∈ X ≡ R 2 × S 1 " x ( t ) ∈ y θ

Probability Density Fn. n 
 l 
 l l 
 
 ol 
 l l ! # l 
 x $ ∈ X ≡ R 2 × S 1 " x ( t ) ∈ y l 
 θ ρ ( x , t ) : X × [ 0, ∞ ) $→ R ≥ 0 ! X ρ d x = 1 for all t ∈ [ 0, ∞ )

Probability Density Fn. Population Density Fn. n 
 n 
 l 
 l 
 l l 
 l l 
 
 ol 
 ol l 
 
 l l l l 
 ! # l 
 x 
 $ ∈ X ≡ R 2 × S 1 " x ( t ) ∈ y l 
 l θ ρ ( x , t ) : X × [ 0, ∞ ) $→ R ≥ 0 ! X ρ d x = 1 for all t ∈ [ 0, ∞ )

Why care about densities?

Prediction Problem Compute joint state PDF Initial Process noise ρ ( x , t ) conditions State density Process Parameters model ρ ( x ( t ) , t ) Trajectory flow: d x ( t ) = f ( x , t ) d t + g ( x , t ) d w ( t ) , d w ( t ) ∼ N ( 0, Q d t ) Density flow: !" $ gQg ⊤ # n ∂ 2 ∂ t = L FP ( ρ ) : = −∇ · ( ρ f ) + 1 ∂ρ ∑ ij ρ 2 ∂ x i ∂ x j i , j = 1

Filtering Problem Filtering Problem Initial Process noise conditions Prior density Process Compute conditional Parameters model ρ � ( x ( t ) , t ) joint state PDF + x t , t ρ + : = ρ ( x , t | z ( s ) , 0 ≤ s ≤ t ) Sensor noise Posterior density Measurement model ρ + ( x ( t ) , t ) Trajectory flow: d x ( t ) = f ( x , t ) d t + g ( x , t ) d w ( t ) , d w ( t ) ∼ N ( 0, Q d t ) d z ( t ) = h ( x , t ) d t + d v ( t ) , d v ( t ) ∼ N ( 0, R d t ) Density flow: % '( & ' ⊤ R − 1 & d ρ + = L FP d t + h ( x , t ) − E ρ + { h ( x , t ) } d z ( t ) − E ρ + { h ( x , t ) } d t ρ +

Control Problem lem Steer joint state PDF via feedback control over finite time horizon ! ! 1 " 0 ! u ! 2 minimize 2 d t E u ∈ U subject to d x = f ( x , u , t ) d t + g ( x , t ) d w , x ( t = 0 ) ∼ ρ 0 , x ( t = 1 ) ∼ ρ 1

PDFs in Mars Entry-Descent-Landing Prediction Problem Filtering Problem Control Problem Navigational uncertainty Heating uncertainty Landing footprint uncertainty Chute deployment uncertainty Image credit: NASA JPL Predict heating rate uncertainty

PDFs in Mars Entry-Descent-Landing Prediction Problem Filtering Problem Control Problem Navigational uncertainty Heating uncertainty Landing footprint uncertainty Chute deployment uncertainty Supersonic parachute Image credit: NASA JPL Predict heating rate uncertainty Estimate state to deploy parachute

PDFs in Mars Entry-Descent-Landing Prediction Problem Filtering Problem Control Problem Navigational uncertainty Heating uncertainty Landing footprint uncertainty Chute deployment uncertainty Supersonic parachute Image credit: NASA JPL Gale Crater (4.49S, 137.42E) Steer state PDF to achieve Predict heating rate uncertainty Estimate state to deploy parachute desired landing footprint accuracy

Solving prediction problem as Wasserstein gradient flow

What’s New? Main idea: Solve Infinite dimensional variational recursion: Proximal operator: Optimal transport cost: Free energy functional:

Geometric Meaning of Gradient Flow Gradient Flow in R n Gradient Flow in P 2 ( R n ) d x ∂ρ ∂ t = �r W Φ ( ρ ) , d t = �r ϕ ( x ) , x (0) = x 0 ρ ( x , 0) = ρ 0 Recursion: Recursion: x k = x k � 1 � h r ϕ ( x k ) ρ k = ρ ( · , t = kh ) ⇢ 1 � ⇢ 1 � 2 k x � x k � 1 k 2 2 W 2 ( ρ , ρ k � 1 ) + h Φ ( ρ ) = arg min 2 + h ϕ ( x ) = arg min x 2 R n ρ 2 P 2 ( R n ) =: prox k · k 2 =: prox W 2 h ϕ ( x k � 1 ) h Φ ( ρ k � 1 ) Convergence: Convergence: x k ! x ( t = kh ) h # 0 as ρ k ! ρ ( · , t = kh ) h # 0 as as Lyapunov function: as Lyapunov functional:

Geometric Meaning of Gradient Flow Gradient Flow in R n Gradient Flow in P 2 ( R n ) ∂ρ d x ∂ t = �r W Φ ( ρ ) , d t = �r ϕ ( x ) , x (0) = x 0 ρ ( x , 0) = ρ 0 D Recursion: Recursion: x k = x k � 1 � h r ϕ ( x k ) ρ k = ρ ( · , t = kh ) ⇢ 1 � ⇢ 1 2 k x � x k � 1 k 2 2 W 2 ( ρ , ρ k � 1 ) + h Φ ( ρ ) = arg min 2 + h ϕ ( x ) = arg min x 2 R n ρ 2 P 2 ( R n ) =: prox k · k 2 =: prox W 2 h ϕ ( x k � 1 ) h Φ ( ρ k � 1 ) Convergence: Convergence: x k ! x ( t = kh ) h # 0 as ρ k ! ρ ( · , t = kh ) h # 0 as

Algorithm: Gradient Ascent on the Dual Space Uncertainty propagation via point clouds No spatial discretization or function approximation

Algorithm: Gradient Ascent on the Dual Space Algorithm: Gradient Ascent on the Du @⇢ @ t = r · ( r ⇢ ) + � � 1 ∆ ⇢ m Proximal Recursion ⇢ 1 � 2 W 2 ( ⇢ , ⇢ k � 1 ) + h Φ ( ⇢ ) ⇢ k = ⇢ ( x , t = kh ) = arg inf ⇢ 2 P 2 ( R n )

Algorithm: Gradient Ascent on the Dual Space Algorithm: Gradient Ascent on the Du @⇢ @ t = r · ( r ⇢ ) + � � 1 ∆ ⇢ m Proximal Recursion ⇢ 1 � ⇢ � 2 W 2 ( ⇢ , ⇢ k � 1 ) + h Φ ( ⇢ ) ⇢ k = ⇢ ( x , t = kh ) = arg inf ⇢ 2 P 2 ( R n ) 2 P + Discrete Primal Formulation ⇢ � 1 2 h C k , M i + h h k � 1 + � � 1 log % , % i % k = arg min min M 2 Π ( % k − 1 , % ) %

Algorithm: Gradient Ascent on the Dual Space Algorithm: Gradient Ascent on the Du @⇢ @ t = r · ( r ⇢ ) + � � 1 ∆ ⇢ m Proximal Recursion ⇢ 1 � ⇢ � 2 W 2 ( ⇢ , ⇢ k � 1 ) + h Φ ( ⇢ ) ⇢ k = ⇢ ( x , t = kh ) = arg inf ⇢ 2 P 2 ( R n ) 2 P + Discrete Primal Formulation ⇢ � 1 2 h C k , M i + h h k � 1 + � � 1 log % , % i % k = arg min ⇢ min � M 2 Π ( % k − 1 , % ) % + Entropic Regularization ⇢ � 1 ⇢ 2 h C k , M i + ✏ H ( M ) + h h k � 1 + � � 1 log % , % i � % k = arg min min M 2 Π ( % k − 1 , % ) % m Dualization ⇢ � opt , � opt h � 0 , % k � 1 i � F ? ( � � 1 ) = arg max 0 1 � 0 , � 1 � 0 ✓ ◆� � ✏ exp( � > 0 h / ✏ ) exp( � C k / 2 ✏ ) exp( � 1 h / ✏ ) h

Recursion on the Cone Theorem: Consider the recursion on the cone R n � 0 ⇥ R n � 0 ⇣ ⌘ = ⇠ k � 1 � z � �✏ > y y � ( Γ k z ) = % k � 1 , Γ k z � h , Then the solution ( y ⇤ , z ⇤ ) gives the proximal update % k = z ⇤ � ( Γ k > y ⇤ )

Algorithmic Setup p Theorem: Block co-ordinate iteration of ( y , z ) recur- sion is contractive on R n > 0 × R n > 0 . the from

Proximal Prediction: 2D Linear Gaussian Proximal Prediction: 2D Linear Gaussian

Proximal Prediction: Nonlinear Non-Gaussian

Computational Time: Nonlinear Non-Gaussian Computational time (seconds) 10 − 6 1 2 3 4 Physical time t k = kh (seconds)

Proximal Prediction: Satellite in Geocentric Orbit 0 1 v x v y B C B C v z 0 1 0 1 0 d x B C B C 0 d y B C B C B C − µ x B C B C B C r 3 + ( f x ) pert − γ v x 0 d z B C p B C B C 2 β − 1 γ = d t + B C , B C B C d v x d w 1 B C B C B C B C B C B C − µ y d v y d w 2 B C @ A @ A r 3 + ( f y ) pert − γ v y B C d v z d w 3 B C B C B C − µ z @ A r 3 + ( f z ) pert − γ v z 0 1 k 3( s θ ) 2 − 1 � � 0 1 0 1 f x s θ c φ c θ c φ − s φ 2 r 4 B C A , k := 3 J 2 R 2 − k = E , µ = constant f y s θ s φ c θ s φ c φ B C @ A @ A r 5 s θ c θ B C 0 f z c θ − s θ @ pert 0

Computational Time: Satellite in Geocentric Orbit Computational time (seconds) 10 − 7 10 − 8 0.005 0.01 Physical time t k = kh (seconds)

Details on Proximal Prediction Extensions: mean-field models for nonlocal interaction, state-dependent diffusions Publications: — K.F. Caluya, and A.H., Proximal Recursion for Solving the Fokker-Planck Equation, ACC 2019 . — K.F. Caluya, and A.H., Gradient Flow Algorithms for Density Propagation in Stochastic Systems, IEEE Trans. Automatic Control 2020, doi: 10.1109/TAC.2019.2951348. Git repo: github.com/kcaluya/UncertaintyPropagation

Solving density control as Wasserstein gradient flow

Finite Horizon Feedback Density Control lem ! ! 1 " 0 ! u ( x , t ) ! 2 minimize 2 d t E u ∈ U subject to # $ √ d x = f ( x , t )+ B ( t ) u ( x , t ) d t + 2 " B ( t ) d w , x ( t = 0 ) ∼ ρ 0 , x ( t = 1 ) ∼ ρ 1

Finite Horizon Feedback Density Control lem ! ! 1 " 0 ! u ( x , t ) ! 2 minimize 2 d t E u ∈ U subject to # $ √ d x = f ( x , t )+ B ( t ) u ( x , t ) d t + 2 " B ( t ) d w , x ( t = 0 ) ∼ ρ 0 , x ( t = 1 ) ∼ ρ 1 Necessary conditions for optimality: coupled nonlinear PDEs (FPK + HJB) % ρ opt % && ∂ρ opt = % 1 ⊤ ' ' ρ opt (( ⊤ ∇ ψ f + B ( t ) D ( t ) ⊙ Hess 1 , + ∇ · ∂ t ∂ t + 1 ∂ψ ⊤ ∇ ψ ) 2 2 ) B ( t ) 2 + 〈∇ ψ , f 〉 = − % 〈 D ( t ) , Hess ( ψ ) 〉 Boundary conditions: Optimal control: ∇ ρ opt ( x , 0 ) = ρ 0 ( x ) , ρ opt ( x , 1 ) = ρ 1 ( x ) u opt ( x , t ) = B ( t ) ⊤ ∇ ψ