Humanoid Robot Co-Design: Coupling Hardware Design with Gait Generation via Hybrid Zero Dynamics

Background and Motivation

Humanoid robots, owing to their human-like morphology, can naturally deploy in environments built for humans and have been a central robotics research direction for decades. With the emergence of cheap high-performance actuators, dynamic humanoid platforms have entered a new development phase. However, these platforms are extremely challenging to develop—mechanical design must balance lightweight construction with impact resistance while maintaining sufficient range of motion; controller synthesis faces the inherent instability and nonlinear hybrid dynamics of legged locomotion.

Due to these challenges, the mainstream approach separates design and control: first complete the mechanical design, then develop controllers on that basis. This serial pipeline leads to a time-consuming and expensive iterative process—design parameters are tightly coupled with controller performance, and changing a single leg-length parameter can completely alter the robot's dynamic behavior, requiring controller redesign. For bipedal robots, the situation is worse: the underlying nonlinear hybrid dynamics make it nearly impossible to intuitively assess whether a design parameter is good or bad.

Existing co-design methods attempt to break this serial barrier. Some use bi-level optimization—first generating motion plans for a simplified model (center of mass + foot positions), then optimizing link lengths and full-body motions. The Vitruvio toolbox uses a single rigid body dynamics (SRBD) model for motion planning combined with a genetic algorithm for leg parameter search. However, these methods either rely on simplified models that lose full hybrid dynamics information, or use gradient-free search with low efficiency. More importantly, they do not formally handle design parameters within trajectory optimization—design variables remain outer-loop search decisions, disconnected from gait generation.

This paper's starting point: can design parameters be directly embedded into the HZD gait generation trajectory optimization, making design and gait variables jointly solvable in a single framework? This avoids model simplification information loss and leverages HZD's rigorous treatment of hybrid dynamics. The answer is yes—through the concept of "design virtual constraints" that encode link lengths as configuration coordinates of virtual joints, placing them on equal footing with standard joint trajectory optimization variables.

Preliminaries

Hybrid Zero Dynamics Framework: Legged robot locomotion consists of multiple continuous phases (single support, double support) and discrete impact events (foot strike), forming a hybrid dynamical system. Given configuration coordinates $\mathbf{q}\in\mathcal{Q}\subset\mathbb{R}^n$ and full state $\mathbf{x}=(\mathbf{q},\dot{\mathbf{q}})$, the second-order mechanical system is:

$$D(\mathbf{q})\ddot{\mathbf{q}} + H(\mathbf{q},\dot{\mathbf{q}}) = B\mathbf{u}$$

where $D(\mathbf{q})\in\mathbb{R}^{n\times n}$ is the inertia matrix, $H(\mathbf{q},\dot{\mathbf{q}})\in\mathbb{R}^n$ is the drift vector (Coriolis and gravity terms), $B\in\mathbb{R}^{n\times m}$ is the actuation matrix, and $\mathbf{u}\in\mathbb{R}^m$ is the input. This equation captures the dynamic coupling between joints—thigh motion affects shin forces and vice versa.

The swing foot height is $z_{\text{nsf}}:\mathcal{Q}\to\mathbb{R}$, and the admissible state set within the domain is:

$$\mathcal{D} := \{(\mathbf{q},\dot{\mathbf{q}})\in\mathcal{X} \mid z_{\text{nsf}}(\mathbf{q}) \geq 0\}$$

The switching surface (foot strike triggers discrete transition) is:

$$\mathcal{S} := \{(\mathbf{q},\dot{\mathbf{q}})\in\mathcal{X} \mid z_{\text{nsf}}(\mathbf{q}) = 0,\;\dot{z}_{\text{nsf}}(\mathbf{q},\dot{\mathbf{q}}) < 0\}$$

The impact causes discrete state jumps via the reset map $\Delta:\mathcal{S}\to\mathcal{X}$:

$$\mathbf{x}^+ = \Delta(\mathbf{x}^-),\quad \mathbf{x}^-\in\mathcal{S}$$

Combining continuous dynamics and impact reset yields the single-domain hybrid control system:

$$\mathcal{HC} = \begin{cases} \dot{\mathbf{x}} = f(\mathbf{x}) + g(\mathbf{x})\mathbf{u} & \mathbf{x}\notin\mathcal{S} \\ \mathbf{x}^+ = \Delta(\mathbf{x}^-) & \mathbf{x}^-\in\mathcal{S} \end{cases}$$

Virtual Constraints and Zero Dynamics: HZD's core idea is dimension reduction via virtual constraints. Define actual outputs $\mathbf{y}_a(\mathbf{q})$ (typically joint angles) and desired outputs $\mathbf{y}_d(\tau(\mathbf{q}),\boldsymbol{\alpha})$ (parameterized by Bézier polynomials), with virtual constraints:

$$\mathbf{y}_\alpha(\mathbf{q}) := \mathbf{y}_a(\mathbf{q}) - \mathbf{y}_d(\tau(\mathbf{q}),\boldsymbol{\alpha})$$

where $\tau:\mathcal{Q}\to\mathbb{R}$ is the phasing variable and $\boldsymbol{\alpha}\in\mathbb{R}^{m\times(b+1)}$ are Bézier coefficients. The zero dynamics surface is:

$$\mathcal{Z}_\alpha := \{\mathbf{x}\in\mathcal{D} \mid \mathbf{y}_\alpha(\mathbf{q})=0,\;\dot{\mathbf{y}}_\alpha(\mathbf{q})=0\}$$

For overall hybrid system stability, outputs must be impact-invariant—the HZD condition:

$$\Delta(\mathcal{S}\cap\mathcal{Z}_\alpha) \subset \mathcal{Z}_\alpha$$

This ensures foot-strike impacts do not push the system off the zero dynamics surface, providing the mathematical guarantee for stable periodic walking.

Figure 1: The HZD co-design framework systematically selects design parameters (leg dimensions) for the ADAM humanoid.

Method

1. Encoding Design Parameters as Virtual Constraints

The core innovation is encoding design parameters (link lengths) as virtual constraints, making them jointly solvable with gait variables in a single HZD optimization. Consider $d$ design parameters representing link lengths. Link length $l_i$ is represented as a configuration coordinate by augmenting the system with $d$ prismatic virtual joints $\tilde{q}_i$—virtual joints are treated as regular joints in optimization but must be fixed in the real world. The virtual link $\tilde{l}_i$ is massless and serves as the parent of the next link $l_{i+1}$.

The augmented system configuration is:

$$\tilde{\mathbf{q}} := (\mathbf{q}^\top, \mathbf{q}_D^\top)^\top \in \tilde{\mathcal{Q}} \subset \mathbb{R}^{n+d}$$

This augmentation is intuitive: insert telescoping "virtual joints" at the hip and knee that the optimizer can freely adjust. Once optimization converges, the virtual joint values are the optimal leg lengths—letting the optimizer "design" the robot's leg dimensions. The augmented continuous dynamics are obtained by applying Euler-Lagrange equations to the augmented configuration coordinates.

2. Design Virtual Constraints

After design parameter augmentation, virtual constraints must be extended. The original constraint $\mathbf{y}_\alpha(\mathbf{q})$ drives joint angles to desired trajectories; design virtual constraints drive actual link lengths to desired design values. The augmented constraint combines both:

$$\mathbf{y}(\tilde{\mathbf{q}}) := \begin{bmatrix} \mathbf{y}_a(\mathbf{q}) - \mathbf{y}_d(\tau(\mathbf{q}),\boldsymbol{\alpha}) \\ \mathbf{q}_D - \boldsymbol{\beta} \end{bmatrix}$$

where $\boldsymbol{\beta}\in\mathbb{R}^d$ are Bézier coefficients for design parameters. The design virtual constraint $\mathbf{q}_D - \boldsymbol{\beta}$ means: the current virtual joint value (i.e., link length) is driven to the optimal value $\boldsymbol{\beta}$ determined by the optimizer. Just as standard virtual constraints are "fictitious" constraints enforced by actuators, design virtual constraints are "fictitious" constraints enforced by mechanical design—once $\boldsymbol{\beta}$ is determined, the robot's legs are manufactured to that length.

3. Extended HZD Optimization Problem

Adding augmented constraints to the standard HZD optimization yields the co-design problem:

$$\arg\min_{\boldsymbol{\alpha},\boldsymbol{\beta},\mathcal{X}} \Phi(\mathcal{X})$$

$$\text{s.t.}\;\dot{\tilde{\mathbf{x}}} = \tilde{f}(\tilde{\mathbf{x}}) + \tilde{g}(\tilde{\mathbf{x}})\mathbf{u}_D(\tilde{\mathbf{x}}) \quad\text{(Closed-Loop Dynamics)}$$

$$\mathbf{y}(\mathbf{v}(\boldsymbol{\alpha},\boldsymbol{\beta})) = \mathbf{0}_{m+d}$$

$$\mathcal{L}_{\tilde{f}_y}(\dot{\mathbf{v}}(\boldsymbol{\alpha},\boldsymbol{\beta})) = \mathbf{0}_{m+d} \quad\text{(Impact-Invariance)}$$

$$\mathcal{X}_{\min} \leq \mathcal{X} \leq \mathcal{X}_{\max} \quad\text{(Decision Variables)}$$

$$\mathbf{c}_{\min} \leq \mathbf{c}(\mathcal{X}) \leq \mathbf{c}_{\max} \quad\text{(Physical Constraints)}$$

$$\mathbf{a}_{\min} \leq \mathbf{a}(\mathcal{X}) \leq \mathbf{a}_{\max} \quad\text{(Essential Constraints)}$$

where $\Phi(\mathcal{X})$ is the cost function (weighted mechanical cost of transport for energy efficiency), physical constraints include friction cones, workspace constraints, and actuator limits, and essential constraints shape gait behavior (step length, duration, speed, foot clearance). The impact-invariance condition now covers both standard outputs and design parameters, ensuring coordination throughout the gait cycle.

Figure 2/3: Link extension method—design parameters encoded as configuration coordinates via virtual joints. Original and augmented ADAM models.

4. Accounting for Changing Inertia

A key challenge: when virtual joint values (link lengths) change, the links' inertial properties also change. Since virtual links are massless, augmented link inertia is entirely determined by design parameters. Solving under fixed inertia assumptions would yield inaccurate results.

The authors employ an iterative outer loop (Algorithm 1): each iteration updates the robot model's inertial properties based on current virtual joint values (via GenerateModel producing a new URDF), then solves the optimization on the updated model. This repeats until virtual joint convergence (change below threshold $\delta$). This step is particularly important when inertia varies significantly with link length.

5. Automatic Robot Model Generation

FROST requires URDF robot models. Normally URDFs are fixed, but the co-design iterative process requires generating URDFs for arbitrary virtual joint values. The GenerateModel($\mathbf{q}_D$) function creates URDFs with accurate kinematic and inertial properties. Since complex link geometry makes analytical inertia calculation intractable, the authors use an empirical approach: uniformly sampling inertial properties for various thigh and shin lengths, then fitting with 3rd-order polynomials. Despite limited samples, the fits achieve near-perfect accuracy—indicating smooth, fittable inertia-length relationships within reasonable ranges.

flowchart TB
    A[Initialize Design Params qD] --> B[GenerateModel: Create URDF]
    B --> C[FROST Trajectory Optimization]
    C --> D{Design Params Converged?}
    D -->|No| E[Update Inertia Properties]
    E --> B
    D -->|Yes| F[Output Optimal Leg Lengths + Gait]
    F --> G[Manufacture ADAM Hardware]
    G --> H[Experimental Gait Validation]

Experimental Results

ADAM Robot Model

ADAM is a 20-DoF humanoid built from a Unitree A1 quadruped's torso and actuators, with legs and arms designed from scratch. Each leg has 4 actuated DoFs: innermost controls yaw, next controls roll, last two control hip pitch and knee pitch. To optimize thigh and shin lengths, 4 virtual joints are added (left/right thigh + left/right shin), yielding a 24-DoF augmented model. Left-right symmetry is enforced via constraints $q_{LT}=q_{RT}$, $q_{LS}=q_{RS}$. For planar walking optimization, non-sagittal joints are fixed, effectively reducing to 15 DoFs. Optimization uses IPOPT in MATLAB with weighted mechanical cost of transport as the cost function.

Optimization Results

Figure 4 / Tables I-II: Link lengths converging from various initial values, with constraint parameters and statistics.

Parameter	Mean	Std
Thigh Length (m)	0.2425	0.0088
Shin Length (m)	0.2515	0.0084
Step Length (m)	0.2281	0.0281
Step Duration (s)	0.7994	0.0165
Walking Speed (m/s)	0.2853	0.0338

Table II: Statistics of optimal gaits from various initial design link lengths.

From 25 runs with different initial values, leg lengths consistently converge to similar values—optimal thigh 0.24m, optimal shin 0.25m, with standard deviations below 0.01m. Only a few cases converge to local optima, but their cost function values are clearly worse (lower energy efficiency), easily identifiable. The shin being slightly longer than the thigh matches human leg proportions, suggesting the optimizer independently "discovered" bio-inspired ratios.

Essential Constraint	Lower	Upper
Step Length (m)	0.15	0.35
Step Duration (s)	0.6	0.9
Thigh Length (m)	0.15	0.35
Shin Length (m)	0.15	0.35

Table I: Essential constraints used in the HZD optimization.

Experimental Validation

Figure 5/6: Phase portrait comparison of simulation vs. hardware, and gait tiles.

The optimal gait was validated on the optimized ADAM hardware. The robot walked on a treadmill with a boom restricting motion to the sagittal plane (note: the boom's additional mass and inertia were not modeled in gait generation). Step length 0.23m, step duration 0.8s, matching simulation. The gait was enforced via joint-level PD control with IMU-based torso pitch estimation. Experimental walking was stable and robust—stability reflected in sustained walking (phase portraits show limit cycle convergence), robustness in maintaining stability across 0.2-0.3 m/s speeds. The high similarity between simulated and hardware gaits validates the HZD co-design framework.

Limitations

Author-stated limitations: Experimental validation is limited to planar walking (physically constrained to sagittal plane by boom), and the boom's additional mass and inertia are not modeled, potentially causing sim-to-real discrepancies. The iterative inertia fitting relies on limited sampling points and 3rd-order polynomials, which may be insufficient when inertia changes drastically with link length. The optimization problem is highly nonlinear and non-convex, requiring good initial guesses (authors use a two-step approach), and is sensitive to initial conditions—occasionally converging to local optima.

Independent assessment: First, the method only optimizes link lengths, not actuator selection, joint configuration, or transmission ratios—incorporating these would require more complex modeling. Second, experiments only validate planar walking; 3D walking involves multi-domain hybrid dynamics (double support, heel-toe rolling), significantly increasing complexity. Third, the cost function uses weighted mechanical cost of transport, but energy-optimal does not equal robustness-optimal—different objectives may yield different design parameters, a sensitivity the paper does not explore.

Conclusion and Future Work

The core contribution is transforming robot design parameters from "a priori fixed values" to "decision variables in trajectory optimization," unifying design and gait within the HZD framework via design virtual constraints. The fundamental advantage: design parameter evaluation no longer relies on intuition or simplified models, but on rigorous trajectory optimization under full hybrid dynamics—the optimal leg length is the value that minimizes energy cost while satisfying all physical constraints. The ADAM experimental validation demonstrates end-to-end feasibility: from mathematical optimization to hardware manufacturing to gait realization.

Future directions include extending to 3D multi-domain walking, incorporating richer design variables (actuator parameters, transmission design), exploring multi-objective optimization (efficiency vs. robustness vs. agility), and transferring the framework to quadrupeds, exoskeletons, and other legged systems. The mathematical rigor of HZD combined with the engineering usability of FROST brings co-design from concept to practice.