3.1. Brushless DC Motor Model
A three-phase six-state inductive (Hall sensor) brushless DC motor is considered here [
26]. The mathematical model for the brushless DC motor is as follows:
In particular, Equation (17) is the voltage balance equation, in which
u is the voltage instantaneous value matrix of each phase stator,
i is the current instantaneous value matrix of each phase stator,
Req is the equivalent resistance value matrix of each phase armature,
L is the inductance matrix, and
Eeq is the instantaneous value matrix of each phase-induced electromotive force.
Equation (18) is the mechanical kinetic equation, in which
Td is the electromagnetic torque of the motor,
n is the actual speed of the motor,
GD2 is the flywheel torque,
Tl is the load torque,
BN is the mechanical damping, and ω is the rotor angular velocity. For simplicity, Equation (19) can be substituted into Equation (18). The inductance
L can be calculated by Equation (20):
The brushless DC motor is a three-phase motor, in which the three-phase winding of the three-phase motor is symmetric. In an electric cycle, the mutual inductance between adjacent two phases is the same, and the self-inductance of each phase is the same [
27]. This can be expressed by Equations (21) and (22), as follows:
where
r is the equivalent resistance of each phase winding. The voltage balance equation can be described by Equation (24):
where
ρ denotes that the current of each phase is differentiated from time. The induced electromotive force of each phase winding can be calculated by Equation (25):
where
N/2
a represents the total number of conductors on a branch,
Bavg represents the average magnetic density which can be calculated by Equation (26),
l is the length of the conductor, and
v is the speed of the conductor cutting the magnetic induction line, which can be calculated by Equation (27):
where
ϕm is the unipolar magnetic flux, and
τ is the polar distance.
P is the number of pole-pairs. The expression of induced electromotive force in Equation (28) can be obtained from Equations (25)–(27):
The back electromotive force constant
Ce can be calculated by Equation (29):
The electromagnetic torque
Td can be calculated by Equation (30):
where
Tavg is the average torque of a conductor, which can be calculated by Equation (31):
where
favg is the average electromagnetic force on a conductor. This can be calculated by Equation (32):
where
i′ is the current through a single conductor. The electromagnetic torque
Td can also be calculated by Equation (33):
Using the formula above, the torque constant
Cd can be obtained. At the same time,
i can be obtained from Equations (18), (19), and (33) as:
Ignoring the disturbed magnetic field
BN and bringing Equation (35) into the voltage balance equation, we obtain:
The Laplace transform and deformation can be obtained by Equation (37):
3.2. Chaotic Random Grey Wolf Proportional Integral Differential (CR-GWO-PID) Control Algorithm
The Grey Wolf Optimization (GWO) algorithm is a meta-heuristic optimization algorithm that simulates the biological behavior of a grey wolf population [
18]. The algorithm is primarily used to simulate and solve the optimization problem through cooperative and competitive biological predation within the grey wolf population. However, the GWO algorithm has the disadvantages of low accuracy, poor global search ability, and falling into local optima. In view of this, a Chaotic Random Grey Wolf Optimization (CR-GWO) algorithm has been developed to address these shortcomings [
28]. The essence of the CR-GWO-PID algorithm also lies in its simulation of swarm intelligence through natural behaviors.
(1) Hybrid Optimization of Population Initialization Improvement Strategy Population initialization significantly influences the convergence speed of global optimization. The purpose of initialization is to generate an initial set of solutions within the search space, typically by random means, designed to cover the entire search space. This allows the algorithm to explore a diverse range of potential solutions. At present, Tent, Chebyshev, Gauss, and other chaotic maps are widely used in group optimization algorithms, due to their randomness, repeatability, and non-recombination [
29]. On the basis of these advantages, the Kent chaotic map uses a mixed sine and cosine random assignment [
30] to initialize the population. The objective is to introduce a degree of randomness into the population initialization process, thereby reducing the risk of convergence to a local optimal. The mathematical model is as follows:
In order to avoid falling into a short period during initialization, it is necessary to ensure that the initial value is not equal to the control parameter
μ in Equation (38); furthermore,
μ ≠ 0.5 and
μ ∈ (0, 1). At the same time, in order to avoid the emergence of repeated populations, a strategy of sine and cosine random allocation can be adopted on the basis of the Kent chaotic map. In Equation (39), the weight coefficient
u = 2, the random coefficient
p ∈ [0, 1], and the random coefficient
v ∈ [0, 1]. Then,
xk is the chaotic sequence value obtained in the
kth iteration. As shown in
Figure 6, first, the region of the sampled value was divided into 10 × 10 grids (the grid is the calculation method program running is the simulation step) to calculate the number of midpoints of each unit grid, following which the density distribution value was calculated by kernel density estimation (KDE). The density distribution value for the CR-GWO algorithm was 0.020737, while the density distribution value for the GWO algorithm was 0.0015117. Therefore, compared with the GWO algorithm, the CR-GWO algorithm had a more uniform distribution, higher randomness, and less repeatability.
(2) Non-Linear Inertia Decreasing Improvement Strategy with Weight Factors. When searching for prey, wolves are affected by a distance weight
A. When |
A| > 1, the grey wolf is far away from the prey, and a large-scale search allows for finding the optimal solution in a more global manner. When |
A| < 1, the search range shrinks and the pursuit of the prey begins [
31]. Compared with the linear decreasing strategy in the GWO algorithm, the improved CR-GWO algorithm introduces a sinusoidal non-linear decreasing factor, which enhances the global search ability in the early stage, improves the local search ability in the later stage, and improves the accuracy of the optimal solution of the algorithm [
32]. At the same time, it helps the algorithm explore a wide range of solution Spaces and increases the possibility of finding promising regions. The expression of the sinusoidal non-linear decreasing factor is:
where
t is the current number of iterations,
tmax is the maximum number of iterations and
n is the decreasing weight;
n = 2 was selected for this study and the maximum number of iterations was set to 100. The GWO algorithm was compared with the CR-GWO algorithm, and the results are shown in
Figure 7.
(3) Random Proportional Displacement Strategy. A random proportional weight was added to the CR-GWO algorithm [
33]. By introducing randomness within the population, the algorithm can converge more swiftly toward potential solutions, reducing unnecessary iteration cycles during the search process and thereby enhancing convergence speed. Simultaneously, the random proportional displacement strategy somewhat reduces the positional disparities among individuals within the population. This aids in accelerating the local search process, enabling the algorithm to converge more rapidly toward the vicinity of local optima, thereby enhancing the algorithm’s precision in local search. The mathematical model expression is:
when the mathematical model is abstracted as a plane image,
X1,
X2, and
X3 are the coordinates of the three vertices of the triangle, where
m1,
m2, and
m3 represent the weight of the distance ratio at any point within the triangle.
X1,
X2, and
X3 correspond to three grey wolves, denoted alpha, beta, and delta, respectively. They determine the location of omega wolves that are constantly close to their prey for attack. The random expansion of the weight factor
ρ ∈ (0, 1) expands the distance ratio of
X2 and
X3 to the internal point, which means that the internal point is closer to
X1. In the algorithm, the position of omega is closer to that of alpha, and alpha’s position is closer to the prey, such that the algorithm needs fewer iterations when seeking the optimal solution and the convergence speed is faster. Based on this improvement, the implementation steps for the CR-GWO algorithm are shown in
Figure 8.
In order to verify its optimization ability six test functions (three unimodal benchmark functions and three multimodal benchmark functions) were selected to test the global solution ability of CR-GWO algorithm and the ability to jump out of local optima, in comparison with the GWO algorithm. The test functions are detailed in
Table 2.
The number of populations was set to 50 and the maximum number of iterations was 100. In order to fully verify the convergence accuracy and evaluation stability of the proposed algorithm to avoid randomness, the CR-GWO and GWO were each run 50 times, and the average and standard deviation of the optimal fitness were taken to measure the performance of the algorithms. The results and data are shown in
Figure 9.
The data results after 50 runs are provided in
Table 3. Compared with the GWO algorithm, the solution accuracy of the CR-GWO algorithm was improved by about 25.13%, 13.68%, 32.65%, 39.20%, 30.82%, and 34.46% on the six test functions, respectively, with an average increase of about 29.323%. In addition, compared with the GWO algorithm, the CR-GWO algorithm consistently obtained the optimal solution with fewer iterations. These results demonstrate that the CR-GWO algorithm has faster convergence speed, higher convergence accuracy, and stronger stability.
On this basis, the parameters for the PID controller of the in-wheel motor were adjusted using the CR-GWO optimization algorithm, and a motor simulation experiment was carried out using the optimized parameters [
34,
35]. The parameters for the brushless DC motor are shown in
Table 4. The preset range of PID parameters was [0, 10], the number of iterations of the algorithm was 100, and the number of populations was 30. In order to test the superiority of the improved CR-GWO-PID control algorithm, compared with other control strategies, open-loop control, self-designed PID parameter value control (PID parameters tuned based on engineering experience), and GWO-PID control were compared in the test. The test was carried out at different speeds (600 rpm, 500 rpm, and 400 rpm) [
25]. The motor simulation process is shown in
Figure 10, and the self-designed PID parameter values were
KP = 0.5,
KI = 0.005, and
KD = 0.001 [
36,
37]. The effects of the four different control methods are shown in
Figure 11.
According to
Table 5,
Table 6 and
Table 7, CR-GWO-PID control demonstrated significant advantages over the other control methods in the simulations at three different speeds. Specifically, when compared with the open-loop control, the rise time of the CR-GWO-PID control was reduced by about 41.25% on average, while the settling time was reduced by 41.92% on average; compared with the self-designed PID control, the rise time of the CR-GWO-PID control was reduced by about 1.587% on average, the settling time was reduced by 20.97% on average, and the peak time was reduced by 27.96% on average. Compared with the GWO-PID control, the rise time of the CR-GWO-PID control was reduced by about 0.957% on average, the settling time was reduced by 9.32% on average, and the peak time was reduced by 1.90% on average. In addition, during simulations at three different speeds, the overshoot of the CR-GWO-PID control did not exceed the range of 10–20%, indicating excellent control accuracy. Open-loop control does not have an effective control effect. While the self-designed PID control exhibits good acceleration performance, it suffers from excessive overshoot. Additionally, manually tuning PID parameters requires extensive testing to achieve the desired results. Compared to the optimized algorithm-calculated PID parameters, it shows significant disadvantages and is therefore not suitable. In contrast to the first two control methods, under GWO-PID control, the system shows relatively good results in terms of overshoot and settling time calculation. However, under CR-GWO-PID control, the system can achieve better acceleration performance and stability with lower overshoot.
The results above demonstrate that the CR-GWO-PID control algorithm presents excellent dynamic performance in the case of simulating different speeds, with a faster rise time, shorter settling time, and lower overshoot, while maintaining excellent acceleration performance and stability.