Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles

Abstract

Due to the short transmission chain, compact structure, and the feature of quick and accurate torque generation, distributed drive electric vehicle (DDEV) has attracted many researchers from academia and industry. The significantly redundant execution characteristic of four independently driven in-wheel motors also provides more potential to guarantee the vehicle dynamics performance. Moreover, the unique torque vector control of DDEV generates the direct yaw moment control mode. It has been proven to be effective to modify the vehicle steering characteristics. Through a reasonable torque vector allocation strategy, the energy-saving can also be realized. This chapter will introduce the distributed drive electric vehicle from the viewpoint of the dynamics modeling, stability performance analysis, and energy-saving strategy.

1. Introduction

Electric vehicles (EVs) have been regarded as one of the effective green transportation in urban traffic due to the zero-emission [1, 2, 3]. As a novel chassis structure, distributed drive electric vehicles (DDEVs) choose in-wheel motors as their actuators [4], which is believed to be a promising electric vehicle architecture [5]. DDEVs with multiple powertrains can provide more control schemes through different torque-vector allocation methods. Such could make full use of the tire force limitation to enhance the vehicle handling stability while improving energy efficiency. However, the limited driving range becomes an important factor that restricts the development of EVs in the industry. Extensive research has focused on how to develop the advanced battery technologies, such as the higher energy-density [6] and superfast charging method [7]. Additionally, improving the work efficiency of in-wheel motors can also be an effective approach to reduce the energy consumption. Thanks to the independently controllable motors of DDEVs, it can be achieved by reasonable torque-vector allocation. It should be noted that the yaw motion control generated by the differential torque inputs between left and right wheels may bring about the vehicle instability. Therefore, it would be an interesting study on how to design a torque-vector allocation framework to realize the energy-saving of DDEVs while enhancing the vehicle handling stability.

The vehicle stability control during lateral motion has been a topic of interest in vehicle dynamics control for many years. One of the intuitive approaches to determining stability is the stability region-based method. This method defines the stability region using various vehicle states as indexes, and then identifies the areas in phase planes that correspond to vehicle instability. While several studies have explored the stability regions of centralized drive vehicles (CDVs), few have explored the same for DDEVs. Wang et al. [8] investigated the impact of driving modes on vehicle stability, and Liu et al. [9] studied the effects of driving-steering coupling on lateral stability. However, both studies failed to consider the effects of DYC on stability. Like steering angle, DYC has the potential to significantly alter the flow pattern of vehicle lateral dynamics, leading to changes in the stability regions for different DYC values. Therefore, further research is needed to better understand the effects of DYC on vehicle stability.

There are several methods available to estimate a vehicle’s stability region, which can be categorized as either numerical or analytical. Numerical methods use a mesh on the phase plane to determine the convergence of grid points, and include cell-to-cell mapping, Lyapunov exponent [10], and ARC length methods [11]. However, these methods tend to be time-consuming despite their high accuracy. Analytical methods, on the other hand, aim to find a function to estimate the stability region, with Lyapunov’s second method being a common approach. Unfortunately, this method is often too conservative. Although some attempts have been made to address this issue, the estimation remains conservative. The Sum of Squares Programming (SOSP) method is a polynomial programming technique that can systematically search for the Lyapunov function. By setting constraints on the SOSP, the optimization of the Lyapunov function can be converted into a convex semi-definite program [12], which ensures both complexity and accuracy. Therefore, this study adopted the SOSP method to find the maximal Lyapunov function for stability region estimation.

Furthermore, a reasonable torque vector control can also realize the energy-saving. The in-wheel motors can work in a high-efficiency zone through the torque allocation and reduce the energy consumption. Related research has been conducted. Reference [13] proposes an offline optimization procedure to replace the traditional motor-efficiency mapping method. The simulation results demonstrate that the proposed controller can reduce the motor power loss under different driving conditions while improving the computational efficiency in real applications. Chen et al. [14] discuss and compare the energy-saving results with different energy-efficient control allocation (EECA) schemes. The simulation and experiment results show that Karush-Kuhn-Tuckert (KKT)-based EECA method consumes the least energy, which also has less computational burden [15]. Analytical solutions are derived in [16] for the torque allocation strategy, which aims to reduce the energy loss on the basis of satisfying the total torque demands. Compared with another two allocation methods, the proposed strategy can achieve both energy-saving and computational efficiency.

The vehicle stability control combined with energy-saving is commonly designed through the hierarchical structure. The upper layer includes the total torque inputs and yaw-moment according to the control objectives of longitudinal speed and lateral stability, respectively. The lower layer allocates the torque considering the energy-efficiency. In [17], the top layer develops a DYC to continuously work and guarantee the cornering stability in extreme conditions. The bottom layer designs a switch-rule based on the friction ellipse constraint to judge the control priority for energy-saving and handling stability during the torque allocation. Hua et al. [18] present a hierarchical structure to realize a trade-off between multi-objectives. The higher motion layer aims to generate the desired total torque and yaw moment based on the sliding mode controller. The lower allocation layer uses model predictive control to optimize the motor efficiency. However, the studies in [19, 20, 21] show that the vehicle yaw moment control has the potential for improving energy-efficiency of EVs. The inappropriate yaw-moment control may lead to extra energy consumption. To this end, this work introduces a relaxation factor in the lower layer (i.e. yaw-moment control layer), which aims to reduce the energy consumption of the excessive yaw motion control when the vehicle has enough safety space. The phase plane [1] is adopted to bound the vehicle stability space and designed as a constraint in the MPC controller. The main contributions of this chapter are shown as follows:

1. This study quantifies the stability region of DDEVs using the SOSP technique, expressing stability regions as analytical Lyapunov functions. The results show that DDEVs has a broader stability region than CDVs. An LMI-based mode decision theorem is developed to determine the boundaries of the drive stability regions. This approach constrains control inputs to a safe region based on the concept of drive stability region.

2. A dual LTV-MPC-based hierarchical control framework is constructed to ensure both energy-saving and stability performance of DDEVs. Specifically, it decouples the torque vector control for the motor efficiency optimization and vehicle yaw motion stability control in different control layers. To reduce the energy consumption caused by the excessive direct yaw motion control, a relaxation factor is introduced to balance different control objectives by evaluating the vehicle stability performance in the $\beta -\gamma$ phase plane.

To further prove the energy-saving performance of the proposed controller (PC), the conventional linear model predictive control method (LMPC), and the average torque (AT) allocation method are set as the comparison.

The comparative results with different performance indices are shown in Figure 1. Specifically, the control indices including the energy consumption ${\mathrm{\Delta }}_1={\sum }_{\eta =1}^{\mathrm{\Gamma }}\frac{1}{E_{\max }\mathrm{\Gamma }}\left|E_{\eta }\right|$, the vehicle speed tracking performance ${\mathrm{\Delta }}_2={\sum }_{\eta =1}^{\mathrm{\Gamma }}\frac{1}{e_{v,max}\mathrm{\Gamma }}\left|e_{v,\eta }\right|$, and the tire slip ratio ${\mathrm{\Delta }}_3={\sum }_{\eta =1}^{\mathrm{\Gamma }}\frac{1}{{\lambda }_{t,max}\mathrm{\Gamma }}\left|{\lambda }_{t,\eta }\right|$ are defined, where $E_{\eta }$, $e_{v,\eta }$, and ${\lambda }_{t,\eta }$ denote the energy consumption, the speed tracking error, and the sum of absolute value of tire longitudinal slip ratio at each sampling time, respectively. $\mathrm{\Gamma }$ is the total test time. $E_{\mathit{max}}$, $e_{v,\mathit{max}}$, and ${\lambda }_{t,\mathit{max}}$ correspond to the maximum absolute value. From the results, the proposed controller behaves better to balance different control objectives compared with other methods. On the basis of guaranteeing the driver’s longitudinal speed control intention, the energy-saving control and longitudinal stability performance can be significantly enhanced.

Figure 1.

The performance indices with different control strategies.

2. The lateral stability region of DDEV for state constraint

2.1 DDEV dynamic model

DDEV with wheel-side or hub motor drive have great advantages in energy saving and emission reduction. Wheel-side or Hub motors operate with low noise, high peak efficiency and high load capacity, and also attract much research attention because of their independently controllable torque and fast and accurate torque response, which can effectively improve the vehicle handling stability and safety. Moreover, the DDEV can also realize the differential steering of the vehicle by independently controlling the drive torque difference between the left and right front wheels. It can serve both as a backup system for steering by wire and as the sole steering system of the vehicle, and the latter can further simplify the vehicle structure. DDEV offers flexible chassis layout options, unconstrained by the design limitations of conventional mechanical transmission, and can leverage the benefits of various drive modes.

A DDEV model with front-wheel steering is established [22, 23]. Ignoring the pitch and roll motions, the vehicle has three planar degrees of freedom for longitudinal motion, lateral motion, and yaw motion. A schematic of the vehicle model is shown in Figure 2. According to the principle of balance of forces and moment, the vehicle model in the longitudinal, lateral, and yaw directions can be expressed as:

Figure 2.

Schematic diagram of a vehicle planar motion model.

$$\left\{\begin{array}{l}{\dot{V}}_x=V_y\gamma +\frac{1}{m}\left[\left(F_{x1}+F_{x2}\right)cos\delta +F_{x3}+F_{x4}-\left(F_{y1}+F_{y2}\right)sin\delta \right]\\ {\dot{V}}_y=-V_x\gamma +\frac{1}{m}\left[\left(F_{y1}+F_{y2}\right)cos{\delta }_f+\left(F_{y3}+F_{y4}\right)\right]\\ \dot{\gamma }=\frac{1}{I_z}\left[\left(F_{y1}+F_{y2}\right)cos{\delta }_f-l_r\left(F_{y3}+F_{y4}\right)\right]+\frac{1}{I_z}{M}_c\end{array}\right.$$E1

The external yaw moment $M_c$ is generated with the longitudinal tire force difference between the left and right wheels.

$$M_c=\left(F_{x1}cos{\delta }_f+F_{x3}\right)t_w-\left(F_{x2}cos{\delta }_f+F_{x4}\right)t_w$$E2

When the tire slip angles are small, the front and rear lateral forces can be modeled as:

$$F_{y1}+F_{y2}=C_f{\alpha }_f,F_{y3}+F_{y4}=C_r{\alpha }_r$$E3

And then the tire slip angles can be expressed as:

$$\left\{\begin{array}{l}{\alpha }_f={\delta }_f-\frac{V_y+\gamma l_f}{V_x}\\ {\alpha }_r=\frac{\gamma l_r-V_y}{V_x}\end{array}\right.$$E4

The rotational dynamics of each wheel can be represented by

$$J_{\mathit{wi}}{\dot{w}}_i=-R_e{F}_{\mathit{xi}}+T_{\mathit{wi}}$$E5

And then the longitudinal tire force at each tire can be rewritten as:

$$F_{\mathit{xi}}=\frac{1}{R_e}\left(T_{\mathit{wi}}-J_{\mathit{wi}}{\dot{w}}_i\right)$$E6

In summary, the overall vehicle model (1) can be rewritten as:

$$\left\{\begin{array}{l}{\dot{V}}_x=V_y\gamma +\frac{1}{{\mathit{mR}}_{\mathrm{e}}}\left(u_1+u_2\right)+d_1\\ {\dot{V}}_y=-\frac{\left(C_f+C_r\right)V_y}{{\mathit{mV}}_x}+\left(\frac{C_r{l}_r-C_f{l}_f}{{\mathit{mV}}_x}-V_x\right)\gamma +\frac{C_f}{m}{u}_3+d_2\\ \dot{\gamma }=\frac{\left(C_r{l}_r-C_f{l}_f\right)V_y-\left(C_f{l}_f^2+C_r{l}_r^2\right)\gamma }{I_z{V}_x}+\frac{C_f{l}_f}{I_z}{u}_3+\frac{l_s}{I_z{R}_{\mathrm{eff}}}\left(u_1-u_2\right)+d_3\end{array}\right.$$E7

where $u_1$ and $u_2$ represent the total motor torque values in the vehicle longitudinal direction on the left and right sides of the vehicle, respectively. $u_3$ represents the front wheel steering angle.

$$\left\{\begin{array}{l}{u}_1=T_{w1}cos\delta +T_{w3}\\ u_2=T_{w2}cos\delta +T_{w4}\\ u_3=\delta \end{array}\right.$$E8

$$\left\{\begin{array}{l}{d}_1=\frac{1}{m}\left(F_{y1}+F_{y2}\right)sin\delta +\frac{1}{{\mathit{mR}}_e}\left(J_{w1}{\dot{w}}_{1}cos\delta +J_{w2}{\dot{w}}_{2}cos\delta +J_{w3}{\dot{w}}_3+J_{w4}{\dot{w}}_4\right)\\ d_2=\frac{1}{m}\left(F_{y1}+F_{y2}\right)\left(cos\delta -1\right)\\ d_3=\frac{l_f}{I_z}\left(F_{y1}+F_{y2}\right)\left(cos\delta -1\right)+\frac{l_s}{I_z{R}_e}\left(J_{w1}{\dot{w}}_{1}cos\delta -J_{w2}{\dot{w}}_{2}cos\delta +J_{w3}{\dot{w}}_3-J_{w4}{\dot{w}}_4\right)\end{array}\right.$$E9

2.2 Rational polynomial based DDEV model

To accurately estimate the stability region, a nonlinear tire model is necessary as the fixed cornering stiffness used in the tire model is insufficient for extreme conditions. Here, we propose to use a rational polynomial function to fit the tire force curve. We define the nominal lateral tire force in Eq. (10) without loss of generality.

$$f_y\left({\alpha }_i\right)=\frac{p_1{{\alpha }_i}^3+p_2{\alpha }_i}{q_1{\alpha }_i^4+q_2{\alpha }_i^2+q_3}$$E10

We choose an odd function for the tire model 10 to cover both positive and negative tire forces. By incorporating the vertical load $F_z$ and friction coefficient $\mu$ into the model, we derive a rational polynomial tire model expressed as Eq. (11).

$$F_{\mathit{yi}}=\frac{F_{\mathit{zi}}{\mu }_i\left(p_1{\left({\alpha }_i/{\mu }_i\right)}^3+p_2\left({\alpha }_i/{\mu }_i\right)\right)}{q_1{\left({\alpha }_i/{\mu }_i\right)}^4+q_2{\left({\alpha }_i/{\mu }_i\right)}^2+q_3}$$E11

After conducting experiments with a 6000 $N$ vertical load, Figure 3 was produced to display the corresponding test data. A least-squares algorithm was then utilized to fit the polynomial coefficients in Eq. (11). Then, the rational polynomial tire model (11) will be employed in the dynamic model of DDEV. For the convenience of deduction, Eq. (11) is recorded as

Figure 3.

Stability region estimation of different degrees.

$$F_{\mathit{yi}}=n_i\left({\alpha }_i\right)/d_i\left({\alpha }_i\right),\left(i=f,r\right)$$E12

According to the kinematic characteristic, the tire slip angle $\alpha$ can be represented by the vehicle state variables.

$$\begin{array}{l}\left\{\begin{array}{c}{\alpha }_f={\delta }_f-\beta -\frac{l_f}{V_x}\gamma \\ {\alpha }_r=\frac{l_r}{V_x}\gamma -\beta \end{array}\right.\\ s.t.\beta =V_y/V_x\end{array}$$E13

The lateral dynamic model of DDEV can be transformed into the standard state space form (14). This results in the rational polynomial dynamic model for stability region estimation.

$$\dot{x}=f\left(x,u\right)=\frac{N\left(x,u\right)}{D\left(x,u\right)}$$E14

2.3 Basic principles of SOSP

The Lyapunov method is commonly utilized to estimate the stability region by seeking the maximal Lyapunov function that approaches the RoA boundary. Although systematic methods for searching the Lyapunov function are scarce, the stability region estimation can be converted into a convex optimization problem by SOSP when dealing with Lyapunov functions in polynomial form. This section introduces the fundamental principles of SOSP.

2.3.1 Sum of Squares Programming

Definition 1 Consider a polynomial function $p\left(x\right)$, $p\left(x\right)$ with $n$ real variables and $m$ degrees. $p\left(x\right)$ is called sum of squares (SOS) if there exist polynomials $f_i\left(x\right)$ such that

$$p\left(x\right)=\sum _{i=1}{f}_i^2\left(x\right)$$E15

Lemma 1 (Quadratic form of polynomial.) For a polynomial $p\left(x\right)$, there definitely exist a symmetric matrix $Q$ such that

$$p\left(x\right)=z^T\left(x\right)\mathit{Qz}\left(x\right)$$E16

where $z\left(x\right)$ is a vector of all monomials of degree less than or equal to $\frac{m}{2}$. Commonly, the matrix $Q$ is not unique. The matrix space of $Q$ could be represented as a function of ${\lambda }_i$ (17).

$$Q\left(\lambda \right)=Q_0+\sum _{i=1}^N{\lambda }_i{M}_i$$E17

Lemma 2 (Sum of Squares Programming.) For a polynomial $p_i\left(x\right)$ is SOS if and only if there exist ${\lambda }_i,i=1,\dots ,N$ such that

$$\begin{array}{l}Q\left(\lambda \right)=Q_0+\sum _{i=1}^N{\lambda }_i{M}_i\ge 0\\ s.t.p\left(x\right)=z^T\left(x\right)Q\left(\lambda \right)z\left(x\right)\end{array}$$E18

It is worth noting that Lemma 2 can be formulated as a linear matrix inequality (LMI) feasibility problem. Moreover, the sum-of-squares (SOS) property of a polynomial is equivalent to the SDP of the corresponding matrix $Q$. Therefore, the SDP approach provides an effective way to solve the SOSP problem.

2.3.2 Generalized S-procedure

The stability region estimation is usually concerned with causality of multiple inequalities. For the convenience of solving, it should be integrated as a single LMI.

Lemma 3 (Generalized S-procedure.) Consider a series of polynomials $p\left(x\right),$ $i=0,\dots ,m$ such that

$$\begin{array}{l}{p}_0\left(x\right)\ge 0\\ s.t.x\in D\end{array}$$E19

where $D$ is domain of $x$ represented as:

$$p_1\left(x\right)\ge 0,\dots ,p_m\left(x\right)\ge 0$$E20

The inequality (19) hold if there exist $q_i\left(x\right),i=1,\dots ,m$ such that

$$\begin{array}{l}{p}_0\left(x\right)-\sum _{i=1}^m{q}_i\left(x\right)p_i\left(x\right)\ge 0\\ s.t.x\in R^n\end{array}$$E21

2.4 Stability region estimation and analysis

Lemma 4 (Invariant subset of RoA.) Consider a function $V$ and $\gamma >0$. Region ${\mathrm{\Omega }}_{V,\gamma }$ is defined as $\left\{x\in R^n:V\left(x\right)\le \gamma \right\}$. If conditions in (22) holds, ${\mathrm{\Omega }}_{V,\gamma }$ is an invariant subset of RoA.

$$\begin{array}{l}{\mathrm{\Omega }}_{V,\gamma }\mathit{\text{is bounded}}\\ V\left(0\right)=0,V\left(x\right)>0\left(\forall x\in R^n\right)\\ {\mathrm{\Omega }}_{V,\gamma }\backslash \left\{0\right\}\subset \left\{x\in R^n:\nabla V\left(x\right)f\left(x\right)\le 0\right\}\end{array}$$E22

To construct an SOSP problem, the Lyapunov function $V$ should be restricted as a polynomial form. Besides, formula ③ can be converted to a SOSP problem according to generalized S-procedure. Except of this, Lemma 4 is not sufficient to find the maximal Lyapunov function. Thus, we set a shape function $s\left(x\right)$ to expand the stability region of DDEV (23). By maximize $\beta$, $V$ tends to approach the maximal Lyapunov function.

$$max\beta :\left\{x\in R|s\left(x\right)<\beta \right\}\subset {\mathrm{\Omega }}_{V,\gamma }$$E23

Theorem 1 For the nonlinear dynamic system, the stability region ${\mathrm{\Omega }}_{V,\gamma }$ can be found by finding $V$, $q_1$, $q_1\in \sum _{n,m}\left[x\right]$ and positive $\gamma >0$ that maximize $\beta$ such that

$$\begin{array}{c}max\beta \\ \left\{\begin{array}{l}V\left(x\right)-{\phi }_1\left(x\right)\in \sum \left[x\right]\\ q_1\left(x\right)\left(s-\beta \right)-\left(V-\gamma \right)\in \sum \left[x\right]\\ -{\nabla }_x\mathit{VN}-D{\phi }_2\left(x\right)+{\mathit{Dq}}_2\left(x\right)\left(V-\gamma \right)\in \sum \left[x\right]\end{array}\right.\end{array}$$E24

The degree of Lyapunov function $V$ has a huge impact on the estimation of stability region. To select the appropriate degree of $V$, we compare the stability regions with different degrees on the phase plane. Figure 3 shows the stability regions of different degrees.

Figure 3 illustrates the phase trajectories near the equilibrium point of straight-running, where the blue lines represent the region of attraction (RoA) and the remaining area indicates instability. Generally, the stability region expands as the degree increases. However, the RoA of the two-degree $V$ only approaches the boundary of yaw rate due to the ellipse shape, which cannot accurately describe the margin of slip angle. In contrast, the six and eight degree $V$ significantly increase the stability region in all directions and better capture the RoA feature around the equilibrium. Given the computational complexity, we opt for the six-degree $V$ to estimate the stability region of DDEV.

2.4.1 Impact of longitudinal velocity and road adhesion

During straight-running conditions, steering angle and DYC control inputs are both set to zero. Longitudinal velocity $V_x$ and road adhesion coefficient $\mu$ are two primary factors that influence lateral stability. To explore their effects on the stability region, we varied $V_x$ and $\mu$ and plotted the estimated stability regions in phase portraits, as shown in Figures 4–7. Our results indicate that as $V_x$ increases, the stability region in the yaw rate direction tends to shrink, which is consistent with the view that high-speed steering can cause vehicle instability. Similarly, when $\mu$ decreases, not only does the available yaw rate reduce, but the stability region in the slip angle direction also sharply decreases due to the restriction of tire adhesion margin from low road friction. As a consequence, the maximum wheel slip angle decreases, leading to a reduction in body slip angle. However, at higher speeds, the tire adhesion margin remains constant, and can still supply sufficient tire force for the vehicle’s lateral motion. These findings suggest that stability control in low road friction conditions may pose a greater challenge than in high-speed running. With these insights, we can better understand the factors influencing lateral stability, and design more effective control strategies for safer vehicle operation.

Figure 4.

$V_x=\mathit{20}m/s\left(\mu =\mathit{1}\right)$.

Figure 5.

$V_x=\mathit{30}m/s\left(\mu =\mathit{1}\right)$.

Figure 6.

$\mu =\mathit{1}\left(V_x=\mathit{20}m/s\right)$.

Figure 7.

$\mu =\mathit{0.5}\left(V_x=\mathit{20}m/s\right)$.

2.4.2 Impact of steering angle and DYC

Compared to straight-running conditions, the estimation of cornering conditions is more complex due to the non-zero control inputs ${\delta }_f$ and $M_z$, which cause the equilibrium to shift away from the origin. To apply Theorem 2, we must first solve for the new equilibrium and then substitute the state variables to transform the equilibrium to the origin, enabling the estimation of the stability region for cornering conditions. For the CDV, the only control input is the steering angle. Figures 8 and 9 illustrate the stability regions for a steering angle of $\pm 12°$, showing that the stability region in the yaw rate direction shrinks in the same direction as the steering angle. This suggests that further increases in steering angle may result in vehicle instability. Furthermore, the stability region for cornering conditions is much narrower than that for straight-running, making the vehicle more susceptible to external disturbances. These findings highlight the importance of understanding the effects of control inputs on vehicle stability and developing effective control strategies to ensure safe and stable operation.

Figure 8.

$M_z=\mathit{0}$ nm $\left(\delta =\mathit{12}°\right)$_.

Figure 9.

$M_z=\mathit{0}$ nm $\left(\delta =-\mathit{12}°\right)$_.

To improve lateral stability, DDEV generates DYC by distributing torque unbalance, affecting the vehicle’s lateral dynamics. However, DYC cannot be arbitrarily imposed. To investigate the impact of DYC on the stability region, we apply different values of $M_z$ in the same or opposite directions of the steering angle, as shown in Figures 10 and 11. For left-turning conditions, a 300 $\mathit{Nm}$ $M_z$ significantly reduces the stability region of the yaw rate. However, an opposite −800 $\mathit{Nm}$ $M_z$ greatly expands the stability region. Compared to 0 and 300 $\mathit{Nm}$ $M_z$, a much higher available yaw rate is achievable with an opposite DYC, allowing the vehicle to withstand greater yaw motion. Thus, we conclude that DDEV can enhance lateral stability by applying an opposite DYC. Simulations in the next section will confirm this hypothesis (Figures 12 and 13).

Figure 10.

$M_z=\mathit{300}$ nm $\left(\delta =\mathit{12}°\right)$_.

Figure 11.

$M_z=-\mathit{300}$ nm $\left(\delta =-\mathit{12}°\right)$_.

Figure 12.

$M_z=-\mathit{800}$ nm $\left(\delta =\mathit{12}°\right)$_.

Figure 13.

$M_z=\mathit{800}$ nm $\left(\delta =-\mathit{12}°\right)$_.

3. The drive stability region of DDEV for control constraint

3.1 The definition of drive stability region

In a hierarchical DYC control scheme, the lateral stability controller at the high level calculates the efficient DYC $M_{\mathit{\text{zopt}}}$, while the driver’s actions provide the total traction force Fxall. To transform $F_{\mathit{\text{xall}}}$ and $M_{\mathit{\text{zopt}}}$ into four-wheel torques, a rational distribution algorithm is necessary. A number of approaches consider $F_{\mathit{\text{xall}}}$ and $M_{\mathit{\text{zopt}}}$ as a strict constraint and employ optimization techniques to solve the problem. However, the excessive values of $F_{\mathit{\text{xopt}}}$ and $M_{\mathit{\text{zopt}}}$ could surpass the adhesion limit during handling limit. The difficulty of distribution control lies in whether to satisfy $F_{\mathit{\text{xall}}}$ and $M_{\mathit{\text{zopt}}}$ or only satisfy one of them. Here we divide them into several situations. By taking the traction force and DYC as the two axes of a plane, we can divide the plane into different areas and make reasonable inferences. First, there is a maximum region $D_{\mathit{max}1}$ in which both $F_{\mathit{\text{xall}}}$ and $M_{\mathit{\text{zopt}}}$ can be fulfilled without sacrificing stability. Second, if we disregard the $F_x$ requirement, we can obtain a larger region $D_{\mathit{max}2}$ where $M_{\mathit{\text{zopt}}}$ can be met without losing stability. For the areas beyond $D_{\mathit{max}2}$, it is not possible to fully comply with either $F_{\mathit{\text{xall}}}$ or $M_{\mathit{\text{zopt}}}$ without compromising stability. As $F_{\mathit{\text{xall}}}$ and $M_{\mathit{\text{zopt}}}$ are linked to the driver’s maneuvers, these feasible regions are referred to as drive stability regions. Figure 14 provides an intuitive depiction of the drive stability region. According to the different drive stability regions, we obtain the definition of three Torque distribution Methods (TDMS).

Figure 14.

The schematic diagram of drive stability region.

3.2 The construction of drive stability region

In order to estimate the drive stability region, a suitable vehicle dynamic model is required. As the drive stability region is related to the longitudinal and lateral dynamics of the vehicle, a three-degree-of-freedom vehicle plane motion model is developed for analysis.

$$\left\{\begin{array}{l}m\left({\dot{v}}_x-\gamma V_y\right)=F_{\mathit{xfl}}+F_{\mathit{xfr}}+F_{\mathit{xrl}}+F_{\mathit{xrr}}\\ m\left({\dot{v}}_y+\gamma V_x\right)=F_{\mathit{yfl}}+F_{\mathit{yfr}}+F_{\mathit{yrl}}+F_{\mathit{yrr}}\\ I_z\dot{\gamma }=l_f\left(F_{\mathit{yfl}}+F_{\mathit{yfr}}\right)-l_r\left(F_{\mathit{yrl}}+F_{\mathit{yrr}}\right)+t_w\left(F_{\mathit{xfr}}-F_{\mathit{xfl}}+F_{\mathit{xrr}}-F_{\mathit{xrl}}\right)/2\end{array}\right.$$E25

where $t_w\left(F_{\mathit{xfr}}-F_{\mathit{xfl}}+F_{\mathit{xrr}}-F_{\mathit{xrl}}\right)/2$ is the actual DYC $M_{\mathit{zc}}$. As the longitudinal tire force increases, the lateral tire force decreases. The maximum longitudinal and lateral tire force envelope forms an ellipse known as the tire adhesion ellipse. This adhesion ellipse is utilized to establish the drive stability region, which is presented in Eq. (26).

$${\left(\frac{F_{\mathit{xij}}}{{\mu }_x}\right)}^2+{\left(\frac{F_{\mathit{yij}}}{{\mu }_y}\right)}^2\le {F_{\mathit{zij}}}^2$$E26

Certainly, the lateral and vertical tire force mentioned in Eq. (25) cannot be directly obtained and must be estimated using the three-freedom vehicle plane model. Accounting for load transfer, the vertical tire force can be calculated using Eq. (27), where $h_g$ represents the height of the center of mass and $W$ is the wheelbase.

$$\begin{array}{l}{F}_{\mathit{zfl}}=\frac{{\mathit{mgl}}_r}{2L}-\frac{m{\dot{v}}_x{h}_g}{2L}-\frac{m{\dot{v}}_y{h}_g{l}_r}{t_{w}L}{F}_{\mathit{zfr}}=\frac{{\mathit{mgl}}_r}{2L}-\frac{m{\dot{v}}_x{h}_g}{2L}+\frac{m{\dot{v}}_y{h}_g{l}_r}{t_{w}L}\\ F_{\mathit{zrl}}=\frac{{\mathit{mgl}}_f}{2L}+\frac{m{\dot{v}}_x{h}_g}{2L}-\frac{m{\dot{v}}_y{h}_g{l}_f}{t_{w}L}{F}_{\mathit{zrr}}=\frac{{\mathit{mgl}}_f}{2L}+\frac{m{\dot{v}}_x{h}_g}{2L}+\frac{m{\dot{v}}_y{h}_g{l}_f}{t_{w}L}\end{array}$$E27

The front and rear lateral tire force is represented as:

$$\left\{\begin{array}{l}{F}_{\mathit{yf}}=F_{\mathit{yfl}}+F_{\mathit{yfr}}=\left(l_{r}m\left({\dot{v}}_y+\gamma V_x\right)+I_z\dot{\gamma }+M_{\mathit{zc}}\right)/L\\ F_{\mathit{yr}}=F_{\mathit{yrl}}+F_{\mathit{yrr}}=\left(l_{f}m\left({\dot{v}}_y+\gamma V_x\right)-I_z\dot{\gamma }-M_{\mathit{zc}}\right)/L\end{array}\right.$$E28

Here, we assume that the lateral force transfer is similar to the vertical force transfer, and define the load transfer coefficients as follows.

$$\begin{array}{l}{k}_{\mathit{fl}}=F_{\mathit{zfl}}/\left(F_{\mathit{zfl}}+F_{\mathit{zfr}}\right)k_{\mathit{fr}}=F_{\mathit{zfr}}/\left(F_{\mathit{zfl}}+F_{\mathit{zfr}}\right)\\ k_{\mathit{rl}}=F_{\mathit{zrl}}/\left(F_{\mathit{zrl}}+F_{\mathit{zrr}}\right)k_{\mathit{rr}}=F_{\mathit{zrr}}/\left(F_{\mathit{zrl}}+F_{\mathit{zrr}}\right)\end{array}$$E29

Combined with Eq. (28), the lateral tire force of each wheel can be calculated, which is represented in Eq. (30).

$$\left\{\begin{array}{l}{F}_{\mathit{yfi}}=k_{\mathit{fi}}\left(l_{r}m\left({\dot{v}}_y+\gamma V_x\right)+I_z\dot{\gamma }+M_{\mathit{zc}}\right)/L\\ F_{\mathit{yri}}=k_{\mathit{ri}}\left(l_{f}m\left({\dot{v}}_y+\gamma V_x\right)-I_z\dot{\gamma }-M_{\mathit{zc}}\right)/L\end{array}\right.$$E30

The longitudinal and lateral tire forces can be estimated using the state variables of the DDEV model. With this information, a linear matrix inequality-based mode decision theorem is formulated below.

Lemma 5 (LMI based conditions of TDM 1) TDM 1 is satisfied for the current vehicle condition, given the total traction force Fxall and optimal DYC $M_{\mathit{\text{zopt}}}$, if and only if the linear matrix inequality (LMI) shown in Eq. (27) is solvable, where $X={\left[F_{\mathit{\text{xall}}}{M}_{\mathit{\text{zopt}}}\right]}^T,j=l,r$.

$$\begin{array}{l}{A}_{m1}X\ge b_{m1}\\ s.t.A_{m1}={\left[\begin{array}{cccc}-\frac{1}{2}& \frac{1}{2}& -\frac{1}{2}& \frac{1}{2}\\ 1& -\frac{1}{t_w}& -1& \frac{1}{t_w}\end{array}\right]}^T\\ b_{m1}={\left[\begin{array}{cccc}-A_{\mathit{fl}}-A_{\mathit{rl}}& -A_{\mathit{rl}}-A_{\mathit{fl}}& -A_{\mathit{fr}}-A_{\mathit{rr}}& -A_{\mathit{rr}}-A_{\mathit{fr}}\end{array}\right]}^T\\ A_{\mathit{fj}}=\sqrt{{\left({\mu }_x{F}_{\mathit{zfj}}\right)}^2-{\left(k_{\mathit{fj}}{\mu }_x\left(l_{r}m\left({\dot{v}}_y+\gamma V_x\right)+I_z\dot{\gamma }+M_{\mathit{zc}}\right)/{\mu }_{y}L\right)}^2}\\ A_{\mathit{rj}}=\sqrt{{\left({\mu }_x{F}_{\mathit{zrj}}\right)}^2-{\left(k_{\mathit{rj}}{\mu }_x\left(l_{f}m\left({\dot{v}}_y+\gamma V_x\right)-I_z\dot{\gamma }-M_{\mathit{zc}}\right)/{\mu }_{y}L\right)}^2}\end{array}$$E31

Lemma 6 (LMI based conditions of TDM 2)

$$\begin{array}{l}{A}_{m2}X\ge b_{m2}\\ s.t.A_{m2}={\left[\begin{array}{cc}-\frac{2}{t_w}& \frac{2}{t_w}\end{array}\right]}^T\\ b_{m2}={\left[\begin{array}{cc}-3A_{\mathit{fl}}-3A_{\mathit{fr}}-A_{\mathit{rl}}-A_{\mathit{rr}}& -3A_{\mathit{fl}}-3A_{\mathit{fr}}-A_{\mathit{rl}}-A_{\mathit{rr}}\end{array}\right]}^T\end{array}$$E32

Lemma 7 (LMI based conditions of TDM 3) For the given total traction force $F_{\mathit{\text{xall}}}$ and optimal DYC $M_{\mathit{\text{zopt}}}$, TDM3 is satisfied, if and only if (24) and (25) is unsolvable (Figures 15–20).

Figure 15.

Static.

Figure 16.

Acceleration.

Figure 17.

Turn left.

Figure 18.

Turn right.

Figure 19.

High friction.

Figure 20.

Low friction.

Thus far, we have developed three TDMs, each with its own set of boundaries (Figures 17–20). These boundaries can be plotted in a two-dimensional plane with $F_x$ (total traction force) and $M_{\mathit{\text{zopt}}}$ (optimal DYC) as the horizontal and vertical axes, respectively. The drive stability region for each TDM is denoted by $D_i\left(i=\mathrm{1,2,3}\right)$. It is important to note that the shape of these regions varies depending on the specific vehicle and road parameters. Figure 8 displays the drive stability regions for different conditions.

It can be observed that TDM 1 is characterized by a quadrilateral shape with curved edges, while TDM 2 is represented by a band shape. Compared to the static condition, the TDM 2 region under acceleration is slightly narrower, indicating a reduction in available DYC. During left turns, TDM 1 tends to tilt towards the left side, and the upper boundary of TDM 2 is significantly reduced due to the saturated lateral tire force. This suggests the need for DYC in the opposite direction to ensure vehicle stability. On low-friction roads, both TDM 1 and TDM 2 regions are much narrower, indicating a decrease in available $M_{\mathit{\text{zopt}}}$ and $F_{\mathit{\text{xall}}}$.

4. Energy-saving oriented torque allocation strategy of distributed drive electric vehicles

4.1 The overall control framework

The four independently driven in-wheel motors endow more potential to enhance the multi-performance control requirement of distributed drive electric vehicles. The hierarchical control scheme can balance the multi-objectives through the layered control methods while simplifying the system complexity. Hence, this section introduces a dual MPC (model predictive control)-based hierarchical scheme to ensure energy conservation and stability control. The control framework details in Figure 21.

Figure 21.

Overall diagram of the proposed control strategy.

In the upper layer, the total torque inputs are generated according to the driver’s speed control requirement. Then using the energy-efficiency map obtained from the dynamometer, the optimal torque inputs are distributed to the front and rear axles. Such design can realize the energy saving through guaranteeing the in-wheel motors work in a high-efficiency zone.

In the lower layer, the additional direct yaw moment control is generated by the differential torque inputs of left and right in-wheel motors, which aims to ensure vehicle handling stability. Considering that the additional torque inputs would degrade the energy-saving performance, a relaxation factor is designed to prevent excessive control inputs based on guaranteeing the vehicle safety. Note that $\beta -\gamma$ phase plane is used to represent the vehicle stability margins.

4.2 The upper layer torque allocation strategy

This section allocates the torque inputs to the front and rear axles according to the motor efficiency map as shown in Figure 22. This is a PD18 RAM in-wheel motor. The design principle is to enable the in-wheel motors to work in a high-efficiency zone. Meanwhile, the system scheme should also take the vehicle longitudinal stability performance into account. Here, we build the wheel dynamics model to represent the rotational motion with the torque control inputs.

Figure 22.

The upper layer torque allocation strategy.

4.2.1 Wheel dynamics model

Through lumping the left and right wheels to the axle, the vehicle longitudinal motion can be represented as follows considering the tire slip ratio.

$$J_w{\dot{w}}_i=T_{\mathit{wi}}-R_e{F}_{\mathit{xi}}$$E33

$$F_{\mathit{xi}}=k_i{\lambda }_{\mathit{wi}},{\lambda }_{\mathit{wi}}=\frac{w_i{R}_e-V_{\mathit{xi},w}}{V_{\mathit{xi},w}},\left(i=f,r\right)$$E34

Combining Eq. (33) and Eq. (34), the tire rotational motion can be expressed by

$${\dot{\lambda }}_{\mathit{wi}}=-\frac{R_e^2}{J_{\mathit{wi}}{V}_{\mathit{xi},w}}{k}_i{\lambda }_{\mathit{wi}}+\frac{R_e}{J_{\mathit{wi}}{V}_{\mathit{xi},w}}{T}_{\mathit{wi}}$$E35

where $V_{\mathit{xi},w}$ and $w_i$ are the longitudinal speed and angular speed of wheel $i$, respectively. $R_e$ and $J_w$ represent rolling radius and inertia moment around $y$ axis of the wheel, respectively. $k_i$ and ${\lambda }_{w}i$ donate the tire longitudinal stiffness and slip ratio, respectively. $T_{\mathit{wi}}$ and $F_{\mathit{xi}}$ represent the torque input and tire longitudinal force, respectively. Then the state space equation of the wheel motion is given by

$$\dot{x}=\mathit{Ax}+\mathit{Bu}$$E36

where $x={\lambda }_{\mathit{wi}}$, $u=T_{\mathit{wi}}$, $A=-\frac{R_e^2}{J_{\mathit{wi}}{V}_{\mathit{xi},w}}{k}_i$, $B=\frac{R_e}{J_{\mathit{wi}}{V}_{\mathit{xi},w}}$.

4.2.2 Energy-saving controller design

The LTV-MPC (linear time varying model predictive control) is employed to handle the uncertain model parameter of longitudinal velocity. The Eq. (36) is required to be discrete first in the predictive controller. The $\mathrm{\Delta }T$ is the sampling time. Then the discrete equation is expressed as

$$x\left(k+1\right)=A^{\prime }x\left(k\right)+B^{\prime }u\left(k\right)$$E37

where $A^{\prime }=e^{A\mathrm{\Delta }T}$, $B^{\prime }={\int }_{k\mathrm{\Delta }T}^{\left(k+1\right)\mathrm{\Delta }T}{e}^{A\left[\left(k+1\right)\mathrm{\Delta }T-t\right]}\mathit{Bdt}$. The vehicle state and torque control input at time $k$ are represented by $x\left(k\right)$ and $u\left(k\right)$, respectively. It should be noted that in the LTV-MPC design, $V_{\mathit{xi},w}$ in the parameter matrices is updating at different sampling time. The LTV-MPC can guarantee the model accuracy, thereby avoiding the invalid direct yaw moment control inputs.

In the upper layer of the torque allocation strategy, the driver’s longitudinal velocity control requirement is satisfied first by the total torque control input. Here, a PI controller is employed to describe the driver longitudinal speed-tracking intention. Hence, the total in-wheel motor torque input is calculated by

$$T_{\mathit{wd}}=K_P{e}_v+K_I\int e_v\mathit{dt}$$E38

where $e_v$ denotes the speed-tracking deviation. $K_i$ and $K_p$ represent the integral and proportional coefficients. Then the following cost function is designed to realize the total torque control.

$$J_1=\sum _{t=1}^{N_p}\rho {\left(2T_{\mathit{wf}}\left(t+k\left|k\right.\right)+2T_{\mathit{wr}}\left(t+k\left|k\right.\right)-T_{\mathit{wd}}\left(k\right)\right)}^2$$E39

where $T_{w,min}\le T_{\mathit{wi}}\le T_{w,max}$. $T_{w,max}$ and $T_{w,min}$ are the admitted maximal and minimal torque control inputs, respectively. $T_{\mathit{wf}}$ and $T_{\mathit{wr}}$ represent the torque control inputs of front and rear in-wheel motors. $\rho$ and $N_p$ denote the weight coefficient and predictive horizon, respectively. Note that the predictive horizon is equal to the control horizon in this paper. Next, the energy-saving control has a priority in the upper layer. The specific method is to guarantee a higher efficiency zone for the motors. Hence, we establish a mapping function between the vehicle speed and energy efficiency. Based on the motor efficiency map in Figure 22, The most energy-efficient torque control at the current speed is selected as the reference value $T_{\mathit{wi},r}\left(i=f,r\right)$ to optimize the control inputs. Then the obtained optimal torque inputs of front and rear axles are evenly distributed to the left and right in-wheel motors. Moreover, the optimization objective of a smaller tire slip ratio is also added to the cost function and expressed by

$$J_2=\sum _{t=1}^{N_p}\left(J_2^1+J_2^2\right)$$E40

$$J_2^1=\sum _{t=1}^{N_p}\left[{\mathrm{\hslash }}_1{\left(T_{\mathit{wf}}\left(k+t\left|k\right.\right)-T_{\mathit{wf},r}\left(k\right)\right)}^2+{\mathrm{\hslash }}_2{\left(T_{\mathit{wr}}\left(k+t\left|k\right.\right)-T_{\mathit{wr},r}\left(k\right)\right)}^2\right]$$E41

$$J_2^2=\sum _{t=1}^{N_p}\left({\alpha }_1{\lambda }_{\mathit{wf}}^2\left(k+t\left|k\right.\right)+{\alpha }_2{\lambda }_{\mathit{wr}}^2\left(k+t\left|k\right.\right)\right)$$E42

where ${\alpha }_1$, ${\alpha }_2$, ${\mathrm{\hslash }}_1$, and ${\mathrm{\hslash }}_2$ are the weighting coefficients. Considering the vehicle longitudinal stability, the front axle has a priority to satisfy the high-efficiency zone. When approaching the tire force limitation, a small slip of the rear wheel would lead to the vehicle instability. Therefore, ${\mathrm{\hslash }}_1$ is endowed with a higher value. Furthermore, to ensure the driver’s longitudinal control intention, a logical judgment is also added. If the optimized torque of the rear wheel is not consistent with the driver’s control intention, the control inputs are set as 0. Through combing Eq. (39)-Eq. (42), the control objective function is represented by

$$J=J_1+J_2$$E43

4.3 The lower layer of direct yaw moment control strategy

The lower layer develops the direct yaw moment control (DYC) to enhance the vehicle handling stability based on the differential torque control inputs of left and right in-wheel motors. To improve the energy efficiency, the relaxation factor is introduced to prevent the excessive yaw moment control inputs. Here, the $\beta -\gamma$ phase plane is used to represent the vehicle stability region.

4.3.1 The vehicle dynamics modeling

A two degree-of-freedom (2-DoF) vehicle model is adopted to describe the vehicle lateral dynamics characteristics. Assuming that the vehicle runs with a small yaw angle and steering input, the vehicle model is expressed as

$$\begin{array}{l}{\mathit{mV}}_x\left(\dot{\beta }+\dot{\phi }\right)=F_{\mathit{yf}}+F_{\mathit{yr}}\\ I_z\dot{\gamma }=l_f{F}_{\mathit{yf}}-l_r{F}_{\mathit{yr}}+M_c\end{array}$$E44

where $I_z$ is the vehicle inertia moment of the yaw motion. $l_f$ and $l_r$ represent the distances from the front and rear axles to the vehicle gravity, respectively. $\beta$ and $\varphi$ denote the vehicle sideslip angle and yaw angle, respectively. $\gamma =\dot{\phi }$, $F_{\mathit{yi}}=2C_i{\alpha }_i$. $\gamma$ represents the vehicle yaw rate. The tire slip angle ${\alpha }_i$ generates the lateral force $F_{\mathit{yi}}$. $M_c$ is the additional direct yaw moment control input.

The tire slip is further written as

$$\left\{\begin{array}{l}{\alpha }_f={\delta }_f-\frac{l_f\gamma }{V_x}-\beta \\ {\alpha }_r=\frac{l_r\gamma }{V_x}-\beta \end{array}\right.$$E45

where ${\delta }_f$ is the driver steering input. Through combing Eq. (44) and Eq. (45), we can obtain

$$\dot{\xi }=\overline{A}\xi +\overline{B}\nu +\overline{C}{\delta }_f$$E46

where $\xi ={\left[\beta ,\gamma \right]}^T$, $\overline{A}=\left[\begin{array}{l}-\frac{C_f+C_r}{{\mathit{mV}}_x}\frac{C_r{l}_r-C_f{l}_f}{{\mathit{mV}}_x^2}-1\\ \frac{C_r{l}_r-C_f{l}_f}{I_z}-\frac{C_f{l}_f^2+C_r{l}_r^2}{I_z{V}_x}\end{array}\right]$, $\overline{B}={\left[01/I_z\right]}^T$, $\overline{C}={\left[\frac{2C_f}{{\mathit{mV}}_x}\frac{2C_f{l}_f}{I_z}\right]}^T$, $\nu =M_c$.

To facilitate the MPC design, the vehicle lateral dynamics model is discrete as follows.

$$\xi \left(k+1\right)={\overline{A}}^{\prime }\xi \left(k\right)+{\overline{B}}^{\prime }\nu \left(k\right)+{\overline{C}}^{\prime }{\delta }_f$$E47

The system matrices are obtained by

$${\overline{A}}^{\prime }=e^{\overline{A}\mathrm{\Delta }T},{\overline{B}}^{\prime }={\int }_{k\mathrm{\Delta }T}^{\left(k+1\right)\mathrm{\Delta }T}{e}^{\overline{A}\left[\left(k+1\right)\mathrm{\Delta }T-t\right]}\overline{B}\mathit{dt},{\overline{C}}^{\prime }={\int }_{k\mathrm{\Delta }T}^{\left(k+1\right)\mathrm{\Delta }T}{e}^{\overline{A}\left[\left(k+1\right)\mathrm{\Delta }T-t\right]}\overline{C}\mathit{dt}$$E48

Due to the uncertain model parameter of vehicle longitudinal velocity, the LTV-MPC is also adopted in the lower layer.

4.3.2 The vehicle yaw motion control design

For the vehicle yaw motion control, the sideslip angle and yaw rate are treated as important indices to represent the handling stability performance. In this paper, the steady yaw rate response and small value of sideslip angle are used as the reference value. Hence, the reference yaw motion can be represented by

$$\left\{\begin{array}{l}{\beta }_{\mathit{ref}}=0\\ {\gamma }_{\mathit{ref}}=\frac{V_x}{l_f+l_r+\frac{{\mathit{mV}}_x^2\left(C_r{l}_r-C_f{l}_f\right)}{2C_f{C}_r\left(l_f+l_r\right)}}{\delta }_f\end{array}\right.$$E49

The cost function for the MPC design can be expressed as

$$\overline{J}=\sum _{t=1}^{N_p}\left[{\lambda }_1{\left(\beta \left(t+k\left|k\right.\right)-{\beta }_{\mathit{ref}}\left(k\right)\right)}^2+{\lambda }_2{\left(\gamma \left(t+k\left|k\right.\right)-{\gamma }_{\mathit{ref}}\left(k\right)\right)}^2+{\lambda }_3{\nu }^2\right]$$E50

$$\left|{\mathrm{\Theta }}_1\gamma -{\mathrm{\Theta }}_2\beta \right|\le {\sigma }_1$$E51

$$\left|{\mathrm{\Phi }}_1\gamma -{\mathrm{\Phi }}_2\beta \right|\le {\sigma }_2$$E52

where $M_{\min }\le M_c\le M_{\mathit{max}}$. $M_{\mathit{max}}$ and $M_{\mathit{min}}$ represent the admitted maximal and minimal yaw moment control input, respectively. ${\lambda }_1$, ${\lambda }_2$, and ${\lambda }_3$ denote the weighting coefficients.

The Eq. (51) and Eq. (52) is widely used as the envelop control to describe the vehicle stability margin [1]. However, as shown in Figure 23, the direct yaw moment control input would also have an effect on the vehicle stability performance. Hence, in this work, the slack factors ${\mathrm{\Phi }}_1$ and ${\mathrm{\Phi }}_2$ in Eq. (53) and Eq. (54) are introduced to permit the vehicle runs out of the traditional stability boundaries to some extent.

$$\left|{\mathrm{\Theta }}_1\gamma -{\mathrm{\Theta }}_2\beta \right|\le {\sigma }_1+{\mathrm{\Psi }}_1$$E53

$$\left|{\mathrm{\Phi }}_1\gamma -{\mathrm{\Phi }}_2\beta \right|\le {\sigma }_2+{\mathrm{\Psi }}_2$$E54

Figure 23.

Effect of yaw-moment control on the vehicle stability region.

Furthermore, considering that the yaw moment control input would also have an effect on the energy saving performance, a small DYC control should be given when the vehicle has enough stability region. Therefore, a relaxation factor $\vartheta$ is adopted to dynamically adjust the weighting coefficients ${\lambda }_1$ and ${\lambda }_2$. Here, the relaxation factor $\vartheta$ can be calculated by

$$\vartheta =\frac{1}{2}\left(\frac{\left|{\overline{\sigma }}_1\right|}{{\sigma }_{1,m}}+\frac{\left|{\overline{\sigma }}_2\right|}{{\sigma }_{2,m}}\right)\times \left(\overline{w}-\underset{¯}{w}\right)+\underset{¯}{w}$$E55

where $\underset{¯}{w}=0.5$, $\overline{w}$, ${\sigma }_{i,m}={\sigma }_i+{\mathrm{\Psi }}_i\left(i=1,2\right)$. Then the weighting coefficients are rewritten as

$${\overline{\lambda }}_1={\mathit{\vartheta \lambda }}_1,{\overline{\lambda }}_2={\mathit{\vartheta \lambda }}_2$$E56

Then the optimal direct yaw moment control inputs are evenly allocated to the left and right in-wheel motors. The total torque inputs for each in-wheel motor are represented by

$$\left\{\begin{array}{l}{T}_{w,\mathit{fl}}=T_{w,f}-\frac{M_c}{t_w}{R}_e\\ T_{w,\mathit{fr}}=T_{w,f}+\frac{M_c}{t_w}{R}_e\\ T_{w,\mathit{rl}}=T_{w,r}\\ T_{w,\mathit{rr}}=T_{w,r}\end{array}\right.$$E57

4.4 Test results

Here, as shown in Figure 24, the hardware-in-the-loop test is conducted to verify the control effect. A high-fidelity distributed drive electric vehicle built by the commercial software Carsim is embedded into the PXI, which provides a real-time simulator. The control strategy is downloaded in the calculation platform dSPACE by the code generation technology. The calculated torque inputs are sent to the PXI through the CAN bus. The DDEV model will execute the control command. Then the vehicle states are regarded as the feedback signal and transmitted to the dSPACE to calculate the optimal control inputs. The driver steering behavior is obtained by the driving simulator. The U-turn maneuver is selected as the test condition in the HIL test. In addition, to demonstrate the proposed method (PC), the traditional torque allocation combined linear quadratic regulator (TLQ) and the proposed controller without the relaxation factor (WRF) are set as the comparative tests. The traditional torque allocation method can be represented by

$$T_{w,\mathit{fl}}=T_{w,\mathit{fr}}=\frac{l_r}{l_f+l_r}\frac{T_{\mathit{wd}}}{2},T_{w,\mathit{rl}}=T_{w,\mathit{rr}}=\frac{l_f}{l_f+l_r}\frac{T_{\mathit{wd}}}{2}$$E58

Figure 24.

HIL bench test.

Considering the limitation of the tire force, Figure 25 shows the reference vehicle speed. The proposed method behaves with good performance to track the desired speed. The vehicle lateral dynamics response is shown in Figures 26 and 27. It can be seen that the proposed method can significantly guarantee the vehicle reference yaw rate tracking performance compared with the TLQ method, while reducing the sideslip angle. Owing to the superiority to handle the uncertain model parameter, the proposed method is effective to enhance vehicle handling stability under the large-curvature road driving condition. However, the tracking error of the TLQ method is a little large. In addition, as observed from the vehicle $\beta -\gamma$ phase plane in Figure 28, the proposed method can guarantee a more safe state by the DYC control. In contrast, the vehicle runs out of the stability margins with the TLQ method. It also proves the proposed method works to balance the multi-performance control.

Figure 25.

Vehicle speed tracking performance.

Figure 26.

Vehicle yaw rate.

Figure 27.

Vehicle sideslip angle.

Figure 28.

Vehicle yaw motion phase plane.

Figures 29 and 30 show the efficiency of front and rear in-wheel motors, respectively. It is clear that the efficiency of the front in-wheel motor with the proposed method is better than that of the rear in-wheel motor. This is because when allocating the torque inputs, the proposed method is first to guarantee the high-efficiency work zone for front in-wheel motors. Hence, the efficiency of the rear in-wheel motor with the proposed method is worse than that of the TLQ method. However, from the power consumption in Figure 31, the proposed method still performs better to ensure the energy-saving performance compared with the TLQ method. However, the test results also demonstrate the proposed controller can be effective to guarantee the prescribed performance. Due to the existence of the disturbance during the HIL tests, there would have some fluctuation in test results. However, the test results also demonstrate the proposed controller can be effective to guarantee the prescribed performance. From Table 1, the proposed method can reduce energy consumption by 3.18% and 10.02% compared to the WRF method and e TLQ method, respectively.

Figure 29.

Vehicle speed tracking performance.

Figure 30.

Vehicle yaw rate.

Figure 31.

Vehicle sideslip angle.

TLQ	WRF	PC
43.1107	40.0642	38.7917

Table 1.

Energy consumption/(kJ).

Furthermore, the feasibility of the proposed energy-saving control method has been proved in the HIL test. This indicates that the proposed controller has great potential to improve the comprehensive performance of the vehicle. In the future, we would concentrate on more efficient ways of energy-saving optimization problems. The proposed strategy can not only be limited to the distributed drive electric vehicles. Meanwhile, the emergency collision avoidance condition would be also considered. In addition, from Figure 32, the tire longitudinal slip ratio with the proposed method is significantly smaller than that of the TLQ method. This demonstrates the proposed method can also improve the vehicle longitudinal stability. The test results verify the proposed method can be effective to improve the energy-saving performance based on guaranteeing the vehicle stability.

Figure 32.

Vehicle yaw motion phase plane.

5. Conclusions

This chapter introduces the distributed drive electric vehicle from the viewpoint of the dynamics modeling, stability performance analysis, and energy-saving strategy. The conventional modeling method of DDEVs is detailed first. Then, the stability region of DDEVs is estimated by establishing a rational polynomial-based DDEV model and adopting the SOSP technique to find the maximal Lyapunov function for estimation. The resulting stability regions with different parameters are presented, and comparison shows that the additional DYC has an expanding effect on the stability region. This suggests that DDEVs have greater potential in terms of stability and safety compared to centralized drive vehicles. Finally, a torque vector control framework for DDEVs is proposed in this paper to reduce the energy-consumption on the basis of maintaining the vehicle stability. The LTV-MPC-based hierarchical strategy is adopted to realize the parallel control of energy-saving and handling stability. A relaxation factor is introduced to reduce the energy consumption caused by additional direct-yaw-moment control input through evaluating the vehicle stability performance.

The proposed stability analysis method also has some issues to solve, in which the developed mode decision theorem and division of drive stability regions are mainly based on the tire adhesion ellipse theorem. However, the nonlinearity of the vehicle dynamics model also has an influence on the stability performance. In future research, theorems of body stability including $\gamma -\beta$ phase diagram and $g-g$ diagram will be considered in the torque distribution method design.

Furthermore, the feasibility of the proposed energy-saving control method has been proved in the HIL test. This indicates that the proposed controller has great potential to improve the comprehensive performance of the vehicle. In the future, we would concentrate on more efficient ways of energy-saving optimization problems. The proposed strategy can not only be limited to the distributed drive electric vehicles. Meanwhile, the emergency collision avoidance condition would be also considered.

Abbreviations

DDEV	distributed drive electric vehicle
EV	electric vehicle
CDV	centralized drive vehicle
DYC	direct yaw moment control
LTV-MPC	linear time varying model predictive control
SOSP	Sum of Squares Programming
EECA	energy-efficient control allocation
KKT	Karush-Kuhn-Tuckert
SOS	sum of squares
LMI	linear matrix inequality
RoA	region of attraction
TDMS	three Torque distribution Methods
MPC	model predictive control
2-DoF	two degree-of-freedom
TLQ	torque allocation combined linear quadratic regulator
WRF	without the relaxation factor
PC	proposed controller
LMPC	linear model predictive control
AT	average torque