Method and system for planning and controlling the balance and locomotion of a humanoid bipedal robot

Abstract

In a method for planning and controlling the balance and locomotion of a humanoid bipedal robot, with at least one control and planning device that performs real-time planning and control of DCM trajectories in 3D using a model of the three-dimensional divergent motion component (3D-DCM), The divergent motion component (3D-DCM) model is designed to include equations of motion for the flight phase.

Description

[0001] The invention relates to a method for planning and controlling the balance and locomotion of a humanoid biped robot according to claim 1, as well as a system for planning and controlling the balance and locomotion of a humanoid biped robot according to claim 11.

[0002] Humanoid robots that walk on two legs are known. Bipedalism, as a form of locomotion, is generally associated with four basic types of movement: standing, walking, running, and jumping. Humans are particularly well adapted to these types of movement, and humanoid robots that walk on two legs are expected to achieve similar or even surpass human performance. Running and jumping, in particular, are locomotion modes that enable legged robots to quickly cover large distances and overcome difficult terrain. Humanoid robots already exist that can run at 10 km / h or even have the ability to traverse parkour-like environments.Despite these successes, however, considerable research efforts are still needed in the field of bipedal locomotion, particularly with regard to the online planning and execution of highly dynamic movements such as running and jumping, crossing uneven and rough terrain, or following a sequence of randomly placed stepping stones.

[0003] In addition to the increasingly popular reinforcement learning approaches [4], [5], a widely used method for implementing bipedal robot locomotion is to focus on the dynamics of the center of mass (CoM). References to publications listed at the end of the description are given below.

[0004] As observed in biological systems [6], CoM dynamics covers the most important aspects of locomotion. Especially in motion planning, whole-body dynamics are often replaced by heuristics [7], simplified models such as the Linear Inverted Pendulum (LIP) [8] or the Spring-Loaded Inverted Pendulum (SLIP) [9], and reduced models such as the three-dimensional Divergent Component of Motion (3D-DCM)

[10] . Of these models, the SLIP model, which is based on biomechanical studies of human walking

[11] , has been widely used to generate walking motions for walking robots

[12] ,

[13] and for bipedal hopping

[14] ,

[15] . However, the dynamics of the SLIP model are nonlinear and cannot be analytically integrated to produce closed-form solutions. Therefore, it is poorly suited to generate walking trajectories for aperiodic locomotion with arbitrary step placement and time parameterization.A method that addresses these limitations was proposed by Englsberger et al. as a Biologically Inspired Deadbeat (BID) controller

[16] , with a more recent work by Egle et al.

[17] implementing the transitions between BID-based running and DCM-based walking.

[0005] An alternative approach is to build on the knowledge and experience gained in bipedal walking research over the past decades. An immediate challenge is that some of the methods developed for bipedal walking, such as the widely used Zero Moment Point (ZMP)

[18] and the Instantaneous Capture Point (iCP)

[19] ,

[20] , are undefined or difficult to apply during the flight phases. The proposed solutions for running [1],

[21] ,

[22] and hopping

[23] combine ZMP-based trajectories during the stance phases with ballistic trajectories during the flight phases. The horizontal and vertical movements are generated independently, and the CoM trajectory is obtained by numerically solving the dynamics.Further development, driven by the requirement to generate vertical CoM movement, led to the VHIP (Variable-Height Inverted Pendulum) model

[24] , an extension of LIP that removes the assumption of a constant CoM height above the ground. The VHIP model has been used in a Model Predictive Controller (MPC) formulation by Ahn and Cho

[25] to generate jumping movements with offline calculated vertical trajectory and by Smaldone et al.

[26] for walking and running movements. An alternative approach to extend the ZMP to three-dimensional space and generate seamless transitions between walking and running was proposed by Sugihara et al.

[27] and used for running

[28] and hopping as a push recovery strategy

[29] .A disadvantage of ZMP-based methods is that the CoM equations of motion are coupled, with the vertical CoM acceleration appearing in the equations of horizontal motion, which is caused by the use of the vertical offset between CoM and ZMP. The 3D-DCM framework

[10] avoids this problem by introducing the Virtual Repellent Point (VRP); the VRP offset above the ground can be used as a free motion parameter instead of the CoM-ZMP offset, resulting in decoupled and identical dynamics on all three axes. The equations of motion can be solved in closed form

[30] , and the 3D-DCM framework has been successfully applied for dynamic walking on mixed terrain

[31] , walking on moving support surfaces

[32] , and multi-contact locomotion

[33] . A time-varying 3D-DCM method was developed by Hopkins et al.

[34] , but the resulting equations of motion can no longer be solved in closed form; a piecewise constant approach for the time-varying DCM with offline calculation of the constant values ​​was recently proposed by Hanasaki et al.

[35] .

[0006] DE102022107587A1 discloses a method for controlling the gait of a robot, comprising the steps of: a) providing a planned task space trajectory, b) providing a planned trajectory of the center of pressure, c) calculating a planned trajectory of the center of mass based on the planned trajectory of the center of pressure, d) receiving a measured trajectory of the robot's center of mass and a measured task space trajectory of the robot, and e) calculating control commands for the robot based on the planned task space trajectory, the measured trajectory of the robot's center of mass, the measured task space trajectory of the robot, the planned trajectory of the center of mass, and learned information about the trajectory of the center of mass angular momentum. Furthermore, a controller, a computer program, and a robot system for implementing this method are disclosed.

[0007] DE102021112485B3 discloses a method for balancing a robot to compensate for external disturbances, the robot comprising: at least one joint with the number of actuated degrees of freedom for the articulated connection of body segments of the robot to at least one actuating device for actively influencing the position of the at least one joint, and a support surface;the method comprising the steps of: a) monitoring and determining the contact moments acting externally on the robot, b) checking whether the determined contact moments exceed at least one predetermined limit value, c) if the at least one predetermined limit value is exceeded in step b) carrying out the steps of: calculating and inducing a required reference center of gravity angular momentum and subsequently reducing the generated reference center of gravity angular momentum by changing the position of at least one joint by means of at least one actuating device of the robot, d) calculating and inducing whole-body movements by means of a whole-body movement optimiser by changing the position of at least one joint by means of at least one actuating device of the robot.;

[0008] EP4321418A1 relates to a method and apparatus for controlling a bipedal robot, and to a readable storage medium. The method includes: detecting a support phase of the bipedal robot in a current state, detecting a current forward velocity and a current lateral velocity of a center of mass of the bipedal robot at an end of the support phase, detecting a torque vector of an ankle joint of a support leg of the bipedal robot according to the current forward velocity and the current lateral velocity, detecting a target torque vector according to the torque vector and a preset quadratic programming model, and controlling the bipedal robot to move according to the target torque vector.EP4164840B1 discloses a method for controlling a kinematically redundant robot to perform multiple tasks, wherein at least one passivity-based first controller module is used, for which at least one first controller module calculates at least one task target description and at least one associated task mapping, for which at least one weighting is calculated, and at least one first controller module is integrated into an overall controller using the at least one weighting, and computer program product comprising instructions which, when the program is executed by means of at least one processor, cause the processor to execute such a method.

[0009] WO2023273107A1 discloses a method and apparatus for controlling the walking of a robot, a robot control apparatus, and a storage medium. The method for controlling the walking of a robot comprises: obtaining a current planned flat-ground walking trajectory of a biped robot and a target inclination angle of a slope on which the biped robot is currently located; adaptively adjusting a planned leg swing trajectory in the planned flat-ground walking trajectory according to the target inclination angle to obtain an expected leg swing trajectory; performing parameter compensation on a planned center of gravity position trajectory in the planned flat-ground walking trajectory according to the target inclination angle to obtain an expected center of gravity position trajectory; and calculating an expected traveling parameter based on a robot algorithm for inverse kinematics.which simultaneously conforms to the expected leg swing trajectory and the expected mass position, and controlling the bipedal robot to move according to the expected driving parameter to enable the bipedal robot to walk on flat ground, so that when the robot encounters a sloping terrain, it runs smoothly on the corresponding slope, thereby improving the walking stability of the bipedal robot.

[0010] US11230001B2 discloses a gait control method for biped robots and a biped robot, the method comprising: obtaining six-dimensional force information and determining a motion state of two legs of the biped robot, calculating a ZMP position of each of the two legs of the biped robot, determining a ZMP expected value of each of the two legs in real time, obtaining a compensation angle of an ankle joint of each of the two legs of the biped robot by inputting the ZMP position, a change rate of the ZMP position, the ZMP expected value, and a change rate of the ZMP expected value to an ankle joint smoothing controller to perform closed-loop ZMP tracking control on each of the two legs, adjusting a current angle of the ankle joint of each of the two legs of the biped robot in real time, and repeating the foregoing steps.

[0011] The object of the present invention is therefore to provide a method and a system for improving the control of humanoid robots that walk on two legs, in particular with regard to the online planning and execution of highly dynamic movements such as running and jumping, crossing uneven and impassable terrain or following a sequence of randomly placed stepping stones.

[0012] The features of claims 1 and 11 serve to solve this problem.

[0013] The invention advantageously provides a method for planning and controlling the balance and locomotion of a humanoid bipedal robot, comprising at least one control and planning device which carries out real-time planning and control of DCM trajectories in 3D by means of a model of the three-dimensional divergent motion component (3D-DCM), wherein the model of the divergent motion component (3D-DCM) contains equations of motion for the flight phase.

[0014] The method can provide that the DCM starting point can be determined as follows: ξφ,0=αφ,ξvφ,0+βφ,ξvφ,T+γφ,ξξφ,T +δφ,ξχφ,0+εφ,ξg. where for the stance phase the coefficients δ φ,ξ and ε φ,ξ are equal to 0, while αφ,ξ =,(1−bTφ+bTφe−Tφb) and βφ,ξ=(bTφ−b+TφTφe−Tφb) and γφ,ξ=e−Tφb and where for the flight phase the coefficients αφ,ξ and b φ,ξ 0 are, while γφ,ξ=e−bb+Tφ and δφ,ξ=Tφb+Tφ and εφ,ξ=bTφ2+2b2Tφ2(b+Tφ) sind.

[0015] Bei dem Verfahren kann preparazent sein, dass die Wegpunkte des Masseschwerpunkts (CoM-Endpunkts) wie folgt bestimmbar ist: xφ,T=αφ,xvφ,0+βφ,xvφ,T+γφ,xξφ,T +δφ,xxφ,0+εφ,xg, bewertung für die Standphase die Koeffizienten d φ,ξ and e φ,ξ are equal to 0, while af,x =,(1−bTφ+bTφe−Tφb) and βφ,ξ=(bTφ−b+TφTφe−Tφb) and γφ,ξ=e−Tφb und behandung für die Flugphase die Koeffizienten a φ,ξ and b φ,ξ 0 are, while γφ,ξ=bb+Tφ δφ,ξ=Tφb+Tφ and εφ,ξ=Tφb+Tφ and εφ,ξ=bTφ2+2b2Tφ2(b+Tφ) correspond.

[0016] The method may include introducing an enhanced centroidal moment pivot point (eCMP) and a virtual repellent point (VRP) into the divergent motion component model (3D-DCM), and using the eCMP and the VRP for real-time planning and control of the DCM trajectories in 3D.

[0017] The method may include a DCM and CoM waypoint calculation Ξ and X in matrix form with respect to the VRP waypoints V, the DCM endpoint ξ f , the CoM starting point x s and the gravitational vector g is as follows: [I−S0TΓξST−S0TΔξS0−STTΓxSTI−STTΔxS0][ΞX]= [S0TAξS0+S0TBξSTSTTAxS0+STTBxST]V+[S0TεξSTTεx]gT +[sf0nω×1]ξfT+[0nω×1ss]xsT, where with V=⌈v1…vnω⌉T∈ℝnω×3 the VRP waypoints, with Ξ=⌈ξ1…ξnω⌉T∈ℝnω×3. the DCM waypoints and with X = [x1 ... x nw ] T ∈ ℝ nw×3the CoM waypoints are shown, where S0 = [I nφ×nφ 0 nφ×1 ] ∈ ℝ nφ×nw selects the starting points of the movement phase, while ST=⌈0nφ×1 Inφ×nφ⌉∈ℝnφ×nω selects the endpoints from an arbitrary waypoint matrix V, Ξ, or X.

[0018] The method can provide for an algorithm for the waypoint calculation with a computational complexity of O(nφ2) is used.

[0019] The method may include using a waypoint calculation algorithm that requires three iterations over all movement phases, wherein in the first two iterations partial solutions for the DCM and CoM waypoints are calculated from the known quantities and the coefficients for the unknown CoM waypoints are tracked, wherein in the third iteration the collected coefficients are used to calculate the complete solution.

[0020] The method may provide for the use of an algorithm I for the waypoint calculation, where the algorithm I looks as follows:

[0021] The procedure may include transition phases between the four types of movement mentioned: standing, walking, running and jumping.

[0022] The method can provide that the control and planning device constructs a sequence of movement phases with alternating stance and flight phases for running, single and double support phases for walking steps, and that gait transition phases are inserted into the phase sequence.

[0023] There may further be provided a system for planning and controlling the balance and locomotion of a humanoid bipedal robot, comprising a humanoid bipedal robot and at least one control and planning device which is designed to carry out a method according to one of claims 1 to 10.

[0024] In the following, embodiments of the present invention are explained in more detail with reference to the drawings.

[0025] They show schematically: Fig. 1 shows the relationships between the CoM and the points defined by the 3D DCM frame, Fig. 2 Overview of the method or movement planner according to the invention, Fig. 3 shows that the 3D-DCM frame is also valid and coherent during the flight phase under the free-falling condition f ext = 0 3×1 Fig. 4 Comparison of the execution times of the two waypoint calculation methods, due to its lower computational complexity O(nφ2) versus O(nφ3) Fig. 5 Example of generating a reference trajectory for two running steps. The gray areas indicate the phases from standstill to flight and from flight to standstill. The dots mark the start and end points of each phase. Fig. 6 Modeling kinematic reach constraints during running. The solid lines illustrate the trajectories during a running step from step p with step vector I using a constant VRP v during the stance phase. The distances between CoM and contact point at touchdown and takeoff are denoted by d td or d to and can be defined as closed functions of the VRP offset above the ground Δz Fig. 7 Example of the transition from standing to walking. The individual movement phases are represented by alternating gray and white backgrounds. The dots indicate the start and end points of the individual phases. Fig. 8 Example of the transition from walking to running. The individual movement phases are represented by alternating gray and white backgrounds. The dots mark the start and end points of the individual phases. Fig. 9 Example of the transition from running to standing. The individual movement phases are represented by alternating gray and white backgrounds, and the dots indicate the start and end points of the individual phases. The beginning and end of the transition from running to standing are highlighted by black vertical lines. Fig. 10 Example of the transition from running to walking. The individual movement phases are represented by alternating gray and white backgrounds, and the dots indicate the start and end points of each phase. The beginning and end of the transition from running to walking are highlighted by black vertical lines. Fig. 11 Jump down a flight of stairs with large steps. Each step is 30 cm long and 40 cm high. Fig. 11 (a) Generated reference CoM trajectories (spatial side view). Solid lines indicate the trajectories during the stance phases; dotted lines are used during the flight phases. The dots mark the location of the start and end points of each movement phase, while the triangles mark the instantaneous variables for the depicted robot. Fig. 11 (b) (b) Required forces on the x- and z-axes and required center of pressure on the x-axis for both feet. Dashed lines indicate the flight phases. Fig. 12 a gait pattern for walking and running, starting with slow walking (0.2 m / s), followed by fast walking (0.4 m / s), slow running (1 m / s) and fast running (2 m / s). Fig. Figure 12 (a) shows simulation images (solid images spaced 1 second apart, translucent images spaced 0.2 s apart during running and 0.5 s apart during walking). 5 s during walking). Fig. 12 (b) Generated reference CoM trajectory (y- and z-axes) showing the start and end points of each movement phase on the y-axis. The gray background represents the dual support phases for walking and the flight phases for running. Fig. 12 (c) Reference and measured data: CoM velocity on the x-axis, center of mass angular momentum (CAM) on all three axes. Fig. (d) Preset forces (z-axis). Fig. 13 Running on a curved track Fig. 13 (a) Simulation images taken at 1 second intervals. Fig. 13 (b) Planned steps and reference CoM trajectory (spatial top-down view). Fig. 13 (c) Planned, reference, and measured waist alignment. The reference trajectory is generated by the angular momentum optimizer. The gray background represents the flight phases. Fig. 14 Walking over a sequence of predefined stepping stones arranged with great variability on all three Cartesian axes. Fig. 14 (a) Images of the simulation. Side view: Fixed images show the robot during the flight phases, translucent images during the stance phases. Top view: Only the flight phases are shown. Fig. 14 (b) Generated reference CoM trajectory. The individual movement phases are shown with alternating gray and white backgrounds. Fig. 14 (c) DCM and CoM tracking errors. The movement phases are labeled with the corresponding contact configuration (DS: double support, LS: left leg support, RS: right leg support, F: flight phase). Fig. 15 BID vs. 3D-DCM comparison of the generated contact forces.

[0026] The present invention describes a method for planning and controlling the balance and locomotion of a humanoid bipedal robot, comprising at least one control and planning device which carries out real-time planning and control of DCM trajectories in 3D by means of a model of the three-dimensional divergent motion component (3D-DCM), wherein the model of the divergent motion component (3D-DCM) contains equations of motion for the flight phase.

[0027] The method is an improvement of the three-dimensional divergent component of motion (3D-DCM) (hereinafter also called three-dimensional divergent component of motion (DCM) or 3D-DCM framework). II Basics

[0028] In the following, the fundamentals of the three-dimensional Divergent Component of Motion (DCM) framework are presented, with a focus on its use for generating reference trajectories. Formally, the motion planning problem for the center of mass x ∈ R 3 of a leg-bearing robot is to find a trajectory that contains the initial state (x s , ẋ s ) with the target state (x f , ẋ f ), taking into account the dynamics of the center of gravity of the second order x¨=1mfext+g,

[0029] Here, m denotes the total mass of the robot, f ext ∈ R 3 is the sum of all external forces acting on the robot and g = (0 0 -g) T is the gravitational acceleration vector.

[0030] The 3D-DCM framework is a reformulation of CoM dynamics (1) without loss of generality. First, the 3D-DCM framework introduces two points separated by a constant vertical offset Δz: the enhanced centroidal moment pivot (eCMP) and the virtual repellent point (VRP). The eCMP e ∈ R 3 encodes the direction and magnitude of the external force f ext above fext=mb2(x−e), where b is a time constant defined as follows b:=Δzg.

[0031] Similarly, the VRP encodes v ∈ R 3 the direction and magnitude of the CoM acceleration ẍ : v=e+(00Δz)T.

[0032] The eCMP point is closely related to the centroidal moment pivot (CMP)

[37] . While the eCMP is a three-dimensional point not necessarily tied to the soil surface, the CMP is located at the intersection of the soil surface with the line connecting the eCMP and the CoM.

[0033] The introduction of the VRP and the VRP offset Δz results in decoupled and identical dynamics on all three Cartesian axes. The VRP offset Δz is a motion parameter that can be chosen arbitrarily as long as it remains strictly positive. For example, if Δz is equal to the vertical distance between CoM and eCMP, Δz = x z - e z , is selected, the vertical CoM acceleration ẍ zto 0, which is typically used when generating walking trajectories on flat terrain. When walking on uneven terrain, Δz can be interpreted as the average CoM height above the ground, assuming that the eCMP is standing on the ground. In this work, we propose an algorithm to calculate Δz for running and jumping such that the generated CoM trajectory satisfies the kinematic reachability conditions (Section IV-G).

[0034] Another point, the three-dimensional divergent motion component (DCM) ξ is defined in

[10] as a linear combination of the CoM position and velocity using the time constant b: ξ:=x+bx˙.

[0035] From (6) and (4) the relationship between DCM and VRP is ξ˙=1b(ξ−v).

[0036] The relationships between the CoM and the points defined by the 3D DCM frame are shown in Fig. 1. By introducing the DCM, the second-order CoM dynamics (4) are split into two first-order dynamics: the stable CoM dynamics relative to the DCM, which results from (6) x˙=−1b(x−ξ), and the unstable DCM dynamics (7) with respect to the VRP. A. Stable DCM dynamics in reverse time

[0037] One of the main ideas presented in

[10] is the generation of the DCM reference trajectory by starting at the end of the movement (t = tf) and calculating the DCM target position ξ (t f ) = ξ f and calculates the DCM trajectory in reverse time. Typically, the target DCM is chosen to be stationary (ξ̇ f = 0), ie, the final DCM and the VRP agree ξ f = v f , according to (7). The most important finding is that the unstable DCM dynamics (7) becomes stable when the time flow is reversed: ξ˙(−t)=−1b(ξ(−t)−v(−t)).

[0038] This property of the DCM dynamics allows the motion planner to generate the CoM trajectory by combining two stable dynamics, the DCM dynamics in backward time and the CoM dynamics in forward time. In the present invention, the reference trajectories generated by integrating the two stable dynamics are referred to as stable trajectories. In the absence of disturbances, the reference trajectories can be tracked by a tracking controller using a pure feedforward approach. However, although the DCM reference trajectory is stable using the backward time method, the actual DCM dynamics remain unstable; any deviation from the reference trajectory requires stabilization, such as that provided by the DCM controller

[10] .

[0039] The time reversal approach offers two major advantages. First, for any VRP reference trajectory, which can be defined as a spatial piecewise linear interpolation over any sequence of Waypoints(vi)i=1n is given, the generated DCM trajectory is stable by design. In fact, the set of all possible DCM trajectories is limited by the convex hull associated with the waypoints v i and ξ f as vertices, as shown in

[30] . Second, the DCM reference trajectory is stable for any duration T of the generated motion, so the duration is an arbitrary parameter of the motion. Lower and upper bounds on the duration are imposed only by the constraints of the physical robot system (kinematic reach limits, friction cone constraints, etc.) and not by trappability concerns

[38] . B. Creation of DCM and CoM reference trajectories

[0040] As shown in

[30] , a complete, C 2 continuous CoM reference trajectory for walking or general multi-contact locomotion by dividing the movement into a sequence of n φPhases are generated. The VRP offset Δz , the VRP trajectory and the phase durations are the free parameters of the complete movement. For each phase φ with duration T φ the reference VRP trajectory v φ (t) as a spatially and temporally linear interpolation between a starting point v φ,0 and an endpoint v φ,T designed: vφ(t)=(1−tTφ)vφ,0+tTφvφ,T, where t ∈ [0,T φ ] is the local time of the motion phase. The DCM trajectory ξ φ (t) is obtained by inserting (10) into (7) and solving the resulting differential equation using the DCM endpoint ξ φ (T φ ) = ξ φ,T as a boundary condition, exploiting the reversal time stability of the DCM dynamics (9): ξφ(t)=(1−t+bTφ+bTφet−Tφb)vφ,0 +(t+bTφ−b+TφTφet−Tφb)vφ,T+et−Tφbξφ,T.

[0041] The DCM starting point ξ φ,0the current movement phase φ is obtained by evaluating (11) for t = 0 : ξφ,0=(1−bTφ+bTφe−tφb)vφ,0+(bTφ−b+TφTφe−Tφb)vφ,T+e−Tφbξφ,T.

[0042] Similarly, the CoM trajectory is derived by substituting (11) into (8) and solving the resulting linear differential equation using the CoM starting point x φ (0) = x φ,0 is solved as a boundary condition: xφ(t)=(1−tTφ+bTφe−Tφbsinh(tb)−e−tb)vφ,0+(tTφ−b+TφTφe−Tφbsinh(tb))vφ,T+e−Tφbsinh(tb)ξφ,T+e−tbxφ,0.

[0043] The CoM endpoint x φ,T of the current movement phase is obtained by evaluating (13) for t=T φ . xφ(t)=(bTφe−Tφbsinh(Tφb)−e−Tφb)vφ,0+(1−b+TφTφe−Tφbsinh(Tφb))vφ,T+e−Tφbsinh(Tφb)ξφ,T+e−Tφbxφ,0.

[0044] The continuity of the complete trajectories is ensured by linking the end points of one movement phase with the starting points of the subsequent phase, ie the following equations apply to all phases: v φ,T = v φ+1,0 , ξ φ,T = ξ φ+1,0 , and x φ,T = x φ+1,0 . The complete reference trajectories for VRP, DCM and CoM can be calculated as piecewise interpolations over sequences of n w = n φ+ 1 waypoints using (10), (11), and (13), respectively, as interpolation functions. Given the VRP waypoints, the DCM and CoM waypoints can be calculated either sequentially or in matrix form using (12) and (14), as shown in

[30] . Once the waypoints are known, the instantaneous VRP, DCM, and CoM positions are calculated using the respective interpolation functions for the current phase φ and the current time t. Finally, the CoM velocity and CoM acceleration are determined by evaluating (8) and (4), respectively. III. OVERVIEW

[0045] An overview of the method and motion planner according to the invention is given in Fig. 2. The input consists of a sequence of steps specified by a step planner, the durations of the movement phases, and the maximum distance between CoM and the contact point during touchdown and takeoff, which is used to calculate the VRP offset Δz. The movement planner constructs a sequence of movement phases with alternating stance and flight phases for running, single and double support phases for walking steps, and inserts gait transition phases into the phase sequence as described in Section V. Based on the movement phase sequence, the instantaneous CoM, the center of mass angular momentum (CAM), and the foot reference trajectories are calculated synchronously with the execution rate of the whole-body controller. Finally, the waist orientation and the reference trajectories of the upper body are generated using a whole-body optimization approach

[39] to determine the whole-body reference trajectory x taskThis is tracked by the passivity-based whole-body torque controller described in

[31] . IV. RUNNING AND JUMPING

[0046] Running and jumping are locomotion patterns characterized by alternating stance and flight phases; in running, the stance leg is changed with each step, whereas in jumping, both legs are used simultaneously for takeoff and landing. In the following, we show how the equations of motion for the flight phase can be described in a way that is compatible with the 3D DCM framework. Furthermore, we establish waypoint equations similar to (12) and (14) for the flight phase, and we present an algorithm for calculating DCM and CoM waypoints for any sequence of movement phases, including flight phases. A. Stance phase

[0047] In bipedal walking, the reference VRP trajectory during the stance phase can be designed either as a fixed point above the stance foot center p with the vertical offset Δz vφ(t)=p+(0 0 Δz)T, or as a heel-toe movement by positioning the VRP start and end points, v φ,0 and v φ,T , with corresponding offsets relative to the foot center. The DCM and CoM reference trajectories are generated using the method presented in Section II-C. B Flight phase

[0048] During the flight phase, the CoM is influenced only by gravity and inevitably follows a parabolic trajectory. The CoM's equations of motion during the free-fall phase can be described as follows: xφ(t)=xφ,0+tx˙φ,0+t22g, x˙φ(t)=x˙φ,0+tg.

[0049] During the flight phase, the VRP trajectory can no longer be freely chosen, as was the case during the stance phase. Nevertheless, the VRP trajectory results from (4), where the free fall condition ẍ = g and the definition of b as vφ(t)xφ(t)−b2g=xφ(t)+(0 0 Δz)T, which shows that the VRP trajectory follows the parabolic CoM trajectory with a constant vertical offset equal to Δz. This is equivalent to the observation that during the flight phase the eCMP trajectory e φ (t) with the CoM trajectory x φ (t) is identical, which can also be confirmed by applying the free fall condition f ext = 0 3×1 on (2). This shows that the 3D-DCM framework is valid and coherent even during the flight phase (see Fig. 3). The DCM trajectory is obtained by multiplying (17) by the time constant b, adding (16) and applying the DCM definition (6) for ξ φ (t) and ξ φ,0. The DCM trajectory with respect to the starting points CoM x φ,0 and DCM ξ φ,0 is ξφ(t)=ξφ,0+tb(ξφ,0−xφ,0)+t2+2tb2g.

[0050] To apply the time reversal approach described above, the DCM starting point ξ φ,0 with respect to the endpoint ξ φ,T ξ φ,T For this purpose, we evaluate (19) for t = T φ and rearrange the terms: ξφ,0=bb+Tφξφ,T+Tφb+Tφxφ,0−bTφ2+2b2Tφ2(b+Tφ)g.

[0051] If the total movement includes a flight phase, the DCM waypoints can no longer be calculated independently of the CoM waypoints, as was the case in

[30] . Here, in (20), the DCM starting point ξ φ,0 the flight phase from CoM starting point x φ.0 However, the DCM and CoM waypoints can be calculated jointly in matrix form or using an algorithmic approach, as shown in the following section. The CoM endpoint xφ,T of the flight phase is obtained by applying (16) for t = T φ evaluates the DCM definition (6) for ξ φ,0 used to calculate the velocity term ẋ φ,0 . to eliminate, and ξ φ,0 replaced by the result from (20): xφ,T=Tφb+Tφξφ,T+bb+Tφxφ,0−bTφ22(b+Tφ)g. C. Unified formulation for waypoint equations

[0052] Recall that the complete movement plan consists of a sequence of heterogeneous movement phases, including stance, flight, and other transition phases, which we introduce in the following sections of this article. To develop a general method for waypoint calculation, we propose a unified formulation of the DCM and CoM waypoint equations. This approach makes the framework extensible with new movement phase types. Combining (12) and (20) in one equation yields the general form for the DCM starting point: ξφ,0=αφ,ξvφ,0+βφ,ξvφ,T+γφ,ξξφ,T +δφ,εxφ,0+εφ,εg.

[0053] For the stance phase, the coefficients δ φ,ξ and ε φ,ξ equal to 0, while α φ,ξ , β φ,ξ and γ φ,ξ correspond to the coefficients in equation (12). Accordingly, for the flight phase, the coefficients α φ,ξ and β φ,ξ are 0, while γ φ,ξ , δ φ,ξ , and ε φ,ξ correspond to the respective coefficients in equation (20). Similarly, combining (14) and (21) in one equation yields the general form for calculating the CoM endpoint: xφ,T=αφ,xvφ,0+βφ,xvφ,T+γφ,xξp,T +δφ,xxφ,0+εφ,xg, with the same considerations regarding the coefficients as discussed above for the equation. D. Waypoint calculation in matrix form

[0054] The complete trajectories for VRP, DCM and CoM can be defined as n φ piecewise interpolations over nw =n φ + 1 waypoints are described. The VRP waypoints are described in V=⌈v1…vnω⌉T∈ℝnω×3, the DCM waypoints Ξ=[ξ1...ξnω]T∈ℝnω×3 and the CoM waypoints X=[x1...xnω]T∈ℝnω×3 The objective is to determine the unknown DCM and CoM waypoints Ξ and X with respect to the VRP waypoints V, the DCM endpoint ξ f ,, the CoM starting point x s and to calculate the gravitational vector g. For this purpose, the derived equations (22) and (23) are used. Since these equations are written in the form of start and end points (v φ,0 , v φ,T , etc.), two selection matrices are introduced: S0 = [I nφ×nφ 0 nφ×1 ] ∈ ℝ nφ×nw selects the starting points of the movement phase, while ST=[0nφ×1 Inφ×nφ]∈ℝnφ×nω selects the endpoints from an arbitrary waypoint matrix V, Ξ, or X. In (22) in matrix form for all n φPhases written as S0Ξ=AξS0V+BξSTV+ΓξSTΞ +ΔξS0X+εξgT, where A ξ , B ξ , Γ ξ and Δξ quadratic (n φ × n φ ), are diagonal matrices containing the coefficients α φ,ξ , β̇ φ,ξ , γ φ,ξ, or δ φ,ξ contain, while ε ξ is a vector containing the coefficients ε φ,ξ contains. (24) consists of n φ equations, while Ξ n w waypoints. To bring (24) into a quadratic form, the DCM terminal constraint ξ nw = ξ f added, which can be expressed in the form of Ξ as [0nφ×nφ0nφ×101×nφ1]︸Sn∈ℝnω×nωΞ=[0nφ×1 1]︸sf∈ℝnωξfT.

[0055] If you multiply (24) on the left side by S0T, adding (25) and grouping the terms Ξ and X on the left side, we get (I−S0TΓξST)Ξ−S0TΔξS0X S0T(AξS0+BξST)V+S0TεξgT+sfξfT, where the equality S0TS0+Sn=Inω×nω=I The CoM waypoints are written (23) in matrix form for all n φ Phases as STX=AxS0V+BxSTV+ΓxSTΞ +ΔxS0X+εxgT, the coefficients α φ,x , β φ,x , γ φ,x , δ φ,x , ε φ,x to the square, diagonal matrices A x , B x , Γ x , Δ x or the vector ε x If you multiply (27) on the left side STT and adds the initial CoM constraint at the beginning of the movement, we get x1 = x s , expressed with respect to X as [101×nφ0nφ×10nφ×nφ]︸S1∈ℝnω×nωX=[ 10nφ×]︸ss∈ℝnωxsT,

[0056] yield (I−STTΔxS0)X−STTΓxSTΞ= STT(AxS0+BxST)V+STTεxgT+ssxsT, where the equality STTST+S1=I Finally, stacking (26) and (29) into one equation leads to the linear system [I−S0TΓξST−S0TΔξS0−STTΓxSTI−STTΔxS0][ΞX]= [S0TAξS0+S0TBξSTSTTAxS0+STTBxST]V+[S0TεξSTTεx]gT +[sf0nω×1]ξfT+[0nω×1 ss]xsT, which can be solved to find the unknowns Ξ and X. E Highly efficient waypoint calculation using an algorithmic approach

[0057] The calculation of the waypoints in matrix form, presented above as equation (30), is a compact and mathematically appealing method for calculating the DCM and CoM waypoints. However, the computational complexity of the algorithm is O(nφ3) due to the matrix multiplications involved and the requirement to create a linear system with 2 nφvariables. From a practical point of view, the high computational complexity is disadvantageous when implementing the algorithm on real-time robot hardware, as it may impose an upper limit on the number of phases n φ Motivated by this inadequacy of the matrix form algorithm, an algorithm for the waypoint calculation with a computational complexity of O(nφ2) provided as Algorithm 1.

[0058] The algorithm requires three iterations across all motion phases. The main idea is to compute partial solutions for the DCM and CoM waypoints from the known quantities in the first two iterations and to track the coefficients for the unknown CoM waypoints. In the third iteration, the collected coefficients are used to compute the complete solution. The input to the algorithm consists of the VRP waypoint matrix V, the DCM endpoint ξ f, the CoM starting point x s , the phase durations, which are listed as T=(Ti)i=1nφ and the DCM time constant b. The term v i denotes the i-th waypoint in the matrix V, with analogous notation for the DCM and CoM waypoint matrices Ξ and X. The algorithm uses two quadratic (n w × n w ) Matrices, C ξ and C x to store the coefficients corresponding to the unknown CoM waypoints. The value in the i-th row and j-th column of the matrix C ξ , written as C ξ [i, j], denotes the coefficient with which the CoM waypoint x j when calculating the DCM waypoint ξ i must be multiplied (Algorithm 1, line 31). The matrix C x is defined in a similar way, where C x [i, j] the CoM waypoint x j with x i(line 27). The partial DCM and CoM waypoints are stored in the Ξp=[ξ1p...ξnωp]T∈ℝnω×3 or Xp=[x1p...xnωp]T∈ℝnω×3 The following describes in detail how the formulas for calculating the partial DCM and CoM waypoints and the coefficients C ξ and C x can be obtained from the previously derived equations. We begin by rewriting (22) in waypoint notation as ξi=αξυi+βξυi+1+γξξi+1+δξxi+εξg.

[0059] It can be proven by induction that the DCM waypoint ξ i can be written as follows ξi=ξip+∑j=inωCξ[i,j]xj.

[0060] For i = n w Equation (32) applies due to the initialization of C ξ and ξnωp in lines 2 and 3. Assuming that equation (32) holds for i +1, ξ i+1in (31) is replaced by the corresponding form from (32) and the terms are rearranged to obtain ξi=αξυi+βξυi+1+γξξi+1p+εξg︷ξip +δξxi+∑j=i+1nωγξCξ[i+1,j]xj.

[0061] It takes advantage of the fact that the first four terms are independent of the unknown CoM waypoints X and p can be summed to obtain the partial DCM waypoint ξip (line 6). The last two terms in (33) yield the corresponding equations for calculating the DCM coefficients C ξ (Lines 7 to 10).

[0062] For the CoM waypoints, the procedure is analogous by rewriting (23) using the waypoint notation as follows xi=αxυi−1+βxυi+γxξi+δxxi−1+εxg, and shows that it can also be written in the following form: xi=xip+∑j=i+1nωCx[i,j]xj.

[0063] The proof is again by induction. For i = 1, equation (35) holds, ie x1 = x s , due to the initialization of C x and x1p in lines 13 and 14. Assuming that equation (35) holds for i - 1, we replace in (34) the DCM waypoint ξ i by the result of (32) and x i -1 by the corresponding equation from (35): xi=αxυi−1+βxυi+γxξip+∑j=inωγxCξ[i,j]xj +δxxi−1p+∑j=inωδxCx[i−1,j]xj+εxg. x i appears on both sides of the equation (on the right-hand side, as part of the sum terms containing x j for j = i). If one removes all terms that contain x i on the left side and rearranging the terms, one can bring (36) into agreement with (35): (1−γxCξ[i,i]−δxCx[i−1,i])︸ηxi= αxυi−1+βxυi+γxξip+δxxi−1p+εxg +∑j=i+1nω(γxCξ[i,j]+δxCx[i−1,j])xj.

[0064] The sum of the first five terms on the right side of the equation divided by η represents the partial CoM waypoint xip (line 18), while the last term is the formula for calculating C x [i, j] (lines 19 to 21). To obtain the complete solution, the algorithm exploits the structure of equations (32) and (35) by iterating backward over all motion phases. During each iteration cycle, the CoM waypoint x i calculated from the known future waypoints.

[0065] x i+1 up to x nw (Lines 25 to 28). As soon as x i is known, the DCM waypoint ξ i be calculated immediately, since it only depends on the CoM waypoints x i up to x nw(lines 29 to 32). To verify the performance of the two waypoint calculation methods, Algorithm 1 and the matrix form (30), both were implemented in C++ and tested on the real-time computer of the humanoid robot TORO. The average execution time for a different number of phases is shown in Fig. 4: The algorithmic approach clearly outperforms the matrix form in all cases due to its lower computational complexity F. Movement phases from flight to stand and from stand to flight

[0066] Discontinuities in the generated reference trajectories should be avoided, as these are undesirable for two reasons. First, a real robot has a limited control bandwidth, meaning that in the event of a discontinuity, the controller is unable to follow the reference trajectory. Second, discontinuities tend to excite unmodeled joint and link elasticities on the real robot, further degrading the controller's tracking performance. The continuity of the DCM and CoM reference trajectories is ensured by linking the start and end points of adjacent phases. In contrast, the VRP and eCMP trajectories presented so far are discontinuous at the beginning and end of the flight phase, leading to discontinuities in the generated forces.For example, the eCMP trajectory during the flight phase is identical to the CoM trajectory, whereas during the stance phase it is conceived as a fixed point on the ground, coinciding with the center of the foot. To eliminate these discontinuities, we introduce two new movement phases: a flight-to-stance and a stance-to-flight phase. For these phases, we attempt to determine the DCM starting point ξ. φ,0 and the CoM endpoint x φ,T in the general form of (22) or (23), so that we can use the waypoint calculation algorithm without changes. During both movement phases, the VRP is calculated from a starting point v φ,0 to an endpoint v φ,Tinterpolated using the interpolation function (10). However, unlike in the stance phase, either the start or the end point is restricted by the free fall condition ẍ = g, the effects of which are explained in more detail below. Flight-to-stance phase: For the flight-to-stance phase, the VRP start point v φ,0 be written as υφ,0=xφ,0−b2g, which corresponds to the previous flight phase. This means that the VRP starting point v φ,0 is an unknown quantity that is determined from the CoM starting point x φ,0 which still has to be calculated using the waypoint calculation algorithm. Inserting (38) into (12) gives ξφ,0=βφ,ξυφ,T+γφ,ξξφ,T+αφ,ξxφ,0−αφ,ξb2g, which corresponds to the general form (22) with α φ,ξ α φ,ξ = 0 and ε φ,ξ = -α φ,ξ b 2 Similarly, inserting (38) into (14) yields xφ,T=βφ,xυφ,T+γφ,xξφ,T +(αφ,x+δφ,x)xφ,0−αφ,xb2g, which corresponds to the general shape of the CoM endpoint (23).

[0067] Stand-to-Flight Phase: For the stand-to-flight phase, the VRP endpoint constraint can be φ,T be written as follows υφ,T=xφ,T−b2g.

[0068] Replacing (41) in (14) and solving for x φ,T leads to a CoM endpoint equation that corresponds to the general form: xφ,T=αφ,x1−βφ,xυφ,0+γφ,x1−βφ,xξφ,T +δφ,x1−βφ,xxφ,0−βφ,xb21−βφ,xg.

[0069] It can easily be verified that β φ,x , the coefficient, the v φ,T in (14), the inequality 0 < β φ,x < 1 for all possible durations T φ > 0 and time constants b > 0. This guarantees that the denominator 1 - β φ,x can never be equal to 0. If we replace (41) in (12) and set x φ,Twith the result from (42), one obtains a DCM starting point equation which corresponds to the general form: ξφ,0=(αφ,ξ+βφ,ξαφ,x1−βφ,x)υφ,0+(γφ,ξ+βφ,ξγφ,x1−βφ,x)ξφ,Tβφ,ξδφ,x1−βφ,xxφ,0−βφ,ξb21−βφ,xg.

[0070] The generated trajectories for two running steps are shown in Fig. 5 and demonstrate all features of the 3D-DCM framework discussed so far. All trajectories are continuous after the introduction of the flight-to-stand and stand-to-flight phases. The eCMP and the VRP are identical on the x- and y-axes and separated by a fixed vertical offset on the z-axis throughout the movement. This explains why the eCMP in Fig. 5 is not visible in the x- and y-plots. During the flight phase, the CoM and the eCMP are identical, while on the x- and y-axes, the offset between the CoM and the DCM is constant, corresponding to a constant CoM velocity. G. Kinematic accessibility restriction

[0071] In general, when walking with humanoid robots, the limits of kinematic reach are avoided by using a constant CoM height above the ground and a limit on the stride length. When walking with longer strides, switching to a toe contact at the end of the stance phase and using a heel contact at touchdown are common solutions to avoid knee singularities and kinematic reach limits. In contrast, in humanoid walking, longer stride lengths are typically used, while the presence of the flight phases induces a significant vertical movement of the DCM and the CoM. In particular, the motion parameter Δz, which indicates the vertical offset between VRP and eCMP, can no longer be chosen and interpreted as the average height of the CoM above the ground, as in bipedal walking. As in Fig. 5 on the z-axis, the VRP during running is significantly lower than the CoM and DCM trajectories during stance phases; in this example, a VRP offset of 0.476 m results in a CoM trajectory with an average height above the ground of 0.87 m. Therefore, the value of Δz during bipedal running and jumping must be chosen so that the kinematic reach constraints during stance phases are not violated. In this section, we propose a method to calculate Δz by specifying the maximum distance between CoM and contact point during touchdown and take-off (see Fig. 6). The choice of this metric is motivated by the fact that the motion planner generates trajectories for a reduced center-of-gravity model without a kinematic body.

[0072] Starting with a step p and a step vector I, we construct a walking cycle consisting of four phases: flight to stance, stance, stance to flight, and flight. During the stance phase, the VRP v is chosen to occupy a fixed position with the vertical offset Δz above the step p. The coordinate system is defined such that p x and p z both are 0, p y corresponds to the offset of the footstep relative to the sagittal plane of the robot and the x- and y-components of the step vector I are positive. The DCM waypoints Ξ=[ξ1...ξ5]T and the CoM waypoints X = [x1 ... x5] T can be determined as closed functions of Δz by solving equation (30) with the following boundary conditions: ξf=ξ5=Trξ1+l, xs=x1=Tr(x5−l), where Tr=[1000−10001].

[0073] Note that both the DCM and CoM waypoints change sign on the y-axis; this results from the alternating sequence of left and right steps when running. When jumping, the procedure is similar, with two differences regarding the y-axis: The VRP is placed in the center of the support polygon, i.e., v y = 0, and the transformation matrix Tr in (44) and (45) is replaced by the identity matrix, since there are no alternating steps when jumping.

[0074] The distance between the CoM waypoint at the touchdown point x1 and the contact point p can be described as follows dtd(Δz)=‖x1(Δz)−p‖ ie as a closed function of Δz . Therefore, we can also take the derivative of d td calculate for Δz in closed form and d td = d max solve iteratively using the Newton method, where d maxis a configuration parameter that specifies the kinematic reachability limit. Similarly, we proceed for the distance between the CoM waypoint and the contact point at the start dto(Δz)=‖x4(Δz)−p‖.

[0075] Finally, the vertical offset Δz is chosen such that the kinematic reachability limit is maintained for both waypoints.

[0076] If the duration of the flight to standstill and the standstill to flight phases are equal, then the CoM movement is symmetrical with respect to the contact point p, and the two distances mentioned above are identical d td = d to V. GAIT TRANSITIONS

[0077] Motivated by the goal of generating continuous reference trajectories, explicit transition phases are provided for connecting the four aforementioned modes: standing, walking, running, and jumping. This section discusses the transitions between standing, walking, and running in detail. They can be applied to jumping with minor modifications, since the differences between running and jumping are only manifested on the y-axis. A. Standing-to-Walking

[0078] In the present invention, the standing state is characterized by a stationary DCM (ξ̇ = 0), which according to (7) means that the DCM is equal to the VRP. In contrast to the commonly used definition for the standing state, the CoM does not have to be stationary, since its stable dynamics means that the CoM converges asymptotically to the DCM position at the end of the movement sequence. The transition from standing to walking is realized by two movement phases: The first phase, the so-called standing-to-moving phase, ensures the continuity of the DCM trajectory during the transition; the second phase, the so-called height-change movement phase, lowers the VRP from the initial CoM height to the relative height above the ground calculated for the walking gait using the kinematic reachability constraint. For both phases, our goal is to determine the DCM starting point ξ φ,0 and the CoM endpoint x φ,Tin their respective general form as equations (22) and (23), which, as shown in Section IV-D, is a prerequisite for applying the waypoint calculation algorithm.

[0079] Stand-to-move phase: During the stand-to-move phase, the VRP is shifted by a calculated amount Δv, which ensures that for a given DCM endpoint ξ φ,T the initial DCM position ξ φ,0 with the initial VRP v φ,0 coincides (see Fig. 7). At the end of the stand-to-move phase, the VRP returns to its starting position, i.e. v φ,T = v φ,0 This phase has two parameters: the duration T shift the VRP shift that occurs at the beginning and end of the phase, and the time T hold , the VRP remains at a standstill; the total duration of the phase is T φ = T hold + 2 T shift. The described VRP movement can also be interpreted as a sequence of three sub-phases: a shift of v φ,0 to an as yet undetermined VRP position φ,1 = v φ,0 + Δv, a holding subphase at v φ,1 with duration T hold and a final displacement back to the starting position. Writing the DCM waypoint equations for the three subphases using (12), the stationary condition at the beginning of the phase (ξ φ,0 = v φ,0 ) and solve for Δv, we obtain the closed form solution Δυ=Tshiftbe−Tφ / b(1−e−Tshift / b)(1−e−(Tshift+Thold) / b)(υφ,0−ξφ,T).

[0080] The CoM endpoint equation for the stand-to-move phase is obtained by rewriting the CoM waypoint equations for the three subphases as (14), combining them with the DCM waypoint equations, and replacing Δv with the result from (48). Remarkably, the resulting equation simplifies to xφ,T=1−e−Tφb2υφ,0+1−e−Tφb2ξφ,T+e−Tφbxφ,0

[0081] In comparison to the transition from standing to walking presented in

[40] , in which an additional VRP waypoint is inserted into the VRP matrix V and obtained by inverting a square matrix (n w × n w), the standing-to-walking phase proposed in this work has two advantages. First, at the end of the phase, the VRP returns to its starting position, with V remaining completely determined, which allows the use of the efficient waypoint calculation algorithm presented above. Furthermore, the calculation of Δv requires only local information such as the VRP starting position v φ,0 , the DCM endpoint ξ φ,T and the durations T shift and T hold Secondly, the proposed stand-to-move phase requires the same phase duration T φ and the same DCM endpoint ξ φ,T a significantly smaller shift of the VRP position than the method from

[40] . Fig. In the example shown in Figure 7, the VRP displacement is reduced by a factor of approximately 2.2. Reducing the VRP displacement increases the feasibility and robustness of the overall movement, as the VRP remains closer to the center of the support surface in the xy plane during this phase. Height change phase: The height change phase is characterized by the vertical movement of the VRP, which can be combined with a stationary VRP or a linearly interpolated movement on the x- and y-axes. The VRP trajectory can be designed independently on the three Cartesian axes. During the transition from standing to walking, the VRP remains stationary on the x-axis, while moving on the y-axis from the center of the support polygon, which corresponds to the double support posture, to the leg that acts as the first stance during the walking movement (see Fig. ). The transitions from running to standing and from walking to running also involve a height change phase, but the VRP movement on the y-axis is different from the one used here, as explained in the relevant sections.

[0082] During the height change phase, the VRP must be lowered along the z-axis from the initial height corresponding to the standing posture to the relative height above the ground calculated for the walking gait cycle. At the same time, the eCMP should maintain its position along the z-axis (see Fig. 7), so that the ground reaction forces are concentrated at the contact point. This relative movement leads to a change in the offset Δz between VRP and eCMP, which corresponds to a change in the DCM time factor b. However, in all previous equations, b was chosen to be constant, since DCM dynamics does not have a closed-form solution for a time-varying b. Therefore, we propose to implement the height change phase with a piecewise constant approach for Δz and b, where the discretization time interval Δt can be chosen arbitrarily small, e.g., equal to the time interval of the controller execution.

[0083] First, we investigate the effects of an instantaneous change in the DCM time constant b on the DCM and VRP positions. For the same CoM position x and velocity ẋ, two DCM positions with different time constants b and b' can be written as DCM definition (6) and ξ'=x+b'x˙.

[0084] If x' in (50) is replaced by the corresponding value from (6), the DCM conversion formula is obtained ξ'=(1−b'b)x+b'bξ.

[0085] Based on the VRP definition (4) and the equation of ẍ in both equations, the VRP conversion formula is as follows υ'=(1−b2b2)x+b'2b2υ, where v' denotes the VRP corresponding to the DCM time constant b'. Intuitively, the conversion formulas can be understood as a decomposition of DCM and VRP into their components, followed by reassembly using a different time constant. In this case, the DCM and VRP calculation becomes dependent on the instantaneous CoM position, while the DCM and VRP trajectories become discontinuous at the moment the time constant changes. While the DCM and VRP trajectories are discontinuous, the generated CoM trajectory remains smooth, meaning the CoM position, velocity, and acceleration are continuous throughout the entire motion.

[0086] In the present embodiment, which consists of several movement phases, the change in the DCM time constant is made at the phase boundaries. For a phase φ, we denote υφ,0' and ξφ,0' the VRP and DCM waypoints v φ,0 and ξ φ,0after converting the local DCM time constant b φ in the previous phase b φ-1 . Neighbouring phases can then be υφ−1,T=υφ,0' and ξφ−1,T=ξφ,0' and apply the algorithm for waypoint calculation. The application of (51) and (52) with b = b φ and b' = b φ-1 to the DCM waypoint equation (22) leads to the following result: ξφ,0'=αφ,ξrφυφ,0'+rφβφ,ξυφ,T+rφγφ,ξξφ,T+(1−rφ(1−αφ,ξ−δφ,ξ)−αφ,ξrφ)xφ,0+rφεφ,ξg where r φ = b φ-1 / b φ gives the relationship between the DCM time constants. The equation for the CoM endpoint is obtained similarly by applying (52) to (23), which leads to the following results xφ,T=αφ,xrφ2υφ,0'+βφ,xυφ,T+γφ,xξφ,T +(αφ,x+δφ,x−αφ,xrφ2)xφ,0+εφ,xg.

[0087] The waypoint equations (53) and (54) can also be used for a sequence of running steps with different time parameterization for each step. According to the approach proposed in Section IV-G, a unique offset Δz is calculated for each step based on the parameters stance and flight duration. The instantaneous change of the DCM time constant can be most easily implemented at the beginning of the flight phase, since the waypoint equations (53) and (54) change due to the fact that the coefficients α φ,ξ , β φ,ξ , α φ,x , and β φ,x for the flight phase are 0. Using the equations derived above, the altitude change phase is defined as a sequence of n h = T φ / Δt subphases implemented, where T φ is the total duration of the phase and Δt is the chosen discretization interval. During each subphase i ∈ {1, ... ,n h} is the vertical offset Δzφ,i constant and equal Δzφ,i=Δzφ,0+inh(Δzφ,T−Δzφ,0), where Δz φ,0 and Δz φ,T denotes the offset at the beginning or end of the altitude change phase. The VRP movement can be described as a linear interpolation over a sequence of VRP waypoints v φ,i be described in relation to the start and end points of φ,0 and v φ,T can be written as follows υφ,i=nh−inhυφ,0+inhυφ,T.

[0088] By applying (53) and (54) to each subphase, the DCM starting point coefficients for the height change phase, α φξ to ε φ,ξ , and the CoM endpoint coefficients, α φ,x to ε φ,x , can be calculated using an algorithm similar to Algorithm 1. The only difference is that instead of DCM and CoM waypoints, the algorithm uses the waypoint coefficients during the backward and forward iterations over the n hThe computational complexity of the algorithm is O(nh2). B. From walking to running

[0089] On the z-axis, the transition from walking to running ( Fig. 8) is realized in a similar way to the transition from standing to walking, whereby the VRP trajectory consists of a stance-to-movement and an elevation change phase. On the x- and y-axes, the VRP is stationary because the transition from walking to running is intended to begin after a double support phase, so that the subsequent single support phase of walking transitions seamlessly into the first stance phase of walking. However, it should be noted that despite its stationary nature, the VRP trajectory consists of three different phases on the x- and y-axes: the single support phase of walking, the elevation change phase, and the standing running phase. In each phase, the DCM and CoM waypoints are calculated using different equations or different DCM time constants b. C. Running to standing

[0090] The transition from running to standing consists of a single elevation change phase in which the VRP is raised from the height of running to the level corresponding to standing (see Fig. 9). On the x-axis, the VRP remains stationary, while on the y-axis it moves from the last step position to the center of the support polygon, where it remains stationary during the standing state. Note that the DCM is stationary at the end of the walking-to-standing transition, while the CoM converges asymptotically to the last common waypoint of the VRP and DCM. D. From running to walking

[0091] The transition from running to walking is implemented similarly to the transition from running to standing and from walking to running. On the z-axis, the VRP trajectory is the same as for the transition from running to standing, while on the x- and y-axes, the last stance phase of the running gait transitions seamlessly into the first single support phase of the walking gait (see Fig. 10). VI. GENERATING AND CONTROLLING WHOLE-BODY MOVEMENTS

[0092] The proposed CoM reference trajectory planner is integrated with the whole-body motion generation and control methods presented in our previous work. The swing leg reference trajectory, including the use of edge contacts (toe-off), is described in

[31] , while the angular momentum generation method from

[39] is used to compensate for the angular momenta induced by the leg movement during running and jumping. The generated reference trajectories are tracked by the passivity-based whole-body controller presented in detail in

[31] and

[41] . In this section, we provide a brief overview of the angular momentum generation method and describe the extension introduced in this work to handle the flight phase. A. Dynamic model

[0093] Following

[42] , the dynamics of the humanoid robot is represented by a floating base model with the CoM position x ∈ ℝ 3 and waist orientation Rb∈SO(3) as base coordinates. The corresponding linear and angular velocities, ẋ ∈ ℝ 3 and ω b ∈ ℝ 3 , become the velocity vector v=(x˙TωbT)T For a humanoid robot with n actuated ωb joints, the equations of motion can be written as follows M(v˙q¨)+C(vq˙)+(−ωg0n×1)=(06×1 τ)+τ¯ext, where M and C represent the inertia and Coriolis matrices, respectively, wg=(mgT 03×1T)T∈ℝ6 the gravitational key, q ∈ ℝ n the joint positions, τ ∈ ℝ n the vector of actuator moments and τ¯ext∈ℝ6+n the generalized external forces. The task space consists of the CoM position, the waist

[0094] Orientation, the Cartesian position and orientation of the feet and the joint positions of the upper body

[31] . The velocity vector of the task can be written as (vvRvLq˙u)︸x˙task=[I6×606×nAdRJRAdLJL0nu×6Su]︸J(vq˙) where v R ∈ ℝ 6 the velocity vector that combines translational and rotational velocities, Ad R ∈ ℝ 6x6 is the adjoint matrix and J R ∈ ℝ 6×n the Jacobian matrix for the right foot. Analogous quantities for the left foot are denoted by the index "L". The velocities of the n u Upper body joints are supported by q̇ u ∈ ℝ nu where S u ∈ ℝ nu×n selects the corresponding joints from the total vector q̇. B. Angular momentum compensation

[0095] The angular momentum of the center of mass (CAM) l c ∈ ℝ 3, which represents the angular momentum around the center of gravity, depends linearly on the velocity vector lc=Mω(vq˙)=MωJ#︸A¯x˙task, where M ω the rotation part (i.e., rows 4 to 6) of the inertia matrix M is

[43] . Here, J# is a generalized damped pseudo-inverse of the task-space Jacobians that is robust to singular configurations J#=JT(JJT+ρ2I)−1 with ρ ∈ ℝ as the damping factor. Whole-body motion optimization

[39] uses a subset of the degrees of freedom (DoF) in the task space to actively generate angular momentum with the goal of providing the system with a reference CAM lcref In this work, we choose the waist orientation R b and the upper body joints q u with the exception of the wrists, since their contribution to the total angular momentum is negligible. a denotes the velocity vector A athe matrix resulting from selecting the corresponding rows from A for this DoF. Starting from a reference velocity in the task space x˙taskref and a reference angular momentum lcref The motion optimization finds optimal speeds x˙aopt, which fulfill the angular momentum task while minimizing the deviation from the reference quantities of the task (for details see

[39] ).

[0096] For running and jumping, we generate the reference CAM individually for the different movement phases presented in this work. During the phases in which the robot touches the ground, such as the stance or transition phases, lcref chosen so that lcref=A¯ax˙aref.

[0097] Intuitively, equation (61) states that the angular momentum induced by the selected DoF is required for the planned movement and should not be compensated by the motion optimizer. More precisely, the reference KAM of the overall system is given only by the selected DoF, while the angular momentum induced by the remaining DoF (i.e., the swing leg movement) is undesirable and must be compensated. In many scenarios, the reference velocities in x˙aref equal to 0, so that the reference angular momentum 0(lcref=03×1) However, this is not always the case, for example when running on a curved track (see Fig. 13).

[0098] During the flight phase, the total angular momentum is constant. To generate consistent trajectories for the whole-body motion, we choose the reference angular momentum equal to the measured angular momentum. lcref=lc

[0099] During the stand-to-flight and flight-to-stand phases, the reference angular momentum is linearly interpolated between the two reference generation methods. In the stand-to-flight phase, lcref for example, created as follows lcref=(1−tTφ)A¯ax˙aref+tTφlc, where T φ the duration of the stand-to-flight phase and t ∈ [0, T φ ] is the local phase time. VII. SIMULATIONS

[0100] We performed extensive simulations of the proposed CoM trajectory generator with the passivity-based whole-body controller in OpenHRP

[44] using the torque-controlled humanoid robot TORO

[36] , a 27-DoF robot with a height of 1.74 m and a total weight of 77.5 kg. Videos of the performed simulations can be found in the Multimedia Appendix. The trajectory generator and the whole-body controller are implemented in Matlab / Simulink and run at a rate of 1 kHz. The real-time capability of the algorithms was verified on the TORO computing hardware, an Intel Core i7 computer; the total computation time ranged from 380 to 540 µs, with the whole-body reference trajectory generator, including the angular momentum optimization, taking 75 to 90 µs.In all simulation scenarios, the robot starts and ends in a standing position to demonstrate the capability of the proposed method to initiate and complete highly dynamic movements such as running and jumping while maintaining balance throughout the movement. Unless otherwise specified, the default duration of the stand-to-fly and flight-to-stand phases is 25 ms, the duration of the height change phase is 100 ms, while the maximum distance between CoM and contact point, which serves as the kinematic reachability constraint for calculating the VRP offset Δz in Section . max The whole-body controller limits the commanded ground reaction forces with a friction coefficient of 0.4 and center of pressure limits of +0.08 m on the x-axis and ±0.035 m on the y-axis relative to the center of the foot, expressed in local foot coordinates. For comparison, the TORO foot is 19 cm long and 9.5 cm wide.

[0101] In the first scenario ( Fig. 11), the robot performs a sequence of three jumps, each 30 cm long and 40 cm high. The duration of the stance phase is set to 300 ms, the duration of the flight phase to 350 ms, while the kinematic reach constraint used to calculate the VRP offset Δz is set to d max = 0.9 m. The high steps in this scenario were deliberately chosen to emphasize the need to use both legs simultaneously (i.e., jumping), since the required force acting on the CoM during landing exceeds 2000 N (see Fig. 11b). Fig. Figure 11b also shows how important the stand-to-flight and flight-to-stand phases are for maintaining the continuity of the commanded contact forces. For example, at t = 1.1 s, during the flight-to-stand phase, the trajectory of the normal force f zsimilar to a linear interpolation from 0 to 1000 N, thereby increasing the feasibility of the movement and reducing the effort required to track the commanded force trajectory.

[0102] The second scenario ( Fig. 12) shows a walking and running gait with increasing forward speed (i.e., on the x -axis). The movement begins with a stand-to-movement phase of a total duration of 250 ms (T shift = 50 ms, T hold = 150 ms, see section VA), which ensures the continuity of the trajectory during the transition from standing to walking. The first three walking steps are 15 cm long and are recorded with one- and two-stage phase durations of T SS = 550 ms and T DS= 200 ms, which corresponds to a walking speed of 0.2 m / s. For the next four steps, the step length is increased to 40 cm by using the toe-distance contacts

[31] , while simultaneously increasing the single-stage phase duration to T SS = 800 ms, resulting in a walking speed of 0.4 m / s. During walking, the angular momentum compensation is deactivated for the x- and y-axis and active only for the z-axis (see Fig. ). With these settings, the waist movement generated by the Motion Optimizer

[39] is significantly reduced, since the required angular momentum on the x- and y-axes can be easily generated by the contact moments. At t = 8 s, the robot transitions from walking to running, using the walking-to-running transition method described in Section VB. The walking gait uses toe-edge contacts, has a stance phase duration of 200 ms and a flight phase duration of 150 ms, resulting in a total stride time of 400 ms, including the stance-to-flight and flight-to-stance phases. During the first 4 seconds of walking, the stride length is fixed at 40 cm, corresponding to a running speed of 1 m / s, while at t = 11.5 s, the stride length is increased to 80 cm, resulting in a maximum running speed of 2 m / s. During walking, the angular momentum is compensated in all three axes ( Fig. 12c). Note the angular momentum peaks caused by the impact of the feet on the ground at the end of the flight phases. At t = 16 s, the robot transitions back to walking using the running-to-walking transition method before finally coming to a stop at t = 19 s.

[0103] In the third scenario ( Fig. 13) The robot runs along a curved path, avoiding large obstacles in its path. The duration of the stance and flight phases is the same as in the second scenario, while the transition from standing to walking is a stance-to-movement phase lasting 350 ms (T shift = 25 ms, T hold= 300 ms). The longer duration of this phase compared to the second scenario ensures that the VRP remains within the support range during the more dynamic type of movement (running vs. walking) (see section VA). During the first four running steps, the step length is gradually increased to 60 cm, which is subsequently used as the nominal running step, corresponding to a running speed of 1.5 m / s. From the fifth step onwards, the steps are rotated in the xy-plane in 8° steps, first for four steps to the left (counterclockwise), then for eight steps to the right, and finally for another four steps to the left, so that the final orientation matches the initial orientation (see Fig. 13b). The yaw angle of the planned hip alignment follows the step orientations using fifth-order polynomial interpolations, while the reference trajectory is generated online by the angular momentum optimizer to compensate for the effects of leg movement. The maximum tracking error of the hip alignment over the entire movement is 1.7° for the roll angle, 2.4° for the pitch angle, and 2.2° for the yaw angle ( Fig. 13c), which demonstrates the feasibility of the reference trajectories and the tracking performance of the passivity-based whole-body controller. In the fourth scenario ( Fig. 14), the robot crosses a sequence of nine stepping stones, modeled as cylindrical columns with a height of 1 m and a diameter of 26 cm. The small contact area provided by each column requires precise tracking of the CoM and foot trajectories to successfully traverse the given terrain. The robot starts on a 1.5 m long and 0.5 m wide platform, from which it accelerates to the nominal running speed of 1.5 m / s, corresponding to a step length of 60 cm, a stance phase duration of 200 ms, and a flight phase duration of 150 ms. The columns are arranged very variably on all three Cartesian axes: on the x-axis, the distance between the columns varies between 45 and 85 cm, on the y-axis between 5 and 50 cm, and on the z-axis, the sixth and eighth columns are lowered by 10 cm, while the seventh column is lowered by 20 cm.The final platform is shorter than the starting platform (1 m long compared to 1.5 m), so the robot must decelerate quickly after entering the platform to come to a stop. The maximum tracking error of the CoM position over the entire movement is 3.7 cm on the x-axis, 0.9 cm on the y-axis, and 1.6 cm on the z-axis (. Fig. 14c). The peaks in the DCM tracking error signal are caused by a delayed contact after the flight phases; for example, at t = 4.4 s, when the robot steps on the seventh pillar, the peak is caused by a 2.5 cm tracking error of the left foot on the z-axis. However, after contact is established, the DCM tracking error is quickly corrected by the whole-body control. VIII. Advantage of the invention

[0104] When discussing the advantages of the proposed invention, one must distinguish between the 3D DCM framework and the specific method for generating trajectories for running, jumping, and the various gait transitions. A 3D DCM frame: advantages

[0105] As a reformulation of CoM dynamics, the 3D-DCM framework maintains its generality and can be used with any legged robot. In contrast to simplified models such as SLIP or LIPM, 3D-DCM is a reduced model of whole-body dynamics. By exploiting the reverse-time stability of DCM dynamics, the framework can be used to generate stable reference trajectories for piecewise linearly interpolated VRP trajectories of arbitrary duration. In particular, the framework generates stable trajectories for arbitrarily slow motions, which poses a particular challenge for methods that only consider the forward-time dynamics of the DCM

[26] . Furthermore, the 3D-DCM framework uses closed-form continuous-time trajectories and does not require discretization of the dynamics, enabling trajectory generation at arbitrary sampling rates.For example, in our implementation, we calculate the instantaneous CoM, DCM, and VRP every millisecond, i.e., at the same rate as the whole-body controller execution rate of 1 kHz.

[0106] Another advantage of computational efficiency is that the framework can generate CoM and DCM trajectories for long-duration movements. For example, for running, we use a preview window of 5 steps, which corresponds to a plan duration of 2.6 seconds, including the transition from running to standing at the end of the movement. The reasons for this choice can be found in the data in Table 1. Table 1 shows discontinuities induced in the instantaneous DCM calculation by a DCM endpoint change of 0.6 m for different lengths of the motion plan preview window. Steps Motion phases Total duration Right x 1 10 1.0 s 0.1438 m 2 14 1.4 s 0.0389 m 3 18 1.8 s 0.0105 m 4 22 2.2 s 0.0028 m 5 26 2.6 s 0.0008 m

[0107] Due to the approach to calculating the backward time, any change in the position of the DCM endpoint ξf to a discontinuity in the calculated instantaneous DCM ξ. The extent of the discontinuity is reduced by longer preview windows; with 5 steps, a change of 0.6 m in the DCM endpoint on the x-axis results in a DCM discontinuity of less than 1 mm. This choice is made possible by the use of a highly efficient waypoint calculation algorithm: The DCM and CoM waypoints of the corresponding 26 motion phases can be calculated by the TORO computing hardware using Algorithm 1 in approximately 25 µs.

[0108] A limitation of the 3D-DCM framework is that it cannot generate a specific contact force profile, such as the biologically inspired forces proposed by BID

[16] . The 3D-DCM framework naturally requires an eCMP trajectory (or an associated VRP trajectory), for which it calculates the corresponding CoM trajectory, while the contact force profile is subsequently given by Equation (2). A comparison of the generated contact forces on the x- and z-axis for a running stance phase is shown in Fig. 15. Note that the required 3D-DCM forces are smaller than the BID forces for the same motion parameters. On the other hand, the 3D-DCM forces, although continuous, are not continuously differentiable, which potentially complicates their generation for robots with elastic elements in the legs. Nevertheless, smoother force profiles can be generated by the 3D-DCM system if required by using higher-order polynomial interpolations in the eCMP reference trajectories. Another advantage of 3D-DCM is that the use of a common framework for walking and running greatly simplifies the implementation of gait transitions, as shown in Section V, whereas combining BID for running with 3D-DCM for walking requires significantly more computational effort

[17] .

[0109] These properties of the 3D DCM framework can be exploited by more computationally intensive methods such as Model Predictive Control (MPC) or Reinforcement Learning (RL) to reduce the overall computational effort. For example, the DCM calculation in the backward direction can be combined with MPC, as shown in

[34] , and RL can be used to learn VRP offsets during the stance phase to generate a contact force profile more similar to BID. We intend to pursue these ideas in our future work. B. Whole-body reference trajectories

[0110] An important factor influencing the robustness of dynamic bipedal locomotion is the controller's ability to maintain the center of pressure close to the center of the foot contact patch during stance phases. This can be achieved by appropriately designed reference trajectories or by using the upper body to generate angular momentum, thereby reducing the required contact moments. In this invention, the eCMP reference trajectory is generated assuming a constant angular momentum; the eCMP waypoints are either placed in the center of the foot contact patch or, alternatively, the eCMP follows a heel-toe trajectory during each stance phase, which has the additional advantage of reducing the magnitude of the required contact forces on the local x-axis. The applied motion optimization

[39] significantly reduces the center of gravity angular momentum during walking, as demonstrated by the simulations performed (see Fig. 12c).

[0111] The leg reference trajectories are generated using simple fifth-order polynomials for each Cartesian translation and rotation axes. Despite their simplicity, the controller is able to track the generated movements without difficulty, as demonstrated by the simulations performed. However, a challenge for highly dynamic movements is avoiding knee singularities. During takeoff and landing, this problem is solved by the VRP offset calculation proposed in Section IV-G, which ensures that the distance between the CoM and the contact point is limited by a configurable maximum value. During the flight phase, however, the independent calculation of the leg and CoM trajectories can lead to violations of the kinematic constraint.

[0112] A limitation of the angular momentum optimization is the immediate calculation of the upper body movement without an explicitly planned angular momentum curve or a preview over a future time horizon. A noticeable side effect is the gradual upward movement of the arms over several running steps, which can be seen during the fast running segment of the second simulation ( Fig. Finally, the full-body simulations demonstrate the feasibility of the generated trajectories in terms of maintaining balance throughout the movement, controlling the angle, and complying with contact constraints (friction cone, unilateral normal force, and limited contact area). However, the required joint velocities and torques exceed the capabilities of the real humanoid robot TORO, preventing us from demonstrating the proposed method in an experimental setup.

[0113] This work proposes a motion planner based on the 3D-DCM framework capable of generating a CoM reference trajectory for movements involving flight phases such as running and jumping. Exploiting the backward time stability of DCM dynamics, the motion planner generates stable reference trajectories for arbitrary contact sequences and time parameterizations. Furthermore, the use of the 3D-DCM framework for running and jumping greatly simplifies the implementation of gait transitions, such as walking to running or running to walking. The proposed motion plan was validated in various simulations using the humanoid robot TORO, demonstrating the feasibility of the generated trajectories.

[0114] Two advantages of the proposed motion planner, the unified approach to describing different motion phases and the high efficiency in calculating the reference trajectories, increase the versatility of the planner and make it suitable for use with a wide variety of bipedal, quadrupedal, or fully humanoid robots.

[0115] The project leading to this application has received funding from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation programme (grant agreement No 819358). Publications referred to [1] T. Takenaka, T. Matsumoto, T. Yoshiike, and S. Shirokura, "Real time motion generation and control for biped robot -2nd report: Running gaitpattern generation," in IEEE / RSJ International Conference on Intelligent Robots and Systems (IROS), 2009, pp. 1092-1099. [2] K. Kojima, Y. Kojio, T. Ishikawa, F. Sugai, Y. Kakiuchi, K. Okada, and M. Inaba, „A robot design method for weight saving aimed at dynamic motions: Design of humanoid JAXON3-P and realization of jump motions,“ in IEEE-RAS International Conference on Humanoid Robots , 2019, pp. 586-593. [3] Boston Dynamics, „Atlas,“ https: / / www.bostondynamics.com / atlas, Accessed: 2023-02-10. [4] J. Siekmann, Y. Godse, A. Fern, and J. Hurst, „Sim-to-real learning of all common bipedal gaits via periodic reward composition,“ in IEEE International Conference on Robotics and Automation (ICRA), 2021, pp. 7309-7315. [5] D. Crowley, J. Dao, H. Duan, K. Green, J. Hurst, and A. Fern, „Optimizing bipedal locomotion for the 100m dash with comparison to human running,“ in IEEE International Conference on Robotics and Automation (ICRA), 2023, pp. 12 205-12 211. [6] R. J. Full and D. E. Koditschek, „Templates and anchors: neuromechan- ical hypotheses of legged locomotion on land,“ Journal of Experimental Biology , vol. 202, no. 23, pp. 3325-3332, 1999. [7] M. H. Raibert, H. B. Brown Jr, and M. Chepponis, „Experiments in balance with a 3D one-legged hopping machine,“ The International Journal of Robotics Research , vol. 3, no. 2, pp. 75-92, 1984. [8] S. Kajita, F. Kanehiro, K. Kaneko, K. Yokoi, and H. Hirukawa, „The 3D linear inverted pendulum mode: A simple modeling for a biped walking pattern generation,“ in IEEE / RSJ International Conference on Intelligent Robots and Systems (IROS), vol. 1, 2001, pp. 239-246 [9] H. Geyer, A. Seyfarth, and R. Blickhan, „Compliant leg behaviour explains basic dynamics of walking and running,“ Proceedings of the Royal Society B: Biological Sciences , vol. 273, no. 1603, pp. 2861-2867, 2006.

[10] J. Englsberger, C. Ott, and A. Albu-Schäffer, „Three-dimensional bipedal walking control based on divergent component of motion,“ IEEE Trans- actions on Robotics , vol. 31, no. 2, pp. 355-368, 2015

[11] R. Blickhan, „The spring-mass model for running and hopping,“ Journal of biomechanics , vol. 22, no. 11-12, pp. 1217-1227, 1989

[12] P. M. Wensing and D. E. Orin, „High-speed humanoid running through control with a 3D-SLIP model,“ in IEEE / RSJ International Conference on Intelligent Robots and Systems (IROS), 2013, pp. 5134-5140.

[13] G. Secer and U. Saranli, „Control of planar spring-mass running through virtual tuning of radial leg damping,“ IEEE Transactions on Robotics , vol. 34, no. 5, pp. 1370-1383, 2018.

[14] X. Xiong and A. D. Ames, „Bipedal hopping: Reduced-order model embedding via optimization-based control,“ in IEEE / RSJ International Conference on Intelligent Robots and Systems (IROS), 2018, pp. 3821- 3828.

[15] S. Li, H. Chen, W. Zhang, and P. M. Wensing, „Quadruped robot hop- ping on two legs,“ in IEEE / RSJ International Conference on Intelligent Robots and Systems (IROS), 2021, pp. 7448-7455.

[16] J. Englsberger, P. Kozłowski, C. Ott, and A. Albu-Schäffer, „Biologically inspired deadbeat control for running: From human analysis to humanoid control and back,“ IEEE Transactions on Robotics , vol. 32, no. 4, pp. 854-867, 2016.

[17] T. Egle, J. Englsberger, and C. Ott, „Analytical center of mass trajectory generation for humanoid walking and running with continuous gait tran- sitions,“ in IEEE-RAS International Conference on Humanoid Robots , 2022, pp. 630-637.

[18] M. Vukobratovi 'c and J. Stepanenko, „On the stability of anthropomor- phic systems,“ Mathematical Biosciences , vol. 15, no. 1-2, pp. 1-37, 1972.

[19] J. Pratt, T. Koolen, T. De Boer, J. Rebula, S. Cotton, J. Carff, M. Johnson, and P. Neuhaus, „Capturability-based analysis and control of legged locomotion, part 2: Application to m2v2, a lower-body humanoid,“ The International Journal of Robotics Research , vol. 31, no. 10, pp. 1117- 1133, 2012.

[20] T. Takenaka, T. Matsumoto, and T. Yoshiike, „Real time motion gen- eration and control for biped robot -1st report: Walking gait pattern generation,“ in IEEE / RSJ International Conference on Intelligent Robots and Systems (IROS), 2009.

[21] S. Kajita, T. Nagasaki, K. Kaneko, and H. Hirukawa, „ZMP-based biped running control,“ IEEE Robotics & Automation Magazine , vol. 14, no. 2, pp. 63-72, 2007.

[22] R. Tajima, D. Honda, and K. Suga, „Fast running experiments involving a humanoid robot,“ in IEEE International Conference on Robotics and Automation (ICRA), 2009, pp. 1571-1576.

[23] B. Ugurlu, E. Sariyildiz, T. Kawasaki, and T. Narikiyo, „Agile and stable running locomotion control for an untethered and one-legged hopping robot,“ Autonomous Robots , vol. 45, pp. 805-819, 2021.

[24] T. Koolen, M. Posa, and R. Tedrake, „Balance control using center of mass height variation: limitations imposed by unilateral contact,“ in IEEE-RAS International Conference on Humanoid Robots , 2016, pp. 8-15.

[25] D. Ahn and B.-K. Cho, „Online jumping motion generation via model predictive control,“ IEEE Transactions on Industrial Electronics , vol. 69, no. 5, pp. 4957-4965, 2021.

[26] F. M. Smaldone, N. Scianca, L. Lanari, and G. Oriolo, „From walking to running: 3D humanoid gait generation via MPC,“ Frontiers in Robotics and AI , vol. 9, 2022.

[27] T. Sugihara, K. Imanishi, T. Yamamoto, and S. Caron, „3D biped locomotion control including seamless transition between walking and running via 3D ZMP manipulation,“ in IEEE International Conference on Robotics and Automation (ICRA), 2021, pp. 6258-6263.

[28] Z. He and K. Yamamoto, „Humanoid running based on 3D COG-ZMP model and resolved centroidal viscoelasticity control,“ in IEEE-RAS International Conference on Humanoid Robots , 2022, pp. 645-650.

[29] K. Yamamoto, N. Kobayashi, T. Ishigaki, and Y. Sakemi, „Integration of variable-height and hopping strategies for humanoid push recovery,“ in IEEE / RSJ International Conference on Intelligent Robots and Systems (IROS), 2022, pp. 12 979-12 985.

[30] G. Mesesan, J. Englsberger, C. Ott, and A. Albu-Schäffer, „Convex properties of center-of-mass trajectories for locomotion based on divergent component of motion,“ IEEE Robotics and Automation Letters , vol. 3, no. 4, pp. 3449-3456, 2018.

[31] G. Mesesan, J. Englsberger, G. Garofalo, C. Ott, and A. Albu-Sch" affer, „Dynamic walking on compliant and uneven terrain using DCM and passivity-based whole-body control,“ in IEEE-RAS International Con- ference on Humanoid Robots , 2019, pp. 25-32.

[32] J. Englsberger, G. Mesesan, C. Ott, and A. Albu-Sch" affer, „DCM-based gait generation for walking on moving support surfaces,“ in IEEE-RAS International Conference on Humanoid Robots , 2018, pp. 1-8.

[33] G. Mesesan, J. Englsberger, B. Henze, and C. Ott, „Dynamic multi- contact transitions for humanoid robots using divergent component of motion,“ in IEEE International Conference on Robotics and Automation (ICRA), 2017, pp. 4108-4115.

[34] M. A. Hopkins, D. W. Hong, and A. Leonessa, „Humanoid locomo- tion on uneven terrain using the time-varying divergent component of motion,“ in IEEE-RAS International Conference on Humanoid Robots , 2014, pp. 266-272.

[35] S. Hanasaki, Y. Tazaki, H. Nagano, and Y. Yokokohji, „Running trajectory generation including gait transition between walking based on the time-varying linear inverted pendulum mode,“ in IEEE-RAS International Conference on Humanoid Robots , 2022, pp. 851-857.

[36] J. Englsberger, A. Werner, C. Ott, B. Henze, M. A. Roa, G. Garofalo, R. Burger, A. Beyer, O. Eiberger, K. Schmid, and A. Albu-Schäffer, „Overview of the torque-controlled humanoid robot TORO,“ in IEEE- RAS International Conference on Humanoid Robots , 2014, pp. 916-923.

[37] M. B. Popovic, A. Goswami, and H. Herr, „Ground reference points in legged locomotion: Definitions, biological trajectories and control implications,“ The International Journal of Robotics Research , vol. 24, no. 12, pp. 1013-1032, 2005.

[38] T. Koolen, T. De Boer, J. Rebula, A. Goswami, and J. Pratt, „Capturability-based analysis and control of legged locomotion, part 1: Theory and application to three simple gait models,“ The International Journal of Robotics Research , vol. 31, no. 9, pp. 1094-1113, 2012.

[39] R. Schuller, G. Mesesan, J. Englsberger, J. Lee, and C. Ott, „Online cen- troidal angular momentum reference generation and motion optimization for humanoid push recovery,“ IEEE Robotics and Automation Letters , vol. 6, no. 3, pp. 5689-5696, 2021.

[40] J. Englsberger, G. Mesesan, and C. Ott, „Smooth trajectory genera- tion and push-recovery based on divergent component of motion,“ in IEEE / RSJ International Conference on Intelligent Robots and Systems (IROS), 2017, pp. 4560-4567.

[41] B. Henze, M. A. Roa, and C. Ott, „Passivity-based whole-body balancing for torque-controlled humanoid robots in multi-contact scenarios,“ The International Journal of Robotics Research , vol. 35, no. 12, pp. 1522- 1543, 2016.

[42] S.-H. Hyon, J. G. Hale, and G. Cheng, „Full-body compliant human- humanoid interaction: balancing in the presence of unknown external forces,“ IEEE Transactions on Robotics , vol. 23, no. 5, pp. 884-898, 2007.

[43] R. Hinata and D. N. Nenchev, „Balance stabilization with angular momentum damping derived from the reaction null-space,“ in IEEE-RAS International Conference on Humanoid Robots , 2018, pp. 188-195.

[44] F. Kanehiro, H. Hirukawa, and S. Kajita, „OpenHRP: Open architecture humanoid robotics platform,“ The International Journal of Robotics Research , vol. 23, no. 2, pp. 155-165, 2004.

[45] R. Schuller, G. Mesesan, J. Englsberger, J. Lee, and C. Ott, „Online learning of centroidal angular momentum towards enhancing dcm- based locomotion,“ in IEEE International Conference on Robotics and Automation (ICRA), 2022, pp. 10 442-10 448.

Claims

[1] Method for planning and controlling the balance and locomotion of a humanoid bipedal robot, with at least one control and planning device which carries out real-time planning and control of DCM trajectories in 3D by means of a model of the three-dimensional divergent movement component (3D-DCM), characterized by that in the divergent motion component model (3D-DCM), equations of motion for the flight phase are included, introducing a unified formulation that can be used for walking, running and jumping movements, allowing to encode both the direction and magnitude of the external (e.g. leg) forces and the total force (i.e. external forces plus gravity) acting on the robot, whereby a method for real-time planning and control of DCM trajectories in 3D is presented based on eCMP, VRP and DCM. [2] Method according to claim 1, characterized by, dass für den DCM-Startpunkt wie folgt bestimmbar ist: ξφ,0=αφ,ξυφ,0+βφ,ξυφ,T+γφ,ξξφ,T +δφ,ξxφ,0+εφ,ξg. bewertung für die Standphase die Koeffizienten d φ,ξ and e φ,ξ are equal to 0, while af,ξ=,(1−bTφ+bTφeTφb) and βφ,ξ=(bTφ−b+TφTφe−Tφb) and γφ,ξ=e−Tφb und behandung für die Flugphase die Koeffizienten a φ,ξ and b φ,ξ 0 are, while γφ,ξ=bb+Tφ and δφ,ξ=Tφb+Tφ and εφ,ξ=bTφ2+2b2Tφ2(b+Tφ) correspond. [3] Verfahren nach Anspruch 1 oder 2, characterized by , dass der Wegpunktes des Masseschwerpunkts (CoM-Endpunkts) wie folgt bestimmbar ist: xφ,T=αφ,xϑφ,0+βφ,xϑφ,T+γφ,xξφ,T +δφ,xxφ,0+εφ,xg, bewertung für die Standphase die Koeffizienten d φ,ξ and e φ,ξ are equal to 0, while af,ξ=,(1−bTφ+bTφe−Tφb) and βφ,ξ=(bTφ−b+TφTφe−Tφb) and γφ,ξ=e−Tφb and where for the flight phase the coefficients α φ,ξ and β φ,ξ 0, while γφ,ξ=bb+Tφ and δφ,ξ=Tφb+Tφ and εφ,ξ=bTφ2+2b2Tφ2(b+Tφ) are equivalent to. [4] Method according to one of claims 1 to 3, characterized by that an “Enhanced Centroidal Moment Pivot point” (eCMP) and a “Virtual Repellent Point” (VRP) are introduced in the divergent motion component model (3D-DCM) and that the eCMP and the VRP are used for real-time planning and control of the DCM trajectories in 3D. [5] Method according to one of claims 1 to 4, characterized by that a DCM and CoM waypoint calculation Ξ and X in matrix form with respect to the VRP waypoints V, the DCM endpoint ξ f , , the CoM starting point x s and the gravitational vector g is as follows: [I−S0TΓξST−S0TΔξS0−STTΓxSTI−STTΔxS0][ΞX]= [S0TAξS0+S0TBξSTSTTAxS0+STTBxST]V+[S0TεξSTTεx]gT +[sf0nw×1]ξfT+[0nw×1 ss]xsT, where with V=⌈v1...vnw⌉T∈ℝnw×3 the VRP waypoints, with Ξ=⌈ξ1...ξnω⌉T∈ℝnω×3, the DCM waypoints and with X = [x1 ... x nw ] T ∈ ℝ nw×3 the CoM waypoints are shown, where S0 = [I nφ×nφ 0 nφ×1 ] ∈ ℝ nφ×nw selects the starting points of the movement phase, while ST=⌈0nφ×1 Inφ×nφ⌉∈ℝnφ×nω selects the endpoints from an arbitrary waypoint matrix V, Ξ, or X. [6] Method according to one of claims 1 to 5, characterized by that an algorithm for the waypoint calculation with a computational complexity of O(nφ2) is used. [7] Method according to one of claims 1 to 6, characterized bythat a waypoint calculation algorithm is used which requires three iterations over all movement phases, where in the first two iterations partial solutions for the DCM and CoM waypoints are calculated from the known quantities and the coefficients for the unknown CoM waypoints are tracked, with the third iteration using the collected coefficients to calculate the complete solution. [8] Method according to one of claims 1 to 7, characterized by that an algorithm I is used for the waypoint calculation, where the algorithm I looks like this: [9] Method according to one of claims 1 to 8, characterized by that transition phases are provided between the four types of movement mentioned: standing, walking, running and jumping. [10] Method according to claim 9, characterized bythat the control and planning device constructs a sequence of movement phases with alternating stance and flight phases for running, single and double support phases for walking steps, and that gait transition phases are inserted into the phase sequence. [11] System for planning and controlling the balance and locomotion of a humanoid bipedal robot, comprising a humanoid bipedal robot and at least one control and planning device which is designed to carry out a method according to one of claims 1 to 10.

Images