Walking control of a biped humanoid robot has always been a challenging and exciting topic in the field of robotics. The complex structure of the humanoid body, with its many degrees of freedom and highly nonlinear dynamics, makes the generation of stable, efficient and natural walking gaits a difficult task. In this article, I present a comprehensive study of multi-objective gait planning for a small-sized biped humanoid robot, named Robonova-1. The work focuses on kinematics modeling, parameterized gait trajectory design, and optimization of gait parameters using a modern multi-objective evolutionary algorithm. The ultimate goal is to obtain a set of Pareto-optimal walking patterns that simultaneously balance stability, energy consumption and walking speed. I first describe the robot’s mechanical structure and derive its forward and inverse kinematics. Then, I decompose the walking cycle into three phases, namely starting, periodic walking and stopping, and use cubic spline interpolation to generate smooth trajectories for the ankle and hip joints. Key gait parameters are identified, and a parameterized gait model is constructed. To optimize these parameters, I formulate three objective functions corresponding to ZMP-based stability, energy consumption and walking speed. Since these objectives are conflicting, I adopt the Strength Pareto Evolutionary Algorithm 2 (SPEA2) to find a set of Pareto-optimal solutions. Simulation and real-robot experiments demonstrate the effectiveness of the proposed approach.
1. Introduction
Biped humanoid robots are designed to mimic human appearance and locomotion. Compared with wheeled robots, they can overcome unstructured terrains and are better suited for human-centered environments. During the past decades, many research groups around the world have developed impressive humanoid platforms. The early WABOT robots in Japan demonstrated the feasibility of biped walking. Later, the Honda P2/P3 and ASIMO robots achieved dynamic walking and even running at moderate speeds. In China, several institutions such as Harbin Institute of Technology, National University of Defense Technology and Beijing Institute of Technology have also made significant progress. These robots are equipped with multiple sensors and actuators, and their control systems are highly complex.
One of the most fundamental abilities of a humanoid robot is to walk in a stable and energy-efficient manner. Gait planning, which determines the joint trajectories of the legs, is the core of biped walking. Researchers have proposed various methods for gait planning, including dynamics-based approaches using inverted pendulum models, human motion capture, evolutionary optimization, passive dynamics and central pattern generators (CPG). Each method has its strengths and limitations. In this study, I focus on optimization-based gait planning because it provides a principled way to handle multiple conflicting performance criteria.
The key issue in optimization-based gait planning is how to define the objective functions and how to search for optimal gait parameters. Many existing works optimize only a single objective, such as stability or energy. However, a single-objective optimum does not necessarily yield a good overall performance. Walking speed, stability and energy consumption are often in conflict. Therefore, I treat the gait parameter optimization problem as a multi-objective optimization problem and use Pareto-based techniques to find a set of trade-off solutions. Among the available evolutionary algorithms, SPEA2 is known for its strength in approximating the Pareto front with good diversity and convergence. I choose SPEA2 to perform the search.
2. Structure and Kinematics Modeling of the Biped Humanoid Robot
2.1 Robot Description
The robot used in this work is the Robonova-1, a small humanoid robot designed for education and research. It is 310 mm tall, 180 mm wide and 90 mm deep, with a mass of about 1.3 kg. The robot has 16 degrees of freedom (DOF): each arm has 3 DOF (shoulder pitch, shoulder roll, elbow), and each leg has 5 DOF (hip pitch, hip roll, knee pitch, ankle pitch, ankle roll). The leg DOFs are essential for walking. The robot is driven by digital servo motors which provide position feedback and overload protection. The overall structure is composed of aluminum frames and plastic joints, making it light and strong.
The technical parameters of the robot are listed in the following table.
| Part | Mass (kg) | Length (mm) |
|---|---|---|
| Upper body | 0.54 | 90 |
| Lower leg | 0.06 | 46 |
| Upper leg | 0.13 | 84 |
| Joint bracket | 0.004 | 24 |
The ranges of motion for the main joints are given in Table 2, which are important constraints for gait planning.
| Joint | Range (deg) |
|---|---|
| Shoulder roll | -20 to 160 |
| Shoulder pitch | -40 to 140 |
| Elbow | 0 to 90 |
| Hip roll | -30 to 90 |
| Hip pitch | -90 to 80 |
| Knee | -90 to 10 |
| Ankle pitch | -90 to 90 |
| Ankle roll | -30 to 90 |
2.2 Forward Kinematics
Since the robot walks mainly in the sagittal plane (the radial plane), I set up a fixed coordinate frame with the origin at the midpoint of the two hip joints projected on the ground. The x-axis points forward, the y-axis points to the left, and the z-axis is vertical. The robot can be simplified as a seven-link model consisting of the upper body, two thighs, two shanks and two feet. Figure 2 (not shown here) depicts the link model. The links are connected by five joints labeled A to E in the sagittal plane.
Let mi be the mass of link i, and let li be its length. The angles between each link and the vertical axis are denoted by θi (i = 1,2,…,7). For the left leg supporting on the ground, the position of the ankle A is fixed at the origin. The forward kinematics equations give the positions of the hip and the swing ankle as functions of θi.
Let (xa, za) be the coordinates of the swing ankle, and (xh, zh) be the coordinates of the hip joint. The constraint equations are:
$$ \begin{aligned} x_h &= x_a + l_2 \sin\theta_2 + l_3 \sin\theta_3, \\ z_h &= z_a + l_2 \cos\theta_2 + l_3 \cos\theta_3, \end{aligned} $$
and for the swing ankle relative to the hip:
$$ \begin{aligned} x_a &= x_h – (l_5 \sin\theta_5 + l_6 \sin\theta_6), \\ z_a &= z_h – (l_5 \cos\theta_5 + l_6 \cos\theta_6). \end{aligned} $$
The center of mass (COM) of each link can be computed using the same transformation. For example, for link 2:
$$ \begin{aligned} x_{c2} &= a_2 \sin\theta_2, \\ z_{c2} &= l_{an} + a_2 \cos\theta_2, \end{aligned} $$
where a2 is the distance from the lower end of the link to its COM. The overall COM of the robot is then:
$$ \begin{aligned} x_G &= \frac{\sum_{i=1}^{7} m_i x_{ci}}{\sum_{i=1}^{7} m_i}, \\ z_G &= \frac{\sum_{i=1}^{7} m_i z_{ci}}{\sum_{i=1}^{7} m_i}. \end{aligned} $$
2.3 Inverse Kinematics
In the gait planning module, given the desired trajectory of the hip and the swing ankle, it is necessary to compute the joint angles. Using the geometric constraints of the seven-link model, I can solve the inverse kinematics uniquely. For the supporting leg, the hip position (xh, zh) and the ankle position (xl, zl) are known. The angles θ2 and θ3 can be obtained from:
$$ \theta_2 = \arccos\left( \frac{(z_h-z_l)^2 + (x_h-x_l)^2 + l_2^2 – l_3^2}{2 l_2 \sqrt{(z_h-z_l)^2 + (x_h-x_l)^2}} \right) – \arctan\left( \frac{x_h-x_l}{z_h-z_l} \right), $$
$$ \theta_3 = -\arccos\left( \frac{(z_h-z_l)^2 + (x_h-x_l)^2 + l_3^2 – l_2^2}{2 l_3 \sqrt{(z_h-z_l)^2 + (x_h-x_l)^2}} \right) – \arctan\left( \frac{x_h-x_l}{z_h-z_l} \right). $$
The trunk is kept upright, so θ4 = 0. For the swing leg, similar expressions give θ5 and θ6. Then the relative joint angles between adjacent links can be converted to servo commands.
3. Parameterized Gait Trajectory Design
3.1 Gait Phase Decomposition
A complete walking cycle consists of double-support phases and single-support phases. In the human gait, the double-support phase takes about 20% of the walking cycle. I decompose the whole walking process into three phases: starting, periodic walking, and stopping. The starting phase moves one leg from the initial upright posture to the beginning of the periodic walking. The periodic walking phase is repeated for each step. The stopping phase brings the legs back together from the periodic walking posture.
Within one walking period Tc, I define five key postures. The periodic walking phase can be described by the following parameters:
- Ds: step length
- Tc: walking cycle time
- Hao: maximum ankle height of the swing foot
- Lao: horizontal displacement of the ankle at the highest point
- Hmax: maximum hip height during the swing
- Hmin: minimum hip height during the swing
- qb: angle between the foot and ground when leaving the ground
- qf: angle between the foot and ground when touching down
- xsd: horizontal distance from the hip to the supporting ankle at the beginning of the single support
- xed: horizontal distance from the hip to the supporting ankle at the end of the single support
These ten parameters form a design vector x = (Ds, Tc, Lao, Hao, Hmax, Hmin, xsd, xed, qb, qf). The task of gait optimization is to find a suitable x that minimizes the desired criteria.
3.2 Spline Interpolation
Given a set of key postures at certain time instants, I need to generate all intermediate postures smoothly. I use cubic spline interpolation because it guarantees continuous first and second derivatives, which is important for preventing jerky motions. Suppose we have n points t1 < t2 < … < tn with function values fj. In each sub-interval [tj, tj+1], the spline S(t) is a cubic polynomial. The second derivative Mj = S”(tj) satisfies the continuity equations. After solving for Mj, the spline can be expressed as:
$$ S(t) = \frac{M_j (t_{j+1}-t)^3 + M_{j+1} (t-t_j)^3}{6h_j} + \left( f_j – \frac{M_j h_j^2}{6} \right)\frac{t_{j+1}-t}{h_j} + \left( f_{j+1} – \frac{M_{j+1} h_j^2}{6} \right)\frac{t-t_j}{h_j}, $$
where hj = tj+1 – tj. The natural boundary conditions set the second derivatives at the endpoints to zero.
3.3 Trajectory Formulations for the Ankle and Hip
For the periodic walking phase, the ankle position xa(t), za(t) and the foot pitch angle θa(t) are defined by key postures. I use a cubic polynomial of the form:
$$ x_a(t) = a_0 + a_1 (t-kT_c) + a_2 (t-kT_c)^2 + a_3 (t-kT_c)^3, $$
where the coefficients are determined by solving a Vandermonde system. Similarly, za(t) and θa(t) are also cubic polynomials. For the hip, I derive:
$$ x_h(t) = b_0 + b_1 (t-kT_c) + b_2 (t-kT_c)^2 + b_3 (t-kT_c)^3, $$
$$ z_h(t) = c_0 + c_1 (t-kT_c) + c_2 (t-kT_c)^2 + c_3 (t-kT_c)^3. $$
The boundary conditions ensure zero velocity at the start and at the end of the step, and continuous acceleration across steps.
For example, at the beginning of a step, the ankle is stationary on the ground, so the velocity is zero. At the end of the double-support phase, the foot leaves the ground with a specified angle qb. These conditions uniquely determine the coefficients.
4. Multi-Objective Optimization of Gait Parameters
4.1 Objective Functions
I now define the three objective functions that evaluate the quality of a given gait parameter vector x.
Stability objective: The most widely used criterion for dynamic walking stability is the Zero Moment Point (ZMP). The ZMP is the point on the ground where the horizontal component of the moment due to gravity and inertia forces is zero. During stable walking, the ZMP must stay within the convex hull of the supporting foot (or feet). To maximize stability, I minimize the difference between the actual ZMP and a desired ZMP trajectory that is kept at the center of the support polygon. Thus:
$$ f_1(x) = \frac{1}{N} \sum_{k=0}^{N-1} \left( x_{zmp}(kT_s) – x_{zmp}^{d}(kT_s) \right)^2 + \left( y_{zmp}(kT_s) – y_{zmp}^{d}(kT_s) \right)^2, $$
where Ts is the sampling period, N = Tc/Ts, and the superscript d denotes the desired value. The actual ZMP is computed from the COM accelerations of all links:
$$ x_{zmp} = \frac{\sum m_i (\ddot{z}_i + g) x_i – \sum m_i \ddot{x}_i z_i – \sum I_{iy} \ddot{\Omega}_{iy}}{\sum m_i (\ddot{z}_i + g)}, $$
$$ y_{zmp} = \frac{\sum m_i (\ddot{z}_i + g) y_i – \sum m_i \ddot{y}_i z_i – \sum I_{ix} \ddot{\Omega}_{ix}}{\sum m_i (\ddot{z}_i + g)}. $$
Energy objective: The energy consumption is related to the torque and angular velocity of each joint. Since the main contribution comes from the swing leg, I use a simplified measure that combines the swing ankle height (which affects the work done against gravity) and the angular acceleration of the swing leg. The angular acceleration is proportional to the required torque. The energy objective is:
$$ f_2(x) = \frac{1}{N} \sum_{k=0}^{N-1} \left( z_a(kT_s) + |\ddot{\theta}_5(kT_s)| \right), $$
where θ5 is the angle of the swing thigh. A smaller value means less energy consumption.
Speed objective: The walking speed is simply the step length divided by the cycle time, v = Ds/Tc. To convert to a minimization problem, I define:
$$ f_3(x) = \frac{T_c}{D_s} = \frac{1}{v}. $$
Clearly, minimizing f3 maximizes the speed.
These three objectives are conflicting. A higher speed generally requires more energy and may reduce stability. Therefore, I seek a set of Pareto-optimal solutions where no improvement in one objective can be achieved without worsening at least one other objective.
4.2 Constraint Handling
The gait parameters must satisfy physical constraints, such as joint angle limits and the requirement that the ZMP stays inside the support polygon during the whole walking cycle. I incorporate these constraints via a penalty function. If a candidate vector x violates any constraint, I add a large penalty to each objective value:
$$ F_i(x) = f_i(x) + C(x), \quad i=1,2,3, $$
where C(x) is a penalty term proportional to the degree of violation. This ensures that infeasible solutions are dominated by feasible ones.
4.3 The SPEA2 Algorithm
Strength Pareto Evolutionary Algorithm 2 (SPEA2) is an elitist multi-objective evolutionary algorithm. It maintains an archive of non-dominated solutions. The main steps are:
- Initialization: Generate an initial population P(0) of size N and an empty archive A(0).
- Fitness assignment: For each individual, compute its strength value (the number of solutions it dominates), then its raw fitness and a density estimate based on the k-th nearest neighbor distance. The final fitness is the sum of raw fitness and density.
- Environmental selection: Copy all non-dominated solutions from P(t) ∪ A(t) into the archive A(t+1). If the archive size exceeds M, use a truncation method to remove the least crowded solutions. If the size is smaller, fill with dominated individuals with low fitness.
- Termination: If the stopping criterion (e.g., maximum generations T) is met, output A(t+1).
- Mating selection: Apply binary tournament selection on A(t+1) to form a mating pool.
- Variation: Apply crossover and mutation operators to create a new population P(t+1).
The parameter encoding uses binary strings of length 10 for each decision variable, resulting in a chromosome length of 100. The algorithm parameters are set as: population size N = 80, archive size M = 80, maximum generations T = 100, crossover probability pc = 0.8, and mutation probability pm = 0.1.
4.4 Experimental Setup and Results
I used MATLAB 7.0 to implement the optimization. The constraints for the gait parameters are defined based on the robot’s dimensions:
$$ \begin{aligned} 2(l_{ab}+l_{af}) &\le D_s \le 3(l_{ab}+l_{af}), \\ l_{an} &\le H_{ao} \le l_{an} + l_{sh}, \\ \frac{l_{an}+l_{ab}+l_{af}}{2} &\le L_{ao} \le \frac{3(l_{ab}+l_{af})}{2}, \\ 0 &\le q_b, q_f \le \frac{\pi}{3}, \\ l_{an}+l_{sh} &\le H_{max} \le l_{an}+l_{sh}+l_{th}+l_{tr}, \\ l_{an}+l_{sh} &\le H_{min} \le H_{max}, \\ 0 &\le x_{sd}, x_{ed} \le l_{ab}+l_{af}, \\ 5 &\le T_c \le 15 \ \text{s}. \end{aligned} $$
To cover the trade-off region uniformly, I ran the algorithm and obtained a set of non-dominated solutions. The resulting Pareto front is shown in Figure 8 (not displayed here). The figure clearly illustrates the conflict among the three objectives. Three representative solutions were selected from the Pareto front and are listed in Table 3.
| Parameter | Gait 1 | Gait 2 | Gait 3 |
|---|---|---|---|
| Ds (mm) | 58.562 | 63.358 | 65.091 |
| Tc (s) | 1.812 | 2.304 | 2.007 |
| Lao (mm) | 48.100 | 51.002 | 52.911 |
| Hao (mm) | 37.960 | 41.018 | 43.815 |
| Hmax (mm) | 138.159 | 136.091 | 132.969 |
| Hmin (mm) | 131.461 | 132.034 | 129.990 |
| xsd (mm) | 17.988 | 15.841 | 17.461 |
| xed (mm) | 26.952 | 30.217 | 31.772 |
| qb (rad) | 0.4498 | 0.3466 | 0.4269 |
| qf (rad) | 0.3914 | 0.7717 | 0.9117 |
| f1(x) | 0.0613 | 0.0519 | 0.0679 |
| f2(x) | 0.2353 | 0.2582 | 0.3481 |
| f3(x) | 0.0309 | 0.0344 | 0.0308 |
Using these optimized parameters, I generated the joint trajectories. Figure 9 (not shown) displays the ankle trajectories for the three gaits, while Figure 10 shows the hip trajectories. The curves are smooth and continuous, confirming the effectiveness of the spline interpolation. The foot angle profile (Figure 11) also varies smoothly, ensuring a comfortable heel-to-toe transition.
To further assess the performance, I compared SPEA2 with a classical weighted-sum method. I ran both algorithms 20 times with the same parameter settings. The statistical results are summarized in Table 4.
| Objective | SPEA2 | Weighted-sum | ||
|---|---|---|---|---|
| Best | Average | Best | Average | |
| f1(x) | 0.0519 | 0.0601 | 0.0608 | 0.0792 |
| f2(x) | 0.2353 | 0.3025 | 0.2991 | 0.3856 |
| f3(x) | 0.0354 | 0.0376 | 0.0392 | 0.0419 |
| Runtime (s) | 30.6 | 34.2 | 39.1 | 45.5 |
The table and the corresponding bar charts in Figure 12 indicate that SPEA2 achieves better or comparable objective values and requires less runtime. The main reason is that the weighted-sum method requires a different run for each weight combination, and the performance strongly depends on the chosen weights. In contrast, SPEA2 simultaneously explores the entire Pareto front and provides multiple solutions in one run.
Finally, I implemented the optimized gait on the actual Robonova-1 robot. The robot successfully walked forward with a stable and smooth motion, as shown in Figure 13. The experiments confirm that the gait parameters obtained by the proposed multi-objective optimization method are feasible and effective in practice.
5. Conclusion and Future Work
In this article, I have presented a systematic framework for multi-objective gait planning of a biped humanoid robot. Starting from the kinematic model of the Robonova-1 robot, I derived the forward and inverse kinematics that are necessary for generating joint trajectories. By decomposing the walking process into three phases and using cubic spline interpolation, I constructed a parameterized gait model with ten key parameters. The gait parameters were optimized using the SPEA2 algorithm, with three conflicting objectives: ZMP-based stability, energy consumption and walking speed. The algorithm successfully produced a well-distributed Pareto front, giving engineers the flexibility to select a gait that best suits the current task. Simulation and real-robot experiments validated that the generated gaits are smooth, stable and energy-efficient.
Future work can extend this approach in several directions. First, the current planning is limited to flat ground walking. To handle stairs, slopes or uneven terrain, the gait model must incorporate terrain information and possibly adaptive control. Second, the objective functions could be refined by including more detailed actuator models and friction effects. Third, I plan to investigate dynamic multi-objective optimization where the robot can re-plan its gait in real time based on sensory feedback. Finally, I would like to combine machine learning techniques, such as reinforcement learning, to automatically discover new gaits from experience.
In conclusion, the multi-objective optimization framework presented here provides a robust and flexible solution to the challenging problem of humanoid robot gait planning. It demonstrates that trade-offs among stability, energy and speed can be effectively managed, thereby enabling more agile and natural locomotion for biped humanoid robots.