Humanoid robot whole body motion control method based on model predictive control and reinforcement learning

Abstract

The invention discloses a control method, and the method comprises the steps: collecting the data of a centroid state and an MPC optimal ground acting force in a plurality of motion states of a robot, designing a full-connection residual network, taking the centroid state information as the input, taking the ground acting force as the output, and designing a loss function based on priority. Upper and lower bounds of a friction cone and a sole reaction force are constrained, and the network is trained offline to fit MPC optimization quantity, so that an efficient mode of online reasoning is realized; the joint angle and the angular velocity are fused with the action at the previous moment, the centroid linear velocity and the angular velocity and the foot-ground acting force fitted through a residual network, a bionic reference trajectory is generated through a human motion public data set, different weights are given to a reward function group according to a dominant function value obtained through calculation, and the control performance is improved; symmetrical loss constraint is added to lower limb movement, and the algorithm convergence speed is increased. And finally, training to obtain a whole-body control strategy based on model priori knowledge and a human motion reference trajectory.

Description

Technical Field

[0001] This invention relates to the field of robot control technology, specifically to a method for controlling the whole-body motion of a humanoid robot based on model predictive control and reinforcement learning. Background Technology

[0002] Humanoid robots are extremely difficult to control due to their high degrees of freedom and nonlinear dynamics. Traditional control methods rely on simplified models and pre-programmed gaits, which lack robustness in complex scenarios and struggle to coordinate full-body movements. While pure deep reinforcement learning (DRL) can handle high-dimensional spaces, it suffers from low sample efficiency, is prone to violating physical constraints, and faces difficulties in transferring simulation to the real-world environment. Model predictive control (MPC), while ensuring stability, relies on accurate models, is computationally complex and difficult to implement in real-time, and traditional methods, which separate "movement" and "operation," are prone to dynamic conflicts, resulting in a triple disconnect between "model and data-driven," "motion and operation," and "simulation and reality." The industry urgently needs a fusion framework, and the biomechanical optimization solution of human motion capture data and large-scale parallel simulation technology offer the possibility of integrating MPC physical guidance, DRL adaptability, human motion priors, and domain randomization. This is the core logical starting point of this invention. Summary of the Invention

[0003] To address the problems existing in the prior art, the purpose of this invention is to provide a method for controlling the whole-body motion of a humanoid robot based on model predictive control and reinforcement learning.

[0004] To solve the above problems, the present invention adopts the following technical solution.

[0005] A method for controlling the whole-body motion of a humanoid robot based on model predictive control and reinforcement learning includes the following steps:

[0006] (1) Construct a single rigid body dynamics model of the humanoid robot;

[0007] (2) Establish and discretize the state-space equation of the humanoid robot with physical constraints, and construct a model predictive control optimization problem to solve the optimal foot-ground force sequence in the rolling time domain;

[0008] (3) Collect centroid information and MPC optimization force data pairs under multiple states, and train a fully connected residual network to fit the MPC output;

[0009] (4) By combining the foot-ground force output by the fitted MPC with the redirected human motion data, a parallel training framework for reinforcement learning is built.

[0010] Taking the Unitree H1 humanoid robot as an example, the single-rigid-body dynamics model of the humanoid robot is derived. The derivation process is as follows: the robot's leg links are considered as massless links, and the total mass is concentrated on the floating base. The position of the humanoid robot's center of mass in the world coordinate system and its Euler angles are defined as p.c =[x,y,z] T and Θ=[φ,θ,ψ] T Where φ, θ, and ψ represent the roll angle, pitch angle, and yaw angle of the center of mass, respectively, and the angular velocity in the center-of-mass coordinate system is ω. Assuming that the roll and pitch angles of the robot's center of mass are almost zero during the movement, the relationship between the robot's angular velocity and the Euler angular velocity is obtained as follows:

[0011]

[0012] Where R z (ψ) represents the rotation matrix of the robot's center-of-mass coordinate system relative to the world coordinate system.

[0013] Assume the contact point between the robot's foot and the ground is a finite value N. c From Newton's second law, we can obtain the translational dynamic equation of the center of mass (Newton's equation):

[0014]

[0015] Where m represents the mass of the robot, g is the acceleration due to gravity, and f i Let represent the ground force at the i-th contact point. From the conservation of angular momentum, the rotational dynamics equations of the center of mass (Euler's equations) can be obtained.

[0016]

[0017] Where I B The moment of inertia r represents the rotation of the robot's center of mass. i This represents the vector from the centroid to the i-th contact point. The inertial torque characterized by the angular acceleration of the center of mass, ω×(I B ω) represents the Coriolis and centrifugal torque terms. This term can be ignored when the humanoid robot walks at low speed and has a small angular velocity. Therefore:

[0018]

[0019] Combining the above formulas, the single-rigid-body dynamic model of the humanoid robot can be obtained as follows:

[0020]

[0021] A model-based predictive control method is constructed: the optimal foot-ground force sequence in the rolling time domain is solved using MPC to provide a force reference for reinforcement learning. The specific steps are as follows:

[0022] Define the ground force on the robot's feet as u = [F1, F2, M1, M2] T , where F i =[Fix ,F iy ,F iz ] T M i =[M iy M iz ] T The subscripts i = 1, 2 represent the left and right legs. Based on the above single rigid body dynamics model, the dynamic model of H1 can be obtained as follows:

[0023]

[0024] Where p1, p2 and p c Let represent the contact points between the left and right legs and the ground, and the position of the center of mass, respectively. Constant matrix:

[0025]

[0026] Establish the system state-space equations and define the state variables. Therefore, the following state-space equations are obtained:

[0027]

[0028] in:

[0029]

[0030] I B =R z (ψ)IR z (ψ) T

[0031] From the formula, we can see that matrix A changes with ψ, and matrix B changes with p1-p c p2-p c The state variable changes along with ψ. Discretizing the above state variable equation yields:

[0032]

[0033] The above formula and Take the average yaw angle ψ and the average foot contact point during the robot's movement. To calculate.

[0034] Next, consider the physical constraints. First is the friction cone constraint: the force exerted by the robot's feet in the horizontal direction must be less than the static friction force; otherwise, slippage will occur. Assuming the ground friction coefficient is μ, the following constraints apply:

[0035]

[0036] 0 <F min ≤Fz ≤F max

[0037] The vertical force exerted by the two feet was also constrained at both the upper and lower limits to ensure that the supporting legs remained in contact with the ground while preventing excessive force from causing the robot to bounce up.

[0038] The mapping relationship between the forces acting on the end effector in Cartesian space and the joint torques can be obtained from the principles of statics:

[0039] τ=J T u

[0040] Where J is the Jacobian matrix, τ is the torque at each joint, and ||τ|| n ||≤τ max This means that the joint torque cannot exceed the upper limit.

[0041] Finally, considering the switching between the swing leg and the supporting leg, a selection matrix S can be used to set the foot force of the swing leg to 0.

[0042] Su = 0

[0043] When the left leg swings:

[0044]

[0045] When the right leg swings

[0046]

[0047] Integrating the above equations, given the reference motion trajectory x ref And prediction step size N p This allows us to construct a complete model predictive control optimization problem:

[0048]

[0049]

[0050] 0 <F min ≤F z ≤F max

[0051] ||τ n ||≤τ max

[0052] τ n [k] = J n u n [k],

[0053] Su[k]=0

[0054] k = i, i+1...i+N pn = 1, 2

[0055] We collect centroid information and MPC optimization force data pairs under multiple states, and train a fully connected residual network to fit the MPC network output. The specific method is as follows:

[0056] (1) Data collection and labeling

[0057] The following data were collected under various robot movement states (standing, walking slowly, walking quickly, turning, etc.) to analyze...

[0058] A. State Input: Kinematic information of the center of mass state and reference trajectory

[0059] B. Label Output: Ground Force Obtained by MPC Optimization

[0060] Data collection covered 10 typical states (standing, walking slowly, walking quickly, turning, etc.), generating 10,000 sets of data for each gait, for a total of 100,000 sets of training samples.

[0061] (2) Design and Training of Deep Residual Networks

[0062] The fully connected residual block consists of two paths: A. Main path (feature transformation): in: Fully connected layer weight matrix; Bias vector; BN(·): Batch normalization layer, which standardizes the features (mean is 0, variance is 1); σ(·): GELU activation function (Gaussian error linear unit), expressed as GELU(x)=x·Φ(x) (Φ(x) is the cumulative distribution function of the standard normal distribution); d in ,d hidden ,d out Input, hidden, and output dimensions (all 256 in this invention).

[0063] B. Shortcut Path (Identity Mapping): When the input and output dimensions are the same, the shortcut path is an identity mapping. If the dimensions do not match (e.g., from the input layer to the first residual block), then a 1×1 convolution Conv is introduced. 1×1 Adjusting dimensions:

[0064] C. Residual block output:

[0065] Total parameters: Approximately 260,000 (Input layer: 15×256+256=4160; Each residual block: 256×256×2+256×2=131328, 5 residual blocks total 656640; Output layer: 256×3+3=771; Total: 4160+656640+771=661,571). Forward propagation per sample requires approximately 1.2×10^5 floating-point operations (FLOPs), supporting real-time inference at the edge (0.1ms on NVIDIA Jetson AGX Orin).

[0066] By combining the foot-ground force output from the fitted MPC network with the redirected human motion data, a parallel training framework is constructed. The specific method is as follows:

[0067] (1) State and motion space design

[0068] The high-dimensional state vector integrates information from joints, center of mass, and foot-ground forces, and also tracks key points in the human body data (108 dimensions in total):

[0069]

[0070]

[0071] The robot's key points are selected from six key points in Cartesian space: the shoulder, elbow, and hand. The coordinate system is the body coordinate system with the center of mass as the origin.

[0072] The motion space comprises the robot's 19 joints, including 4 joints in the arms, 1 joint in the waist, and 10 joints in the legs. The policy network output is the increment Δq of the joint position. i The motor force is calculated by the PD controller to form a position closed loop:

[0073]

[0074] Where K p For proportional gain, K d For differential gain, q default To initialize joint positions, This represents the current joint velocity.

[0075] (2) Reward Function Design

[0076] The reward function r(s,a) consists of three weighted parts: r = ω reg ·r reg +ω track ·r track +ω mpc ·r mpc

[0077] A. Regular reward rreg Encourage standing and punish abnormal joint movements.

[0078] • Standing stability bonus (based on center of gravity height and tilt angle):

[0079]

[0080] Where z com For the actual center of mass height, z com,ref =1.05 is the reference centroid height. k1 and k2 are the sum of the squares of the roll and pitch angles, respectively, and are hyperparameters.

[0081] • Joint speed penalty (limits joint overspeed):

[0082] Where N=19 is the number of driving joints. The joint angular velocity, denoted as , where is the maximum angular velocity of each joint, and k3 is a hyperparameter.

[0083] • Joint acceleration penalty (inhibiting joint impact): in Joint angular acceleration, K represents the maximum angular acceleration of each joint, and k4 is a hyperparameter.

[0084] B. Key Point Tracking Rewards track Encourage robots to track actions in datasets.

[0085]

[0086] Where: M=6 is the number of tracking key points (left and right shoulders, elbows, and hands); p j,act p represents the actual key point location of the robot. j,ref k is the location of the reference keypoint in the dataset; k5 is a hyperparameter.

[0087] C. MPC Reference Value Reward mpc Encourage robots to contact the ground with appropriate force.

[0088]

[0089] in: The reference ground force output by the residual network; The actual ground force calculated for the simulation environment; F weight =500 is the robot's gravity. As an indicator function, this reward only applies upon ground contact.

[0090] (3) Training process of the improved proximal policy optimization algorithm (PPO)

[0091] An improved PPO algorithm employing a multi-critic network, generative adversarial networks, and symmetric loss is described below.

[0092] A. Multi-critic network advantage estimation: Mean and variance of the output value functions of three parallel evaluation networks: A multi-critic architecture (each critic corresponds to a set of reward functions) is adopted, combining multiple independent advantage functions. The overall advantage A is aggregated as follows:

[0093]

[0094] in, This involves standardizing each dominance function (subtracting the batch mean μ and dividing by the standard deviation σ) to eliminate differences in numerical scales between different dominance functions and prevent the dominance of a particular critic from being overemphasized. i The weights assigned to each standardized advantage are used to adjust the importance of different reward groups in the overall strategy update.

[0095] B. Symmetrical Loss Function: Utilizing the robot's left-right symmetry, a symmetric loss term is introduced:

[0096] By constructing a left-right symmetry transformation matrix for the robot, the mirror mapping of the robot's actions is realized. The policy network is forced to output symmetrical actions for left-right symmetry states, and the generalization of the policy is improved by utilizing the structural symmetry of the robot.

[0097]

[0098] Where: λ = 0.5 is the loss weight coefficient; π(o) is the policy output action corresponding to the original state o; π rec (o) is the mirror recovery action (the action of restoring the image after symmetric transformation).

[0099] C. Generative Adversarial Networks (GANs) are used to train a discriminator to distinguish between state transitions generated by the policy and those of experts in the dataset. Adversarial learning constrains the state transition distribution output by the policy to approximate the expert data distribution, improving the realism of the policy's actions. The total discriminator loss consists of two parts: adversarial loss and gradient penalty loss.

[0100]

[0101] For real samples in the reference dataset D, the discriminator is encouraged to output values ​​close to 1 (labeled "real"), achieved by minimizing the squared error:

[0102]

[0103] For agent policy π θ The generated samples encourage the discriminator to output values ​​close to -1 (marked as "fake"), also achieved by minimizing the squared error:

[0104]

[0105] To limit the magnitude of the discriminator's gradient on real samples, if the gradient of the discriminator on real data samples is non-zero, it will cause the agent to generate samples that "deviation" from the real data distribution. Gradient penalty can alleviate this problem, reduce oscillations in GAN training, and improve stability. The expression is as follows:

[0106]

[0107] Where w gp and λ gp This is a hyperparameter.

[0108] The first two items combine the design of least squares GAN (LSGAN) with the aim of minimizing the Pearson divergence between the reference data distribution and the agent generation distribution, so that the agent behavior distribution is as close as possible to the reference data.

[0109] D. Global randomization strategy

[0110] During training, a global randomization strategy was employed to enhance the robustness of the strategy. Specific measures included: randomly perturbing the center of mass weight within a range of ±10% (e.g., the base mass of the Yushu H1 robot is 50kg, and the randomization range is 45-55kg) to simulate load changes or minor deviations in mass distribution; adding a random offset of ±5° to the initial joint position to restore the small error in the initial joint posture during actual machine startup; and randomly adjusting the proportional gain Kp (base value 200) and differential gain Kd (base value 5) of the PD controller within a range of ±15% to simulate the aging or calibration deviation of the actuator control parameters, ensuring that the strategy adapts to the joint response characteristics under different control parameters.

[0111] Beneficial effects of the present invention

[0112] Compared with the prior art, the advantages of this invention are:

[0113] 1. The innovative "single rigid body dynamics model + MPC optimization with physical constraints" mechanism effectively avoids problems such as foot slippage and joint overload by using constraints such as friction cones and joint torques, while maintaining physical rationality and improving the stability control capability of the robot's whole-body motion.

[0114] 2. Based on the strategy of "fitting MPC output with fully connected residual network", this paper addresses the pain points of high computational cost of online MPC solution and difficulty in adapting to parallel training of reinforcement learning. It achieves highly efficient reference generation of "offline training-online inference" and greatly improves the real-time performance of control.

[0115] 3. This method is applicable to walking, balancing and whole-body coordinated operation tasks of high-degree-of-freedom humanoid robots. In complex scenarios such as model errors, sensor noise and rough terrain, it shows better stability and environmental adaptability compared to pure DRL or pure MPC methods.

[0116] 4. By using a parallel training framework to fuse and fit the foot-ground force and redirected human motion data of MPC, the stability of model-driven and the adaptability of data-driven approaches are significantly integrated, greatly improving the robot's walking robustness under dynamic disturbances and its ability to transfer simulation data to real objects. Attached Figure Description

[0117] Figure 1 This is the overall control flowchart of the method of the present invention;

[0118] Figure 2 This is a schematic diagram of the model;

[0119] Figure 3 The torque curves of the two legs when walking using the above controller;

[0120] Figure 4 The above controller is used to generate curves for center of mass position tracking and velocity tracking during walking.

[0121] Figure 5 This refers to the foot placement point when walking using the aforementioned controller;

[0122] Figure 6 Select a schematic diagram for the key points. Detailed Implementation

[0123] The technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments of the present invention. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without creative effort are within the scope of protection of the present invention.

[0124] Please see Figures 1 to 6A method for controlling the whole-body motion of a humanoid robot based on model predictive control and reinforcement learning, comprising the following steps:

[0125] 1. Construct a single rigid body dynamic model for a humanoid robot.

[0126] The single-rigid-body dynamic model of the humanoid robot is as follows:

[0127]

[0128] 2. Constructing a model for predictive control optimization problems

[0129]

[0130] 0 <F min ≤F z ≤F max

[0131] ||τ n ||≤τ max

[0132] τ n [k] = J n u n [k],

[0133] Su[k]=0

[0134] k = i, i+1...i+N p n = 1, 2

[0135] 3. Collect centroid information and MPC optimization force data pairs, and train a fully connected residual network to fit the MPC output. (1) Data collection and labeling

[0136] The following data were collected under various robot movement states (standing, walking slowly, walking quickly, turning, etc.) to analyze...

[0137] A. State Input: Kinematic information of the center of mass state and reference trajectory

[0138] B. Label Output: Ground Force Obtained by MPC Optimization

[0139] Data collection covered 10 typical states (standing, walking slowly, walking quickly, turning, etc.), generating 10,000 sets of data for each gait, for a total of 100,000 sets of training samples.

[0140] (2) Design and Training of Deep Residual Networks

[0141] The fully connected residual block consists of two paths: A. Main path (feature transformation): in: Fully connected layer weight matrix; Bias vector; BN(·): Batch normalization layer, which standardizes the features (mean is 0, variance is 1); σ(·): GELU activation function (Gaussian error linear unit), expressed as GELU(x)=x·Φ(x) (Φ(x) is the cumulative distribution function of the standard normal distribution); d in ,d hidden ,d out Input, hidden, and output dimensions (all 256 in this invention).

[0142] B. Shortcut Path (Identity Mapping): When the input and output dimensions are the same, the shortcut path is an identity mapping. If the dimensions do not match, a 1×1 convolution Conv is introduced. 1×1 Adjusting dimensions:

[0143] C. Residual block output:

[0144] D. Loss Function Design

[0145] D1 core loss: The weighted mean squared error (MSE) is expressed as follows:

[0146]

[0147] D2 Constraint Penalty: Physical feasibility constraint of force, targeting the fundamental physical constraints of ground action forces (avoiding...).

[0148]

[0149] Friction cone penalty: To prevent excessive horizontal force from causing foot slippage; when the predicted force exceeds the friction cone's range. When μ = 0.6 is the coefficient of ground friction, a second penalty is applied:

[0150]

[0151] Vertical force upper and lower bound penalty: Ensures the support force is within a safe range; when the predicted vertical force exceeds the range of [50N, 1500N], a penalty is applied.

[0152]

[0153] Penalty coefficient: λ f=10 (priority of enhanced friction constraints), λ z =5 Ensures the vertical force safety range). D3 total loss function is the sum of core fitting loss and force constraint penalty term:

[0154]

[0155] 4. By combining the foot-ground force output from the fitted MPC with the redirected human motion data, a parallel training framework for reinforcement learning is built.

[0156] The above description is merely a preferred embodiment of the present invention; however, the scope of protection of the present invention is not limited thereto. Any equivalent substitutions or modifications made by those skilled in the art within the scope of the technology disclosed in the present invention, based on the technical solution and its improved concepts, should be covered within the scope of protection of the present invention.

Claims

1. A method for controlling the whole-body motion of a humanoid robot based on model predictive control and reinforcement learning, characterized by the following: step: (1) Construct a single rigid body dynamics model of the humanoid robot; (2) Establish and discretize the state-space equation of the humanoid robot with physical constraints, and construct a model predictive control optimization problem to solve the optimal foot-ground force sequence in the rolling time domain; (3) Collect centroid information and MPC optimization force data pairs under multiple states, and use a composite loss function combining weighted mean square error with friction cone and vertical force upper and lower bound constraint penalty terms to train a fully connected residual network to accurately fit the ground action force output by MPC, while ensuring the physical feasibility of the network's predicted force. (4) By combining the foot-ground force output by the fitted MPC with the human motion dataset, a parallel training framework for reinforcement learning is built.

2. The method for controlling the whole-body motion of a humanoid robot based on model predictive control and reinforcement learning according to claim 1, characterized in that: In the construction of the humanoid robot's single rigid body dynamic model in step (1), the robot's leg links are considered as massless links, and the entire body mass is concentrated on the floating base. The position of the humanoid robot's center of mass in the world coordinate system and its Euler angles are defined as p. c =[x,y,z] T and Θ=[φ,θ,ψ] T Where φ, θ, and ψ represent the roll, pitch, and yaw angles of the center of mass, respectively, and the angular velocity in the center-of-mass coordinate system is ω. Assuming that the roll and pitch angles of the robot's center of mass are almost zero during the walking process, the relationship between the robot's angular velocity and the Euler angular velocity is obtained as follows: Where R z (ψ) represents the rotation matrix of the robot's center-of-mass coordinate system relative to the world coordinate system. Assume the contact point between the robot's foot and the ground is a finite value N. c From Newton's second law, we can obtain the translational dynamic equation of the center of mass (Newton's equation): Where R z (ψ) represents the rotation matrix of the robot's center-of-mass coordinate system relative to the world coordinate system. Assume the contact point between the robot's foot and the ground is a finite value N. c From Newton's second law, we can obtain the translational dynamic equation of the center of mass: Where m represents the mass of the robot, g is the acceleration due to gravity, and f i Let the ground force at the i-th contact point be represented. By conservation of angular momentum, the rotational dynamic equation of the center of mass can be obtained: Where I B The moment of inertia r represents the rotation of the robot's center of mass. i I represents the vector from the centroid to the i-th contact point. B ω represents the inertial torque caused by the angular acceleration of the center of mass, ω×(I B ω) represents the Coriolis and centrifugal torque terms. This term can be ignored when the humanoid robot walks at low speed and has a small angular velocity. Therefore: Combining the above formulas, the single-rigid-body dynamic model of the humanoid robot can be obtained as follows:

3. The method for controlling the whole-body motion of a humanoid robot based on model predictive control and reinforcement learning according to claim 2, characterized in that: In step (2), the robot's foot-ground force u is defined as [F1, F2, M1, M2]. T , where F i =[F ix ,F iy ,F iz ] T M i =[M iy M iz ] T The subscripts i = 1, 2 represent the left and right legs, and the dynamic model of H1 can be obtained as follows: Where p1, p2 and p c Let represent the contact points between the left and right legs and the ground, and the position of the center of mass, respectively. Constant matrix: Establish the system state-space equations and define the state variables. Therefore, the following state-space equations are obtained: in: I B =R z (ψ)AND z (ψ) T In the formula, matrix A changes with ψ, and matrix B changes with p1-p c p2-p c The state variable equation changes with ψ. Discretizing the above state variable equation yields: The above formula and Take the average yaw angle ψ and the average foot contact point during the robot's movement. To calculate, Friction cone constraint: The force exerted by the robot's feet in the horizontal direction must be less than the static friction force; otherwise, slippage will occur. Assuming the ground friction coefficient is μ, the following constraints apply: 0<F min ≤F z ≤F max The vertical force exerted by the legs was also constrained, both vertically and horizontally, to ensure the supporting legs remained in contact with the ground while preventing excessive force from causing the robot to bounce up. The mapping relationship between the forces acting on the end effector in Cartesian space and the joint torques can be obtained from the principles of statics: τ=J T you Where J is the Jacobian matrix, τ is the torque at each joint, and ||τ|| n ||≤τ max This means that the joint torque must not exceed the upper limit. Finally, considering the switching between the swing leg and the supporting leg, a selection matrix S can be used to set the foot-ground force of the swing leg to 0: Su = 0 When the left leg swings: When the right leg swings Integrating the above equations, given the reference motion trajectory x ref And prediction step size N p This allows us to construct a complete model predictive control optimization problem: 0<F min ≤F z ≤F max ||t n ||≤τ max τ n [k]=J n u n [k], Su[k]=0 k=i,i+1...i+N p ,n=1,2。 4. The method for controlling the whole-body motion of a humanoid robot based on model predictive control and reinforcement learning according to claim 3, characterized in that: In step (3), the following data pairs (s, ): A. State Input: Kinematic information of the center of mass state and reference trajectory B. Label Output: Ground Force Obtained by MPC Optimization The data collection covers the following typical states, with 10,000 sets of data generated for each state, for a total of 90,000 training samples, which are divided into training and validation sets in an 8:2 ratio. Static standing test: center of gravity 1m high, maintained for 10s, sampling frequency 1kHz, generating 10,000 sets of data. Slow walking: stride length 0.3m, cadence 0.8Hz, speed 0.24m / s, 500 gait cycles were collected, and 10,000 valid samples were obtained after cropping. Brisk walking: stride length 0.5m, cadence 1.2Hz, speed 0.6m / s, 400 gait cycles were collected, for a total of 10,000 sets. Turning in place: angular velocities of 5° / s, 10° / s, and 15° / s were collected, with approximately 3,333 sets each, for a total of 10,000 sets. Lateral movement: stride length 0.15m, stride frequency 0.8Hz, speed 0.12m / s, 500 gait cycles were collected, and 10,000 sets were obtained after clipping. Standing under load: ~3,333 sets of data were collected each at a center of gravity of 5kg, 10kg, and 15kg, for a total of 10,000 sets. Walking with a 5kg load: stride length 0.4m, cadence 1.0Hz, speed 0.4m / s, 400 gait cycles were collected, for a total of 10,000 sets. Disturbance resistance: Apply horizontal pulse forces of 5N and 10N for 0.2s while standing or walking slowly. Collect approximately 5,000 data sets for each condition, for a total of 10,000 data sets. Dynamic start-stop: 5,000 data sets were collected each time from acceleration from a standstill to 0.5 m / s, 2 seconds, and deceleration from 0.5 m / s to a standstill, for a total of 10,000 data sets; (2) Design and Training of Deep Residual Networks Using fully connected residual blocks with two paths composition: A. Main path: in: : Weight matrix of fully connected layer; : Bias vector; BN(·): Batch normalization layer, which standardizes the features to a mean of 0 and a variance of 1; σ(·): GELU activation function, expressed as GELU(x)=x·Φ(x) (Φ(x) is the cumulative distribution function of the standard normal distribution); d in ,d hidden ,d out Input, hidden, and output dimensions. B. Shortcut Path: When the input and output dimensions are the same, the shortcut path is an identity mapping. If the dimensions do not match, a 1×1 convolution Conv is introduced. 1×1 Adjusting dimensions: C. Residual block output: Total parameters: approximately 260,000; Input layer: 15×256+256=4160; Each residual block: 256×256×2+256×2=131328, 5 residual blocks total 656640; Output layer: 256×3+3=771; Total: 4160+656640+771=661,571. Forward propagation of a single sample requires approximately 1.2×10^5 floating-point operations, supporting real-time inference at the edge. D. Loss function design: In the training of fully connected residual networks, a composite loss function combining weighted mean square error and force constraint penalty term is used to constrain the physical feasibility of ground action forces. The specific definition is as follows: D1 Core Loss: The weighted mean square error assigns weights based on the physical importance differences of the 6-dimensional foot-ground forces (i.e., the three-dimensional forces of the left and right feet), focusing on improving the fitting accuracy of the vertical forces that have the greatest impact on support stability. The expression is: Wherein: F pred The predicted foot-ground forces L / R represent the left and right feet respectively, (x / y / z) corresponding to the three-dimensional directions, and F mpc The target force for MPC optimization; the weighting coefficients are set according to the physical action of the force: w z =1.5, w x =w y =1.0; mean square error is defined as N is the batch size. D2 constraint penalty term: The physical feasibility constraint of the force is designed with a penalty term based on the fundamental physical constraints of the ground action force to ensure that the network's predicted force conforms to the actual interaction law. Friction cone penalty: To prevent excessive horizontal force from causing foot slippage; when the predicted force exceeds the friction cone's range. When μ = 0.6 is the coefficient of ground friction, a second penalty is applied: Vertical force upper and lower bound penalty: Ensures the support force is within a safe range; when the predicted vertical force exceeds the range of [50N, 1500N], a penalty is applied. Penalty coefficient: λ f =10, λ z =5 Ensures the vertical force safety range; D3 total loss function is the sum of the core fitting loss and the force constraint penalty term:

5. The method for controlling the whole-body motion of a humanoid robot based on model predictive control and reinforcement learning according to claim 4, characterized in that: In step (4), (1) the key points of the robot body are selected from the Cartesian space positions of 6 key points: shoulder, elbow and hand. The coordinate system is the body coordinate system with the center of mass as the origin. The motion space consists of the robot's 19 joints, including 4 joints in the arms, 1 joint in the waist, and 10 joints in the legs. The policy network output is the increment Δq of the joint position. i The motor force is calculated by the PD controller to form a position closed loop: Where K p For proportional gain, K d For differential gain, q default To initialize joint positions, This represents the current joint velocity. (2) Reward function design: The reward function r(s,a) consists of three weighted parts: r = ω reg ·r reg +ω track ·r track +ω mpc ·r mpc , A. Regular reward r reg Encourage standing and punish abnormal joint movements. • Standing stability bonus (based on center of gravity height and tilt angle): Where z com For the actual center of mass height, z com,ref =1.05 is the reference centroid height. k1 and k2 are the sum of the squares of the roll and pitch angles, respectively, and are hyperparameters. • Joint speed penalty (limits joint overspeed): Where N=19 is the number of driving joints. The joint angular velocity, k is the maximum angular velocity of each joint, and k3 is a hyperparameter. • Joint acceleration penalty (inhibiting joint impact): in Joint angular acceleration, k is the maximum angular acceleration of each joint, and k4 is a hyperparameter; B. Key Point Tracking Rewards track Encourage robots to track actions in datasets. Where: M = 6 is the number of tracking key points; p j,act p represents the actual key point location of the robot. j,ref k5 represents the reference keypoint locations in the dataset; k5 is a hyperparameter. C. MPC Reference Value Reward mpc Encourage robots to contact the ground with appropriate force. in: The reference ground force output by the residual network; The actual ground force calculated for the simulation environment; F weight =500 is the robot's gravity. As an indicator function, this reward only applies upon ground contact. (3) Training process of the improved proximal policy optimization algorithm (PPO): An improved PPO algorithm employing a multi-critic network and symmetric loss is described below, with the core formula as follows: A. Multi-evaluation network dominance estimation: Mean and variance of the output value functions of three parallel evaluation networks: Employing a multi-critic architecture, multiple independent advantage functions Aφ are used. i The overall advantage A is aggregated as follows: in, This involves standardizing each dominance function by subtracting the batch mean μ and dividing by the standard deviation σ. The aim is to eliminate differences in numerical scale between different dominance functions and prevent the dominance of a particular critic from being overemphasized. i The weights for each standardized advantage are used to adjust the importance of different reward groups in the overall strategy update. B. Symmetrical Loss Function: Utilizing the robot's left-right symmetry, a symmetric loss term is introduced: By constructing a left-right symmetry transformation matrix to achieve mirror mapping of robot actions, the policy network is forced to output symmetrical actions for left-right symmetry states, thereby improving the generalization of the policy by leveraging the structural symmetry of the robot. Where: λ = 0.5 is the loss weight coefficient; π(o) is the policy output action corresponding to the original state o; π rec (o) represents the image restoration action. C. Generative Adversarial Networks (GANs): These networks train a discriminator to distinguish between state transitions generated by the policy and those of experts in the dataset. Adversarial learning constrains the state transition distribution output by the policy to approximate the expert data distribution, improving the realism of the policy's actions. The total discriminator loss consists of two parts: adversarial loss and gradient penalty loss. For real samples in the reference dataset D, the discriminator is encouraged to output values ​​close to 1, which are then labeled as "real". This is achieved by minimizing the squared error. For agent policy π θ The generated samples encourage the discriminator to output values ​​close to -1 as "fakes," also achieved by minimizing the squared error: To limit the magnitude of the discriminator's gradient on real samples, if the gradient of the discriminator on real data samples is non-zero, it will cause the agent to generate samples that "deviation" from the real data distribution. Gradient penalty can alleviate this problem, reduce oscillations in GAN training, and improve stability. The expression is as follows: Where w gp and λ gp For hyperparameters, The first two items combine the design of least squares GAN, aiming to minimize the Pearson divergence between the reference data distribution and the agent generation distribution, so that the agent behavior distribution is as close as possible to the reference data. D. Global randomization strategy During training, a global randomization strategy was employed to enhance the robustness of the strategy. Specific measures included: randomly perturbing the center of mass weight within a range of ±10% to simulate load changes or minor deviations in mass distribution; adding a random offset of ±5° to the initial joint position to reproduce the small errors in the initial joint posture during actual machine startup; and randomly adjusting the proportional gain Kp and derivative gain Kd of the PD controller within a range of ±15% to simulate the aging or calibration deviations of the actuator control parameters, ensuring that the strategy adapts to the joint response characteristics under different control parameters.

Images