Distributed drive heavy truck flexible support gravity center coupling stability control method and system

Abstract

The application provides a distributed driving heavy truck flexible support gravity center coupling stability control method and system. The method is based on vehicle posture, dynamic wheel load and suspension spring-damping relationship fusion optimization to solve each suspension deformation variable, reconstructs the flexible support plane, obtains three-dimensional dynamic gravity center and its estimated covariance in combination with wheel load moment balance; and recursively suspensions deformation, wheel load, gravity center trajectory and dynamic support domain in the prediction time domain, and makes a graded stability determination according to the robust stability margin, lateral load transfer rate, minimum wheel load ratio and instability crossing time; through the limited prediction time domain optimization calculation target additional yaw moment, the wheel end motor torque and brake force joint constraint optimization is adopted to complete the execution distribution, so as to improve the vehicle instability risk identification and early control ability under heavy load, partial load and complex road conditions.

Description

Technical Field

[0001] This invention relates to the field of active safety and chassis cooperative control technology for heavy vehicles, specifically to a distributed drive heavy truck flexible support center of gravity coupling stability control method and system. Background Technology

[0002] Distributed drive heavy-duty trucks, relying on the independent and controllable drive and braking torque output characteristics of each wheel, possess outstanding advantages such as fast wheel-end response speed, high degree of freedom in vehicle control, and flexible construction of additional yaw control torque. This represents a significant development direction for active safety and intelligent chassis control in heavy-duty vehicles. Compared to traditional centralized drive heavy-duty trucks, the distributed drive configuration can achieve refined control of vehicle dynamics through independent adjustment of the force at each wheel end, exhibiting excellent stability control potential. However, heavy-duty trucks have structural and operational characteristics such as large vehicle mass, high center of gravity installation position, wide dynamic fluctuation range of transport load, and complex and variable driving conditions. Under typical conditions such as cornering, emergency lane changes, hill starts, off-center loading, and driving on uneven roads, vehicles are prone to significant lateral load transfer and substantial attenuation of the inner wheel load, which can induce dangerous conditions such as vehicle yaw instability, wheel lift-off, and rollover, seriously affecting vehicle driving safety and adaptability to various operating conditions.

[0003] Currently, mainstream technologies for vehicle lateral stability control are generally based on the rigid body assumption and fixed center of gravity modeling theory for controller design. They primarily use single or combined state parameters such as yaw rate, center of gravity sideslip angle, lateral acceleration, and lateral load transfer rate as criteria for vehicle stability. While this modeling and control approach achieves good control results under light vehicle conditions, paved roads, and uniform loads, it struggles to adapt to the complex dynamic characteristics of heavy trucks. During actual heavy truck operation, dynamic factors such as asymmetric suspension compression, flexible elastic deformation of the chassis, random vertical excitation from the road surface, and uneven load distribution caused by cargo offset continuously alter the vehicle's actual support posture, causing the projected position of the vehicle's center of gravity relative to the ground support area to shift dynamically in real time. Continuing to use traditional modeling and control logic with a fixed center of gravity and fixed wheelbase support boundaries fails to accurately represent the true dynamic characteristics of heavy trucks, significantly underestimating the degree of center of gravity shift and the tendency for wheels to lift off the ground. This results in problems such as delayed intervention of stability control strategies, overly optimistic stability margin judgments, and insufficient control accuracy under extreme conditions, making it difficult to meet the safety control requirements under heavy loads, complex road conditions, and dynamic off-center loading.

[0004] To address the issue that fixed center of gravity and rigid support modeling cannot accurately reflect the dynamic support posture of heavy-duty trucks, existing publicly available technologies have focused on distributed vehicle stability control and wheel-end force distribution. For example, CN110239621B delineates the stability region based on the road adhesion coefficient and the real-time vehicle state, and achieves distributed electric vehicle yaw stability control by solving for the desired additional yaw moment; CN112644455B optimizes the yaw control moment through multi-objective joint decision-making based on the differences in operating conditions inside and outside the vehicle's stability domain; CN112668093A proposes a distributed drive vehicle all-wheel longitudinal force optimization distribution strategy. However, existing solutions mainly focus on vehicle yaw moment regulation and wheel-end drive force distribution optimization, without considering the changes in vehicle support posture caused by the flexible deformation of the heavy-duty truck frame and the dynamic asymmetric deformation of the multi-axle suspension. The modeling process ignores the dynamic offset characteristics of the center of gravity, resulting in significant shortcomings in model accuracy and adaptability to complex operating conditions of heavy-duty trucks.

[0005] Beyond the aforementioned yaw control and wheel-end force distribution schemes, there are also technical deficiencies in suspension state identification and stability boundary construction. Most current heavy-duty truck stability control schemes rely solely on macroscopic state quantities such as vehicle roll angle and left / right wheel load difference to infer suspension deformation states. This fails to effectively isolate multi-dimensional coupled interference factors such as vehicle attitude changes, road surface undulations, static compression of suspension springs, and dynamic damping force of shock absorbers, resulting in low accuracy in suspension deformation identification and poor reliability in state estimation. Furthermore, existing stability boundaries are all statically defined based on a fixed wheelbase, neglecting the dynamic attenuation characteristics of the effective ground contact area as the wheels approach ground clearance, and ignoring the physical laws governing the dynamic changes in the effective support domain of multi-axle heavy-duty trucks with wheel load. Furthermore, existing technologies do not incorporate estimation errors caused by vehicle state observation and parameter identification into the safety margin design, resulting in insufficient robustness in stability determination and risks of missed or misjudged instability. For example, CN110267832A only uses road surface preview information to achieve active adjustment of suspension damping force, and only optimizes the performance of a single dimension of the suspension. It does not relate to the flexible support posture of the whole vehicle, dynamic center of gravity offset and dynamic evolution law of the support domain, and cannot solve the dynamic instability problem of heavy trucks at the whole vehicle level.

[0006] In summary, existing distributed drive heavy truck stability control technologies mostly focus independently on single modules such as yaw stability control, wheel-end force distribution, and active suspension control. They have not built an integrated and collaborative control system that integrates suspension deformation identification, flexible support plane reconstruction, dynamic three-dimensional center of gravity prediction, and dynamic support domain correction for the operating characteristics of multi-axle heavy trucks. As a result, they cannot solve the core problem of vehicle dynamic instability under heavy load, complex road surface, and off-center load conditions. Summary of the Invention

[0007] To address the shortcomings of existing technologies, this invention provides a distributed drive heavy truck flexible support center of gravity coupling stability control method and system.

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

[0009] A distributed drive heavy-duty truck flexible support center of gravity coupling stability control method includes the following steps:

[0010] Step S1: Unify controller input, collect road height sequence and road adhesion coefficient, and calculate road vertical excitation;

[0011] Step S2: Estimate the vehicle's longitudinal velocity, lateral velocity, center of gravity sideslip angle, and vehicle body heave; calculate the tire longitudinal force and tire lateral force; and reconstruct the dynamic wheel load of each wheel.

[0012] Step S3: Each suspension connection point is used as a flexible support node. The suspension geometric deformation is obtained based on the vehicle body heave, roll, pitch and road vertical excitation. The suspension geometric deformation, suspension spring-damping dynamics and dynamic wheel load consistency relationship are constrained and fused to obtain the optimized suspension deformation and wheel load consistency residual of each suspension connection point.

[0013] Step S4: Based on the optimized suspension deformation, road vertical excitation and initial installation height of suspension connection points, restore the spatial height of each support point, and use a robust weighted method to fit the flexible support plane to obtain the parameter estimates, unit normal vector, equivalent roll angle and equivalent pitch angle of the flexible support plane.

[0014] Step S5: The static center of gravity is geometrically transformed according to the pose change of the flexible support plane to form a prior geometric center of gravity. The prior geometric center of gravity is constrained and corrected by the longitudinal and lateral moment balance relationship formed by the dynamic wheel load of each wheel to obtain the current dynamic three-dimensional center of gravity of the vehicle and its estimated covariance.

[0015] Step S6: Considering the driver's future driving input, based on steps S2 to S5, recalculate the flexible support plane, dynamic three-dimensional center of gravity, and dynamic support domain formed by the effective contact wheel in each prediction step in the future prediction time domain.

[0016] Step S7: Use robust stability margin, lateral load transfer rate, minimum wheel load ratio, and instability crossing time to construct a predictive time-domain comprehensive risk index to judge vehicle stability. The vehicle judgment results include stable state, warning state, instability state, or rollover risk state.

[0017] Step S8: When the vehicle is in a warning state, an unstable state, or a rollover risk state, the target additional yaw moment is calculated by finite prediction time domain optimization, and the execution allocation is completed by joint constraint optimization of wheel-end motor torque and braking force.

[0018] Furthermore, in step S1, the controller input is uniformly represented as:

[0019]

[0020] In the formula, For the first The controller input vector at each sampling time; The longitudinal speed of the vehicle; , , These are the longitudinal, lateral, and vertical accelerations in the vehicle coordinate system, respectively. , , These are the roll, pitch, and yaw angular velocities, respectively. , These are the roll angle and the pitch angle, respectively. This refers to the steering angle; , and The first Wheel speed, wheel-end motor torque, and braking pressure of each wheel; This represents the total number of wheel ends. Indicates the first The longitudinal distance of each wheel relative to the road surface sensing device;

[0021] The formula for vertical excitation of the road surface is expressed as:

[0022]

[0023]

[0024] In the formula, For the first The road surface excitation propagates to the first sampling time. The time it takes for a wheel to pass For the first The first wheel in the Vertical excitation of the road surface at each sampling time, For the first The longitudinal speed of the vehicle at each sampling time. For the first Suspension deformation at each sampling time, To minimize the converted vehicle speed, This is a sequence of road surface heights.

[0025] Furthermore, in step S2, the dynamic wheel load reconstruction formula for each wheel is:

[0026]

[0027] In the formula, For the first Estimated dynamic wheel load vectors for each wheel at each sampling time. For dynamic wheel load, For the first The prior load at each sampling time, This represents the vertical resultant force equilibrium constraint coefficient matrix; The coefficient matrix representing the balance constraints of longitudinal and lateral moments. and All are positive definite weight matrices.

[0028] Furthermore, in step S3, the suspension geometry deformation... for:

[0029]

[0030] In the formula, This is the suspension geometry mapping matrix; Based on the first step S1, the first... The set of vertical excitations of the road surface at each sampling time; For the first Generalized coordinates of vehicle body heave, roll, and pitch at each sampling moment;

[0031] The suspension spring-damping dynamics relationship is as follows:

[0032]

[0033] In the formula, The equivalent mass matrix for each support node; and These are the equivalent damping matrix and stiffness matrix of the suspension, respectively; and These are deformation velocity and deformation acceleration, respectively. This is the dynamic wheel load estimation vector obtained based on step S2. This is a static wheel load; For equivalent perturbation;

[0034] The formula for the optimized suspension deformation at each suspension connection point is expressed as follows:

[0035]

[0036] In the formula, For the first The suspension deformation variables after fusion and optimization at each sampling time point This represents the suspension deformation vector to be solved; and Both are weight matrices; For spatial or interaxial structural continuity matrices; The structural continuity regularization weight coefficient; No. The equivalent input of wheel load dynamics at each sampling time, For coefficients, For the first Geometric variables at each sampling time.

[0037] Furthermore, in step S4, the spatial height of the support point is:

[0038]

[0039] In the formula, For the first The sampling time of the first sampling moment Current spatial height of each support point; For the first Initial installation height of each support point; For vertical excitation of the road surface; The optimized suspension deformation;

[0040] The estimated parameters of the flexible support plane are:

[0041]

[0042] In the formula, For the estimation of the plane parameters of the flexible support, For the first Parameters of the flexible support plane at each sampling time. For the first The sampling time of the first sampling moment Fitting weights for each support point; The Huber robust loss function; For the first Planar regression vector of each support point; For planar time regularization weights; This represents the total number of wheel ends.

[0043] Furthermore, in step S5, the geometric centroid is prior. for:

[0044]

[0045] In the formula, For the first Translation vector of the flexible support plane at each sampling time; For the first The support attitude rotation matrix at each sampling time; As the center of gravity; This is the compliance mapping matrix from suspension deformation to equivalent center of gravity shift; The compliance mapping matrix caused by acceleration; The optimized suspension deformation; For static suspension deformation; and The first The longitudinal and lateral accelerations of the vehicle at each sampling time;

[0046] The dynamic three-dimensional center of gravity is:

[0047]

[0048] In the formula, For the first The dynamic three-dimensional centroid estimate at each sampling time; Center of gravity Feasible region; For the first The coefficient matrix of wheel load torque balance relationship at each sampling time; For the first The torque balance constraint vector at each sampling time; This is the weight matrix for the wheel load torque balance terms; Let be the weight matrix of the prior terms of the geometric centroid.

[0049] Furthermore, in step S6, the dynamic support domain is defined as:

[0050]

[0051] In the formula, To ensure effective contact with the wheel assembly; The wheel load threshold at which the wheel approaches the ground; This represents the predicted dynamic wheel load value for the i-th wheel in the future; For dynamic support domain, (·) represents the convex hull operation of the projection point on the horizontal plane. , The first The longitudinal and transverse coordinates of each support point.

[0052] Furthermore, in step S7, the predicted time-domain comprehensive risk index is:

[0053]

[0054] In the formula, To predict the comprehensive risk index in the time domain; , , and Weights for four types of risk items; , , and The first The robust stability margin, lateral load transfer rate, minimum wheel load, and normalized risk terms of instability travel time for each prediction step; For prediction in the time domain.

[0055] Furthermore, in step S8, the formula for calculating the target additional yaw moment is:

[0056]

[0057] In the formula, For the first The target for each prediction step includes an additional yaw moment. For prediction steps; and The first The yaw rate of each prediction step and its reference value; , The first The centroid sideslip angle and roll angular velocity of each prediction step; , The first Robust stability margin and lateral load transfer rate for each prediction step; For the first Slack variables for each prediction step; , , , , , and All are non-negative objective weights.

[0058] A distributed drive heavy-duty truck flexible support center of gravity coupling stability control system, used to implement the distributed drive heavy-duty truck flexible support center of gravity coupling stability control method, including:

[0059] The signal acquisition layer is used to execute step S1, unifying the controller input and calculating the road surface vertical excitation of the wheels;

[0060] The data and basic estimation layer is used to perform step S2, estimate the vehicle's longitudinal velocity, lateral velocity, center of gravity sideslip angle, and vehicle body heave, calculate the tire longitudinal force and tire lateral force, and reconstruct the dynamic wheel load of each wheel.

[0061] The suspension deformation derivation module is used to execute step S3 and estimate the deformation of each suspension connection point and the wheel load consistency residual.

[0062] The dynamic center of gravity reconstruction module is used to execute steps S4 and S5 to estimate the flexible support plane, the dynamic three-dimensional center of gravity and its covariance.

[0063] The predictive stability control module is used to execute steps S6 and S7 and output the future center of gravity trajectory, dynamic support domain and stability state.

[0064] The execution feedback layer is used to execute step S8, output the target additional yaw moment, send commands to the wheel-end motor and brake actuator, and feed back the actual executed torque, braking force and wheel load to the next control cycle.

[0065] Compared with the prior art, the present invention has the following beneficial effects:

[0066] (1) By integrating vehicle body attitude geometry, suspension dynamics and wheel load consistency constraints to solve suspension deformation, the limitations of traditional single detection methods are eliminated. The deformation state of each support node can be obtained with high precision without the need for full coverage suspension travel sensors, which effectively reduces the cost of vehicle sensing hardware and wiring complexity, and improves the accuracy of suspension state identification and adaptability to working conditions.

[0067] (2) By robustly fitting and reconstructing the flexible support plane, and combining wheel load torque balance to achieve accurate dynamic three-dimensional center of gravity solution, it can truly reflect the flexible deformation of the frame and dynamic offset of the center of gravity under heavy load, off-center load and uneven road surface. It overcomes the defects of traditional rigid body and fixed center of gravity model modeling distortion, and significantly improves the accuracy of vehicle dynamics modeling.

[0068] (3) Based on the real-time wheel load update dynamic support domain, the stability margin is robustly corrected by combining the state estimation covariance, and the instability crossover time is solved by predicting the time domain to realize early warning of instability risk. It avoids the optimistic judgment error of the fixed boundary model and the lag problem of traditional instantaneous control intervention, and greatly improves the robustness of instability judgment and control lead.

[0069] (4) The target additional yaw moment is solved by finite-time domain optimization, and a joint distribution strategy of motor and brake is constructed to strictly meet the tire adhesion constraints and the driver's longitudinal power requirements. Through the coordinated logic of motor priority control and brake redundancy compensation, actuator saturation failure is avoided, and the lateral stability and control reliability of the vehicle under extreme conditions are effectively improved. Attached Figure Description

[0070] Figure 1 This is a flowchart of the distributed drive heavy truck flexible support center of gravity coupling stability control method described in this invention.

[0071] Figure 2 This is a schematic diagram of the distributed drive heavy truck flexible support center of gravity coupling stability control system described in this invention. Detailed Implementation

[0072] 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.

[0073] This embodiment defines a unified vehicle coordinate system:

[0074] (1) Longitudinal: The normal straight-ahead direction of a vehicle is longitudinal ( (Axle), with the front of the vehicle pointing in the longitudinal direction; along the longitudinal direction from front to back, it is divided into front, middle, and rear areas, with the front side of the vehicle being the front, the middle section of the vehicle being the middle, and the rear side of the vehicle being the rear.

[0075] (2) Lateral: The left and right directions parallel to the ground and perpendicular to the longitudinal direction are called lateral ( (Axle), with the left side of the vehicle as the lateral direction;

[0076] (3) Vertical: Vertical means perpendicular to the ground and pointing upwards. axis).

[0077] This embodiment uses a 6×6 distributed drive heavy truck as the implementation object. Each of the six wheel ends is equipped with a wheel-end drive motor and a brake actuator, and each wheel has a corresponding suspension connection point. The vehicle is also equipped with sensors (including an inertial measurement unit (IMU), steering angle sensor, wheel speed sensor, brake pressure sensor, and load sensor or air suspension pressure sensor), and a road surface prediction device. The aforementioned sensors, drive motors, and brake actuators are connected to the vehicle controller via an onboard communication network.

[0078] A distributed drive heavy-duty truck flexible support center of gravity coupling stability control method, such as Figure 1 As shown, it includes the following steps:

[0079] In step S1, the unified controller input is used to collect the road surface height sequence and road surface adhesion coefficient, and to calculate the road surface vertical excitation, thereby realizing multi-source signal acquisition.

[0080] In this embodiment, a vehicle coordinate system is established. The axis points forward of the vehicle. The axis points to the left side of the vehicle. The axis is vertically upward. For the aforementioned Wheel arrangement, with six wheels numbered sequentially as follows: , , , , and These represent the front left, front right, middle left, middle right, rear left, and rear right wheel ends, respectively; The planar coordinates of the suspension connection points are: , In the first At each sampling time, the controller input is uniformly represented as:

[0081] (1)

[0082] In the formula, For the first The controller input vector at each sampling time; The longitudinal speed of the vehicle; , , These are the longitudinal acceleration, lateral acceleration, and vertical acceleration in the vehicle coordinate system, respectively. , , These are the roll rate, pitch rate, and yaw rate, respectively. , These are the roll angle and the pitch angle, respectively. This refers to the steering angle; , and The first Wheel speed, wheel-end motor torque, and braking pressure of each wheel; This represents the total number of wheel ends. For sampling time index, Indicates the first The longitudinal distance of each wheel relative to the road surface sensing device.

[0083] Obtain road height sequence through road surface pre-aiming device and road surface adhesion coefficient Based on the vehicle's longitudinal speed and the first The longitudinal distance of each wheel relative to the road surface sensing device is used to convert the road surface height sequence into the corresponding vertical road surface excitation for each wheel. The formula is as follows:

[0084] (2)

[0085] (3)

[0086] In the formula, For the first The road surface excitation propagates to the first sampling time. The time it takes for a wheel to pass For the first The first wheel in the Vertical excitation of the road surface at each sampling time, For the first The longitudinal velocity of the vehicle at each sampling time can be directly obtained using the vehicle state observer. For the first The suspension deformation at each sampling time is as follows: the initial suspension deformation is the original deformation, and the subsequent suspension deformation at each sampling time is derived from the feedback in the subsequent step S3. This is the minimum converted vehicle speed.

[0087] In step S2, based on the multi-source signals from step S1, the basic vehicle state variables such as longitudinal speed, lateral speed, center of gravity sideslip angle, and vehicle body heave are estimated, the tire longitudinal force and tire lateral force are calculated, and the dynamic wheel load of each wheel is reconstructed.

[0088] Vehicle longitudinal velocity, lateral velocity, sideslip angle, vehicle heave, and tire lateral force can be obtained using extended Kalman filtering, unscented Kalman filtering, sliding mode observer, or other mature methods. To avoid excessive complexity in the basic calculations, this implementation will... The longitudinal force of a tire can be simplified as follows:

[0089] (4)

[0090] In the formula, For the first The sampling time of the first sampling moment The longitudinal force of each wheel's tire; For the first The sampling time of the first sampling moment The torque of the motor at the wheel end of each wheel; For the first The effective rolling radius of each wheel; For the first The sampling time of the first sampling moment Braking pressure of each wheel; For the first Conversion factor from braking pressure to braking force for each wheel; For the first The sampling time of the first sampling moment The rolling resistance of each wheel.

[0091] When the lateral force of the tire is directly obtained from the vehicle state observer, it is subject to the tire-road contact adhesion model. Restrictions. Among them, , , They represent the first The estimated values ​​of the longitudinal force, lateral force, and wheel load of the i-th tire at each sampling time.

[0092] The dynamic wheel load of each wheel is not distributed using a single algebraic method. Instead, it is reconstructed using the static wheel load or the wheel load from the previous moment as prior information, and constrained by vertical resultant force, lateral moment, and longitudinal moment. The formula is expressed as:

[0093] (5)

[0094] In the formula, For the first Estimated dynamic wheel load vectors for each wheel at each sampling time; For dynamic wheel load; For the first The prior vector of wheel load at each sampling time; This represents the vertical resultant force equilibrium constraint coefficient matrix; The coefficient matrix representing the balance constraints of longitudinal and lateral moments; and All are positive definite weight matrices.

[0095] in:

[0096] (6)

[0097] In the formula, and These represent the longitudinal and lateral coordinate vectors of each wheel end, respectively. For the first Vehicle mass estimation at each sampling time; It is the acceleration due to gravity; , and The first The longitudinal acceleration, lateral acceleration, and vertical acceleration of the vehicle at each sampling time point; , , They represent the first Estimated values ​​of the longitudinal, lateral, and vertical coordinates of the dynamic three-dimensional centroid at each sampling time. and These are the moments of inertia about the longitudinal axis of the vehicle body and the moments of inertia about the lateral axis, respectively. and The first The vehicle body roll angle acceleration and pitch angle acceleration at each sampling time.

[0098] Equation (5) uses the dynamic three-dimensional centroid of the previous moment (see step S5 for specific calculation) to construct the wheel load constraint at the current moment. Then, based on the wheel load obtained by the solution, the dynamic three-dimensional centroid is updated, and one or more iterations are performed in the current time step until convergence, thereby avoiding the direct algebraic closed loop caused by the mutual dependence of wheel load and centroid.

[0099] In step S3, the suspension deformation variables are derived by fusion. Each suspension connection point is taken as a flexible support node, and the suspension geometric variables are obtained based on the vehicle body heave, roll, pitch and vertical road surface excitation. The suspension geometric variables, suspension spring-damping dynamics and dynamic wheel load consistency relationship are constrained and fused to obtain the optimized suspension deformation variables and wheel load consistency residuals for each suspension connection point.

[0100] Let the generalized coordinates of vehicle body heave, roll, and pitch be:

[0101] (7)

[0102] In the formula, For the first Generalized coordinates of vehicle body heave, roll, and pitch at each sampling moment; For the first Vehicle body lift displacement at each sampling time; and The first The vehicle body roll angle and pitch angle at each sampling time.

[0103] Constructing the suspension geometry mapping matrix , its first Behavior The suspension geometry at each suspension connection point is as follows:

[0104] (8)

[0105] In the formula, The suspension geometric variables are obtained from the vehicle body attitude and the vertical excitation of the road surface; Based on the first step S1, the first... The set of vertical excitations of the road surface at each sampling time.

[0106] Equation (8) reflects the geometric contribution of the overall vehicle body heave, pitch, roll and road surface elevation to the deformation of each suspension.

[0107] On the other hand, if each suspension connection point is equivalent to a spring-damping support, the suspension spring-damping dynamic relationship is as follows:

[0108] (9)

[0109] In the formula, The equivalent mass matrix for each support node; and These are the equivalent damping matrix and stiffness matrix of the suspension, respectively; and These are deformation velocity and deformation acceleration, respectively. This is the dynamic wheel load estimation vector obtained based on step S2. This is a static wheel load; This represents the equivalent disturbance caused by road impact, chassis torsion, and unmodeled factors.

[0110] suspension deformation Using backward differential, the deformation rate can be obtained. And further, the deformation acceleration is obtained. .in, and The first The sampling time and the first sampling time Geometric variables at each sampling time, The sampling period.

[0111] After backward difference and Substituting into equation (9) and further simplifying, we get: Let: coefficient ,

[0112] (10)

[0113] In the formula, No. The wheel load dynamic equivalent input at each sampling time. In Equation (10), the sampling period is used to convert the continuous-time derivative into the discrete-time difference; the continuous-time equivalent stiffness, damping and mass matrix are combined to form the discrete equivalent stiffness matrix after backward difference, and the wheel load and disturbance term are combined to form the wheel load dynamic equivalent input.

[0114] Therefore, the geometric deformation and wheel load dynamics are unified into the following constraint optimization:

[0115] (11)

[0116] st

[0117] In the formula, For the first The suspension deformation variables after fusion and optimization at each sampling time point This represents the suspension deformation to be solved; and These are the geometric prior weight matrix and the wheel load dynamic consistency term weight matrix, respectively; For spatial or interaxial structural continuity matrices; The structural continuity regularization weight coefficient; , These represent the minimum and maximum values ​​of the suspension deformation, respectively. , These represent the minimum and maximum values ​​of the deformation rate, respectively. The typical value range is 0.01 to 10. and The diagonal elements can be 0.1 to 10.

[0118] When the boundary constraints of equation (11) are not activated, its closed-form solution is:

[0119] (12)

[0120] In the formula, To achieve the optimized equivalent information matrix, It is the inverse of the equivalent information matrix.

[0121] (13)

[0122] Equation (13) shows that the suspension deformation is not directly determined by a single formula, but is determined by a weighted sum of geometric information, wheel load information and structural continuity at the previous moment. and The inverse matrix of the geometric estimation covariance and the inverse matrix of the wheel load estimation covariance can be taken separately. When the wheel load is close to zero, the sensor is abnormal, or the road impact is too large, the corresponding weight can be reduced by scaling the diagonal elements of the corresponding weight matrix to avoid abnormal support points dominating the results.

[0123] Wheel load consistency residual The deviation between dynamic wheel load and suspension spring-damping dynamics calculation results is defined as:

[0124] (14)

[0125] In step S4, the flexible support plane is reconstructed. Based on the optimized suspension deformation, road vertical excitation, and initial installation height of the suspension connection points, the spatial height of each support point is restored. A robust weighted method is used to fit the flexible support plane to obtain the estimated parameters of the support plane, the unit normal vector, the equivalent roll angle, and the equivalent pitch angle.

[0126] According to the Initial installation height of each suspension connection point (support point) Suspension deformation and road vertical excitation are used to restore the current support point spatial height.

[0127] (15)

[0128] In the formula, For the first The sampling time of the first sampling moment Current spatial height of each support point; For the first Initial installation height of each support point.

[0129] Let the flexible support plane be:

[0130] (16)

[0131] In the formula, The height of the flexible support plane, For the first Parameters of the flexible support plane at each sampling time; For the first The plane intercept at each sampling time, and The first The slope of the plane along the longitudinal coordinate of the vehicle body and the slope along the lateral coordinate of the vehicle body at each sampling time. , The first The longitudinal and transverse coordinates of each support point.

[0132] Considering that some support points may experience outlier errors due to insufficient wheel load, sensor malfunction, or localized road impact, a robust weighted estimation method is adopted for the first support point. Parameters of the flexible support plane at each sampling time:

[0133] (17)

[0134] In the formula, For the estimation of the plane parameters of the flexible support, For the first The sampling time of the first sampling moment Fitting weights for each support point; The Huber robust loss function; For the first Planar regression vector of each support point; These are the planar time regularized weights. By the The dynamic wheel load, suspension deformation estimation variance and fitting residual at each sampling time are jointly determined and normalized to be between 0 and 1. Typical values ​​range from 0.01 to 10, and are calibrated jointly by plane fitting residuals and attitude change continuity.

[0135] Through iterative reweighted least squares, the th The estimated values ​​of the flexible support plane parameters at each sampling time can be written as:

[0136] (18)

[0137] In the formula, For the first The iterative reweighting matrix at each sampling time point; Regression matrix for each support point It is the identity matrix. For the first The plane intercept estimate at each sampling time. and The first The slope estimates of the plane along the longitudinal direction of the vehicle body and the slope estimates of the lateral direction at each sampling time point.

[0138] The unit normal vector of the flexible support plane is:

[0139] (19)

[0140] In the formula, For the first The unit normal vector of the flexible support plane at each sampling time; , and These represent the unit longitudinal axis component, unit transverse axis component, and unit vertical axis component of the flexible support plane, respectively.

[0141] According to the positive coordinate direction of this embodiment, the equivalent pitch angle and equivalent roll angle can be calculated from the normal vector:

[0142] (20)

[0143] In the formula, For the first The equivalent pitch angle at each sampling time; For the first The equivalent tilt angle at each sampling time; It is the arctangent function with quadrant determination.

[0144] Compared to calculating the parameters of the flexible support plane by dividing the average difference between the left and right suspensions by the wheel track, equations (17) to (20) can simultaneously utilize the six effective support points of a 6×6 wheeled vehicle and take into account the inter-axle attitude differences of the front, middle, and rear axles. When a local wheel is close to lifting off the ground or is subjected to abnormal road impact, the robustness weight can reduce the influence of the corresponding support point on the plane fitting, thereby avoiding the abnormality of a single wheel end dominating the overall support attitude.

[0145] In step S5, the dynamic three-dimensional center of gravity is reconstructed. The static center of gravity is geometrically transformed according to the pose change of the flexible support plane to form a prior geometric center of gravity. The prior geometric center of gravity is then constrained and corrected using the longitudinal and lateral moment balance relationship formed by the dynamic wheel loads of each wheel, thus obtaining the current dynamic three-dimensional center of gravity of the vehicle and its estimated covariance.

[0146] Perform dynamic 3D center of gravity reconstruction. Let the center of gravity obtained from the vehicle's static or low-dynamic calibration be... ,in, , , These are the longitudinal, lateral, and vertical coordinates of the vehicle's static center of gravity, respectively. The translation vector of the flexible support plane is... .

[0147] The support attitude rotation matrix is:

[0148] (twenty one)

[0149] In the formula, For the first The support attitude rotation matrix at each sampling time; and These are the rotation matrices around the transverse axis and the longitudinal axis of the vehicle body, respectively. and For the first The equivalent pitch angle and equivalent roll angle at each sampling time.

[0150] By performing rigid body pose transformations on the flexible support plane using static or low-dynamic calibration of the center of gravity, and adding correction terms caused by the flexible compliance of the frame, suspension, and load, the prior geometric center of gravity can be obtained. :

[0151] (twenty two)

[0152] In the formula, This is the compliance mapping matrix from suspension deformation to equivalent center of gravity shift; The compliance mapping matrix caused by acceleration; For static suspension deformation; and The first The longitudinal and lateral accelerations of the vehicle at each sampling time.

[0153] Based on the longitudinal and lateral moment balance formed by dynamic wheel loads, construct the equivalent gravity term and moment balance constraint vector:

[0154] (twenty three)

[0155] In the formula, For the first The equivalent gravity term after considering vertical inertial acceleration at each sampling time. It is the acceleration due to gravity; For the first The coefficient matrix of wheel load torque balance relationship at each sampling time; for The estimated value of the vehicle mass at each sampling time.

[0156] (twenty four)

[0157] In the formula, For the first The torque balance constraint vector at each sampling moment is composed of the dynamic wheel load of each wheel and the angular acceleration of the vehicle. Based on the result obtained in step 3 The sampling time of the first sampling moment Estimated dynamic wheel load for each wheel; and Let be the moment of inertia of the vehicle body. , These correspond to pitch acceleration and roll acceleration, respectively.

[0158] The final dynamic 3D center of gravity is obtained by fusing the following constraints:

[0159] (25)

[0160] In the formula, For the first The dynamic three-dimensional centroid estimate at each sampling time; The center of gravity is determined by the vehicle's structural dimensions, cargo box range, maximum load capacity, and historical estimates. Feasible region; This is the weight matrix for the wheel load torque balance terms; This is the weight matrix for the prior terms of the geometric centroid. The diagonal elements of the two weight matrices can take values ​​from 0.1 to 10 after input normalization, and the accuracy is jointly calibrated by judging the centroid error, wheel load residuals, and stability.

[0161] If the centroid feasible region constraint in equation (25) is not activated, that is, if the optimal dynamic three-dimensional centroid obtained without considering the boundary is located inside the vehicle boundary, the closed-form solution of equation (25) is:

[0162] (26)

[0163] In the formula, To integrate the information matrix, its inverse matrix gives the uncertainty of the dynamic centroid estimation.

[0164] The covariance of dynamic three-dimensional centroid estimation can be approximated as:

[0165] (27)

[0166] In the formula, For the first The dynamic three-dimensional centroid estimation covariance matrix at each sampling time is used for subsequent zero-moment point covariance propagation and stability boundary confidence contraction.

[0167] In step S6, the future center of gravity trajectory and dynamic support domain are predicted. Considering the driver's future driving input, based on steps S2 to S5, the flexible support plane, dynamic three-dimensional center of gravity, and dynamic support domain formed by the effectively contacting wheels are recalculated in each prediction step within the future prediction time domain.

[0168] Considering the driver's future driving input in the prediction time domain Since the internal dynamics are generally unknown, this embodiment divides future inputs into two categories: driving inputs and external disturbances, and predicts them separately. Future driving inputs include driver steering inputs and driver longitudinal driving force requirements; road vertical excitations are external disturbances and are not included in the driving inputs, but are predicted separately.

[0169] For predictive processing of driver steering input: a first-order hold algorithm is preferentially used to predict the steering angle, and the prediction formula is as follows: ,in, For the first The measured steering angle at each sampling time. Indicates the speed of change of steering angle. This is the prediction step number. Simultaneously, the predicted steering angle is... Upper and lower boundary constraints are set for the steering angle transformation speed; if the system determines that the reliability of the steering angle transformation speed observation result is insufficient, it switches to zero-order hold mode, and the predicted value is fixed to the current measured steering angle. .

[0170] For predicting and processing the driver's longitudinal driving force demand: a zero-order hold mode can be selected, or a first-order linear prediction can be used based on the current rate of change of driving force. The longitudinal driving force amplitude boundary limit is applied throughout the process to avoid the predicted output from exceeding the physical limit of the actuator.

[0171] For predicting the vertical excitation of the road surface in external disturbances When equipped with a road surface pre-aiming device, the road surface height sequence is indexed by longitudinal displacement, and the vertical excitation of the road surface at each wheel time is mapped by the wheelbase time delay.

[0172] Through the above processing, clear data sources of future driving inputs and external disturbances can be provided, eliminating the problem of undefined inputs during prediction and ensuring that the recursive calculation is complete and solvable. Based on steps S2 to S5, the estimation of suspension deformation, support plane parameters, and dynamic three-dimensional center of gravity is re-executed in each prediction step, resulting in:

[0173] (28)

[0174] In the formula, , and These represent the mapping relationships of suspension deformation (obtained from step S3), support plane parameters (obtained from step S4), and dynamic three-dimensional center of gravity (obtained from step S5), respectively. This represents the prediction result of the vehicle state variables at time k for time k+j; This represents the prediction of the dynamic wheel load at time k to time k+j; This represents the prediction result of the suspension deformation at time k for time k+j; This represents the prediction result of the support plane parameters at time k for time k+j; This represents the prediction result of the dynamic three-dimensional centroid at time k for time k+j.

[0175] The dynamic support domain is defined as:

[0176] (29)

[0177] In the formula, To ensure effective contact with the wheel assembly; The wheel load threshold at which the wheel approaches the ground; This represents the predicted dynamic wheel load value for the i-th wheel in the future; For dynamic support domain, (·) represents the convex hull operation of the projection points onto the horizontal plane. For For wheeled vehicles, the dynamic support domain is formed by the effective contact points that meet the wheel load threshold among the six wheel ends: front left, front right, middle left, middle right, rear left, and rear right. When any wheel loses effective support, that wheel end withdraws from the set of wheel ends with effective contact, and the convex hulls of the remaining contact points automatically form a new support boundary.

[0178] Compared with the fixed rectangular support boundary, Equation (29) can more realistically reflect the actual rollover resistance boundary of a three-axle six-wheel vehicle under different wheel ground clearance combinations.

[0179] In step S7, robust stability margin, lateral load transfer rate (LTR), minimum wheel load ratio, and instability crossing time are used to construct a comprehensive risk index for predicting time domain to robustly determine vehicle stability. The vehicle determination results include stable state, warning state, instability state, or rollover risk state.

[0180] The dynamic zero torque point is represented as:

[0181] (30)

[0182] In the formula, The horizontal coordinates of the dynamic zero torque point; For dynamic three-dimensional horizontal projection of the center of gravity; This is the vehicle's horizontal acceleration vector; The rate of change of horizontal angular momentum about the center of gravity; Equivalent mass of the whole vehicle; It is the acceleration due to gravity; The height of the center of gravity; It is a two-dimensional antisymmetric matrix.

[0183] Equation (30) simultaneously considers the position of the center of gravity, translational inertial force, and changes in roll and pitch angular momentum. If the change in angular momentum is ignored, it can degenerate into the commonly used dynamic three-dimensional center of gravity projection expression.

[0184] The dynamic support domain Edges are represented as ,in, Let represent the coordinates of any candidate point within the horizontal support plane. Then, the directed distance from the zero-moment point to this boundary is:

[0185] (31)

[0186] In the formula, and All are dynamic support domains Boundary coefficient of the strip edge; From the zero torque point to the th The directed distance of an edge; This is for nominal stability margin.

[0187] This indicates that the zero moment point is located within the support region. This indicates that the boundary has been reached. This indicates an out-of-bounds error. To avoid overly optimistic judgments due to estimation errors, the dynamic three-dimensional center of gravity and vehicle state covariance are propagated to the zero-moment point covariance. And construct a robust stability margin:

[0188] (32)

[0189] In the formula, To account for the robust stability margin after considering the estimated uncertainty; The covariance matrix at the zero moment point; This is the confidence coefficient, typically ranging from 1.96 to 3.0.

[0190] Lateral load transfer rate and minimum wheel load ratio are defined as follows:

[0191] (33)

[0192] In the formula, For lateral load transfer rate; Minimum wheel load ratio; To prevent positive numbers with a denominator of zero; and The first The dynamic wheel load and the calibrated static wheel load of each wheel.

[0193] Instability travel time is defined as:

[0194] (34)

[0195] In the formula, To predict the robust stability margin in the time domain, the first time it should not be greater than zero or the minimum wheel load. First time not greater than the ground threshold The time; j represents the prediction step number; This represents the robust stability margin obtained by predicting the future time (k+j) at time k. If the above conditions do not occur within the prediction time domain, then the upper limit of the prediction time domain is taken. If the conditions of equation (34) do not occur within the prediction time domain, then... The upper limit of the prediction time domain is uniformly taken.

[0196] To utilize multiple indicators mentioned above simultaneously, a normalized risk term is constructed:

[0197] (35)

[0198] In the formula, , , and These are the normalized risk terms for robust stability margin, lateral load transfer rate, minimum wheel load, and instability travel time, respectively. Width of the support area; For robust stability margin; This is the limit for the lateral load transfer rate; For reference wheel load; Minimum wheel load; This is a time-scale parameter.

[0199] The predicted time-domain comprehensive risk index is:

[0200] (36)

[0201] In the formula, To predict the comprehensive risk index in the time domain; , , and All are risk item weights, all are non-negative and their sum is 1. Initially, they can be 0.4, 0.3, 0.2 and 0.1 respectively. Then, grid search or Bayesian calibration is performed based on the early warning lead time and false alarm rate under typical working conditions. , , and The first Normalized risk terms for robust stability margin, lateral load transfer rate, minimum wheel load, and instability travel time for each prediction step.

[0202] set up :when When it is determined to be a stable state; when When it is determined to be a warning state; when When it is determined to be in an unstable state; when or Or the minimum wheel load is not greater than It was determined to be in a state of risk of rollover.

[0203] In step S8, the target additional yaw moment is calculated and joint allocation is performed. When the vehicle is in a warning state, an instability state, or a rollover risk state, the target additional yaw moment is calculated through finite prediction time domain optimization, and the allocation is completed through joint constraint optimization of wheel-end motor torque and braking force.

[0204] Under steady-state conditions, the controller maintains the driver's original wheel-end torque demand; under warning, instability, or rollover risk conditions, it uses finite prediction time-domain optimization to generate a target additional yaw moment. For the first The additional yaw moment for each control step, the optimization problem is:

[0205] (37)

[0206] st

[0207] In the formula, For the first The target for each prediction step includes an additional yaw moment. For prediction steps; and The first The yaw rate of each prediction step and its reference value; , The first The centroid sideslip angle and roll angular velocity of each prediction step; For the first Slack variables for each prediction step; , , , , , and All are non-negative objective weights. , These are the state vector and the external disturbance, respectively. This is a mapping relationship. Minimum wheel load; , and These represent the minimum, maximum, and maximum change of the additional yaw moment, respectively.

[0208] Equation (37) directly penalizes insufficient robust stability margin and lateral load transfer, in addition to yaw rate tracking and centroid sideslip angle control, and allows temporary infeasibility under extreme conditions to be handled through relaxation variables. Only the first additional yaw moment in the optimal sequence is executed in each control cycle, and the solution is recalculated in the next cycle based on the new state.

[0209] For the combined distribution of the wheel-end motor and brake, let the first... The torque increment of the wheel-end motor of each wheel is Braking force is Then the longitudinal force increment of the wheel is:

[0210] (38)

[0211] In the formula, For the first The longitudinal force increment of each wheel tire; For the first The torque increment of the motor at the wheel end of each wheel; For the first Braking force of each wheel; For the first Effective rolling radius of each wheel.

[0212] Using the wheel-end motor torque increment and braking force as common decision variables, an electromechanical joint allocation optimization model is constructed:

[0213] (39)

[0214] st

[0215] In the formula, This represents the torque increment of the motor at the wheel end. For braking force on the wheels, and These are the cost weight matrices for the motor torque increment and the braking force, respectively. and These are the slack variables for maintaining the yaw moment and longitudinal force constraints, respectively. and All are slack penalty weights; The optimal solution for the additional yaw moment; , These are the lower limit of the wheel-end motor torque and the upper limit of the wheel-end motor torque, respectively. This represents the original wheel-end motor torque; This is the upper limit of the wheel braking force; This refers to the tire adhesion coefficient; , and The first Estimates of longitudinal force, lateral force and dynamic wheel load for each wheel. , After normalization, the diagonal elements can be taken as 0.1 to 10, and are jointly calibrated by actuator energy consumption, additional yaw moment error, and driver longitudinal force holding error; to ensure motor priority, a setting can be made. The corresponding cost is greater than .

[0216] The final motor torque and braking force commands are as follows:

[0217] (40)

[0218] In the formula, and The first Motor torque command and braking pressure command for each wheel. For the first The driver's original motor torque requirement for each wheel; and These are the optimal solutions for the wheel end motor torque increment and the optimal solutions for the wheel braking force obtained from equation (39); For the first The conversion factor from braking pressure to braking force for each wheel.

[0219] like Figure 2 As shown, a distributed drive heavy truck flexible support center of gravity coupling stability control system includes a signal acquisition layer, a data and basic estimation layer, a suspension deformation derivation module, a dynamic center of gravity reconstruction module, a predictive stability control module, and an execution feedback layer.

[0220] The signal acquisition layer is used to execute step S1, which unifies the controller input and calculates the road surface vertical excitation of the wheels.

[0221] The data and basic estimation layer is used to perform step S2, estimating basic vehicle state variables such as vehicle longitudinal velocity, lateral velocity, center of gravity sideslip angle, and vehicle body heave, calculating tire longitudinal force and tire lateral force, and reconstructing the dynamic wheel load of each wheel.

[0222] The suspension deformation derivation module is used to execute step S3, which estimates the deformation of each suspension connection point and the wheel load consistency residual.

[0223] The dynamic center of gravity reconstruction module is used to perform steps S4 and S5 to estimate the flexible support plane, the dynamic three-dimensional center of gravity, and their covariance.

[0224] The predictive stability control module is used to execute steps S6 and S7, and outputs the future center of gravity trajectory, dynamic support domain, and stability state.

[0225] The execution feedback layer is used to execute step S8, output the target additional yaw moment, send commands to the wheel-end motor and brake actuator, and feed back the actual executed torque, braking force and wheel load to the next control cycle.

[0226] 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 distributed drive heavy-duty truck flexible support center of gravity coupling stability control method, characterized in that, Includes the following steps: Step S1: Unify controller input, collect road height sequence and road adhesion coefficient, and calculate road vertical excitation; Step S2: Estimate the vehicle's longitudinal velocity, lateral velocity, center of gravity sideslip angle, and vehicle body heave; calculate the tire longitudinal force and tire lateral force; and reconstruct the dynamic wheel load of each wheel. Step S3: Each suspension connection point is used as a flexible support node. The suspension geometric deformation is obtained based on the vehicle body heave, roll, pitch and road vertical excitation. The suspension geometric deformation, suspension spring-damping dynamics and dynamic wheel load consistency relationship are constrained and fused to obtain the optimized suspension deformation and wheel load consistency residual of each suspension connection point. Step S4: Based on the optimized suspension deformation, road vertical excitation and initial installation height of suspension connection points, restore the spatial height of each support point, and use a robust weighted method to fit the flexible support plane to obtain the parameter estimates, unit normal vector, equivalent roll angle and equivalent pitch angle of the flexible support plane. Step S5: The static center of gravity is geometrically transformed according to the pose change of the flexible support plane to form a prior geometric center of gravity. The prior geometric center of gravity is constrained and corrected by the longitudinal and lateral moment balance relationship formed by the dynamic wheel load of each wheel to obtain the current dynamic three-dimensional center of gravity of the vehicle and its estimated covariance. Step S6: Considering the driver's future driving input, based on steps S2 to S5, recalculate the flexible support plane, dynamic three-dimensional center of gravity, and dynamic support domain formed by the effective contact wheel in each prediction step in the future prediction time domain. Step S7: Use robust stability margin, lateral load transfer rate, minimum wheel load ratio, and instability crossing time to construct a predictive time-domain comprehensive risk index to judge vehicle stability. The vehicle judgment results include stable state, warning state, instability state, or rollover risk state. Step S8: When the vehicle is in a warning state, an unstable state, or a rollover risk state, the target additional yaw moment is calculated by finite prediction time domain optimization, and the execution allocation is completed by joint constraint optimization of wheel-end motor torque and braking force.

2. The method according to claim 1, characterized in that, In step S1, the controller input is uniformly represented as follows: ， In the formula, For the first The controller input vector at each sampling time; The longitudinal speed of the vehicle; , , These are the longitudinal, lateral, and vertical accelerations in the vehicle coordinate system, respectively. , , These are the roll, pitch, and yaw angular velocities, respectively. , These are the roll angle and the pitch angle, respectively. This refers to the steering angle; , and The first Wheel speed, wheel-end motor torque, and braking pressure of each wheel; This represents the total number of wheel ends. Indicates the first The longitudinal distance of each wheel relative to the road surface sensing device; The formula for vertical excitation of the road surface is expressed as: ， ， In the formula, For the first The road surface excitation propagates to the first sampling time. The time it takes for a wheel to pass For the first The first wheel in the Vertical excitation of the road surface at each sampling time, For the first The longitudinal speed of the vehicle at each sampling time. For the first Suspension deformation at each sampling time, To minimize the converted vehicle speed, This is a sequence of road surface heights.

3. The method according to claim 1, characterized in that, In step S2, the dynamic wheel load reconstruction formula for each wheel is as follows: ， In the formula, For the first Estimated dynamic wheel load vectors for each wheel at each sampling time. For dynamic wheel load, For the first The prior load at each sampling time, This represents the vertical resultant force equilibrium constraint coefficient matrix; The coefficient matrix representing the balance constraints of longitudinal and lateral moments. and All are positive definite weight matrices.

4. The method according to claim 1, characterized in that, In step S3, the suspension geometry deformation for: ， In the formula, This is the suspension geometry mapping matrix; Based on the first step S1, the first... The set of vertical excitations of the road surface at each sampling time; For the first Generalized coordinates of vehicle body heave, roll, and pitch at each sampling moment; The suspension spring-damping dynamics relationship is as follows: ， In the formula, The equivalent mass matrix for each support node; and These are the equivalent damping matrix and stiffness matrix of the suspension, respectively; and These are deformation velocity and deformation acceleration, respectively. This is the dynamic wheel load estimation vector obtained based on step S2. This is a static wheel load; For equivalent perturbation; The formula for the optimized suspension deformation at each suspension connection point is expressed as follows: ， In the formula, For the first The suspension deformation variables after fusion and optimization at each sampling time point This represents the suspension deformation vector to be solved; and Both are weight matrices; For spatial or interaxial structural continuity matrices; The structural continuity regularization weight coefficient; No. The equivalent input of wheel load dynamics at each sampling time, For coefficients, For the first Geometric variables at each sampling time.

5. The method according to claim 1, characterized in that, In step S4, the spatial height of the support point is: ， In the formula, For the first The sampling time of the first sampling moment Current spatial height of each support point; For the first Initial installation height of each support point; For vertical excitation of the road surface; The optimized suspension deformation; The estimated parameters of the flexible support plane are: ， In the formula, For the estimation of the plane parameters of the flexible support, For the first Parameters of the flexible support plane at each sampling time. For the first The sampling time of the first sampling moment Fitting weights for each support point; The Huber robust loss function; For the first Planar regression vector of each support point; For planar time regularization weights; This represents the total number of wheel ends.

6. The method according to claim 1, characterized in that, In step S5, the geometric centroid is prior. for: ， In the formula, For the first Translation vector of the flexible support plane at each sampling time; For the first The support attitude rotation matrix at each sampling time; As the center of gravity; This is the compliance mapping matrix from suspension deformation to equivalent center of gravity shift; The compliance mapping matrix caused by acceleration; The optimized suspension deformation; For static suspension deformation; and The first The longitudinal and lateral accelerations of the vehicle at each sampling time; The dynamic three-dimensional center of gravity is: ， In the formula, For the first The dynamic three-dimensional centroid estimate at each sampling time; Center of gravity Feasible region; For the first The coefficient matrix of wheel load torque balance relationship at each sampling time; For the first The torque balance constraint vector at each sampling time; This is the weight matrix for the wheel load torque balance terms; Let be the weight matrix of the prior terms of the geometric centroid.

7. The method according to claim 1, characterized in that, In step S6, the dynamic support domain is defined as: ， In the formula, To ensure effective contact with the wheel assembly; The wheel load threshold at which the wheel approaches the ground; This represents the predicted dynamic wheel load value for the i-th wheel in the future; For dynamic support domain, (·) represents the convex hull operation of the projection point on the horizontal plane. , The first The longitudinal and transverse coordinates of each support point.

8. The method according to claim 1, characterized in that, In step S7, the predicted time-domain comprehensive risk index is: ， In the formula, To predict the comprehensive risk index in the time domain; , , and Weights for four types of risk items; , , and The first The robust stability margin, lateral load transfer rate, minimum wheel load, and normalized risk terms of instability travel time for each prediction step; For prediction in the time domain.

9. The method according to claim 1, characterized in that, In step S8, the formula for calculating the target additional yaw moment is: ， In the formula, For the first The target for each prediction step includes an additional yaw moment. For prediction steps; and The first The yaw rate of each prediction step and its reference value; , The first The centroid sideslip angle and roll angular velocity of each prediction step; , The first Robust stability margin and lateral load transfer rate for each prediction step; For the first Slack variables for each prediction step; , , , , , and All are non-negative objective weights.

10. A distributed drive heavy-duty truck flexible support center of gravity coupling stability control system, used to implement the method according to any one of claims 1 to 9, characterized in that, include: The signal acquisition layer is used to execute step S1, unifying the controller input and calculating the road surface vertical excitation of the wheels; The data and basic estimation layer is used to perform step S2, estimate the vehicle's longitudinal velocity, lateral velocity, center of gravity sideslip angle, and vehicle body heave, calculate the tire longitudinal force and tire lateral force, and reconstruct the dynamic wheel load of each wheel. The suspension deformation derivation module is used to execute step S3 and estimate the deformation of each suspension connection point and the wheel load consistency residual. The dynamic center of gravity reconstruction module is used to execute steps S4 and S5 to estimate the flexible support plane, the dynamic three-dimensional center of gravity and its covariance. The predictive stability control module is used to execute steps S6 and S7 and output the future center of gravity trajectory, dynamic support domain and stability state. The execution feedback layer is used to execute step S8, output the target additional yaw moment, send commands to the wheel-end motor and brake actuator, and feed back the actual executed torque, braking force and wheel load to the next control cycle.

Images