GWO-Based Multi-Stage Algorithm for PMDC Motor Parameter Estimation

During the design of a wheeled mobile robot, the problem of the proper selection of the parameters of its motor controllers was encountered. Knowing the parameters of the robot’s Permanent Magnet Direct Current (PMDC) motors, precise tuning of the controllers can be performed, which then results in improved robot dynamics. Among many methods of parametric model identification, optimization-based techniques, particularly genetic algorithms, have gained more and more interest recently. The articles on this topic present the results of parameter identification, but they do not refer to the search ranges for individual parameters. With too wide a range, genetic algorithms do not find solutions or are time-inefficient. This article introduces a method for determining the parameters of a PMDC motor. The proposed method performs an initial estimation of the range of searched parameters to shorten the estimation time of the bioinspired optimization algorithm.


Introduction
PMDC motors are widely utilized, e.g., as drive units of wheeled mobile robots, or as DC-DC converters [1]. High-speed accurate performance of a robot can be provided by setting the proper parameters of the motor speed/torque controllers. The most frequently used controllers are digital PID controllers [2], and their modifications [3]. Other algorithms for PMDC control can also be found. The authors in [4] introduced the model reference adaptive system (MRAS)-based method to control the motor, and compared the results with a PI controller. Different simulation methods can be used to design the controller of the PMDC motor, yet the parameters of the motor model need to be known in advance [5]. In the literature, a simple linear viscous friction model is mostly implemented. The friction model is more complex and takes into account Coulomb, Stribeck, and sticky friction as presented by the authors in paper [6]. The shape of the friction characteristics is a nonlinear curve and is often approximated by an exponential or polynomial function, as demonstrated by the authors in paper [7].
There are many methods of parametric model identification [8]. The most popular of the optimization methods used for identification is the Recursive Least-Squares (RLS) method [9]. It is recommended for its simplicity and applicability to real systems in addition to its advantage of lower computational complexity. An interesting solution for PMDC model estimation using binary-valued measurements is presented in [10]. In the paper, the suitable regression model is used to model the relationship between the input voltage and the output speed of the motor, to support online estimation utilizing cheap and lowcomplex binary sensors. A similar approach was utilized by the authors in [5,11], where the authors presented a deeper description and a wider scope of examples. A sensorless method for speed estimation of a PMDC motor is shown in [12]. The authors utilized the RLS method and the first-order mathematical model of the PMDC to reduce the number of samples required for online implementation. The authors in [13] presented a method for estimating a PMDC motor speed considering the effect of the armature teeth-slots and commutation. Different motor responses are used to gather the data required for parameter estimation. In [14], authors used a Pseudo-Random Binary Sequence (PRBS) for the system identification process. The authors in [15] utilized multiparametric programming to estimate the parameters of the DC motor. Multiparametric programming is a mathematical optimization technique that can be used to solve parameter estimation problems. In this method, the parameters of a system are estimated by minimizing a cost function subject to a set of constraints. The authors in paper [16] suggested a novel approach for parameter estimation based on the Rao-1 algorithm, which is a modified version of the Kalman filter. Another approach to motor parameter estimation is proposed in [17]. The article presents the method for identifying the parameters of a DC motor using a compound least-squares (CLS) approach. Articles [18,19] propose novel methods for identifying the parameters of a DC motor based on algebraic equations. Motor parameter estimation methods utilizing Swarm Intelligence algorithms are presented in [20,21].
Optimization-based methods, in particular the techniques based on genetic algorithms, have gained broad interest [22]. Particularly interesting is the work presented by the authors in [23], which uses the Gray Wolf Optimization (GWO) genetic algorithm to determine the parameters of a brushed electric motor. However, the downside of the method is a very wide range of searched parameters, which results in a long estimation period and uncertainty of the obtained results. Furthermore, the fitness/cost function presented by the authors can be questioned due to the lack of weighting for current and speed components. Due to the significant order-of-magnitude difference between the current and speed, the above-mentioned fitness/cost function focuses on matching the speed waveform, almost completely ignoring the current response of the motor.
All the above-mentioned methods can be characterized in terms of the accuracy, complexity, and time of parametric model identification. The time of identifying the parameters is directly related to the range in which the solution is sought. To increase the efficiency of optimization methods, it is necessary to roughly define the values of the parameters, to establish search ranges for individual parameters.
In the following paper, a new method for quick estimation of the PMDC motor parameters is proposed. The method consists of two stages. In the former one, rough values of parameters are calculated, while in the latter one, they are carefully evaluated to estimate the fine values. Since the second stage requires a little more time to finish, the method intends to shorten the whole procedure by starting fine parameter estimation with rough estimates of the values from the first stage.

Methods
This section presents the methods of modeling and estimating PMDC motor parameters.

Modeling of the PMDC Motor
The principle of operation of a PMDC motor is related to basic electrical and mechanical laws [24]. The flow of current through the rotor winding creates a force as a result of the interaction of the stator's magnetic field with the rotor's magnetic field. The torque of the force acting on the rotor is proportional to the current flowing through the rotor winding, and the mechanical constant K m . A DC motor with a permanent magnet can be modeled in the form of two analogous circuits: electrical and mechanical (see Figure 1). The parameters of the electric circuit are rotor winding resistance R, rotor winding inductance L, and electromotive constant K φ . The parameters of the mechanical circuit are friction coefficient b, rotor inertia J, mechanical constant K m , and loading torque T L . Due to the analogy between the two circuits, the mechanical circuit is modeled using electrical symbols. The mathematical model of the motor can be described by two differential equations. Using the second Kirchhoff law, the following relationships can be obtained: The total power of the motor (or any other electrical device) is expressed as P = u(t)i(t). According to Equation (1), the total power of the motor is: The motor power that is converted into mechanical power is described by the third term of Equation (3) (the first term describes the power associated with Joule-Lenz heat, while the second describes the magnetization of the coil). Since the mechanical power of a rotating motor is also defined as P m = ωK m i(t), it can be assumed that in all the previous and following Equations, there is an equality: Equations (1) and (2) form a system of differential equations, where the unknown functions are i(t) and ω(t). The most common approach to solve this differential equations system is to use numerical methods (e.g., using Simulink). However, this system can be solved analytically. The advantage of such an approach is that it speeds up further calculations performed during numerical optimization. Analytical solutions of Equations (1) and (2) take the form: where

Characteristic Points Method
To simplify and speed up calculations, the characteristic points method (CPM) of parametric model identification is introduced. The method assumes that the values of model parameters are calculated based on current and speed values at characteristic points for which the derivative is equal to zero.
Equation (1) for characteristic points t 0 and t 1 (see Figure 2) can be written as: Equation (2) in time t 1 , on the other hand, can be written in the following form: Based on Equations (7)-(9) with condition u(t 0 ) = u(t 1 ) = U, the following parameters of the PMDC motor model can be obtained: where τ m and τ e are mechanical and electrical time constants, respectively. Observing a certain dependence for most DC motors that τ e «τ m , we can simply assume that it is a first-order system. We can approximately assume that moments after forcing the voltage on the motor terminals, the system behaves like a regular RL circuit. Therefore, the electrical time constant can be estimated as the time it takes for the current to reach 63.2% of its maximum value. The error in this estimation will increase as the τ m /τ e ratio increases. Due to the fact that the CPM method is used only to determine the ranges, such a simplification is justified. It should be noted that the described CPM method is not precise, because it is strongly sensitive to unavoidable measurement noise at characteristic points t 0 and t 1 and an error in the graphical determination of T e and T m . However, it can roughly estimate the parameter values and determine the scope of their search for the optimization method.

Optimization Method
The purpose of the below analysis is to find a nonlinear fit between experimental and theoretical data using an optimization method. The experimental data are the angular speed ω(t i ) and currentî(t i ) of the motor (for example, such as those shown in Figure 2), while the theoretical data are the results ω(t i ) and current i(t i ) obtained from Equations (5) and (6), which depend on motor parameters R, K, b, J, and L.
The goal of the optimization algorithm is to find such R, K, b, J, L that minimize the following objective function: where i is the sample number and N is the number of step response samples measured in the experiment. Thus, the minimization problem can be presented in the following form: To normalize the values of motor speed and current, a weighting factor in Equation (15) is introduced. This coefficient takes the form: Equation (15) clearly shows the advantage of analytical solutions of differential Equations (1) and (2)-the values of ω(t i ) and i(t i ) can be calculated immediately from Equations (5) and (6) instead of numerical simulation.
Optimization algorithms search the solution space in different ways. Some are better suited to given issues, and some are worse. Optimization methods differ in effectiveness and time of finding solutions. To select the best bioinspired algorithm to solve the presented problem, preliminary simulation tests were performed. A set of test data was prepared, and parameters were identified using five different algorithms. The value of the objective function and the number of iterations served as a qualitative indicator for each of the algorithms. Algorithms that were selected for the comparison group were GWO algorithm [25], whale optimization algorithm (WOA) [26], particle swarm optimization (PSO) algorithm [27], chimp optimization algorithm (CHIMP) [28], and bat algorithm (BAT) [29]. The obtained results are presented in Figure 3. In solving the presented identification problem, the GWO algorithm turned out to be the most effective, which is why it was chosen.
GWO is a metaheuristic algorithm inspired by the hunting behavior of gray wolves. The algorithm mimics the social hierarchy and hunting mechanism of gray wolves to solve optimization problems.
The algorithm starts with the initialization of a population of candidate solutions, which are called wolves. The position of each wolf in the search space represents a potential solution to the problem. The population of wolves is then divided into four groups, namely alpha, beta, delta, and omega wolves Figure 4, based on their fitness values. The alpha wolf has the best fitness value, while the omega wolf has the worst fitness value. In each iteration of the algorithm, the position of each wolf is updated based on its position and the positions of the alpha, beta, and delta wolves. This update is performed using three different equations, each representing one of the three types of wolf behavior: searching, encircling, and attacking. The search behavior is used by the omega wolf to explore the search space, while the encircling and attacking behaviors are used by the alpha, beta, and delta wolves to converge toward the optimal solution.
The algorithm terminates when a stopping criterion is met, such as a maximum number of iterations or a desired level of fitness is reached. The final solution is the position of the alpha wolf, which represents the best solution found by the algorithm.
Overall, the Gray Wolf Optimizer is a simple yet effective algorithm that has shown promising results in solving a wide range of optimization problems, including engineering, finance, and machine learning. Its ability to balance exploration and exploitation makes it a popular choice for solving complex optimization problems.
The following equations are proposed to mathematically model the encircling behavior of gray wolves during the hunt: where X p is the victim position vector, X is the gray wolf position vector, t indicates the current iteration, A and D are coefficient vectors. Vectors A and C are computed as follows: where vector a linearly decreases during iteration from 2 to 0, r 1 and r 2 , are random vectors with a value between 1 and 0.
To mathematically simulate the hunting behavior of gray wolves, it is assumed that the alpha (the best candidate solution), beta, and delta have a better understanding of the potential location of the prey. Thus, the first three best solutions obtained thus far are saved, and the other search agents (including the omegas) are required to update their positions according to the position of the best search agent. In this regard, the following formulas are proposed: A more detailed description of the algorithm can be found in [25].

Combined CPM-GWO Method
A common problem that occurs when solving optimization tasks is the possibility of the objective function becoming stuck in a local minimum. To prevent this, constraints can be imposed on the parameters. This is generally problematic because it requires a priori knowledge of the range of solutions being sought.
To overcome this difficulty, a new method called CPM-GWO is proposed, which combines the two algorithms: the characteristic points method (CPM) described in Section 2.2 and the GWO introduced in Section 2.3. The estimation process in this case is as follows: first, make a rough identification of model parameters using the characteristic points method, then use model parameters p obtained during the first stage (i.e., CPM) to define the solutions range according to the following notation: The range is defined as the area of searching for individual parameters in the optimization method. Figure 5 shows an example comparison of the GWO method with the searching area taken as [0.0000001, 10] (which appears to be a reasonable area for the PMDC motor parameters) and the CPM-GWO method where the search area was narrowed using the CPM solution. It can be seen that despite using the same objective function given by Equation (15), the GWO solution is stuck in the local minimum.
The discrepancy in the GWO results can be seen in detail in Table 1, where the methods used are compared.
Step responses i(t) and ω(t) of the selected motor were measured. Motor parameters were then identified at a nominal voltage of 12 V using CPM, GWO and CPM-GWO methods.

Test Rig
A block diagram of the test rig is shown in Figure 6 and the photo of the rig in Figure 7. The motor was connected to a supply voltage source KA3005D set for the rated voltage of 12 V (1). An oscilloscope (2) was used to record the voltage and current of the motor (5). A current probe (4) was placed at the power cord to measure the current value. Motor rotational speed was measured using an optical encoder and a USB frequency recorder Teensy 3.2 (3) triggered by an oscilloscope. The measurement of the step response of the motor was forced by switching on the key controlled from the Teensy board.

Parameters Estimation Using the CPM
Based on the recorded step responses of the Buhler gear motor, the points for which the di dt = 0 condition is met were determined. Then, using Equations (10)- (14), parameters R, L, K, J, b were calculated. In the next stage, using the graphical method, time constants τ m and τ e were determined. Next, parameters J and L were calculated according to Equations (13) and (14). Step responses for model parameters obtained with the proposed characteristic points method were generated and are presented in Figure 8. Parameter values obtained through the characteristic points method are listed in Table 1. Step response for parameters calculated by the CPM.

Motor Model Identification Using CPM-GWO
Using the results obtained from the CPM method, search ranges for individual parameters were determined in accordance with Equation (24). The objective function was assumed in the form given by Equation (15). A comparison of experimental and simulation step responses obtained for the identified model parameters is presented in Figure 9. Step response for parameters calculated by the CPM-GWO.
The obtained values of the PMDC model parameters are collated in Table 1, together with the reference values given by the manufacturer of the tested motor. To improve the overall credibility of the obtained results, parameter estimation was repeated 100 times for each of the CPM-GWO and GWO methods. Table 1 describes the average values of the parameters obtained throughout this process, and the minimum, maximum and standard deviation of the estimated values.
The results shown in Table 1 clearly demonstrate that the GWO method with an a priori assumed searching range [0.0000001, 10] produces unreliable results. The average values of parameters K, J, b are not close to the factory reference values. In addition, these results are characterized by a large scatter, which can be seen in the standard deviation, min and max values. In contrast, it can be seen that the results for the CPM-GWO method are comparable to the reference factory data. According to the data in Table 1, the CPM-GWO method produces a little scatter in the results, as demonstrated by the value of the standard deviation, and in the worst cases (see min and max) the results differ from the average by a few percent at most. The exception is the value of the friction coefficient for which the standard deviation is a percentage point larger (σ = 7.81 × 10 −7 to avg = 6.11 × 10 −6 ) compared to the rest of the parameters. The problem with the friction coefficient is discussed further.
Step responses of the motor for the calculated parameters and the reference values provided by the manufacturer are shown in Figure 10.

Motor Identification for a Wide Range of Supply Voltages
The identification of motor parameters described in Sections 3.2 and 3.3 was performed for a nominal voltage of 12 V. Thus, the results obtained can be compared with the parameters declared by the motor manufacturer. However, for purposes such as simulation, it is required to know the parameters of the motor for any voltage and not just the nominal one. To determine the parameters for a wide range of supply voltages, the measurements of motor step responses were carried out for voltages from 5 V to 16 V. In the next step, the CPM-GWO method was applied to the obtained data and the motor parameters were determined. The obtained results are shown in Figure 11 (for a few selected supply voltages), in Figure 12 and in Table 2.   As can be seen in Figure 12 and Table 2, the values of almost all estimated parameters for different voltages are constant and they fluctuate around the average value. However, in the case of friction coefficient b, its value depends on the angular speed of the shaft. This problem with the friction coefficient is discussed further in the next section.
Similar results as presented in Table 2 can be also obtained using a single optimization calculation, using the measured parameters for different voltages to obtain a single set of parameters for the motor, valid for different voltages. Such techniques of parameter optimization have been recently proposed by the authors in [31].

Discussion
To confirm the dependence between friction coefficient b and angular speed ω, the friction characteristics of the motor as a function of its angular speed were measured. The measurement was carried out at the same test stand as explained in Section 3. The obtained results are shown in Figure 13b. The obtained friction force moment bω (Figure 13a,b) corresponds to the Stribeck curve [32], which is characteristic for DC motors [33]. For lower rotational speeds (less than 200 rad/s) the dependency between the friction moment and angular speed is highly nonlinear. For higher rotational speeds (over 200 rad/s) the dependency is almost linear (both for positive and negative speeds). These observations agree with Stribeck friction theory, according to which at lower speeds dry friction (boundary lubrication) is present, and for higher rotational speeds viscous friction (hydrodynamic lubrication) takes advantage. In the viscous friction area, the viscous friction coefficient is usually assumed to be constant, and these almost constant values of b can be observed in Figure 13b for rotational speeds greater than 200 rad/s. For lower speeds (less than 200 rad/s) the relationship between b and ω is highly nonlinear, because, in this speed range, dry friction takes advantage. Thus, the obtained results confirm that the commonly assumed linear viscous friction model is not correct for all rotational speed ranges, and Stribeck or other more accurate friction models [6,7] should be included. To compare the results of the two methods, the root mean square error (RMSE) qualitative indicator was proposed: where p is the step response sample for the parameters calculated with a given method,p is the experimental step response sample, N is the number of the step response samples. The indicator has been calculated for measurement data obtained from step responses at supply voltages ranging from 5 to 16 V. The more accurately the step response is reproduced, the smaller the indicator. The value of the RMSE for the current and speed is shown in Figure 14. As can be seen, the application of the CPM-GWO method reduces the error by about 50% when compared to the CPM method. As shown in Figure 5 the GWO method used solely has great difficulties in estimating parameter values of the PMDC model. The GWO achieves significantly worse accuracy than the CPM-GWO, and stops improving just after the 15th iteration.

Conclusions
In the current paper, the CPM-GWO method has been proposed and evaluated in experimental tests. It has been demonstrated that the proposed method offers significantly lower matching error when compared to the CPM and GWO methods used solely. More-over, the proposed method provides much lower computational time required to estimate the values of model parameters when compared to the GWO method.
The GWO method requires a priori information to provide correct and accurate estimates of the parameters. In the CPM-GWO method, this information comes from the initial estimates of the parameters achieved by the CPM, making the proposed method present correct, accurate, and reliable estimates.
The manufacturer provides values of nominal parameters that average the characteristics of the electric motors produced. Differences in the parameter values obtained using the presented CPM-GWO method in relation to the values declared by the manufacturer may result from differences arising at the production stage, wear of motor components (e.g., brushes, bearings).
The conducted tests have shown that friction coefficient b depends on the angular speed ω of the motor, and that this relation is in good agreement with Stribeck friction theory. Due to the observed variability of parameter b, it is worth including a more accurate friction model in future research related to this subject. In addition, future research could focus on more complex applications, such as when the motor is running under a constant or variable load.