A Convex Optimization Approach to Compute Trapping Regions for Lossless Quadratic Systems

TL;DR

This work addresses the challenge of certifying and characterizing boundedness for lossless quadratic systems by introducing a convex semidefinite programming framework that yields a necessary and sufficient trapping-region condition. If a trapping region exists, a second SDP computes the least-conservative spherical region, with strong duality enabling exact radius determination and identification of farthest points on the boundary with non-decreasing energy. The approach outperforms prior non-convex methods by providing guarantees, tight region estimates, and boundary-energy insights, and scales to systems with up to roughly $O(100)$ states. Numerical examples, including a two-dimensional system, the Lorenz attractor, and high-dimensional stacked Lorenz subsystems, demonstrate efficiency, accuracy, and scalability. The results offer a practical tool for modeling and control of lossless quadratic dynamical systems and motivate extensions to robust and data-informed settings.

Abstract

Quadratic systems with lossless quadratic terms arise in many applications, including models of atmosphere and incompressible fluid flows. Such systems have a trapping region if all trajectories eventually converge to and stay within a bounded set. Conditions for the existence and characterization of trapping regions have been established in prior works for boundedness analysis. However, prior solutions have used non-convex optimization methods, resulting in conservative estimates. In this paper, we build on this prior work and provide a convex semidefinite programming condition for the existence of a trapping region. The condition allows precise verification or falsification of the existence of a trapping region. If a trapping region exists, then we provide a second semidefinite program to compute the least conservative trapping region in the form of a ball. Two low-dimensional systems are provided as examples to illustrate the results. A third high-dimensional example is also included to demonstrate that the computation required for the analysis can be scaled to systems of up to $\sim O(100)$ states. The proposed method provides a precise and computationally efficient numerical approach for computing trapping regions. We anticipate this work will benefit future studies on modeling and control of lossless quadratic dynamical systems.

Figures (4)

• Figure 1: An illustration of trapping region for the two-state system \ref{['eq:TwoStateSys']} with kinetic energy $K_0(x)=\frac{1}{2}x^\top x$. All trajectories with random initial conditions (gray) converge to the equilibrium point at $(0,0.25)$ (magenta square). The set of states with non-decreasing energy $E$ (green solid ellipsoid) and a conservative trapping region $B(0,R_0)$ (blue-dashed circle) are shown. The tightest trapping region centered at the origin can be computed by our proposed method (Section \ref{['sect:Verif_Computing_beta']}) and is shown in the red-dashed circle.
• Figure 2: Illustrations of trapping regions of the bounded two-state system \ref{['eq:TwoStateSys_example']}. The left plot shows the ellipsoid with non-decreasing energy $E$, valid trapping regions $B(0,R_0)$ and $B(0,R_0^*)$, and critical points $x^*$. Note that the radius $R_0$ (in blue) is an over-estimation discussed in Section \ref{['sect:prel_LorenzCond']}, and $R_0^*$ (in red) is the tightest radius computed by the SDP-based method proposed in Section \ref{['sect:Verif_Computing_beta']}. The critical points $x^*$ (purple diamonds) are computed by the process proposed in Section \ref{['sect:Verif_ystar']}. The right plot shows the energy versus time evaluated along the state trajectories shown in the left plot. All energy trajectories decrease below $R_0^*$ in this system and well below $R_0$.
• Figure 3: Trapping regions of the Lorenz attractor computed by Theorem \ref{['theorem:TrappingRegion']} (blue) and the proposed SDP-based analysis (red). The left plot visualizes trajectories in the $x_2-x_3$ plane with $x_1=0$. The green ellipsoid is the region where energy is instantaneously non-decreasing. The right plot shows the energy versus time evaluated along the state trajectories shown in the left plot. All trajectories are ultimately bounded in the region computed by the proposed SDP-based analysis (red), which captures the boundedness behavior more accurately than Theorem \ref{['theorem:TrappingRegion']} (blue).
• Figure 4: Computation times $T$ for the proposed SDP-based analysis versus the number of states $n=3K$ in the high-dimensional system \ref{['eq:HighDimSys']}. The computation includes solving the SDP \ref{['SDP:min_lambda']} and SDP \ref{['eq:QCQPdual']}. The mean (blue) increases as number of states grows, and the standard deviation (gray) is negligible for large state dimensions. For large state dimensions, the execution time scales with $O(n^{3.98})$ and is close to the theoretical complexity $O(n^4)$.

Theorems & Definitions (4)

• definition 1: Boundedness
• definition 2: Trapping region
• theorem 1: Theorem 1 of Schlegel and Noack schlegel2015long
• theorem 2