Full-Order Safety for Nonlinear Reduced-Order MPC
Reduced-order models (ROMs) make model predictive control (MPC) of high-dimensional nonlinear systems tractable, but the resulting modeling error can lead to constraint violations on the real system. We introduce CoRE-MPC, a robust MPC formulation that plans with a nonlinear ROM and guarantees recursive feasibility and constraint satisfaction for the full-order system.
Aus dem Abstract. Unsere Zusammenfassung folgt.
Vollständige Studie
Full-Order Safety for Nonlinear Reduced-Order MPC
Veröffentlicht unter CC BY 4.0. Wiedergabe mit Namensnennung, Original unter arXiv:2610.19250 (PDF).
Abstract
Reduced-order models (ROMs) make model predictive control (MPC) of high-dimensional nonlinear systems tractable, but the resulting modeling error can lead to constraint violations on the real system. We introduce CoRE-MPC, a robust MPC formulation that plans with a nonlinear ROM and guarantees recursive feasibility and constraint satisfaction for the full-order system. It can be applied with a broad class of nonlinear model reduction techniques and the resulting MPC does not require evaluating the full-order model online, thus ensuring computational efficiency. The key insight is in how the ROM affects the dynamics of the prediction error. We account for this with two error radii that are propagated alongside the ROM and provide a guaranteed bound on the prediction error. The approach is validated in simulation on an 80-state oscillator chain, a Burgers equation, and a 9{,}768-state finite-element soft robot, where a two-state ROM with data-driven bounds keeps all 1{,}628 nodes clear of a wall at a median solve time below 1\,\mathrm{ms} for a 20\,\mathrm{ms} control period. Across these examples, CoRE-MPC is 13 to 61 times faster than full-order robust MPC at only 27 to 35\% higher closed-loop cost, whereas existing certified schemes are either infeasible or at least twice as conservative. Code: [released upon publication].
9{,}768-state finite-element simulation of a soft robot, and the wire shells show a two-state ROM prediction of it at the same instant (t=1.64\,\mathrm{s}). With the nominal model, ROM-MPC (gray) predicts 0.1\,\mathrm{mm} clearance but actually penetrates the wall by 13.4\,\mathrm{mm} (red) as it fails to account for the prediction error of the ROM. By separately bounding the error along the approximation manifold and off it, including their coupling, CoRE-MPC (blue) keeps its predicted error envelope (shaded), and hence every node, clear of the wall with the same ROM, horizon, cost function, and solver. The dashed ellipse represents the common XY reference trajectory at the tip height and the star marks the current target along the trajectory.I Introduction
Reduced-order models (ROMs) enable computationally tractable predictive control for high-dimensional nonlinear systems such as fluid flows, flexible structures, soft robots, and humanoids [1, 2, 3, 4]. Safety, however, is a property of the real system, whose trajectory the ROM only approximates, and the two can differ significantly. For example, a ROM can predict a collision-free trajectory while the physical robot still hits an obstacle (Fig. 1). In practice, this mismatch is absorbed by hand-tuned constraint margins or by a larger ROM, which trade performance and computation for safety without guaranteeing it.
Robust MPC accounts for model error by bounding it as a function of the measured state and input and propagating the bound as a tube along the plan [5]. The prediction error of a ROM is different: it depends on the full state, of which only the reduced coordinates are measured, so its unmeasured part evolves with its own dynamics and feeds back into the reduced prediction.
For linear FOMs, the error dynamics are linear and reduced-order MPC with full-order guarantees is well developed, with uniform [6, 7, 8, 9] or input-dependent [10, 11] error bounds. Treating the nonlinearities as additional bounded disturbances extends these approaches to nonlinear dynamics, which is, however, prohibitively conservative (cf. Section V-A). Nonlinear ROMs, such as proper orthogonal decomposition (POD)–Galerkin projections [12, 2], quadratic manifolds [13], and spectral submanifolds [14], can provide a more accurate ROM for nonlinear dynamics. However, using these ROMs directly in an MPC formulation can only guarantee constraints on the ROM prediction, not the real system. Applying a standard nonlinear robust MPC [5] to ROMs treats the unmodelled dynamics as an independent model mismatch, resulting in overly conservative or even infeasible operation. To the best of our knowledge, the robust nonlinear reduced-order MPC (RN-ROMPC) in [15] is the only scheme that certifies full-order constraints using only ROM prediction. However, its application requires an invariant manifold, which is available only for specific systems near an equilibrium. Thus, there exists no MPC scheme that can ensure constraint satisfaction using the nonlinear ROMs that are common in high-dimensional soft robots.
Contribution: We present CoRE-MPC, a reduced-order MPC that enforces constraints on the full-order system for a general class of nonlinear ROMs. The prediction error within the reduced coordinates is measured, and can thus be reset at every update, but the error outside them is not. CoRE-MPC therefore bounds the two separately, as a tangential and a transverse error, together with their coupling.
- Error bound: A coupled reduction-error certificate (CoRE, Thm. 1) bounds the prediction error for general ROMs using a two-dimensional comparison system.
- Theoretical guarantees: Online, CoRE-MPC propagates two error radii with the ROM, at a computational cost comparable to nominal ROM-MPC with two extra states, and guarantees recursive feasibility and constraint satisfaction for the original full-order system (Thm. 5).
- Numerical experiments: We demonstrate CoRE-MPC on an oscillator chain, a Burgers equation, and the
9{,}768-state finite-element soft robot of Fig. 1, where nominal ROM-MPC violates the constraints. On the oscillator chain, all assumptions of Theorem 5 are satisfied. Compared to full-order robust MPC, CoRE-MPC is13to61times faster at27to35\%higher closed-loop cost. Even when invariant-manifold assumption of RN-ROMPC [15] holds, CoRE-MPC incurs less than a quarter of its excess closed-loop cost. On the soft robot, where full-order robust MPC is computationally prohibitive, CoRE-MPC keeps all1{,}628nodes clear of the wall at a median solve time below1ms.
Outline: Section II states the problem, Section III derives the error certificate, Section IV presents CoRE-MPC and its closed-loop guarantees, and Section V reports the experiments.
Notation: \|\cdot\| denotes the Euclidean norm and its induced matrix norm, and \langle\cdot,\cdot\rangle the Euclidean inner product. Inequalities between vectors or matrices hold componentwise. A map is C^{k} if it is k times continuously differentiable. D denotes the Jacobian and D_{x} the Jacobian with respect to x. D^{+} denotes the upper right Dini derivative, D^{+}v(t):=\limsup_{h\downarrow 0}\,[v(t{+}h)-v(t)]/h, which reduces to \dot{v} wherever v is differentiable.
II Problem Formulation
Consider the full-order model (FOM)
\dot{x}=F(x,u)+d,with state x\in\mathbb{R}^{N}, input u\in\mathbb{R}^{m}, C^{1} dynamics F, and an unknown piecewise continuous disturbance d with \|d(t)\|\leq\bar{d}. We impose the input and state constraints
u\in\mathcal{U},\qquad h(x)\leq 0,at all times, with a compact set \mathcal{U}\subset\mathbb{R}^{m} and a C^{1} constraint function h=(h_{1},\ldots,h_{n_{h}}). The input constraints encode, for instance, actuator limits, and the state constraints collision avoidance and workspace bounds over the entire body of a soft robot (Fig. 1). Real-time control relies on the reduced-order model (ROM)
\dot{z}=R(z,u),\qquad\widehat{x}=\Phi(z),with reduced state z\in\mathbb{R}^{r}, r\ll N, locally Lipschitz dynamics R, and a C^{2} lifting \Phi:\mathbb{R}^{r}\to\mathbb{R}^{N} that reconstructs the full state from z. Its image \Phi(\mathbb{R}^{r}), the states the ROM can represent, is the approximation manifold. A matrix P\in\mathbb{R}^{r\times N} extracts the reduced coordinates Px of the full state and is a left inverse of the lifting,
P\Phi(z)=z.This setup includes linear projections, \Phi(z)=Vz and P=W^{\top} with W^{\top}V=I, such as POD–Galerkin [12], Petrov–Galerkin [1], and modal truncation, as well as nonlinear liftings such as quadratic manifolds [13] and spectral submanifolds [14], with R projected from F or learned from data. Beyond the left inverse (4), the results below require neither invariance of the approximation manifold nor stability of the ROM, although both influence the growth of the error bounds (Remark 3).
Online, the controller measures only the reduced coordinates Px(t_{k}) at the update times t_{k}=k\delta, \delta>0, at which the MPC re-plans. We further assume a known bound \sigma_{0}\geq\|x(0)-\Phi(Px(0))\| on the initial distance to the approximation manifold (\sigma_{0}=0 if the system starts on it).
The problem is to design a controller that (i) plans with the ROM (3) without evaluating the FOM online, (ii) guarantees the constraints (2) for the closed-loop FOM under every admissible disturbance, and (iii) minimizes a task cost.
III Coupled Reduction-Error Certificates
Guaranteeing (2) for the FOM requires bounding the prediction error x-\widehat{x}=x-\Phi(z) of the ROM. We do so with radii \rho=(\rho_{\parallel},\rho_{\perp}) for its two components (5). Section III-A derives the exact error dynamics, and Section III-B bounds their growth and coupling rates and constructs the bounds offline.
III-A Tangential and Transverse Error Dynamics
Analogous to RN-ROMPC [15], we split the prediction error into two components via the nonlinear projection \Phi\circ P onto the approximation manifold \Phi(\mathbb{R}^{r}),
e_{\parallel}:=Px-z\in\mathbb{R}^{r},\qquad e_{\perp}:=x-\Phi(Px)\in\mathbb{R}^{N}.Here the tangential error e_{\parallel} measures prediction error in the reduced coordinates, whereas the transverse error e_{\perp} measures departure from the approximation manifold. By (4),
Pe_{\perp}=0,\qquad x=\Phi(z+e_{\parallel})+e_{\perp}.As Pe_{\perp}=0, the measurement Px only determines e_{\parallel} and leaves the high-dimensional e_{\perp} unmeasured.
The dynamics of the transverse error are given by
\dot{e}_{\perp}=\dot{x}-D\Phi(Px)P\dot{x}=Q(Px)\dot{x},where we define Q(\zeta):=I-D\Phi(\zeta)P. Differentiating (4) gives PD\Phi(\zeta)=I, so Q(\zeta)^{2}=Q(\zeta), \operatorname{range}Q(\zeta)=\ker P, and Q(\zeta)D\Phi(\zeta)=0, so Q(\zeta) is the oblique projection onto \ker P along the tangent space of the manifold. For a linear lifting \Phi(z)=Vz, Q=I-VP is constant.
Substituting (1), (3), and (6) into \dot{e}_{\parallel}=P\dot{x}-\dot{z} and \dot{e}_{\perp}=Q(Px)\dot{x} gives the exact coupled dynamics
\displaystyle\dot{e}_{\parallel} \\ \displaystyle=\underbrace{PF(\Phi(z+e_{\parallel})+e_{\perp},u)-R(z,u)}_{=:\,g_{\parallel}(z,u,e_{\parallel},e_{\perp})}+Pd, \\ \displaystyle\dot{e}_{\perp} \\ \displaystyle=\underbrace{Q(z+e_{\parallel})F(\Phi(z+e_{\parallel})+e_{\perp},u)}_{=:\,g_{\perp}(z,u,e_{\parallel},e_{\perp})}+Q(z+e_{\parallel})d.We call the disturbance-free part g:=(g_{\parallel},g_{\perp}) the error drift and its value at zero error the residual. The tangential residual g_{\parallel}(z,u,0,0)=PF(\Phi(z),u)-R(z,u) measures how closely R reproduces the on-manifold flow and vanishes when the ROM is projected as R(z,u)=PF(\Phi(z),u). The transverse residual g_{\perp}(z,u,0,0) is the rate at which the FOM flow leaves the manifold.
III-B Comparison System for the Error Norms
For radii \rho=(\rho_{\parallel},\rho_{\perp})\geq 0, the errors within the radii that satisfy Pe_{\perp}=0 form the compact convex set1
\mathcal{E}_{\rho}:=\{(e_{\parallel},e_{\perp}):Pe_{\perp}=0,\|e_{\parallel}\|\leq\rho_{\parallel},\|e_{\perp}\|\leq\rho_{\perp}\}.Since F is C^{1} and \Phi is C^{2}, the drift g is C^{1} in the errors, hence Lipschitz and one-sided Lipschitz on every compact convex set of errors [17, Lemma 1]. We compute these bounds on user-chosen compact sets \mathcal{E}_{\bar{\rho}} with a fixed \bar{\rho}=(\bar{\rho}_{\parallel},\bar{\rho}_{\perp})>0 and a compact region \mathcal{Z}\subset\mathbb{R}^{r} of reduced states containing those the task visits, trading domain of validity against conservatism. In particular, at each (z,u)\in\mathcal{Z}\times\mathcal{U}, there exist continuous one-sided Lipschitz bounds \mu_{\parallel}(z,u),\mu_{\perp}(z,u) and Lipschitz bounds \gamma_{\parallel\perp}(z,u),\gamma_{\perp\parallel}(z,u)\geq 0 such that, for every (e_{\parallel},e_{\perp})\in\mathcal{E}_{\bar{\rho}},
\displaystyle\langle e_{\parallel},\,g_{\parallel}(z,u,e_{\parallel},0)-g_{\parallel}(z,u,0,0)\rangle \\ \displaystyle\leq\mu_{\parallel}\|e_{\parallel}\|^{2}, \\ \displaystyle\|g_{\perp}(z,u,e_{\parallel},0)-g_{\perp}(z,u,0,0)\| \\ \displaystyle\leq\gamma_{\perp\parallel}\|e_{\parallel}\|, \\ \displaystyle\|g_{\parallel}(z,u,e_{\parallel},e_{\perp})-g_{\parallel}(z,u,e_{\parallel},0)\| \\ \displaystyle\leq\gamma_{\parallel\perp}\|e_{\perp}\|, \\ \displaystyle\langle e_{\perp},\,g_{\perp}(z,u,e_{\parallel},e_{\perp})-g_{\perp}(z,u,e_{\parallel},0)\rangle \\ \displaystyle\leq\mu_{\perp}\|e_{\perp}\|^{2}. \\ b_{\parallel}(z,u),b_{\perp}(z,u)\geq 0 \\ \|e_{\parallel}\|\leq\bar{\rho}_{\parallel} \\ \displaystyle b_{\parallel}(z,u) \\ \displaystyle\geq\|g_{\parallel}(z,u,0,0)\|+\|P\|\,\bar{d}, \\ \displaystyle b_{\perp}(z,u) \\ \displaystyle\geq\|g_{\perp}(z,u,0,0)\|+\|Q(z+e_{\parallel})\|\,\bar{d}.The four rates in (11a) assemble into the matrix
M(z,u):=\begin{bmatrix}\mu_{\parallel}(z,u)&\gamma_{\parallel\perp}(z,u)\\
\gamma_{\perp\parallel}(z,u)&\mu_{\perp}(z,u)\end{bmatrix}.The growth rates \mu_{\parallel},\mu_{\perp} on the diagonal are negative if the respective error dynamics contract. The transverse error perturbs the reduced prediction through \gamma_{\parallel\perp}, and the tangential error drives the state off the manifold through \gamma_{\perp\parallel}. The map b bounds the forcings that persist at zero error: the residual g(z,u,0,0) and the disturbance d.
For a planned ROM trajectory (z,u), we predict the radii \rho alongside the ROM with the comparison system
\dot{\rho}=M(z,u)\rho+b(z,u).Theorem 1 (Coupled reduction-error certificate (CoRE)).
Let M and b satisfy (11), let (z,u) be a planned ROM trajectory–input pair with values in \mathcal{Z}\times\mathcal{U} on [t_{0},t_{1}], and let x solve the FOM (1) on [t_{0},t_{1}] under the same input u and any disturbance with \|d(t)\|\leq\bar{d}, with the errors e_{\parallel},e_{\perp} defined by (5). Let \rho solve the comparison system (13) with \rho_{\parallel}(t_{0})\geq\|e_{\parallel}(t_{0})\| and \rho_{\perp}(t_{0})\geq\|e_{\perp}(t_{0})\|, and suppose that \rho<\bar{\rho} on [t_{0},t_{1}]. Then \|e_{\parallel}\|\leq\rho_{\parallel} and \|e_{\perp}\|\leq\rho_{\perp} on [t_{0},t_{1}].
Proof.
Throughout, we omit the time argument and the (z,u)-dependence. We first show that, for (e_{\parallel},e_{\perp})\in\mathcal{E}_{\bar{\rho}},
D^{+}\!\begin{bmatrix}\|e_{\parallel}\|\\
\|e_{\perp}\|\end{bmatrix}\leq M\begin{bmatrix}\|e_{\parallel}\|\\
\|e_{\perp}\|\end{bmatrix}+b,with D^{+} acting componentwise. For e_{\parallel}\neq 0, (8) and the upper Dini derivative of the Euclidean norm give
\displaystyle D^{+}\|e_{\parallel}\| \\ \displaystyle=\frac{e_{\parallel}^{\top}}{\|e_{\parallel}\|}\bigl[g_{\parallel}(e_{\parallel},e_{\perp})+Pd\bigr] \\ \displaystyle\leq\frac{e_{\parallel}^{\top}}{\|e_{\parallel}\|}\bigl[g_{\parallel}(e_{\parallel},0)-g_{\parallel}(0,0)\bigr]+\|g_{\parallel}(0,0)\| \\ \displaystyle\quad+\|g_{\parallel}(e_{\parallel},e_{\perp})-g_{\parallel}(e_{\parallel},0)\|+\|Pd\| \\ \displaystyle\leq\mu_{\parallel}\|e_{\parallel}\|+\gamma_{\parallel\perp}\|e_{\perp}\|+b_{\parallel},where the first inequality follows from telescoping g_{\parallel}(e_{\parallel},e_{\perp}) around (e_{\parallel},0) and (0,0) and from Cauchy–Schwarz (except the first term, which keeps its sign) and the second from (11) and \|d\|\leq\bar{d}. The case e_{\parallel}=0 follows from D^{+}\|e_{\parallel}\|\leq\|\dot{e}_{\parallel}\|. A similar argument applied to (9) gives the second row of (14). Since M is Metzler, the right-hand side of (13) is quasimonotone. Hence, the comparison theorem for systems of differential inequalities [18, Cor. 1.7.1] gives \|e_{\parallel}\|\leq\rho_{\parallel} and \|e_{\perp}\|\leq\rho_{\perp} as long as (e_{\parallel},e_{\perp})\in\mathcal{E}_{\bar{\rho}}, which Pe_{\perp}=0 from (6), these bounds, and \rho<\bar{\rho} ensure on [t_{0},t_{1}] by a continuation argument. ∎
Corollary 2 (Radius-dependent bounds).
Let M(z,u,\rho) and b(z,u,\rho) be continuous, locally Lipschitz and componentwise nondecreasing in \rho, and satisfy (11) on \mathcal{E}_{\rho} for every \rho\geq 0. Then Theorem 1 holds for \dot{\rho}=M(z,u,\rho)\rho+b(z,u,\rho) without the condition \rho<\bar{\rho}.
Proof.
Let v:=(\|e_{\parallel}\|,\|e_{\perp}\|). Since (e_{\parallel},e_{\perp})\in\mathcal{E}_{v}, the derivation of (14) with the bounds evaluated at \rho=v gives D^{+}v\leq G(t,v):=M(z,u,v)v+b(z,u,v) at every t. For 0\leq w\leq w^{\prime} with w_{i}=w_{i}^{\prime}, the monotonicity of M and b and \gamma\geq 0 give G_{i}(t,w)\leq G_{i}(t,w^{\prime}), so G is quasimonotone nondecreasing and locally Lipschitz in w, and the comparison theorem [18, Cor. 1.7.1] yields v\leq\rho on [t_{0},t_{1}]. ∎
Offline construction: Any offline bounding technique can be used to get (11), e.g., composing Lipschitz constants of the FOM nonlinearities [15]. Sharper values follow from the Jacobian of the drift with respect to the errors. Abbreviating \zeta:=z+e_{\parallel}, x_{e}:=\Phi(\zeta)+e_{\perp}, F:=F(x_{e},u), F_{x}:=D_{x}F(x_{e},u), Q:=Q(\zeta), and D\Phi:=D\Phi(\zeta), its blocks at any (z,u,e_{\parallel},e_{\perp}), ordered as in (12), are
\begin{bmatrix}D_{e_{\parallel}}g_{\parallel}&D_{e_{\perp}}g_{\parallel}\\
D_{e_{\parallel}}g_{\perp}&D_{e_{\perp}}g_{\perp}\end{bmatrix}=\begin{bmatrix}PF_{x}D\Phi&PF_{x}\\
QF_{x}D\Phi-D^{2}\Phi[PF]&QF_{x}\end{bmatrix},where D^{2}\Phi[v] denotes the derivative of D\Phi at \zeta in direction v. On the convex set \mathcal{E}_{\bar{\rho}}, one may take \gamma_{\perp\parallel},\gamma_{\parallel\perp} as the suprema of the norms of the off-diagonal blocks and \mu_{\parallel},\mu_{\perp} as the suprema of the logarithmic norms \mu_{2}(A):=\lambda_{\max}((A{+}A^{\top})/2) of the diagonal blocks [17, Lemma 1]. The blocks with respect to e_{\parallel} are evaluated at e_{\perp}=0, as in the first two inequalities of (11a). Those with respect to e_{\perp} act on \ker P only.
These suprema can be upper-bounded in many ways. We use three, one in each example of Section V. For polynomial dynamics with a linear lift, closed-form bounds in (z,u,\rho) hold on the current tube (Corollary 2, Section V-B). Otherwise, we fit upper envelopes to the block norms at samples over \mathcal{Z}\times\mathcal{U}\times\mathcal{E}_{\bar{\rho}} (Section V-A).2 When the FOM is too large to sample the block Jacobian, we directly fit M and b to prediction errors recorded along closed-loop trajectories [15]. The bounds (11) are then checked on the data rather than certified (Section V-C).
Remark 3 (Stability and structure of M).
While Theorem 1 applies to general matrices M constructed according to (11), the stability of the growth of the radii depends on whether M admits a Hurwitz constant bound. For M(z,u)\leq\bar{M} and b(z,u)\leq\bar{b} on \mathcal{Z}\times\mathcal{U}, \bar{M} is Hurwitz iff \bar{\mu}_{\parallel},\bar{\mu}_{\perp}<0 and \bar{\mu}_{\parallel}\bar{\mu}_{\perp}>\bar{\gamma}_{\parallel\perp}\bar{\gamma}_{\perp\parallel}. When \bar{M} is Hurwitz, \rho is ultimately bounded by -\bar{M}^{-1}\bar{b}. Negative \bar{\mu}_{\parallel},\bar{\mu}_{\perp} require stable error dynamics and a suitable norm or change of coordinates (footnote 1). As we see in the numerical example (Sec. V-B), CoRE-MPC can also be applied when M is not Hurwitz, since re-centering on the measurement can reset \rho_{\parallel} at every update. The carried transverse radius then stays bounded if, over one sampling interval, its own contraction outweighs the growth it receives through the coupling (thus, a negative \mu_{\perp} is necessary but not sufficient). For invariant manifolds [14], Q(\zeta)F(\Phi(\zeta),0)=0, the autonomous flow never leaves the manifold, so \gamma_{\perp\parallel}(z,0)=0 and M(z,0) is upper triangular. Thus, M(z,0) is Hurwitz iff \mu_{\parallel},\mu_{\perp}<0. This is the setting of RN-ROMPC [15]. In contrast, for general POD projections (Section V-B), the general structure of M proposed in CoRE-MPC is required.
IV CoRE-MPC: Predictive Control with Coupled Error Certificates
We now use the radii of Theorem 1 to tighten the FOM constraints inside a reduced-order MPC and establish its closed-loop guarantees (Theorem 5).
Constraint tightening: By (6), x=\Phi(z+e_{\parallel})+e_{\perp} with (e_{\parallel},e_{\perp})\in\mathcal{E}_{\rho} whenever the errors lie within the radii \rho. As in Section III-B, let c_{j}^{\parallel}(z),c_{j}^{\perp}(z)\geq 0 satisfy, for every z\in\mathcal{Z} and (e_{\parallel},e_{\perp})\in\mathcal{E}_{\bar{\rho}},
\displaystyle h_{j}(\Phi(z+e_{\parallel}))-h_{j}(\Phi(z)) \\ \displaystyle\leq c_{j}^{\parallel}(z)\|e_{\parallel}\|, \\ \displaystyle h_{j}(\Phi(z+e_{\parallel})+e_{\perp})-h_{j}(\Phi(z+e_{\parallel})) \\ \displaystyle\leq c_{j}^{\perp}(z)\|e_{\perp}\|,e.g., bounds over \mathcal{E}_{\bar{\rho}} on \|Dh_{j}\,D\Phi\| and on \|Dh_{j}\|_{\ker P}, the maximum of \|Dh_{j}\,v\| over unit v\in\ker P. Adding the two lines, (e_{\parallel},e_{\perp})\in\mathcal{E}_{\rho} implies h_{j}(x)\leq\widehat{h}_{j}(z,\rho) for the tightened constraint
\widehat{h}_{j}(z,\rho):=h_{j}(\Phi(z))+c_{j}^{\parallel}(z)\rho_{\parallel}+c_{j}^{\perp}(z)\rho_{\perp},which CoRE-MPC enforces in place of h_{j}\leq 0. For affine constraints h_{j}=a_{j}^{\top}x-\bar{h}_{j}, c_{j}^{\perp}=\|a_{j}^{\top}\|_{\ker P} is exact.
Let T\geq\delta be the prediction horizon. Given nonnegative stage and terminal costs \ell and V_{f} and writing \chi:=(z,\rho) for the augmented prediction state, CoRE-MPC solves at each t_{k} the following optimization problem:
\displaystyle V_{T}(Px(t_{k}),\sigma_{k}):=\min_{z,\rho,u}\int_{0}^{T}\!\ell(z(\tau),u(\tau))\,\mathrm{d}\tau+V_{f}(z(T)) \\ \displaystyle\text{s.t. }\ \dot{z}(\tau)=R(z(\tau),u(\tau)),\quad\forall\tau\in[0,T], \\ \displaystyle\hphantom{\text{s.t. }\ }\dot{\rho}(\tau)=M(z(\tau),u(\tau))\rho(\tau)+b(z(\tau),u(\tau)), \\ \displaystyle\hphantom{\text{s.t. }\ }\|Px(t_{k})-z(0)\|\leq\rho_{\parallel}(0),\quad\rho_{\perp}(0)=\sigma_{k}, \\ \displaystyle\hphantom{\text{s.t. }\ }\widehat{h}_{j}(z(\tau),\rho(\tau))\leq 0,\quad j=1,\ldots,n_{h}, \\ \displaystyle\hphantom{\text{s.t. }\ }(z(\tau),u(\tau))\in\mathcal{Z}\times\mathcal{U},\quad 0\leq\rho(\tau)<\bar{\rho}, \\ \displaystyle\hphantom{\text{s.t. }\ }\chi(T)\in\mathcal{X}_{f}.Here the initial condition (17d) leaves the nominal state z(0) free within the initial tube, as is standard in tube MPC. In particular, it admits any z(0) between the measurement Px(t_{k}), with \rho_{\parallel}(0)=0, and the shifted previous plan. The constraint (17f) keeps the plan inside the domain \mathcal{Z}\times\mathcal{U} of the bounds (11), the scalar \sigma_{k}\geq\|x(t_{k})-\Phi(Px(t_{k}))\| is a transverse-error bound maintained by Algorithm 1, and the terminal set \mathcal{X}_{f} is specified in Assumption 4. Let (z_{k}^{\star},\rho_{k}^{\star},u_{k}^{\star}) denote a minimizer of (17), assumed to exist. Algorithm 1 summarizes offline design and online operation.
Compared to a nominal ROM-MPC scheme, (17) adds the two scalar states \rho and replaces the state constraints by the precomputed tightened constraints \widehat{h}_{j}, so its computational demand is comparable to that of a nominal MPC for an (r{+}2)-dimensional system and the FOM vector field F is never evaluated online. Through M(z,u) and b(z,u), the tube grows where the local reduction geometry is unfavorable, so the optimizer trades performance against tube growth.
Closed-loop guarantees: The following assumption captures the required terminal ingredients for the augmented state \chi, whose dynamics we write as \dot{\chi}=f_{\rm C}(\chi,u):=(R(z,u),\,M(z,u)\rho+b(z,u)).
Assumption 4 (Terminal ingredients).
There exist a closed terminal set \mathcal{X}_{f}\subset\mathcal{Z}\times[0,\bar{\rho}) and a continuous terminal law \kappa_{f}:\mathcal{X}_{f}\to\mathcal{U} such that, for every \chi=(z,\rho)\in\mathcal{X}_{f},
\displaystyle\widehat{h}_{j}(z,\rho)\leq 0,\quad j=1,\ldots,n_{h}, \\ \displaystyle\mathcal{X}_{f}\text{ is positively invariant under }\dot{\chi}=f_{\rm C}(\chi,\kappa_{f}(\chi)).These are standard MPC terminal conditions (admissibility, viability) [20], posed for the augmented state \chi=(z,\rho) as in [5]. Viability now requires the terminal region to tolerate the growth of \rho, which stability of the ROM alone does not ensure, since nothing is reset within a plan. Under a constant Hurwitz bound \bar{M}, \bar{b} as in Remark 3, valid on \mathcal{X}_{f}, the box [0,\rho_{\rm f}] with \bar{M}\rho_{\rm f}+\bar{b}\leq 0, \rho_{\rm f}<\bar{\rho}, is positively invariant for (17c) and serves as the radius component of a product terminal set. Online, (17c) retains the state- and input-dependent M(z,u) and b(z,u).
Theorem 5 (CoRE-MPC closed-loop guarantees).
Let M and b satisfy (11), let c_{j}^{\parallel},c_{j}^{\perp} satisfy (15), let Assumption 4 hold, let \sigma_{0}\geq\|e_{\perp}(0)\|, and suppose that the CoRE-MPC problem (17) is feasible at t_{0}=0. Then, under Algorithm 1, (17) is feasible at every t_{k} and the closed-loop FOM satisfies
u(t)\in\mathcal{U},\qquad h(x(t))\leq 0,\qquad t\geq 0.Proof.
The argument is that of [15, Thm. 1], with the coupled comparison system (13) in place of two one-way coupled scalar tubes. Constraint satisfaction: If (17) is feasible at t_{k} with \sigma_{k}\geq\|e_{\perp}(t_{k})\|, then the initial-tube constraint (17d) gives (e_{\parallel}(t_{k}),e_{\perp}(t_{k}))\in\mathcal{E}_{\rho_{k}^{\star}(0)}, the domain constraint (17f) keeps the certificate valid, and Theorem 1 gives (e_{\parallel}(t_{k}+\tau),e_{\perp}(t_{k}+\tau))\in\mathcal{E}_{\rho_{k}^{\star}(\tau)} for \tau\in[0,\delta]. The tightened constraints (17e) then imply h(x)\leq 0 there by (16), and the input constraint in (17f) gives u\in\mathcal{U}. Recursive feasibility: The same containment at \tau=\delta is the initial condition (17d) for the shifted plan (\chi_{k}^{\star}(\cdot+\delta),u_{k}^{\star}(\cdot+\delta)) at t_{k+1}, and appending \kappa_{f} on [T-\delta,T] keeps \chi in \mathcal{X}_{f} by invariance (19). Induction over k completes the proof. ∎
If, in addition, \ell is positive definite in (z,u) and V_{f} is C^{1} with DV_{f}(z)\,R(z,\kappa_{f}(\chi))\leq-\ell(z,\kappa_{f}(\chi)) on \mathcal{X}_{f}, the value-decrease argument of [20] applies to the shifted candidate solution of the proof, so the planned nominal state z_{k}^{\star}(0) converges to the origin. The measured state Px(t_{k}) stays within \bar{\rho}_{\parallel} of it by the initial-tube constraint (17d).
48 cases relative to FOM-MPC; nominal ROM-MPC is excluded because it violates the constraint. (d) Solve times over 495 warm updates per controller from five cases run serially in balanced order on one thread. Boxes show quartiles, lines medians, whiskers the most extreme samples within 1.5 interquartile ranges, and circles outliers.V Experiments
We evaluate CoRE-MPC on three systems of increasing dimension:
- the oscillator chain (Section V-A) quantifies the trade-off between full-order safety, closed-loop cost, and computation on a reduction where the invariant-manifold assumption of RN-ROMPC [15] holds (Remark 3).
- the Burgers equation (Section V-B) tests the certificate on a POD reduction whose manifold is not invariant, with closed-form bounds that hold for every state and radius.
- the finite-element soft robot (Section V-C) tests whole-body constraint satisfaction, on all
1{,}628nodes of a9{,}768-state robot, with a two-state ROM and bounds fitted from data.
We compare the following controllers:
- CoRE-MPC: the proposed scheme (17).
- Constant CoRE: CoRE-MPC with
Mandbreplaced by constants, their domain maxima on the oscillator chain and separately fitted constants on the soft robot. - Nominal ROM-MPC: the same ROM, horizon, and cost function without accounting for model mismatch.
- FOM-MPC [5]: robust MPC on the FOM with a single tube propagated by a sampled logarithmic-norm bound and re-centered on the measured full state at every update. As it uses the full model and full state, it serves as the benchmark for closed-loop cost and computation time.
- RN-ROMPC [15]: the special case of (13) with constant
Mand\gamma_{\perp\parallel}=0by invariance of the manifold (Remark 3), and an input-dependent forcingb(u)=b_{0}+b_{u}\|u\|with constant gains, which absorbs the input-induced departure from the manifold. Its bounds are sampled from the model on the oscillator chain and fitted from rollouts on the soft robot. - Linear ROMPC [7, 8]: MPC on the linearized ROM with the nonlinearity bounded as a disturbance.
Within each example, all controllers share the task, horizon, and cost function, and use a weighted Euclidean norm (footnote 1). FOM-MPC measures the full state, whereas the reduced-order controllers measure only Px. All controllers use a quadratic terminal cost. On the oscillator chain, CoRE-MPC also imposes the terminal set (17g) (Section V-A), so Theorem 5 applies. The baselines and the remaining examples omit it for ease of implementation, which relaxes recursive feasibility and so does not disadvantage the baselines. Constraint satisfaction on each sampling interval still follows from Theorem 1 whenever the problem is feasible and the bounds hold. Solve times are measured using the SQP solver in acados [21] with HPIPM on an Apple M3 Max. All controllers shown fully converge at every control update, so no fallback to the shifted previous plan of the proof of Theorem 5 was needed. Implementation details and certificate parameters are documented with the accompanying code.
t^{*}=2.2 s with the reference (dotted), the limit (dashed), and the actuator locations (triangles). (b) Largest normalized field value over all nodes. (c) CoRE-MPC radii (solid) and realized errors (dashed) over the first second; the tangential radius, chosen by the optimizer, restarts at zero after the first update, whereas the transverse radius is carried over. (d) Solve times over 316 updates per controller from four cases, with box conventions as in Fig. 2.V-A Oscillator Chain
Setup: We consider the oscillator chain from [22] with 40 masses connected by nonlinear springs and linear dampers, resulting in N=80 states (Fig. 2(a)). A single actuator with the spatial shape of the first linear mode drives the chain, and a bounded force disturbance acts on the last mass. Every displacement is constrained, |q_{i}|\leq q_{\max}, and the reference requests a first-mode peak one third above this limit, so the constraint is active. The ROM has r=2 states on a spectral submanifold [14], which is invariant for the autonomous FOM, so the assumption of RN-ROMPC holds. The horizon is 2 s and the sampling time 100 ms. The bounds of CoRE-MPC and RN-ROMPC are obtained from samples on the same domain. Those of CoRE-MPC include the covering margins of footnote 2, so (11) holds on all of \mathcal{Z}\times\mathcal{U}\times\mathcal{E}_{\bar{\rho}}. CoRE-MPC further imposes the terminal set (17g). The augmented terminal region uses a continuous-time, two-radius adaptation of the construction in [23, Prop. 6]. Its numerical definition and its admissibility and invariance certificate, (18) and (19), are provided with the accompanying code. We test 48 initial conditions for 10 s each, 16 of which were generated to stress the constraint handling of the controllers.
Results: Fig. 2 shows the case with the largest nominal violation and the statistics over all cases. Nominal ROM-MPC violates the displacement limit in 23 of 48 cases, by up to 4\% of q_{\max} (Fig. 2(b)). In contrast, all robust controllers satisfy the constraints in all cases, and the realized errors remain inside the predicted tubes throughout, in accordance with Theorem 1. The excess closed-loop cost over FOM-MPC (Fig. 2(c)) increases with the conservatism of the bounds, from 35\% for CoRE-MPC to 49\% for Constant CoRE and 149\% for RN-ROMPC. Constant CoRE and RN-ROMPC replace the state- and input-dependent bounds by constants, so their tubes grow at the worst-case rate along the entire plan. Consequently, RN-ROMPC requires a large conservative backoff from the constraints, while CoRE-MPC can operate close to them (Fig. 2(b)). Linear ROMPC is infeasible in all 48 cases, since bounding the nonlinear spring force as a disturbance is too conservative. In computation (Fig. 2(d)), CoRE-MPC solves 61 times faster than FOM-MPC at 2.9 times the solve time of nominal ROM-MPC, with a 95th percentile of 10.6 ms against the 100 ms sampling time. Constant CoRE and RN-ROMPC, which also propagate two radii but omit the terminal set and use simpler bounds, solve in about 80\% of the time of CoRE-MPC. Thus, even where the assumption of RN-ROMPC holds, CoRE-MPC incurs less than a quarter of its excess closed-loop cost at comparable computation.
0.3 to 3.3 s). (c) Minimum clearance over all 1{,}628 nodes; red marks wall penetration. (d) Controller wall times over 190 updates, including the check of all nodes, with box conventions as in Fig. 2; times exclude FEM integration.V-B Burgers Equation
Setup: We consider the viscous Burgers equation on the unit interval with viscosity \nu=0.06, discretized on 128 nodes and driven by three Gaussian actuators (Fig. 3(a)). A six-mode POD–Galerkin ROM [24] captures 99.998\% of the snapshot energy, and its approximation manifold is not invariant. The magnitude of the field is bounded by w_{\max} at every node, while the reference peaks 54\% above this limit. The horizon is 0.3 s and the sampling time 100 ms. A bounded, spatially localized sinusoidal disturbance acts on the field, and its 24 equally spaced phases define the 24 cases, each run for 8 s from rest. The FOM dynamics are quadratic in the state and the lift is linear, so the Jacobian blocks and residuals are polynomial in (z,u,e) and the norm inequalities give closed-form upper bounds satisfying (11) on \mathcal{E}_{\rho} for every \rho\geq 0 instead of on a fixed \mathcal{E}_{\bar{\rho}}, affine in \rho with nonnegative coefficients. The comparison system (13) becomes \dot{\rho}=M(z,u,\rho)\rho+b(z,u,\rho), to which Theorem 1 extends by Corollary 2.
Results: Fig. 3 shows the case with the largest nominal violation and the solve times. Despite the small discarded snapshot energy, nominal ROM-MPC violates the limit in all 24 cases, by up to 0.9\% of w_{\max} and for most of the run (Fig. 3(b)). CoRE-MPC and FOM-MPC complete all cases without violation and converge to a similar steady state at about 0.90\,w_{\max}. Along the reference, both couplings and the transverse residual are nonzero, so the coupled structure of (13) is required, and M is not Hurwitz. CoRE-MPC nevertheless remains feasible in every case: with the free initialization (17d), the optimizer selects \rho_{\parallel}(0)=0 at every update after the first, so the tangential radius restarts, and the carried transverse radius stays bounded (the gain from \sigma_{k} to \sigma_{k+1} is at most 0.73 over all intervals, cf. Remark 3). The certificate is valid and tight: every tube contains the realized errors, which reach 63\% of the transverse and 98\% of the tangential radius within the applied intervals over all cases (Fig. 3(c)), and along the evaluated trajectories the coupling bounds exceed the exact Jacobian block norms by factors of at most 1.15 and 1.38. CoRE-MPC incurs 27\% higher closed-loop cost than FOM-MPC and solves 13 times faster, at 1.5 times the solve time of nominal ROM-MPC, with a 95th percentile of 12.3 ms against the 100 ms sampling time (Fig. 3(d)).
V-C Finite-Element Soft Robot
Setup: We consider the cable-driven diamond robot of [2], simulated in SOFA [25] with 9{,}768 states, 1{,}628 nodes, and four cable forces between 0 and 3 N. The task is to track an elliptical tip reference that extends beyond a wall (Fig. 4(a)), so it cannot be tracked exactly on its upper arc, while every node must stay clear of the wall. The ROM has r=2 states on a spectral submanifold of the FEM [14], whose lifting reconstructs all nodes. The horizon is 60 ms and the sampling time 20 ms. Each controller runs the task for 3.8 s without disturbance. Since the FOM is too large to sample the block Jacobian, an input-dependent M and a state- and input-dependent b are fitted to the prediction errors of ten earlier rollouts on the same robot and wall with nearby elliptical references (Section III-B), and the bounds of RN-ROMPC and Constant CoRE are fitted to the same data, so all three are empirical. The weighting of the error norm is optimized on the residuals of the calibration rollouts to minimize the wall tightening over all nodes (footnote 1). Only the nodes with the largest predicted wall envelope enter the optimization, and every plan is checked against all 1{,}628 nodes, which CoRE-MPC never violates.
Results: The clearance in Fig. 4(c) summarizes the outcome: nominal ROM-MPC penetrates the wall by 13 mm, whereas the three robust controllers keep every node clear, with minimum clearances of 3.7 mm for CoRE-MPC, 5.4 mm for Constant CoRE, and 0.5 mm for RN-ROMPC. RN-ROMPC [15], with its constant upper-triangular M, backs away from the wall along the upper arc yet comes closest to it. Averaged over matched prediction stages with identical initial radii, the wall tightening of RN-ROMPC is 2.3 times that of CoRE-MPC (15.3 against 6.6 mm). CoRE-MPC, whose M and b follow the plan, tracks with 29\% less error than Constant CoRE and 38\% less than RN-ROMPC (Fig. 4(b)), and the lower error of nominal ROM-MPC is obtained on a trajectory that penetrates the wall. In computation (Fig. 4(d)), all four controllers solve in under 2 ms at the median, CoRE-MPC in 0.9 ms against 0.55 ms for nominal ROM-MPC and never above 5 ms, whereas two RN-ROMPC updates exceed the 20 ms period. Full-order MPC is excluded, since one FEM horizon prediction alone takes 0.56 s, or 28 sampling periods.
Validation: To test the fitted bounds beyond the calibration references, we fit M and b of the three robust controllers above to the prediction errors of 202 rollouts on 24 elliptical references with a 10\% margin over the recorded errors. We then apply CoRE-MPC and RN-ROMPC on five held-out elliptical references, one of which is repeated with an initial transverse error of 0.5 mm and again, as a stress test, with a 50 mN force toward the wall exceeding \bar{d}. In all 14 held-out runs, every node stays clear of the wall and both radii contain the realized errors at every update, reaching at most 90\% of the tangential and 77\% of the transverse radius. On the task of Fig. 4, the larger calibration set increases the RMS tracking error of all controllers. CoRE-MPC increases to 17.4 mm, which is still 20\% less than Constant CoRE (21.7 mm) and 50\% less than RN-ROMPC (34.7 mm).
VI Conclusion
We have presented CoRE-MPC, a reduced-order MPC that enforces full-order constraints for a general class of nonlinear ROMs, including non-invariant approximation manifolds and ROMs learned from data. The error within the reduced coordinates can be reset at every update, but the error outside them cannot, so CoRE-MPC bounds the two separately and propagates both radii with the ROM at the cost of two extra states. In the numerical studies, this yields a fraction of the excess cost of existing certified schemes where they apply, and real-time whole-body collision avoidance for a soft robot where full-order robust MPC is prohibitive. Open issues include the certification of comparably tight bounds at the finite-element scale, where the bounds are currently fitted to data, and additional sensing of the transverse error, which is propagated open-loop.
Acknowledgment
GPT-6 Astra (OpenAI) was used to assist in writing the scripts for Figs. 1–4, and Claude Fable 5.1 (Anthropic) was used for editing and grammar. The authors independently verified every derivation, implementation, and result.