Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles

*Jinhao Liang,Tong Shen, Ruiqi Fang and Faan Wang*

## **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.

**Keywords:** distributed drive electric vehicle, vehicle system dynamics, torque vector control, stability performance analysis, energy-saving

### **1. Introduction**

Electric vehicles (EVs) have been regarded as one of the effective green transportation in urban traffic due to the zero-emission [1–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 torquevector 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 motorefficiency 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 energyefficiency. 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

### *Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles DOI: http://dx.doi.org/10.5772/intechopen.111908*

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–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:


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 <sup>Δ</sup><sup>1</sup> <sup>¼</sup> <sup>P</sup><sup>Γ</sup> *η*¼1 1 *<sup>E</sup>*max<sup>Γ</sup> *E<sup>η</sup>* � � � �, the vehicle speed tracking performance <sup>Δ</sup><sup>2</sup> <sup>¼</sup> <sup>P</sup><sup>Γ</sup> *η*¼1 1 *ev*, max <sup>Γ</sup> *ev*,*<sup>η</sup>* � � � �, and the tire slip ratio <sup>Δ</sup><sup>3</sup> <sup>¼</sup> <sup>P</sup><sup>Γ</sup> *η*¼1 1 *<sup>λ</sup>t*, max <sup>Γ</sup> *λ<sup>t</sup>*,*<sup>η</sup>* � � � � are defined, where *Eη*, *ev*,*<sup>η</sup>*, and *λ<sup>t</sup>*,*<sup>η</sup>* denote

**Figure 1.** *The performance indices with different control strategies.*

the energy consumption, the speed tracking error, and the sum of absolute value of tire longitudinal slip ratio at each sampling time, respectively. Γ is the total test time. *Emax*, *ev*, *max* , and *λt*, *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.

### **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:

$$\begin{cases} \dot{V}\_{\text{x}} = V\_{\text{y}}\chi + \frac{1}{m} \left[ (F\_{\text{x1}} + F\_{\text{x2}}) \cos \delta + F\_{\text{x3}} + F\_{\text{x4}} - \left( F\_{\text{y1}} + F\_{\text{y2}} \right) \sin \delta \right] \\\ \dot{V}\_{\text{y}} = -V\_{\text{x}}\chi + \frac{1}{m} \left[ \left( F\_{\text{y1}} + F\_{\text{y2}} \right) \cos \delta\_{\text{f}} + \left( F\_{\text{y3}} + F\_{\text{y4}} \right) \right] \\\ \dot{\gamma} = \frac{1}{I\_{x}} \left[ \left( F\_{\text{y1}} + F\_{\text{y2}} \right) \cos \delta\_{\text{f}} - l\_{r} \left( F\_{\text{y3}} + F\_{\text{y4}} \right) \right] + \frac{1}{I\_{x}} M\_{\text{c}} \end{cases} \tag{1}$$

**Figure 2.** *Schematic diagram of a vehicle planar motion model.* *Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles DOI: http://dx.doi.org/10.5772/intechopen.111908*

The external yaw moment *Mc* 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\tag{2}$$

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

$$F\_{\mathbf{y}1} + F\_{\mathbf{y}2} = \mathbf{C}\_{\mathbf{f}} a\_{\mathbf{f}},\\ F\_{\mathbf{y}3} + F\_{\mathbf{y}4} = \mathbf{C}\_{r} a\_{r} \tag{3}$$

And then the tire slip angles can be expressed as:

$$\begin{cases} a\_f = \delta\_f - \frac{V\_\mathcal{V} + \eta l\_f}{V\_\mathcal{X}} \\\\ a\_r = \frac{\eta l\_r - V\_\mathcal{V}}{V\_\mathcal{X}} \end{cases} \tag{4}$$

The rotational dynamics of each wheel can be represented by

$$J\_{wi}\dot{w}\_i = -R\_c F\_{xi} + T\_{wi} \tag{5}$$

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

$$F\_{xi} = \frac{1}{R\_{\varepsilon}} (T\_{wi} - J\_{wi}\dot{w}\_i) \tag{6}$$

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

$$\begin{cases} \dot{V}\_x = V\_\gamma \chi + \frac{1}{mR\_e}(u\_1 + u\_2) + d\_1\\ \dot{V}\_\gamma = -\frac{(\mathbf{C}\_f + \mathbf{C}\_r)V\_\gamma}{mV\_x} + \left(\frac{\mathbf{C}\_r l\_r - \mathbf{C}\_f l\_f}{mV\_x} - V\_x\right)\chi + \frac{\mathbf{C}\_f}{m}u\_3 + d\_2\\ \dot{\gamma} = \frac{\left(\mathbf{C}\_r l\_r - \mathbf{C}\_f l\_f\right)V\_\gamma - \left(\mathbf{C}\_f l\_f^2 + \mathbf{C}\_r l\_r^2\right)\chi}{I\_x V\_x} + \frac{\mathbf{C}\_f l\_f}{I\_x} u\_3 + \frac{l\_s}{I\_x R\_{\text{eff}}}(u\_1 - u\_2) + d\_3 \end{cases} (7)$$

where *u*<sup>1</sup> and *u*<sup>2</sup> represent the total motor torque values in the vehicle longitudinal direction on the left and right sides of the vehicle, respectively. *u*<sup>3</sup> represents the front wheel steering angle.

$$\begin{cases} u\_1 = T\_{w1} \cos \delta + T\_{w3} \\ u\_2 = T\_{w2} \cos \delta + T\_{w4} \\ u\_3 = \delta \end{cases} \tag{8}$$

$$\begin{cases} d\_1 = \frac{1}{m} \left( F\_{\mathcal{Y}1} + F\_{\mathcal{Y}2} \right) \sin \delta + \frac{1}{m R\_\epsilon} \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\_{\mathcal{Y}1} + F\_{\mathcal{Y}2} \right) (\cos \delta - 1) \\\ d\_3 = \frac{l\_f}{I\_x} \left( F\_{\mathcal{Y}1} + F\_{\mathcal{Y}2} \right) (\cos \delta - 1) + \frac{l\_s}{I\_x R\_\epsilon} \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{cases} \tag{9}$$

### **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\_{\chi}(a\_i) = \frac{p\_1 a\_i^3 + p\_2 a\_i}{q\_1 a\_i^4 + q\_2 a\_i^2 + q\_3} \tag{10}$$

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

$$F\_{yi} = \frac{F\_{xi}\mu\_i \left(p\_1(a\_i/\mu\_i)^3 + p\_2(a\_i/\mu\_i)\right)}{q\_1(a\_i/\mu\_i)^4 + q\_2(a\_i/\mu\_i)^2 + q\_3} \tag{11}$$

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

$$F\_{yi} = n\_i(a\_i) / d\_i(a\_i), (i = f, r) \tag{12}$$

According to the kinematic characteristic, the tire slip angle *α* can be represented by the vehicle state variables.

**Figure 3.** *Stability region estimation of different degrees.*

*Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles DOI: http://dx.doi.org/10.5772/intechopen.111908*

$$\begin{cases} a\_f = \delta\_f - \beta - \frac{l\_f}{V\_\chi} \gamma \\\\ a\_r = \frac{l\_r}{V\_\chi} \gamma - \beta \\\\ s.t. \quad \beta = V\_\gamma / V\_\chi \end{cases} \tag{13}$$

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{\mathbf{x}} = f(\mathbf{x}, \boldsymbol{\mu}) = \frac{N(\mathbf{x}, \boldsymbol{\mu})}{D(\mathbf{x}, \boldsymbol{\mu})} \tag{14}$$

### **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 x*ð Þ, *p x*ð Þ with *n* real variables and *m* degrees. *p x*ð Þ is called sum of squares (SOS) if there exist polynomials *fi* ð Þ *x* such that

$$p(\mathbf{x}) = \sum\_{i=1} f\_i^2(\mathbf{x}) \tag{15}$$

Lemma 1 (Quadratic form of polynomial.) For a polynomial *p x*ð Þ, there definitely exist a symmetric matrix *Q* such that

$$p(\mathbf{x}) = \mathbf{z}^T(\mathbf{x})\mathbf{Q}\mathbf{z}(\mathbf{x})\tag{16}$$

where *z x*ð Þ is a vector of all monomials of degree less than or equal to *<sup>m</sup>* <sup>2</sup> . Commonly, the matrix *Q* is not unique. The matrix space of *Q* could be represented as a function of *λ<sup>i</sup>* (17).

$$\mathbf{Q}(\lambda) = \mathbf{Q}\_0 + \sum\_{i=1}^{N} \lambda\_i \mathbf{M}\_i \tag{17}$$

Lemma 2 (Sum of Squares Programming.) For a polynomial *pi* ð Þ *x* is SOS if and only if there exist *λi*, *i* ¼ 1, … , *N* such that

$$\begin{aligned} \mathbf{Q}(\lambda) &= \mathbf{Q}\_0 + \sum\_{i=1}^{N} \lambda\_i \mathbf{M}\_i \ge \mathbf{0} \\ \text{s.t.} &p(\mathbf{x}) = \mathbf{z}^T(\mathbf{x}) \mathbf{Q}(\lambda) \mathbf{z}(\mathbf{x}) \end{aligned} \tag{18}$$

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 x*ð Þ, *i* ¼ 0, … , *m* such that

$$\begin{array}{l} p\_0(\mathbf{x}) \ge \mathbf{0} \\ \text{s.t.} \mathbf{x} \in D \end{array} \tag{19}$$

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

$$p\_1(\mathbf{x}) \ge \mathbf{0}, \dots, p\_m(\mathbf{x}) \ge \mathbf{0} \tag{20}$$

The inequality (19) hold if there exist *qi* ð Þ *x* , *i* ¼ 1, … , *m* such that

$$\begin{aligned} p\_0(\mathbf{x}) - \sum\_{i=1}^m q\_i(\mathbf{x}) p\_i(\mathbf{x}) &\ge \mathbf{0} \\ \text{s.t.} &\mathbf{x} \in \mathbb{R}^n \end{aligned} \tag{21}$$

### **2.4 Stability region estimation and analysis**

Lemma 4 (Invariant subset of RoA.) Consider a function *V* and *γ* >0. Region Ω*<sup>V</sup>*,*<sup>γ</sup>* is defined as *<sup>x</sup>*∈*R<sup>n</sup>* f g : *V x*ð Þ≤*<sup>γ</sup>* . If conditions in (22) holds, <sup>Ω</sup>*<sup>V</sup>*,*<sup>γ</sup>* is an invariant subset of RoA.

$$\begin{aligned} \Omega\_{V, \chi} &\text{ is bounded} \\ V(\mathbf{0}) &= \mathbf{0}, V(\mathbf{x}) > \mathbf{0} (\forall \mathbf{x} \in \mathbf{R}^n) \\ \Omega\_{V, \chi} &\{ \mathbf{0} \} \subset \{ \mathbf{x} \in \mathbf{R}^n : \nabla V(\mathbf{x}) f(\mathbf{x}) \le \mathbf{0} \} \end{aligned} \tag{22}$$

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 x*ð Þ to expand the stability region of DDEV (23). By maximize *β*, *V* tends to approach the maximal Lyapunov function.

$$\max \beta : \{ \mathfrak{x} \in \mathbb{R} | \mathfrak{s}(\mathfrak{x}) < \beta \} \subset \Omega\_{V, \mathfrak{y}} \tag{23}$$

Theorem 1 For the nonlinear dynamic system, the stability region Ω*<sup>V</sup>*,*<sup>γ</sup>* can be found by finding *<sup>V</sup>*, *<sup>q</sup>*1, *<sup>q</sup>*<sup>1</sup> <sup>∈</sup> <sup>P</sup> *<sup>n</sup>*, *<sup>m</sup>* ½ � *x* and positive *γ* >0 that maximize *β* such that

$$\max \beta$$

$$\begin{cases} V(\mathbf{x}) - \varrho\_1(\mathbf{x}) \in \sum [\mathbf{x}] \\ q\_1(\mathbf{x})(\mathbf{s} - \beta) - (V - \gamma) \in \sum [\mathbf{x}] \\ -\nabla\_{\mathbf{x}} \text{VN} - D\varrho\_2(\mathbf{x}) + Dq\_2(\mathbf{x})(V - \gamma) \in \sum [\mathbf{x}] \end{cases} \tag{24}$$

### *Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles DOI: http://dx.doi.org/10.5772/intechopen.111908*

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 straightrunning, 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 *Vx* and road adhesion coefficient *μ* are two primary factors that influence lateral stability. To explore their effects on the stability region, we varied *Vx* and *μ* and plotted the estimated stability regions in phase portraits, as shown in **Figures 4**–**7**. Our results indicate that as *Vx* 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 *μ* 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

**Figure 4.** *Vx* ¼ *20m=s* ð Þ *μ* ¼ *1 .*

*Vx* ¼ *30m=s* ð Þ *μ* ¼ *1 .*

**Figure 6.** *μ* ¼ *1 V*ð Þ *<sup>x</sup>* ¼ *20m=s .*

these insights, we can better understand the factors influencing lateral stability, and design more effective control strategies for safer vehicle operation.

### *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 *δ<sup>f</sup>* and *Mz*, 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.

*Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles DOI: http://dx.doi.org/10.5772/intechopen.111908*

**Figure 7.** *μ* ¼ *0:5 V*ð Þ *<sup>x</sup>* ¼ *20m=s .*

**Figure 8.** *Mz* ¼ *0 nm* ð Þ *δ* ¼ *12*° *.*

For the CDV, the only control input is the steering angle. **Figures 8** and **9** illustrate the stability regions for a steering angle of �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 straightrunning, 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 9.** *Mz* ¼ *0 nm* ð Þ *δ* ¼ �*12*° *.*

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 *Mz* in the same or opposite directions of the steering angle, as shown in **Figures 10** and **11**. For left-turning conditions, a 300 *Nm Mz* significantly reduces the stability region of the yaw rate. However, an opposite �800 *Nm Mz* greatly expands the stability region. Compared to 0 and 300 *Nm Mz*, 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.** *Mz* ¼ *300 nm* ð Þ *δ* ¼ *12*° *.*

*Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles DOI: http://dx.doi.org/10.5772/intechopen.111908*

**Figure 11.** *Mz* ¼ �*300 nm* ð Þ *δ* ¼ �*12*° *.*

**Figure 12.** *Mz* ¼ �*800 nm* ð Þ *δ* ¼ *12*° *.*

### **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 *Mzopt*, while the driver's actions provide the total traction force Fxall. To transform *Fxall* and *Mzopt* into four-wheel torques, a rational distribution algorithm is necessary. A number of approaches consider *Fxall* and *Mzopt* as a strict constraint and employ optimization techniques to solve the problem. However, the excessive values of *Fxopt* and *Mzopt* could surpass the adhesion limit during

**Figure 13.** *Mz* ¼ *800 nm* ð Þ *δ* ¼ �*12*° *.*

handling limit. The difficulty of distribution control lies in whether to satisfy *Fxall* and *Mzopt* 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 max* <sup>1</sup> in which both *Fxall* and *Mzopt* can be fulfilled without sacrificing stability. Second, if we disregard the *Fx* requirement, we can obtain a larger region *D max* <sup>2</sup> where *Mzopt* can be met without losing stability. For the areas beyond *D max* 2, it is not possible to fully comply with either *Fxall* or *Mzopt* without compromising stability. As *Fxall* and *Mzopt* 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.*

*Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles DOI: http://dx.doi.org/10.5772/intechopen.111908*

### **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.

$$\begin{cases} m\left(\dot{\nu}\_{\rm x} - \gamma V\_{\rm y}\right) = F\_{\rm xfl} + F\_{\rm xfr} + F\_{\rm xrl} + F\_{\rm xrr} \\ m\left(\dot{\nu}\_{\rm y} + \gamma V\_{\rm x}\right) = F\_{\rm yfl} + F\_{\rm yfr} + F\_{\rm yrl} + F\_{\rm yrr} \\ I\_{x}\dot{\nu} = l\_{f}\left(F\_{\rm yfl} + F\_{\rm yfr}\right) - l\_{r}\left(F\_{\rm yrl} + F\_{\rm yrr}\right) + t\_{w}\left(F\_{\rm xfr} - F\_{\rm xgl} + F\_{\rm xrr} - F\_{\rm xrl}\right) / 2 \end{cases} \tag{25}$$

where *tw Fxfr* � *Fxfl* <sup>þ</sup> *Fxrr* � *Fxrl* � �*=*2 is the actual DYC *Mzc*. 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\_{x\circ j}}{\mu\_x}\right)^2 + \left(\frac{F\_{y\circ j}}{\mu\_y}\right)^2 \le F\_{x\circ j}\,^2\tag{26}$$

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 *hg* represents the height of the center of mass and *W* is the wheelbase.

$$\begin{aligned} F\_{xfl} &= \frac{mgl\_r}{2L} - \frac{m\dot{v}\_x h\_\text{g}}{2L} - \frac{m\dot{v}\_y h\_\text{g} l\_r}{t\_w L} \\ F\_{xrl} &= \frac{mgl\_f}{2L} + \frac{m\dot{v}\_x h\_\text{g}}{2L} - \frac{m\dot{v}\_y h\_\text{g} l\_f}{t\_w L} \end{aligned} \qquad F\_{xrr} = \frac{mgl\_r}{2L} + \frac{m\dot{v}\_x h\_\text{g} l\_r}{t\_w L} \tag{27}$$

The front and rear lateral tire force is represented as:

$$\begin{cases} F\_{\mathcal{Y}} = F\_{\mathcal{Y}l} + F\_{\mathcal{Y}\dot{r}} = \left( l\_r m \left( \dot{\nu}\_{\mathcal{Y}} + \gamma V\_x \right) + I\_x \dot{\gamma} + M\_{xc} \right) / L\\ F\_{\mathcal{Y}r} = F\_{\mathcal{Y}rl} + F\_{\mathcal{Y}r} = \left( l\_f m \left( \dot{\nu}\_{\mathcal{Y}} + \gamma V\_x \right) - I\_x \dot{\gamma} - M\_{xc} \right) / L \end{cases} \tag{28}$$

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

$$\begin{aligned} k\_{\rm fl} &= F\_{\rm xfl} / \left( F\_{\rm xfl} + F\_{\rm xfr} \right) \\ k\_{\rm rl} &= F\_{\rm xrl} / \left( F\_{\rm xrl} + F\_{\rm xrr} \right) \end{aligned} \qquad \begin{aligned} k\_{\rm fr} &= F\_{\rm xfr} / \left( F\_{\rm xfl} + F\_{\rm xfr} \right) \\ k\_{\rm rr} &= F\_{\rm xrr} / \left( F\_{\rm xrl} + F\_{\rm xrr} \right) \end{aligned} \tag{29}$$

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

$$\begin{cases} F\_{\rm yfi} = k\_{\rm fi} \left( l\_r m \left( \dot{\upsilon}\_\gamma + \chi V\_\chi \right) + I\_z \dot{\chi} + M\_\infty \right) / L \\\ F\_{\rm yri} = k\_{\rm ri} \left( l\_f m \left( \dot{\upsilon}\_\gamma + \chi V\_\chi \right) - I\_z \dot{\chi} - M\_\infty \right) / L \end{cases} \tag{30}$$

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 *Mzopt*, if and only if the linear matrix inequality (LMI) shown in Eq. (27) is solvable, where *<sup>X</sup>* <sup>¼</sup> *Fxall Mzopt* � �*<sup>T</sup>* , *j* ¼ *l*,*r.*

$$\begin{aligned} &A\_{m1}X \ge b\_{m1} \\ &s.t.A\_{m1} = \begin{bmatrix} -\frac{1}{2} & \frac{1}{2} & -\frac{1}{2} & \frac{1}{2} \\ 1 & -\frac{1}{t\_w} & -1 & \frac{1}{t\_w} \end{bmatrix}^T \\ &b\_{m1} = \begin{bmatrix} -A\_{fl} - A\_{rl} & -A\_{rl} - A\_{fl} & -A\_{fr} - A\_{rr} & -A\_{rr} - A\_{fr} \end{bmatrix}^T \\ &A\_{f\bar{\eta}} = \sqrt{\left(\mu\_x F\_{\sharp\bar{\eta}}\right)^2 - \left(k\_{f\bar{\eta}}\mu\_x \left(l\_r m \left(\dot{\bar{\upsilon}}\_y + \gamma V\_x\right) + I\_z \dot{\bar{\eta}} + M\_{\text{xc}}\right) / \mu\_y L\right)^2} \\ &A\_{r\bar{\jmath}} = \sqrt{\left(\mu\_x F\_{\text{z}\eta}\right)^2 - \left(k\_{\bar{\eta}}\mu\_x \left(l\_y m \left(\dot{\bar{\upsilon}}\_y + \gamma V\_x\right) - I\_z \dot{\bar{\nu}} - M\_{\text{xc}}\right) / \mu\_{\text{y}} L\right)^2} \end{aligned} \tag{31}$$

Lemma 6 (LMI based conditions of TDM 2)

$$\begin{aligned} A\_{m2}X &\ge b\_{m2} \\ \text{s.t.} \ A\_{m2} &= \left[ -\frac{2}{t\_w} \quad \frac{2}{t\_w} \right]^T \\ b\_{m2} &= \left[ -3A\_{fl} - 3A\_{fr} - A\_{rl} - A\_{rr} \quad -3A\_{fl} - 3A\_{fr} - A\_{rl} - A\_{rr} \right]^T \end{aligned} \tag{32}$$

Lemma 7 (LMI based conditions of TDM 3) For the given total traction force *Fxall* and optimal DYC *Mzopt*, TDM3 is satisfied, if and only if (24) and (25) is unsolvable (**Figures 15**–**20**).

**Figure 15.** *Static.*

*Dynamics Modeling and Characteristics Analysis of Distributed Drive Electric Vehicles DOI: http://dx.doi.org/10.5772/intechopen.111908*

**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 *Fx* (total traction force) and *Mzopt* (optimal DYC) as the horizontal and vertical axes, respectively. The drive stability region for each TDM is denoted by *Di*ð Þ *i* ¼ 1,2,3 . 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 *Mzopt* and *Fxall*.
