Robust Walking for Humanoid Robot Based on Divergent Component of Motion

In order to perform various complex tasks in place of humans, humanoid robots should walk robustly in the presence of interference. In the paper, an improved model predictive control (MPC) method based on the divergent components of motion (DCM) is proposed. Firstly, the humanoid robot model is simplified to a finite-sized foot-pendulum model. Then, the gait of the humanoid robot in the single-support phase (SSP) and double-support phase (DSP) is planned based on DCM. The center of mass (CoM) of the robot will converge to the DCM, which simplifies the feedback control process. Finally, an MPC controller incorporating an extended Kalman filter (EKF) is proposed to realize the tracking of the desired DCM trajectory. By adjusting the step duration, the controller can compensate for CoM trajectory errors caused by disturbances. Simulation results show that—compared with the traditional method—the method we propose achieves improvements in both disturbed walking and uneven-terrain walking.


Introduction
In recent years, research on humanoid robot has received high priority. As the first step for humanoid robots to replace humans, it is necessary to ensure robust walking for humanoid robots, especially when walking on uneven terrain [1]. The research on bipedal robust walking of humanoid robots has become one of the most challenging research topics [2].
When walking on uneven ground, humanoid robots need the ability to adjust their gait to avoid falling caused by disturbance. Due to the high dimension of contact force between the sole of the humanoid robot and the ground, the robot's leg-motion control is very difficult [3,4]. Nowadays, there are many theories for the balance control of humanoid robots, such as the zero-moment point (ZMP) [5] and DCM [6] criteria of stability, etc. In addition, an excellent controller is necessary to ensure robust walking of robots. Since MPC controllers can solve optimization problems with multiple constraints, they are widely used in the field of bipedal robot walking.
In this paper, an improved MPC controller is proposed to track the desired DCM trajectory, which enhances robustness of the robot suffering from disturbance. Desired DCM trajectories for the SSP and the DSP are generated. An inverted pendulum model with a finite-sized foot [7] is used in gait planning. It takes into account the ankle-joint moment to overcome the low robot dexterity caused by insufficient ankle drive. In addition, the MPC controller expanded by an EFK is designed, which can compensate for a certain range of disturbances by adjusting step duration during walking. The simulation results show that the hybrid controller we propose can change step length and width by adjusting step duration. large disturbances, experts in humanoid robotics have started using the DCM criterion for gait planning and controller design.
Englsberger et al. [31] proposed two control loop design methods, CPS and CPT, based on the 2D-DCM stability criterion. CPS was a real-time control method, which had high requirements on the hardware of the robot. CPT needed to generate the gait offline. A trajectory-tracking controller was designed in the paper and took the error between the desired and actual DCM as input. Then, Englsberger et al. [32] proposed a trajectorytracking controller based on ZMP, which compensated for changes in CoM height and angular momentum during the walking process. Sedye et al. [33] analyzed the influence of the angular momentum generated by the swinging legs during high-speed walking, which predicted the angular momentum of CoM in the planning stage. The generator produced smooth and continuous trajectories based on the instantaneous DCM.
On the one hand, the traditional controller based on ZMP has difficulty dealing with complex nonlinear constraints caused by the change in CoM height [34]. On the other hand, to improve stability, many constraints, such as structural size and stability margin, must be considered [35]. MPC can solve nonlinear multiple constraint optimization problems, so it is widely used in humanoids [36][37][38]. Krause et al. [39] were the first to design an MPC controller based on the DCM. They used double closed-loop control, which takes ZMP as input, the outer is an MPC controller, and the inner is a ZMP controller. Griffin et al. [40] analyzed the variation in DCM during walking and presented a time-varying DCM. To generate a robust gait, an MPC frame was formulated based on time-varying DCM that took step position and rotation into account. Simulation results verified the effectiveness of the method. Shafiee-Ashtiani et al. [41] designed a robust stepping MPC controller that could adjust step location and CoP. Experimental results showed that the robot could walk in various environments by taking time-varying DCM into account. To decrease the reaction time, Kasaei et al. [42] decoupled the control and replanning procedures. They formulated an offline optimal controller based on Linear-Quadratic-Gaussian (LQG) control [43] to achieve a reference trajectory-tracking process. Simulation in Matlab showed the capability of this algorithm. However, none of these methods consider the effect of noise when optimizing. We used EKF [44] to estimate the cost function to improve the robustness of bipedal walking.

Gait Planning for Humanoids
To adapt to real-world environments, bipedal robots should have human-like walking patterns [45]. The reproduction of a walking cycle for a walking robot is based on the generation of a step sequence [46]. Each step is composed of two successive phases: the single-support phase (SSP) and double-support phase (DSP) [47]. In the SSP, the foot of the swing leg does not directly affect any force, and the ZMP is in the support area formed by the support foot. In the DSP, both feet are subject to the reaction force of the ground at the same time. The support area at this time is a quadrilateral formed by two feet. In this phase, the ZMP moves from the support area formed by the left (right) foot to the support area formed by the right (left) foot [48]. The transfer of the ZMP between the feet allows the robot to keep moving forward. In addition, since the robot has a larger support area in the DSP, ZMP has a larger stability margin in the DSP [49]. Overall, in order to make the robot walk stably, the SSP and DSP are indispensable. Next, we introduce the inverted pendulum model used and the gait planning of the SSP and DSP.

LIPM with Finite-Sized Feet
The simplified model of bipedal walking is shown in Figure 1. It takes the ankle torque into account, which makes the CoP controllable [50]. As shown in Figure 1, m represents the total mass of the robot. We assume that the mass of the robot is concentrated in one point r = (x y z) T , which moves on the same level and z = z 0 . f = ( f x f y f z ) T is the reaction force of the floor acting on the ankle joint. r ankle = (x ankle y ankle 0) T is the projection of the ankle joint onto the ground. τ ankle = (τ ankle,x τ ankle,y τ ankle,z ) T is the ankle torque. Pr is the projection of the CoM position on the ground. and g = (0 0 −g) T is the gravitational acceleration vector.
Micromachines 2022, 13, x FOR PEER REVIEW 4 of 18 is the reaction force of the floor acting on the ankle joint.  The force and moment equations of LIPM with finite-sized feet are as follows [51]: mm =+ r f g (1) ankle ankle Due to the constraints of z = z0, z  and z  , both are zero, and it can be known from Formula (1) In summary, the dynamic equation of the bipedal walking model is:  The force and moment equations of LIPM with finite-sized feet are as follows [51]: Due to the constraints of z = z 0 , . z and ..
z, both are zero, and it can be known from Formula (1) that f z = mg. Taking it into Formula (2) to get f x , f y and τ ankle,z : In summary, the dynamic equation of the bipedal walking model is: where ω = g/z 0 . r cop = r ankle + ∆r cop Since the vertical motion of the CoM is fixed, Formula (4) can be written as Formula (5).
After homogeneous transformation, the unstable component of the state-space formula, which is also called DCM, can be obtained:

Single-Support Phase
The bipedal walking of the humanoid robot includes an SSP and a DSP. This paper assumes the transition of DCM and CoP between the feet is instantaneous when designing the SSP. It is shown in Figure 2.
After homogeneous transformation, the unstable component of the state-space formula, which is also called DCM, can be obtained:

Single-Support Phase
The bipedal walking of the humanoid robot includes an SSP and a DSP. This paper assumes the transition of DCM and CoP between the feet is instantaneous when designing the SSP. It is shown in Figure 2. Combining Formulas (4) and (6) leads to Formula (7): The solution to Formula (7) is:  Combining Formulas (4) and (6) leads to Formula (7): The solution to Formula (7) is: where ξ 0 is the initial position of the DCM, time t ∈ 0, t step i , and t step i is the duration of step i.
In the last step, the DCM overlap with the CoP and the bipedal motion stops. The desired DCM position at the end of each step can be obtained by Formula (8): Combining Formulas (8) and (9), the desired trajectory of DCM can be obtained: According to Formula (10), the expected DCM trajectory of each step can be obtained recursively.

Double-Support Phase
Considering only two SSPs in one walking cycle make the desired CoP trajectory discontinuous, the joint trajectory calculated by inverse kinematics is not continuous either, so, we plan the DSP trajectory. From Formula (7), we can get: From Formula (11), it can be seen that the continuity of the desired CoP trajectory requires continuity of the first derivative of the DCM trajectory, We smooth the desired DCM trajectory using a cubic interpolation function [52]: At the same time, Equation (12) must satisfy the boundary conditions of position and velocity: The polynomial coefficients can be solved by Equation (13): where T i is the duration of the i-th DSP. Thus, the DCM trajectory of the double-support stage is: The desired DCM trajectory is shown in Figure 3, the SSP trajectory is an exponential interpolation function, and the DSP is a cubic interpolation function.
After the above planning, the first derivative of the DCM trajectory is continuous everywhere and the CoP trajectory is continuous too. The desired DCM trajectory is shown in Figure 3, the SSP trajectory is an exponential interpolation function, and the DSP is a cubic interpolation function. After the above planning, the first derivative of the DCM trajectory is continuous everywhere and the CoP trajectory is continuous too.

Model Predictive Controller
Model predictive control is a control method that combines control and optimization [53,54]. It can solve optimization problems with multiple constraints. MPC relies on dynamic models of the process, most often linear empirical models obtained by system identification. The main advantage of MPC is the fact that it allows the current timeslot to be optimized, while keeping future timeslots in account [55].

Predictive Model
Given a continuous linear-time invariant system, such as Formula (7), it can be transformed as: However, the MPC model must be discrete. Thus, its discrete counterpart, evaluated using the zero-order hold (ZOH) technique, is: where T is the sampling time and Fd and Gd are: The prediction model of model predictive control can be obtained:

Model Predictive Controller
Model predictive control is a control method that combines control and optimization [53,54]. It can solve optimization problems with multiple constraints. MPC relies on dynamic models of the process, most often linear empirical models obtained by system identification. The main advantage of MPC is the fact that it allows the current timeslot to be optimized, while keeping future timeslots in account [55].

Predictive Model
Given a continuous linear-time invariant system, such as Formula (7), it can be transformed as: . ξ = ωξ − ωr cop (15) However, the MPC model must be discrete. Thus, its discrete counterpart, evaluated using the zero-order hold (ZOH) technique, is: where T is the sampling time and F d and G d are: The prediction model of model predictive control can be obtained: So, the prediction model based on state-space variables is: where P x is the matrix of coefficients formed by F. H x is the matrix of coefficients formed by F and G, and y k+1 , ξ k , r cop k are deviations of output, variable, and control signal respectively.

Constraints of Model Predictive Control
To ensure the viability of the gait, it is desired that the CoP cannot leave the support area created by the foot. The blue polygon in Figure 4 is the support polygon area of the DSP, and the red polygon is the support area of the SSP.
So, the prediction model based on state-space variables is: where Px is the matrix of coefficients formed by F. Hx is the matrix of coefficients formed by F and G, and r are deviations of output, variable, and control signal respectively.

Constraints of Model Predictive Control
To ensure the viability of the gait, it is desired that the CoP cannot leave the support area created by the foot. The blue polygon in Figure 4 is the support polygon area of the DSP, and the red polygon is the support area of the SSP. The equation of the straight line ab is: The next step (taking the left foot as an example) expects the CoP needs to be satisfied: which is: Then, the constraint to be satisfied by CoP at time k is: The equation of the straight line ab is: The next step (taking the left foot as an example) expects the CoP needs to be satisfied: Then, the constraint to be satisfied by CoP at time k is:

Cost Function
The cost function is used to ensure that the controller tracks the desired trajectory. This paper adopts the following cost function: where, e k+1 = ξ k+1 − ξ re f k+1 , Q, and R are positive symmetric matrices to correct variables that deviate from the reference trajectory. The cost function can be divided into two parts: the part related to the state space and the part related to the control input. The former is used to reduce the error between the expected and predicted DCM positions, and the latter is used to minimize the CoP rate of change, resulting in a smooth ZMP trajectory.
By merging Formulas (19) and (23) we get the optimization function: In summary, we can find the control input:

Extended Kalman Filter
To effectively deal with the nonlinearities caused by changes in the height of CoM during walking, we use the EKF, which reduces the impact of measurement noise and external disturbances on the system [56,57]. The equations used in the EKF are the same as the ones in [58], where the stages of update and measure are expressed as follows: Measure: whereξ − k andξ − k−1 are the estimated value at now and the last moment. w k is the system noise. y k and r cop k are the output and input respectively, q and s are the weights, v k is the measurement noise, and these values determine the output. P − k is the covariance matrix which depends on the system matrix A k , the covariance matrix of the previous moment P k−1 , the noise covariance matrix Q k−1 , and the system noise matrix W k . K k represents the Kalman gain, which depends on H k , V k , and P − k , where H k and V k are matrices of measurement and measuring.
By predicting one step ahead, the EKF updates or revises the forecast using system measurements to estimate the value in the process [59]. A detailed explanation of the above equations can be found in the literature [60]. The matrix H determines the magnitude of the Kalman filter gain, which is shown in Formula (28). In summary, the LIPM with finite-sized foot and the optimal algorithm determine the procedure model. In order to obtain the optimal solution and iteratively update the weight matrices, Q, R, and H must be linearized: It is necessary to obtain general expressions for MPC and EKF. Combining Formulas (18) and (25), we can get linear equivalence: ξ 1 , ξ 2 are the penalizations q and s, respectively, which are also used in MPC. The state-space variables of the robot model are the DCM in the x and y directions, respectively. We assign the same weight to both directions. In other words, we consider the DCM in the x and y directions to be equally important. The overall control block diagram is shown in Figure 5.
1  , 2  are the penalizations q and s, respectively, which are also used in MPC.
The state-space variables of the robot model are the DCM in the x and y directions, respectively. We assign the same weight to both directions. In other words, we consider the DCM in the x and y directions to be equally important. The overall control block diagram is shown in Figure 5.

Experiments
We used matlab2020b as a simulation experiment platform to build a bipedal model. Figure 6 shows the link structure (left) and simulation model (right) of a humanoid robot. Where, Li and mi are the length and mass of the links, respectively,

Experiments
We used matlab2020b as a simulation experiment platform to build a bipedal model.  In order to cope with complex environments, humanoid robots should have the ability to walk on uneven ground and cope with external disturbances. Three simulation experiments were designed to demonstrate the effectiveness of the DCM planning and control method.
Simulation 1 The robot walks at a constant speed on a flat surface, and there is no disturbance during the walking.
Simulation 2 The robot walks at a constant speed on flat ground. During the walking, an impact disturbance is applied. The disturbance is applied to the CoM, and it changes the speed and direction of the CoM. In this simulation, the robot stands with a single support on its right leg at 9 s. To obtain the maximum performance of the controller, a disturbance will be applied to the right hip joint at 9 s and last 10 ms, because applying disturbance to the support leg will impact the robot greatly.
Simulation 3 When a robot is walking outdoors, the height of the center of mass of the robot does not remain on a horizontal plane. However, in this case, the calculated DCM trajectory will not guarantee the robustness of the bipedal walking. To ensure that the controller can compensate for the effects of CoM height changes, we built an uneven In order to cope with complex environments, humanoid robots should have the ability to walk on uneven ground and cope with external disturbances. Three simulation experiments were designed to demonstrate the effectiveness of the DCM planning and control method.
Simulation 1 The robot walks at a constant speed on a flat surface, and there is no disturbance during the walking.
Simulation 2 The robot walks at a constant speed on flat ground. During the walking, an impact disturbance is applied. The disturbance is applied to the CoM, and it changes the speed and direction of the CoM. In this simulation, the robot stands with a single support on its right leg at 9 s. To obtain the maximum performance of the controller, a disturbance will be applied to the right hip joint at 9 s and last 10 ms, because applying disturbance to the support leg will impact the robot greatly. Simulation 3 When a robot is walking outdoors, the height of the center of mass of the robot does not remain on a horizontal plane. However, in this case, the calculated DCM trajectory will not guarantee the robustness of the bipedal walking. To ensure that the controller can compensate for the effects of CoM height changes, we built an uneven ground in Simulink and the ground height change was ∆h ∈ [−0.05m, +0.05m]. In this simulation, the robot walks on uneven terrain at a constant speed.
The results of experiments are compared in Figures 7-9. Figure 7 indicates the real CoM position in the coronal plane in three simulations. The CoM position in the sagittal plane is shown in Figure 8. Figure 9 illustrates the height variation of the robot's CoM in three simulations. The red curve represents the result of simulation 1, the blue curve represents the result of simulation 2, and the green curve describes the result of simulation 3.  Figure 7 shows that the coronal motion of the CoM can be considered as a periodic motion similar to a sinusoid wave when the robot is walking at a constant speed without disturbance. During 0 s to 5 s, the robot is lowering the CoM by bending legs and shifting the CoM to the left leg, so three experiments have the same results. The first two cycles of the robot are walking on flat ground. Simulation 2 (blue curve) shows that when being disturbed by a large shock at 9 s, after 3.5 cycles of correction by the controller, the motion of the coronal plane will balance in a new position and continue to maintain a periodic motion similar to a sinusoid. Thus, the robot takes 7 steps (3.5 cycles) during 9 s to 13 s, one more step than walking on flat ground. Taking one more step in the same amount of time shortens the step length and prevents it from exceeding the maximum limit. Besides this, the robot has one more double-support phase that helps the ZMP stabilize in the center of the support area. Simulation 3 (green curve) shows that when walking on uneven terrain, although the coronal motion of the robot is different from that of undisturbed walking in each cycle, its equilibrium position is always near the equilibrium position of the coronal motion during undisturbed walking. Similar to our description in simulation 2, the controller will adjust the gait of the robot by shortening or lengthening the step time Δt to compensate for the velocity and position errors of CoM. Over all, from Figure 7, we draw the following conclusion: the controller compensates the trajectory error of the CoM by adjusting the step width and increasing the number of DSP.   Figure 7 shows that the coronal motion of the CoM can be considered as a periodic motion similar to a sinusoid wave when the robot is walking at a constant speed without disturbance. During 0 s to 5 s, the robot is lowering the CoM by bending legs and shifting the CoM to the left leg, so three experiments have the same results. The first two cycles of the robot are walking on flat ground. Simulation 2 (blue curve) shows that when being disturbed by a large shock at 9 s, after 3.5 cycles of correction by the controller, the motion of the coronal plane will balance in a new position and continue to maintain a periodic motion similar to a sinusoid. Thus, the robot takes 7 steps (3.5 cycles) during 9 s to 13 s, one more step than walking on flat ground. Taking one more step in the same amount of time shortens the step length and prevents it from exceeding the maximum limit. Besides this, the robot has one more double-support phase that helps the ZMP stabilize in the center of the support area. Simulation 3 (green curve) shows that when walking on uneven terrain, although the coronal motion of the robot is different from that of undisturbed walking in each cycle, its equilibrium position is always near the equilibrium position of the coronal motion during undisturbed walking. Similar to our description in simulation 2, the controller will adjust the gait of the robot by shortening or lengthening the step time Δt to compensate for the velocity and position errors of CoM. Over all, from Figure 7, we draw the following conclusion: the controller compensates the trajectory error of the CoM by adjusting the step width and increasing the number of DSP. From Figure 8, it can be seen that when there is no external interference, the sagittal displacement of the CoM is proportional to the time after the robot starts to walk. The red  Figure 7 shows that the coronal motion of the CoM can be considered as a periodic motion similar to a sinusoid wave when the robot is walking at a constant speed without disturbance. During 0 s to 5 s, the robot is lowering the CoM by bending legs and shifting the CoM to the left leg, so three experiments have the same results. The first two cycles of the robot are walking on flat ground. Simulation 2 (blue curve) shows that when being disturbed by a large shock at 9 s, after 3.5 cycles of correction by the controller, the motion of the coronal plane will balance in a new position and continue to maintain a periodic motion similar to a sinusoid. Thus, the robot takes 7 steps (3.5 cycles) during 9 s to 13 s, one more step than walking on flat ground. Taking one more step in the same amount of time shortens the step length and prevents it from exceeding the maximum limit. Besides this, the robot has one more double-support phase that helps the ZMP stabilize in the center of the support area. Simulation 3 (green curve) shows that when walking on uneven terrain, although the coronal motion of the robot is different from that of undisturbed walking in each cycle, its equilibrium position is always near the equilibrium position of the coronal motion during undisturbed walking. Similar to our description in simulation 2, the controller will adjust the gait of the robot by shortening or lengthening the step time ∆t to compensate for the velocity and position errors of CoM. Over all, from Figure 7, we draw the following conclusion: the controller compensates the trajectory error of the CoM by adjusting the step width and increasing the number of DSP. fluctuation disappears after the equilibrium is restored. Compared to the desired DCM, the real DCM will produce a larger displacement for the same time interval Δt due to the increased CoM velocity. This will make the step distance of the robot exceed the maximum limit, and even if there is still a margin for the step, a longer step will have a great impact on the next cycle of walking. Thus, the controller shortens the Δt to make the robot's feet land more frequently. In this way, the robot can adjust the ZMP more quickly, which stabilizes the robot. When walking on uneven terrain (green curve), the sagittal motion of the CoM fluctuates continuously compared with normal walking, but the robot can still walk forward. In this simulation, the velocity of the CoM will increase or decrease. When the velocity of the CoM is increased, the controller will shorten Δt, as we described above. However, the velocity of the CoM may be decreased sometimes when walking on uneven terrain. To solve this problem, the controller lengthens Δt. As shown in Figure 8, during 9.5 s to 10.1s and 12 s to 13 s in simulation 3 (green curve), we can see a significant increase in the step distance of the robot. Overall, from Figure 8, we can see that the controller compensates for the trajectory error of the CoM by adjusting the step length. Figure 9. Variation in the height of the CoM. Figure 9 indicates that when the robot walks without disturbance, its CoM is always at the same horizontal height. When it is disturbed by impact, its CoM height changes, and gradually returns to the original horizontal height under the action of the controller. The results of the robot walking on uneven terrain show that even if the height of the CoM of the robot changes continuously, it will never fall out of balance. Figure 10 is the comparison of ZMP and DCM in simulation 2 and simulation 3, respectively. Compared to the DCM method, a robot controlled using the ZMP method falls to the left when it suffers an impact disturbance and cannot walk on uneven terrain. The ZMP-based MPC, which uses CoM as the controller input, cannot compensate for displacement and the change of acceleration caused by the external force to the CoM. Walking on uneven terrain causes the height of the robot's CoM to change constantly; similarly, the acceleration in the z-direction of the CoM can also be difficult to compensate for. From Figure 8, it can be seen that when there is no external interference, the sagittal displacement of the CoM is proportional to the time after the robot starts to walk. The red curve shows that when the robot walks on flat ground, it is a uniform linear motion. When the robot is disturbed, there is a small fluctuation in the CoM trajectory, and the fluctuation disappears after the equilibrium is restored. Compared to the desired DCM, the real DCM will produce a larger displacement for the same time interval ∆t due to the increased CoM velocity. This will make the step distance of the robot exceed the maximum limit, and even if there is still a margin for the step, a longer step will have a great impact on the next cycle of walking. Thus, the controller shortens the ∆t to make the robot's feet land more frequently. In this way, the robot can adjust the ZMP more quickly, which stabilizes the robot. When walking on uneven terrain (green curve), the sagittal motion of the CoM fluctuates continuously compared with normal walking, but the robot can still walk forward. In this simulation, the velocity of the CoM will increase or decrease. When the velocity of the CoM is increased, the controller will shorten ∆t, as we described above. However, the velocity of the CoM may be decreased sometimes when walking on uneven terrain. To solve this problem, the controller lengthens ∆t. As shown in Figure 8, during 9.5 s to 10.1s and 12 s to 13 s in simulation 3 (green curve), we can see a significant increase in the step distance of the robot. Overall, from Figure 8, we can see that the controller compensates for the trajectory error of the CoM by adjusting the step length. Figure 9 indicates that when the robot walks without disturbance, its CoM is always at the same horizontal height. When it is disturbed by impact, its CoM height changes, and gradually returns to the original horizontal height under the action of the controller. The results of the robot walking on uneven terrain show that even if the height of the CoM of the robot changes continuously, it will never fall out of balance. Figure 10 is the comparison of ZMP and DCM in simulation 2 and simulation 3, respectively. Compared to the DCM method, a robot controlled using the ZMP method falls to the left when it suffers an impact disturbance and cannot walk on uneven terrain. The ZMP-based MPC, which uses CoM as the controller input, cannot compensate for displace-ment and the change of acceleration caused by the external force to the CoM. Walking on uneven terrain causes the height of the robot's CoM to change constantly; similarly, the acceleration in the z-direction of the CoM can also be difficult to compensate for.  Figure 11 shows the real position of ZMP for the MPC and the MPC+EFK method. Figure 11a,b show the real ZMP position in the coronal plane and sagittal plane, respectively. In the coronal plane, the control method we proposed can ensure the ZMP is always positioned on the central axis of the support area and the ZMP trajectory is smooth when it switches between feet. However, the position of the ZMP trajectory is offset from the central axis of the support area and closer to the inside of the foot. Besides, the ZMP trajectory will suddenly change when it switches between two feet. In the sagittal plane, both methods have almost the same result. Overall, the method we proposed generates smoother ZMP trajectories with larger stability margins.

Conclusions
We propose a robust biped robot gait-generation and control method. We use LIPM with finite-sized feet as the mathematical model of the robot to plan the expected DCM  Figure 11 shows the real position of ZMP for the MPC and the MPC+EFK method. Figure 11a,b show the real ZMP position in the coronal plane and sagittal plane, respectively. In the coronal plane, the control method we proposed can ensure the ZMP is always positioned on the central axis of the support area and the ZMP trajectory is smooth when it switches between feet. However, the position of the ZMP trajectory is offset from the central axis of the support area and closer to the inside of the foot. Besides, the ZMP trajectory will suddenly change when it switches between two feet. In the sagittal plane, both methods have almost the same result. Overall, the method we proposed generates smoother ZMP trajectories with larger stability margins.  Figure 11 shows the real position of ZMP for the MPC and the MPC+EFK method. Figure 11a,b show the real ZMP position in the coronal plane and sagittal plane, respectively. In the coronal plane, the control method we proposed can ensure the ZMP is always positioned on the central axis of the support area and the ZMP trajectory is smooth when it switches between feet. However, the position of the ZMP trajectory is offset from the central axis of the support area and closer to the inside of the foot. Besides, the ZMP trajectory will suddenly change when it switches between two feet. In the sagittal plane, both methods have almost the same result. Overall, the method we proposed generates smoother ZMP trajectories with larger stability margins.

Conclusions
We propose a robust biped robot gait-generation and control method. We use LIPM with finite-sized feet as the mathematical model of the robot to plan the expected DCM trajectory of the robot in the SSP and DSP. Since the CoM will automatically converge to

Conclusions
We propose a robust biped robot gait-generation and control method. We use LIPM with finite-sized feet as the mathematical model of the robot to plan the expected DCM trajectory of the robot in the SSP and DSP. Since the CoM will automatically converge to the DCM, this paper proposes a DCM-based MPC controller to achieve indirect control of the CoM and ZMP by tracking the desired DCM trajectory. The combination of EFK and MPC reduces interference from the outside world. After being disturbed, the robot no longer moves in the desired trajectory. The controller will change the landing point of the actual DCM by optimizing the step time. The change in DCM implies an adjustment in step size and width so that the robot can compensate for CoM velocity and position errors caused by disturbances. Particularly, when the velocity of the CoM is increasing, the controller will shorten the step time to increase the number of DSPs and decrease the step size. It can adjust the ZMP and compensate for the velocity error of the CoM in real time. The experimental results show that the proposed method can effectively deal with the impact disturbance during the walking process of the robot, and at the same time the robot can walk on uneven ground.
In the future, we plan to give different weights to state variables in EKF. This idea is to realize the anthropomorphism of robots. For example, the DCM in the y-direction has a more significant impact on the stability of the robot when a humanoid robot is crossing a bridge or walking on a narrow road. We will let the robot walk on a suspended stone road. It should focus on the lateral external forces when suffering from disturbance. Even if the disturbance cannot be compensated for, it should fall in the forward direction. This can avoid irreparable damage caused by the robot falling from a swell.
Author Contributions: Z.Z. and L.Z. conceived and designed the experiments; S.X., N.X. and X.W. made suggestions for revision of the paper; and Z.Z. wrote the paper. All authors have read and agreed to the published version of the manuscript.