Stability Criterion for Dynamic Gaits of Quadruped Robot

Dynamic-stability criteria are crucial for robot’s motion planning and balance recovery. Nevertheless, few studies focus on the motion stability of quadruped robots with dynamic gait, none of which have accurately evaluated the robots’ stability. To fill the gaps in this field, this paper presents a new stability criterion for the motion of quadruped robots with dynamic gaits running over irregular terrain. The traditional zero-moment point (ZMP) is improved to analyze the motion on irregular terrain precisely for dynamic gaits. A dynamic-stability criterion and measurement are proposed to determine the stability state of the robot and to evaluate its stability. The simulation results show the limitations of the existing stability criteria for dynamic gaits and indicate that the criterion proposed in this paper can accurately and efficiently evaluate the stability of a quadruped robot using such gaits.


Introduction
Legged robots that imitate animals have flexible joints and interact with the ground intermittently through their feet, giving them excellent environmental adaptability, which makes them suitable for working in complex environments such as mountains, disaster sites, and warehouses.Compared with other legged robots, quadruped robots have better maneuverability and load capacity.Over the past several decades, there has been an increasing interest in quadruped robots [1][2][3][4].
When considering the dynamic motion of legged robots, one must pay attention to their stability.It is very important for robots to complete the assigned task while maintaining balance and moving on an even keel.Researchers have proposed a variety of static-and dynamic-stability-analysis methods to evaluate the stability of legged robot.The static-stability criteria and margins are only applicable for analyzing low-speed static gaits rather than the dynamic motion because they do not consider inertial forces or external impacts [5].Fukuoka and Kimura proposed a stability criterion named wide stability margin (WSM) [6].They hypothesized that a robot can maintain its balance if the projection of the center of mass (CoM) of the body is within the support polygon formed by the projection of the current support and swing feet on the horizontal plane; the shortest distance from this projection point to the boundaries of the support polygon is used to evaluate the robot's stability.This stability criterion, in essence, is a static-stability criterion, and it is not suitable for analyzing the stability of high-speed motion.The rationality of the definition of the support polygon has not been explained; there is no experiment to verify that this stability criterion and margin can accurately evaluate the robot's stability.Many dynamic-stability criteria have been presented to address the limitations of static-stability criteria.
The main idea behind the dynamic-stability criterion is that, at a certain state, the motion of the robot is considered to be stable if the torque caused by the ground-reaction force can prevent the robot from tumbling around any support boundary.
The most common dynamic-stability criteria are based on the zero-moment point (ZMP) or the center of pressure (CoP).ZMP refers to a point on the ground where the net torque in the direction parallel to the ground caused by the gravity and inertial forces of each part of the robot is zero.ZMP shows the tendency of the robot to tumble.Researchers believe that the robot is stable if ZMP is located inside the support polygon, which is connected to the support feet.The minimum distance between ZMP and the support boundaries is used as the stability margin.The larger the minimum distance is, the better the stability that the robot can achieve [7][8][9][10].CoP refers to the point where the net torque caused by the ground-reaction forces is zero; it coincides exactly with ZMP when the robot is dynamically stable [11,12].The traditional ZMP-based-stability criteria can only be applied to planar motion, and some researchers have made improvement to analyze the motion on irregular terrain by adjusting the support plane [11,13,14].The stability criteria based on ZMP do not consider the influence of the robot's current speed on its stability, but speed is very important for dynamic stability [15].
Methods using the leg-end-supporting moment (LSM) and tumble stability margin (TSM) [13,[16][17][18][19], among others, directly consider the resultant moment around each support boundary.The robot will be stable if the resultant moments can counteract the torque which makes the robot tumble.A simpler stability criterion is the force-angle-stability measure (FASM) [20], which considers the angle between the resultant force of gravity, inertial force, and other external forces at the CoM and the perpendicular line from CoM to the support boundary.Other researchers have proposed some energy-based stability criteria [21,22], e.g., the normalized dynamic energy stability margin (NDESM).These methods evaluate the stability of the robot by quantifying the maximum impact energy that the robot can absorb without losing its stability.
The limitation of the above stability criteria is that they are all suitable to the situation in which a support polygon exists when the legged robots are moving.They can be used to analyze the stability of walking gait (when the supporting triangle exists), but are not suitable for dynamic gaits, which only have a support line in most cases, such as trotting, pacing, and bounding.
There have been few studies concerning the stability criteria for quadruped robots with dynamic gaits, all of which focus on trotting.Most researchers believe that a robot with a trot gait is stable only when the current ZMP is located on the support line, because there is no tumbling moment around the support line.The distance between ZMP and the support line (substantially equivalent to the tumbling moment around the support line) is used to guide the motion control of the robot [14,23,24].Such analytical methods do not consider the effect of the dynamic process of support-leg alternation on the robot's motion stability.In practice, even if ZMP is not located on the current support line, the robot will not necessarily lose its balance and fall.Hence, it is difficult to correctly evaluate the stability of the robot by directly using the distance between ZMP and the current support line.Therefore, the dynamic-stability criterion requires further consideration.
Lang et al. [25] proposed a stability analysis method based on Lyapunov theory of switching system.They built the dynamic model of the trotting robot as a nonlinear switching system, and designed a Lyapunov function to prove that the switching system is uniform asymptotic stable.This method is used to choose the landing position of the swing feet to make the position tracking errors of the robot's body converge to zero.The computational complexity of this method is too high, and it cannot evaluate the stability of the robot's motion.Another stability measure is called Landing Accordance Ratio (LAR) [26].This method considers the time difference between the moments of the feet in the same phase touch down.If the discordance amount is zero, there is no unexpected moment, and the robot moves stably.There are three major shortcomings of this method: (1) it cannot evaluate the stability of the robot in time, because LAR is updated when the two feet in the same phase both touch the ground; (2) it is not suitable for the situation where the ground-reaction forces at current support feet cannot prevent the robot from tumbling; and (3) it is hard to decide when to update the value of LAR if the robot tips over, because one of the feet will never touch the ground.
The good stability criterion and evaluation method should have the following characteristics: (1) it should determine the stability state of the robot efficiently and correctly; (2) it should evaluate the robot's stability accurately; and (3) it should provide guidance for the robot's motion planning and balance recovery.The purpose of this paper is to present a new stability criterion for the quadruped robots with dynamic gaits over irregular terrain.
When animals experience impacts in nature, they quickly adjust the position of their next set of support feet to maintain stability.This indicates that dynamic stability depends not only on their current motion and support state but also on the positions of the next set of support points.Inspired by quadruped animal moving gesture and referring to the previous dynamic-stability criteria, this paper assumes that a robot is stable if the current or next set of support feet can provide the moment necessary to prevent the robot from tipping in any direction.The less likely the robot is to tip over, the smaller the difference among ground-reaction forces at each support foot are, meaning that the robot has higher capability to stay stable.
Based on the above ideas, by defining and adjusting the virtual-support plane, the traditional ZMP can be extended to three-dimensional (3D) space so that it can be used to analyze the motion state on irregular terrain.An extra term containing the velocity of the robot's body is added to the equation for the calculation of ZMP to better reflect the robot's motion state.Here, we define the virtual-support quadrilateral and propose a stability criterion for quadruped robots with dynamic gaits.The robot can stay stable if the expanded ZMP is located within this quadrilateral, because the current or next set of support feet can provide the moment necessary to prevent the robot from tipping over.The stability is quantified using the distance between ZMP and the boundary of the virtual-support quadrilateral, the angle between the vector from CoM to ZMP and the normal vector of the support plane, and the distance between ZMP and the support line.The simulation results show that the trotting robot can still run stably even if ZMP is not on the current support line, and the distance between ZMP and the current support line cannot accurately assess the robot's stability.The results also indicate that the proposed stability criterion and measure can efficiently and accurately evaluate the stability of a quadruped robot using dynamic gaits.
This paper is outlined as follows.The method for expanding ZMP to 3D space is described in Section 2. The stability criterion for the motion of a quadruped robot with a dynamic gait over irregular terrain is introduced in Section 3. The method of calculating the dynamic-stability measurement is presented in Section 4. The simulation results and discussion are reported in Section 5. Finally, conclusions are presented in Section 6.

ZMP for Irregular Terrain
In this section, the traditional ZMP is modified to deal with the motion on irregular terrain.Because most of the mass is concentrated in the body of the robot, the mass of the legs is ignored to simplify our stability analysis.Traditional ZMP refers to a point on the ground at which the net torque caused by gravity, the inertial force, and the inertial moment of each part of the robot is zero; it is only applicable when analyzing the motion on the horizontal plane.In this paper, ZMP is extended to 3D space by defining and adjusting the plane of the support feet, and the velocity term is introduced to reflect the motion more efficiently and accurately.
On irregular terrain, ZMP should be considered as a point in 3D space.
As shown in Figure 1, the position of the body's CoM (C) in the inertial coordinate system (i) . The gravity of the body is , and the inertial moment is . The position of ZMP in the inertial coordinate system is assumed to be P i ZMP = x i ZMP y i ZMP z i ZMP T ; hence, according to the definition of ZMP, Only two of the three equations are linearly independent [14], but there are three unknowns.An additional condition must be given to solve this problem.
Referring to the assumptions made in traditional ZMP, it is stipulated in this paper that ZMP is always located on the current support plane (which can be nonhorizontal).For the convenience of calculation, the support coordinate system () (Figure 2) is defined.The origin of this coordinate system is located at the midpoint  of the current support line  ̅̅̅̅ , the direction of   is from the rear foot to the front foot of the current support line, and   is perpendicular to the current support plane and points downward.Another coordinate system used in this paper is the principal axes coordinate system of the body (), whose origin is located at the body's CoM.The following can be obtained from Equation (1): Only two of the three equations are linearly independent [14], but there are three unknowns.An additional condition must be given to solve this problem.
Referring to the assumptions made in traditional ZMP, it is stipulated in this paper that ZMP is always located on the current support plane (which can be nonhorizontal).For the convenience of calculation, the support coordinate system (s) (Figure 2) is defined.The origin of this coordinate system is located at the midpoint S of the current support line AB, the direction of z s is from the rear foot to the front foot of the current support line, and x s is perpendicular to the current support plane and points downward.Another coordinate system used in this paper is the principal axes coordinate system of the body (c), whose origin is located at the body's CoM.
Only two of the three equations are linearly independent [14], but there are three unknowns.An additional condition must be given to solve this problem.
Referring to the assumptions made in traditional ZMP, it is stipulated in this paper that ZMP is always located on the current support plane (which can be nonhorizontal).For the convenience of calculation, the support coordinate system () (Figure 2) is defined.The origin of this coordinate system is located at the midpoint  of the current support line  ̅̅̅̅ , the direction of   is from the rear foot to the front foot of the current support line, and   is perpendicular to the current support plane and points downward.Another coordinate system used in this paper is the principal axes coordinate system of the body (), whose origin is located at the body's CoM.The position of ZMP in the support coordinate system P s ZMP = 0 y s ZMP z s ZMP T is calculated by considering the force and moment balance in this system: x s ZMP = 0; (5) It is difficult to ascertain the actual support plane for some dynamic gaits (trotting, pacing, bounding, etc.) because only two support feet are used most of the time.To solve this problem, the positions of the previous support line, A F B F , and the next support line, A L B L , are comprehensively considered (for gaits with only one support point, such as gallops, more previous and later support points should be used).The definition of the virtual-support plane is proposed.
The following assumptions are made for the current virtual-support plane: (1) this plane always contains the current support line AB; and (2) from the disappearance of the previous support line to the disappearance of the current support line, the virtual-support plane gradually changes from Σ 1 (determined by the current and previous support lines) to Σ 2 (determined by the current and next support lines).In practice, for uneven terrain, AB and A F B F (or A L B L ) may not lie on the same plane.This paper presents a method to determine the virtual-support plane.
The unit normal of the plane on which A, B, and A F lie is n 11 (the angle between this normal and x i is acute), and that of the plane on which A, B, and B F lie is n 12 (it also has an acute angle with x i ).The unit normal vector n 1 of the plane Σ 1 lies on the bisector of the angle formed by n 11 and n 12 (Figure 3): This paper presents a method to determine the virtual-support plane.
The unit normal of the plane on which , , and   lie is  11 (the angle between this normal and   is acute), and that of the plane on which , , and   lie is  12 (it also has an acute angle with   ).The unit normal vector  1 of the plane  1 lies on the bisector of the angle formed by  11 and  12 (Figure 3):  Set AA F is equal to the zero vector, and n 11 = −n 11 , if the x i component of n 11 is negative.A similar method is used to deal with Equation (9).Hence, n 21 is the unit normal of the plane on which A, B, and A L lie and n 22 is the unit normal of the plane on which A, B, and B L lie. n 2 is the unit normal vector of Σ 2 : Then, from the disappearance of the previous support line to that of the current support line, the unit normal vector of the virtual-support plane (i.e., x s in the inertial coordinate system) is where T b is the duration between two adjacent steps (half the duration of a gait cycle for trot gait), and t is the time gap between the current running time and the time point when the previous support line disappeared.The z s axis of the support coordinate system is and the y s axis is The rotation matrix from the support to the inertial coordinate system is The current velocity of the body, which can reflect the robot's current motion state, is not taken into account in the traditional ZMP-based stability criteria, however, the velocity is very important for maintaining balance [15], especially under some special condition.For example, when the current ground-reaction force cannot just offset the torque, the robot is made to tumble around a certain direction; however, if the robot's velocity is high enough to allow it to quickly enter a stable area, the robot can maintain its balance.A more accurate assessment of the robot's stability can be proposed by taking the velocity into account.In this paper, the velocity terms are added to the calculation of ZMP: The modified ZMP is denoted by ZMP 0 .ZMP 0 still lies on the current virtual-support plane, and the differences between the current actual velocity, v s r , and the expected velocity, v s d , in the y s and z s directions are considered in the calculation of ZMP.η is a constant-value coefficient considering the magnitude of the size and velocity of the robot: L is the length of the robot's body and W is the width.The expected benefits of doing so are the following: (1) the current stability of the robot can be assessed more accurately and efficiently; and (2) a reference can be provided to eliminate undesired velocity during motion planning.

Dynamic-Stability Criterion
The first problem to be solved when considering the stability of the robot is determining whether the robot can maintain balance at the current state.The stability of a robot with a dynamic gait is closely related to two factors: (1) the ground-reaction forces that act on the current support feet; and (2) the ground-reaction forces that will acti on the next support feet.The robot is unstable if neither of these forces can prevent it from tumbling.In addition, the elimination of unexpected speed is also taken into account when evaluating the robot's stability by adding the velocity term into the calculation of ZMP.
There are always two legs moving in the same phase for most dynamic gaits; hence, these two legs can be equivalent to one virtual leg from the midpoint of the support line to the robot's CoM.The stability criterion is considered in the current virtual-support plane.A L and B L are the projection points on the support plane of the next set of support feet, A L and B L , respectively, as shown in Figure 4. taken into account when evaluating the robot's stability by adding the velocity term into the calculation of ZMP.
There are always two legs moving in the same phase for most dynamic gaits; hence, these two legs can be equivalent to one virtual leg from the midpoint of the support line to the robot's CoM.The stability criterion is considered in the current virtual-support plane.  ′ and   ′ are the projection points on the support plane of the next set of support feet,   and   , respectively, as shown in Figure 4.It should be noted that the direction of the support force,   , is perpendicular to the support plane and points outward, and its main role is to support the robot standing and to maintain the stability of the robot's rotation around the   and   axes.The frictional force   is located on the support plane, which is opposite to the direction of the feet's motion trend with respect to this plane; its function is to push the robot forward and to prevent the robot from rotating around the   axis.
The contact forces act on the two support feet,  and .Hence, the current support line cannot provide the torque around the   axis; its main role is to prevent the robot from tumbling around the   axis.The torque that prevents the robot from tipping around the   axis can only be provided by the next set of support feet.The contact forces at   ′ and   ′ also play a supplementary role in maintaining rotational stability around   .The rotation of the robot around   is stopped by the friction at the current support and next support lines, but this is not the main reason the robot loses stability.ZMP0 can be considered to be on the extension line of the equivalent force,   , acting at the CoM of the body (Figure 5).In the     plane,  1 represents the point with the largest   -coordinate among the four support points (, ,   ′ , and   ′ ), whereas the   -coordinate of  4 is the smallest one, as shown in Figure 6.It should be noted that the direction of the support force, F n , is perpendicular to the support plane and points outward, and its main role is to support the robot standing and to maintain the stability of the robot's rotation around the y s and z s axes.The frictional force F f is located on the support plane, which is opposite to the direction of the feet's motion trend with respect to this plane; its function is to push the robot forward and to prevent the robot from rotating around the x s axis.
The contact forces act on the two support feet, A and B. Hence, the current support line cannot provide the torque around the z s axis; its main role is to prevent the robot from tumbling around the y s axis.The torque that prevents the robot from tipping around the z s axis can only be provided by the next set of support feet.The contact forces at A L and B L also play a supplementary role in maintaining rotational stability around y s .The rotation of the robot around x s is stopped by the friction at the current support and next support lines, but this is not the main reason the robot loses stability.
ZMP 0 can be considered to be on the extension line of the equivalent force, F e , acting at the CoM of the body (Figure 5).There are always two legs moving in the same phase for most dynamic gaits; hence, these two legs can be equivalent to one virtual leg from the midpoint of the support line to the robot's CoM.The stability criterion is considered in the current virtual-support plane.  ′ and   ′ are the projection points on the support plane of the next set of support feet,   and   , respectively, as shown in Figure 4.It should be noted that the direction of the support force,   , is perpendicular to the support plane and points outward, and its main role is to support the robot standing and to maintain the stability of the robot's rotation around the   and   axes.The frictional force   is located on the support plane, which is opposite to the direction of the feet's motion trend with respect to this plane; its function is to push the robot forward and to prevent the robot from rotating around the   axis.
The contact forces act on the two support feet,  and .Hence, the current support line cannot provide the torque around the   axis; its main role is to prevent the robot from tumbling around the   axis.The torque that prevents the robot from tipping around the   axis can only be provided by the next set of support feet.The contact forces at   ′ and   ′ also play a supplementary role in maintaining rotational stability around   .The rotation of the robot around   is stopped by the friction at the current support and next support lines, but this is not the main reason the robot loses stability.ZMP0 can be considered to be on the extension line of the equivalent force,   , acting at the CoM of the body (Figure 5).In the     plane,  1 represents the point with the largest   -coordinate among the four support points (, ,   ′ , and   ′ ), whereas the   -coordinate of  4 is the smallest one, as shown in Figure 6.In the x s Sz s plane, D 1 represents the point with the largest z s -coordinate among the four support points (A, B, A L , and B L ), whereas the z s -coordinate of D 4 is the smallest one, as shown in Figure 6.For any moment, is the angle between   and the support plane; hence, The ground-reaction force at the current or the next set of support feet is unable to maintain the rotational balance of the robot around the   axis if   0  is less than the minimum   -coordinate of the four support points, or larger than the maximum value, because the direction of the support forces can only be perpendicular to the support plane and point outward.Similarly, the ground-reaction force at the next set of support feet cannot prevent the robot from tumbling around the   axis if   0  is less than the minimum or larger than the maximum   -coordinate of the four support points.However, it is insufficient to consider only the rotation around these two directions for stability analysis.The direction of friction is determined, and the magnitude is also limited.Hence, in some specific cases, it is difficult for the ground to provide sufficient friction to push the robot forward and to maintain the rotational balance around the   axis.Besides, it is necessary to consider whether the robot can maintain stability easily when the next support line is formed and the current one disappears because the two support lines do not always exist simultaneously.Therefore, referring to the traditional stability criteria based on ZMP, this paper holds that the robot stays stable if ZMP0 is located within the convex quadrilateral formed by , ,   ′ , and   ′ .This quadrilateral is called the virtual-support quadrilateral.The quadrilaterals   ′   ′ or   ′   ′ are formed depending on whether  intersects   ′   ′ , as shown in Figure 7.For any moment, τ xz is the angle between F exz and the support plane; hence, The ground-reaction force at the current or the next set of support feet is unable to maintain the rotational balance of the robot around the y s axis if z s ZMP 0 is less than the minimum z s -coordinate of the four support points, or larger than the maximum value, because the direction of the support forces can only be perpendicular to the support plane and point outward.Similarly, the ground-reaction force at the next set of support feet cannot prevent the robot from tumbling around the z s axis if y s ZMP 0 is less than the minimum or larger than the maximum y s -coordinate of the four support points.
However, it is insufficient to consider only the rotation around these two directions for stability analysis.The direction of friction is determined, and the magnitude is also limited.Hence, in some specific cases, it is difficult for the ground to provide sufficient friction to push the robot forward and to maintain the rotational balance around the x s axis.Besides, it is necessary to consider whether the robot can maintain stability easily when the next support line is formed and the current one disappears because the two support lines do not always exist simultaneously.Therefore, referring to the traditional stability criteria based on ZMP, this paper holds that the robot stays stable if ZMP 0 is located within the convex quadrilateral formed by A, B, A L , and B L .This quadrilateral is called the virtual-support quadrilateral.The quadrilaterals AA L BB L or ABB L A L are formed depending on whether AB intersects A L B L , as shown in Figure 7.For any moment, is the angle between   and the support plane; hence, The ground-reaction force at the current or the next set of support feet is unable to maintain the rotational balance of the robot around the   axis if   0  is less than the minimum   -coordinate of the four support points, or larger than the maximum value, because the direction of the support forces can only be perpendicular to the support plane and point outward.Similarly, the ground-reaction force at the next set of support feet cannot prevent the robot from tumbling around the   axis if   0  is less than the minimum or larger than the maximum   -coordinate of the four support points.However, it is insufficient to consider only the rotation around these two directions for stability analysis.The direction of friction is determined, and the magnitude is also limited.Hence, in some specific cases, it is difficult for the ground to provide sufficient friction to push the robot forward and to maintain the rotational balance around the   axis.Besides, it is necessary to consider whether the robot can maintain stability easily when the next support line is formed and the current one disappears because the two support lines do not always exist simultaneously.Therefore, referring to the traditional stability criteria based on ZMP, this paper holds that the robot stays stable if ZMP0 is located within the convex quadrilateral formed by , ,   ′ , and   ′ .This quadrilateral is called the virtual-support quadrilateral.The quadrilaterals   ′   ′ or   ′   ′ are formed depending on whether  intersects   ′   ′ , as shown in Figure 7.

Dynamic-Stability Measurement
The next problem to be solved is how to evaluate stability if the dynamic motion of the robot is stable.It is clear that the robot rarely tips over if the difference between ground-reaction forces at the support feet is small.
The following three quantities are taken into account to quantify the stability of the robot: (1) Distances between ZMP 0 and the boundaries of the virtual-support quadrilateral in the virtual-support plane.
These distances not only show the stability of the robot with its current support line but also consider the stability after the next support line appears.The vertical distance between ZMP 0 and each boundary is a comprehensive representation of the distances between ZMP 0 and the two support feet forming the boundary.From Equation ( 25), it can be seen that, in a certain direction, larger minimum distances between ZMP 0 and the two support feet located at both ends result in smaller differences between the ground-reaction forces at these support feet.The instability of the robot is mainly caused by unbalanced rotation in two mutually perpendicular directions.Therefore, the distances between ZMP 0 and the two sets of opposite boundaries (d 11 and d 12 , and d 21 and d 22 ) should be considered when evaluating the stability of the robot (Figure 8).

Dynamic-Stability Measurement
The next problem to be solved is how to evaluate stability if the dynamic motion of the robot is stable.It is clear that the robot rarely tips over if the difference between ground-reaction forces at the support feet is small.
The following three quantities are taken into account to quantify the stability of the robot: (1) Distances between ZMP0 and the boundaries of the virtual-support quadrilateral in the virtualsupport plane.
These distances not only show the stability of the robot with its current support line but also consider the stability after the next support line appears.The vertical distance between ZMP0 and each boundary is a comprehensive representation of the distances between ZMP0 and the two support feet forming the boundary.From Equation ( 25), it can be seen that, in a certain direction, larger minimum distances between ZMP0 and the two support feet located at both ends result in smaller differences between the ground-reaction forces at these support feet.The instability of the robot is mainly caused by unbalanced rotation in two mutually perpendicular directions.Therefore, the distances between ZMP0 and the two sets of opposite boundaries ( 11 and  12 , and  21 and  22 ) should be considered when evaluating the stability of the robot (Figure 8).The two support points forming a certain boundary of the quadrilateral do not necessarily exist at present.Hence, the distances between ZMP0 and boundaries should be multiplied by their corresponding coefficients   , which are defined as follows: The integrated minimum distance between ZMP0 and the boundaries is defined as is expected to be as large as possible, and it is the most important quantity for evaluating the robot's stability.
(2) Angle between the vector pointing from CoM to ZMP0 and the normal vector of the virtualsupport plane.
The ground-reaction force can be divided into two parts: support and friction forces.The support force is perpendicular to the support plane, and, in theory, it can be provided sufficiently to guarantee the force balance in this direction.The direction of the friction force is opposite to that of the trend of the foot's motion with respect to the support plane if the foot does not slip; the largest static-friction force that the ground can provide is proportional to the magnitude of the support force.
Under a constant equivalent force   , larger angles, , between the vector pointing from CoM to ZMP0 (the direction is the same as   ) and the normal vector of the virtual-support plane (Figure 9) result in lower support forces and greater required friction forces; thus, the feet are more likely to slip, and the stability of the robot becomes worse.|sin | is used as a key quantity to measure the robot's stability, in which  represents the vector pointing from CoM to ZMP0: The two support points forming a certain boundary of the quadrilateral do not necessarily exist at present.Hence, the distances between ZMP 0 and boundaries should be multiplied by their corresponding coefficients ω i , which are defined as follows: two support f eet both exist 0.5, only one support f oot exists The integrated minimum distance between ZMP 0 and the boundaries is defined as d edge is expected to be as large as possible, and it is the most important quantity for evaluating the robot's stability.
(2) Angle between the vector pointing from CoM to ZMP 0 and the normal vector of the virtual-support plane.
The ground-reaction force can be divided into two parts: support and friction forces.The support force is perpendicular to the support plane, and, in theory, it can be provided sufficiently to guarantee the force balance in this direction.The direction of the friction force is opposite to that of the trend of the foot's motion with respect to the support plane if the foot does not slip; the largest static-friction force that the ground can provide is proportional to the magnitude of the support force.
Under a constant equivalent force F e , larger angles, τ, between the vector pointing from CoM to ZMP 0 (the direction is the same as F e ) and the normal vector of the virtual-support plane (Figure 9) result in lower support forces and greater required friction forces; thus, the feet are more likely to slip, and the stability of the robot becomes worse.|sin τ| is used as a key quantity to measure the robot's stability, in which u represents the vector pointing from CoM to ZMP 0 : (3) Distance between ZMP0 and the support line.
For dynamic gaits, the support lines are alternately formed so that the ground-reaction force can push the robot to move forward and prevent it from falling.What really interacts with the ground are the support legs; therefore, the distance between ZMP0 and the support lines is also very important for the evaluation of robot stability.The larger is the distance between ZMP0 and a certain support line, the easier it is for the robot to tip over around this support line.
The stability measurement may change suddenly when the current support line disappears and the next support line forms if only the vertical distance between ZMP0 and the current support line (denoted by   ) is considered.The vertical distance between ZMP0 and the next support line is represented by   , as shown in Figure 10.In practice,   should be as short as possible during the first half of the support phase, and   should be as short as possible during the second half.The distance between ZMP0 and the support line is defined as where   and   are time-varying weights: The weights are illustrated in Figure 11.(3) Distance between ZMP 0 and the support line.
For dynamic gaits, the support lines are alternately formed so that the ground-reaction force can push the robot to move forward and prevent it from falling.What really interacts with the ground are the support legs; therefore, the distance between ZMP 0 and the support lines is also very important for the evaluation of robot stability.The larger is the distance between ZMP 0 and a certain support line, the easier it is for the robot to tip over around this support line.
The stability measurement may change suddenly when the current support line disappears and the next support line forms if only the vertical distance between ZMP 0 and the current support line (denoted by d C ) is considered.The vertical distance between ZMP 0 and the next support line is represented by d L , as shown in Figure 10.In practice, d C should be as short as possible during the first half of the support phase, and d L should be as short as possible during the second half.
(3) Distance between ZMP0 and the support line.
For dynamic gaits, the support lines are alternately formed so that the ground-reaction force can push the robot to move forward and prevent it from falling.What really interacts with the ground are the support legs; therefore, the distance between ZMP0 and the support lines is also very important for the evaluation of robot stability.The larger is the distance between ZMP0 and a certain support line, the easier it is for the robot to tip over around this support line.
The stability measurement may change suddenly when the current support line disappears and the next support line forms if only the vertical distance between ZMP0 and the current support line (denoted by   ) is considered.The vertical distance between ZMP0 and the next support line is represented by   , as shown in Figure 10.In practice,   should be as short as possible during the first half of the support phase, and   should be as short as possible during the second half.The distance between ZMP0 and the support line is defined as where   and   are time-varying weights: The weights are illustrated in Figure 11.The distance between ZMP 0 and the support line is defined as where ω C and ω L are time-varying weights: The weights are illustrated in Figure 11.Using the above three quantities, this paper proposes the stability measurement   : A larger   value indicates better motion stability.The value of the coefficients   should consider: (1) unification of the order of magnitude of the three terms on the right-hand side of Equation (33) (which are related to the size of the robot); and (2) the importance of each parameter.According to the definitions of   , |sin |, and   , for the dynamically stable gait: where   can be calculated by: According to the Equations ( 34)-(39), which indicates the possible minimum value of   , but this does not mean that the motion is stable if   > 0.35( + ).

Results and Discussion
The simplified and mechanical models of the quadruped robot used in this paper are shown in Figure 12.Each leg of the robot has three degrees of freedom.Using the above three quantities, this paper proposes the stability measurement M stability : A larger M stability value indicates better motion stability.The value of the coefficients λ i should consider: (1) unification of the order of magnitude of the three terms on the right-hand side of Equation (33) (which are related to the size of the robot); and (2) the importance of each parameter.According to the definitions of d edge , |sin τ|, and d spt , for the dynamically stable gait: where λ i can be calculated by: According to the Equations ( 34)-(39), which indicates the possible minimum value of M stability , but this does not mean that the motion is stable if M stability > 0.35(L + W).

Results and Discussion
The simplified and mechanical models of the quadruped robot used in this paper are shown in Figure 12.Each leg of the robot has three degrees of freedom.
which indicates the possible minimum value of   , but this does not mean that the motion is stable if   > 0.35( + ).

Results and Discussion
The simplified and mechanical models of the quadruped robot used in this paper are shown in Figure 12.Each leg of the robot has three degrees of freedom.The geometric and weight parameters are shown in Table 1.The moments of inertia of the body around the principal axes are J xx = 0.74, J yy = 0.56, and J zz = 0.22 (kg•m 2 ).
A simulation model of the quadruped robot was built in the dynamic-simulation software ADAMS.The contact parameters between each foot and ground are shown in Table 2 (referring to the contact parameters between dry rubber and dry asphalt).Trotting is the most commonly used dynamic gait and has been the focus of previous studies; thus, this study used the trot gait as the basic experimental gait.The trot-gait parameters for all simulations are shown in Table 3.

Gait Parameters (Unit) Motion on Horizontal Plane
Step length (m) 0.2 Maximum height for the tip of swing legs (m) 0.03 Duration of a gait cycle (s) 1 Velocity (m/s) 0.2 Duty factor 0.6 Three quadruped-robot offline-gait simulations were designed as follow: • Simulation 1 verified that the quadruped robot can maintain its stable state during dynamic motion when ZMP does not lie on the current support line.

•
Simulation 2 verified that the stability criterion and measurement proposed in this paper can determine the robot's stability state efficiently and evaluate the stability accurately compared with the stability margins that only consider the distance between ZMP and the current support line and the stability analysis method called LAR.

•
Simulation 3 showed that the stability-analysis method proposed in this paper is applicable to motion in 3D space.

Simulation 1
Most existing stability criteria for trotting gaits assume that the robot is in stable motion only when ZMP is located on the support line formed by the current support feet.However, this is inconsistent with the movement of animals in nature.To illustrate that these stability criteria cannot accurately determine the stability state of the robot, the following simulation was designed.
The projection of CoM in the direction of gravity onto the virtual-support plane is denoted by CoG.In the trot-gait-motion planning used in this simulation, the CoM of the body moves forward at a constant speed and, at each phase-changing point, CoG is located at the center of the quadrilateral formed by the just-disappearing support line and the newly formed line (see Figure 13).

Simulation 1
Most existing stability criteria for trotting gaits assume that the robot is in stable motion only when ZMP is located on the support line formed by the current support feet.However, this is inconsistent with the movement of animals in nature.To illustrate that these stability criteria cannot accurately determine the stability state of the robot, the following simulation was designed.
The projection of CoM in the direction of gravity onto the virtual-support plane is denoted by CoG.In the trot-gait-motion planning used in this simulation, the CoM of the body moves forward at a constant speed and, at each phase-changing point, CoG is located at the center of the quadrilateral formed by the just-disappearing support line and the newly formed line (see Figure 13).The gray lines in Figure 14 are the support lines and are shown as several groups of line clusters for two reasons.First, the active or passive compliance was not considered in this simulation process, therefore the feet were bouncing for a short period after landing.Second, the feet were slipping during the support phases.The blue line is the actual motion trajectory of CoG on the horizontal plane.The motion trajectory of CoG was not completely consistent with the planned trajectory under the influence of the unbalanced moment around the current support line, the legs' movement, and the impact force at the feet.The red line is the trajectory of traditional ZMP.It can be clearly seen that, during this simulation, ZMP was not always or mostly located on the current support line.
As shown in Figure 15, the running direction of the body was slightly skewed (see yaw angle), and the absolute values of the pitch and roll angles were less than 0.11 rad (6°).The curves of these The simulation animation is in Video S1.The position changes of the support lines, CoG, and ZMP on the support plane (z i Oy i plane) during the entire simulation process are shown in Figure 14.

Simulation 1
Most existing stability criteria for trotting gaits assume that the robot is in stable motion only when ZMP is located on the support line formed by the current support feet.However, this is inconsistent with the movement of animals in nature.To illustrate that these stability criteria cannot accurately determine the stability state of the robot, the following simulation was designed.
The projection of CoM in the direction of gravity onto the virtual-support plane is denoted by CoG.In the trot-gait-motion planning used in this simulation, the CoM of the body moves forward at a constant speed and, at each phase-changing point, CoG is located at the center of the quadrilateral formed by the just-disappearing support line and the newly formed line (see Figure 13).The gray lines in Figure 14 are the support lines and are shown as several groups of line clusters for two reasons.First, the active or passive compliance was not considered in this simulation process, therefore the feet were bouncing for a short period after landing.Second, the feet were slipping during the support phases.The blue line is the actual motion trajectory of CoG on the horizontal plane.The motion trajectory of CoG was not completely consistent with the planned trajectory under the influence of the unbalanced moment around the current support line, the legs' movement, and the impact force at the feet.The red line is the trajectory of traditional ZMP.It can be clearly seen that, during this simulation, ZMP was not always or mostly located on the current support line.
As shown in Figure 15, the running direction of the body was slightly skewed (see yaw angle), and the absolute values of the pitch and roll angles were less than 0.11 rad (6°).The curves of these The gray lines in Figure 14 are the support lines and are shown as several groups of line clusters for two reasons.First, the active or passive compliance was not considered in this simulation process, therefore the feet were bouncing for a short period after landing.Second, the feet were slipping during the support phases.The blue line is the actual motion trajectory of CoG on the horizontal plane.The motion trajectory of CoG was not completely consistent with the planned trajectory under the influence of the unbalanced moment around the current support line, the legs' movement, and the impact force at the feet.The red line is the trajectory of traditional ZMP.It can be clearly seen that, during this simulation, ZMP was not always or mostly located on the current support line.
As shown in Figure 15, the running direction of the body was slightly skewed (see yaw angle), and the absolute values of the pitch and roll angles were less than 0.11 rad (6 • ).The curves of these attitude angles show that the robot moved stably and had no obvious tendency to tumble.This indicates that the existing stability criteria for the trot gait cannot accurately determine the stability state of the robot.

Simulation 2
A second set of comparative simulations was designed to verify the following: (1) it is difficult to correctly evaluate the stability by only considering the distance between ZMP and the current support line; (2) the accuracy and effectivity of LAR are worse than the stability-evaluation method proposed in this paper; and (3) the stability-evaluation method proposed in this paper can satisfy the first two requirements of a good-stability criterion and evaluation method presented in Section 1.
In addition to the motion trajectory of the CoG used in Simulation 1 (denoted as Test 1), two other CoG trajectories were used for off-line planning of the quadruped robot's movement (the actual trajectories of CoG did not necessarily coincide with the planned trajectory); refer to the animation (Video S2).


Test 2: CoG moves around the current support line (the velocity parallel to   is 0.2 m/s, and the velocity and acceleration parallel to   at the two ends of this direction are zero), as shown in Figure 16.

Simulation 2
A second set of comparative simulations was designed to verify the following: (1) it is difficult to correctly evaluate the stability by only considering the distance between ZMP and the current support line; (2) the accuracy and effectivity of LAR are worse than the method proposed in this paper; and (3) the stability-evaluation method proposed in this paper can satisfy the first two requirements of a good-stability criterion and evaluation method presented in Section 1.
In addition to the motion trajectory of the CoG used in Simulation 1 (denoted as Test 1), two other CoG trajectories were used for off-line planning of the quadruped robot's movement (the actual trajectories of CoG did not necessarily coincide with the planned trajectory); refer to the animation (Video S2).

•
Test 2: CoG moves around the current support line (the velocity parallel to z i is 0.2 m/s, and the velocity and acceleration parallel to y i at the two ends of this direction are zero), as shown in Figure 17.

•
Test 3: CoG starts to move from the midpoint of the current support line to that of the next support line at a constant speed when the previous support line has just disappeared, as shown in Figure 16.

Simulation 2
A second set of comparative simulations was designed to verify the following: (1) it is difficult to correctly evaluate the stability by only considering the distance between ZMP and the current support line; (2) the accuracy and effectivity of LAR are worse than the stability-evaluation method proposed in this paper; and (3) the stability-evaluation method proposed in this paper can satisfy the first two requirements of a good-stability criterion and evaluation method presented in Section 1.
In addition to the motion trajectory of the CoG used in Simulation 1 (denoted as Test 1), two other CoG trajectories were used for off-line planning of the quadruped robot's movement (the actual trajectories of CoG did not necessarily coincide with the planned trajectory); refer to the animation (Video S2).


Test 2: CoG moves around the current support line (the velocity parallel to   is 0.2 m/s, and the velocity and acceleration parallel to   at the two ends of this direction are zero), as shown in Figure 16.The posture angles, positions, and velocities of the body for the three tests are shown in Figure 18a-c, respectively.The driving moments of the three joints of the left hind leg are shown in Figure 18d.trajectories of CoG did not necessarily coincide with the planned trajectory); refer to the animation (Video S2).


Test 2: CoG moves around the current support line (the velocity parallel to   is 0.2 m/s, and the velocity and acceleration parallel to   at the two ends of this direction are zero), as shown in Figure 16.The wave heights of the data curves referring to the difference between adjacent peaks and troughs show the severity of data fluctuations if the fluctuation period is similar.The wave-heightvariation range of each dataset and the driving-torque-variation ranges of the three joints of the left hind leg are shown in Figure 19.As shown in Figure 19, the motion of Test 3 was the steadiest and required the minimum jointdriving torque, and the stabilities of the other two tests were about the same.In fact, stability cannot be evaluated accurately using individual-sensor data, which is one of the reasons the stability criterion and measurement must be proposed.
The generally accepted stability-evaluation method (denoted as Method 1), i.e., the one that measures the distance between the traditional ZMP and the current support line; LAR (denoted as Method 2); and the proposed stability measurement for the trot gait of a quadruped robot (denoted as Method 3) were used to evaluate the stability of the three tests (Figure 20a-c).As shown in Figure 19, the motion of Test 3 was the steadiest and required the minimum joint-driving torque, and the stabilities of the other two tests were about the same.In fact, stability cannot be evaluated accurately using individual-sensor data, which is one of the reasons the stability criterion and measurement must be proposed.
The generally accepted stability-evaluation method (denoted as Method 1), i.e., the one that measures the distance between the traditional ZMP and the current support line; LAR (denoted as Method 2); and the proposed stability measurement for the trot gait of a quadruped robot (denoted as Method 3) were used to evaluate the stability of the three tests (Figure 20a-c).As shown in Figure 19, the motion of Test 3 was the steadiest and required the minimum jointdriving torque, and the stabilities of the other two tests were about the same.In fact, stability cannot be evaluated accurately using individual-sensor data, which is one of the reasons the stability criterion and measurement must be proposed.
The generally accepted stability-evaluation method (denoted as Method 1), i.e., the one that measures the distance between the traditional ZMP and the current support line; LAR (denoted as Method 2); and the proposed stability measurement for the trot gait of a quadruped robot (denoted as Method 3) were used to evaluate the stability of the three tests (Figure 20a-c).The average values (over 15 s) of the stability evaluations obtained using these three methods are shown in Table 4. Test 2 had the lowest stability, followed by Test 3 and then Test 1 according to Method 1.However, this result was not consistent with the actual movement of the robot.Method 2 indicated that the stability of the robot in Test 3 was better than Test 2, followed by Test 1.Using the second method, the stability of Test 1 was found to be slightly better than that of Test 2, with Test 3 being the most stable case.The last two methods could both indicate that the stability of Test 3 was best, but one more simulation should be designed to compare the stability of Test 1 and Test 2.
Two impact forces parallel to   in the opposite direction were both applied to the CoM of the body during the trotting in Test 1 and Test 2. The magnitudes, timings, and durations of these forces are shown in Table 5.For example, in the case shown in Figure 21, the impact forces in the positive and negative directions of   can both drive ZMP0 closer to the boundaries of the virtual-support quadrilateral.However, the force in the positive direction can reduce the distance between ZMP and the support The average values (over 15 s) of the stability evaluations obtained using these three methods are shown in Table 4. Test 2 had the lowest stability, followed by Test 3 and then Test 1 according to Method 1.However, this result was not consistent with the actual movement of the robot.Method 2 indicated that the stability of the robot in Test 3 was better than Test 2, followed by Test 1.Using the second method, the stability of Test 1 was found to be slightly better than that of Test 2, with Test 3 being the most stable case.The last two methods could both indicate that the stability of Test 3 was best, but one more simulation should be designed to compare the stability of Test 1 and Test 2.
Two impact forces parallel to y c in the opposite direction were both applied to the CoM of the body during the trotting in Test 1 and Test 2. The magnitudes, timings, and durations of these forces are shown in Table 5.The robot is supported by the same set of diagonal legs during the first half of each gait cycle of the trotting gait.The magnitude of the forces that the robot can bear differs in the positive and negative directions of y c if the positions of ZMP 0 and the next support line are fixed.
For example, in the case shown in Figure 21, the impact forces in the positive and negative directions of y c can both drive ZMP 0 closer to the boundaries of the virtual-support quadrilateral.However, the force in the positive direction can reduce the distance between ZMP and the support lines, whereas the force in the negative direction increases this distance.Thus, the robot can endure a greater positive impact force along y c .The average values (over 15 s) of the stability evaluations obtained using these three methods are shown in Table 4. Test 2 had the lowest stability, followed by Test 3 and then Test 1 according to Method 1.However, this result was not consistent with the actual movement of the robot.Method 2 indicated that the stability of the robot in Test 3 was better than Test 2, followed by Test 1.Using the second method, the stability of Test 1 was found to be slightly better than that of Test 2, with Test 3 being the most stable case.The last two methods could both indicate that the stability of Test 3 was best, but one more simulation should be designed to compare the stability of Test 1 and Test 2.
Two impact forces parallel to   in the opposite direction were both applied to the CoM of the body during the trotting in Test 1 and Test 2. The magnitudes, timings, and durations of these forces are shown in Table 5.For example, in the case shown in Figure 21, the impact forces in the positive and negative directions of   can both drive ZMP0 closer to the boundaries of the virtual-support quadrilateral.However, the force in the positive direction can reduce the distance between ZMP and the support The robot in Test 2 tumbled to the ground after those two impact forces acting in different directions, while the robot in Test 3 could sustain the impacts without losing its stability.This means that the stability of the motion in Test 1 was better than Test 2, thus the method proposed in this paper could evaluate the stability of the quadruped robot with trotting gait more accurately than LAR.
Three impact forces parallel to   were applied to the CoM of the body in the three trotting tests, as described above, to further illustrate the accuracy of the proposed stability criterion and measurement.The magnitudes, timings, and durations of these forces are shown in Table 6.Table 6.Impact forces applied to the body.Test 2 was the least able to bear the impact force in the positive direction of   , since the positions of the next set of support feet differed with respect to ZMP0 (as can be seen simply from the positions The robot in Test 2 tumbled to the ground after those two impact forces acting in different directions, while the robot in Test 3 could sustain the impacts without losing its stability.This means that the stability of the motion in Test 1 was better than Test 2, thus the method proposed in this paper could evaluate the stability of the quadruped robot with trotting gait more accurately than LAR.
Three impact forces parallel to y c were applied to the CoM of the body in the three trotting tests, as described above, to further illustrate the accuracy of the proposed stability criterion and measurement.The magnitudes, timings, and durations of these forces are shown in Table 6.The simulation results of these three tests after impacts are shown in Figure 23.The robot in Test 2 tumbled to the ground after those two impact forces acting in different directions, while the robot in Test 3 could sustain the impacts without losing its stability.This means that the stability of the motion in Test 1 was better than Test 2, thus the method proposed in this paper could evaluate the stability of the quadruped robot with trotting gait more accurately than LAR.
Three impact forces parallel to   were applied to the CoM of the body in the three trotting tests, as described above, to further illustrate the accuracy of the proposed stability criterion and measurement.The magnitudes, timings, and durations of these forces are shown in Table 6.Table 6.Impact forces applied to the body.Test 2 was the least able to bear the impact force in the positive direction of   , since the positions of the next set of support feet differed with respect to ZMP0 (as can be seen simply from the positions Test 2 was the least able to bear the impact force in the positive direction of y c , since the positions of the next set of support feet differed with respect to ZMP 0 (as can be seen simply from the positions of the next support lines with respect to CoG), tumbling after the first impact (Figure 23b).The robot walked stably in Test 1 until the second impact was applied to the body (Figure 23a), and the robot in Test 3 lost its balance after the last impact (Figure 23c).
The three stability-evaluation methods were also used for the analysis of the three tests with external impacts.To illustrate this intuitively, for Method 2, the stability measurement was set to −0.3 if ZMP 0 was outside the virtual-support quadrilateral and the robot lost its balance according to Equation (40).
As shown in Figure 24, the distance between ZMP and the current support line was not zero for most of the time, making it very difficult to determine when the robot would lose its balance.Using Method 2, it was hard to decide when to update the value of LAR if the robot tips over, because one of the feet would never touch the ground.From the results of Method 3, the robot was considered to lose its balance at 9.753 s (Test 1), 5.109 s (Test 2), and 14.2 s (Test 3), which was consistent with the actual motions of the robot. of the next support lines with respect to CoG), tumbling after the first impact (Figure 23b).The robot walked stably in Test 1 until the second impact was applied to the body (Figure 23a), and the robot in Test 3 lost its balance after the last impact (Figure 23c).The three stability-evaluation methods were also used for the analysis of the three tests with external impacts.To illustrate this intuitively, for Method 2, the stability measurement was set to −0.3 if ZMP0 was outside the virtual-support quadrilateral and the robot lost its balance according to Equation (40).
As shown in Figure 24, the distance between ZMP and the current support line was not zero for most of the time, making it very difficult to determine when the robot would lose its balance.Using Method 2, it was hard to decide when to update the value of LAR if the robot tips over, because one of the feet would never touch the ground.From the results of Method 3, the robot was considered to lose its balance at 9.753 s (Test 1), 5.109 s (Test 2), and 14.2 s (Test 3), which was consistent with the actual motions of the robot.The simulation data from Test 3 were analyzed separately.Here, the three posture angles of the principal-axes coordinate system and the velocities of the body in the y i and z i directions are shown in Figures 25 and 26, respectively. of the next support lines with respect to CoG), tumbling after the first impact (Figure 23b).The robot walked stably in Test 1 until the second impact was applied to the body (Figure 23a), and the robot in Test 3 lost its balance after the last impact (Figure 23c).The three stability-evaluation methods were also used for the analysis of the three tests with external impacts.To illustrate this intuitively, for Method 2, the stability measurement was set to −0.3 if ZMP0 was outside the virtual-support quadrilateral and the robot lost its balance according to Equation (40).
As shown in Figure 24, the distance between ZMP and the current support line was not zero for most of the time, making it very difficult to determine when the robot would lose its balance.Using Method 2, it was hard to decide when to update the value of LAR if the robot tips over, because one of the feet would never touch the ground.From the results of Method 3, the robot was considered to lose its balance at 9.753 s (Test 1), 5.109 s (Test 2), and 14.2 s (Test 3), which was consistent with the actual motions of the robot.The impact forces caused large changes in the velocities along the   and   directions since the current support lines and these forces were not vertical; thus, the posture angles of the body exceeded the angle-fluctuation range of normal stable motion, reducing the stability of the robot.After the first two impacts, the current ground-reaction force prevented the continued increase of undesired speeds in these two directions to a certain extent since ZMP0 was located inside the virtual-support quadrilateral, but did not affect the speed in the direction perpendicular to the support line.This part of speed could only be eliminated when the next set of support feet hit the ground, causing a great impact force at the next set of support feet, such that the stability of the robot was still poor and it had to take several steps for the robot to restore normal stable state.Following the last impact, the speed in the   direction, the pitch angle, and the roll angle of the body changed sharply, and the ground-reaction force at the support feet could not prevent the robot from tipping over, and the robot lost its balance.
It can be seen in Figure 24 that, after the first impact, the distance between ZMP and the current support line did not exceed the range of the normal stable trotting, and, after the second impact, it was not until 9.8 s that the evaluated stability of Method 1 decreased markedly.Method 2 could show the reduction of the robot's stability correctly after impacts were applied in this test, but it could not evaluate the stability in time, because LAR was updated when the two feet in the same phase both touch the ground.However, from the curve of Method 3, the stability measurement of the robot decreased very quickly and markedly after the first two impacts and then fluctuated significantly within a few steps after impact.After the last impact, the robot lost its balance.This was consistent with the movement previously analyzed.
The positions of ZMP and ZMP0 were both calculated in this simulation to illustrate that the modified ZMP with a velocity term can reflect the motion state better than the traditional ZMP. Figure 27 shows the curves of ZMP, ZMP0, and CoG in the   direction.After each impact,   of ZMP0 could show the movement trend faster than that of ZMP and CoG, and there was not much difference between the peaks (or troughs) of ZMP and ZMP0 during the whole movement, meaning that the measurement using ZMP0 could decrease response time without misjudging the stability state of the robot.As shown in Figure 28, after the final impact, the stability criterion using ZMP0 determined that the robot lost its balance at 14.2 s, whereas the criterion using ZMP showed instability at 14.5 s.The time difference between the two adjacent steps of the trot The impact forces caused large changes in the velocities along the y i and z i directions since the current support lines and these forces were not vertical; thus, the posture angles of the body exceeded the angle-fluctuation range of normal stable motion, reducing the stability of the robot.After the first two impacts, the current ground-reaction force prevented the continued increase of undesired speeds in these two directions to a certain extent since ZMP 0 was located inside the virtual-support quadrilateral, but did not affect the speed in the direction perpendicular to the support line.This part of speed could only be eliminated when the next set of support feet hit the ground, causing a great impact force at the next set of support feet, such that the stability of the robot was still poor and it had to take several steps for the robot to restore normal stable state.Following the last impact, the speed in the y i direction, the pitch angle, and the roll angle of the body changed sharply, and the ground-reaction force at the support feet could not prevent the robot from tipping over, and the robot lost its balance.
It can be seen in Figure 24 that, after the first impact, the distance between ZMP and the current support line did not exceed the range of the normal stable trotting, and, after the second impact, it was not until 9.8 s that the evaluated stability of Method 1 decreased markedly.Method 2 could show the reduction of the robot's stability correctly after impacts were applied in this test, but it could not evaluate the stability in time, because LAR was updated when the two feet in the same phase both touch the ground.However, from the curve of Method 3, the stability measurement of the robot decreased very quickly and markedly after the first two impacts and then fluctuated significantly within a few steps after impact.After the last impact, the robot lost its balance.This was consistent with the movement previously analyzed.
The positions of ZMP and ZMP 0 were both calculated in this simulation to illustrate that the modified ZMP with a velocity term can reflect the motion state better than the traditional ZMP. Figure 27 shows the curves of ZMP, ZMP 0 , and CoG in the y i direction.The impact forces caused large changes in the velocities along the   and   directions since the current support lines and these forces were not vertical; thus, the posture angles of the body exceeded the angle-fluctuation range of normal stable motion, reducing the stability of the robot.After the first two impacts, the current ground-reaction force prevented the continued increase of undesired speeds in these two directions to a certain extent since ZMP0 was located inside the virtual-support quadrilateral, but did not affect the speed in the direction perpendicular to the support line.This part of speed could only be eliminated when the next set of support feet hit the ground, causing a great impact force at the next set of support feet, such that the stability of the robot was still poor and it had to take several steps for the robot to restore normal stable state.Following the last impact, the speed in the   direction, the pitch angle, and the roll angle of the body changed sharply, and the ground-reaction force at the support feet could not prevent the robot from tipping over, and the robot lost its balance.
It can be seen in Figure 24 that, after the first impact, the distance between ZMP and the current support line did not exceed the range of the normal stable trotting, and, after the second impact, it was not until 9.8 s that the evaluated stability of Method 1 decreased markedly.Method 2 could show the reduction of the robot's stability correctly after impacts were applied in this test, but it could not evaluate the stability in time, because LAR was updated when the two feet in the same phase both touch the ground.However, from the curve of Method 3, the stability measurement of the robot decreased very quickly and markedly after the first two impacts and then fluctuated significantly within a few steps after impact.After the last impact, the robot lost its balance.This was consistent with the movement previously analyzed.
The positions of ZMP and ZMP0 were both calculated in this simulation to illustrate that the modified ZMP with a velocity term can reflect the motion state better than the traditional ZMP. Figure 27 shows the curves of ZMP, ZMP0, and CoG in the   direction.After each impact,   of ZMP0 could show the movement trend faster than that of ZMP and CoG, and there was not much difference between the peaks (or troughs) of ZMP and ZMP0 during the whole movement, meaning that the measurement using ZMP0 could decrease response time without misjudging the stability state of the robot.As shown in Figure 28, after the final impact, the stability criterion using ZMP0 determined that the robot lost its balance at 14.2 s, whereas the criterion using ZMP showed instability at 14.5 s.The time difference between the two adjacent steps of the trot After each impact, y i of ZMP 0 could show the movement trend faster than that of ZMP and CoG, and there was not much difference between the peaks (or troughs) of ZMP and ZMP 0 during the whole movement, meaning that the measurement using ZMP 0 could decrease response time without misjudging the stability state of the robot.As shown in Figure 28, after the final impact, the stability criterion using ZMP 0 determined that the robot lost its balance at 14.2 s, whereas the criterion using ZMP showed instability at 14.5 s.The time difference between the two adjacent steps of the trot gait used in this study was 0.5 s.Such a time difference is crucial for the robot to adjust the landing positions of the next support line in time to maintain its balance during robot-motion control.

Simulation 3
Simulation 3 was designed to prove that the stability-evaluation method proposed in this paper is also applicable to the motion of a robot in 3D space (Figure 29).The robot moved on a 12° slope as shown in the animation (Video S3), and the body's trajectory was similar to that used in Test 3 of the previous simulation: CoG started to move from the midpoint of the current support line to that of the next support line at a constant speed when the previous support line had just disappeared.The contact parameters between each foot and the ground and the trot-gait parameters of this simulation were set to be identical to Test 3.

Simulation 3
Simulation 3 was designed to prove that the stability-evaluation method proposed in this paper is also applicable to the motion of a robot in 3D space (Figure 29).The robot moved on a 12 • slope as shown in the animation (Video S3), and the body's trajectory was similar to that used in Test 3 of the previous simulation: CoG started to move from the midpoint of the current support line to that of the next support line at a constant speed when the previous support line had just disappeared.The contact parameters between each foot and the ground and the trot-gait parameters of this simulation were set to be identical to Test 3.

Simulation 3
Simulation 3 was designed to prove that the stability-evaluation method proposed in this paper is also applicable to the motion of a robot in 3D space (Figure 29).The robot moved on a 12° slope as shown in the animation (Video S3), and the body's trajectory was similar to that used in Test 3 of the previous simulation: CoG started to move from the midpoint of the current support line to that of the next support line at a constant speed when the previous support line had just disappeared.The contact parameters between each foot and the ground and the trot-gait parameters of this simulation were set to be identical to Test 3. The attitude angles of the principal-axes coordinate system, and the variation of the stability measurements calculated using above three methods with time are shown in Figure 30a  The attitude angles of the principal-axes coordinate system, and the variation of the stability measurements calculated using above three methods with time are shown in Figure 30a

Simulation 3
Simulation 3 was designed to prove that the stability-evaluation method proposed in this paper is also applicable to the motion of a robot in 3D space (Figure 29).The robot moved on a 12° slope as shown in the animation (Video S3), and the body's trajectory was similar to that used in Test 3 of the previous simulation: CoG started to move from the midpoint of the current support line to that of the next support line at a constant speed when the previous support line had just disappeared.The contact parameters between each foot and the ground and the trot-gait parameters of this simulation were set to be identical to Test 3.  The movement of the robot trotting on the slope was steady.The average values (over 15 s) of the stability evaluations obtained using these three methods of Simulation 3 and Test 3 are listed in Table 7. Method 1 showed that those two motions had almost the same stability.Method 2 indicated that the stability of the motion on slope was much higher than that on horizontal plane, which was contrary to the result of Method 3. Three same impact forces (as shown in Table 6) were applied to the CoM of the body during Simulation 3. The posture angles of the principal-axes coordinate system, the positions and velocities of the body are shown in Figure 31a-c, respectively.The movement of the robot trotting on the slope was steady.The average values (over 15 s) of the stability evaluations obtained using these three methods of Simulation 3 and Test 3 are listed in Table 7.
Table 7.Average values of stability evaluations when the robot trotting on the 12° slope and horizontal plane.Method 1 showed that those two motions had almost the same stability.Method 2 indicated that the stability of the motion on slope was much higher than that on horizontal plane, which was contrary to the result of Method 3. Three same impact forces (as shown in Table 6) were applied to the CoM of the body during Simulation 3. The posture angles of the principal-axes coordinate system, the positions and velocities of the body are shown in Figure 31a-c The robot did not lose its balance before the last impact, but its feet skidded on the slope after the first impact because it was hard for the ground to provide enough friction to push the robot upward.The wave heights of the posture angles were not able to recover to the states before impact.Therefore, the robot was more likely to lose its balance on the slope, and the stability of this motion was worse than in Test 3. Three Methods were used to analyze the motion on slope with impacts (Figure 32).The stability-evaluation method proposed in this paper could measure the stability of motion in 3D space more accurately and efficiently than Methods 1 and 2. The robot did not lose its balance before the last impact, but its feet skidded on the slope after the first impact because it was hard for the ground to provide enough friction to push the robot upward.The wave heights of the posture angles were not able to recover to the states before impact.Therefore, the robot was more likely to lose its balance on the slope, and the stability of this motion was worse than in Test 3. Three Methods were used to analyze the motion on slope with impacts (Figure 32).The stability-evaluation method proposed in this paper could measure the stability of motion in 3D space more accurately and efficiently than Methods 1 and 2.

)Figure 1 .
Figure 1.Force-schematic diagram of the robot body.

Figure 2 .Figure 1 .
Figure 2. Settings of the coordinate systems.The position of ZMP in the support coordinate system    = [0

Figure 1 .
Figure 1.Force-schematic diagram of the robot body.

Figure 2 .
Figure 2. Settings of the coordinate systems.The position of ZMP in the support coordinate system    = [0

Figure 2 .
Figure 2. Settings of the coordinate systems.

Figure 3 .
Figure 3. Determination of the virtual-support plane.

Figure 4 .
Figure 4. Projection of the next set of support feet and contact force at the current support feet.

Figure 4 .
Figure 4. Projection of the next set of support feet and contact force at the current support feet.
Appl.Sci.2018, 8, x FOR PEER REVIEW 7 of 25 taken into account when evaluating the robot's stability by adding the velocity term into the calculation of ZMP.

Figure 4 .
Figure 4. Projection of the next set of support feet and contact force at the current support feet.

25 Figure 6 .
Figure 6.Force diagram in the     plane.

Figure 6 .
Figure 6.Force diagram in the x s Sz s plane.

25 Figure 6 .
Figure 6.Force diagram in the     plane.

Figure 8 .
Figure 8. Distances between ZMP0 and each boundary of the virtual-support quadrilateral.

Figure 8 .
Figure 8. Distances between ZMP 0 and each boundary of the virtual-support quadrilateral.

Figure 10 .
Figure 10.Distances between  0 and the support lines (  and   ).

Figure 9 .
Figure 9. Angle between F e and x s .

Figure 10 .
Figure 10.Distances between  0 and the support lines (  and   ).

Figure 10 .
Figure 10.Distances between ZMP 0 and the support lines (d C and d L ).

25 Figure 11 .
Figure 11.Weights for the distance between  0 and the two support lines.

Figure 11 .
Figure 11.Weights for the distance between ZMP 0 and the two support lines.

Figure 12 .
Figure 12.(a) Simplified; and (b) mechanical models of the quadruped robot.Figure 12. (a) Simplified; and (b) mechanical models of the quadruped robot.

Figure 12 .
Figure 12.(a) Simplified; and (b) mechanical models of the quadruped robot.Figure 12. (a) Simplified; and (b) mechanical models of the quadruped robot.

Figure 14 .
Figure 14.Positions of the support line and ZMP.

Figure 14 .
Figure 14.Positions of the support line and ZMP.

Figure 14 .
Figure 14.Positions of the support line and ZMP.

25 Figure 15 .
Figure 15.Attitude angles of the body.

Figure 15 .
Figure 15.Attitude angles of the body.

25 Figure 15 .
Figure 15.Attitude angles of the body.

Figure 16 .
Figure 16.Trajectory of the CoG in Test 2.

Figure 17 .
Figure 17.Trajectory of the CoG in Test 3.

Figure 16 .
Figure 16.Trajectory of the CoG in Test 3.

Figure 17 .
Figure 17.Trajectory of the CoG in Test 3.

Figure 18 .
Figure 18.Sensor data of the motion in Simulation 2: (a) posture angles; (b) positions; (c) velocities of the body; and (d) driving moments of the left hind leg.The wave heights of the data curves referring to the difference between adjacent peaks and troughs show the severity of data fluctuations if the fluctuation period is similar.The wave-height-variation range of each dataset and the driving-torque-variation ranges of the three joints of the left hind leg are shown in Figure19.

Figure 19 .
Figure 19.Variation ranges of: (a) wave heights of posture angles; (b) wave heights of positions; (c) wave heights of velocities; and (d) joint torques.

Figure 19 .
Figure 19.Variation ranges of: (a) wave heights of posture angles; (b) wave heights of positions; (c) wave heights of velocities; and (d) joint torques.

Figure 19 .
Figure 19.Variation ranges of: (a) wave heights of posture angles; (b) wave heights of positions; (c) wave heights of velocities; and (d) joint torques.
supported by the same set of diagonal legs during the first half of each gait cycle of the trotting gait.The magnitude of the forces that the robot can bear differs in the positive and negative directions of   if the positions of  0 and the next support line are fixed.

Figure 21 .
Figure 21.Position of  0 and the support lines
supported by the same set of diagonal legs during the first half of each gait cycle of the trotting gait.The magnitude of the forces that the robot can bear differs in the positive and negative directions of   if the positions of  0 and the next support line are fixed.

Figure 21 .
Figure 21.Position of  0 and the support lines

Figure 21 .Figure 22 .
Figure 21.Position of ZMP 0 and the support lines

Figure 22 .
Figure 22.Posture angles of the two tests with impact forces: (a) 100 N; and (b) −80 N.

Figure 22 .
Figure 22.Posture angles of the two tests with impact forces: (a) 100 N; and (b) −80 N.

Figure 25 .
Figure 25.Posture angles of the principal-axes coordinate system.

Figure 25 .
Figure 25.Posture angles of the principal-axes coordinate system.Figure 25.Posture angles of the principal-axes coordinate system.

Figure 25 .
Figure 25.Posture angles of the principal-axes coordinate system.Figure 25.Posture angles of the principal-axes coordinate system.

Figure 26 .
Figure 26.Velocities of the body in the   and   directions.

Figure 27 .
Figure 27.Positions of ZMP, ZMP0, and CoG in the   direction.

Figure 26 .
Figure 26.Velocities of the body in the y i and z i directions.

25 Figure 26 .
Figure 26.Velocities of the body in the   and   directions.

Figure 27 .
Figure 27.Positions of ZMP, ZMP0, and CoG in the   direction.

Figure 27 .
Figure 27.Positions of ZMP, ZMP 0 , and CoG in the y i direction.
Appl.Sci.2018, 8, x FOR PEER REVIEW 21 of 25gait used in this study was 0.5 s.Such a time difference is crucial for the robot to adjust the landing positions of the next support line in time to maintain its balance during robot-motion control.

Figure 29 .
Figure 29.Simulation of Simulation 3.The attitude angles of the principal-axes coordinate system, and the variation of the stability measurements calculated using above three methods with time are shown in Figure30a,b.
Appl.Sci.2018, 8, x FOR PEER REVIEW 21 of 25gait used in this study was 0.5 s.Such a time difference is crucial for the robot to adjust the landing positions of the next support line in time to maintain its balance during robot-motion control.
Appl.Sci.2018, 8, x FOR PEER REVIEW 21 of 25gait used in this study was 0.5 s.Such a time difference is crucial for the robot to adjust the landing positions of the next support line in time to maintain its balance during robot-motion control.

Figure 29 .
Figure 29.Simulation of Simulation 3.The attitude angles of the principal-axes coordinate system, and the variation of the stability measurements calculated using above three methods with time are shown in Figure30a,b.

Figure 31 .
Figure 31.Sensor data of the motion with impact in Simulation 3: (a) posture angles of the principalaxes coordinate system; (b) positions; and (c) velocities of the body.

Figure 31 .
Figure 31.Sensor data of the motion with impact in Simulation 3: (a) posture angles of the principal-axes coordinate system; (b) positions; and (c) velocities of the body. ) Determination of the virtual-support plane.Set  11 =  12 , if  ⃗⃗⃗⃗⃗ ×   ⃗⃗⃗⃗⃗⃗⃗⃗ is equal to the zero vector, and  11 = − 11 , if the   component of  11 is negative.A similar method is used to deal with Equation (9).Hence,  21 is the unit normal of the plane on which , , and   lie and  22 is the unit normal of the plane on which , , and   lie. 2 is the unit normal vector of  2 : min{ 11  11 ,  12  12 } + min{ 21  21 ,  22  22 }.(

Table 1 .
Geometric and weight parameters of the robot.

Table 2 .
Contact parameters between the feet and the ground.

Table 4 .
Average values of stability evaluations.

Table 5 .
Impact forces applied in Test 1 and Test 2.

Table 4 .
Average values of stability evaluations.

Table 5 .
Impact forces applied in Test 1 and Test 2.

Table 4 .
Average values of stability evaluations.

Table 5 .
Impact forces applied in Test 1 and Test 2.

Table 6 .
Impact forces applied to the body.

Table 7 .
Average values of stability evaluations when the robot trotting on the 12 • slope and horizontal plane.