Indirect Optimal Control for Constructing Libraries of Periodic Trajectories: Application to Bipedal Walking

Maximilian Raff; Kathrin Flaßkamp; C. David Remy

cs.ROOct 12, 2024

Indirect Optimal Control for Constructing Libraries of Periodic Trajectories: Application to Bipedal Walking

Authors: Maximilian Raff, Kathrin Flaßkamp, C. David Remy

Organizations: Institute for Adaptive Mechanical Systems, University of Stuttgart, Germany. · Systems Modeling and Simulation, Saarland University, Germany

Abstract

Trajectory optimization is a fundamental tool for generating efficient, dynamically consistent periodic motions, commonly referred to as gaits, in legged locomotion. This paper explores the indirect method of periodic trajectory optimization, emphasizing its application in creating optimal gaits for legged systems and contrasting it with the more common direct method. While the direct method offers implementation flexibility, it is limited by its need for an input space parameterization. In contrast, the indirect method achieves higher accuracy by computing the control input from states and costates along the optimal trajectory. In this work, we tackle the convergence challenges associated with indirect shooting methods by utilizing numerical continuation methods. This is particularly useful for the systematic development of gait libraries. Our contributions include: (1) the formalization of a general periodic trajectory optimization problem that extends existing first-order necessary conditions to a broader range of cost functions and operating conditions; (2) a methodology for efficiently generating optimal gait libraries using single shooting and numerical continuation; (3) a novel approach for reconstructing Lagrange multipliers and costates from passive gaits; (4) a comparative analysis of the indirect and direct shooting methods using a compass-gait walker as a case study, demonstrating the higher accuracy of the indirect method; and (5) demonstrating applicability to the more complex legged robot RABBIT, with ten dynamic states and four inputs. These results establish a systematic indirect optimal control framework for synthesizing families of optimal periodic trajectories in legged locomotion.

Figures & tables

Figure 1: The compass-gait walker on a slope

\gamma

, with parameters normalized by total mass

m

, gravity

g

, and leg length

l_{\mathrm{o}}

, and proportions

\nicefrac{{a}}{{b}}=1

and

\nicefrac{{m_{\mathrm{h}}}}{{m_{\mathrm{l}}}}=2

Panel of Figure 2 — Figure 2: Projections of the implicit curve $\boldsymbol{r}^{-1}(\boldsymbol{0})$ for the compass-gait walker (CGW), highlighted in varying shades of red and computed using Algorithm 1 with continuation parameter $\sigma:=\gamma$ and fixed average speed $v_{\text{avg}}=0.1~\sqrt{gl_{\text{o}}}$ . Towards level ground $\gamma=0^{\circ}$ (marked by $\ast$ ), the curve $\boldsymbol{r}^{-1}(\boldsymbol{0})$ connects the passive gait families with longer ( $\Gamma_{\text{pas.-l.}}$ ) and shorter ( $\Gamma_{\text{pas.-s.}}$ ) periods at slopes $\gamma_{\text{pas.-l.}}$ and $\gamma_{\text{pas.-s.}}$ , respectively. Filled circles mark where passive and optimal gaits coincide.

Panel of Figure 3 — Figure 3: Visualization of an optimal walking gait of the CGW at an average speed of $v_{\text{avg}}=0.1\sqrt{gl_{\text{o}}}$ on level ground ( $\gamma=0^{\circ}$ ). Keyframes illustrate the walking motion, while input and state trajectories (swing leg ”sw”, stance leg ”st”) compare indirect (red) and direct (black) shooting methods. The direct method parameterizes the input with a single cubic B-spline ( $n_{\xi}=4$ ).

Panel of Figure 4 — Figure 4: Comparison of direct ( $\hat{\bar{\boldsymbol{r}}}$ ) and indirect ( $\bar{\boldsymbol{r}}$ ) approaches for the CGW, implemented in MATLAB on an Intel Core i5-8500 (3.00 GHz, 8 GB RAM). The direct method parameterizes the input using Bézier polynomials or cubic B-splines with $n_{\xi}$ variables. Results in (a) and (b) correspond to level-ground walking ( $\gamma=0^{\circ}$ ) at $v_{\text{avg}}=0.1\sqrt{gl_{\text{o}}}$ . The corresponding (locally) optimal gaits were generated via continuation from a passive gait using Algorithm 1 , as shown in Figure 2 .

Figure 5: Projections of the implicit curves from the indirect ( $\boldsymbol{r}^{-1}(\boldsymbol{0})$ ) and direct ( $\hat{\boldsymbol{r}}^{-1}(\boldsymbol{0})$ ) methods, computed using Algorithm 1 with continuation parameter $\sigma=v_{\text{avg}}$ . Both curves describe level-ground walking ( $\gamma=0^{\circ}$ ), starting from an optimal gait at $v_{\text{avg}}=0.1~\sqrt{gl_{\text{o}}}$ (marked by $\ast$ ). The indirect method consistently achieves lower costs, while the direct method (utilizing cubic B-splines with $n_{\xi}=4$ ) distinguishes local minima (solid black) from saddle points (dotted black) using second-order optimality conditions.

Figure 6: The RABBIT robot on a slope $\gamma$ , with parameters taken from ( Chevallereau et al., 2003 ; Westervelt et al., 2007 ) normalized by total mass $m$ , gravity $g$ , and leg length $l_{\mathrm{o}}$ .

Panel of Figure 7 — Figure 7: Projections of the implicit curve $\boldsymbol{r}^{-1}(\boldsymbol{0})$ for the robot RABBIT, computed using Algorithm 1 with the continuation parameter $\sigma:=\gamma$ . The curve corresponds to a fixed average speed of $v_{\text{avg}}=0.05~\sqrt{gl_{\text{o}}}$ and originates from a passive gait family at slope $\gamma_{\text{pas}}$ . The filled circle indicates the point where the passive and optimal gait coincide.

Panel of Figure 8 — Figure 8: Optimal walking gait of RABBIT at an average speed of $v_{\text{avg}}=0.05\sqrt{gl_{\text{o}}}$ on level ground ( $\gamma=0^{\circ}$ ). Keyframes show successive configurations of RABBIT throughout the gait cycle, while input trajectories compare indirect (red) and direct (black) shooting methods. The direct method parameterizes each input $\hat{u}_{i}$ with a cubic B-spline ( $n_{\xi}=4$ ).

Appendix figures & tables1 asset

Supplementary material from the paper’s appendix.

Appendix


$\displaystyle\frac{\mathrm{d}L_{1}}{\mathrm{d}\varepsilon}\bigg\|_{\varepsilon=0}=$	$\displaystyle~\frac{\mathrm{d}c}{\mathrm{d}\varepsilon}\bigg(\tilde{T},\boldsymbol{z}^{\ast}(\tilde{T})+\varepsilon\cdot\delta\boldsymbol{z}(\tilde{T})\bigg)\bigg\|_{\varepsilon=0}+\frac{\mathrm{d}}{\mathrm{d}\varepsilon}\bigg(\tilde{\boldsymbol{\lambda}}^{\mathop{\mathrm{T}}}\boldsymbol{h}\big(\tilde{T},~\boldsymbol{x}^{\ast}(\tilde{T})+\varepsilon\delta\boldsymbol{x}(\tilde{T}),~\boldsymbol{x}^{\ast}(0)+\varepsilon\delta\boldsymbol{x}(0)\big)\bigg)\bigg\|_{\varepsilon=0}$	(28a)
$\displaystyle=$	$\displaystyle~\left[\frac{\partial c}{\partial T}\bigg\|_{\ast}+\frac{\partial c}{\partial\boldsymbol{z}_{T}}\bigg\|_{\ast}\dot{\boldsymbol{z}}^{\ast}{(T^{\ast})}+{\boldsymbol{\lambda}^{\ast}}^{\mathop{\mathrm{T}}}\left(\frac{\partial\boldsymbol{h}}{\partial T}\bigg\|_{\ast}+\frac{\partial\boldsymbol{h}}{\partial\boldsymbol{x}_{T}}\bigg\|_{\ast}\dot{\boldsymbol{x}}^{\ast}{(T^{\ast})}\right)\right]\delta T+\delta\boldsymbol{\lambda}^{\mathop{\mathrm{T}}}\boldsymbol{h}\|_{\ast}$
	$\displaystyle~+\frac{\partial c}{\partial y_{T}}\bigg\|_{\ast}\delta y(T^{\ast})+\left[\frac{\partial c}{\partial\boldsymbol{x}_{T}}\bigg\|_{\ast}+{\boldsymbol{\lambda}^{\ast}}^{\mathop{\mathrm{T}}}\frac{\partial\boldsymbol{h}}{\partial\boldsymbol{x}_{T}}\bigg\|_{\ast}\right]\delta\boldsymbol{x}(T^{\ast})+{\boldsymbol{\lambda}^{\ast}}^{\mathop{\mathrm{T}}}\frac{\partial\boldsymbol{h}}{\partial\boldsymbol{x}_{0}}\bigg\|_{\ast}\delta\boldsymbol{x}(0),$	(28b)

Appendix

Table 9

Bipedal Locomotion

See all 22 papers

Oct 14, 2026

LEAP: Gradient Flow Augmented Lagrangian Trajectory Optimization for Contact Constrained Robots

Generating dynamically feasible motions for legged and manipulating robots requires solving a large, nonlinear, non-convex optimal control problem (OCP) in which contact wrenches must satisfy friction-cone constraints throughout the motion, not only at isolated instants. Existing solvers are sensitive to the initial guess and scale poorly with the number of joints. This paper introduces LEAP, which transcribes the OCP with pseudospectral collocation and, instead of a Newton or SQP solve, integrates a primal--dual gradient flow of a Powell--Hestenes--Rockafellar augmented Lagrangian to compute a solution. Because the flow is stiff, this paper designs a linearly implicit integrator with a Gauss--Newton metric, momentum, and a merit-based step controller that keeps the augmented Lagrangian monotone. Across locomotion benchmarks on five robot platforms, from planar hoppers to 32-joint humanoids, LEAP is the only method to solve every tested problem, achieves substantial speedups over Crocoddyl and Aligator on several benchmarks, and reaches trackable references at smaller compute budgets. The resulting trajectories are validated by tracking them in a high-fidelity simulator and on hardware. We will make the code public upon final acceptance of this paper.

César E. Ramos Chuquiure, Ryan Tsai, Yutong Zhou +8

Department of Robotics, University of Michigan, Ann Arbor, MI, USA. · Computer Science and Engineering, University of Michigan, Ann Arbor, MI, USA. · Department of Electrical Engineering and Computer Science, University of Michigan, Ann Arbor, MI, USA.

Hypergradient Optimization Trajectory Optimization Robotic Optimization Augmented Lagrangian Method Robot Locomotion

Ask PDF

Sep 18, 2026

Duty Factor Predicts Robust Constrained Quadrupedal Locomotion Across Gait Types

Quadrupedal robots are increasingly deployed in environments where locomotion must remain robust to disturbances and constrained terrain. Gait type, such as walking or trotting, is commonly used to characterize quadrupedal locomotion. However, gait type does not uniquely define locomotion, as parameters such as duty factor, speed, and stance width can vary within a single gait type. In this work, we investigate the relationship between these gait parameters using three distinct quadrupedal locomotion control approaches. First, using whole body trajectory optimization with LQR feedback, we show that duty factor is a stronger predictor of local error convergence than nominal gait type. Second, we investigate duty factor selection with a learned locomotion controller, suggesting how duty factor may serve as a low-dimensional parameter for adapting locomotion robustness in narrow-terrain environments. Finally, we show that these trends persist under a centroidal model predictive control framework and validate them through narrow-terrain experiments on a physical quadruped. These results show that duty factor provides a simple and effective basis for understanding and selecting robust quadrupedal locomotion across gait types and control architectures.

James Zhu, David Ologan, George Ortiz +5

Department of Electrical and Computer Engineering, University of Miami, Coral Gables, FL, USA · Department of Mechanical Engineering, 3 Robotics Institute, Carnegie Mellon University, Pittsburgh, PA, USA

Compliant Legged Locomotion Quadruped Locomotion Legged Robotics Gait Analysis

Ask PDF

May 7, 2026

On the Emergence of Pendular Structure in Multi-Contact Locomotion

LIPM is everywhere in legged-locomotion control, but almost always as a modeling choice rather than as something the controller's cost actually prefers. This note tries to make that link more explicit. Working from a small centroidal OCP that penalizes the rate of angular momentum, we look at what its optimum tends to look like. Three things come out. With full-rank stance, the optimum drifts toward a pendular force pattern at a rate determined by the SVD of the moment Jacobian; the constant is set by foot-span geometry and matches the experiments to within 16%. With N=2 stance, as in trot, the friction cone introduces a lower bound on

\|\dot{H}_G\|

that no amount of weight tuning fixes; we also see a non-smooth feasibility kink at a critical horizontal acceleration that we can write in closed form. Adding a task term that asks for a nonzero

\dot{H}_G

moves the optimum off the pendular set in a predictable way. None of this is far from the classical ZMP/DCM picture. We test these claims on a point-mass quadruped and on the Unitree Go1 in MuJoCo (open-loop QP and a torque-level closed-loop controller), and we note where the asymptotic story stops being a good description of what the closed loop actually does.

Lingxue Lyu, Zihui Liu

School of Engineering and Applied Science, University of Pennsylvania, Philadelphia, PA 19104, USA · Department of Aeronautics & Astronautics, Stanford University, Stanford, CA 94305, USA

Legged Robotics Trajectory Optimization Momentum Methods Optimal Control Quadruped Locomotion

Ask PDF


$\displaystyle\frac{\mathrm{d}L_{1}}{\mathrm{d}\varepsilon}\bigg\|_{\varepsilon=0}=$	$\displaystyle~\frac{\mathrm{d}c}{\mathrm{d}\varepsilon}\bigg(\tilde{T},\boldsymbol{z}^{\ast}(\tilde{T})+\varepsilon\cdot\delta\boldsymbol{z}(\tilde{T})\bigg)\bigg\|_{\varepsilon=0}+\frac{\mathrm{d}}{\mathrm{d}\varepsilon}\bigg(\tilde{\boldsymbol{\lambda}}^{\mathop{\mathrm{T}}}\boldsymbol{h}\big(\tilde{T},~\boldsymbol{x}^{\ast}(\tilde{T})+\varepsilon\delta\boldsymbol{x}(\tilde{T}),~\boldsymbol{x}^{\ast}(0)+\varepsilon\delta\boldsymbol{x}(0)\big)\bigg)\bigg\|_{\varepsilon=0}$	(28a)
$\displaystyle=$	$\displaystyle~\left[\frac{\partial c}{\partial T}\bigg\|_{\ast}+\frac{\partial c}{\partial\boldsymbol{z}_{T}}\bigg\|_{\ast}\dot{\boldsymbol{z}}^{\ast}{(T^{\ast})}+{\boldsymbol{\lambda}^{\ast}}^{\mathop{\mathrm{T}}}\left(\frac{\partial\boldsymbol{h}}{\partial T}\bigg\|_{\ast}+\frac{\partial\boldsymbol{h}}{\partial\boldsymbol{x}_{T}}\bigg\|_{\ast}\dot{\boldsymbol{x}}^{\ast}{(T^{\ast})}\right)\right]\delta T+\delta\boldsymbol{\lambda}^{\mathop{\mathrm{T}}}\boldsymbol{h}\|_{\ast}$
	$\displaystyle~+\frac{\partial c}{\partial y_{T}}\bigg\|_{\ast}\delta y(T^{\ast})+\left[\frac{\partial c}{\partial\boldsymbol{x}_{T}}\bigg\|_{\ast}+{\boldsymbol{\lambda}^{\ast}}^{\mathop{\mathrm{T}}}\frac{\partial\boldsymbol{h}}{\partial\boldsymbol{x}_{T}}\bigg\|_{\ast}\right]\delta\boldsymbol{x}(T^{\ast})+{\boldsymbol{\lambda}^{\ast}}^{\mathop{\mathrm{T}}}\frac{\partial\boldsymbol{h}}{\partial\boldsymbol{x}_{0}}\bigg\|_{\ast}\delta\boldsymbol{x}(0),$	(28b)

Abstract

Bipedal Locomotion

Explore similar work

LEAP: Gradient Flow Augmented Lagrangian Trajectory Optimization for Contact Constrained Robots

Duty Factor Predicts Robust Constrained Quadrupedal Locomotion Across Gait Types

On the Emergence of Pendular Structure in Multi-Contact Locomotion