# Discretized Mid-Value CLVI-PDNN Based Redundancy Resolution for Single Leg of Quadruped Robot

Authors: Jinhui Rao; Taihui Zhang; Honglei An; Hongxu Ma
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

The two most important performance indicators of quadruped robot are load capacity and walking speed, and these performance indicators of the whole robot finally reflect on the joint torques and angular velocities. To satisfy different requirements of walking speed and load capacity when quadruped robots implement different tasks, the joint torques and angular velocities need to be balanced with physical constraints of the joints. A single leg with redundant DOF (degree of freedom) could optimize the distribution of joint torques or angular velocities based on different performance requirements. This paper presents a kind of new recurrent neural networks taking joint torques and angular velocities simultaneously into consideration and proposes mid-value CLVI-PDNN to achieve the optimal joint torques and angular velocities with physical constraints of the mechanism as described in our previous paper. Because the continuous mid-value CLVI-PDNN has difficulty in real-time operation because of too much calculation workload, two kinds of methods are proposed to discretize the mid-value CLVI-PDNN for application on computer or digital circuit. The simulation results demonstrate the efficacy of the algorithm proposed in this paper.


---

## Body

## 1. Introduction

Legged robots could be applied for rescue or service in unstructured environment due to its agility on rough terrain. Because of much more motion stability than biped robots and less complexity compared to hexapod robots or octopod robots, quadruped robots have a broad application prospect. Figure1 shows the hydraulic actuated quadruped robot developed by our laboratory, which has four legs assembled symmetrically, four revolute joints on each leg, a torso with computer, and other electronic equipment on board.Figure 1
Hydraulic actuated quadruped robot [31].The schematic diagram of the single leg is shown in Figure2, which has a redundant joint in the sagittal plane. The moment arms of actuating forces are marked in Figure 2 using red lines. If the performance of the hydraulic actuator is finitude, the maximum joint torques and angular velocity are determined by the length of these red lines. Dimension design of the moment arms of actuating forces directly affects the performance of walking speed and load capacity of the robot.Figure 2
Schematic diagram of the single leg [31].Designing the mechanism dimensions of the legs especially the moment arms needs to optimize redundancy resolution of the single leg depending on performance demands then distribute joint torques and angular velocity on the condition of joint physical constraints. Simultaneously, optimizing the distribution of joint torques and angular velocity could be used for designing the control algorithm of the robot. Generally, two kinds of methods are widely used to deal with inverse kinematics of redundant robots, i.e., analytical method and numerical analysis.In some analytical methods, inverse kinematics of redundant robot could be described as a parameterized function with parameters including joint angles [1, 2], the angle between the reference plane and robot [3, 4]. Manually adding constraints is another analytical method to calculate the inverse kinematics of redundant robot [5]. The inverse kinematics problem could be analytically solved while the redundant DOFs matches number of added constraints.Because redundant robot systems are generally strongly nonlinear and highly coupled, it is not a simple task to resolve inverse kinematics problem. Redundancy resolution is not only beneficial to obstacles avoidance and mechanism physical constraints avoidance [6–8] but also is used for improving motion performance based on diverse criteria [9–12]. Pseudoinverse type solution is a kind of widely used redundancy resolution [13, 14]. Some researchers solved inverse kinematics problem using different methods including interpolation [15, 16] and gradient projection [17–19], but these methods were restricted to computation speed and error. A number of novel or so called intelligent numerical resolving methods are used for redundant inverse kinematics problem, including quadratic programming (QP), artificial neural networks, quaternion method, online learning algorithm, and genetic algorithm [20–24]. Among these methods, quadratic programming is tractable and has well expansibility to real-time application. Minimizing sum of squares of joint velocity is the most widely utilized optimization criterion which is computationally tractable but ignores dynamics at all [25]. Besides, the optimization criterion that minimizes sum of squares of joint torques is also exploited to resolve redundant inverse kinematics which may lead to instability.To acquire a comprehensive performance combining kinematics and dynamics, Ma proposed an optimization criterion simultaneously minimizing joint velocity and torques [26]. Yunong Zhang also proposed an optimization criterion that consists of minimizing sum of squares of joint velocity and torques simultaneously [27]. The optimization problem on different level (i.e., joint velocity and joint torques) was converted into a standard QP problem with equality and inequality constraints on joint-acceleration level. Differing from pseudoinverse type solution used by Ma, Zhang numerically solved the constrained QP problem using a kind of recurrent neural network which is named linear variational inequality-based primal-dual neural network (LVI-PDNN) [28–30].It is different from manipulators that legged robots alter between support phase and flight (swing) phase when they are walking or running. Therefore, the legs need to swing forward and support body in different phases alternatively. Accordingly, the legs need to provide sufficient joint velocity and joint torques in different phases alternatively. On the basis of the above analysis, a new kind of optimization criterion consisting of support phase and flight (swing) phase is proposed in our previous paper and converted onto angular acceleration level based on ZNN. A feasible distribution of the joint torques and angular velocity is obtained according to continuous mid-value CLVI-PDNN, but the algorithm used in our previous paper has a drawback, i.e., hard to real-time application because of too large computation cost [31].For applying recurrent neural networks based real-time optimization on computer or digital circuit, discrete neural networks gradually attract attention from researchers. M.J. Perez-Ilzarbe proposed a kind of discrete recurrent neural network to solve bounded strictly convex quadratic programming problem [32]. Then he applied Wolf duality theory to form a kind of simple discrete recurrent neural network for solving linear constraints quadratic programming problem [33]. H. Tang et al. proposed a general recurrent neural network to solve hybrid constraints quadratic programming problem [34]. Q. Liu et al. used discrete recurrent neural network to solve equality constraints quadratic programming problem and proved the global exponential stability of its solution [35]. Y. Zhang, B. Liao, L. Jin et al. researched on a series of problems such as time-variant matrix inversion, time-variant polynomial root searching, and time-variant nonlinear optimization. They compared several discretization methods and proved the convergency of discrete recurrent neural network [36–43].The main contributions of this paper are as follows: firstly, two kinds of discretization methods are proposed to discretize the mid-value CLVI-PDNN for computer control of the quadruped robot, i.e., bilinear transform-type and Taylor-type discretization methods, while the recursion formulas of the two kinds of discrete mid-value CLVI-PDNN are also presented; secondly, inherent nature of the two kinds of discretization algorithms is analyzed theoretically, including the order of truncation error and stability; thirdly, simulation results of the proposed methods indicate the efficacy for redundancy resolution of single leg.


## 2. The Optimization Criterion and QP Formulation

In this paper, the redundant single leg of the quadruped robot is illustrated to show how to deal with its inverse kinematics problem using optimization method. Firstly, the kinematics and dynamics of the single leg system in sagittal plane with three joints are presented in this section. The kinematics of single leg is as follows:(1)r=fθ,where θ=θ1θ2θ3T is the position of joints in joint space, while r=r1r2T is the position of foot in Cartesian space. The differential kinematics is as follows:(2)r˙=Jθθ˙,where θ˙=θ˙1θ˙2θ˙3T is the velocity of joints in joint space, while r˙=r˙1r˙2T is the velocity of foot in Cartesian space. J=∂f/∂θ is Jacobian matrix. Second order differential equation is (3)r¨=Jθθ¨+J˙θθ˙,where θ¨=θ¨1θ¨2θ¨3T is the acceleration of joint in joint space, while r¨=r¨1r¨2T is the acceleration of foot in Cartesian space. J˙=∂J/∂t is differentiation ofJwith respect to timet. To simplify the description, Equation (3) is rewritten as follows:(4)Jθθ¨=r¨a,in which r¨a=r¨-J˙θ˙.As the general form of robot dynamics, the single leg also is presented as follows:(5)Mθθ¨+Cθ,θ˙+Gθ=τ,in which M(θ), C(θ,θ˙), and G(θ) denote inertia matrix, centripetal force and Coriolis force, and the gravitational term, respectively. τ=τ1τ2τ3T is joint torques of the single leg.It has been analyzed that a single optimization criterion is not suitable for legged robots that alter between support phase and flight (or swing) phase when walking or running. Then a modified optimization criterion that minimizes two norms including joint velocity and torques was proposed in our previous paper as follows:(6)min12αsθ˙TVθ˙+βsτTWτ,supportphase12αfθ˙TVθ˙+βfτTWτ,flightphase,with equality constraints Equations (4) and (5) and inequality constraints, i.e., physical constraints of the joint as follows:(7)θ-≤θ≤θ+θ˙-≤θ˙≤θ˙+θ¨-≤θ¨≤θ¨+,αs and βs are weight factors in support phase, while αf and βf are weight factors in flight (or swing) phase. The values of these weight factors fall in the interval 01 and satisfy αs+βs=1, αf+βf=1. V and W are diagonal weight matrices reflecting effects of different joints.The optimization problem formed by Equations (4)-(7) is converted into a standard QP problem on the joint angular acceleration level using the method base on ZNN in [27] as follows:(8)min12θ¨TQsθ¨+psTθ¨,supportphase12θ¨TQfθ¨+pfTθ¨,flightphase,subjected to equality and inequality constraints(9)r¨a=Jθ¨,θ¨-≤θ¨≤θ¨+,in which(10)Qs=αsV+1-αsMTWMpsT=αsλVθ˙+1-αsMTWTC+G,(11)Qf=αfV+1-αfMTWMpfT=αfλVθ˙+1-αfMTWTC+G,λ is a positive constant corresponding to convergence rate of ZNN.The standard QP problem formulated by Equations (8) and (9) is equivalent to a kind of neural network named LVI-PDNN (LVI-based primal-dual neural network) as follows:(12)y˙=ηI+HsTPΩy-Hsy+qs-y,supportphaseηI+HfTPΩy-Hfy+qf-y,flightphase,where y=[θ¨T,uT]T,qs=[ps,-r¨aT]T,qf=[pf,-r¨aT]T,Hs=Qs-JLTJL0,Hf=Qf-JLTJL0, and u denotes the dual decision variable vector consisting of coefficients of equality constraints in Equation (9). η is a positive constant that determines convergence rate of LVI-PDNN and should be selected large enough. PΩ(·) is a piecewise linear projection operator as follows:(13)PΩ·:Rn+m+l→Ω=y∣ς-≤y≤ς+PΩyi=ςi-,yi<ςi-yi,ςi-≤yi≤ςi+ςi+,yi>ςi+,i=1,⋯,n+m+l.


## 3. CLVI-PDNN and Mid-Value CLVI-PDNN

In our previous paper [31] two kinds of neural network are proposed based on LVI-PDNN, i.e., CLVI-PDNN and mid-value CLVI-PDNN. CLVI-PDNN improves the drawback of open-loop LVI-PDNN such as divergency of the optimized results from the desired value. Mid-value CLVI-PDNN further solves the problem that the optimized joint angles violate the physical constraints of the mechanism and deviate from the desired configuration especially for a legged robot.CLVI-PDNN means closed-loop LVI-PDNN that introduces the output to the input of LVI-PDNN as feedback shown in Figure3. In Figure 3, the block diagram in the blue dotted line shows LVI-PDNN proposed by Zhang [27]. The block diagram outside of the blue dotted line is a closed loop introducing feedback to the input of the neural network system as follows:(14)r¨a=r¨d+kprrd-fθ+kir∫0trd-fθdt+kpvr˙d-Jθθ˙+kiv∫0tr˙d-Jθθ˙dt-J˙θθ˙,in which rd, r˙d, and r¨d denote the desired motion of the foot in Cartesian space, respectively, while kpr, kir, kpv, and kiv are coefficient matrices.Figure 3
Block diagram of CLVI-PDNN [31].CLVI-PDNN has a critical defect that the optimized joint angles exceed the joint physical limitations leading to its uselessness for motion control of the quadruped robot. To improve this weakness of CLVI-PDNN, it could be modified as follows: defining a slack variableε=θ-θd, where θd denotes the desired joint configuration of the single leg when walking or running. ε could keep stable using a simple feedback control law(15)θ=θ0+kpθd-θ+∫0tθ˙dt,in which θ0 means initial values and kp means positive coefficient matrix. To prevent joint angles from exceeding the joint physical constraints, value of θd should be selected close to the mid-value of the joint motion range. That is also why this algorithm is named mid-value CLVI-PDNN. Figure 4 shows the block diagram of mid-value CLVI-PDNN.Figure 4
The block diagram of mid-value CLVI-PDNN [31].As analyzed in our previous paper [31], the mid-value CLVI-PDNN just compromise between the inequality and equality constrains. That is to say, the output of the QP solver maintains in motion range of the joints ensuring a feasible rather than optimal solution on condition that the error between the desired and the optimized trajectory is not enough to destroy normal motion control. To some extent, in many situations especially when the quadruped robots walk or run in the field, the output of mid-value CLVI-PDNN would be useful for motion control.For some circumstances, such as high dimensions of the augmented matrixMs and Mf, restricted performance of the onboard controller or computer, the discrete form of the QP solvers mid-value CLVI-PDNN may be more applicable for real-time control. We will deduce the discrete form of mid-value CLVI-PDNN and analyze the characteristics in next section.


## 4. Discretize Mid-Value CLVI-PDNN

To solve the problem that continuous mid-value CLVI-PDNN is difficult to be applied in real-time control system, two kinds of discretization methods are presented in this section. For a general unconstraint optimization problem(16)minxt∈Rnfxt,t.Suppose the extreme value of equation as follows:(17)Ω∗t=x∗t,t∣gx∗t,t=∂fx∗t,t∂x∗t=0.Obviously, condition of the extreme value is g(x(t),t)=0, so that the error equation is selected as e(t)=g(x(t),t). Then the error dynamics is designed on the basis of ZNN as follows:(18)e˙t=-λet,where λ is a positive constant to scale the rate of convergency of e(t). Substituting e(t)=g(x(t),t) into Equation (18) could obtain(19)Hxt,tx˙t=-λgxt,t-g′xt,t,where(20)Hxt,t=∂gxt,t∂xT=∂2fxt,t∂xt∂xTt,g′xt,t=∂gxt,t∂t=∂2fxt,t∂xt∂t.If the Hessian matrix H(x(t),t) is nonsingular, multiplying H-1(x(t),t) both sides of Equation (19) results in(21)x˙t=-H-1xt,tλgxt,t+g′xt,t.Given initial value x(0), the optimal state x∗(t) will be achieved on condition that Hessian matrix H(x(t),t) is positive definite. Actually, for our quadruped robot, Hessian matrix always maintains positive definite.Taking forward Euler method as example, the discretization process could be presented as follows: rewriting Equation (21) as x˙=G(x(t),t), the forward Euler equation is as follows:(22)xk+1=xk+δGxk,tk,where xk=x(tk) means the value of x(t) at moment tk. tk+1=tk+δ, in which δ denotes the sample interval. Therefore, discretization equation of x(t) is as follows:(23)xk+1=xk-H-1xk,tkδλgxk,tk+δg′xk,tk.Next the above-mentioned discretization process is applied for redundancy resolution of single leg. For a time-variant QP problem, it could be described at moment tk as follows:(24)min12θ¨kTQkθ¨k+pkTθ¨ks.t.JLkθ¨k=r¨akξ-≤θ¨k≤ξ+.Moreover, the dynamics of LVI-PDNN at moment tk is as follows:(25)y˙k=ηI+HkTPΩyk-Hkyk+qk-yk,where yk=θ¨kuk,qk=pk-r¨ak,Hk=Qk-JLkTJLk0.For different discretization methods, accuracy and rapidity differ from each other. Taylor-type and bilinear transform-type discretization methods are used for mid-value CLVI-PDNN in this paper. The discretized state of neural network using bilinear transform is as follows:(26)yk+1=yk+δ2y˙k+1+y˙k.Substituting Equation (25) into Equation (26) obtains(27)yk+1=yk+δ2ηI+HTPΩyk+1-Hk+1yk+1+qk+1-yk+1+PΩyk-Hkyk+qk-yk.Then, the recursion formula of yk could be achieved:(28)yk+1=2-δ2+δyk+δ2+δηI+HTPΩyk+1-Hk+1yk+1+qk+1+PΩyk-Hkyk+qk.yk+1 in the right-hand side of Equation (28) could be estimated according to forward difference(29)yk+1=yk+δy˙k=yk+δηI+HTPΩyk-Hkyk+qk-yk,Hk+1=Hk+δH˙k=Hk+ΔHk,qk+1=qk+δq˙k=qk+Δqk,where ΔHk=Hk-Hk-1 and Δq=qk-qk-1.Taylor-type discretization method exploits high order terms of the Taylor series and then utilizes multiple points to discretize mid-value CLVI-PDNN. This kind of discretization method possesses the advantage of higher computational accuracy. Higher order corresponds higher accuracy while consuming more computation. Eliminating third-order derivative terms of Taylor series obtains(30)y˙k=2yk+1-3yk+2yk-1-yk-22δ+Oδ2.Substituting Equation (25) into Equation (30) achieves the recursion formula as follows:(31)yk+1=1.5yk-yk-1+0.5yk-2+δy˙k=1.5yk-yk-1+0.5yk-2+δηI+HkTPΩyk-Hkyk+qk-yk.


## 5. Theory Analysis

In this section, convergency and steady-state residual error of the discretization methods used in this paper is analyzed.


### 5.1. Bilinear Transform-Type Discretization Method

Theorem 1.
Bilinear transform-type recursion formula as Equation (28) is 0-stable.Proof.
The characteristic polynomial of bilinear transform-type recursion formula as Equation (28) is(32)Pz=z-2-δ2+δ=0,which has a unique root on/in the unit disk, i.e., z=(2-δ)/(2+δ) with 0<δ<1. Therefore, bilinear transform-type recursion formula as Equation (28) is 0-stable. The proof is thus completed.Theorem 2.
Bilinear transform-type recursion formula as Equation (28) is consistent and convergent, which converges with the order of truncation error O(δ2).Proof.
Discretization formula ofyk based on forward difference formula is as follows:(33)y˙k=yk+1-ykδ+Oδ.While discretization formula of yk based on backward difference formula is as follows:(34)y˙k+1=yk+1-ykδ+Oδ.Adding Equations (33) and (34) could obtain the bilinear transform-type discretization formula(35)yk+1=yk+δ2y˙k+1+y˙k+Oδ2.Note that dropping O(δ2) from Equation (35) yields exactly bilinear transform-type recursion formula as Equation (26), and, thus, the truncation error of Equation (26) is O(δ2). Consequently, bilinear transform-type recursion formula as Equation (28) is consistent with order 2. Combining withTheorem1, we know that the model in Equation (28) is both 0-stable and consistent, so that it can be conclusively concluded that bilinear transform-type recursion formula as Equation (28) is consistent and convergent, which converges with the order of its truncation error being O(δ2). The proof is now completed.Theorem 3.
Considering the QP solver mid-value CLVI-PDNN based on ZNN, the steady-state residual errorlimk→∞⁡y˙k2 of bilinear transform-type discretization recursion formula as Equation (28) is O(δ2) in the sense of Euclidean norm.Proof.
Letyk∗ denote the theoretical optimal solution of the QP problem described by Equations (8) and (9), we know that PΩ[yk∗-(Hkyk∗+qk)]-yk∗=0. On the basis ofTheorems1 and 2, it could be concluded that yk=yk∗+O(δ2) with k large enough. Therefore, we obtain the limit error (the definition of error follows [42]):(36)limk→∞y˙k2=limk→∞ηI+HkTPΩyk-Hkyk+qk-yk2=limk→∞ηI+HkTPΩyk∗+Oδ2-Hkyk∗+Oδ2+qk-yk∗-Oδ22=limk→∞ηI+HkTPΩyk∗-Hkyk∗+qk-yk∗+ηI+HkTI-Hk-IOδ22=limk→∞ηI+HkTI-Hk-IOδ22.Considering that ((I+HkT)(I-Hk)-I) is a constant matrix at a certain time instant, we further obtain(37)limk→∞y˙k2=limk→∞ηI+HkTI-Hk-IOδ22=Oδ2.The proof is completed.


### 5.2. Taylor-Type Discretization Method

Taylor-type recursion formula as Equation (30) is 0-stable. Taylor-type recursion formula as Equation (30) is consistent and convergent, which converges with the order of truncation error O(δ3). These conclusions are proved in [42]. Therefore, the steady-state residual error could be concluded as follows.Theorem 4.
Considering the QP solver mid-value CLVI-PDNN based on ZNN, the steady-state residual errorlimk→∞⁡y˙k2 of Taylor-type discretization recursion formula as Equation (31) is O(δ3) in the sense of Euclidean norm.Proof.
Letyk∗ denote the theoretical optimal solution of the QP problem described by Equations (8) and (9), so that PΩ[yk∗-(Hkyk∗+qk)]-yk∗=0 and yk=yk∗+O(δ3) withklarge enough. Therefore, we obtain the limit error as follows:(38)limk→∞y˙k2=limk→∞ηI+HkTPΩyk-Hkyk+qk-yk2=limk→∞ηI+HkTPΩyk∗+Oδ3-Hkyk∗+Oδ3+qk-yk∗-Oδ32=limk→∞ηI+HkTPΩyk∗-Hkyk∗+qk-yk∗+ηI+HkTI-Hk-IOδ32=limk→∞ηI+HkTI-Hk-IOδ32=Oδ3.The proof is completed.


## 5.1. Bilinear Transform-Type Discretization Method

Theorem 1.
Bilinear transform-type recursion formula as Equation (28) is 0-stable.Proof.
The characteristic polynomial of bilinear transform-type recursion formula as Equation (28) is(32)Pz=z-2-δ2+δ=0,which has a unique root on/in the unit disk, i.e., z=(2-δ)/(2+δ) with 0<δ<1. Therefore, bilinear transform-type recursion formula as Equation (28) is 0-stable. The proof is thus completed.Theorem 2.
Bilinear transform-type recursion formula as Equation (28) is consistent and convergent, which converges with the order of truncation error O(δ2).Proof.
Discretization formula ofyk based on forward difference formula is as follows:(33)y˙k=yk+1-ykδ+Oδ.While discretization formula of yk based on backward difference formula is as follows:(34)y˙k+1=yk+1-ykδ+Oδ.Adding Equations (33) and (34) could obtain the bilinear transform-type discretization formula(35)yk+1=yk+δ2y˙k+1+y˙k+Oδ2.Note that dropping O(δ2) from Equation (35) yields exactly bilinear transform-type recursion formula as Equation (26), and, thus, the truncation error of Equation (26) is O(δ2). Consequently, bilinear transform-type recursion formula as Equation (28) is consistent with order 2. Combining withTheorem1, we know that the model in Equation (28) is both 0-stable and consistent, so that it can be conclusively concluded that bilinear transform-type recursion formula as Equation (28) is consistent and convergent, which converges with the order of its truncation error being O(δ2). The proof is now completed.Theorem 3.
Considering the QP solver mid-value CLVI-PDNN based on ZNN, the steady-state residual errorlimk→∞⁡y˙k2 of bilinear transform-type discretization recursion formula as Equation (28) is O(δ2) in the sense of Euclidean norm.Proof.
Letyk∗ denote the theoretical optimal solution of the QP problem described by Equations (8) and (9), we know that PΩ[yk∗-(Hkyk∗+qk)]-yk∗=0. On the basis ofTheorems1 and 2, it could be concluded that yk=yk∗+O(δ2) with k large enough. Therefore, we obtain the limit error (the definition of error follows [42]):(36)limk→∞y˙k2=limk→∞ηI+HkTPΩyk-Hkyk+qk-yk2=limk→∞ηI+HkTPΩyk∗+Oδ2-Hkyk∗+Oδ2+qk-yk∗-Oδ22=limk→∞ηI+HkTPΩyk∗-Hkyk∗+qk-yk∗+ηI+HkTI-Hk-IOδ22=limk→∞ηI+HkTI-Hk-IOδ22.Considering that ((I+HkT)(I-Hk)-I) is a constant matrix at a certain time instant, we further obtain(37)limk→∞y˙k2=limk→∞ηI+HkTI-Hk-IOδ22=Oδ2.The proof is completed.


## 5.2. Taylor-Type Discretization Method

Taylor-type recursion formula as Equation (30) is 0-stable. Taylor-type recursion formula as Equation (30) is consistent and convergent, which converges with the order of truncation error O(δ3). These conclusions are proved in [42]. Therefore, the steady-state residual error could be concluded as follows.Theorem 4.
Considering the QP solver mid-value CLVI-PDNN based on ZNN, the steady-state residual errorlimk→∞⁡y˙k2 of Taylor-type discretization recursion formula as Equation (31) is O(δ3) in the sense of Euclidean norm.Proof.
Letyk∗ denote the theoretical optimal solution of the QP problem described by Equations (8) and (9), so that PΩ[yk∗-(Hkyk∗+qk)]-yk∗=0 and yk=yk∗+O(δ3) withklarge enough. Therefore, we obtain the limit error as follows:(38)limk→∞y˙k2=limk→∞ηI+HkTPΩyk-Hkyk+qk-yk2=limk→∞ηI+HkTPΩyk∗+Oδ3-Hkyk∗+Oδ3+qk-yk∗-Oδ32=limk→∞ηI+HkTPΩyk∗-Hkyk∗+qk-yk∗+ηI+HkTI-Hk-IOδ32=limk→∞ηI+HkTI-Hk-IOδ32=Oδ3.The proof is completed.


## 6. Results and Discussion

The proposed discrete QP algorithms are operated in Matlab to solve the optimization problem formulated by Equations (8) and (9). In the simulation environment, the quadruped robot moves along horizontal direction with trotting gait in sagittal plane. The body of the quadruped robot moves at a constant speed on horizontal direction as well as constant height. The parameters of trotting gait are listed in Table 1.h is a constant meaning the height of the hip. hapex is the maximum height of the foot when the quadruped robot is walking. v is a constant meaning the desired velocity of the body. T is the gait cycle period of trotting gait. g is the gravitational acceleration. Parameters of the single leg system with half mass of the whole body in the sagittal plane as shown in Figure 2 are listed in Table 2. M is mass of body of the single system, i.e., half mass of the whole body. The others mean the masses and lengths of the three links of the single leg in sagittal plane.Table 1
Parameters of trotting gait.
parameter value unit h 0.7 m h a p e x 0.15 m v 6.0 Km/h T 0.6 sTable 2
Physical parameters of single leg.
parameter value unit M 50 kg M 1 3 kg M 2 3 kg M 3 3 kg L 1 0.25 m L 2 0.30 m L 3 0.30 mThe joint physical constraints, i.e., the upper and lower bound of joint angles, angular velocity and acceleration are listed in Table3. The parameters of the discretized mid-value CLVI-PDNN are listed in Table 4.Table 3
Physical constraints of joints.
parameter value unit θ + 17π/36,5π/6,5π/6 rad θ - π/9,π/9,π/9 rad θ ˙ + 12,12,20 rad/s θ ˙ - -12,12,20 rad/s θ ¨ + 1700 ,  2200 ,  2800 rad/s θ ¨ - - 1700 ,  2200 ,  2800 rad/sTable 4
Parameters of mid-value CLVI-PDNN.
parameter value λ 3600 η 107 k p r 1000 k i r 100000 k p v 100000 k i v 100000In addition, motion trajectory of foot is designed as a curve respect to timet on the horizontal and vertical direction respectively. The expected motion curves in horizontal direction and vertical direction are described as follows:(39)x=-π2sin20π3t+53t-14,0≤t≤0.314-53t-0.3,0.3<t≤0.6,(40)y=-0.7+237t3-2370.4t4+7901.2t5-8779.1t6,0≤t≤0.3-0.7,0.3<t≤0.6.The step length is set the same as the control cycle, i.e., δ=1ms. The simulation results of the single leg under the aforementioned conditions are illustrated in Figures 5–10, in which Figures 5–7 present simulation results acquired by using Taylor-type method to discretize the mid-value CLVI-PDNN, while Figures 8–10 present simulation results acquired with bilinear transform-type discretization method.Figure 5
Taylor-type discrete mid-value CLVI-PDNN with parameter set 1.Figure 6
Taylor-type discrete mid-value CLVI-PDNN with parameter set 2.Figure 7
Taylor-type discrete mid-value CLVI-PDNN with parameter set 3.Figure 8
Bilinear transform-type discrete mid-value CLVI-PDNN with parameter set 1.Figure 9
Bilinear transform-type discrete mid-value CLVI-PDNN with parameter set 2.Figure 10
Bilinear transform-type discrete mid-value CLVI-PDNN with parameter set 3.The parameter set 1 in Figure5 are V=W=diag([1,1,1]) (diag(X) means a diagonal matrix whose diagonal elements are vector X), αs=0.9,αf=0.1,θd=(0.6379,1.8059,1.0472); the parameter set 2 in Figure 6 are V=W=diag([1,1,1]),αs=0.1,αf=0.9, θd=(0.6379,1.8059,1.0472); the parameter set 3 in Figure 7 are V=diag([4,2,1]),W=diag([1,2,4]), αs=0.1, αf=0.9, θd=0.63791.80591.0472. It is concluded that the optimized trajectory of the single leg in operation space could track the desired value well when using Taylor-type discretization method with position error about ±0.5cm in horizontal direction and ±2.5cm in vertical direction. For quadruped robot walking in the field that requires less position precision, the position error shown in Figures 5–7 could satisfy the design and control requirements relative to the motion range of the single leg. On the other hand, different parameter sets result in different optimized motion curves in joint space. As illustrated in Figures 5 and 6, the joint angular velocities in Figure 6 with αs=0.1 are smoother than the curves in Figure 5 with αs=0.9 demonstrating that norm of minimizing joint torques would lead to instability. The angular velocity of hip joint with V=diag([1,1,1]) shown in Figure 6 is larger than the value with V=diag([4,2,1]) shown in Figure 7. While angular velocity of ankle joint shown in Figure 6 is smaller than the value shown in Figure 7. Moreover, the motion ranges of joint angles in Figures 6 and 7 are quite distinct correspondingly. The comparison between Figures 6 and 7 indicates the effect of different velocity weights among the three joints.The parameter set 1 in Figure8 are V=W=diag([1,1,1]), αs=0.9,αf=0.1,θd=(0.6379,1.8059,1.0472); the parameters set 2 in Figure 9 are V=W=diag([1,1,1]),αs=0.1, αf=0.9, θd=(0.6379,1.8059,1.0472); the parameters set 3 in Figure 10 are V=diag([4,2,1]),W=diag([1,2,4]),αs=0.1, αf=0.9,θd=0.63791.80591.0472. It is concluded that the optimized trajectory of the single leg in operation space could track the desired value well using bilinear transform-type discretization method with position error about ±0.5cm in horizontal direction and ±2.5cm in vertical direction. But it is disappointing that the motion trajectory in joint space is unsmooth and fluctuant so that perhaps it is difficult to be utilized in control of the single leg on the real system.According to the simulation results of the two kinds of discretization methods, i.e., Taylor-type and bilinear transform-type, a conclusion could be achieved easily that the optimized motion trajectory in operation space tracks the desired value well and satisfies the control requirements. But the control of the single leg executes in joint space actually, the optimized motion trajectory in joint space by using Taylor-type discretization method is better than the other one obviously and easier to be applied in the control of the real mechanism.


## 7. Conclusion

As the discussion in our previous paper [31], we try to find a practicable QP solver to solve the inverse kinematics of the redundant single leg of the quadruped robot. In this paper we adopt two kinds of discretization methods to discretize the mid-value CLVI-PDNN proposed in our previous paper. According to the simulation results and analysis we could conclude that the optimized motion trajectories in operation space obtained by both Taylor-type and bilinear transform-type discretization methods track the desired values well if we do not overly concern about the position control precision of the single leg in the wild. Meanwhile, it is obvious that the optimized motion trajectory in joint space obtained by Taylor-type discretization method is superior to the other one. Of course, higher order Taylor-type discretization methods will achieve more accurate inverse kinematics of the redundant single leg.It was analyzed in our previous paper that the mid-value CLVI-PDNN QP solver just outputs feasible solutions rather than intuitional optimal solutions corresponding to the optimization criterion formulated in Equation (6). In addition, the deviation between the optimized value and the desired value is much larger than theory analysis. The error may be caused by the constraints such as integration threshold, the deviation between the criteria in Equations (6) and (8), and probably the crucial reason, i.e., the introduction of expected joint configuration. Although the optimized value deviates from the desired value to some extent, it is satisfactory to achieve a useable answer in the end.



---

*Source: 5071254-2019-06-24.xml*