<figure><img src="image_1.png" /></figure>

HAL
open science

# CPC: Cascaded Predictive Control

Nicola Scianca, Enrico Mingo Hoffman

► To cite this version:

Nicola Scianca, Enrico Mingo Hoffman. CPC: Cascaded Predictive Control. 2025 IEEE-RAS 24th International Conference on Humanoid Robots, Sep 2025, Seoul, South Korea. hal-05070222

HAL Id: hal-05070222

https://hal.science/hal-05070222v1

Submitted on 16 May 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
CPC: Cascaded Predictive Control

Nicola Scianca<fnref n="1" /> and Enrico Mingo Hoffman<fnref n="2" />

**Abstract**— We present *Cascaded Predictive Control* (CPC), a framework that decomposes a model-predictive control horizon into multiple stages, each solved with the model and method best suited to its accuracy requirements. CPC transfers information across stages through the value function: the terminal cost produced by a downstream stage acts as a boundary condition for the preceding one, enabling heterogeneous combinations of different techniques and solvers. Two proof-of-concept studies illustrate the approach. First, a humanoid walking controller uses a double-DDP cascade that schedules single-rigid-body and linear-inverted-pendulum models to extend the prediction horizon at a lower computational cost. Second, a cart-pendulum swing-up couples a DDP stage with a value function learned by an actor-critic policy; the hybrid controller succeeds where the policy alone fails. Although the examples are preliminary, they indicate that CPC can provide a flexible framework for achieving effective results with reduced computational effort by combining the strengths of different techniques.

## I. INTRODUCTION

Recent successes in agile control of legged robots are, for the most part, obtained using either Model Predictive Control (MPC) [1] or Reinforcement Learning (RL) methods [2]. Whereas RL excels at contact-agnostic motion (e.g., walking on rough terrains) and imitation of human movements. These policies are usually not easy to repurpose and do not offer any guarantees beyond what can be ascertained with extensive testing campaigns. MPC, on the other hand, offers more general frameworks but is often limited by the quality of the contact models, and long horizons are still computationally prohibitive.

An effective and efficient technique for online MPC is Differential Dynamic Programming (DDP) [3], [4], which exploits the recursive structure of the control problem to offer a solution that is linear in the length of the horizon. In fact, most modern techniques for optimization-based robot control are based on DDP-like techniques [5], [6], [7], [8]. However, despite the efficient formulations and implementations, accurate models still do not allow for long prediction horizons [9]. To reduce the computational load, previous research has proposed to use different models throughout the prediction horizon [8].

RL-based approaches do not suffer from problems related to horizon length, as all computations are done when training

the network. However, a policy capable of effective bipedal locomotion is not easily usable for related tasks, such as moving across obstacles or stepping on stairs, and usually requires extensive retraining. Furthermore, there is no way to ascertain that the policy is effective across the entire state space, if not through extensive testing.

These complementary strengths and weaknesses have sparked research on hybrid optimization–learning control. Examples include [10], which embeds a differentiable MPC layer inside an actor–critic agent, and [11], where the authors propose to use a trajectory optimization algorithm in place of the critic network during training of a DDPG policy. In [12], MPC is used to enforce constraints on the output of a learned policy. [13] uses MPC in parallel with a policy to achieve effective footstep placement. [14] adds compensations to an MPC controller using a learned policy.

In this paper, we propose the idea of Cascaded Predictive Control (CPC). In our formulation, the horizon is constituted by a sequence of stages, each optimizing a suitable representation of the robot trajectory. While this concept has been previously explored in a less general fashion [8], we extend the idea not just to using a different model of the robot, but entirely different optimization techniques, which we demonstrate by combining DDP with an RL actor-critic method.

## II. GENERAL FORMULATION

Our approach works by solving a cascade of problems in series over a prediction horizon. Each problem is referred to as *stage*, and is identified by an upper index. For example, the variables in stage *j* will be

$$ 
\begin{aligned}
\boldsymbol{u}^j &= (\boldsymbol{u}_0^j, \dots, \boldsymbol{u}_{N_j-1}^j), \\
\boldsymbol{x}^j &= (\boldsymbol{x}_0^j, \dots, \boldsymbol{x}_{N_j-1}^j, \boldsymbol{x}_{N_j}^j).
\end{aligned}
\qquad (1) 
$$

i.e., $N_j$ input variables and $N_j + 1$ state variables Note that input and state symbols without a lower index, e.g., $\boldsymbol{x}^j$ or $\boldsymbol{u}^j$, represent all variables in stage *j*, whereas $\boldsymbol{x}$ and $\boldsymbol{u}$ represent the total sets of variables in all stages. The problem is constituted by *M* stages plus a terminal stage. The latter always has a length of 0, and is thus just constitute by the terminal state of the overall trajectory $\boldsymbol{x}_0^M$.

These sequences must satisfy the discrete-time dynamics of the prediction model used during stage *j*, i.e.,

$$ \boldsymbol{x}_{k+1}^j = f(\boldsymbol{x}_k^j, \boldsymbol{u}_k^j). \qquad (2) $$

Furthermore, the terminal state of each stage must be consistent with the initial state of the next stage, i.e.,

$$ g(\boldsymbol{x}_{N_j}^j) = \boldsymbol{x}_0^{j+1}, \qquad (3) $$

<footnote n="1">¹Dipartimento di Ingegneria Informatica, Automatica e Gestionale, Sapienza Università di Roma, Italy; E-mail: scianca@diag.uniromal.it</footnote>

<footnote n="2">²Université de Lorraine, CNRS, Inria LARSEN Team, LORIA, France. E-mail: enrico.mingo-hoffman@inria.fr</footnote>

<footnote>This work was granted by the French National Research Agency (ANR) under the project ANR-24-CE33-0753-01 (MeRLin); Nicola Scianca has been fully supported by PNRR MUR project PE0000013-FAIR.</footnote>
where $g(x)$ is a function that maps the state between consecutive stages<fnref n="1" />.

The cost function to be minimized can be formulated as a sum of the running costs across all timesteps in all stages, plus a terminal cost. The overall optimization problem is then formulated as

$$ \min_{\mathbf{x},\mathbf{u}} \sum_{j=0}^{M-1} \sum_{i=0}^{N_j-1} \ell_i^j(\mathbf{x}_i^j, \mathbf{u}_i^j) + \ell_f(\mathbf{x}_0^M) $$
subject to:
$$ \begin{aligned} \mathbf{x}_{k+1}^j &= f(\mathbf{x}_k^j, \mathbf{u}_k^j) \\ \mathbf{x}_0^0 &= \bar{\mathbf{x}}_{\text{curr}} \\ \mathbf{x}_0^{j+1} &= g(\mathbf{x}_{N_j}^j). \end{aligned} \quad (4) $$

Note that the individual stages do not include a terminal cost, as there is only an overall terminal cost given by the final stage.

The value function at the beginning of the k-th stage can be written as

$$ V(\mathbf{x}_0^k) = \min_{\mathbf{u}^k} \sum_{j=k}^{N_j-1} \ell_i^j(\mathbf{x}_i^j, \mathbf{u}_i^j) + V(\mathbf{x}_0^{k+1}). \quad (5) $$

This is essentially the Bellman equation applied to the entire stage, and shows that given knowledge of the value function at stage k+1 we can solve the stage k as an independent optimization problem, and possibly use a different optimization technique in each stage.

While there are several interesting combinations of techniques that could be explored, in this preliminary study, we chose to focus on two in particular: DDP and Actor-Critic reinforcement learning.

## A. DDP

DDP/iLQR is an optimization technique<fnref n="2" /> tailored to optimal control problems. The basic formulation is unconstrained, except for the constraints on dynamics and initial state, but Augmented Lagrangian approaches allow for more general constraint integration [16]. Here we describe a basic unconstrained formulation. For ease of notation, we drop the j index relative to the stage, as all computations here are relative to the variables within the generic j-th stage.

We recall Bellman's optimality principle for the finite-horizon problem

$$ \begin{aligned} V_k(\mathbf{x}_k) &= \min_{\mathbf{u}_k} \left[ \ell_k(\mathbf{x}_k, \mathbf{u}_k) + V'(f(\mathbf{x}_k, \mathbf{u}_k)) \right] \\ &= \min_{\mathbf{u}_k} \left[ Q(\mathbf{x}_k, \mathbf{u}_k) \right] \end{aligned} \quad (6) $$

where $V'$ represents the value function at $k+1$. DDP approximates the next state value function $V'(f(\mathbf{x}, \mathbf{u}_k))$

using a second-order Taylor expansion around the current guess $\bar{\mathbf{x}}_{k+1}$ of the next state, i.e.,

$$ V'(\bar{\mathbf{x}}_{k+1} + \delta \mathbf{x}_{k+1}) \approx \bar{V}' + V_x' \delta \mathbf{x}_{k+1} + \frac{1}{2} \delta \mathbf{x}_{k+1}^\top V_{xx}' \delta \mathbf{x}_{k+1}, $$

where $\bar{V}', V_x'$ and $V_{xx}'$ represents, respectively, the zero-order term, gradient, and Hessian of the currently available approximation of the value function of the next state $\mathbf{x}_{k+1}$.

In basic DDP, the initial state guess is propagated forward using the dynamics. We instead use a modification to this approach that computes a dynamic defect [5]

$$ d_k = \mathbf{x}_{k+1} - f(\mathbf{x}_k, \mathbf{u}_k). \quad (7) $$

In DDP, we approximate the Q-function $Q(\mathbf{x}_k, \mathbf{u}_k)$ using a second-order Taylor expansion

$$ \bar{Q}_k + \begin{bmatrix} Q_{x,k} \\ Q_{u,k} \end{bmatrix}^T \begin{bmatrix} \delta \mathbf{x}_k \\ \delta \mathbf{u}_k \end{bmatrix} \frac{1}{2} + \begin{bmatrix} \delta \mathbf{x}_k \\ \delta \mathbf{u}_k \end{bmatrix}^T \begin{bmatrix} Q_{xx,k} & Q_{xu,k} \\ Q_{ux,k} & Q_{uu,k} \end{bmatrix} \begin{bmatrix} \delta \mathbf{x}_k \\ \delta \mathbf{u}_k \end{bmatrix}, $$

where $\bar{Q}_k$ is the zero-order term, and the first and second order terms can be identified by the appropriate subscripts. These are found to be

$$ Q_{x,k} = \ell_{x,k} + A_k^T (V_x' + V_{xx}' d_{k+1}) $$

$$ Q_{u,k} = \ell_{u,k} + B_k^T (V_x' + V_{xx}' d_{k+1}) $$

$$ Q_{xx,k} = \ell_{uu,k} + A_k^T V_{xx}' A_k \quad (8) $$

$$ Q_{uu,k} = \ell_{uu,k} + B_k^T V_{xx}' B_k $$

$$ Q_{ux,k} = \ell_{ux,k} + B_k^T V_{xx}' A_k, $$

with

$$ A_k = \left. \frac{\partial f}{\partial \mathbf{x}} \right|_{\bar{\mathbf{x}}_k}, \quad B_k = \left. \frac{\partial f}{\partial \mathbf{u}} \right|_{\bar{\mathbf{u}}_k}, \quad (9) $$

Solving the minimization problem in (6) yields the well-known local Riccati feedback law

$$ \delta \mathbf{u}_k^* = \mathbf{k}_k + \mathbf{K}_k \delta \mathbf{x}_k \quad (10) $$

with

$$ \mathbf{k}_k = -Q_{uu,k}^{-1} Q_{u,k} \quad (11) $$

$$ \mathbf{K}_k = -Q_{uu,k}^{-1} Q_{ux,k}. $$

The value function approximation is then updated as

$$ V_x = Q_{x,k} - \mathbf{K}_k^T Q_{u,k} \mathbf{k}_k $$

$$ V_{xx} = Q_{xx,k} - \mathbf{K}_k^T Q_{ux,k} \mathbf{K}_k $$

Given the gains $(\mathbf{k}_k, \mathbf{K}_k)$ we compute a new guess for the trajectory as

$$ \mathbf{u}_k = \bar{\mathbf{u}}_k + \alpha \mathbf{k}_k + \mathbf{K}_k (\mathbf{x}_k - \bar{\mathbf{x}}_k) \quad (12) $$

$$ \mathbf{x}_{k+1} = f(\mathbf{x}_k, \mathbf{u}_k) - (1 - \alpha) d_{k+1} $$

with a line-search parameter $\alpha \in (0, 1]$. Then, the defects $d$ are correspondingly updated.

<footnote n="1">¹For example, if stage j uses joint variables to represent the state and stage j+1 uses an end-effector pose, the mapping function will be the forward kinematics of said end-effector.</footnote>

<footnote n="2">²Technically, DDP and iLQR refer to the same algorithm, but iLQR neglects second-order dynamic derivatives whereas DDP does not. However, since this is a very common approximation, the terms are mostly used interchangeably. For this reason, we prefer to use the term DDP. The interested reader can refer to [15] for more information.</footnote>
## B. Actor-Critic RL stage

We learn an Actor policy $\pi_{\theta}^{j}(u|x)$ that minimizes<fnref n="3" /> the expected discounted trajectory cost

$$ \mathbb{E}_{\pi_{\theta}^{j}} \left[ \sum_{k=0}^{N_j-1} \gamma^k \ell_k^j \right], \quad \gamma \in (0, 1), $$

and a Critic $V_{\phi}^{j}(x)$ that approximates the value function. Both the Actor and the Critic are neural networks with a stacked parameter vector $\theta$. We train the networks using Proximal Policy Optimization (PPO). The Actor loss is

$$ -\frac{1}{K} \sum_{k=0}^{K-1} \min \left( \rho_k \hat{A}_k, \text{clip} \left( \rho_k, 1-\varepsilon, 1+\varepsilon \right) \hat{A}_k \right). $$

where $K$ is the batch size, $\rho_k = \pi_{\theta}(u_k|x_k)/\pi_{\theta_{\text{old}}}(u_k|x_k)$ is the likelihood ratio, i.e., the ratio between the probability of taking the action under new policy over the old policy, and $\epsilon$ is the clipping parameter. The advantage $\hat{A}_k$ is computed with Generalized Advantage Estimation

$$ \hat{A}_k = \sum_{l=0}^{N_j-k-1} (\gamma \lambda)^l (\ell_{k+l}^j + \gamma V_{\phi}^j(x_{k+l}^j) - V_{\phi}^j(x_{k+l}^j)), $$

with trace parameter $\lambda \in [0, 1]$. The critic is trained according to the loss

$$ \frac{1}{K} \sum_{k=0}^{K-1} \text{SmoothL1}(V_{\phi}(x_k) - \hat{A}_k - V_{\phi}(x_k)), $$

where SmoothL1($z$) is the smooth-$L_1$ penalty.

In this work we only use this as a terminal stage, and those only use the critic network to provide the initial value function approximation. This is computed using automatic differentiation on the critic network itself. Note that, in order to be able to compute the Hessian, activation functions must be smooth, thus ReLU is not a viable option.

## C. Connecting consecutive stages

The terminal state of stage $j$, $x_{N_j}^j$, is mapped to the initial state of the next stage through the map (3), which is enforced during the forward pass.

During the backward pass, we must project the value function approximation to the new state. Let

$$ G_j = \frac{\partial g_j}{\partial x_{N_j}^j}, \qquad (13) $$

the Jacobian of the map $g$ w.r.t. the state in stage $j$. With the quadratic value expansion at the start of stage $j+1$,

$$ V^{j+1}(x) = \bar{V}^{j+1} + V_x^{j+1} \Delta x + \frac{1}{2} \Delta x^T V_{xx}^{j+1} \Delta x, $$

the derivatives at the end of stage $j$ become

$$ V_x^j = G_j^T V_x^{j+1}, \qquad (14) $$

$$ V_{xx}^j = G_j^T V_{xx}^{j+1} G_j. \qquad (15) $$

Equation (15) serves as the boundary condition for the backward sweep of stage $j$.

## Algorithm 1: CPC algorithm

1. prepare initial guess ($x, u$)

2. initialize value function approximation $V_x, V_{xx}$ of the terminal stage

3. repeat

4. for $j \in [M-1, \dots, 0]$ do

5. project value function approximation to the state of the $j$-th stage;

6. compute backward pass;

7. update value function approximation $V_x, V_{xx}$

8. for $j \in [0, \dots, M-1]$ do

9. compute forward pass;

10. compute initial state for next stage

11. $x_0^{j+1} = g(x_{N_j}^j)$;

11. until convergence or max iterations;

## D. Cascaded Predictive Control Algorithm

Algorithm 1 shows the high-level operations that must be executed in our CPC architecture. The first operation is to initialize the value function approximation for the terminal stage. In our examples, this is either done using a terminal cost term or using a learned value function approximation.

The value function approximation is then passed backwards to stage $M-1$ and becomes the value function approximation at the end of this stage. If the state is different, it is mapped to the new state using (15). Then, a backward pass is executed, giving us a value function approximation at the beginning of stage $M-1$. This is passed to the preceding stage, which performs its own backward pass, until we get to stage 0.

Now the algorithm performs a forward pass starting from the initial (measured) state $\bar{x}_{\text{curr}}$. This means propagating the dynamics forward using the newly computed optimal inputs, possibly including a line search for determining an optimal step size. After stage 0 has finished its forward pass, its terminal stage is used to determine the initial state of the next stage using (3). This is repeated until we reach the final stage, at which point we have a new guess for the trajectory across all stages.

This process is repeated until we reach either a convergence condition or a maximum number of iterations. In our case, we used the latter, as this is much more common in MPC, to contain the computation time to a predictable value.

## III. EXAMPLE 1: DOUBLE DDP FOR HUMANOID WALKING

Our first example focuses on humanoid locomotion. We use a three stage optimization (two prediction stages plus the final stage) that employ different approximations of the humanoid dynamics. The centroidal dynamics of a humanoid can be approximated using several simplified models, such as:

* the Single Rigid-Body Dynamics (SRBD) approximation, in which the robot is approximated as having a

<footnote n="3">³We use the cost $\ell$ in place of a reward for consistency with the optimization problem, so signs are flipped with respect to conventional RL language.</footnote>
<figure><img src="image_3.png" /><figcaption>Fig. 1: Two-stage MPC using SRBD as the first stage dynamic model and a LIP as the second stage dynamic model. The predicted contact forces at the beginning of each stage are displayed, along with CoM and feet trajectories. The ellipsoid visualizes the rigid body at select predicted nodes, and the sphere is centered around the CoM during LIP predicted nodes.</figcaption></figure>

fixed inertia matrix, simplifying the computation of the angular momentum;

* the Linear Inverted Pendulum (LIP) approximation, in which rotational components of the dynamics are neglected, and contact forces are encoded using the Zero Moment Point (ZMP).

The first model is more accurate, but since it is nonlinear and has more states and inputs, it demands more computational resources. The second model is considerably simpler as it is linear and has fewer variables, but completely neglects angular contributions. The use of the two models in conjunction is aimed at combining the strengths of each model to obtain the best of both worlds. The more complex model is used first in the horizon, where higher fidelity is more important, and the simpler model is used last, as approximations can be well tolerated in the later stages of the prediction. The final stage is constituted by the final cost terms on the LIP.

The first stage uses the SRBD, in which the state and input vectors are constituted by

$$x^1 = \begin{pmatrix} c \\ r \\ p \\ \dot{c} \\ \omega \\ \dot{p} \end{pmatrix}, \quad u^1 = \begin{pmatrix} f_i \\ \vdots \\ f_{n_c} \\ \ddot{p} \end{pmatrix}, \qquad (16)$$

where $c \in \mathbb{R}^3$ is the position of the Center of Mass (CoM), $r \in SO(3)$ is the orientation of the rigid body represented as a quaternion, $\omega \in \mathbb{R}^3$ is its angular velocity, and $p \in \mathbb{R}^{3 \times n_c}$ is the stack of all the $n_c = 4$ contact positions$^4$, and $f_i$, with $i = 0, \dots, n_c - 1$ are the contact forces. The dynamics in

the first stage are

$$\dot{x}^1 = \frac{d}{dt} \begin{pmatrix} c \\ r \\ p \\ \dot{c} \\ \omega \\ \dot{p} \end{pmatrix} = \begin{pmatrix} \dot{c} \\ \frac{1}{2} r \otimes \begin{bmatrix} 0 \\ \omega \end{bmatrix} \\ \dot{p} \\ \frac{1}{m} \sum_{i=0}^{n_c-1} f_i + g \\ I^{-1} \left( \sum_{i=0}^{n_c-1} (p_i - c) \times f_i + \omega \times (I\omega) \right) \\ \ddot{p} \end{pmatrix},$$

where $m$ is its mass of the rigid body, $g = (0, 0, -9, 81)$ m/s$^2$ is the gravity acceleration, $I$ is the inertia matrix, and $\otimes$ denoting the quaternion multiplication. The discrete model is obtained with Euler integration:

$$x_{k+1}^1 = f_{\text{SRB}}(x_k^1, u_k^1). \qquad (17)$$

The second stage uses a 3D LIP approximation of the dynamics, in which rotational contributions to the dynamics are neglected. Different from the classical LIP approximation, in which the CoM height is kept constant, we use a 3D extension that still results in a linear model but allows vertical CoM motion [17].

The state and input vector in the second stage are constituted by

$$x^2 = \begin{pmatrix} c \\ p \\ \dot{c} \\ \dot{p} \end{pmatrix}, \quad u^2 = \begin{pmatrix} z \\ \ddot{p} \end{pmatrix}. \qquad (18)$$

Note that the state vector is a subset of the state vector used in the previous stage, whereas the input vector is changed, as now we use the position of the 3D ZMP $z \in \mathbb{R}^3$ in place of the contact forces. The dynamics in the second stage are

$$\dot{x}^2 = \frac{d}{dt} \begin{pmatrix} c \\ p \\ \dot{c} \\ \dot{p} \end{pmatrix} = \begin{pmatrix} \dot{c} \\ \dot{p} \\ \eta^2 (c-z) + g \\ \ddot{p} \end{pmatrix}, \qquad (19)$$

with $g = (0, 0, -9.81)$ m/s$^2$ and $\eta^2 = g/\bar{h}$, where $\bar{h}$ is a target CoM height.

The dynamics in the second stage are also discretized with Euler integration:

$$x_{k+1}^2 = f_{\text{LIP}}(x_k^2, u_k^2). \qquad (20)$$

The cost terms for the two stages and the terminal stage

<footnote n="4">⁴Here we are using line feet assumption.</footnote>
are respectively:

$$
\begin{align*}
\ell^1(\mathbf{x}, \mathbf{u}) = & \underbrace{\alpha_z (c^z - \bar{h})^2}_{\text{CoM z tracking}} + \underbrace{\alpha_\omega \|\boldsymbol{\omega} - \boldsymbol{\omega}^{\text{ref}}\|^2}_{\text{orientation tracking}} + \underbrace{\alpha_{\dot{r}} \|\dot{\boldsymbol{c}} - \dot{\boldsymbol{c}}^{\text{ref}}\|^2}_{\text{CoM velocity}} \\
& + \underbrace{\alpha_{\Delta p_1} \|\boldsymbol{p}_0^{xy} - \boldsymbol{p}_2^{xy} - \boldsymbol{d}_L\|^2 + \alpha_{\Delta p_2} \|\boldsymbol{p}_1^{xy} - \boldsymbol{p}_3^{xy} - \boldsymbol{d}_R\|^2}_{\text{contact spacing}} \\
& + \underbrace{\alpha_{\ddot{c}} \|\ddot{\boldsymbol{c}}\|^2 + \alpha_{\dot{\omega}} \|\dot{\boldsymbol{\omega}}\|^2}_{\text{acceleration}} \\
& + \sum_{i=0}^{n_c-1} \left( \underbrace{\alpha_f \|\boldsymbol{f}_i\|^2}_{\text{force regularization}} + \underbrace{\alpha_{f,\text{sw}} (1-\sigma_i)^2 \|\boldsymbol{f}_i\|^2}_{\text{contact switching}} \right),
\end{align*}
$$

$$
\begin{align*}
\ell^2(\mathbf{x}, \mathbf{u}) = & \underbrace{\alpha_z (c^z - \bar{h})^2}_{\text{CoM z tracking}} + \underbrace{\alpha_{xy} \|\mathbf{c}^{xy} - \frac{1}{4} \sum_{i=0}^{n_c-1} \mathbf{p}_i^{xy}\|^2}_{\substack{\text{CoM xy to support center} \\ \text{ZMP xy tracking}}} + \underbrace{\alpha_{\dot{c}} \|\dot{\mathbf{c}} - \dot{\mathbf{c}}^{\text{ref}}\|^2}_{\text{CoM velocity}} \\
& + \underbrace{\alpha_{\text{zmp},xy} \|\mathbf{z}^{xy} - \frac{\sum_{i=0}^{n_c-1} \sigma_i \mathbf{p}_i^{xy}}{\sum_{i=0}^{n_c-1} \sigma_i}\|^2}_{\text{ZMP xy tracking}} + \underbrace{\alpha_{\text{zmp},z} \|\mathbf{z}^z - \mathbf{z}_{\text{ground}}\|^2}_{\text{ZMP z tracking}} \\
& + \underbrace{\alpha_{\Delta p_1} \|\mathbf{p}_0^{xy} - \mathbf{p}_2^{xy} - \mathbf{d}_L\|^2 + \alpha_{\Delta p_2} \|\mathbf{p}_1^{xy} - \mathbf{p}_3^{xy} - \mathbf{d}_R\|^2}_{\text{contact spacing}} \\
& + \underbrace{\alpha_{\dot{p}} \|\ddot{\mathbf{p}}\|^2}_{\text{contact acceleration}}
\end{align*}
$$

$$
\ell_f(\mathbf{x}) = \underbrace{\alpha_z (c^z - \bar{h})^2}_{\text{CoM z tracking}} + \underbrace{\alpha_{xy} \|c^{xy} - \frac{1}{4} \sum_{i=0}^{3} c_i^{xy}\|^2}_{\text{CoM xy tracking}} + \underbrace{\alpha_{\dot{r}} \|\dot{c} - \dot{c}^{\text{ref}}\|^2}_{\text{CoM velocity tracking}} \\
+ \underbrace{\alpha_{\Delta p_1} \|\mathbf{p}_0^{xy} - \mathbf{p}_2^{xy} - \mathbf{d}_L\|^2 + \alpha_{\Delta p_2} \|\mathbf{p}_1^{xy} - \mathbf{p}_3^{xy} - \mathbf{d}_R\|^2}_{\text{contact spacing}}
$$

where $\dot{\mathbf{c}}^{\text{ref}}$ and $\omega^{\text{ref}}$ are reference velocities that we generate by reading keyboard input from the user, $\mathbf{d}_R \in \mathbb{R}^2$ and $\mathbf{d}_L \in \mathbb{R}^2$ are vectors dependent on the foot sizes that keep the contact positions aligned with the robot feet, $\sigma_i$ is a boolean that indicates whether contact $i$ is currently active.

The mapping function simply selects the CoM and contact positions to be carried over to the next stage, i.e.,

$$
g(\mathbf{x}_{N_0}^1) = \begin{pmatrix} c_{N_0}^1 \\ \mathbf{p}_{N_0}^1 \\ \dot{\mathbf{c}}_{N_0}^1 \\ \dot{\mathbf{p}}_{N_0}^1 \end{pmatrix}. \qquad (21)
$$

A. *Humanoid simulations*

The first simulations show the model of the humanoid
with two DDP prediction stages, which have $N_1 = 10$ and
$N_2 = 10$, with a total time $T = 1$ s. We used the Horizon
framework [18] to write the MPC optimal control problem.

Figure 1 reports the computed trajectories during walking and jumping examples, visualized open-loop in RVIZ.

Figure 2 shows the trajectories of the states and control predictions at a specific time during the simulation. State trajectories are continuous during the transition. The inputs are not continuous, nor are they required to be; in fact, they may represent entirely different quantities. Note that the first stage also predicts angular quantities not represented in the plot.

Figure 3 presents a comparison of the solution times when using only the SRBD, only the LIP, and the proposed model scheduling approach, all with the same number of nodes. It can be observed that the proposed method requires a computation time comparable to that of the LIP. In fact, it appears to be even lower, although the difference is very minor and likely due to some implementation quirk that we were unable to identify.

IV. EXAMPLE 2: DDP WITH TRAINED POLICY FOR CART-PENDULUM

In the second example, we also use a two-stage approach,
but the stages use different techniques, rather than different
models. In particular, the first stage is an iLQR on the actual
system model, while the second stage is an actor-critic pair
trained with an RL technique.

We consider a cart-pendulum system with state

$$
\mathbf{x}^{1,2} = \begin{pmatrix} p \\ \theta \\ \dot{p} \\ \dot{\theta} \end{pmatrix} \in \mathbb{R}^4, \quad \mathbf{u}^{1,2} = (f) \in \mathbb{R}. \qquad (22)
$$

where $p \in \mathbb{R}$ is the cart position, $\theta \in \text{SO}(2)$ is the pendulum angle, and $f \in \mathbb{R}$ is the horizontal force on the cart. Note that the state variables are the same in both stages. The dynamics of the system are

$$
\dot{\mathbf{x}} = \begin{pmatrix} \dot{p} \\ \dot{\theta} \\ M(\theta)^{-1} \left( \begin{pmatrix} f \\ 0 \end{pmatrix} - \mathbf{n}(\mathbf{x}) \right) \end{pmatrix}, \quad (23)
$$

where $\mathbf{n}(\mathbf{x})$ collects the nonlinear terms. The first stage
minimizes only the sum of square inputs, i.e.

$$
\ell_i^1(\mathbf{x}, \mathbf{u}) = \mathbf{u}^T \mathbf{u}, \quad i = 0, \dots, N_1. \qquad (24)
$$

The horizon of the first stage is $N_1 = 30$, and the discretized
dynamics has a sampling time of 0.01 s.

The second stage is a trained actor-critic pair trained with PPO. Both networks have 2 hidden layers of size 64, with tanh activation functions. This is particularly important, as using ReLU would not allow computation of the Hessian, because it can only produce piecewise linear approximations. The networks were trained over a number of episodes, interrupted when the total number of samples reached 10⁶. All episodes are initialized near the upright configuration $\mathbf{x}_{\text{up}} = (0, \pi, 0, 0)$, shifted on each component by a random amount sampled uniformly from [-0.1, 0.1]. Episodes were terminated under any of the following conditions:
<figure><img src="image_4.png" /><figcaption>Fig. 2: Predicted trajectories for states and controls at one point in time, using the double DDP approach for humanoid control. The plots are split to highlight the different stages using, respectively, the SRBD (left, continuous lines) and LIP (right, dashed lines) models. States are continuous because they are common to both models, whereas the inputs are different (forces in the first stage, ZMP in the second). Angular components are left out for clarity.</figcaption></figure>

* $|p| \ge 2$ m;

* $|\theta - \pi| \ge \pi/4$ rad;

* $T \ge 10$ s;

where $T$ denotes the total episode duration.

The reward function is computed as:

$$r(\mathbf{x}, \mathbf{u}) = -\ell^2(\mathbf{x}, \mathbf{u}) = 10^3 \left( 1 - 0.1(\mathbf{x} - \mathbf{x}_{\text{up}})^T (\mathbf{x} - \mathbf{x}_{\text{up}}) \right).$$

Note that, since the RL is solving a maximization problem, whereas the DDP stage expects a value function associated with a minimization problem, the critic approximation and its derivative should be multiplied by $-1$ when fed to the preceding stage. Table I reports the hyperparameters used during training.

### A. Cart-Pendulum simulations

This subsection presents simulations on the cart-pendulum model. The MPC optimal control problem was written using Horizon, and the training of the networks was performed

<table><caption>TABLE I: Hyperparameters used during training.</caption><thead><tr><td>Parameter</td><td>Value</td></tr></thead><tbody><tr><td>γ</td><td>0.99</td></tr><tr><td>λ<sub>GAE</sub></td><td>0.95</td></tr><tr><td>n<sub>step</sub></td><td>2048</td></tr><tr><td>Batch size</td><td>64</td></tr><tr><td>Number of epochs</td><td>16</td></tr><tr><td>Clipping ϵ</td><td>0.1</td></tr><tr><td>Value function coefficient</td><td>1.0</td></tr><tr><td>Entropy coefficient</td><td>0.0</td></tr><tr><td>Max gradient norm</td><td>1.0</td></tr><tr><td>Evaluation frequency (episodes)</td><td>4</td></tr><tr><td>Total timesteps</td><td>100,000 × 10</td></tr></tbody></table>

using PyTorch. Figure 4 shows results obtained on the cart-pendulum. Since the training did not cover the entire state space, the policy alone is not capable of swinging up the pendulum from the downward configuration (left plots in the figure). However, chaining the two stages together achieves a correct swing-up maneuver (right plots in the figure).
<figure><img src="image_5.png" /><figcaption>Fig. 3: MPC solution time comparison using: LIP model only (blue), SRBD model only (orange), our two-stage approach that combines both models (green).</figcaption></figure>

The initial guess used to linearize the iLQR during the first iteration is chosen as a linear interpolation of the state from the start to the goal. This is fundamental because the network approximating the value function must be evaluated in states within the basin considered during training.

## V. DISCUSSION

While the results presented here are preliminary, we believe they show the potential of this approach. Using different models throughout the optimization horizon provides great advantages in terms of foresight capability, at a relatively low computational cost, as demonstrated in the double DDP example of Section III.

The cart-pendulum example presented in Section IV, on the other hand, demonstrates that it can be much more general than just a model change, and that interesting results can also be obtained by mixing radically different techniques. Specifically, we showed that good performance can be obtained even using low-quality training, as the policy by itself cannot execute the swing up, whereas the combination of DDP and the trained policy can.

In this second example, the reader might question whether the second stage is necessary at all. After all, trajectory optimization by itself can certainly execute a successful cart-pendulum swing-up motion. This is certainly true, but it is only because it is very easy to select a terminal condition for the optimization that executes the job, i.e., the distance between the terminal state and the upright configuration must be heavily penalized.

More complex systems, such as legged robots, do not allow for such simplistic terminal conditions, so it is often necessary to resort to longer horizons, and possibly heuristics and/or pre-optimized reference trajectories. However, if one could obtain a network value function that correctly encodes the quality of the final state, one could greatly reduce the necessity of excessively tuning the optimization, and possibly

get away with a much simpler objective and a shorter horizon.

Using the policy by itself, which would be the other extreme, is certainly viable, but we reckon that would not necessarily be the best solution. The optimization stage typically offers a higher quality trajectory, more consistency, and also, using an appropriate solution technique, the possibility to include constraints.

<figure><img src="image_6.png" /><figcaption>Fig. 4: Failed attempt at swing-up from a downward configuration using policy only (left). Successful swing-up from a downward configuration using a combination of DDP and the value function approximation provided by the critic network (right).</figcaption></figure>

## VI. CONCLUSION AND FUTURE WORK

In this work, we introduced the concept of Cascaded Predictive Control, where the prediction horizon is structured as a sequence of stages, each optimizing a different representation of the robot's trajectory. Building upon ideas from [8], we extend the approach beyond switching models by integrating entirely different optimization techniques, specifically demonstrating the combination of DDP with an actor-critic RL. The proposed method has been validated through numerical examples, including model scheduling of templates for humanoid bipedal locomotion and a cart-pendulum swing-up task, where DDP is seamlessly combined with a value function approximation network trained using PPO. This second example, in particular, highlights an elegant and effective way to merge model-based and learning-based control techniques.

Future work will focus on showcasing the CPC algorithm by staging multiple solver techniques and models together. In particular, the method could involve an initial stage using the robot's full-body model, followed by a simplified template model, with each stage solved using different techniques, for example, Sequential Quadratic Programming (SQP) for the initial stage to handle complex constraints, DDP for the intermediate stage, and finally a learned value function approximation for the last stage.

## REFERENCES

[1] C. Khazoom, S. Hong, M. Chignoli, E. Stanger-Jones, and S. Kim, "Tailoring solution accuracy for fast whole-body model predictive control of legged robots," *IEEE Robotics and Automation Letters*, 2024.
[2] T. He, J. Gao, W. Xiao, Y. Zhang, Z. Wang, J. Wang, Z. Luo, G. He, N. Sobanbabu, C. Pan, Z. Yi, G. Qu, K. Kitani, J. Hodgins, L. J. Fan, Y. Zhu, C. Liu, and G. Shi, “Asap: Aligning simulation and real-world physics for learning agile humanoid whole-body skills,” *arXiv*, 2025.

[3] Y. Tassa, N. Mansard, and E. Todorov, “Control-limited differential dynamic programming,” in *IEEE International Conference on Robotics and Automation (ICRA)*, 2014, pp. 1168–1175.

[4] P. M. Wensing, M. Posa, Y. Hu, A. Escande, N. Mansard, and A. D. Prete, “Optimization-based control for dynamic legged robots,” *IEEE Transactions on Robotics*, vol. 39, no. 4, pp. 2810–2832, 2023.

[5] C. Mastalli, R. Budhiraja, W. Merkt, G. Saurel, B. Hammoud, M. Naveau, J. Carpentier, L. Righetti, S. Vijayakumar, and N. Mansard, “Crocoddyl: An efficient and versatile framework for multi-contact optimal control,” in *IEEE International Conference on Robotics and Automation (ICRA)*, 2020, pp. 2536–2542.

[6] T. A. Howell, B. E. Jackson, and Z. Manchester, “Altro: A fast solver for constrained trajectory optimization,” in *2019 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)*. IEEE, 2019, pp. 7674–7679.

[7] W. Jallet, A. Bambade, E. Arlaud, S. El-Kazdadi, N. Mansard, and J. Carpentier, “Proxxddp: Proximal constrained trajectory optimization,” *IEEE Transactions on Robotics*, 2025.

[8] H. Li and P. M. Wensing, “CaFe-mpc: A cascaded-fidelity model predictive control framework with tuning-free whole-body control,” *IEEE Transactions on Robotics*, 2024.

[9] E. Dantec, W. Jallet, and J. Carpentier, “From centroidal to whole-body models for legged locomotion: a comparative analysis,” in *IEEE-RAS International Conference on Humanoid Robots (Humanoids)*, 2024, pp. 335–342.

[10] A. Romero, E. Aljalbout, Y. Song, and D. Scaramuzza, “Actor-critic model predictive control: Differentiable optimization meets reinforce-

ment learning,” *arXiv preprint*, 2024.

[11] G. Grandesso, E. Alboni, G. P. R. Papini, P. M. Wensing, and A. D. Prete, “Cacto: Continuous actor-critic with trajectory optimization—towards global optimality,” *IEEE Robotics and Automation Letters*, vol. 8, no. 6, pp. 3318–3325, 2023.

[12] G. Turrisi, B. B. Carlos, M. Cefalo, V. Modugno, L. Lanari, and G. Oriolo, “Enforcing constraints over learned policies via nonlinear mpc: Application to the pendubot,” *IFAC-PapersOnLine*, vol. 53, no. 2, pp. 9502–9507, 2020.

[13] S. H. Bang, C. A. Jové, and L. Sentis, “Rl-augmented mpc framework for agile and robust bipedal footstep locomotion planning and control,” in *IEEE-RAS International Conference on Humanoid Robots (Humanoids)*, 2024, pp. 607–614.

[14] Y. Chen and Q. Nguyen, “Learning agile locomotion and adaptive behaviors via rl-augmented mpc,” in *IEEE International Conference on Robotics and Automation (ICRA)*, 2024, pp. 11436–11442.

[15] S. Singh, R. P. Russell, and P. M. Wensing, “Analytical second-order derivatives of rigid-body contact dynamics: Application to multi-shooting ddp,” in *IEEE-RAS International Conference on Humanoid Robots (Humanoids)*, 2023, pp. 1–8.

[16] W. Jallet, A. Bambade, N. Mansard, and J. Carpentier, “Constrained differential dynamic programming: A primal-dual augmented la-grangian approach,” in *IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)*, 2022, pp. 13371–13378.

[17] M. Cipriano, P. Ferrari, N. Scianca, L. Lanari, and G. Oriolo, “Humanoid motion generation in a world of stairs,” *Robotics and Autonomous Systems*, vol. 168, p. 104495, 2023.

[18] F. Ruscelli, A. Laurenzi, N. G. Tsagarakis, and E. Mingo Hoffman, “Horizon: A trajectory optimization framework for robotic systems,” *Frontiers in Robotics and AI*, vol. 9, 2022.