Mechanism Analysis and Control of Lateral Instability of 4WID Vehicle Based on Phase Plane Analysis Considering Front Wheel Angle

The current vehicle stability judgment methods in the stability control of distributed drive electric vehicles focuses on using the phase plane method, but the existing phase plane stability boundary delineation methods ignore the influence of the front wheel rotation angle. Therefore, a new phase plane partitioning method (NPPPM) is proposed in this paper, and a new phase plane stability domain library is established. Based on this, the direct yaw moment control of distributed drive electric vehicle through phase plane analysis is established. Firstly, the influence of the front wheel turning angle on the phase plane is considered in the design of the phase plane stability boundary, and a new stability boundary function is fitted by the least squares method. This approach avoids the lateral instability caused by too large a front-wheel turning angle and improves the accuracy of vehicle stability judgment. Next, a hierarchical direct transverse moment controller is constructed using Matlab/Simulink. Among them, the upper layer controller adopts sliding mode control with the sideslip angle error of the center of mass as the tracking target. The lower controller adopts the optimal configuration algorithm with the optimal road adhesion utilization as the objective function. Finally, the effectiveness and superiority of the proposed phase plane partitioning method (NPPPM) for distributed drive electric vehicles with sinusoidal and angular order inputs for the front wheel rotation angle are verified by simulation.


Introduction
Electric vehicles are currently the focus of scholarly research due to their advantages in terms of energy efficiency and environmental protection. In the development process of distributed drive electric vehicle, the vehicle driven by a four-wheel hub motor plays a milestone role. It can make good use of the drive and braking torque from the wheel hub motors to enhance the control of the vehicle. Meanwhile, the distributed drive electric vehicle eliminates the drive shaft and transmission components, simplifies the chassis structure of the vehicle, improves the space utilization rate of the vehicle, achieves the purpose of being lightweight, and, thus, reduces energy consumption [1][2][3]. At the same time, the distributed hub motor drive vehicle eliminates the complex transmission mechanism, and the driver's operation can be quickly transferred to the hub motor, which improves the control speed and accuracy. However, the distributed drive electric vehicles driven by hub motors, due to the larger mass of hub motor, will lead to increased difficulty in lateral stability control when the vehicle occurs yaw and sideslip, so the distributed drive electric vehicle lateral stability control research has become indispensable in the research of new energy vehicles [4][5][6][7][8].
vehicle. Gong Tianyang et al. [29] designed a fuzzy control method, which takes the errors of yaw rate and sideslip angle as tracking targets, outputs additional yaw moment and rear wheel angle, and uses a genetic algorithm to optimize the membership function of fuzzy control. Yitong Song [30] proposed the estimation method of tire lateral force and vehicle sideslip angle based on an odorless Kalman filter and proposed a hierarchical control strategy where the upper controller is the sliding mode control and the lower controller is the optimal distribution controller, which improves the vehicle handling stability and energy efficiency. However, the nonlinear characteristics of vehicle motion cannot be avoided by using simple sliding mode control. Hongwei Wang et al. [31] proposed a vehicle lateral stability controller based on integral sliding mode control and optimal allocation. By adding an exponential term to the traditional approach rate, the smooth transition of the sliding mode surface is realized, and the effect of buffeting is reduced. Xiaoqiang Sun et al. [32] proposed a method based on non-singular terminal sliding mode control for direct yaw moment control of distributed drive electric vehicles. This method can effectively avoid the nonlinear characteristics of vehicle motion. However, the stability region division of the phase plane does not take into account the influence of the front wheel angle.
In summary, a new phase plane partitioning method (NPPPM) is proposed, and the phase plane stability domain library is established on this basis. In the process of designing the phase plane stable boundary, the stable boundary function about the front wheel angle is designed, which avoids the situation that the phase plane stable boundary changes little when the front wheel angle changes in the past. The lateral stability of the vehicle with sinusoidal input and angular step input is verified by the layered direct yaw moment control. The upper layer is controlled by sliding mode, and the error of the centroid side deflection angle is the tracking target. The lower controller is the optimal allocation and takes the optimal road adhesion utilization rate as the objective function. Finally, the effectiveness of the proposed phase plane division method (NPPPM) is verified by simulation under sinusoidal and angular step input conditions of the front wheel angle.
This paper consists of six parts. The first part introduces the research background and motivation for this paper. The Section 2 introduces the dynamic model of distributed drive electric vehicle. The Section 3 introduces the mechanism of lateral instability of vehicle based on phase plane analysis and the establishment of a phase plane stability domain library. The fourth part is the construction of the layered direct yaw control system. In Section 5, the stability performance of the vehicle under the same working conditions is evaluated by using the NPPPM method and OPPPM method, respectively, and the simulation results are given. Finally, the conclusion is given in the last section.

Distributed Drive Vehicle Seven-Degree-of-Freedom Dynamics Model
In this paper, the body dynamics of the vehicle in three directions, longitudinal, lateral, and yaw, were established using Newton-Euler theory, as shown in Figure 1.
To comply with the assumption of consistent front wheel angles in the ideal model of a two-degree-of-freedom vehicle below, we assume in our study that the left and right front wheel angles are equal in a seven-degree-of-freedom vehicle model. This assumption, known as the front wheel cell assumption, is commonly used in simplified vehicle dynamic models to reduce complexity and computational effort. When the speed of the vehicle in the longitudinal direction is constant, then the equations for the longitudinal, lateral, and yaw movements of the vehicle are as follows: ∑ F x = m . V x − V y γ = F x_ f l + F x_ f r cosδ + F x_rl + F x_rr − F y_ f l + F y_ f r sinδ (1) ∑ F y = m . V y + V x γ = F y_ f l + F y_ f r cosδ + F y_rl + F y_rr + F x_ f l + F x_ f r sinδ (2) In the above equation, m is the overall vehicle mass, I z is the vehicle moment of inertia, γ is the vehicle yaw rate, and δ represents the front wheel angle. In order to conform to the ideal model of a two-degree-of-freedom vehicle presented below, we make the assumption that the turning angles of the left and right front wheels are identical. V x is the vehicle longitudinal velocity, and V y is the vehicle lateral velocity. V is the vehicle center of mass velocity, L is the wheelbase, a is the distance from the centre of mass to the front axle, b is the distance from the center of mass to the rear axle, F x_ij is the longitudinal force at each wheel, and F y_ij is the lateral force per wheel. Vehicle parameters as shown in Table 1.  The mechanical analysis of the wheels during vehicle travel is shown in Figure 2 where the differential equation for the four-wheel moments of a distributed drive electric vehicle can be expressed as follows: where is the rotational moment of inertia of the wheel, is the angular acceleration of the vehicle rotation, is the wheel braking torque, is the wheel driving torque, is the longitudinal force of the tire, and is the tire radius.

Magic Formula Tyre Model
In this paper, the Magic Formula tire model is used as an empirical model to represent the tire forces as a combination of trigonometric functions. It can directly reflect the mathematical relationship between the input variables: side slip angle of the tire, slip rate of the tire, and the output variables: tire longitudinal force; tire lateral force; and return torque. When is the side slip angle, ( ) is sought as the lateral force; when is the slip rate, ( ) is sought as the longitudinal force. The general expressions are as follows: where D is the peak tire curve factor, C is the curve shape factor of the tire, B is the stiffness factor of the tire, and E is the curve curvature factor of the tire. All of the above parameters relate to the vertical load on the tire. Under the effect of lateral acceleration, axial load transfer occurs between the wheels, and the vertical loads on each wheel are as follows:

Magic Formula Tyre Model
In this paper, the Magic Formula tire model is used as an empirical model to represent the tire forces as a combination of trigonometric functions. It can directly reflect the mathematical relationship between the input variables: side slip angle of the tire, slip rate of the tire, and the output variables: tire longitudinal force; tire lateral force; and return torque. When x is the side slip angle, F(x) is sought as the lateral force; when x is the slip rate, F(x) is sought as the longitudinal force. The general expressions are as follows: where D is the peak tire curve factor, C is the curve shape factor of the tire, B is the stiffness factor of the tire, and E is the curve curvature factor of the tire. All of the above parameters relate to the vertical load on the tire. Under the effect of lateral acceleration, axial load transfer occurs between the wheels, and the vertical loads on each wheel are as follows: where h is the distance from the center of mass to the ground, and a y is the lateral acceleration. Side slip angles α i (i = 1, 2, 3, 4) of all four tires are as follows: where L is the distance between the front and rear wheels, respectively.

In-Wheel Motor Model
The parameters of the In-Wheel motor used in this paper are shown in Table 2. The power, torque, and speed of the motor are functionally related as follows, and the external characteristic curve is shown in Figure 3.

T =
T max 0 n n d 9550 P max n n d n n max (8) where T max is the Maximum torque of the motor, P max is the maximum power of the motor, n is the speed of the motor, n d is Motor base speed, and n max is a maximum speed of the motor.

2-DOF Vehicle Reference Model
The vehicle two-degree-of-freedom model is shown in Figure 4. Vehicle two-degreeof-freedom force analysis diagram. including lateral motion and yaw motion. Through the two-degree-of-freedom model, the ideal yaw rate and the sideslip angle of the vehicle can be calculated, and the sideslip angle, the yaw rate phase plane, can be obtained. The dynamic equation of lateral motion and yaw motion can be expressed as follows: when the vehicle is in a stable running state, its centroid side declination angle and yaw angle speed should meet the following formula:

2-DOF Vehicle Reference Model
The vehicle two-degree-of-freedom model is shown in Figure 4. Vehicle two-degreeof-freedom force analysis diagram. including lateral motion and yaw motion. Through the two-degree-of-freedom model, the ideal yaw rate and the sideslip angle of the vehicle can be calculated, and the sideslip angle, the yaw rate phase plane, can be obtained.

2-DOF Vehicle Reference Model
The vehicle two-degree-of-freedom model is shown in Figure 4. Vehicle two-degreeof-freedom force analysis diagram. including lateral motion and yaw motion. Through the two-degree-of-freedom model, the ideal yaw rate and the sideslip angle of the vehicle can be calculated, and the sideslip angle, the yaw rate phase plane, can be obtained. The dynamic equation of lateral motion and yaw motion can be expressed as follows: The dynamic equation of lateral motion and yaw motion can be expressed as follows: when the vehicle is in a stable running state, its centroid side declination angle and yaw angle speed should meet the following formula: .
Thus, the desired yaw rate and sideslip angle can be derived as follows: where δ is the front wheel angle, L is the wheelbase, K is the coefficient of stability, , k f , k r are the lateral stiffness of the front and rear wheels.
Due to the limitations of ground adhesion, the desired yaw rate needs to meet the following: where α is the safety factor (α = 0.85), and g is the acceleration of gravity (g = 9.8).
The final desired yaw angular velocity is: After the sideslip angle of the center of mass is expected to be deformed: Then it is expected that the sideslip angle of the center of mass is limited to Therefore, the final expected sideslip angle of the center of mass is

Establishment of Phase Plane
The phase plane analysis method is based on a two-degree-of-freedom model of the vehicle, as shown in Figure 5. The phase plane of the vehicle under different driving conditions is obtained through multiple cycles of the assignment. The phase plane mainly consists of the sideslip angle-sideslip angular velocity phase plane (β − . β) and sideslip angle-Yaw rate phase plane (β − γ). The differential equation for the system in continuous time is: From Equation (17), we can obtain the following: From Equation (18), we can know that when ( , ) = 0 and ( , ) = 0 are established simultaneously, then ( , ) is the stable point or saddle point of the system, when ( , ) = 0 and ( , ) = 0 are not established simultaneously. As can be seen from Equation (18), the slope of the passing point ( , ) is unchanged at this time, so the phase trajectory through this point is unique. In a continuous time domain system, state initial value, as the time tends to infinity, tends to a constant value, namely, lim → ( ) = , then the system in the local scope gradually reaching stable state points in the system domain. In the phase plane, the phase trajectory of the stable state point will converge at the stable point, while the phase trajectory of the unstable state point will eventually diverge.

Analysis of Influencing Factors of Phase Plane Stability Region
The phase trajectory of distributed drive EV is mainly affected by longitudinal velocity, front wheel angle, and road adhesion coefficient. In this paper, the longitudinal vehicle speed is set at a fixed interval of 10 km/h, ranging from 10 km/h to 50 km/h. Turn the front wheel angle at a fixed interval of 1°, from 0° to 5°. The road adhesion coefficient is fixed at 0.1 intervals, from 0.1 to 1. According to the sample points, 300 groups of sample conditions are obtained. By analyzing the phase plane corresponding to these 300 groups of sample conditions, we can obtain the change rule of the phase plane category and the trend of phase plane change. The differential equation for the system in continuous time is: .
From Equation (17), we can obtain the following: From Equation (18), we can know that when f 1 (x, y) = 0 and f 2 (x, y) = 0 are established simultaneously, then (x, y) is the stable point or saddle point of the system, when f 1 (x, y) = 0 and f 2 (x, y) = 0 are not established simultaneously. As can be seen from Equation (18), the slope of the passing point (x, y) is unchanged at this time, so the phase trajectory through this point is unique.
In a continuous time domain system, state initial value, as the time tends to infinity, tends to a constant value, namely, lim t→∞ x(t) = x e , then the system in the local scope gradually reaching stable state points in the system domain. In the phase plane, the phase trajectory of the stable state point will converge at the stable point, while the phase trajectory of the unstable state point will eventually diverge.

Analysis of Influencing Factors of Phase Plane Stability Region
The phase trajectory of distributed drive EV is mainly affected by longitudinal velocity, front wheel angle, and road adhesion coefficient. In this paper, the longitudinal vehicle speed is set at a fixed interval of 10 km/h, ranging from 10 km/h to 50 km/h. Turn the front wheel angle at a fixed interval of 1 • , from 0 • to 5 • . The road adhesion coefficient is fixed at 0.1 intervals, from 0.1 to 1. According to the sample points, 300 groups of sample conditions are obtained. By analyzing the phase plane corresponding to these 300 groups of sample conditions, we can obtain the change rule of the phase plane category and the trend of phase plane change.

Influence of Longitudinal Speed on the Phase Plane
Longitudinal vehicle speeds were selected as 20, 30, and 40 km/h, respectively, and the front wheel angle was guaranteed to be 0 • , and the road adhesion coefficient was unchanged at 0.8. As shown in Figure 6a-c, the trajectory of the phase tending to stability mainly distributes symmetrically according to the origin in the second and fourth quadrants, while the trajectory of the first and third quadrants mainly diverges. It can also be seen from Figure 6 that when the front wheel angle is 0 • , the type of phase plane will not change. Longitudinal vehicle speeds were selected as 20, 30, and 40 km/h, respectively, and the front wheel angle was guaranteed to be 0°, and the road adhesion coefficient was unchanged at 0.8. As shown in Figure 6a-c, the trajectory of the phase tending to stability mainly distributes symmetrically according to the origin in the second and fourth quadrants, while the trajectory of the first and third quadrants mainly diverges. It can also be seen from Figure 6 that when the front wheel angle is 0°, the type of phase plane will not change.

The Influence of the Road Adhesion Coefficient on the Phase Plane
The road adhesion coefficients were selected as 0.3, 0.6, and 0.9, respectively; the front wheel angle was guaranteed to be 0°, and the longitudinal speed was unchanged at 40 km/h. As shown in Figure 7a-c, the trajectory of the phase tending to stability is mainly distributed symmetrically according to the origin of the second and fourth quadrants, while the trajectory of the first and third quadrants mainly diverges. With the increase in

The Influence of the Road Adhesion Coefficient on the Phase Plane
The road adhesion coefficients were selected as 0.3, 0.6, and 0.9, respectively; the front wheel angle was guaranteed to be 0 • , and the longitudinal speed was unchanged at 40 km/h. As shown in Figure 7a-c, the trajectory of the phase tending to stability is mainly distributed symmetrically according to the origin of the second and fourth quadrants, while the trajectory of the first and third quadrants mainly diverges. With the increase in road adhesion coefficient, the slope of the stable boundary and the stable region increase gradually. Combined with Figure 7, it can also be seen that when the front wheel angle is 0, the type of phase plane will not change.
Actuators 2023, 12, x FOR PEER REVIEW 11 of 25 road adhesion coefficient, the slope of the stable boundary and the stable region increase gradually. Combined with Figure 7, it can also be seen that when the front wheel angle is 0, the type of phase plane will not change.

The Influence of the Front Wheel Angle on the Phase Plane
The front wheel angles were selected as 3°, 4°, and 5°, respectively; the road adhesion coefficient was guaranteed to be 0.5, and the longitudinal speed was constant at 50 km/h. As shown in Figure 8 the stability region of the phase plane becomes the second type of stability region, and most of the stability region is in the second quadrant. With the increase in the front wheel angle, the length BD between the intersection points of the boundary of the stability region gradually increases, but the distance AC between two points, A and C, gradually decreases, and the area of the stability region also gradually decreases.

The Influence of the Front Wheel Angle on the Phase Plane
The front wheel angles were selected as 3 • , 4 • , and 5 • , respectively; the road adhesion coefficient was guaranteed to be 0.5, and the longitudinal speed was constant at 50 km/h. As shown in Figure 8 the stability region of the phase plane becomes the second type of stability region, and most of the stability region is in the second quadrant. With the increase in the front wheel angle, the length BD between the intersection points of the boundary of the stability region gradually increases, but the distance AC between two points, A and C, gradually decreases, and the area of the stability region also gradually decreases.

Judgment of Boundary Type of Stable Region
According to the analysis in Section 3.2, there is a close relationship between the type of stability domain and driving conditions, including front wheel steering angle, speed, and road friction coefficient. When the front wheel angle is 0°, the phase plane category remains unchanged in the linear stability domain. When the front wheel angle changes, the phase plane category will also change. Therefore, this paper divides the phase plane into the linear stability domain and quadrilateral stability domain, as shown in Figure 9.

Judgment of Boundary Type of Stable Region
According to the analysis in Section 3.2, there is a close relationship between the type of stability domain and driving conditions, including front wheel steering angle, speed, and road friction coefficient. When the front wheel angle is 0 • , the phase plane category remains unchanged in the linear stability domain. When the front wheel angle changes, the phase plane category will also change. Therefore, this paper divides the phase plane into the linear stability domain and quadrilateral stability domain, as shown in Figure 9.
Before judging the stability of the current state of the vehicle, it is necessary to determine the corresponding phase plane type of the vehicle at the moment. Through the analysis of the 300 groups of data, the front wheel steering angle and vehicle speed as independent variables, and the road adhesion coefficient as a dependent variable, the three-dimensional graph was drawn. As shown in Figure 9c, when the evaluation value corresponding to the formal condition is located at the top of the three-dimensional surface, the stability domain of the vehicle corresponds to the first type, a linear stability domain. When the corresponding point is located at the bottom of the three-dimensional surface, the vehicle stability domain corresponds to the second type, the quadrilateral stability domain. Before judging the stability of the current state of the vehicle, it is necessary to determine the corresponding phase plane type of the vehicle at the moment. Through the analysis of the 300 groups of data, the front wheel steering angle and vehicle speed as independent variables, and the road adhesion coefficient as a dependent variable, the threedimensional graph was drawn. As shown in Figure 9c, when the evaluation value corresponding to the formal condition is located at the top of the three-dimensional surface, the stability domain of the vehicle corresponds to the first type, a linear stability domain. When the corresponding point is located at the bottom of the three-dimensional surface, the vehicle stability domain corresponds to the second type, the quadrilateral stability domain.

The Boundary Design of the First Type of Stability Region
In the linear stability domain, in order to ensure the accuracy of the design boundary in the control process and the processing of data will not be too complicated. In this paper, linear expression is adopted. As shown in Figure 9a, it can be seen that the two lines are symmetric, with the origin as the center. The right line is L 1 , and the left line is L 2 . The formula is as follows, and the Coefficient table of fitting function is shown in Table 3. The relationship between vehicle state and fitting coefficient is shown in Figure 10: In the linear stability domain, in order to ensure the accuracy of the design boundary in the control process and the processing of data will not be too complicated. In this paper, linear expression is adopted. As shown in Figure 9a, it can be seen that the two lines are symmetric, with the origin as the center. The right line is , and the left line is . The formula is as follows, and the Coefficient table of fitting function is shown in Table 3. The relationship between vehicle state and fitting coefficient is shown in Figure 10:   β B . Let the boundaries be as follows: When the longitudinal speed remains constant, the trend of β A , β D ,  Figure 11: where . β A ,β B are shown as a function of the angle of rotation of the front wheels and the coefficient of road adhesion as follows: . β A = f 5 (δ, µ) = 0.342 − 0.166δ + 0.684µ + 0.017δ 2 − 0.059δµ + 0.904µ 2 β B = f 6 (δ, µ) = −0.94 + 0.207δ − 0.387µ − 0.016δ 2 − 0.06δµ + 0.189µ 2 (22) Actuators 2023, 12, x FOR PEER REVIEW 15 of 25

The Boundary Design of the Second Type of Stability Region
When the front wheel angle changes, the − phase plane type will change too. When the vehicle speed exceeds 30 km/h, and the front wheel angle exceeds 3°, the phase plane type will become quadrilateral stability domain. As shown in Figure 9b, A quadrilateral is constructed with four points of A, B, C, and D, where the horizontal coordinates of A and C are equal, and the vertical coordinates of B and D are equal. Assuming that the coordinates of point A are ( , ), of point C are ( , ), of point B are ( , ), and of point D are ( , ). Therefore, it is sufficient to determine , , , , , . Let the boundaries be as follows: When the longitudinal speed remains constant, the trend of , , , is not significant. However, the value of , change significantly as the front wheel angle and road adhesion coefficient change. So, design the function ( ) to represent , , , . Design ( , ) to represent , , where , , , are related as a function of longitudinal vehicle speed as follows. The relationship between vehicle longitudinal and Coordinate point is shown in Figure 11:

Lateral Stability Controller Design
The front-wheel active steering system (AFS) and the direct yaw moment control system (DYC) are mainly used for vehicle yaw stability control strategies. The AFS uses the lateral force of the tires to generate additional yaw moments, and the system control function fails when the lateral force of the tires is saturated [33][34][35]. Therefore, in this paper, we use DYC to control the stability of the vehicle by using additional transverse moments.
The distributed drive electric vehicle lateral stability controller in this paper uses a hierarchical direct yaw moment control method. First, in order to obtain the optimal additional transverse sway moment, the upper layer controller uses the SMC control algorithm. The additional transverse moment is obtained by tracking the vehicle sideslip angle error. Then, the lower controller uses the optimal distribution algorithm to distribute the torque of the four wheels with the objective of optimal tire slip rate. The control framework is shown in Figure 12. The control framework.

Lateral Stability Controller Design
The front-wheel active steering system (AFS) and the direct yaw moment control system (DYC) are mainly used for vehicle yaw stability control strategies. The AFS uses the lateral force of the tires to generate additional yaw moments, and the system control function fails when the lateral force of the tires is saturated [33][34][35]. Therefore, in this paper, we use DYC to control the stability of the vehicle by using additional transverse moments.
The distributed drive electric vehicle lateral stability controller in this paper uses a hierarchical direct yaw moment control method. First, in order to obtain the optimal additional transverse sway moment, the upper layer controller uses the SMC control algorithm. The additional transverse moment is obtained by tracking the vehicle sideslip angle error. Then, the lower controller uses the optimal distribution algorithm to distribute the torque of the four wheels with the objective of optimal tire slip rate. The control framework is shown in Figure 12. The control framework.

Lateral Stability Controller Design
The front-wheel active steering system (AFS) and the direct yaw moment control system (DYC) are mainly used for vehicle yaw stability control strategies. The AFS uses the lateral force of the tires to generate additional yaw moments, and the system control function fails when the lateral force of the tires is saturated [33][34][35]. Therefore, in this paper, we use DYC to control the stability of the vehicle by using additional transverse moments.
The distributed drive electric vehicle lateral stability controller in this paper uses a hierarchical direct yaw moment control method. First, in order to obtain the optimal additional transverse sway moment, the upper layer controller uses the SMC control algorithm. The additional transverse moment is obtained by tracking the vehicle sideslip angle error. Then, the lower controller uses the optimal distribution algorithm to distribute the torque of the four wheels with the objective of optimal tire slip rate. The control framework is shown in Figure 12. The control framework.

Upper Controller Design
In this paper, the sideslip angle is selected to participate in the design of the sliding mode controller. Due to the addition of an additional yaw moment in the control process, the original two-degree-of-freedom model differential equation can be rewritten as follows: Here ∆M z is the additional yaw moment.

Yaw Moment Calculation Based on Sideslip Angle Error
The tracking error of the vehicle sideslip angle and its derivative are defined as follows: On this basis, the time derivative of sliding mode surface s β can be defined as follows: where c β is the relative weight coefficient between the error and the rate of error change, and its value is greater than 0. The sliding mode approach mode selects the isokinetic approach law can be defined as follows: where K β is the control convergence velocity to the sliding mode surface and is a positive constant. Then: .
The additional yaw moment obtained by tracking the sideslip angle error can be defined as follows:

Controller Stability Verification
When the system state is outside the boundary layer, the control rate and switching logic remain unchanged, but if it is inside the boundary layer, the symbolic function is replaced by the saturation function sat(s).
where H is the boundary layer thickness. In order to reduce the chattering phenomenon of the control system, the function sgn(s) is replaced by the function sat(s) in this paper, so the calculated yaw moment is expressed as follows: Define the Lyapunov function is as follows: Actuators 2023, 12, 121 17 of 23 The first derivative is as follows: .
Substituting ∆M β into the above formula, we can obtain the following: .

Lower Controller Design
In the previous paper, the sideslip angle is selected as the control variable, and the SMC algorithm is applied to control it so as to follow the ideal sideslip angle and solve for the additional yaw moment ∆M that can keep the vehicle laterally stable. The task of the lower controller is to distribute the solved ∆M to the four-wheel hub motors according to a reasonable algorithm in order to ensure the lateral stability of the vehicle.
In this paper, the torque optimal allocation algorithm is used, considering the tire road optimal adhesion utilization as the objective function. The road surface adhesion utilization rate can be used to describe the stability margin of the vehicle. When the road surface adhesion utilization is higher, it means that the distance between the tire force and the adhesion limit is smaller at this time, and the stability is lower. In the distributed drive electric vehicle, all four wheels can generate driving force, so the adhesion utilization rate of all four wheels needs to be considered. Therefore, the objective function used in this paper, which considers the optimal adhesion utilization of the tire road, is as follows: In the optimization process, the objective needs to be achieved by satisfying the equation constraint requirements as the upper force and torque requirements, i.e., satisfying the cross-swing torque and the whole vehicle longitudinal force requirements, while being limited by the road attachment and motor performance. The inequality constraint requirements are the torque constraint as well as the adhesion constraint. The torque constraint is influenced by the maximum driving torque and T max , and the adhesion constraint is mainly influenced by the tire droop load and the ground adhesion coefficient.
According to the above analysis, the constraint function can be written as follows: The solution is as follows:

Simulation Verification
In order to verify the effectiveness of the distributed drive electric vehicle direct transverse moment control system under the phase plane analysis method based on the front wheel rotation angle proposed in this paper. This section conducts Matlab/Simulinkbased simulation to verify it. In this paper, a new phase plane partition method (NPPPM) for phase plane stability domain partitioning and a new stability boundary expression proposed considering the effect of front wheel rotation angle on the phase plane stability boundary are presented. In this section, the proposed NPPPM is compared with the old phase plane partition method (OPPPM) and verified under DYC control at the same time.

Sine Input of the Front Wheel Angle
Set the front wheel angle as a sine input with an amplitude of 0.08 rad, as shown in Figure 13 The friction coefficient of the road surface is set at 0.4 to simulate the wet road surface. Longitudinal speed is set at 50 km/h.

Simulation Verification
In order to verify the effectiveness of the distributed drive electric vehicle direct transverse moment control system under the phase plane analysis method based on the front wheel rotation angle proposed in this paper. This section conducts Matlab/Simulinkbased simulation to verify it. In this paper, a new phase plane partition method (NPPPM) for phase plane stability domain partitioning and a new stability boundary expression proposed considering the effect of front wheel rotation angle on the phase plane stability boundary are presented. In this section, the proposed NPPPM is compared with the old phase plane partition method (OPPPM) and verified under DYC control at the same time.

Sine Input of the Front Wheel Angle
Set the front wheel angle as a sine input with an amplitude of 0.08 rad, as shown in Figure 13 The friction coefficient of the road surface is set at 0.4 to simulate the wet road surface. Longitudinal speed is set at 50 km/h.  As can be seen from Figure 14, the wheel turning angle changes at 3 s. Before that, the vehicle keeps driving straight, and both the sideslip angle and the yaw rate can follow the ideal value. After 3 s, the actual state of the vehicle starts to deviate from the ideal value. The controller with NPPPM responded immediately after 3 s, indicating that the state of the vehicle was outside the stability domain of that phase plane at that time, the vehicle handling performance deteriorated, and the control system intervened, while the controller with OPPPM started to control only after about 3.3 s. It can be inferred that NPPPM can control the vehicle faster when the front wheel angle is changed. Moreover, we can see from the graph that the amplitude of the sideslip angle and yaw rate are 0.0079 (rad) and 0.0307 (rad/s) lower compared with those of the OPPPM method. In the maneuvering stability control system, the magnitude of the sideslip angle should be limited as much as possible to facilitate the stability control of the vehicle. As can be seen from Figure 14, the wheel turning angle changes at 3 s. Before that, the vehicle keeps driving straight, and both the sideslip angle and the yaw rate can follow the ideal value. After 3 s, the actual state of the vehicle starts to deviate from the ideal value. The controller with NPPPM responded immediately after 3 s, indicating that the state of the vehicle was outside the stability domain of that phase plane at that time, the vehicle handling performance deteriorated, and the control system intervened, while the controller with OPPPM started to control only after about 3.3 s. It can be inferred that NPPPM can control the vehicle faster when the front wheel angle is changed. Moreover, we can see from the graph that the amplitude of the sideslip angle and yaw rate are 0.0079 (rad) and 0.0307 (rad/s) lower compared with those of the OPPPM method. In the maneuvering stability control system, the magnitude of the sideslip angle should be limited as much as possible to facilitate the stability control of the vehicle.  Figure 15 shows the simulation comparison of NPPPM and OPPPM under the same working conditions. As can be seen from Figure 15a, the improved phase-plane division method has better lateral control performance than the control strategy without considering the front wheel angle. Figure 15b shows the comparison of additional yaw moment based on the sideslip angle under the two division methods. The additional yaw moment  Figure 15 shows the simulation comparison of NPPPM and OPPPM under the same working conditions. As can be seen from Figure 15a, the improved phase-plane division method has better lateral control performance than the control strategy without considering the front wheel angle. Figure 15b shows the comparison of additional yaw moment based on the sideslip angle under the two division methods. The additional yaw moment under the NPPPM control method is higher than that under the OPPPM control method. The greater the additional yaw torque is applied, which ensures the lateral stability of the vehicle under the ultimate working conditions but reduces the driver's comfort to a certain extent. Under the NPPPM control method, the torque distribution trend of the same side wheel is the same, but the amplitude is different. Front-wheel torque is greater than rear-wheel torque because the vehicle has front-wheel steering.

Step Input of the Front Wheel Angle
Set the front wheel angle as a step input with an amplitude of 0.08 rad, as shown in Figure 16. The friction coefficient of the road surface is set at 0.4 to simulate the wet road surface. Longitudinal speed is set at 50 km/h.

Step Input of the Front Wheel Angle
Set the front wheel angle as a step input with an amplitude of 0.08 rad, as shown in Figure 16. The friction coefficient of the road surface is set at 0.4 to simulate the wet road surface. Longitudinal speed is set at 50 km/h.

Step Input of the Front Wheel Angle
Set the front wheel angle as a step input with an amplitude of 0.08 rad, as shown in Figure 16. The friction coefficient of the road surface is set at 0.4 to simulate the wet road surface. Longitudinal speed is set at 50 km/h.

Step Input of the Front Wheel Angle
Set the front wheel angle as a step input with an amplitude of 0.08 rad, as shown in Figure 16. The friction coefficient of the road surface is set at 0.4 to simulate the wet road surface. Longitudinal speed is set at 50 km/h.

Figure 16.
Step input of the front wheel angle.
At this time, the sideslip angle and yaw rate curve of the distributed drive electric vehicle are shown in Figure 17. As can be seen from the figure, the peak value of the At this time, the sideslip angle and yaw rate curve of the distributed drive electric vehicle are shown in Figure 17. As can be seen from the figure, the peak value of the sideslip angle and yaw rate curve of the controlled centroid is lower than that of the uncontrolled one, and the curve is more stable. At the same time, it can be seen that NPPPM, which takes into account the influence of the front wheel angle on the phase plane stable boundary, has a better control effect than OPPPM under the same conditions. sideslip angle and yaw rate curve of the controlled centroid is lower than that of the uncontrolled one, and the curve is more stable. At the same time, it can be seen that NPPPM, which takes into account the influence of the front wheel angle on the phase plane stable boundary, has a better control effect than OPPPM under the same conditions.  Figure 18 shows the simulation comparison of NPPPM and OPPPM under the same working conditions. Figure 18a shows that the lateral displacement of the trajectory based on the NPPPM control method is small, so the lateral control of the vehicle is slightly better than the OPPPM control method. Figure 18b shows the comparison of the additional yaw moment based on the lateral yaw angle error of the center of mass under the two division methods. The additional yaw moment under the NPPPM control method is higher than that under the OPPPM control method, which leaves less operating space for the driver and reduces the driver's handling comfort to a certain extent but ensures the lateral stability performance of the vehicle under extreme working conditions. Under the NPPPM control method, the torque distribution trend of the same side wheel is the same, but the amplitude is different. Front-wheel torque is greater than rear-wheel torque because the vehicle is front-wheel steering. The phase trajectory based on the NPPPM control method can be closed, and the phase plane is closer to the origin than that under the OPPPM control method, and the area is smaller, which well reflects the superiority of this control method.  Figure 18 shows the simulation comparison of NPPPM and OPPPM under the same working conditions. Figure 18a shows that the lateral displacement of the trajectory based on the NPPPM control method is small, so the lateral control of the vehicle is slightly better than the OPPPM control method. Figure 18b shows the comparison of the additional yaw moment based on the lateral yaw angle error of the center of mass under the two division methods. The additional yaw moment under the NPPPM control method is higher than that under the OPPPM control method, which leaves less operating space for the driver and reduces the driver's handling comfort to a certain extent but ensures the lateral stability performance of the vehicle under extreme working conditions. Under the NPPPM control method, the torque distribution trend of the same side wheel is the same, but the amplitude is different. Front-wheel torque is greater than rear-wheel torque because the vehicle is front-wheel steering. The phase trajectory based on the NPPPM control method can be closed, and the phase plane is closer to the origin than that under the OPPPM control method, and the area is smaller, which well reflects the superiority of this control method. division methods. The additional yaw moment under the NPPPM control method is higher than that under the OPPPM control method, which leaves less operating space for the driver and reduces the driver's handling comfort to a certain extent but ensures the lateral stability performance of the vehicle under extreme working conditions. Under the NPPPM control method, the torque distribution trend of the same side wheel is the same, but the amplitude is different. Front-wheel torque is greater than rear-wheel torque because the vehicle is front-wheel steering. The phase trajectory based on the NPPPM control method can be closed, and the phase plane is closer to the origin than that under the OPPPM control method, and the area is smaller, which well reflects the superiority of this control method.

Conclusions
The innovation of this paper is to propose a new method to divide and represent the phase plane stability region and establish a new phase plane stability region library based on this method. The stability boundary is designed by taking longitudinal vehicle speed, road adhesion coefficient, and front wheel angle as influencing factors. In the process of fitting the stable boundary, the least square method is used to fit 300 groups of stable points. Then, the designed phase plane stability domain library is substituted into the direct yaw moment control of the vehicle to determine the stability of the vehicle. In DYC control, a hierarchical control method is adopted. The upper layer is sliding-mode control,

Conclusions
The innovation of this paper is to propose a new method to divide and represent the phase plane stability region and establish a new phase plane stability region library based on this method. The stability boundary is designed by taking longitudinal vehicle speed, road adhesion coefficient, and front wheel angle as influencing factors. In the process of fitting the stable boundary, the least square method is used to fit 300 groups of stable points. Then, the designed phase plane stability domain library is substituted into the direct yaw moment control of the vehicle to determine the stability of the vehicle. In DYC control, a hierarchical control method is adopted. The upper layer is sliding-mode control, and the error of the centroid side deflection angle is the tracking target. The lower controller is the optimal allocation and takes the optimal road adhesion utilization rate as the objective function.
The NPPPM, OPPPM, and no control conditions with sine input and step input of the front wheel angle are compared. The simulation results show that compared with OPPPM, NPPPM has better stability judgment ability, reduces the stability range of the vehicle, improves the stability performance of the vehicle, and the four-wheel torque distribution is reasonable. Although the control comfort of the driver is reduced to a certain extent, the lateral stability of the vehicle is improved to ensure the safety of the vehicle.