Source-linked AI summary
Algorithms for Overcoming the Curse of Dimensionality for Certain Hamilton-Jacobi Equations Arising in Control Theory and Elsewhere
Jérôme Darbon, Stanley Osher
TL;DR
Grid-based HJ-PDE methods face exponential complexity in dimension, motivating alternatives for high-dimensional control and level-set problems. The paper uses Hopf formulas and convex optimization, reports grid-free polynomial-appearing algorithms with very fast evaluations, and notes derivative discontinuities at certain interfaces.
Problem
Grid-based HJ-PDE approximations are essentially impossible for n ≥4 because their computational complexity is exponential in dimension.
Method
The paper develops grid-free Hopf-formula algorithms, borrowing convex-optimization ideas, to compute HJ-PDE solutions, gradients, and closest points for convex-set geometries.
Results
The algorithms appear polynomial in dimension, use very low memory, are totally parallelizable, and evaluate solutions in 10^-4–10^-8 seconds on a laptop.
Takeaways & Limitations
The approach supports high-dimensional level-set, fast marching, fast sweeping, and closest-point computations without grids.
Abstract
from arXiv · showhide
It is well known that time dependent Hamilton-Jacobi-Isaacs partial differential equations (HJ PDE), play an important role in analyzing continuous dynamic games and control theory problems. An important tool for such problems when they involve geometric motion is the level set method. This was first used for reachability problems. The cost of these algorithms, and, in fact, all PDE numerical approximations is exponential in the space dimension and time. In this work we propose and test methods for solving a large class of the HJ PDE relevant to optimal control problems without the use of grids or numerical approximations. Rather we use the classical Hopf formulas for solving initial value problems for HJ PDE. We have noticed that if the Hamiltonian is convex and positively homogeneous of degree one that very fast methods exist to solve the resulting optimization problem. This is very much related to fast methods for solving problems in compressive sensing, based on $\ell_1$ optimization. We seem to obtain methods which are polynomial in the dimension. Our algorithm is very fast, requires very low memory and is totally parallelizable. We can evaluate the solution and its gradient in very high dimensions at $10^{-4}$ to $10^{-8}$ seconds per evaluation on a laptop. We carefully explain how to compute numerically the optimal control from the numerical solution of the associated initial valued HJ-PDE for a class of optimal control problems. We show that our algorithms compute all the quantities we need to obtain easily the controller. The term curse of dimensionality, was coined by Richard Bellman in 1957 when considering problems in dynamic optimization.
1. Introduction to Hopf Formulas, HJ PDEs and Level Set Evolutions
The paper addresses the exponential cost of grid-based HJ-PDE approximations by using Hopf formulas and convex optimization to compute solutions, gradients, and level-set evolutions without discretization. Examples cover distance functions, projections, and high-dimensional computations for convex and related nonconvex settings.
- Motivation: Grid-based numerical approximations become essentially impossible in dimensions n ≥4 because their complexity is exponential in n.Higher-order methods do not remove this curse of dimensionality.
- Hopf-formula approach: The initial-value HJ problem is solved without discretization using the Hopf formula and Fenchel-Legendre convex analysis.The resulting optimization problem can also provide the spatial gradient under additional convexity, differentiability, and coercivity conditions.
- Level-set evolutions: The zero level set represents points at prescribed Manhattan or Euclidean distance from a convex set when H is respectively ∥·∥1 or ∥·∥2.For the Euclidean case, this interpretation supports projection computations onto compact convex sets.
- Computational performance: The authors report computing selected HJ-PDE solutions in fairly high dimensions at below a millisecond per evaluation on a standard laptop.The examples are illustrative, with timing results presented separately.
- Convex optimization structure: For H = ∥·∥1, the optimization decouples into scalar problems whose unique minimizer is the soft thresholding operator.This yields an analytical solution in which the dimension n plays no significant role for the special case considered.
- Extensions: The approach generalizes to signed-distance computations for unions of compact convex sets and to more elaborate positively 1-homogeneous Hamiltonians.The paper also discusses merging fronts, closest points, and proximal points for certain nonconvex unions.
2. Introduction to Optimal Control
The paper connects optimal control problems with Hamilton-Jacobi PDEs and uses Hopf formulas to recover solutions and controls. It extends this framework to positively homogeneous Hamiltonians, nonconvex initial data, and level-set problems involving unions of convex sets.
- Optimal control and HJ-PDE: The value function for the control problem satisfies a Hamilton-Jacobi-Bellman equation with terminal value.
- Optimal control and HJ-PDE: The transformed value function is the viscosity solution of an HJ PDE with Hamiltonian H(∇xϕ).
- Optimal control and HJ-PDE: For differentiable solutions, the optimal control is recovered from the gradient through β(t) = ∇H(∇xϕ(x, T −t)).
- Optimal control and HJ-PDE: The algorithm computes both ϕ(x, t) and ∇xϕ(x, t), providing the quantities needed to obtain the optimal control without numerical gradient approximations.
- Extensions: For convex initial data and convex 1-coercive Hamiltonians, Hopf-Lax formulas solve the initial-value problem.
- Extensions: Taking pointwise minima of solutions extends the method to unions of convex sets, enabling level-set motion with merging fronts and closest-point computations.
3. Overcoming the Curse of Dimensionality for Convex Initial Data and Convex Homogeneous Degree One Hamiltonians – Optimal Control
The paper develops a grid-free algorithm for convex initial data and convex, positively 1-homogeneous Hamiltonians. It combines Hopf formulas with proximal and splitting methods to evaluate the solution, its gradient, and projections efficiently in high dimensions.
- HJ-PDE evaluation: Hopf formulas replace grid-based numerical approximations for solving the HJ-PDE in high dimensions.
- HJ-PDE evaluation: The method combines Hopf formulas, Moreau’s identity, and split Bregman iterations to evaluate the solution and its gradient.
- Performance: 10^-8 to 10^-4 seconds per evaluation is reported for the paper’s examples on a standard laptop, with apparent polynomial time complexity in n.
- HJ-PDE evaluation: When the minimizer is unique, the split Bregman iterates converge to ∇xϕ(x, t).
- Proximal computation: For convex positively 1-homogeneous H, the proximal map can be computed through projection onto its associated Wulff shape C.
- Projection and level sets: The approach also computes projections onto compact convex sets and supports closest-point calculations for level-set evolutions.
4. Numerical Results
The paper evaluates Hopf-formula-based methods for several Hamiltonians and initial data, combining explicit proximal-map formulas with projection and search procedures. Experiments report rapid solution evaluation on a single core and strong scaling across multiple cores.
- Algorithms: The ℓ_1 and ℓ_2 proximal-map cases are computed in linear time with respect to the dimension n.For the ℓ_∞ case, breakpoint sorting and bitonic search each take O(n log n).
- Experimental setup: The experiments cover ℓ_p Hamiltonians, positive-definite quadratic forms, and corresponding norm-based initial data.The tested matrices include diagonal positive-definite and symmetric positive-definite cases.
- Experimental setup: The numerical evaluation uses uniformly sampled (x, t) values in [-10, 10]^n × [0, 10] and averages 1,000,000 runs.Computations use IEEE double precision with denormalized-number mode disabled.
- Parallel performance: The multi-core implementation scales well because the algorithm requires little memory and fits in each processor core's L1 cache.The authors suggest this makes the approach suitable for low-energy embedded systems.
5. Conclusion
The paper concludes that Hopf-formula algorithms can rapidly solve certain convex, positively 1-homogeneous Hamilton-Jacobi equations while also computing gradients and related geometric quantities. The approach extends to level-set evolution, distance computation, and closest-point problems, with several nonconvex and polyhedral extensions identified for future work.
- Conclusion: The algorithms solve certain Hamilton-Jacobi equations rapidly, including cases with convex, positively 1-homogeneous Hamiltonians, while computing solution gradients.The gradient also provides the optimal control in the considered control problems.
- Conclusion: 10^-4–10^-8 seconds per evaluation is reported for fairly high-dimensional problems on a laptop, with low memory use and total parallelizability.The approach uses no grids, and its apparent complexity is polynomial in dimension with small constants.
- Projection and distance computation: An equally fast method finds a closest point in a finite union of compact convex sets with nonempty interiors.The paper also states that distances to these sets can be computed faster than Dijkstra-type fast methods.
- Level-set applications: The same framework can solve certain fast marching and fast sweeping problems in high dimensions through an HJ-PDE with prescribed boundary data.The construction uses a function ψ that is zero on the boundary and satisfies ∥∇xψ(x)∥2 = 1.
- Future work: Future work includes controls with x- and t-dependent Hamiltonians, polyhedral regions, nonconvex positively 1-homogeneous Hamiltonians, and differential-game problems.The proposed extensions include feasibility and closest-boundary-point computations.
Appendix A. Gauge
The appendix characterizes nonnegative, convex, positively 1-homogeneous Hamiltonians through gauges, support functions, and polar convex sets. Under symmetry, the gauge becomes a seminorm, while the polar set can be recovered from the Hamiltonian.
- A nonnegative Hamiltonian H is both the support function of C and the gauge of a characterized convex set.The gauge is defined for a closed convex set containing the origin.
- The gauge γ_D is lower semicontinuous, convex, and positively 1-homogeneous; if 0 lies in the interior of D, it is finite everywhere.
- If D is symmetric, γ_D is a seminorm.
- The polar-set relations give (Ω◦)◦ = Ω and identify the gauge γ_D with the support function of D◦.
- Taking D = C◦ expresses H as a gauge, with C◦ = {p ∈ Rn | H(p) ≤ 1}.