Mechanical Design of a Biped Robot FORREST and an Extended Capture-Point-Based Walking Pattern Generator

Abstract

In recent years, many studies have shown that soft robots with elastic actuators enable robust interaction with the environment. Compliant joints can protect mechanical systems and provide better dynamic performance, thus offering huge potential for further developments of humanoid robots. This paper proposes a new biped robot. The new robot combines a torque sensor-based active elastic hip and a spring-based passive elastic knee/ankle. In the first part, the mechanical design is introduced, and in the second part, the kinematics and dynamics capabilities are described. Furthermore, we introduce a new extended capture-point-based walking pattern generator that calculates footstep positions, which are used as input for the controller of our new biped robot. The main contribution of this article is the novel mechanical design and an extended walking pattern generator. The new design offers a unique solution for cable-driven bipeds to achieve both balancing and walking. Meanwhile, the new walking pattern generator can generate smooth desired curves, which is an improvement over traditional generators that use a constant zero-moment-point (ZMP). A simple cartesian controller is applied to test the performance of the walking pattern generator. Although the robot has been built, all experiments regarding the pattern generator are still simulated using MATLAB/Simulink. The focus of this work is to analyze the mechanical design and show the capabilities of the robot by applying a new pattern generator.

1. Introduction

While humans can easily walk with two legs, it remains very challenging to realize two-legged humanoid robots. Over the years, several biped robots have been developed. Compared to wheeled mobile robots, two-legged humanoid robots offer significant advantages, as they can move in difficult or uneven terrains, climb stairs, walk over obstacles, and reach spaces as humans do, even in challenging environments. Therefore, when designing biped robots, researchers need to ensure that they can move like humans with minimal energy consumption and that they are highly compliant when interacting with their environment or with humans.

In the past, the joints of most humanoid robots have been directly driven by actuators, which consist of motors and gears. One of the earliest famous humanoid robots was the ASIMO [1], developed by Honda and first introduced in 2000. It can walk, run, and jump. Additionally, AIST developed the humanoid robot series HRP-(1–5) [2,3,4,5,6]. Unlike the others, the HRP-4 is a smart version with a slim lightweight body and a female shape. The latest generation, HRP-5, is designed to perform heavy tasks or to operate in hazardous environments. Toyota developed the robot T-HR3 [7], which enables the entire robot body to be controlled via a wearable device that maps the movements of the user’s hands, arms, and feet to the robot. Takanishi Laboratory developed the robot WABIAN-2R [8], which can walk with straight knees, heel contact, and toe-off. Another excellent humanoid robot is the robot TORO [9], developed by the German Aerospace Center (DLR). The robot TORO is built based on the legs of the former DLR biped [10] and is completely impedance-controlled based on the torque controller drive technology, similar to the DLR lightweight arms. It can adapt to an uneven environment similar to a foam mattress. Researchers at TU Munich developed two robots, Johnnie [11] and Lola [12]. The robot Johnnie can walk on flat and uneven ground and around curves. The next-generation robot, Lola, can recognize obstacles and achieve dynamic real-time movement through planning and control. Researchers at Oregon State University invented the biped robot Cassie [13], which has two bird-like legs. The researchers at PAL Robotics developed the biped robot REEM-C in 2013 [14]. In 2017, PAL built their second biped robot TALOS [15], which uses torque control. The robot can walk on uneven terrain. The newer biped robot, called Kangaroo, was released in 2022 [16]. All the robot’s joints in Kangaroo are driven by ball-screw linear actuators. The company Tesla also released their biped robot called Optimus in 2022 [17].

All of the abovementioned robots use motor-gear drive systems and rigid joints. However, their dynamic behavior is limited by motor performance. Thus, researchers have developed bipeds that are driven by elastic actuators to improve the dynamic performance and safe interaction with environments. Atlas [18] is one of the world’s most famous humanoid robots developed by Boston Dynamics. Atlas uses a hydraulic actuator system to drive all the joints. Compared to motor-gear drivers, the dynamic performance of a hydraulic driver system is more powerful. At Ritsumeikan University, researchers developed a torque-controlled hydraulic humanoid robot called TaeMu [19]. It is capable of full-body compliant balancing. In addition to hydraulic bipeds, another possibility to improve bipeds is to combine motor-gearbox systems with elastic springs. The robot cCub [20] is such a compliant robot that uses compact units based on series elastic actuators. A serial elastic actuator combines the common motor with an elastic spring and thus can provide an elasticity with constant stiffness [21]. Furthermore, the adjustable stiffness of a leg may improve the performance in walking on a variety of terrains [22]. Although variable stiffness actuators that can adjust their stiffness are available [23], they have not yet been applied to biped robots. The humanoid robot Valkyrie [24], which was developed by NASA, applies serial elastic actuators to drive every joint. The robot Valkyrie aims to solve challenging tasks in space. Another elastic biped robot is MABEL [25]. Its hip is constrained to revolute motion in the sagittal plane. It can only walk in a circle around a tower. In [26], the Technical University of Darmstadt presented an elastic biped robot called BioBiped. The robot BioBiped has two three degrees of freedom (DOF) legs and uses musculoskeletal technology. This robot can jump while its trunk is externally constrained to vertical translation. Another similar elastic biped robot, the C-Runner [27], was developed by the DLR. The elastic elements give the legs high impact robustness and are able to store energy in mechanical springs during dynamic movements. The bipeds MABEL, BioBiped, and C-Runner use elastic cable-driven systems. The springs can absorb shocks and store energy for running or jumping. They show good dynamic performance and huge potential for running and jumping. Due to the design difficulties of elastic cable-driven systems for multi-DoF joints (hips and ankles), such bipeds lack of the necessary number of DoFs for balancing. They can only walk around a central tower or in a sagittal plane. In order to achieve a fully compliant leg, we add a spring in the knee and design a new elastic cable-driven ankle with two DoFs to overcome these disadvantages. By applying the new ankle, higher velocities and torques can be reached. With our new leg, our robot is able to walk in a human-like way and save energy. Moreover, jumping and running will become possible in the near future.

In the past, many researchers have developed different control strategies to achieve stable balancing. A comprehensive robot control strategy for bipeds comprises a walking pattern generator and a balance controller. A walking pattern generator plans the desired trajectory for the biped. A simple generator based on the linear inverted pendulum (LIP) and ZMP concepts was developed by AIST [28]. The walking biped robot model is commonly simplified as a LIP, which assumes that the robot’s center of mass (CoM) can be approximated as a point mass located above the stance foot. This simplification facilitates the development of control strategies for maintaining dynamic balance during walking, and this technique forms the basis of most walking pattern generators used today. Moreover, a walking pattern generator using preview control of ZMP was described in [29]. A virtual ZMP tracking control was utilized to generate the actual trajectory from the previewed reference ZMP in future steps. In [30], a third-order polynomial interpolation curve was employed for online gait generation, with the aim of facilitating its implementation in real systems. The real-time walking and running gait pattern generation for the bipedal robot ASIMO was presented in [31,32], where the generator utilized a model consisting of three point masses. One point mass was located at the end of the LIP, another was situated at the ankle of the supporting foot, and the third was positioned at the ankle of the swing foot. The compensation of the dynamics error between the approximate dynamic model and the real ASIMO robot was discussed in [33]. A walking pattern generator based on the capture point (CP) was described in [34]. The capture point was determined from the orbital energy of the LIP model. Specifically, the CP represents the point at which the velocity of the CoM is zero. This property enables the CP to be used for calculating a suitable foothold, thus preventing the robot from falling. The walking robot was stabilized by CP and ZMP controllers. Most generators for walking rely on a constant CoM height to simplify the planning process. However, this approach limits the robot to bending its knees during locomotion. In contrast, human-like walking involves a stretched knee, heel strike, and toe push-off, resulting in variations in the CoM height throughout the gait cycle. An enhanced centroidal moment pivot and virtual repellent point were presented in [35] to expand the 2D capture point concept into a 3D divergent component of motion (DCM). This approach facilitates the planning of the robot’s CoM height trajectory, thereby enhancing the robot’s ability to adapt to uneven terrain. The method, as presented in [36], accomplished the CoM height trajectory by adjusting the natural frequency of the DCM for locomotion on uneven terrain. A concept of spatially quantized dynamics (SQD) was proposed in [37] to achieve a more human-like walking pattern with a stretched knee. It involved discretizing the trajectory into constant distance intervals and transforming the spatial walking pattern into the time domain. In [38], a pattern generator for walking with variable height was presented and enabled 3D walking over uneven terrains based on capture inputs. The paper [39] proposed a planning algorithm capable of generating continuous-time walking patterns, including seamless transitions between flat-contact and heel-to-toe walking gaits. In addition to the conventional linear inverted pendulum (LIP) model, a spring-loaded inverted pendulum (SLIP) model was proposed in [40]. This model incorporated a virtual spring along the inverted pendulum, which can be compressed to absorb energy upon the landing of the swing foot and subsequently accelerate the CoM when the swing foot takes off. The SLIP model proved to be a valuable tool in controlling and analyzing the running and hopping of bipeds. In [41], a 3D SLIP model was employed to achieve high-speed running for a whole-body humanoid robot via simulation. Furthermore, in [42], a task-decomposed energy-exchange dynamics learning method was proposed, which combined model-based reinforcement learning to capture the simplified SLIP biped dynamics and utilize them for control. The difference between the LIP and SLIP is explained in Appendix A.

The balance controller is utilized to stabilize bipedal robots and prevent them from falling over. Controllers for ASIMO, including ground reaction force control and model ZMP control, were presented in [43]. In [44], an approach for balancing a humanoid robot with multiple contacts was presented and implemented in the robot TORO. A passivity-based controller was applied to the robot TORO in [45] to achieve balancing on soft terrain, such as a mattress. The controllers for multicontact and divergent components of motion were combined in [46]. In [47], the passivity controller was extended so that the robot TORO could balance itself on an unstable ground surface. An energy-efficient controller based on optimization was presented in [48]. The researcher minimized the energy cost of walking and calculated the target value of each step incrementally. An online foot position compensator was proposed in [49] to improve the robustness of walking. The control strategy of the biped robot DURUS was introduced in [50]. It divided walking into a sequence of distinct events and realized multicontact walking. The control strategy for the biped robot HRP-4 was presented in [51]. A 3D LIP method was implemented, and feedback linearization was used for joint tracking. A two-level variable horizon predictive controller was proposed in [52]. The two levels calculated the landing location and generated trajectories for landing in the desired time. A model-predictive-based control law, which uses extended centroidal dynamics, was presented in [53] to consider heavy limbs. An online nonlinear model predictive control approach was developed in [54] to realize the desired walking behaviors. In [55], the authors provided a robust controller using model predictive control to compensate for the gap between simulation and reality. The authors of [56] proposed a full-body predictive model control scheme based on differential dynamic programming that took into account the full dynamics of the system and determined the optimal actuation for the robot’s lower body. A passivity-based inverse dynamics controller using a global energy tank was introduced in [57]. The approach used a task space inverse dynamics quadratic programming to calculate the desired torque for satisfying a set of tasks. Currently, most controllers are designed for rigid biped robots, and only a few studies have introduced control approaches for passive elastic bipeds [20,24]. The control strategy with passive elasticity is a challenge for our future work.

In this paper, a new design of a biped robot, called FORREST, is described, which was developed at the TU Chemnitz and is shown in Figure 1. The new biped robot combines a software-based elastic hip and spring-based passive elastic knees/ankles for walking and running. Each leg provides six DoFs. Currently, most cable-driven biped robots cannot maintain balance. We designed a novel elastic cable-driven ankle to fill this gap, by providing better dynamics. A CP-based extended walking pattern generator is introduced in this work. A simple CP controller and a torque-based Cartesian PD controller are used to test the walking pattern generator. This article is structured as follows: First, Section 2 provides an overview of the robot. In in Section 3, the mechanical design of each part of the robot is described. Section 4 analyzes the kinematic and dynamic performance of the knee and ankle, and the dynamic model is introduced. The extended walking pattern generator is presented in Section 5. Finally, Section 6 describes the control strategy and shows the experiments of locomotion with FORREST conducted in simulation.

2. Overview of FORREST’s Design

This section provides an overview of the design of the biped robot FORREST. Figure 2 shows the kinematic scheme of the biped robot and the CAD model of one leg. Each leg consists of six joints: a three DoF hip, a one DoF knee, and a two DoF ankle. Most biped robots have six DoFs per leg. The total weight of the robot is 35 $\mathrm{k}$$\mathrm{g}$.

Table 1. Overview of weight and height.

Segment	Weight [kg]	Height [mm]
base	17.14	250
thigh	5.83	405
calf	2.8	405
foot	0.35	105
total	35.1	1165

lists the weight and the size of each segment.

The joints of the hip are driven by actuators from SENSODRIVE, which consist of a BLDC motor, a harmonic gear, and an integrated torque sensor. This allows for software-based compliance at the hip. The knee has one degree of freedom and is driven by a motor with a ball screw spindle with a pitch of 4 $\mathrm{m}$$\mathrm{m}$. The ankle has two degrees of freedom and is driven by two MAXON motors using a parallel mechanism. The performance of the joints, considering the gear and drive system, is listed in

Table 2. Torque and speed performance of joints.

Joint	Max Toque [Nm]	Max Speed [rpm]
hip 1	120	31.9
hip 2	107	19.9
hip 3	315	29
knee	195	15.8
ankle 1	212	35
ankle 2	138	47

.

Considering the human anatomy, the muscles in the thigh and calf regions play an important role while humans are running. These muscles can dampen impact, store energy, and release it during bursts, exhibiting behavior similar to that of springs, which we will utilize in our work. However, due to the design difficulty of elastic cable-drive systems for multi-DoF joints such as the hip and ankle, most cable-driven bipeds cannot balance [25,26,27]. Our goal is to overcome this disadvantage by designing a cable-driven ankle for our prototype. Only the knee and ankle joints use springs to store impact energy and achieve human-like running. FORREST’s main purpose is to serve as an experimental platform for combining a torque-controlled software-based elastic hip with mechanical elastic knees/ankles. Further advantages of this design lie in the high reactivity of the new ankle joint, improving the possibility for balancing, and the combination of spring-based deformable joints. Compared to the other bipeds shown in

Table 3. Comparison with other bipeds.

	FORREST	Toro	Lola	TALOS	C-Runner	MABEL	BioBiped
balance	◯	◯	◯	◯	×	×	×
compliance	active + passive	active	active	active	passive	passive	passive

, a typical biped can balance itself but only has active compliance, while most cable-driven bipeds have passive compliance but cannot balance themselves. However, our biped not only includes passive compliance but can also balance itself.

3. Design of the New Biped Robot

In this section, the mechanical design of each part of the leg is discussed. The hip joints are driven by a motor–gearbox system, while the knee is driven by an elastic ball screw spindle system. The ankle, on the other hand, is driven by an elastic parallel mechanism. The new biped robot, FORREST, aims to combine software-based compliance, which means joints where the output torque is measured at the link side and fed back into the control, and mechanical compliance, which enables it to be more responsive than previous solutions.

3.1. Hip

The structure of the hip is similar to that of most other biped robots to simplify the design. Figure 3 shows the front and back views of the hip. All the joints of the hip are rigid joints without springs and are driven by SENSORDRIVE actuators. The compliance of the joints is achieved by the software-based impedance controller, which offers compliance in one DOF. A high-precision inertia measurement unit is mounted on the hip to estimate the orientation and angular velocity.

3.2. Knee

As is often the case, the knee joint of the biped uses one actuator, which is mounted on the knee axis. To achieve an elastic knee joint, the knee of our robot is driven by a ball screw linear system, which is presented in Figure 4. The linear system consists of a Maxon EC Flat Motor, a ball screw spindle with a 4 mm pitch, a nut, and a compression spring. In the human body, muscles in the thigh can absorb shocks when jumping. The elastic linear actuator in our robot provides the same function as these muscles. Another advantage is that the motors, which drive the ankle, can also be mounted in the thigh when the knee is driven by a ball screw. Compared to other robots with actuators mounted on joint axes, our knee and ankle motors are mounted on the thigh to raise the CoM of the biped and reduce the weight of the calf. The stiffness of the spring in the knee drive system is 74 $\mathrm{N}$/$\mathrm{m}$$\mathrm{m}$, and it can produce a maximum force of 1800 $\mathrm{N}$. Assuming that the weight of the upper body is 30 kg and using the formula for kinetic energy, this spring has a maximum capacity to absorb 33 $\mathrm{J}$ of energy and can convert the required kinetic energy for the upper body to move at a velocity of approximately $1.5$ $\mathrm{m}$/$\mathrm{s}$ into elastic potential energy, which is sufficient for the purpose of studying walking. The choice of spring must also satisfy the limitations of the motor. The stiffness will be adjusted to correspond to the muscle stiffness during running and jumping in the future, but at this stage, we have chosen the highest possible stiffness to reduce the difficulty of control.

3.3. Ankle

A parallel mechanism is often used to drive the ankle of a biped robot, as it can provide more power and greater stability than serial kinematics. Our new biped robot uses a novel cable-driven system to drive the ankle, which is based on a parallel mechanism. Two motors, mounted in the thigh, as shown in Figure 4, are used to drive the two DoFs of the ankle. Each motor output connects two steel wires, where one wire drives the positive rotation, and the other drives the negative rotation, as shown in Figure 4. These wires are marked with red and black lines, respectively, and are transmitted through the knees by several fixed pulleys. One end of each wire is connected to a motor, while the other end is connected to an extension spring. The two wires drive a linear slider via a moving pulley. The stiffness of each extension spring in the ankle drive system is 64 $\mathrm{N}$/$\mathrm{m}$$\mathrm{m}$, and it can produce a maximum force of 760 $\mathrm{N}$. The selection of the ankle spring stiffness is based on the same methodology utilized for the selection of the knee spring stiffness.

Figure 5 shows the working principle of this cable-driven system. One wire connects the motor with the one side of the spring and passes through the movable pulley. The other side of the spring is fixed on the calf. The movable pulley is mounted on the slider. When the motor rotates in positive direction, the slider can move up. According to the characteristics of the moving pulley, the output force can be doubled. Because the wire can transmit only pulling forces, a second wire is necessary to drive the slider in a negative direction. The right hand side of Figure 5 shows the path of the second wire. This wire must pass through several fixed pulleys and a movable pulley of the slider. When the motor rotates in the negative direction, the slider moves down.

Two linear sliders on the calf are used to drive the ankle. Each linear slider is connected to a link by a rod end. The other end of the link connects to the foot via a universal joint. With the cooperation of two sliders, one ankle can be driven. This design provides a solution that combines a wire-driven system and a parallel mechanism to actuate the ankle. We can mount the motor that drives the ankle joint on the thigh using the wire-driven system. The output force of the motor on the parallel mechanism can be doubled by using movable pulleys, and furthermore, a single motor with a cable can drive the positive and negative directions of each input of the parallel mechanism.

3.4. Foot

Most biped robots use flat feet. When such a biped robot walks, the floating foot must be parallel to the ground. Our proposed biped robot uses a movable toe and heel. With an additional degree of freedom of the toes and heels, the landing and lifting of the floating foot can be improved. Springs are used on the toe and heel to reduce the shock effect when the floating foot contacts the ground. Figure 6 shows the foot of the FORREST. In addition to the elastic toe and heel, rubber mats were installed under the toe and heel. The rubber mat can dampen the shock when a foot touches the ground. A 6-axis force/torque sensor will be installed between the foot and ankle. With the measurement of this sensor, the ZMP of the robot can be calculated.

4. Kinematic/Dynamic Analysis and Dynamic Model

In this section, we first introduce the kinematic parameters of FORREST. The hip joints are driven directly by actuators, and their kinematic and dynamic performance is equivalent to that of the motor. The knee and ankle joints are actuated by lead screws and parallel mechanisms, and we modeled and analyzed their kinematic and dynamic performance. Additionally, we introduce the dynamic model of FORREST.

4.1. Kinematic and Dynamic Parameter of FORREST

We used a floating-base tree-structure model to define the kinematic and dynamic parameters of FORREST.

Table 4. Kinematic and dynamic parameters of FORREST.

Link Id	Child Id	Parent Id	${\mathit{T}}_{\mathit{i}}^{\mathit{p}}$	m [kg]	r [m]	${\mathit{I}}_{\mathit{d}}$ [kg m2]
1	2	0	${\mathit{T}}_{roty}(\pi /2)$
2	3	1	${\mathit{T}}_{rotx}(-\pi /2)$
3	4	2	${\mathit{T}}_{roty}(-\pi /2)$
4	5	3	${\mathit{T}}_{rotx}(\pi /2)$
5	6	4	${\mathit{T}}_{roty}(-\pi /2)$
6	7, 13	5	${\mathit{T}}_{rotx}(-\pi /2){\mathit{T}}_{rotz}\left(\pi \right)$	$4.45$	${\left[\begin{array}{ccc}0& 0& 0.44\end{array}\right]}^T$	$\left[\begin{array}{ccc}0.23& 0.18& 0.06\end{array}\right]$
7	8	6	${\mathit{T}}_{trans}\left([0,0.12,0]\right)$	$2.26$	${\left[\begin{array}{ccc}-0.13& 0& 0.02\end{array}\right]}^T$	$\left[\begin{array}{ccc}0.01& 0.06& 0.05\end{array}\right]$
8	9	7	${\mathit{T}}_{roty}(-\pi /2)$	$4.22$	${\left[\begin{array}{ccc}0& 0& 0\end{array}\right]}^T$	$\left[\begin{array}{ccc}0.01& 0.01& 0.01\end{array}\right]$
9	10	8	${\mathit{T}}_{rotx}(\pi /2)$	$4.06$	${\left[\begin{array}{ccc}-0.21& 0& 0\end{array}\right]}^T$	$\left[\begin{array}{ccc}0.01& 0.22& 0.22\end{array}\right]$
10	11	9	${\mathit{T}}_{trans}\left([0,0,-0.405]\right)$	$1.95$	${\left[\begin{array}{ccc}-0.17& 0& 0\end{array}\right]}^T$	$\left[\begin{array}{ccc}0& 0.1& 0.1\end{array}\right]$
11	12	10	${\mathit{T}}_{trans}\left([0,0,-0.405]\right)$	$0.08$	${\left[\begin{array}{ccc}0& 0& 0\end{array}\right]}^T$	$\left[\begin{array}{ccc}0& 0& 0\end{array}\right]$
12		11	${\mathit{T}}_{rotx}(\pi /2)$	$1.46$	${\left[\begin{array}{ccc}-0.07& 0& 0.01\end{array}\right]}^T$	$\left[\begin{array}{ccc}0.01& 0.02& 0.01\end{array}\right]$
13	14	6	${\mathit{T}}_{trans}\left([0,-0.12,0]\right)$	$2.26$	${\left[\begin{array}{ccc}-0.13& 0& 0.02\end{array}\right]}^T$	$\left[\begin{array}{ccc}0.01& 0.06& 0.05\end{array}\right]$
14	15	13	${\mathit{T}}_{roty}(-\pi /2)$	$4.22$	${\left[\begin{array}{ccc}0& 0& 0\end{array}\right]}^T$	$\left[\begin{array}{ccc}0.01& 0.01& 0.01\end{array}\right]$
15	16	14	${\mathit{T}}_{rotx}(\pi /2)$	$4.06$	${\left[\begin{array}{ccc}-0.21& 0& 0\end{array}\right]}^T$	$\left[\begin{array}{ccc}0.01& 0.22& 0.22\end{array}\right]$
16	17	15	${\mathit{T}}_{trans}\left([0,0,-0.405]\right)$	$1.95$	${\left[\begin{array}{ccc}-0.17& 0& 0\end{array}\right]}^T$	$\left[\begin{array}{ccc}0& 0.01& 0.01\end{array}\right]$
17	18	16	${\mathit{T}}_{trans}\left([0,0,-0.405]\right)$	$0.08$	${\left[\begin{array}{ccc}0& 0& 0\end{array}\right]}^T$	$\left[\begin{array}{ccc}0& 0& 0\end{array}\right]$
18		17	${\mathit{T}}_{rotx}(\pi /2)$	$1.46$	${\left[\begin{array}{ccc}-0.07& 0& 0.01\end{array}\right]}^T$	$\left[\begin{array}{ccc}0.01& 0.02& 0.01\end{array}\right]$

describes these parameters of FORREST. Each link had its own ID, child ID, and parent ID to define the tree structure. The links of the driver system of the knee and ankle are not included. Links 1–6 are the translational (links 1–3) and rotational (links 4–6) joints of the floating base. Links 7–12 build the left leg, and links 13–18 build the right leg. ${\mathit{T}}_i^p$ defines the homogeneous transformation matrix of link i with respect to the coordinate frame of the parent link p. Compared to the traditional Denavit–Hartenberg (DH) parameter, we can define any kinematics chain by a sequence of homogeneous transformation matrices. In addition, link i rotates about axis i, which is identical to the modified DH parameter. The vector $\mathit{r}{c_i}^i$ is the CoM of the link expressed in its own coordinate frame. ${\mathit{I}}_d$ contains the diagonal parameters of the inertia matrix. The off-diagonal elements of the inertia matrix are not described here. The kinematic and dynamic parameters were used for the simulation and control design, as described in the next section.

4.2. Knee

Figure 7 shows the dimensions of the knee drive system. Points A and B are two free joints that connect the ball screw unit with the thigh and calf, respectively. It was assumed that the origin of the coordinate system was fixed at the knee (point O), and the thigh was fixed in the coordinate system along the y-axis. When the knee joint $q_4$ rotated, the calf could move in the coordinate system. The picture shows the original position of points A and B.

$${\mathit{x}}_A=\left[\begin{array}{c}0.235\\ 0.065\end{array}\right],{\mathit{x}}_B=\left[\begin{array}{c}-0.08\\ 0.065\end{array}\right],$$

(1)

when $q_4=0$. When the knee joint rotated, the point B moved and was calculated by

$${\mathit{x}}_{B^{*}}=\mathit{R}\left(q_4\right){\mathit{x}}_B=\left[\begin{array}{cc}cos\left(q_4\right)& -sin\left(q_4\right)\\ sin\left(q_4\right)& cos\left(q_4\right)\end{array}\right]{\mathit{x}}_B.$$

(2)

According to the structure of the knee in Figure 7, we obtained the angle of the two passive joints $q_{19}$ and $q_{20}$ at both ends of the linear drive system, with

$$q_{19}=atan2(x_A-x_{B^{*}},y_A-y_{B^{*}}),$$

(3)

$$q_{20}=q_{19}-q_4.$$

(4)

By using Heron’s formula, the height h of the triangle $AO{B}^{*}$ was computed with

$$\begin{array}{cc}\hfill h=& {\frac{2\eta }{\parallel \overrightarrow{A{B}^{*}}\parallel }}_2,\hfill \\ \hfill \eta =& \sqrt{s(s-{\parallel \overrightarrow{OA}\parallel }_2)(s-{\parallel \overrightarrow{O{B}^{*}}\parallel }_2)(s-{\parallel \overrightarrow{A{B}^{*}}\parallel }_2)},\hfill \\ \hfill s=& \frac{{\parallel \overrightarrow{OA}\parallel }_2+{\parallel \overrightarrow{O{B}^{*}}\parallel }_2+{\parallel \overrightarrow{A{B}^{*}}\parallel }_2}{2}.\hfill \end{array}$$

(5)

The relationship of the speed ${\dot{q}}_k$ and the torque ${\tau }_k$ between the motor and the knee joint and the compliance $k_k$ of the knee were:

$${\dot{q}}_k=\frac{n_{m}p}{60h},$$

(6)

$${\tau }_k=\frac{2\pi {\tau }_{m}h}{p}.$$

(7)

$$k_k=\frac{\Delta {\tau }_{ks}}{\Delta q_4}=\frac{k_s·\Delta {\parallel \overrightarrow{A{B}^{*}}\parallel }_{2}h}{J_{AB}^{-1}·\Delta {\parallel \overrightarrow{A{B}^{*}}\parallel }_2}=\frac{k_s·h}{J_{AB}^{-1}},$$

(8)

$$J_{AB}=\frac{d{\parallel \overrightarrow{A{B}^{*}}\parallel }_2}{d{q}_4}.$$

(9)

Herein, $n_m$ and ${\tau }_m$ represent the motor speed and torque, respectively, and p is the screw pitch. By using these values, we determined the kinematic and dynamic performance of the knee joint. ${\tau }_{ks}$ denotes the torque caused by the spring, $k_s$ is the spring stiffness, and $J_{AB}$ is the Jacobian. Figure 8 illustrates the relationship between the maximum joint speed/torque/compliance and the knee joint angle. It is obvious that the joint torque and stiffness were highest at ${60}^{\circ }$, where the joint velocity was minimal.

4.3. Ankle

The ankle joint has two degrees of freedom and is driven by a parallel mechanism. Figure 9 illustrates the dimensions of the parallel mechanism. Points A and B represent the ends of a rod that connects a slider with point A on the calf and the foot with point B. $A_1$ and $B_1$ are the endpoints of the left rod, while $A_2$ and $B_2$ correspond to the right rod. The length of the rod was $0.185$ $\mathrm{m}$. We assumed that the origin of a coordinate system was fixed on the ankle (point O) to obtain the kinematic and dynamic relationship between the ankle and the motors. The calf was fixed in this coordinate system, while the foot was movable. The original positions of $A_{1/2}$ and $B_{1/2}$ were

$$A_1=\left[\begin{array}{c}-0.015\\ 0.04\\ -0.06\end{array}\right],B_1=\left[\begin{array}{c}{x}_{B_1}\\ 0.04\\ -0.0345\end{array}\right],A_2=\left[\begin{array}{c}-0.015\\ -0.04\\ -0.06\end{array}\right],B_2=\left[\begin{array}{c}{x}_{B_2}\\ -0.04\\ -0.0345\end{array}\right].$$

(10)

When the ankle joints $q_{11}$ and $q_{12}$ rotated, the new position of points $A_1$ were calculated by

$${\mathit{A}}_1^{*}={\mathit{R}}_y\left(q_{11}\right){\mathit{R}}_z\left(q_{12}\right){\mathit{A}}_1=\left[\begin{array}{ccc}cos\left(q_{11}\right)& 0& sin\left(q_{11}\right)\\ 0& 1& 0\\ -sin\left(q_{11}\right)& 0& cos\left(q_{11}\right)\end{array}\right]\left[\begin{array}{ccc}cos\left(q_{12}\right)& -sin\left(q_{12}\right)& 0\\ sin\left(q_{12}\right)& cos\left(q_{12}\right)& 0\\ 0& 0& 1\end{array}\right]{\mathit{A}}_1.$$

(11)

As point $B_1$ can only slide along the x-axis, the y and z positions of $B_1$ are constant. Using the Euclidean norm, the new x position of $B_1$ was determined as follows:

$$x_{B_1^{*}}=x_{A_1^{*}}+\sqrt{L^2-{(y_{A_1^{*}}-y_{B_1^{*}})}^2-{(z_{A_1^{*}}-z_{B_1^{*}})}^2}.$$

(12)

The calculation of $z_{B_2^{*}}$ was identical. All the passive joints of the parallel mechanism were obtained with

$$\begin{array}{cccccc}\hfill q_{25}& =0.405-x_{b_1}^{*},\hfill & \hfill q_{26}& =asin\left(\frac{y_{b1}^{*}-y_{a1}}{L·cos\left(q_{27}\right)}\right),\hfill & \hfill q_{27}& =asin\left(\frac{x_{b1}^{*}-x_{a1}}{L}\right),\hfill \\ \hfill q_{28}& =0.405-x_{b_2}^{*},\hfill & \hfill q_{29}& =asin\left(\frac{y_{b2}^{*}-y_{a2}}{L·cos\left(q_{30}\right)}\right),\hfill & \hfill q_{30}& =asin\left(\frac{x_{b2}^{*}-x_{a2}}{L}\right).\hfill \end{array}$$

(13)

Then, we calculated the Jacobian matrix, which describes the relationship between the ankle joints and sliders.

$$\mathit{J}=\left[\begin{array}{cc}\frac{\partial q_{25}}{\partial q_{11}}& \frac{\partial q_{25}}{\partial q_{12}}\\ \frac{\partial q_{28}}{\partial q_{11}}& \frac{\partial q_{28}}{\partial q_{12}}\end{array}\right].$$

(14)

The relationship of the speed and torque between the motor and knee joint was

$$\left[\begin{array}{c}{\dot{q}}_{11}\\ {\dot{q}}_{12}\end{array}\right]={\mathit{J}}^{-1}\left[\begin{array}{c}2\pi r_1{n}_{m_1}\\ 2\pi r_2{n}_{m_2}\end{array}\right],$$

(15)

$$\left[\begin{array}{c}{\tau }_{11}\\ {\tau }_{12}\end{array}\right]={\mathit{J}}^T\left[\begin{array}{c}2{\tau }_{m_1}/r_1\\ 2{\tau }_{m_2}/r_2\end{array}\right],$$

(16)

$$k_{11}=\frac{\Delta {\tau }_{11s}}{\Delta q_{11}}=\frac{{\mathit{J}}_{[1,-]}^T{\mathit{k}}_s\mathit{J}{\left[\begin{array}{cc}\Delta q_{11}& 0\end{array}\right]}^T}{\Delta q_{11}}.$$

(17)

Herein, $n_{m_1}$, $n_{m_2}$, ${\tau }_{m_1}$, and ${\tau }_{m_2}$ are the speeds and torques of the motors, respectively. ${\dot{q}}_{11}$, ${\dot{q}}_{12}$, ${\tau }_{11}$, and ${\tau }_{12}$ are the ankle joint speeds and torques, respectively. $r_1$ and $r_2$ are the radii of the cable drums, which are driven by motors. $k_{11}$ is the stiffness of the ankle joint 11. The subscript of the Jacobi matrix ${\mathit{J}}_{[1,-]}^T$ represents the elements of the first row. $k_s$ is the stiffness of the spring matrix. Assuming that joint 12 does not move when the stiffness of joint 11 is calculated, the method of calculating the joint 12 stiffness is the same. Figure 10, Figure 11 and Figure 12 show the relationship between the maximum joint speed/torque and the angle of the ankle joints.

4.4. Dynamic Model

The dynamic model of FORREST is a floating base model with closed chains and elastic joints. We neglected the cables and pulleys of the cable-driven parallel mechanism to simplify the modeling. The drive system of the ankle joints was simplified as two virtual linear drivers mounted on the calf, and the motors and sliders were connected by virtual springs. The virtual drivers and springs can be found in Figure 9. We used the Newton–Euler method to calculate the dynamics of our biped. To solve the inverse dynamics, we first calculated the velocity and acceleration of all the passive joints that were not equipped with encoders. Using (18) and (19), we obtained the Jacobian matrices of the passive joints of the left leg. The following formulas with the joint numbers in this section were used to calculate the dynamics of the left leg. The calculation of the right leg was identical.

$$J_i=\frac{\partial q_i}{\partial q_{10}},{\dot{J}}_i=\frac{\partial q_i}{\partial q_{10}\partial q_{10}}{\dot{q}}_{10}andi=19,20,$$

(18)

$$\begin{array}{cccc}\hfill J_i& =\left[\begin{array}{cc}\frac{\partial q_i}{\partial q_{11}}& \frac{\partial q_i}{\partial q_{12}}\end{array}\right]\hfill & \hfill ,i& =25-30,\hfill \\ \hfill {\dot{J}}_i& =\left[\begin{array}{cc}\frac{\partial q_i}{\partial q_{11}\partial q_{11}}\dot{q_{11}}+\frac{\partial q_i}{\partial q_{11}\partial q_{12}}\dot{q_{12}}& \frac{\partial q_i}{\partial q_{12}\partial q_{11}}\dot{q_{11}}+\frac{\partial q_i}{\partial q_{12}\partial q_{12}}\dot{q_{12}}\end{array}\right]\hfill & \hfill ,i& =25-30.\hfill \end{array}$$

(19)

The velocity and acceleration of all the passive joints were calculated by

$$\begin{array}{cccc}\hfill {\dot{q}}_i& =J_i\left(q_{10}\right)\dot{q_{10}}\hfill & \hfill ,i& =19,20,\hfill \\ \hfill {\ddot{q}}_i& ={\dot{J}}_i(q_{10},{\dot{q}}_{10}){\dot{q}}_{10}+J_i\left(q_{10}\right){\ddot{q}}_{10}\hfill & \hfill ,i& =19,20,\hfill \end{array}$$

(20)

$$\begin{array}{cccc}\hfill {\dot{q}}_i& =J_i(q_{11},q_{12})\left[\begin{array}{c}{\dot{q}}_{11}\\ {\dot{q}}_{12}\end{array}\right]\hfill & \hfill ,i& =25-30,\hfill \\ \hfill {\ddot{q}}_i& ={\dot{J}}_i(q_{11},q_{12})\left[\begin{array}{c}{\dot{q}}_{11}\\ {\dot{q}}_{12}\end{array}\right]+J_i(q_{11},q_{12})\left[\begin{array}{c}{\ddot{q}}_{11}\\ {\ddot{q}}_{12}\end{array}\right]\hfill & \hfill ,i& =25-30.\hfill \end{array}$$

(21)

Then, we calculated the angular and linear velocity and acceleration of each link coordinate and CoM in the world coordinate system by (22)–(28). Herein, ${\mathit{\omega }}_i$, ${\dot{\mathit{\omega }}}_i$, ${\mathit{v}}_i$, and ${\mathit{a}}_i$ are the angular velocity/acceleration and linear velocity/acceleration of link coordinate i, respectively. ${\dot{\mathit{\omega }}}_{p_i}$ and ${\mathit{a}}_{p_i}$ are the accelerations of the parent link of link i. For a floating base biped model, the base acceleration is related to the foot contact force/moment. The contact force/moment was measured by sensors that were mounted on the ankle. The acceleration of the base was still unknown. Therefore, the equations of the link acceleration needed to be divided into a term with base accelerations and a term without base accelerations (${\dot{\mathit{\omega }}}_i^{*}$ and ${\mathit{a}}_i^{*}$). ${\dot{\mathit{\omega }}}_b$ and ${\mathit{a}}_b$ are the accelerations of the floating base coordinate (link 6). The operator $S\left(\right)$ converts a vector into a skew-symmetric matrix and simplifies the calculation of a cross product.

$${\mathit{\omega }}_i=\left\{\begin{array}{cc}{\mathit{\omega }}_{p_i}+{\dot{q}}_i{\mathit{z}}_i\hfill & ,\mathrm{revolute},\hfill \\ {\mathit{\omega }}_{p_i}\hfill & ,\mathrm{prismatic},\hfill \end{array}\right.$$

(22)

$${\dot{\mathit{\omega }}}_i={\dot{\mathit{\omega }}}_b+{\dot{\mathit{\omega }}}_i^{*},$$

(23)

$${\dot{\mathit{\omega }}}_i^{*}=\left\{\begin{array}{cc}{\dot{\mathit{\omega }}}_{p_i}+{\ddot{q}}_i{\mathit{z}}_i+\mathit{S}\left({\mathit{\omega }}_{p_i}\right)\left({\dot{q}}_i{\mathit{z}}_i\right)\hfill & ,\mathrm{revolute},\hfill \\ {\dot{\mathit{\omega }}}_{p_i}^{*}\hfill & ,\mathrm{prismatic},\hfill \end{array}\right.$$

(24)

$${\mathit{v}}_i=\left\{\begin{array}{cc}{\mathit{v}}_{p_i}\hfill & ,\mathrm{revolute},\hfill \\ {\mathit{v}}_{p_i}+{\dot{q}}_i{\mathit{z}}_i\hfill & ,\mathrm{prismatic},\hfill \end{array}\right.$$

(25)

$${\mathit{a}}_i={\mathit{a}}_b-\mathit{S}\left(r_{i,b}\right){\dot{\mathit{\omega }}}_b+{\mathit{a}}_i^{*},$$

(26)

$${\mathit{a}}_i^{*}=\left\{\begin{array}{cc}{\mathit{a}}_{p_i}^{*}-\mathit{S}\left({\mathit{r}}_{i,p_i}\right){\dot{\mathit{\omega }}}_{p_i}^{*}+\mathit{S}\left({\mathit{\omega }}_{p_i}\right)\mathit{S}\left({\mathit{\omega }}_{p_i}\right){\mathit{r}}_{i,p_i}\hfill & ,\mathrm{revolute},\hfill \\ {\mathit{a}}_{p_i}^{*}-\mathit{S}\left({\mathit{r}}_{i,p_i}\right){\dot{\mathit{\omega }}}_{p_i}^{*}+\mathit{S}\left({\mathit{\omega }}_{p_i}\right)\mathit{S}\left({\mathit{\omega }}_{p_i}\right){\mathit{r}}_{i,p_i}+{\ddot{q}}_i{\mathit{z}}_i+2\mathit{S}\left({\mathit{\omega }}_{p_i}\right)\left({\dot{q}}_i{\mathit{z}}_i\right)\hfill & ,\mathrm{prismatic},\hfill \end{array}\right.$$

(27)

$${\dot{\mathit{\omega }}}_{c_i}={\dot{\mathit{\omega }}}_i.$$

(28)

We needed to calculate the acceleration of the CoM of all links to obtain the Newton–Euler equations. The angular acceleration of the CoM was the same as that of the link coordinate. ${\mathit{a}}_{c_i}$ represents the linear acceleration of the CoM of link i.

$$\begin{array}{cc}\hfill {\mathit{a}}_{c_i}& ={\mathit{a}}_i-\mathit{S}\left({\mathit{r}}_{c_i,i}\right){\dot{\mathit{\omega }}}_i+\mathit{S}\left({\mathit{\omega }}_i\right)\mathit{S}\left({\mathit{\omega }}_i\right){\mathit{r}}_{c_i,i}\hfill \\ \hfill & ={\mathit{a}}_b-\mathit{S}\left({\mathit{r}}_{c_i,b}\right){\dot{\mathit{\omega }}}_b+{\mathit{a}}_i^{*}-\mathit{S}\left({\mathit{r}}_{c_i,i}\right){\dot{\mathit{\omega }}}_i^{*}+\mathit{S}\left({\mathit{\omega }}_i\right)\mathit{S}\left({\mathit{\omega }}_i\right){\mathit{r}}_{c_i,i}.\hfill \end{array}$$

(29)

When the accelerations were calculated, we obtained the Newton–Euler equations of each link. ${\mathit{F}}_i^k$ and ${\mathit{M}}_i^k$ are the forces and moments exerted on link i by the parent and child of link i, respectively. Here, because of the elastic drive system, the springs connected link $2, 0$ and $2, 1$, link $2, 5$ and $3, 7$, and link $2, 8$ and $3, 8$. The spring force can be considered an interaction force between two links and was included in the Newton–Euler equations. Usually, when the robot kinematics are calculated, the deformation of the spring can be obtained, and the spring force can be calculated. On the other hand, we can also compute it with the Newton–Euler method.

$$m_i({\mathit{a}}_b-\mathit{S}\left({\mathit{r}}_{c_i,b}\right){\dot{\mathit{\omega }}}_b+{\mathit{a}}_i^{*}-\mathit{S}\left({\mathit{r}}_{c_i,i}\right){\dot{\mathit{\omega }}}_i^{*}+\mathit{S}\left({\mathit{\omega }}_i\right)\mathit{S}\left({\mathit{\omega }}_i\right){\mathit{r}}_{c_i,i})=\sum _k{\mathit{F}}_i^k+m_i\mathit{g},$$

(30)

$${\mathit{I}}_i{\dot{\mathit{\omega }}}_b+{\mathit{I}}_i{\dot{\mathit{\omega }}}_i^{*}+\mathit{S}\left({\mathit{\omega }}_i\right){\mathit{I}}_i{\mathit{\omega }}_i=-\sum _k\mathit{S}\left({\mathit{r}}_{c_i,k}\right){\mathit{F}}_i^k+\sum _k{\mathit{M}}_i^k.$$

(31)

The Newton–Euler Equation of each link was transformed into a matrix form (32). The acceleration of the floating base, forces, and torques on the link were unknown variables. $\mathit{E}$ is a $3\times 3$ unit matrix. ${\mathit{F}}_i$ and ${\mathit{M}}_i$ represent the vector of all applied forces ${\mathit{F}}_i^k$ and moments ${\mathit{M}}_i^k$ by the parent and child k. ${\mathit{S}}_{{\mathit{r}}_{c_i,k}}$ is a matrix, which contains all the skew-symmetric matrices $\mathit{S}\left({\mathit{r}}_{c_i,k}\right)$ for all the applied forces ${\mathit{F}}_i^k$ in (31).

$$\left[\begin{array}{cccc}{m}_i\mathit{E}& -m_i\mathit{S}\left({\mathit{r}}_{c_i,b}\right)& -\mathit{E}& \mathbf{0}\\ \mathbf{0}& {\mathit{I}}_i& {\mathit{S}}_{{\mathit{r}}_{c_i,k}}& -\mathit{E}\end{array}\right]\left[\begin{array}{c}{\mathit{a}}_b\\ {\mathit{\omega }}_b\\ {\mathit{F}}_i\\ {\mathit{M}}_i\end{array}\right]=\left[\begin{array}{c}{m}_i(-{\mathit{a}}_i^{*}+\mathit{S}\left({\mathit{r}}_{c_i,i}\right){\dot{\mathit{\omega }}}_i^{*}-\mathit{S}\left({\mathit{\omega }}_i\right)\mathit{S}\left({\mathit{\omega }}_i\right){\mathit{r}}_{c_i,i}+\mathit{g})\\ -{\mathit{I}}_i{\dot{\mathit{\omega }}}_i^{*}-\mathit{S}\left({\mathit{\omega }}_i\right){\mathit{I}}_i{\mathit{\omega }}_i\end{array}\right].$$

(32)

Furthermore, the constraints of the passive joints needed to be considered. The moment along the passive joint axis was zero. It can also be described that the interaction moment between link i and its parent link $p_i$ along the joint axis was zero. All moments on the two spherical joints, which connected link 25 (28) and 27 (30), were zero. The force on the two sliders (link 25, 28) along the linear axis can only be produced by the spring, and the other force source (link 10) along the linear axis must be zero according to (35).

$${\mathit{z}}_i^{p_i,T}{\mathit{M}}_i^{p_i}=0,i=10,11,12,19,20,25–30$$

(33)

$${\mathit{M}}_{25}^{27}={\left[\begin{array}{ccc}0& 0& 0\end{array}\right]}^T,{\mathit{M}}_{28}^{30}={\left[\begin{array}{ccc}0& 0& 0\end{array}\right]}^T$$

(34)

$${\mathit{z}}_i^{10,T}{\mathit{F}}_i^{10}=0,i=25,28.$$

(35)

By combining the Newton–Euler Equation (32) and the constraints (33)–(35) of all links, we obtained the matrix form of the system of linear equations for the biped (36). Here, $\mathit{F}$ and $\mathit{M}$ are the vectors of all the interaction forces and moments, including the spring forces. The acceleration of the base and all forces/moments were calculated by solving (36). However, due to the presence of the elastic knee and ankle joints, these joints were under-actuated and could not be directly driven by motors. The angular acceleration of the joint was related to the spring force, not the motor output force. Therefore, this inverse dynamics method could not be directly used for inverse dynamics-based control. When setting the desired joint acceleration, the spring force did not match. To solve this problem, it was necessary to extend the dynamics to a fourth order system (including jerk and snap), which will be covered in further work.

$$\mathit{H}\left[\begin{array}{c}{\mathit{a}}_b\\ {\mathit{\omega }}_b\\ \mathit{F}\\ \mathit{M}\end{array}\right]=\mathit{Y}$$

(36)

5. Extended CP-Based Walking Pattern Generator

We used an LIP-based walking pattern generator to generate the walking gait trajectory. A common LIP model has a resultant external force by the ground ${\mathit{F}}_{ext}$, which must be along the LIP. The vertical component of the external force is cancelled out by gravity, which keeps the height $z_c$ of the CoM constant. In our work, we aimed to extend the walking pattern generator by adding a variable height $z_c$ to the CoM. A similar result using a different method was introduced in [35]. We assumed that the LIP was subject to an external force $F_{az}$ that could produce vertical acceleration in addition to the external forces $F_{ax}$, $F_{ay}$, and $F_{ag}$, which were the components along the LIP of the resultant external force. The new dynamic equation of the LIP was

$$\begin{array}{cc}\hfill m{\ddot{x}}_c& =\frac{x_c-p_x}{r}f,\hfill \\ \hfill m{\ddot{y}}_c& =\frac{y_c-p_y}{r}f,\hfill \\ \hfill m{\ddot{z}}_c& =\frac{z_c-p_z}{r}f+\frac{z_c}{r}f-mg,\hfill \end{array}$$

(37)

and the scalar of resultant force f satisfied

$$\frac{z_c}{r}f=mg.$$

(38)

Herein, the three components of $\mathit{x}=\left[x_c{y}_c{z}_c\right]$ were the position of the CoM. $\mathit{p}={\left[p_x{p}_y{p}_z\right]}^T$ was the position of the extended ZMP (eZMP), which had the same x–y position as the common 2D ZMP of the LIP and an additional z position. Its position on the z-axis was symmetrical to the endpoint of the vector $F_{az}$ about a plane passing through the CoM and perpendicular to the z-axis. Figure 13 shows the forces on the LIP Model. The scalar f represents the resultant force of $F_{ax}$, $F_{ay}$, and $F_{ag}$ and should be along the LIP. The scalar r represents the length of the vector from ZMP to the CoM. With the additional force $F_{az}$, the resultant force on the LIP moves the CoM upward.

By combining (37) and (38), we obtained the new dynamic equation of the LIP

$$\ddot{\mathit{x}}=\frac{1}{m}({\mathit{F}}_{ext}+{\mathit{F}}_g)=\frac{1}{m}\left(\left[\begin{array}{c}{F}_{ax}\\ F_{ay}\\ F_{az}+F_{ag}\end{array}\right]+\left[\begin{array}{c}0\\ 0\\ -mg\end{array}\right]\right)={\omega }^2(\mathit{x}-\mathit{p})=\frac{g}{z_c}\left[\begin{array}{c}{x}_c-p_x\\ y_c-p_y\\ z_c-p_z\end{array}\right],$$

(39)

and the components of external force satisfied

$$\frac{F_{ax}}{x_c-p_x}=\frac{F_{ay}}{y_c-p_y}=\frac{F_{ag}}{z_c}=\frac{mg}{z_c}.$$

(40)

The time constant $\omega$ was $\sqrt{g/z_c}$. In the following, we describe the use of the 3D capture point $\mathit{\xi }$ dynamics to generate our desired trajectory. The 3D CP is also called the divergent component of motion (DCM). The introduction of the CP and DCM can be found in [34,35].

$$\dot{\mathit{x}}=-\omega (\mathit{x}-\mathit{\xi })$$

(41)

$$\dot{\mathit{\xi }}=\omega (\mathit{\xi }-\mathit{p}).$$

(42)

The dynamics of the CoM (41) constitute a stable first-order open loop system, which ensures that the position of the CoM will always converge to the CP. On the other hand, the dynamics of the CP (42) constitute an unstable first-order open loop system. Based on the behavior of the CoM, CP, and ZMP, it is sufficient to control the ZMP to track the desired CP trajectory. However, a common walking pattern generator utilizes a constant ZMP position, leading to instantaneous ZMP exchanges and sudden changes in CoM acceleration, which can result in an unsmooth velocity trajectory. To overcome this disadvantage, we introduced three additional phases between two ZMPs. Figure 14 illustrates an example of the curve of the ZMP and CP. Each step comprised three phases: two exchange phases, lasting $T_1$ and $T_3$, and one constant ZMP phase, lasting $T_2$. We assumed that the initial and final ZMPs, $P_0$ and $P_2$, were the same as the initial and final CPs, ${\xi }_0$ and ${\xi }_e$, in each step. Our objective was to search for a desired ZMP $P_1$ to satisfy the above assumptions.

The function of the ZMP is:

$$P\left(t\right)=\left\{\begin{array}{cc}\hfill k_{p_1}t+p_0& t_0\le t<t_1\hfill \\ \hfill p_1& t_1\le t\le t_2\hfill \\ \hfill k_{p_2}t+p_1& t_2<t\le t_3\hfill \end{array}\right.$$

(43)

$$k_{p_1}=\frac{p_1-p_0}{T_1},k_{p_2}=\frac{p_2-p_1}{T_3}.$$

(44)

By solving the ordinary differential equation, we obtained the position equation of the CoM about time

$$x\left(t\right)=\left\{\begin{array}{c}(x_{0_1}-p_0)cosh\left(\omega t\right)+\frac{{\dot{x}}_{0_1}-k_{p_1}}{\omega }sinh\left(\omega t\right)+k_{p_1}t+P_0{t}_0\le t\le t_1\hfill \\ (x_{0_2}-p_1)cosh\left(\omega t\right)+\frac{{\dot{x}}_{0_2}}{\omega }sinh\left(\omega t\right)+p_1{t}_1<t\le t_2\hfill \\ (x_{0_3}-p_1)cosh\left(\omega t\right)+\frac{{\dot{x}}_{0_3}-k_{p_2}}{\omega }sinh\left(\omega t\right)+k_{p_2}t+p_1{t}_2<t\le t_3,\hfill \end{array}\right.$$

(45)

and the CP equation

$$\xi \left(t\right)=\left\{\begin{array}{c}{\xi }_{0_1}{e}^{\omega t}+\left(P_0+\frac{P_1-P_0}{\omega T_1}\right)\left(1-e^{\omega t}\right)+\frac{P_1-P_0}{T_1}t{t}_0\le t\le t_1\hfill \\ {\xi }_{0_2}{e}^{\omega t}+P_1(1-e^{\omega t})t_1<t\le t_2\hfill \\ {\xi }_{0_3}{e}^{\omega t}+\left(P_1+\frac{P_2-P_1}{\omega T_3}\right)\left(1-e^{\omega t}\right)+\frac{P_2-P_1}{T_3}t{t}_2<t\le t_3.\hfill \end{array}\right.$$

(46)

Herein, $x_{0_1}$, $x_{0_2}$, and $x_{0_3}$ are the initial positions of the CoM for each phase. ${\dot{x}}_{0_1}$, ${\dot{x}}_{0_2}$ and ${\dot{x}}_{0_3}$ are its initial velocity. ${\xi }_{0_1}$, ${\xi }_{0_2}$ and ${\xi }_{0_3}$ are the initial position of the CP. It is known that ${\xi }_{0_1}={\xi }_0$ and $\xi \left(t_3\right)={\xi }_e=P_2$. Because the final CP of first phase is the initial CP of the second phase, and the final CP of the second phase is the initial CP of the third phase, we combined the three formulas of (46) and obtained

$$\begin{array}{c}\hfill {\xi }_e\left(t_3\right)=P_2={\xi }_0{e}^{\omega (T_1+T_2+T_3)}+\frac{P_0}{\omega T_1}\left(e^{\omega (T_1+T_2+T_3)}-e^{\omega (T_2+T_3)}-\omega T_1{e}^{\omega (T_1+T_2+T_3)}\right)\\ \hfill +\frac{P_1}{\omega T_1{T}_3}\left(T_1{e}^{\omega T_3}+T_3{e}^{\omega (T_2+T_3)}-T_3{e}^{\omega (T_1+T_2+T_3)}-T_1\right)+\frac{P_2}{\omega T_3}\left(1-e^{\omega T_3}+\omega T_3\right).\end{array}$$

(47)

Once $P_0$ and $P_2$ are known, we can use (48) to calculate the desired $P_1$.

$$\begin{array}{cc}\hfill P_1=& -\frac{{\xi }_0\omega T_1{T}_3{e}^{\omega (T_1+T_2+T_3)}+P_0{T}_3\left(e^{\omega (T_1+T_2+T_3)}-e^{\omega (T_2+T_3)}-\omega T_1{e}^{\omega (T_1+T_2+T_3)}\right)}{T_1{e}^{\omega T_3}+T_3{e}^{\omega (T_2+T_3)}-T_3{e}^{\omega (T_1+T_2+T_3)}-T_1}\hfill \\ \hfill & -\frac{P_2{T}_1\left(1-e^{\omega T_3}\right)}{T_1{e}^{\omega T_3}+T_3{e}^{\omega (T_2+T_3)}-T_3{e}^{\omega (T_1+T_2+T_3)}-T_1}.\hfill \end{array}$$

(48)

Then, we can obtain the trajectory of the CoM and CP by (45) and (46). Figure 15 shows the curves of the ZMP, CP, CoM, and the velocity curve of the CoM. By using this new method, we obtained a position curve of the CoM with a smooth velocity curve. The variable height of the CoM also became possible. However, for this article, the CoM height only decreased from its initial position to a constant value. We will discuss in future articles the application of variable height in efficient walking, such as walking with extended knees. Furthermore, the ZMP and the CP were the same at the end of each step, ensuring CP stability. Figure 16 shows the ZMP, CP, and CoM curve in the $xy$ plane. It is evident that the trajectories of the CP and ZMP overlapped because the CP and ZMP shared the same starting and ending points for each footstep. At the beginning of each step, the ZMP moved away from the CP towards the $P_1$, generating a suitable repulsive force to accelerate the CP away from ZMP. As the end of each step neared, the ZMP quickly approached the CP and coincided with it at the conclusion.

6. Control Strategy

In this section, we introduce the control strategy for the biped FORREST. Currently, the control strategy for balancing the robot is under development. We used MATLAB/Simulink to test and validate the performance of the controller, which can be implemented in our biped in the future. Figure 17 shows the control scheme for our biped, which consists of an LIP-based walking pattern generator, a capture point controller, a cartesian PD controller, an optimizer for contact force distribution, and an inverse dynamics module. We defined the step distance and step time, and the walking pattern generator calculated the desired trajectories of the CoM, CP, and feet. The capture point controller calculated the feedforward acceleration of the CoM, while the PD controller obtained the desired force on the CoM according to the desired trajectories. An optimization algorithm was implemented to obtain the optimized force distribution from the desired force. Finally, the inverse dynamics module calculated the desired torque of each joint and sent it to the robot. Firstly, we introduce the walking pattern generator of the robot.

6.1. Capture Point Control

We used the method from [35] to track the CP trajectory. We defined a stable dynamic equation of the error of CP for $k>0$.

$${\dot{\mathit{\xi }}}_d-\dot{\mathit{\xi }}=-k({\mathit{\xi }}_d-\mathit{\xi }).$$

(49)

By substituting (42) into (49), we obtained the controlling ZMP ${\mathit{p}}^c$

$${\mathit{p}}^c={\mathit{p}}^d+(1+k/\omega )(\mathit{\xi }-{\mathit{\xi }}_d).$$

(50)

Finally, we used ${\mathit{p}}^c$ to calculate the desired acceleration of the CoM. The desired acceleration was used as a feedforward input in the Cartesian PD control.

$${\ddot{\mathit{x}}}_c^f={\omega }^2({\mathit{x}}_c^d-{\mathit{p}}^c).$$

(51)

6.2. Cartesian PD Controller

In this section, we introduce the Cartesian PD controller, which is based on [44], the optimization of contact force distribution, and the inverse dynamics that we implemented in our simulation. We used the Cartesian PD controller with feed forward from (51) to keep the robot in balance. The control law is

$$\begin{array}{cc}\hfill \left[\begin{array}{c}{\mathit{a}}_c^d\\ {\mathit{\omega }}_c^d\end{array}\right]& ={\mathit{K}}_{pc}{\tilde{\mathit{x}}}_c+{\mathit{K}}_{dc}{\dot{\tilde{\mathit{x}}}}_c+{\ddot{\mathit{x}}}_c^f\hfill \\ \hfill \left[\begin{array}{c}{\mathit{a}}_f^d\\ {\mathit{\omega }}_f^d\end{array}\right]& ={\mathit{K}}_{pf}{\tilde{\mathit{x}}}_f+{\mathit{K}}_{di}{\dot{\tilde{\mathit{x}}}}_f.\hfill \end{array}$$

(52)

${\tilde{\mathit{x}}}_c$ contains the vector of the position and orientation error of the CoM. ${\tilde{\mathit{x}}}_f$ are the vectors of the errors of the two feet. By using (53), the desired wrench ${\mathit{F}}_c^d$ on the CoM was obtained.

$${\mathit{F}}_c^d=\left[\begin{array}{c}{m}_c{\mathit{a}}_c^d\\ {\mathit{I}}_c{\dot{\mathit{\omega }}}_c^d+{\mathit{\omega }}_c\times \left({\mathit{I}}_c{\mathit{\omega }}_c\right)\end{array}\right]$$

(53)

In order to calculate the force distribution, we defined the optimized wrench on the feet, $\mathit{F}{f}^{opt}$, as shown in (54). This wrench consists of the vertical and horizontal contact forces, $\mathit{f}{l}^{opt}$ and $\mathit{f}{r}^{opt}$, of the left and right feet. The horizontal elements of the contact moments, $\mathit{M}{l}^{opt}$ and $\mathit{M}{r}^{opt}$, were zero. Additionally, $\mathit{p}{l}^{opt}$ and ${\mathit{p}}_r^{opt}$ represent the center of pressure (CoP) of the left and right feet, respectively.

$${\mathit{F}}_f^{opt}=\left[\begin{array}{ccc}{\mathit{F}}_l^{opt}& {\mathit{M}}_l^{opt}& {\mathit{p}}_l^{opt}\\ {\mathit{F}}_r^{opt}& {\mathit{M}}_r^{opt}& {\mathit{p}}_r^{opt}\end{array}\right]=\left[\begin{array}{cccccccc}{f}_{l,x}^{opt}& f_{l,y}^{opt}& f_{l,z}^{opt}& 0& 0& {\tau }_{l,z}^{opt}& p_{l,x}^{opt}& p_{l,y}^{opt}\\ f_{r,x}^{opt}& f_{r,y}^{opt}& f_{r,z}^{opt}& 0& 0& {\tau }_{r,z}^{opt}& p_{r,x}^{opt}& p_{r,y}^{opt}\end{array}\right]$$

(54)

We used the following constrained quadratic optimization problem:

$$\underset{{\mathit{F}}_f^{opt}}{min}({\mathit{F}}_c^d-{\mathit{F}}_c^{opt})\mathit{W}({\mathit{F}}_c^d-{\mathit{F}}_c^{opt})$$

(55)

with

$${\mathit{F}}_c^{opt}=\left[\begin{array}{c}{\mathit{F}}_l^{opt}+{\mathit{F}}_r^{opt}\\ {\mathit{F}}_l^{opt}\times {\mathit{p}}_{l,c}+{\mathit{F}}_r^{opt}\times {\mathit{p}}_{r,c}+{\mathit{M}}_l^{opt}+{\mathit{M}}_r^{opt}\end{array}\right]$$

(56)

and with the constraints for the support foot. The optimized vertical force must be larger than the minimum vertical force. The horizontal forces must be smaller than the friction force, where $\mu$ is the friction factor. The position of the CoP must be within the range of the support polygon S of the feet. The optimized wrench was substituted into (53) to calculate the optimized acceleration of the CoM. By combining the optimized acceleration of the CoM and the desired acceleration of the two feet, we calculated the desired acceleration of the joints using (57).

$$\left[\begin{array}{c}{\mathit{a}}_c^{opt}\\ {\mathit{a}}_l^d\\ {\mathit{a}}_r^d\end{array}\right]=\dot{\mathit{J}}\dot{\mathit{q}}+\mathit{J}{\ddot{\mathit{q}}}^d.$$

(57)

Finally, we substituted the current joint position, velocity, desired acceleration, and optimized wrench of the feet into the floating base inverse dynamics of the biped, which was described in the previous section, to calculate the desired torque ${\mathit{\tau }}^d$ of the joints.

$${\mathit{\tau }}^d=invdyn(\mathit{q},\dot{\mathit{q}},{\ddot{\mathit{q}}}^d,{\mathit{F}}_l^{opt},{\mathit{F}}_r^{opt})$$

(58)

6.3. Results

We implemented our new walking pattern generator in a simulation using the Simscape Multibody Toolbox of Simulink to validate its performance. The dynamic model of the robot was exported from the CAD model. Four contact points were set to the four corners of the foot bottom plate to simulate the contact between the robot and the ground. The step length was set to 20 $\mathrm{c}$$\mathrm{m}$, and the time per step was 1 $\mathrm{s}$. The sample time of the simulation was $0.1$ $\mathrm{m}$$\mathrm{s}$. The cycle time of the controller was 1 $\mathrm{m}$$\mathrm{s}$. The process of walking is displayed in Figure 18. The robot took two steps from 1 $\mathrm{s}$ to 3 $\mathrm{s}$. The left foot was the support foot between 1 $\mathrm{s}$ and 2 $\mathrm{s}$, and the right foot was the support foot between 2 $\mathrm{s}$ and 3 $\mathrm{s}$. The red circle is the CoM of the biped. The green circle is the CoP. Figure 19 shows the trajectory of the CoM and CoP in the $xy$-plane. The red curve is the desired trajectory, and the blue curve is the measured trajectory. The biped tracked the trajectory accurately in simulation. We compared the position error of the extended generator with the traditional CP-based generator using the same controller to verify the performance of the extended CP-based walking pattern generator. Figure 20 shows the position error of the extended CP and traditional CP. As seen in Figure 20, the difference was not very significant. By using both walking pattern generators, the biped could walk and keep balance. In order to confirm the performance, we compared the absolute value of both errors by (59). If f is larger than zero, it means that the error of extended CP is smaller. The result was that $67.5\%$ of the time, the x position error of the extended CP was smaller, and $57.4\%$ of the time the y position error of the extended CP was also smaller.

$$f=\left|{\tilde{x}}_{trad}\right|-\left|{\tilde{x}}_{ext}\right|$$

(59)

7. Conclusions

In this paper, we proposed a new biped robot called FORREST. The robot consists of a hip with torque-controlled joints and mechanical elastic knee/ankle joints. The massive compliance in the robot will enable human–robot cooperation and bring humanoid robots into daily life. FORREST serves as an experimental platform for combining active elastic joints and passive elastic joints. The knee is actuated by an elastic ball screw system, and we use a novel elastic cable-driven parallel mechanism to drive the ankle joint. Unlike most cable-driven bipeds, our robot can maintain balance. Thus, in the future we will show how to use the robot in our daily surroundings. In this paper, we first described the structure and working principles of the knee and ankle, analyzed the kinematic and dynamic performance of the joints, and introduced the dynamic model of the entire robot. We proposed an extended CP-based walking pattern generator for our robot to achieve a smoother walking pattern including smooth trajectories. We validated the performance of the new generator using a simple control strategy in simulation, and the results show that the new walking pattern generator performs better than traditional generators. In the near future, we plan to realize the real walking of our biped and develop a new control strategy that takes into account its elasticity.