Kernel-Based LMI Approaches to Solving the Hamilton-Jacobi-Bellman Equation and Nonlinear Optimal Control

1. Problem Statement

The paper addresses the fundamental challenge of solving the Hamilton–Jacobi–Bellman (HJB) equation for nonlinear optimal control problems.

• System: Control-affine systems of the form $\dot{x} = f(x) + g(x)u$ with an infinite-horizon cost functional.
• Challenge: The HJB equation is a nonlinear partial differential inequality (PDI). Direct solution is computationally intractable for general nonlinear systems due to the "curse of dimensionality."
• Limitations of Existing Methods:
  • Finite Difference/DP: Suffer from exponential computational growth with state dimension.
  • Neural Networks: Lack rigorous theoretical guarantees on optimality and stability.
  • Sum-of-Squares (SOS): Limited to polynomial systems and scale poorly with degree.
  • Operator Methods (Koopman/Perron-Frobenius): Often rely on RKHS frameworks for implementation but lack a direct convex formulation for the HJB inequality itself.

The goal is to develop a method that provides rigorous theoretical guarantees (optimality, stability, convergence) while remaining computationally tractable for nonlinear systems, including those with unstable equilibria.

2. Methodology

The authors propose a Kernel-Based Linear Matrix Inequality (LMI) approach that transforms the non-convex HJB problem into a convex Semidefinite Program (SDP). The methodology consists of four core components:

A. LMI Reformulation via Schur Complement

The nonlinear HJB inequality is converted into a linear matrix inequality (LMI) using the Schur Complement Lemma.

• The original HJB inequality involves a quadratic term in the gradient of the value function ($\nabla V$).
• By defining a block matrix $M(x)$, the authors show that the HJB inequality is equivalent to $M(x) \succeq 0$.
• Significance: This transforms the problem from a nonlinear PDI into a constraint that is linear in the optimization variables (the coefficients of the value function approximation), enabling convex optimization.

B. Reproducing Kernel Hilbert Space (RKHS) Approximation

To handle the infinite-dimensional nature of the HJB equation, the value function gradient $\nabla V(x)$ is approximated within a finite-dimensional subspace of an RKHS.

• Representation: The value function is approximated as $V(x) = \sum p_i \kappa(x, x_i)$, where $\kappa$ is a kernel function and $p_i$ are coefficients.
• Derivatives: The gradient and Hessian are computed analytically from the kernel expansion, allowing the LMI constraints to be expressed as linear functions of the coefficient vector $p$.
• Discretization: The infinite-dimensional LMI is reduced to a finite set of constraints via collocation at specific points in the state space.

C. The Riccati Hessian Constraint (Critical Innovation)

A major hurdle in HJB optimization is the existence of trivial solutions (e.g., $V(x) \equiv 0$), which satisfy the inequality but yield no control.

• Solution: The authors impose a hard constraint at the equilibrium point ($x=0$):
  $$\nabla^2 V(0) = P$$
  where $P$ is the solution to the Algebraic Riccati Equation (ARE) of the linearized system.
• Purpose:
  1. Prevents Triviality: Since $P \succ 0$ for stabilizable systems, $V \equiv 0$ is mathematically impossible.
  2. Ensures Consistency: Guarantees that the nonlinear controller behaves like the optimal Linear Quadratic Regulator (LQR) near the equilibrium.
  3. Stability Guarantee: Ensures local exponential stability with convergence rates predicted by linear theory.

D. Convex Semidefinite Programming (SDP)

The final optimization problem minimizes a norm of the coefficients (e.g., $\|p\|^2$) subject to:

1. Boundary conditions ($V(0)=0, \nabla V(0)=0$).
2. The Riccati Hessian constraint ($\nabla^2 V(0) = P$).
3. The LMI constraints at collocation points ($M(x_j) \succeq 0$).
   This results in a convex SDP that can be solved efficiently using standard solvers (e.g., SCS).

3. Key Contributions

1. Convex Reformulation: Established the equivalence between the nonlinear HJB inequality and a linear matrix inequality via Schur complement, enabling the use of convex optimization.
2. RKHS Framework: Demonstrated that value function gradients admit finite-dimensional kernel representations with rigorous approximation-theoretic guarantees.
3. Trivial Solution Elimination: Introduced the Riccati Hessian constraint, a novel mechanism to prevent trivial solutions while enforcing consistency with linearized optimal control theory.
4. Theoretical Guarantees: Provided complete proofs for:
   • Suboptimality Bounds: Bounding the cost difference between the approximate and optimal control.
   • Exponential Stability: Proving that closed-loop trajectories converge exponentially to the origin.
   • Convergence Rates: Establishing error bounds based on the fill distance of kernel centers.
5. Unified Framework: Showed that classical Lyapunov function construction (e.g., Giesl's method) is a special case of this framework when the control input is zero and the system is already stable.

4. Numerical Results

The method was validated on three distinct systems, demonstrating robustness and scalability:

• 1D Scalar Polynomial System:
  • A system with an unstable cubic drift ($\dot{x} = x + x^3 + u$).
  • Result: Successfully stabilized the unstable equilibrium from 6 different initial conditions. All trajectories converged exponentially to the origin ($|x| < 10^{-6}$) despite moderate global approximation errors.
• 2D Radially Symmetric System:
  • A 2D system with rotational symmetry and open-loop instability.
  • Result: Stabilized from 8 initial conditions distributed on a unit circle. Trajectories spiraled inward to the origin, confirming omnidirectional stability.
• Van der Pol Oscillator:
  • A classic non-polynomial, nonlinear oscillator with an unstable equilibrium.
  • Result: The kernel-based LMI approach successfully stabilized the system from 8 initial conditions. The method handled non-polynomial dynamics while maintaining convexity, with all trajectories converging to the origin ($|x| < 10^{-7}$).

Key Observation: In all cases, global approximation accuracy was moderate, but local accuracy (enforced by the Riccati constraint) was sufficient to guarantee global exponential stability.

5. Significance

This work represents a significant advancement in nonlinear control theory by bridging the gap between dynamic programming and modern convex optimization.

• Theoretical Rigor: Unlike many data-driven or learning-based control methods, this approach provides mathematical certificates for stability and suboptimality.
• Handling Instability: It uniquely addresses the problem of stabilizing unstable equilibria (control synthesis) rather than just verifying stability for already stable systems.
• Computational Tractability: By leveraging RKHS and SDP, it avoids the curse of dimensionality associated with grid-based methods and offers a scalable alternative to SOS programming for non-polynomial systems.
• Practical Applicability: The method is computationally efficient (solving 2D problems in seconds) and robust to various initial conditions, making it a viable candidate for real-world nonlinear control applications.