Research on operating domain optimization of power split hybrid electric vehicle based on global bifurcation and chaos threshold

Abstract

The power-split hybrid electric vehicle (PS-HEV) has multiple working modes to maintain high operation efficiency according to different conditions. The main modes involved in the vehicle driving process are pure electric mode and the hybrid driving mode. Because the electromechanical coupling problem is involved in the above two working modes, the transmission system exhibits strong non-linear characteristics. If the operation range of the engine and motor are unreasonable, the rotor system will vibrate and become instability. In this paper, the non-linear dynamic equations of the electromechanical coupling of the transmission system are established for electric driving mode and hybrid driving mode. The closed-homoclinic phase trajectory equation at the center point of the disturbance-free Hamilton system is determined. The chaotic thresholds for the pure electric and hybrid driving modes are derived through the Melnikov’s method to obtain the optimal working domain of the engine and motor. Finally, numerical simulation analysis is conducted to verify the feasibility of the work domain optimization scheme. Simulation results show that the proposed engine and motor working area optimization scheme can effectively avoid the homoclinic bifurcation in the PS-HEV during the driving process and prevent the vehicle from entering the chaotic state.

Keywords

Electromechanical coupling global bifurcation chaos Melnikov analysis

Introduction

The configuration of the hybrid vehicle powertrain differs significantly from that of the traditional automobile. For the hybrid vehicle, there is electromechanical coupling between the motor and the load, as well as the motor and the engine during the power transmission process, which demonstrates strong nonlinear dynamics. One of the key challenges of hybrid electric vehicles is the electromechanical coupling due to multiple energy sources.^1,2 If the operation range of the engine and motor are unreasonable, the rotor system will vibrate and become instability, even resulting in the shaft fracture of the power system. The nonlinear electromechanical coupling characteristics directly affect the performance, work efficiency, and safety of hybrid electric vehicles. Therefore, research on the dynamic characteristics and instability mechanism of hybrid systems has important theoretical significance on the system design and application.

With the development of the non-linear dynamics theory, non-linear characteristics of non-linear systems, such as vibration, instability, periodic solutions, bifurcations, and chaos, have become research hotspots. Hu et al. proposed an irregular internal excitation (engine excitation and motor excitation) and an external excitation (road excitation) would cause torsional vibrations of the transmission system. Besides, a simplified two-mass nonlinear dynamic model of a hybrid vehicle power system was established, and then applied to the prediction of the non-linear dynamics of the torsional instability range of the hybrid electric vehicle (HEV) powertrain. The analysis results were used to optimize the hybrid control strategy in order to improve the torsional stability of the hybrid vehicle power system under three typical driving modes.³ The research on the torsional vibration of a hybrid vehicle power system mainly focused on the influence of the engine excitation, natural frequency analysis, and the corresponding modal analysis. Further, Hu et al. proposed a non-linear dynamic model of a hybrid power transmission system was established considering the multi-frequency excitation, and then the model was transformed into a fast-slow model with two-time scales. The De Moivre equation was introduced to unify the slow variables into a single variable, and the dynamic equations at different excitation frequencies and amplitudes were derived. Also, the bifurcation theory was used to study the bifurcation and stability analysis of the equilibrium point. Through the numerical analysis, the effects of the excitation frequency and amplitude on a dynamic behavior were studied using the equilibrium point curve, phase change phase diagram, and time history curve. The simulation results showed that the bifurcation would cause the jumping of the system trajectory and the corresponding blasting oscillation. With the continuous power and speed improvements, the electromechanical coupling effect became more and more prominent.⁴ Zhou et al. studied the nonlinear dynamic characteristics of a system under the internal and external excitations were studied. The effects of the speed and load fluctuations on the system dynamic response were analyzed. The results showed that the increase in the excitation amplitude would cause the system to enter the chaotic motion through frequency doubling, and the speed change could effectively control system vibrations.⁵ Liu et al. studied the dynamic characteristics of the gear transmission system under a multi-frequency excitation based on the Hopf bifurcation theory and revealed the unstable mathematical mechanism. The largest Lyapunov index revealed the influence of the system parameters, including the friction coefficient and modal damping, on the dynamic response from the table period to the unstable, chaotic motion.⁶ Zhao et al. studied in-depth research on the torsional vibration of wind turbine transmission systems was conducted. By using the maximum Lyapunov exponent and Poincare interface, the effects of the internal and external incentives on the system dynamic response were studied. The results showed that the external excitation has the greatest influence on the system dynamic characteristics, and changes in the excitation amplitude led to the system chaotic motion. In addition, it was also found that various internal stimuli, such as stiffness and damping, could significantly change the system dynamic response.^7–9 Bai et al. proposed a more comprehensive system model, including a nonlinear permeance network model (PNM), of a motor with a lateral torsion dynamic planetary gear system was established. The non-linear PNM considered the winding distribution of the motor, slotting of the stator and rotor, and non-linearity of the material. The lateral torsional dynamic planetary gear system considered the time-varying gear meshing stiffness, and the time-varying meshing and damping in the planetary gear system were studied in the proposed model. The simulation results were compared with the results predicted by the electromechanical model of a dynamic motor model with constant inductance. The results showed that the electromechanical coupling effect could cause additional severe vibration of the gear. In addition, it was found that external conditions, especially voltage transients, could have a significant impact on the dynamic characteristics of electromechanical systems.¹⁰

Although the electromechanical coupling characteristics of the hybrid electric vehicle have been studied by many scholars, most of the studies focus on the influence of motor parameters, motor torque fluctuations, or engine torque fluctuations on the electromechanical coupling characteristics, or have only considered the electromechanical coupling process between the motor and the load. Namely, different working modes have been rarely included in studying the electromechanical coupling process between the motor and engine, and the motor and load. The chaotic phenomenon of the electromechanical coupling system was not considered during the mode switching coordination process. The main innovation of this paper is to analyze the system chaotic and bifurcation characteristics in pure electric mode and hybrid drive mode by establishing the electromechanical coupling torsional vibration model, so as to provide the basis for the torque distribution scheme during the mode switching of the PS-HEV. In the chaos analysis of the electromechanical coupling characteristics, the effects of different excitations, such as the engine, motor and the load, on the system are fully considered.

The remainder of this paper is organized as follows. Section I analyzes the configuration of the HEV powertrain and the corresponding working mode. Section II establishes the motor torque, engine torque models, and electromechanical coupling torsional vibration models of the power system in the pure electric and hybrid driving modes. Section IV analyzes the bifurcation and chaos thresholds of the two working modes, and Section IV determines the optimization scheme of the engine and motor operation area based on the chaos threshold. Section V conducts the numerical simulation analysis. Section VI concludes the paper.

The main contributions of this paper can be summarized as follows.

The nonlinear dynamic equations of the power split system for the pure electric and hybrid driving modes are established to demonstrate the system electromechanical coupling characteristics. The equilibrium point of the disturbance-free Hamilton system is determined and analyzed. The closeness of the disturbance-free Hamilton system at the center Homoclinic orbit phase trajectory equation is obtained.

By adding the disturbance term to the Hamilton system, the chaos thresholds for the pure electric and hybrid driving modes are determined with Melnikov method to obtain the optimization scheme of the working domains of the engine and motor.

PS-HEV structure

PS-HEV powertrain structure

The schematic diagram of the PS-HEV powertrain system studied in this paper is shown in Figure 1. The schematic diagram of the PS-HEV powertrain¹¹ studied in this paper is shown in Figure 2, where it can be seen that the transmission system consists of the engine, motor MG1, motor MG2, brakes CB1, CB2, and two planetary rows PG1 and PG2. The engine is connected to carrier C1, MG1 is connected to the sun gear S1, and MG2 is connected to the sun gear S2. In addition, the ring gear R1, carrier C2, and the output shaft are connected, and the ring gear R2 is fixed. Since the torque of MG1 can regulate the speed of the engine, the transmission in the power system is eliminated, leading to more compact configuration. The brakes CB1 and CB2 can control the engine and motor MG1, making them rotate freely or be fixed, and can also reduce the torque ripple of the engine and motor MG1.

Figure 1.

The schematic diagram of a hybrid powertrains system.

Figure 2.

The schematic diagram of hybrid power transmission system.

By controlling the working status of each power source and brake, the hybrid vehicle can switch between different working modes. According to the working status of the power source and the brake, eight working modes of the hybrid vehicle can be realized: electric driving mode, compound pure electric driving mode with MG1 and MG2, engine starting mode, hybrid driving mode, charging while standstill mode, compound braking mode, mechanical braking mode and stopping mode. Among them, the switching process from electric driving mode to hybrid driving mode involves various intermediate states such as engine startup and speed regulation. Compared with other driving mode switching, it is more complicated, as shown in Table 1.

Table 1.

Mode division of power split hybrid electric vehicles.

Mode	Eng	MG1	MG2	CB1	CB2
electric driving	Off	Off	On	Off	Off
hybrid driving mode	On	On	On	Off	Off

The electromechanical coupling of the transmission system in the pure electric mode includes mainly the coupling between motor MG2 and load. On the other hand, the electromechanical coupling of the transmission system in the hybrid driving mode includes mainly the coupling between motor MG2, engine, and load.

Therefore, the nonlinear dynamics of an electromechanical coupling system that are studied in this paper are as follows: (1) the nonlinear torsional vibrations caused by the electromechanical coupling between motor MG2 and load; (2) the nonlinear torsional vibrations caused by the electromechanical coupling between motor MG2, engine, and load.

Electromechanical coupling torsional vibration model of power system

Electromagnetic torque model of permanent magnet synchronous motor

In order to simplify the analysis, the following assumptions are made^12,13:

(1) The winding current changes with the time according to the cosine law, and the windings are star-connected; (2) the rotor is cylindrical, and the air gap is uniform; (3) the rotating air gap magnetic field ignores the effects of higher harmonics; and (4) the motor loss is neglected.

When the coil has the alternating current that changes with the cosine law, and the influence of harmonics is not considered, according to the theory of electrical machinery, the magnetomotive force of the three-phase symmetrical stator winding is expressed as^12,14:

F_{s} (α, t) = F_{sm} \cos (ω_{1} t - α)

(1)

where $F_{sm}$ denotes the fundamental amplitude of the stator magnetomotive force, and F_s_m=1.35NI_mk_N1 /p; N represents the number of per phase windings, I_m denotes the current, pre presents the pole pairs of a PMSG; k_N1 is the coefficient of the fundamental wave magnetization, and $ω_{1}$ is electrical angular frequency of the motor.

The magnetomotive force of the rotor is given by:

F_{r} (α, t) = F_{rm} \cos (ω_{1} t - α)

(2)

where $F_{rm}$ represents the fundamental amplitude of the magnetomotive rotor force, and $F_{rm} = B_{r} h_{m} \sin (\frac{α_{p} π}{2}) / π / μ_{0}$ ; $B_{r}$ denotes the material remanence of the permanent magnet, h_m is the permanent magnet thickness, $μ_{0}$ represents the air permeability coefficient, and $μ_{0} = 4 π \times 10^{- 7}$ ; and $α_{p}$ denotes the arc coefficient of the permanent magnet.

The air-gap magnetic field energy of a permanent magnet synchronous motor is given by:

\begin{matrix} W = \frac{Rl}{2} \int_{0}^{2 p π} Λ (α, t) F^{2} (α, t) d α \\ = \frac{pRl Λ π}{2} [F_{sm}^{2} + F_{rm}^{2} - 2 F_{sm} F_{rm} \sin (p φ + Ψ)] \end{matrix}

(3)

where R denotes an average radius of the air gap, and $l$ denotes the effective length of a rotor.

By considering the torsional vibration angle of the PMSG, the rotor synthesis fundamental wave magnetomotive force can be expressed as:

\begin{matrix} F (α, t) = F_{sm} \cos (ω_{1} t - α) \\ + F_{rm} \cos [ω_{1} t - α + (p φ + Ψ + \frac{π}{2})] \end{matrix}

(4)

If the rotor torsional angle is denoted as $φ$ , then the electromagnetic torque is expressed as:

\begin{matrix} T_{e} = \frac{dW}{d (p φ)} = F_{m} \cos (p φ + Ψ) \\ = [\sin Ψ \sin (p φ) - \cos Ψ \cos (p φ)] \end{matrix}

(5)

where $F_{m} = p π Rl Λ F_{sm} F_{rm}$ , $R$ denotes the inner radius of a stator, $l$ represents the effective length of a rotor, and $Λ$ is the air gap permeability of a PMSG.

The simplified equation of the electromagnetic torque for the PMSG can be obtained by Taylor approximation.

T_{e} = k_{0} + k_{1} φ - k_{2} φ^{2} - k_{3} φ^{3}

(6)

where $k_{0} = F_{m} \cos φ$ , $k_{1} = p F_{m} \sin φ$ , $k_{2} = p^{2} F_{m} \cos φ / 2$ , $k_{3} = p^{3} F_{m} \sin φ / 6$ , and $k_{0}$ , $k_{1}$ , $k_{2}$ , and $k_{3}$ represent the electromagnetic parameters. When the rotor torsional vibration angle is equal to zero, that is, $φ$ = 0, the motor is in a stable running state, and $(k_{1} φ - k_{2} φ^{2} - k_{3} φ^{3})$ denotes a disturbance term.

Engine shaft torque model

The unbalanced torque the engine is mainly composed of the reciprocating moment of inertia $T_{j}$ and the combustion gas moment $T_{g}$ .^15,16 The variation law of the combustion gas torque is complicated, which can be expressed by Fourier series theory as follows:

T_{g} = T_{0} + \sum_{r = 1}^{\infty} (a_{r} sinr ω_{2} t + b_{r} cosr ω_{2} t)

(7)

where $T_{0}$ denotes an average driving torque, and $ω_{2}$ denotes the angular velocity of the shafting.

Since the fourth- and higher-order unbalanced moments account for only 2.5% of that of the second-order, only the second-order unbalanced torque is considered and given by:

{\begin{matrix} T_{j} = 2 λ m_{j} r^{2} ω_{2}^{2} \sin 2 α \\ T_{g} = \frac{π r D^{2}}{4} (a_{2} \sin 2 α + b_{2} \cos 2 α) \end{matrix}

(8)

Where $m_{j}$ represents the reciprocating equivalent mass of the connecting rod and piston group, $r$ denotes the crank radius, $D$ is the piston diameter, and $λ$ represents the connecting rod ratio is expressed as $λ = r / l$ .

The cosine component of the combustion gas moment $T_{g}$ is small compared with its sine component. When the torsional vibration angle is considered, the simplified torque model of the engine is given by:

\begin{matrix} T_{E} = T_{E 0} + (2 λ m_{j} r^{2} ω_{2}^{2} - \frac{π r D^{2}}{4} a_{2}) \sin 2 ω_{2} t \\ = T_{E 0} + F_{2} \sin 2 ω_{2} \end{matrix}

(9)

where $T_{E 0}$ denotes a constant part of the engine torque.

Torsional vibration model

Gear pair dynamics model

In order to study the electromechanical coupling characteristics of the PS-HEV powertrain, the dynamical model of the planetary gear is introduced and shown in Figure 3.

Figure 3.

The planetary gear dynamics model.

The general equation of the pure torsional dynamics of the planetary gears dynamics model that is presented in Figure 3 is as follows.¹⁷

For the sun wheel:

J_{s} {\overset{\cdot\cdot}{θ}}_{s} + \sum_{i = 1}^{n} (C_{spi} {\overset{\cdot}{δ}}_{spi} + k_{spi} δ_{spi}) R_{s} = T_{1}

(10)

where $θ_{s}$ denotes the rotational angle of the sun gear, $J_{s}$ represents the rotational inertia of the system equivalent to the sun gear, $R_{s}$ denotes the radius of the base circle for the sun gear, $T_{1}$ denotes the driving torque; $C_{spi}$ denotes the meshing damping between the sun gear and the ith planetary gear, $k_{spi}$ denotes the meshing stiffness between the sun gear and the ith planetary gear, $i$ is the number of planetary gears, and $δ_{spi}$ is the relative displacement between the sun gear to the ith planet gear.

For planetary wheels:

J_{p} {\overset{\cdot\cdot}{θ}}_{pi} + [\begin{matrix} (C_{spi} {\overset{\cdot}{δ}}_{spi} + k_{spi} δ_{spi}) + \\ (C_{rpi} {\overset{\cdot}{δ}}_{rpi} + k_{rpi} δ_{rpi}) \end{matrix}] R_{P} = 0

(11)

where $θ_{pi}$ denotes the rotation angle of the planet gear, $J_{p}$ denotes the rotational inertia of the system equivalent to the planetary wheel, $R_{P}$ denotes the radius of the base circle for the planetary wheel radius, $C_{rpi}$ denotes the meshing damping between the ring gear and the ith planetary gear, $k_{rpi}$ denotes the meshing stiffness between the ring gear and the ith planetary gear, and $δ_{rpi}$ is the relative displacement of the ring gear to the ith planetary gear.

For the planet carrier:

J_{C}^{'} {\overset{\cdot\cdot}{θ}}_{C} - \sum_{i = 1}^{n} [\begin{matrix} (C_{spi} {\overset{\cdot}{δ}}_{spi} + k_{spi} δ_{spi}) + \\ (C_{rpi} {\overset{\cdot}{δ}}_{rpi} + k_{rpi} δ_{rpi}) \end{matrix}] R_{C} = - T_{2}

(12)

where $θ_{C}$ denotes the rotation angle of the planet carrier, $J_{C}^{'}$ denotes the rotational of inertia of the system equivalent to the planet carrier, and $J_{C}^{'} = (J_{C} + \sum_{i = 1}^{n} m_{pi} R_{c}^{2})$ ; $m_{pi}$ is the mass of the ith planet wheel, and $R_{C}$ denotes the radius of the base circle of the planet carrier.

Also, $δ_{spi}$ and $δ_{rpi}$ can be respectively expressed as:

{\begin{matrix} δ_{spi} = R_{s} (θ_{s} - θ_{C}) + R_{p} (θ_{p} - θ_{C}) \\ δ_{rpi} = R_{p} (θ_{p} - θ_{C}) + R_{r} θ_{C} \end{matrix}

(13)

where $R_{p}$ denotes the radius of the base circle of the planetary gear, $θ_{p}$ denotes the rotation angle of the planetary gear, and $R_{r}$ denotes the radius of the base circle of the ring gear.

Nonlinear model of electromechanical coupling rotor in electric driving mode

According to the torsional vibration model of the gear pair given by equations (10)–(12) the electromechanical coupled torsional vibration model at the pure electric mode is established,¹⁸ and it is given by:

{\begin{matrix} T_{1} = T_{e} = T_{e 0} + Δ T_{e} \\ T_{2} = T_{L} = T_{L 0} + Δ T_{L} \end{matrix}

(14)

where $T_{e 0}$ denotes the constant part of the electromagnetic torque, $T_{L 0}$ denotes the constant part of the load torque, $Δ T_{e}$ denotes the electromagnetic torque disturbance, and $Δ T_{L}$ is the load torque disturbance.

According to the torsional vibration equation of the gear pair, equation (14) can be substituted into equations (10)–(12), we have:

{\begin{matrix} \begin{matrix} J_{s} {\overset{\cdot\cdot}{θ}}_{s} + \sum_{i = 1}^{n} (C_{spi} {\overset{\cdot}{δ}}_{spi} + k_{spi} δ_{spi}) R_{s} = \\ k_{1} θ - k_{2} θ^{2} - k_{3} θ^{3} \end{matrix} \\ J_{p} {\overset{\cdot\cdot}{θ}}_{pi} + [(\begin{matrix} C_{spi} {\overset{\cdot}{δ}}_{spi} + \\ k_{spi} δ_{spi} \end{matrix}) + (\begin{matrix} C_{rpi} {\overset{\cdot}{δ}}_{rpi} + \\ k_{rpi} δ_{rpi} \end{matrix})] R_{P} = 0 \\ \begin{matrix} J_{C}^{'} {\overset{\cdot\cdot}{θ}}_{C} - \sum_{i = 1}^{n} [\begin{matrix} (C_{spi} {\overset{\cdot}{δ}}_{spi} + k_{spi} δ_{spi}) + \\ (C_{rpi} {\overset{\cdot}{δ}}_{rpi} + k_{rpi} δ_{rpi}) \end{matrix}] R_{C} = \\ - F_{1} \cos ω_{1} t \end{matrix} \end{matrix}

(15)

New variables are set to $x = R_{s} (θ_{s} - θ_{C}) + R_{p} (θ_{p} - θ_{C})$ . Then, the nonlinear torsional vibration model of the hybrid powertrain is given by:

\overset{\cdot\cdot}{x} + μ \overset{\cdot}{x} + (1 - κ) x + β x^{2} + γ x^{3} = f \cos ω_{1} t + Δ

(16)

where $μ = \frac{R_{s}^{2} C_{spi}}{J_{s}} - \frac{R_{s} R_{C} C_{spi}}{J_{C}^{'}} - \frac{R_{p} R_{C} C_{spi}}{J_{C}^{'}} + \frac{R_{p}^{2} C_{spi}}{J_{p}}$ ; $(1 - κ) = (\frac{R_{s}^{2} C_{spi}}{J_{s}} - \frac{R_{s} R_{C} C_{spi}}{J_{C}^{'}} - \frac{R_{p} R_{C} C_{spi}}{J_{C}^{'}} + \frac{R_{p}^{2} C_{spi}}{J_{p}}) + \frac{R_{s} r}{J_{s}} k_{1}$ ; $β = \frac{R_{s}}{J_{s}} r^{2} k_{2}$ ; $γ = \frac{R_{s}}{J_{s}} r^{3} k_{3}$ ; $r = \frac{1}{1 + v}$ ; $v = \frac{J_{s}}{J_{C}^{'}}$ ; $f = \frac{(R_{s} + R_{p})}{J_{C}^{'}} F_{1}$ .

Also, $Δ$ represents the ring gear disturbance item, and

\begin{matrix} Δ = \frac{(C_{rpi} {\overset{\cdot}{δ}}_{rpi} + k_{rpi} δ_{rpi}) R_{P}}{J_{p}} \\ - \frac{(R_{s} + R_{p}) \sum_{i = 1}^{n} (C_{rpi} {\overset{\cdot}{δ}}_{rpi} + k_{rpi} δ_{rpi}) R_{C}}{J_{C}^{'}} . \end{matrix}

Nonlinear model of electromechanical coupling rotor in hybrid driving mode

According to the torsional vibration model of the gear pair given by equations (10)–(12) the electromechanical coupled torsional vibration model at the hybrid driving mode is established,³ and it is expressed as:

{\begin{matrix} T_{1} = T_{e} = T_{e 0} + Δ T_{e} \\ T_{2} = T_{L} + T_{E} = T_{L 0} + Δ T_{L} + T_{E 0} + Δ T_{E} \end{matrix}

(17)

where $T_{E 0}$ denotes the constant part of the engine shaft torque, and $Δ T_{E}$ denotes the engine torque disturbance.

According to the torsional vibration equation of the gear pair, equation (17) can be substituted into equations (10)–(12) which results in the following relations:

{\begin{matrix} \begin{matrix} J_{s} {\overset{\cdot\cdot}{θ}}_{s} + \sum_{i = 1}^{n} (C_{spi} {\overset{\cdot}{δ}}_{spi} + k_{spi} δ_{spi}) R_{s} = \\ k_{1} θ - k_{2} θ^{2} - k_{3} θ^{3} \end{matrix} \\ J_{p} {\overset{\cdot\cdot}{θ}}_{pi} + [(\begin{matrix} C_{spi} {\overset{\cdot}{δ}}_{spi} + \\ k_{spi} δ_{spi} \end{matrix}) + (\begin{matrix} C_{rpi} {\overset{\cdot}{δ}}_{rpi} + \\ k_{rpi} δ_{rpi} \end{matrix})] R_{P} = 0 \\ \begin{matrix} J_{C}^{'} {\overset{\cdot\cdot}{θ}}_{C} - \sum_{i = 1}^{n} [\begin{matrix} (C_{spi} {\overset{\cdot}{δ}}_{spi} + k_{spi} δ_{spi}) + \\ (C_{rpi} {\overset{\cdot}{δ}}_{rpi} + k_{rpi} δ_{rpi}) \end{matrix}] R_{C} = \\ - F_{1} \cos ω_{1} t - F_{2} \sin 2 ω_{2} t \end{matrix} \end{matrix}

(18)

New variables are set to $x = R_{s} (θ_{s} - θ_{C}) + R_{p} (θ_{p} - θ_{C})$ . Then, the nonlinear torsional vibration model of the hybrid powertrain is given by:

\begin{matrix} \overset{\cdot\cdot}{x} + μ \overset{\cdot}{x} + (1 - κ) x + β x^{2} + γ x^{3} \\ = f_{1} \cos ω_{1} t + f_{2} \sin 2 ω_{2} t + Δ \end{matrix}

(19)

Work area optimization based on global bifurcation and chaos threshold

Equilibrium point analysis of undisturbed Hamilton system

The excitations of the engine and load input, and the engine damping terms are regarded as disturbance. Equations (16) and (19) are respectively rewritten as the form of the state equations^19,20:

{\begin{matrix} \overset{\cdot}{x} = y \\ \overset{\cdot}{y} = - (1 - κ) x - β x^{2} - γ x^{3} - ε μ \overset{\cdot}{x} + \\ ε (f \cos ω_{1} t + Δ) \end{matrix}

(20)

{\begin{matrix} \overset{\cdot}{x} = y \\ \overset{\cdot}{y} = - (1 - κ) x - β x^{2} - γ x^{3} - ε μ \overset{\cdot}{x} + \\ ε (\begin{matrix} f_{1} \cos ω_{1} t + f_{2} \\ \sin 2 ω_{2} t + Δ \end{matrix}) \end{matrix}

(21)

where $ε$ denotes the small parameter, and $0 < ε < 1$ ; also, the value of $ε Δ$ is small enough and can be ignored here.

The potential function of the undisturbed system has an important on studying different nonlinear dynamic behaviors. Therefore, it is necessary to study the function of the undisturbed system in conducting electromechanical coupled nonlinear dynamics. When $ε = 0$ , the damping and disturbance terms are excluded, and equations (20) and (21) can be transformed into disturbance-free system, which can be expressed as follows:

{\begin{matrix} \overset{\cdot}{x} = y \\ \overset{\cdot}{y} = - (1 - κ) x - β x^{2} - γ x^{3} \end{matrix}

(22)

It should be noted that by discarding the damping and harmonic excitation terms, the transmission system can be converted into a Hamiltonian system. The two working modes have the same Hamilton form, and the properties of their equilibrium points are also the same. The Hamiltonian function is given by:

H (x, y) = (1 - κ) \frac{x^{2}}{2} + β \frac{x^{3}}{3} + γ \frac{x^{4}}{4} + \frac{y^{2}}{2}

(23)

And the asymmetric potential energy is expressed as:

V (x) = (1 - κ) \frac{x^{2}}{2} + β \frac{x^{3}}{3} + γ \frac{x^{4}}{4}

(24)

Let $y = \overset{\cdot}{y} = 0$ in equation (24), we obtain the equilibrium points $P_{0} (0, 0)$ , $P_{1} (x_{1}, 0)$ , $P_{2} (x_{2}, 0)$ .

Where:

x_{1, 2} = \frac{- β \pm \sqrt{β^{2} - 4 γ (1 - κ)}}{2 γ}

Then, the equilibrium point and its stability are analyzed according to the conditions of the system parameters $γ$ , $β$ , and $(1 - κ)$ . Based on the theory of electrical engineering, the properties of a fixed point $(u_{s}, 0)$ can be determined as follows.

If $V ″ (u_{s}) < 0$ , where $V (u_{s})$ is the maximum, the fixed point is the saddle point.

If $V ″ (u_{s}) > 0$ , where $V (u_{s})$ is the minimum, the fixed point is the center point.

The stability of the equilibrium point of the undisturbed Hamilton system is as follows:

(A) $β^{2} \leq 4 γ (1 - κ)$ , there is only one equilibrium point $P_{0} (0, 0)$ :

Case A₁: if $κ < 1$ , $P_{0}$ is the center.

Case A₂: if $κ > 1$ , $P_{0}$ is the saddle.

There exit no homoclinic or heteroclinic orbit.

(B) $β^{2} > 4 γ (1 - κ) and 4 γ (1 - κ) > 0$ , there are three equilibrium points $P_{0} (0, 0), P_{1} (x_{1}, 0), and P_{2} (x_{2}, 0)$ :

(1) Case B₁: if $γ > 0$ , $P_{1} (x_{1}, 0)$ is the saddle, and $P_{0} (0, 0) and P_{2} (x_{2}, 0)$ are the centers.

There exit the trajectory to intersect the homoclinic orbit.

(2) Case B₂: if $γ < 0$ , $P_{1} (x_{1}, 0)$ is the center, and $P_{0} (0, 0) and P_{2} (x_{2}, 0)$ are the saddles.

There exit heteroclinic orbit in the response of an unperturbed system.

(C) $β^{2} > 4 γ (1 - κ) and 4 γ (1 - κ) < 0$ , there are three equilibrium points $P_{0} (0, 0), P_{1} (x_{1}, 0), and P_{2} (x_{2}, 0)$ :

(1) Case C₁: if $γ > 0$ , $P_{0} (0, 0)$ is the saddle, $P_{1} (x_{1}, 0) and P_{2} (x_{2}, 0)$ are the centers.

There exit the trajectory to intersect the homoclinic orbit.

(2) Case C₂: if $γ < 0$ , $P_{0} (0, 0)$ is the center, $P_{1} (x_{1}, 0) and P_{2} (x_{2}, 0)$ are the saddles.

There exit heteroclinic orbit in the response of an unperturbed system.

Based on the parameters given in Table 2, $β = 0.03$ , $γ = 0.25$ and $κ = 1.5$ can be calculated. Then, the equilibrium point of the disturbance-free Hamiltonian system satisfies the condition C1, that is, $P_{0} (0, 0)$ is the saddle point, and $P_{1} (x_{1}, 0)$ and $P_{2} (x_{2}, 0)$ are the center points. The potential function and phase trajectory of the undisturbed Hamilton system are shown in Figure 4. It can be seen that the undisturbed Hamilton system has a closed orbit at the center point and a homoclinic orbit. When the Hamilton system is disturbed the external interference, the homoclinic orbit is split by disturbance; thus, the system produces homoclinic bifurcation, which causes chaos.

Table 2.

Main parameters of a hybrid powertrain.

Parameter	Value	Unit
Phase, m	3
Pole pairs, p	4
Each phase series windings, N	200
Inner radius of rotor, R	55	mm
Magnet thickness, hm	8	mm
Inner power factor angle, ψ	π/3
Magnetic saturation coefficient, $k_{μ}$	1.5
Length of stator core, l	85	mm
Permanent magnet remanence, $B_{r}$	1.07	T
Length of air gap, $δ_{0}$	1.5	mm
Engine shaft moment of inertia	20	$kg \cdot m^{2}$
PMSG moment of inertia	0.8	$kg \cdot m^{2}$

Figure 4.

The potential function and phase trajectory of an undisturbed Hamilton system.

In this case, $P_{0} (0, 0)$ is the saddle point, $P_{1} (x_{1}, 0)$ and $P_{2} (x_{2}, 0)$ are the center points. Then, the Hamiltonian of the system is given by:

H (x, y) = - \frac{x^{2}}{4} + \frac{x^{4}}{16} + \frac{y^{2}}{2}

(25)

When $H (0, 0) = 0$ , there are two homoclinic orbits ${q^{0} (τ)}$ connecting the saddle point $(0, 0)$ ; thus, the homoclinic orbits satisfy the following conditions:

{\begin{matrix} H (x, y) = - \frac{x^{2}}{4} + \frac{x^{4}}{16} + \frac{y^{2}}{2} = 0 \\ \overset{\cdot}{x} = y \end{matrix}

(26)

Thus, the equations of x and y can be obtained, that is, the phase trajectory equations of the homoclinic orbit of the Hamilton system, and they are as follows:

{\begin{matrix} x^{\pm} = \pm \frac{\sqrt{2}}{\cosh (\frac{1}{2} t)} \\ y^{\pm} = \mp \frac{\sqrt{2} \sinh (\frac{1}{2} t)}{2 \cosh {(\frac{1}{2} t)}^{2}} \end{matrix}

(27)

Chaos threshold determination

According to the above analysis, the undisturbed Hamilton system transformed by the PS-HEV powertrain has the homoclinic orbits, and the PS-HEV transmission system has uncertain disturbances. When the disturbance terms is added to the Hamilton system, the closed homoclinic orbit is destroyed, which will lead to the same-direction unfolding bifurcation The existence of global homoclinic bifurcation often indicates the chaotic motion. Therefore, Melnikov method can be used to estimate the parameter space that the homoclinic bifurcation may be transformed into the chaos. In this way, the threshold of chaotic motion of the system can be determined.

The essence of the Melnikov method is to calculate the distance between the stable and unstable manifolds after the disturbance of the homoclinic orbit of the undisturbed Hamilton system. Determine whether the stable and unstable orbits meet the conditions of cross-section intersection. Based on the given threshold conditions for the chaotic vibration of electromechanical coupling to determine whether there is chaos in the system.^21,22

Since the disturbance terms for the pure electric and hybrid driving modes are different, the Melnikov analysis needs to be performed separately for these two working modes.

Melnikov analysis of electric driving mode

equation (20) can be transformed into:

\overset{\cdot}{X} = [\begin{matrix} p_{1} (x, y) \\ p_{2} (x, y) \end{matrix}] + ε [\begin{matrix} q_{1} (x, y, τ) \\ q_{1} (x, y, τ) \end{matrix}]

(28)

where $X = [\begin{matrix} x \\ y \end{matrix}], p_{1} (x, y) = y$ , and $q_{1} (x, y, τ) = 0$ , $p_{2} (x, y) = - (1 - κ) x - β x^{2} - γ x^{3}$ , $q_{2} (x, y, τ) = - μ \overset{\cdot}{x} + fcos ω_{1} t$ .

According to Melnikov’s principle, the distance between stable and unstable manifolds can be expressed by the Melnikov equation,^23,24 as follows:

M^{\pm} (τ_{0}) = \int_{- \infty}^{+ \infty} p (X_{h} (τ)) \land q (X_{h} (τ), τ + τ_{0}) d τ

(29)

where $X_{h} = (x^{\pm}, y^{\pm})$ represents the homoclinic orbit of the Hamilton system, and $p \land q = p_{1} q_{2} - p_{2} q_{1}$ .

By substituting equation (28) into equation (29), the Melnikov integral is obtained as:

\begin{matrix} M^{\pm} (τ_{0}) = \\ \int_{- \infty}^{+ \infty} [- μ y^{\pm} + f \cos ω_{1} (τ + τ_{0})] y^{\pm} d τ = - μ I_{1} + f I_{2} \sin ω_{1} τ_{0} \end{matrix}

(30)

where $I_{1} = \int_{- \infty}^{+ \infty} {[y^{\pm} (τ)]}^{2} d τ$ and $I_{2} = \int_{- \infty}^{+ \infty} y^{\pm} (τ) \sin (ω_{1} τ) d τ$ .

The first integral $I_{1}$ can be directly calculated using the mathematical software, whereas the second integral $I_{2}$ is difficult to solve. Consequently, the numerical integration method is used in the paper to derive the integral $I_{2}$ . and the chaotic threshold curve are drawn by using the numerical integration method. Based on the threshold curve, the system can be qualitatively and quantitatively analyzed.

According to the Melnikov theory, if $M (τ_{0}) = 0$ , then $| μ I_{1} | = | f I_{2} \sin ω_{1} τ_{0} |$ , and if $M (τ_{0}) \neq 0$ , then $| \sin ω_{1} τ_{0} | < 1$ . The condition of producing chaos with a fixed frequency is expressed as:

- 1 < \frac{μ I_{1}}{f I_{2}} < 1

Assuming $f$ is the control parameter and $μ = 0.08$ , the Melnikov threshold curve is shown in Figure 5. The three-point $f$ -valuesin Figure 5 are $f = 1.4$ , $f = 1.6,$ and $f = 1.8$ . $ω_{1} = 610$ corresponds to $f =$ 1.58 on the threshold curve, and the corresponding amplitude of the load excitation is 120 $N \cdot M$ . In Figure 5, it can be seen that the change in $f$ determines whether the system appears chaotic motion. When $f < 1.58$ , the load excitation is smaller than 120 $N \cdot M$ , the system is in the stable state with periodic vibration, and there is no chaos; on the other hand, when $f \geq 1.58$ , the load excitation is $F_{1} \geq$ 120 $N \cdot M$ , the system passes through the stable region and then enters into the chaotic area, and the chaos phenomenon occurs.

Figure 5.

Threshold curve of the homoclinic bifurcation ( $f - ω_{1}$ ).

Melnikov analysis of hybrid driving mode

Assume $ω_{1} = 2 ω_{2}$ , and express equation (21) as:

\overset{\cdot}{X} = [\begin{matrix} p_{1} (x, y) \\ p_{2} (x, y) \end{matrix}] + ε [\begin{matrix} q_{1} (x, y, τ) \\ q_{1} (x, y, τ) \end{matrix}]

(31)

where $X = [\begin{matrix} x \\ y \end{matrix}], p_{1} (x, y) = y$ , and $q_{1} (x, y, τ) = 0$ , $p_{2} (x, y) = - (1 - κ) x - β x^{2} - γ x^{3}$ , $q_{2} (x, y, τ) = - μ \overset{\cdot}{x} + (f_{1} \cos ω_{1} t + f_{2} \sin ω_{1} t)$ .

According to Melnikov’s principle, the distance between stable and unstable manifolds can be expressed by the Melnikov equation as follows:

M^{\pm} (τ_{0}) = \int_{- \infty}^{+ \infty} p (X_{h} (τ)) \land q (X_{h} (τ), τ + τ_{0}) d τ

(32)

where $X_{h} = (x^{\pm}, y^{\pm})$ represents the homoclinic orbit of the Hamilton system, and $p \land q = p_{1} q_{2} - p_{2} q_{1}$ .

By substituting equation (29) into equation (30), the Melnikov integral is given by:

\begin{matrix} M^{\pm} (τ_{0}) = \int_{- \infty}^{+ \infty} [\begin{matrix} - μ y^{\pm} + f_{1} \cos ω_{1} (τ + τ_{0}) \\ + f_{2} \sin ω_{1} (τ + τ_{0}) \end{matrix}] \\ y^{\pm} d τ = - μ I_{1} + I_{2} \sin ω_{1} τ_{0} + f_{2} I_{3} \sin ω_{1} τ_{0} \end{matrix}

(33)

Where $I_{1} = \int_{- \infty}^{+ \infty} {[y^{\pm} (τ)]}^{2} d τ$ , $I_{2} = \int_{- \infty}^{+ \infty} y^{\pm} (τ) \sin (ω_{1} τ) d τ$ , and $I_{3} = - \int_{- \infty}^{+ \infty} y^{\pm} (τ) \cos (ω_{1} τ) d τ$

Similarly, the first integral $I_{1}$ can be directly calculated using the mathematical software, whereas the second and third integral $I_{2}$ and $I_{3}$ are difficult to solve, so in this work, the numerical integration method is used to solve them. Although the corresponding analytical formula cannot be obtained, the chaotic threshold curve can be drawn using the numerical integration method. According to the threshold curve, the system can also be qualitatively and quantitatively analyzed.

According to the Melnikov theory, if $M (τ_{0}) = 0$ , then $| μ I_{1} | = | f I_{2} \sin ω_{1} τ_{0} + f_{2} I_{3} \sin ω_{1} τ_{0} |$ , and if $M (τ_{0}) \neq 0$ , then $| \sin ω_{1} τ_{0} | < 1$ . The condition of producing chaos with fixed frequency is expressed as:

$- 1 < \frac{μ I_{1}}{f I_{2} + f_{2} I_{3}} < 1$

Assume $f_{1}$ is a control parameter and $μ = 0.08$ , then the Melnikov threshold curve is as shown in Figure 6. The three-point $f$ -values in Figure 6 are $f_{1} = 1.5$ , $f_{1} = 1.7,$ and $f_{1} = 1.8$ ; where $ω_{1} = 715$ corresponds to $f_{1} =$ 1.62 on the threshold curve, and the corresponding load excitation amplitude is 125 $N \cdot M$ . In Figure 6, it can be seen that the change in $f_{1}$ determines whether the system appears chaotic motion. When $f_{1} < 1.62$ , the load excitation is $F_{1} <$ 125 $N \cdot M,$ the system is in a stable state of periodic vibration, and there is no chaos; on the other hand, when $f_{1} \geq 1.62$ , the load excitation $F_{1} \geq$ 125 $N \cdot M$ , the system passes through the stable region, and then enters into the chaotic area; thus, the chaos phenomenon arises.

Figure 6.

The threshold curve of homoclinic bifurcation ( $f_{1} - ω_{1}$ ).

Assume $f_{2}$ is a control parameter and $μ = 0.08$ , then the Melnikov threshold curve is as shown in Figure 7. The blue line represents the relationship between the engine excitation amplitude and the angular velocity $ω_{1}$ , and the red line represents the stable boundary. When $ω_{1} < 887$ or $ω_{1} > 1079$ , $f_{2}$ is in the stable region at this time, and the system runs stably. When $887 < ω_{1} < 1079$ , the system runs in an unstable area, and chaos occurs.

Figure 7.

The threshold curve of homoclinic bifurcation ( $f_{2} - ω_{1}$ ).

Engine and motor working area optimization scheme

According to the Melnikov analysis of the PS-HEV powertrain in pure electric and hybrid driving modes, the optimization scheme of the engine and motor working domain based on the global bifurcation and chaos threshold is designed. It is assumed that the load excitation is 120 $N \cdot M$ and that it always smaller than 120 $N \cdot M$ during the vehicle running process.

In the pure electric mode, the speed of motor MG2 and the corresponding working torque are 5830 $r / \min (ω = 610 rad / s)$ and 74 $N \cdot M$ to avoid chaos. If the speed is too small, and the system will cross the chaos threshold curve, the system will enter the chaotic region. Therefore, in the pure electric mode, the torque of motor MG2 should be controlled to be not more than $74 N \cdot M,$ that is $T_{MG 2} \leq 74 N \cdot M$ .

In the hybrid driving mode, the speed of motor MG2 and the corresponding working torque needed to ensure that the chaos does not appear are 6800 $r / \min (ω = 715 rad / s)$ and $T_{MG 2} =$ 63 $N \cdot M$ . Therefore, under the condition that the required torque of the vehicle is satisfied, the torque of the control motor MG2 is $T_{MG 2} \leq 63 N \cdot M$ . If the torque $T_{MG 2}$ provided by motor MG2 is greater than 63 $N \cdot M$ , the system will enter the chaotic region, and chaos will occur. At the same time, according to the established control strategy, the engine rotational speed is $1500 r / \min \leq n_{E} \leq 6000 r / \min$ , and the angular velocity is $157 rad / s \leq ω_{2} \leq 628 rad / s$ . It has been analyzed that $887 rad / s < ω_{1} < 1079 rad / s$ , that is, when $443.5 rad / s \leq ω_{2} \leq 539.5 rad / s$ , the system runs in an unstable area. Therefore, the optimal angular speed range of the engine to ensure that the system does not appear chaotic is $157 rad / s \leq ω_{2} \leq 443.5 rad / s$ , and the corresponding optimal speed range is $1500 r / \min \leq n_{E} \leq 4237 r / \min$ .

In summary, when the vehicle is operating in the pure electric mode, MG2 provides the driving torque to drive the vehicle. The torque of MG2 is smaller than 74 $N \cdot M$ . Under this condition, if the torque provided by MG2 is insufficient, the mode should be changed and switched to the hybrid driving mode. When the vehicle is operating in the hybrid driving mode, MG2 and the engine provide the driving torque simultaneously, the torque of MG2 is $T_{MG 2} < 63 N \cdot M$ , and the remaining driving torque is provided by the engine. The optimal speed is $n_{E} =$ 1500 $r / \min$ . Combined with the fuel economy operation area of the engine, the speed range of the control engine is $1500 r / \min < n_{E} < 4237 r / \min$ .

Numerical simulation verification

In order to verify the accuracy of the results obtained by the above-mentioned Melnikov analysis, Runge-Kutta method is used to verify the torsional vibration equations (20) and (21) of the transmission system. The calculated time response history curve (excluding the transient instability part at the beginning of the simulation), phase trajectory diagram, Poincare section diagram and amplitude spectrum are analyzed.

The purpose of the numerical simulation is to verify the effectiveness of the proposed work area optimization scheme.

Numerical analysis of electric driving mode

When $T_{MG 2} = 65 N \cdot M$ , the time history, phase plan, Poincare map, and Fourier spectra obtained by the numerical analysis are shown in Figure 8. It can be seen that the time response history curve has an orderly cycle, the phase trajectory is a circular circle, the Poincare map has only one point, and only one energy pulse appears in the Fourier spectra. This indicates that the transmission system is in a periodic motion state, and the maximum dynamic tangential displacement of the system is relatively small. At this time, the system has not lost its stability, and no bifurcation and chaotic motion have occurred.

Figure 8.

The numerical simulation results at $T_{MG 2} = 65 N \cdot M$ . (a) Time history; (b) phase plane; (c) Poincare map; (d) Fourier spectra.

When $T_{MG 2} = 85 N \cdot M$ , the time history, phase plan, Poincare map, and Fourier spectra obtained by the numerical analysis of the transmission system are shown in Figure 9. As shown in Figure 9, the time history curve starts to be chaotic and disordered. The phase trajectory has multiple intersecting loops. There are multiple points on the Poincare map. The Fourier spectra has regions of energy concentration. This indicates that the time-domain response of the transmission system is non-periodic, and the maximum dynamic tangential displacement of the system is relatively large, which indicates that the system loses stability, a global homoclinic bifurcation occurs, and the system is in chaos.

Figure 9.

The numerical simulation results at $T_{MG 2} = 85 N \cdot M$ . (a) Time history; (b) phase plane; (c) Poincare map; (d) Fourier spectra.

The numerical simulation analysis results show that in pure electric mode, the torque of MG2 was $T_{MG 2} \leq 74 N \cdot M$ , and it could ensure the system stability without bifurcation and chaotic motion. Simulation results are consistent with the results obtained by the Melnikov analysis method.

Numerical analysis of hybrid driving mode

$When T_{MG 2} = 55 N \cdot M$ , $ω_{2} = 400 rad / s$ , the time history, phase plan, Poincare map, and Fourier spectra obtained by the numerical analysis of the transmission system are shown in Figure 10. It can be seen that the time response history curve has an orderly cycle, the phase trajectory is a circular circle, the Poincare map has only one point, and only one energy pulse appears in the Fourier spectra. This indicates that the transmission system is in a periodic motion state, and the maximum dynamic tangential displacement of the system is relatively small. At this time, the system has not lost its stability, and no bifurcation and chaotic motion have occurred.

Figure 10.

The numerical simulation results at $T_{MG 2} = 55 N \cdot M$ , $ω_{2} = 400 rad / s$ . (a) Time history; (b) phase plane; (c) Poincare map; (d) Fourier spectra.

When $T_{MG 2} = 70 N \cdot M$ , $ω_{2} = 500 rad / s$ , the time history, phase plan, Poincare map, and Fourier spectra obtained by the numerical analysis of the transmission system are shown in Figure 11. As shown in Figure 11, the time history curve starts to be chaotic and disordered. The phase trajectory has multiple intersecting loops. There are multiple points on the Poincare map. The Fourier spectra has regions of energy concentration. This indicates that the time-domain response of the transmission system is non-periodic, and the maximum dynamic tangential displacement of the system is relatively large, which indicates that the system loses stability, a global homoclinic bifurcation occurs, and the system is in chaos.

Figure 11.

The numerical simulation results at $T_{MG 2} = 70 N \cdot M$ ÿ¥ $ω_{2} = 500 rad / s$ . (a) Time history; (b) phase plane; (c) Poincare map; (d) Fourier spectra.

The numerical simulation analysis results showed that in the hybrid driving mode, the torque of motor MG2 was $T_{MG 2} \leq 63 N \cdot M$ , and it was controlled. Because the engine was involved in the driving process, in order to ensure the stability of the transmission system, no bifurcation and chaotic motion occurred, and the engine could work in a low fuel consumption area, achieving good fuel economy, while controlling the engine output torque to ensure meeting the driving requirements, and the engine speed range was $1500 r / \min < n_{E} < 4237 r / \min$ . These results are consistent with those obtained by the Melnikov analysis method.

Conclusion

In this paper, a PS-HEV is studied, and its two working modes that are involved in the mode switching process are analyzed. Also, the bifurcation and chaos characteristics of the hybrid vehicle transmission system are studied in different working modes, and the electric driving mode and hybrid driving mode are analyzed by Melnikov’s method. In addition, based on the chaos threshold, the working domain optimization scheme of the engine and motor is obtained. Finally, numerical simulation analysis is used to verify the feasibility of the work area optimization scheme.

By studying an undisturbed Hamilton system corresponding to the electromechanical coupling nonlinear equation of a transmission system, by analyzing the equilibrium point characteristics of a system, and by solving the undisturbed Hamilton system potential function and phase trajectory, it is concluded that the system has a closed orbit at the center point, and homoclinic exists. In orbit, when the Hamilton system experiences the external interference, the homoclinic orbit is split by disturbance, and the system produces homoclinic bifurcation, which causes the chaos.

Due to the uncertain disturbances of a PS-HEV transmission system, when the disturbance term is added to the Hamilton system, the closed homoclinic orbit is destroyed, resulting in the global homoclinic bifurcation and chaotic motion. The homoclinic can be estimated by Melnikov’s method, and the bifurcation may be transformed into a parameter space of chaos so that the threshold of the chaotic motion of the torsional vibration of electromechanical coupling can be determined. In the pure electric mode, in order to ensure the stability of the transmission system without bifurcation and chaotic motion, the torque of motor MG2 should be $T_{MG 2} \leq 74 N \cdot M$ ; and, in the hybrid drive mode, it should be $T_{MG 2} \leq 63 N \cdot M$ . Because the engine is involved in the driving process, no bifurcation and chaotic motion occurred, and the engine could work in the low fuel consumption area, achieving better fuel economy. By controlling the engine output torque to ensure the driving requirements, the engine speed range was $1500 r / \min < n_{E} < 4237 r / \min$ .

However, there are still some shortcomings that need to be addressed. For instance, a simplified construction of the transmission system mode should be obtained, and the internal motor parameters, the number of pole pairs of the motor, the magnetic circuit saturation coefficient, and other electromagnetic characteristics should also be considered.

Footnotes

Handling Editor: James Baldwin

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This work has been supported by the National Key Research and Development Program of China (2018YFB0105000 and 2017YFB0102603), National Natural Science Foundation of China (51875255, 61601203, 61773184, U1564201, U1664258, U1764257 and U1762264), Natural Science Foundation of Jiangsu Province (BK20180100), Six Talent Peaks Project of Jiangsu Province (2018-TD-GDZB-022), Key Project for the Development of Strategic Emerging Industries of Jiangsu Province (2016-1094), Key Research and Development Program of Zhenjiang City (GY2017006), National Natural Science Foundation of China (Grant No. 51905219), Natural Science Founding of Jiangsu Province (Grant No. BK20190844), and Natural Science Research Project of Jiangsu Higher Education Institutions (Grant No. 19KJB580001).

ORCID iD

Dou Lei

References

Yong-bin

Zeng-nian

Wen-xin

, et al. Design and development of few-tooth-deference planetary gear Powertrain Coupling system with biased bent shaft in electric and HEV for all-terrain vehicles. In: IEEE vehicle power and propulsion conference, Harbin, 3–5 September 2008, pp. 1–5. New York: IEEE.

Sun

Shi

Lei

, et al. Analysis and design optimization of a permanent magnet synchronous motor for a campus patrol electric vehicle. IEEE Trans Veh Tech 2019; 68: 10535–10544.

Yan

. Optimization methodology for control strategy of parallel hybrid electric vehicle based on chaos prediction. AIP Adv 2018; 8: 115305.

Yan

, et al. Dynamics analysis of the hybrid powertrain under multi-frequency excitations with two-time scales. AIP Adv 2018; 8: 065212 .

Shihua

Song

Zhaohui

, et al. Nonlinear dynamic analysis of coupled gear-rotor-bearing system with the effect of internal and external excitations. Chin J Mech Eng 2016; 29: 281–292.

Liu

Parker

. Nonlinear parametrically excited dynamics of two-stage spur gear trains with mesh stiffness fluctuation. Proc IMechE, Part C: J Mech Eng Sci 2012; 226: 1939–1957.

Zhao

. Nonlinear torsional vibrations of a wind turbine gearbox. Appl Math Model 2016; 39: 4928–4950.

Zhao

Qin

Wang

. Yaw and lateral stability control for four-wheel steer-by-wire system. IEEE/ASME Trans Mechatron 2018; 23: 2628–2637.

Gao

. Bifurcation and chaotic vibration in electromechanical integrated toroidal drive. J Vib Contr 2013; 21: 1556–1565.

10.

Bai

Qin

Wang

, et al. Dynamic characteristics of motor-gear system under load saltations and voltage transients. Mech Syst Signal Process 2018; 100: 1–16.

11.

Dehua

. Research on hybrid dynamic modeling and optimization control of energy distribution system for power split-type hybrid electric vehicles. Jiangsu: JiangSu University, 2017.

12.

. Optimization of operating domain for generator set in series hybrid electric with power follower control strategy. AIP Adv 2019; 9: 1–9.

13.

Sun

Lei

, et al. An improved model predictive current control for PMSM drives based on current track circle. IEEE Trans Ind Electron. Epub ahead of print April 2020. DOI: 10.1109/TIE.2020.2984433.

14.

Sun

Zhu

, et al. MPTC for PMSMs of EVs with multi-motor driven system considering optimal energy allocation. IEEE Trans Magn 2019; 55: 1–6.

15.

Ying

Meng

Jing

. Turbocharger rotor dynamics with foundation excitation. Arch Appl Mech 2009; 79: 287–299.

16.

Shi

Sun

Cai

, et al. Torque analysis and dynamic performance improvement of a PMSM for EVs by skew angle optimization. IEEE Trans Applied Supercond 2019; 29: 0600305.

17.

Al-shyyab

Kahraman

. A non-linear dynamic model for planetary gear sets. Proc IMechE, Part K: J Multibody Dyn 2007; 221: 567–576.

18.

Meintz

Ferdowsi

. Control strategy optimization for a parallel hybrid electric vehicle. In: Proc IEEE vehicle power and propulsion conference, Harbin, 3–5 September 2008, pp. 1–5. IEEE.

19.

Tchoukuegno

Woafo

. Dynamics and active control of motion of a particle in a

ϕ

6 potential with a parametric forcing. Physica D 2002; 167: 86–100.

20.

Chen

Peng

, et al. Bifurcation and chaos analysis of torsional vibration in a PMSM-based driven system considering electromechanically coupled effect. Nonlinear Dyn 2017; 88: 277–292.

21.

Yang

Jing

. Inhibition of chaos in a pendulum equation. Chaos Solitons Fractals 2008; 35: 726–737.

22.

Zhang

Yao

Zhang

. Using the extended Melnikov method to study the multi-pulse global bifurcations and chaos of a cantilever beam. J Sound Vibr 2009; 319: 541–569.

23.

Zhou

Chen

. Global bifurcation analysis and chaos of an arch structure with parametric and forced excitation. Mech Res Commun 2010; 37: 67–71.

24.

Najafifi

Ghazavi

Jafari

. Stability and hamiltonian hopf bifurcation for a nonlinear symmetric bladed rotor. Nonlinear Dyn 2014; 78: 1049–1064.