A hierarchical optimization control method for dynamic walking of wheeled biped robot

Abstract

The application discloses a layered optimization control method for dynamic marching of a wheeled biped robot, and particularly relates to the following steps: obtaining optimal reference forward position, optimal reference forward velocity, optimal reference forward acceleration and optimal future trajectory gamma by model prediction from target velocity and current velocity, target position and current position * ; obtaining optimal reference center of mass inclination theta by fitting an inner dynamic inverse mapping function through a nonlinear least square method, using gamma * and the inner dynamic of a first-order cart-pole model * ref ; establishing a third-order cart-pole model, solving optimal reference joint angles based on and constructing a QP optimization problem related to and by using multi-joint dynamics, and solving the QP optimization problem to obtain optimal control force and torque, which are used for controlling the biped robot. The application has small calculation requirement, so that the robot can maintain self balance and track a reference trajectory.

Description

Technical Field

[0001] This invention belongs to the field of bipedal robot technology, specifically relating to a hierarchical optimization control method for the dynamic movement of wheeled bipedal robots. Background Technology

[0002] Wheeled robots and bipedal robots are the two main forms of terrestrial mobile robots today, and these two types of mobile robots have significant complementary advantages. Bipedal robots can handle more complex environments than traditional wheeled robots. For example, bipedal robots can traverse obstacles and unstructured surfaces, and with bionic arms, they can perform manipulative tasks such as grasping, carrying, and driving. However, unlike wheeled robots, high-speed, long-distance movement is more difficult for bipedal robots, and walking consumes significantly more energy than rolling. Therefore, to combine the complementary advantages of both, wheeled bipedal robots have received increasing attention.

[0003] Due to the significant advantages of wheeled bipedal robots, researchers from various countries have proposed numerous related robots. Boston Dynamics and Harbin Institute of Technology have each proposed two full-size hydraulically driven tandem wheeled bipedal robots. Tencent and ETHZ have each designed a small-sized motor-driven parallel wheeled bipedal robot. While wheeled bipedal robots combine many advantages of wheels and bipedalism, the combination of wheels and feet also brings some significant problems: First, without a foot plate system, it is difficult to maintain balance because the contact between the wheeled bipedal robot and the ground degenerates from surface contact to point contact, causing the dynamic system to change from a locally statically stable system to a critically stable system; second, because the robot's longitudinal motion is coupled with its pitch motion, the robot needs a certain tilt angle to accelerate. As the tilt angle of the center of gravity increases, the nonlinear characteristics of the system will be significantly enhanced, thus disrupting the balance. This contradiction makes it difficult for the robot to accurately track the longitudinal motion reference value while maintaining balance.

[0004] To address the aforementioned issues, the traditional approach is to simplify the robot as a first-order inverted pendulum system and constrain its center of mass to a vertical position during movement, which prevents the robot from responding promptly to changes in the target speed. Summary of the Invention

[0005] To address the shortcomings of existing technologies, this invention provides a hierarchical optimization control method for the dynamic movement of wheeled bipedal robots.

[0006] The present invention achieves the above-mentioned technical objectives through the following technical means.

[0007] A hierarchical optimization control method for dynamic movement of wheeled bipedal robots:

[0008] From the target velocity vtar and current velocity v, target position x ar Given the current position x, the optimal reference forward position is obtained through model prediction. Optimal reference forward velocity Optimal reference forward acceleration and the optimal future trajectory γ * ;

[0009] Using the optimal future trajectory γ * The internal dynamics of the first-order inverted pendulum model of the cart were obtained by fitting the inverse mapping function of the internal dynamics using the nonlinear least squares method, thus obtaining the optimal reference acceleration. The optimal reference centroid tilt angle θ * ref ;

[0010] Establish a third-order inverted pendulum model of a cart, based on and Solving for the optimal reference joint angle Constructing with multi-joint dynamics and The relevant QP optimization problem is solved to obtain the optimal control force and torque, which are then applied to the bipedal robot for control.

[0011] A further technical solution is that the constraint equations predicted by the model are:

[0012]

[0013]

[0014] v ref (k+i)≤|v max |

[0015] |a ref (k+i)|≤|Ψ max |

[0016]

[0017]

[0018] Where: X ref (k) represents the generalized state variable of the k-th prediction, x ref (k) represents the reference forward position of the k-th prediction, v ref (k) represents the reference forward velocity predicted in the k-th prediction, x(k), v(k), θ(k), a ref (k) represents the robot's actual forward position, actual forward velocity, actual center of mass tilt angle, and actual forward acceleration observed during a single prediction, respectively. initThis represents the initial value of Ψ in each prediction. max This represents the maximum value limit of Ψ, where Ψ represents the robot's internal dynamic positive mapping function F. PD The output value, v ref (k+i), a ref (k+i) represent the reference forward velocity and acceleration of the robot from the k-th prediction onwards to the i-th prediction, respectively. max ΔT represents the rated speed of the robot wheel motor. mpc This indicates the prediction step size.

[0019] A further technical solution involves the model prediction using the following method:

[0020]

[0021] Among them: J N (k) represents the cost function, Υ(k) is an intermediate quantity, and U ref (k) represents the generalized input state variable of the k-th prediction, and N is the total number of prediction steps per step; matrix H α =2(Φ+Γ T ΩΓ),F α =2Γ T ΩΘ, Ω、Ф、Г、Θ、 Λ represents a matrix; Q is the error weight matrix in the prediction process.

[0022] A further technical solution is to obtain the optimal reference centroid tilt angle in the following manner:

[0023] Take the first term of the optimal solution for the current prediction period. Combination and Solve the problem. The optimal fitting matrix is... For the fitting matrix, Ψ -1 This is the output value of the inner dynamic inverse mapping function. For Ψ -1 An approximation, where g is the acceleration due to gravity.

[0024] A further technical solution is that the internal dynamic inverse mapping function is specifically obtained by constructing a cost function for the error between the internal dynamic inverse mapping function and its approximate inverse mapping function based on the optimal forward acceleration and the internal dynamic forward mapping function of the previous prediction period, and solving it using the nonlinear least squares method.

[0025] A further technical solution is that the cost function is:

[0026]

[0027] Where: k represents the k-th model prediction, i represents the i-th subsequent prediction, J ID (k) is the cost function.

[0028] Further technical solutions, and and The relevant QP optimization problem includes an objective function and constraints, wherein the objective function is:

[0029]

[0030] Among them: A qp It is the weight coefficient matrix, b qp It is a feedforward matrix, and:

[0031]

[0032] The reference adjustment acceleration is the feedback control output between the generalized state variable and its optimal reference value in QP optimization. The joint torque τ at the reference position leg feedforward reference value, and K represents the reference values ​​for the ground horizontal force F1 and ground support force F2 converted from wheel joint torque. p K d The proportional and differential coefficients of this feedback, χ ref* (k) is the optimal reference value of the generalized state variable. H represents the optimal reference position in the z-direction, χ is the generalized state variable of the model, and H -1 It is the pseudo-inverse matrix of the joint driving matrix H, where M, C, and G represent the inertia matrix, Coriolis force matrix, and gravity matrix of the third-order inverted pendulum model of the trolley, respectively, and J leg F is the Jacobian matrix of the leg end relative to the world coordinate system. ext K represents the external force acting on the robot's chassis. pv ,K dv ,K pl ,K dl M represents the proportional and derivative adjustment coefficients during the tracking process. w q1 represents the joint angle of the robot's thigh, q2 represents the joint angle of the robot's torso, g is the acceleration due to gravity, and x and v are the robot's current velocity and current position, respectively. This represents the optimal reference joint angle for the thigh joint. This represents the optimal reference joint angle for the lower leg joint.

[0033] The constraints include equality constraints and inequality constraints, wherein the equality constraints are: The inequality constraints are:

[0034] Where: τ max F corresponds to the joint torque limit values ​​of the lower leg and thigh. max The limit value of the external input force is given by μ, which is the coefficient of friction between the tire and the ground.

[0035] The beneficial effects of this invention are as follows:

[0036] (1) This invention uses different simplified models and optimization methods to optimize the planning and control at the kinematic level and the control at the dynamic level separately, which requires less computing power.

[0037] (2) The present invention sets constraint equations in the kinematic planning and control process based on model prediction. The obtained reference forward position, velocity and acceleration can track the target command, meet the requirements of robot dynamics, and ensure the stability of robot motion.

[0038] (3) Based on the internal dynamics and least squares fitting method of the first-order inverted pendulum of the trolley, the present invention obtains the reference centroid tilt angle trajectory that matches the reference forward acceleration, and solves the problem of large speed tracking error caused by the mutual interference between the reference forward acceleration and the reference centroid tilt angle trajectory.

[0039] (4) Based on the design of the multi-joint dynamics optimization control framework and constraint equations, this invention obtains the optimal input force and torque of the robot, enabling the robot to maintain its own balance and track the reference trajectory, thus solving the problem of difficult multi-joint cooperative control of wheeled bipedal robots. Attached Figure Description

[0040] Figure 1 This is a block diagram of the hierarchical optimization control for dynamic movement of a wheeled bipedal robot as described in this invention;

[0041] Figure 2 This is a schematic diagram of the two-dimensional dynamic model described in this invention;

[0042] Figure 3 This is a schematic diagram of the three-order inverted pendulum model of the trolley described in this invention. Detailed Implementation

[0043] The present invention will be further described below with reference to the accompanying drawings and specific embodiments, but the scope of protection of the present invention is not limited thereto.

[0044] This invention discloses a hierarchical optimization control method for the dynamic movement of a wheeled bipedal robot. Specifically, it is implemented in three parts: kinematic planning and control based on model prediction, trajectory planning based on the inverse internal dynamic solution of the center of mass tilt angle, and motion control based on multi-joint dynamics optimization. This method does not consider the robot's turning motion, but only performs analysis in the robot's sagittal plane. The kinematic planning and control based on model prediction and the trajectory planning based on the inverse internal dynamic solution of the center of mass tilt angle belong to the upper-level kinematic planning part, while the motion control based on multi-joint dynamics optimization belongs to the lower-level dynamics control part. Figure 1 As shown, the hierarchical optimization control method is as follows: First, the industrial control computer determines the target speed v. tar and current velocity v, target position x ar Given the current position x, the optimal reference forward position is obtained through model prediction. Optimal reference forward velocity Optimal reference forward acceleration and the optimal future trajectory γ * Then, using the optimal future trajectory γ * The internal dynamics of the first-order inverted pendulum model of the cart were obtained by fitting the inverse mapping function of the internal dynamics using the nonlinear least squares method, thus obtaining the optimal reference acceleration. The optimal reference centroid tilt trajectory θ * ref Finally, a third-order inverted pendulum model of the robot is established. Multi-joint dynamics is used to construct a solution framework that satisfies the QP optimization problem. The multi-task optimization problem that can maintain the robot's posture and track the reference trajectory is solved to obtain the optimal control force and torque, thereby controlling the robot to quickly and stably track the reference trajectory.

[0045] The following is a detailed introduction to each part:

[0046] (1) Kinematic planning and control based on model prediction

[0047] The main objective of kinematic planning and control for wheeled bipedal robots is to generate the optimal reference forward position that can track the target velocity. Optimal reference forward velocity and optimal forward acceleration The model prediction is performed once in each control cycle of the robot's industrial computer, and each prediction can predict the next N prediction steps ΔT. mpc The speed following error is analyzed and optimized to obtain the optimal predicted trajectory that satisfies the constraints and minimizes the cost. To reduce the computational difficulty of optimization, the prediction model used in this invention is a pure kinematic model, but internal dynamic constraints describing the coupling relationship between the tilt angle and forward acceleration of the wheeled bipedal robot are added at the kinematic constraint level. The total number of prediction steps for each step in this prediction method is N=10, and the prediction step length ΔT for each step is...mpc Control cycle ΔT of the robot industrial computer ctrl The same, i.e., ΔT ctrl =ΔT mpc To reduce computational power requirements, the model prediction frequency is four times that of the robot control system, with model prediction performed only once every four control cycles.

[0048] The kinematic prediction iterative equations used are as follows:

[0049] X ref (k+1)=AX ref (k)+BU ref (k)

[0050]

[0051] X ref (k+i)=AX ref (k+i-1)+BU ref (k+i-1)

[0052]

[0053] X ref (k+N)=AX ref (k+N-1)+BU ref (k+N-1),

[0054]

[0055] X ref (k+i)=[x ref (k+i)v ref (k+i)] T U ref (k+i)=a ref (k+i);

[0056] Where: k represents the k-th model prediction, a positive integer; i represents the i-th prediction step forward, a non-negative integer; N is the total number of prediction steps per step, a positive integer; x ref (k+i),v ref (k+i),a ref (k+i) represent the reference forward position, velocity, and acceleration of the robot from the k-th prediction onwards to the i-th prediction; X ref (k+i) represents the generalized state variable from the k-th prediction onwards, up to the i-th prediction step, U ref(k+i) represents the generalized input state variable from the k-th prediction to the i-th prediction step; A is the state matrix of the discrete system of the prediction model, and B is the input matrix of the discrete system of the prediction model. This model is a discrete linear time-invariant equation, applicable to the standard form of linear model prediction.

[0057] The cost function J for model prediction is defined as follows:

[0058]

[0059] Among them: J N X is the total cost function in the prediction time domain; P is the error weight matrix at the end of the prediction time domain, Q is the error weight matrix during the prediction process, and R is the control weight matrix. All three matrices are symmetric positive definite matrices; tar =[x tar v tar ] T This refers to the input commands from the robot's external control terminal. In this invention, it is assumed that the external input command X is used during a single prediction process. tar constant.

[0060] The constraint equations for model prediction are defined as follows:

[0061]

[0062]

[0063]

[0064]

[0065]

[0066] The internal dynamics primarily describe the relationship between the robot's forward acceleration and its center of mass tilt angle. To maintain balance, the reference center of mass tilt angle cannot experience significant abrupt changes to ensure that the robot's upright loop control input does not change abruptly. Therefore, three saturation constraints are applied to the reference forward acceleration in the constraint equations. During a single prediction process, the reference forward acceleration at each step is required to not exceed Ψ. max The first-order difference of the reference forward acceleration at each step is required to not exceed [a certain value]. The second-order difference of the reference forward acceleration at each step is required to not exceed At the start of each prediction, a is required ref (k) and the previous time a ref Compared to (k-1), no mutation occurs, and its value is required to be equal to Ψ. init For the state variable x ref (k), v ref(k) only has its initial value and maximum value constrained.

[0067] In the constraint equations: x(k), v(k), and θ(k) represent the robot's actual forward position, actual forward velocity, and actual center of mass tilt angle observed during a single prediction; v max Ψ represents the rated speed of the robot's wheel motors; F represents the robot's internal dynamic positive mapping function. PD The output value, whose specific form is shown in the next section, Ψ init This describes the initial value of Ψ for each prediction. max This describes the maximum value limit for Ψ.

[0068] Since both the iterative equations and the constraint equations are linear equations, they can be simplified and written in the following form:

[0069] X(k)=ΘX ref (k)+ΓΥ(k)

[0070]

[0071] Substitute the cost function J N This yields the following form:

[0072] J N (k)=X ref T (k)QX ref (k)+X T (k)ΩX T (k)+Υ T (k)ΦΥ T (k),

[0073]

[0074] Furthermore, since matrices R, P, and Q are symmetric positive definite matrices, matrices Ω and Φ are also positive definite matrices. Therefore, the cost function J... N (k) can be further written in the following form:

[0075]

[0076] H α =2(Φ+Γ T ΩΓ),F α =2Γ T ΩΘ

[0077] Since matrices Ω and Φ are positive definite, therefore H α It is a positive definite matrix; since X(k) is the known observation of the robot at the initial time in the prediction time domain, J N (k) is a quadratic convex function of Υ(k).

[0078] Similarly, since the inequality constraints in the constraint equations are linear constraints, they can also be transformed into:

[0079]

[0080]

[0081]

[0082]

[0083] Since X(k) = ΘX ref Substituting (k)+ΓΥ(k) into the above inequality constraint and further combining like terms, the above inequality constraint can be further simplified to about X. ref The linear inequalities of (k) and Υ(k) are:

[0084]

[0085] Therefore, this prediction process can be transformed into a standard QP optimization form, namely:

[0086]

[0087] The cost function J is obtained through QP optimization. N (k) takes a local minimum value, that is, the optimal solution Υ(k) * (k):

[0088]

[0089] It describes the sequence of optimal control variables obtained from each prediction; according to the rolling optimization principle, only the first term is retained. As the control variable at the current prediction time, the control variable at the next time is generated as the prediction time domain rolls. ref * The physical meaning of (k) is the reference optimal forward acceleration at the current prediction time. The optimal reference forward velocity at the current moment can be obtained by integration. Optimal reference forward position at the current moment These three optimal references are important reference quantities in the lower-level dynamic control; at the same time, this method also utilizes the optimal solution Y. * (k) Solve for the internal dynamic inverse mapping function at the next time step.

[0090] (2) Centroid tilt angle trajectory planning based on internal dynamic inverse solution

[0091] In trajectory planning, to reduce computational burden, the model is usually simplified. The error between the simplified model and the actual model is typically treated as a disturbance term in the lower-level dynamic control section. A wheeled bipedal robot has 14 degrees of freedom and is a highly complex dynamic system. Therefore, in the trajectory planning for the center of mass tilt angle, the robot is simplified to a first-order inverted pendulum model. Since the motion involved in this invention is mainly the robot's forward motion, the model is simplified to a two-dimensional dynamic model built on the sagittal plane. Specifically, the robot chassis is simplified to a cart, and the robot torso is simplified to a point mass. For example... Figure 2 As shown, where M w It is the mass of the car, m b I is the mass of the center of mass, l is the distance to the hinge point of the trolley, and θ is the angle of inclination of the line connecting the center of mass and the hinge point of the trolley relative to the vertical direction. b It is the rotational inertia of the robot's upper body relative to the wheel axle, and F is the external driving force on the robot, which is the external input of the system.

[0092] The dynamic equations obtained by modeling the two-dimensional dynamic model using the Newton-Euler method are shown below:

[0093]

[0094] By further reorganizing, combining like terms, and eliminating F, we can obtain the following subsystem:

[0095]

[0096] This subsystem, lacking external input and describing only the coupling relationships between the system's state variables, is called the internal dynamic subsystem. Specifically, it describes the relationship between acceleration and the tilt angle of the body's center of mass during acceleration in a wheeled bipedal robot. The corresponding equation for this subsystem is a multivariable coupled equation with differential components. The forward mapping function of this equation is relatively easy to solve. It describes the mapping relationship from the tilt angle of the center of mass to the forward acceleration; however, the inverse mapping function of this equation... Explicit analytical solutions are difficult to obtain, making it challenging to find an inverse mapping function describing the acceleration to the center-of-mass tilt angle. Therefore, to obtain a center-of-mass tilt angle trajectory that matches the planned optimal acceleration trajectory, a nonlinear least squares method is used to fit the inverse mapping function of the internal dynamics.

[0097] First, at the initial time of each prediction time domain, three fitting functions are selected to fit the internal dynamic inverse mapping function, in the following form:

[0098]

[0099] Among them: Ψ -1 For the internal dynamic inverse mapping function FID The output value; For Ψ -1 approximation; To fit the matrix, a 3x3 diagonal matrix is ​​chosen, with the diagonal elements denoted as . and The optimal gain that minimizes the fitting error The solution is obtained using nonlinear least squares.

[0100] Using the optimal forward acceleration Υ obtained from the previous prediction period * (k-1), combined with the internal dynamic positive mapping function, construct the cost function J. ID (k):

[0101]

[0102] The cost function J ID (k) constrains the error between the internal dynamic inverse mapping function and its approximate inverse mapping function. Simultaneously, since the reference tilt angle and angular velocity of the robot's center of mass are important control reference quantities in robot balance control, the error value between the reference tilt angle difference and angular velocity obtained from the inverse function is also constrained.

[0103] The cost function J ID (k) The solution is obtained using a nonlinear least squares method. Since the least squares method introduces errors in fitting the nonlinear curve, the optimal fitting matrix is ​​recalculated at the beginning of each prediction period to minimize these errors. The calculation is performed iteratively, with the initial value being the best-fit matrix from the previous calculation. Then, the optimal fitting matrix that currently satisfies the minimum error is... Substituting into the inverse mapping function equation, we can obtain the approximate inverse mapping function expression that satisfies the minimum error, namely:

[0104]

[0105] Similarly, based on the approximate inverse mapping function expression, Ψ in the model prediction inequality constraints can also be obtained. max , The approximate value, that is:

[0106]

[0107]

[0108]

[0109] Where: θ max This represents the maximum angle of inclination of the line connecting the center of mass and the hinge point of the trolley relative to the vertical direction. The optimal fitting matrix The diagonal elements in.

[0110] Finally, based on the prediction process, the optimal solution Y for the current prediction period is selected. * The first term of (k) Combination The current time and Matching optimal reference tilt angle With the optimal reference angular velocity

[0111] (3) Motion control based on multi-joint dynamics optimization

[0112] When considering multi-joint dynamics, a wheeled bipedal robot can be viewed in the sagittal plane as a third-order inverted pendulum system with a cart. The robot's wheels are simplified to a cart, and its upper body is simplified to a third-order inverted pendulum model, as shown below. Figure 3 As shown, the system has three joints: the torso, the thigh, and the lower leg. The lower leg joint is a passive joint, while the torso and thigh joints are active joints and generate two joint torques. The three-stage inverted pendulum model of the trolley also has two external input forces in the horizontal and vertical directions.

[0113] Let the generalized state variables of this model be... Where x and z represent the horizontal and vertical displacements of the trolley; q leg =[q0 q1 q2] represents the joint angles of the lower leg, thigh, and torso; τ leg =[τ1 τ2] T Representing joint torques, τ1 represents the lower leg joint torque, τ2 represents the thigh joint torque; F ext =[F1 F2] T This represents the external force acting on the car, where F1 is the torque τ at the wheel joint. w The converted horizontal force on the ground, i.e., F1 = τ w R w R w F1 is the radius of the wheel, and F2 is the ground support force; the model's dynamic equations are constructed using Lagrange dynamics:

[0114]

[0115] Where M(χ) is the inertia matrix of the third-order inverted pendulum system. G(χ) is the Coriolis force matrix of the third-order inverted pendulum system, G(χ) is the gravity matrix of the third-order inverted pendulum system, H is the joint drive matrix, and J is the joint drive matrix. leg It is the Jacobian matrix of the leg end relative to the world coordinate system.

[0116] Motion control for multi-joint dynamics optimization involves two problems: tracking the robot's reference forward velocity and reference center of mass tilt angle; and obtaining motion control through inverse kinematics combined with the upper-level trajectory. and Capable of solving for optimal reference joint angles and angular velocities (The specific solution process is based on existing technology). Because the frequency predicted by the model is lower than the control frequency of the robot control system, therefore... and It is not necessary to update every control cycle.

[0117] Construct a QP optimization problem, with the optimization variables being: The objective function is A qp It is a weight coefficient matrix, and It is a diagonal matrix, with diagonal elements being the weight coefficients of each optimization variable, b qp It is a feedforward matrix:

[0118]

[0119] in: The optimal reference value of the generalized state variable obtained by trajectory planning The optimal trajectory obtained from the previous trajectory planning is used, but since the longitudinal reference position of the robot's off-ground motion is not considered, then... Since model predictions are not executed in every control cycle, therefore χ ref* (k) represents discrete values; K represents the reference adjustment acceleration of the feedback control output between the generalized state variable and its optimal reference value in QP optimization. p K d It is the ratio and differential coefficient of the feedback; τ is at the reference position leg Feedforward reference value; H -1 It is the pseudo-inverse matrix of H; and These are reference values ​​for F1 and F2, and their main function is to track the forward velocity v of the vehicle. ref With passive degrees of freedom of the legs K pv ,K dv ,K pl ,K dl These are the proportional and derivative adjustment coefficients for the tracking process.

[0120] The equality constraint is:

[0121] The inequality constraints are:

[0122] The purpose of equality constraints is to ensure that the QP optimization problem conforms to the physical model of the actual robot, while the purpose of inequalities is to limit the input force and torque and provide slippage-free constraints; in the above inequality constraints: τ max F corresponds to the joint torque limit values ​​of the lower leg and thigh. max The limit value of the external input force is given by μ, which is the coefficient of friction between the tire and the ground.

[0123] By solving the above QP optimization problem, the optimal joint torque and the optimal external force on the robot wheels are obtained and applied to the bipedal robot for control.

[0124] The embodiments described above are preferred embodiments of the present invention, but the present invention is not limited to the above embodiments. Any obvious improvements, substitutions or modifications that can be made by those skilled in the art without departing from the essence of the present invention shall fall within the protection scope of the present invention.

Claims

1. A hierarchical optimization control method for dynamic walking of a wheeled biped robot, characterized in that: By target speed v tar and current speed v, target position x tar and current position x, the optimal reference forward position optimal reference forward speed optimal reference forward acceleration and optimal future trajectory γ * ; Utilizing the optimal future trajectory γ * and the inner dynamics of the first-order cart-inverted pendulum model, the inner dynamics inverse mapping function is fitted by the nonlinear least squares method, and the optimal reference acceleration matching the optimal reference center of mass inclination θ * ref ; A third-order cart-pole model is established, and the optimal reference joint angles are solved based on and The QP optimization problem related to and is constructed and solved by using multi-joint dynamics, and the optimal control force and torque are obtained to control the biped robot.​ With and The QP optimization problem associated with includes an objective function and constraints, which are: where: P qp is the optimization variable, A qp is the weight coefficient matrix, b qp is the feedforward matrix, and: is the reference adjustment acceleration of the feedback control output between the generalized state variable in QP optimization and its optimal reference value, is the feedforward reference value of the joint torque τ leg under the reference position, F1 ref and is the reference value of the ground horizontal force F1 and the ground support force F2 converted from the wheel joint torque, K p , K d is the proportional and differential coefficient of the feedback, χ ref* (k) is the optimal reference value of the generalized state variable represents the optimal reference position in the z direction, χ is the generalized state variable of the model, H -1 is the pseudo-inverse matrix of the joint driving matrix H, M, C and G represent the inertia matrix, the Coriolis force matrix and the gravity matrix of the three-order car inverted pendulum model, J leg is the Jacobian matrix of the leg end relative to the world coordinate system, F ext represents the external force received by the robot chassis, K pv , K dv , K pl , K dl is the proportional and differential adjustment coefficient in the process of tracking the forward velocity and the passive degree of freedom of the leg, M w is the mass of the robot chassis, q1 represents the joint angle of the robot thigh, q2 represents the joint angle of the robot torso, g is the acceleration of gravity, x and v are the current speed and current position of the robot, represents the optimal reference joint angle of the thigh joint, represents the optimal reference joint angular velocity of the shank joint; the constraint conditions include equality constraints and inequality constraints, and the equality constraints are: The inequality constraints are: where: τ max corresponding to the joint torque limit value of the lower leg and the upper leg, F max is the limit value of the external input force, and μ is the friction coefficient between the tire and the ground.

2. The hierarchical optimization control method for dynamic walking of a wheeled biped robot according to claim 1, wherein, the constraint equation of the model prediction is: v ref (k+i)≤|v max | |a ref (k+i)|≤|Ψ max | where: X ref (k) represents the kth predicted generalized state quantity, x ref (k) represents the kth predicted reference forward position, v ref (k) represents the kth predicted reference forward velocity, x(k), v(k), θ(k), a ref (k) are the actual forward position, actual forward velocity, actual center of mass pitch angle and actual forward acceleration of the robot observed at each single prediction, respectively, Ψ init represents the initial value of Ψ at each prediction, Ψ max represents the maximum value limit of Ψ, Ψ represents the output value of the forward mapping function F of the internal dynamics of the robot, v PD (k+i) represents the kth predicted reference forward velocity and acceleration of the robot, v ref (k+i) respectively represent the reference forward velocity and acceleration of the robot from the kth prediction to the (k+i)th prediction, v ref (k+i) respectively represent the reference forward velocity and acceleration of the robot from the kth prediction to the (k+i)th prediction, v max represents the rated speed of the robot wheel motor, ΔT mpc represents the prediction step length.

3. The hierarchical optimization control method for dynamic walking of a wheeled biped robot according to claim 2, wherein, the model prediction adopts the following manner: where: J N (k) is the cost function predicted by the model, Y(k) is an intermediate quantity, and U ref (k) is the generalized input state quantity predicted at the kth step, N is the total number of steps of prediction; the matrix H α = 2(Φ+Γ T ΩΓ), F α = 2Γ T ΩΘ, Ω, Ф, Г, Θ, Λ are matrices; Q is the error weight matrix during prediction.

4. The hierarchical optimization control method for dynamic walking of a wheeled biped robot according to claim 3, wherein, the optimal reference center of mass tilt angle is obtained by the following manner: the first term of the optimal solution of the current prediction cycle combining and solved, is the optimal fitting matrix, is the fitting matrix, Ψ -1 is the output value of the inner dynamic inverse mapping function, is Ψ -1 approximation, g is the acceleration of gravity.

5. The hierarchical optimization control method for dynamic walking of a wheeled biped robot according to claim 4, wherein, the inner dynamic inverse mapping function is specifically constructed according to the optimal forward acceleration of the last prediction period and the inner dynamic forward mapping function, a cost function about the inner dynamic inverse mapping function and the error of its approximate inverse mapping function is constructed, and the nonlinear least squares method is used to solve and obtain.

6. The hierarchical optimization control method for dynamic walking of a wheeled biped robot according to claim 5, wherein, the cost function is: where: k denotes the kth model prediction, i denotes the prediction i steps ahead, J ID (k) is a cost function.

Images