The Use of Evolutionary Methods for the Determination of a DC Motor and Drive Parameters Based on the Current and Angular Speed Response

Abstract

Determination of the seven parameters of a Direct Current (DC) motor and drive is presented, based on the speed and current step responses. The method is extended for the motor and drive parameter determination in the case of a controlled drive. The influence of a speed controller on the responses is considered in the motor model with the use of the measured voltage. Current limitation of the supply unit is also considered in the DC motor model. For parameter determination, a motor model is used, which is determined with two coupled differential equations. Euler’s first-order and Runge–Kutta fourth-order methods are used for the motor model simulations. For parameter determination, evolutionary methods are used and compared to each other. Methods used are Genetic Algorithm, Differential Evolutions with two strategies, Teaching–Learning-Based Optimization, and Artificial Bee Colony. To improve results, deviation of the motor model simulation time is used and Memory Assistance with three different approaches is analyzed to shorten the calculation time. The tests showed that Differential Evolution (DE)/rand/1/exp is the most appropriate for the presented problem. The division of the motor model simulation time improves the results. For the presented problem, short-term memory assistance can be suggested for calculation time reduction.

1. Introduction

Drives with Direct Current motors (DC motors) are used widely in industrial applications. DC motors are easy to control and model, and therefore they are often used for industry control systems. The problem is that often, DC motor and a drive parameters are not known. Without known parameters, a precise model of the DC motor cannot be made, although it may be needed, among other things, for analytical control system design and optimization. Motor parameters, given in the motor specifications provided by the motor manufacturer, could have relatively large tolerances in their electric and mechanical parameters, especially for cheaper DC motors, and often not all parameters are given. Some of them can be measured, but appropriate equipment is needed, and for some measurements the coils should be disconnected, which can be a problem.

Also, drive parameters, such as inertia of the drive and friction, are important parameters in dimensioning of the drive control system. Inertia of the drive influences the start, stop times, and times needed for speed changes. Also, energies needed for transient phenomena are strongly dependent on the inertia. Friction means losses, and it is welcome to evaluate the friction of the motor or of the drive. Usually the inertia and friction of all the drive parts are not known, and because of that, they cannot be determined based on the manufacturer data.

In the past, many authors tried to determine DC motor parameters using different methods. Some of these approaches are: The use of the curve fitting approach [1], the use of the constraint optimization technique [2], the pattern search and last square algorithm [3], parameter identification via the moments method [4], parameter identification via inverse problem methodology using the conjugate gradients method and Tikhonov’s regularization method [5], parameter determination using the regression method [6], etc. Also, evolutionary methods are used for parameter estimation [7,8,9,10].

Determination of seven parameters of the DC motor and drive, based on the speed and current step responses, is presented in our work. These parameters are: Three DC motor parameters, inertia of the motor or drive, and three parameters describing the friction of the motor or drive.

The input data for the parameter determination measured speed and current responses. Measured responses are compared to calculated responses, obtained using the DC motor simulation. Simulation is made based on numerical solving of a differential equations system which describes the DC motor drive. Two methods are tested for numerical solving of the differential equations system, which are the simple first-order Euler method and fourth-order Runge–Kutta method. Parameters used as input for simulation should be selected in such a way that simulated responses are as similar as possible to measured responses. We are dealing with an inverse problem, which is also an optimization problem; we are searching for the minimum difference between measured and simulated responses. We decided to use evolutionary methods. Different methods behave differently on different problems, therefore, we tested different methods with the aim to find the most appropriate for the presented problem. We tested Genetic Algorithm (GA) [11,12,13,14,15,16], Differential Evolution (DE) [17,18,19,20,21,22,23,24,25,26] with two different strategies, which are DE/rand/1/exp and DE/best/1/bin, Teaching–Learning-Based Optimization (TLBO) [27,28,29,30,31,32,33,34,35] and Artificial Bee Colony (ABC) [36,37,38,39,40,41,42,43,44,45]. To evaluate the quality of the methods better, tests were also made for simulated input data (all parameters are known), and not only for measured input (mostly exact values of the parameters are not known).

Our contributions in this work are:

• Determination of seven DC motor and drive parameters only based on speed and current responses. We tested different evolutionary methods, while other authors used only one standard or evolutionary method. The most appropriate evolutionary method between those selected is proposed based on tests.
• Different parameters can be determined with the presented approach: Only the motor (parameters of the motor and inertia and friction of the motor), DC motor drive without load (parameters of the motor and inertia and friction of the drive), DC motor drive with the load (parameters of the motor and inertia of the drive and load characteristic of the drive).
• The method is also extended for the motor and drive parameter determination in the case of a controlled drive. The influence of a speed controller on the current and speed responses is considered in the DC motor model with the use of the measured voltage. Also, current limitation of the supply unit is considered in the DC motor model. Such an approach has not been found in the literature, as mostly authors use step responses of speed and current, or, in the case of controlled drive, they use a model of the controller.
• The accuracy of the DC motor model, used for the Objective Function calculation, is improved by dividing the time step of the measurement, which cannot be reduced, into smaller steps in the Objective Function calculation. Such approach to improve the accuracy of the calculation was not found in the literature.
• To reduce calculation time short-term memory assistance (STMA) is used, which records not only the current population, but also population from the previous calculation step. Also presented is long-term memory assistance (LTMA), which records the entire search history. In the case of LTMA two different search strategies are used for the search in history.

Some other authors tried to determine DC motor parameters with some of the evolutionary methods. In [7] Sankardos et al. used a GA for parameter determination of the permanent magnet DC motor. Lankarani et al. also made DC motor parameter estimation using a GA [8]. In [9] Dupois et al. proposed a multi-objective elitist GA for DC motor parameter estimation. Udomsuk used in [10] the adaptive tabu search technique for separately excited DC motor parameter identification. None of them used measured voltage for simulation of the controlled drive, as we did.

The paper consists of five sections. The principle of parameter determination is described in the second section. Methods are also presented which were used for simulation of the motor model and the evolutionary methods used for parameter determination. Results obtained using simulated data as input are presented in the third section. Results obtained using measured data as input are presented in the fourth section. Cross-validation, testing of result correctness at stationary points, the influence of different weights in Objective Function, improvement of the motor model simulations using time division and analysis of Memory Assistance are also shown. Conclusions are given in the last, fifth section.

2. Principle of Parameter Determination

The aim of the presented procedure is determination of the DC motor and drive parameters (which are often not given, or they are not precise enough) based on the input data, which are current and speed step responses. The input data are simulated current and speed responses that are used for analysis of the solving procedure and solving methods. Or, measured current and speed responses that are used for determination of the real DC motor and drive parameters. The parameter determination principle is based on the comparison of the input data (current and speed responses) and simulated data, obtained with DC motor model. In the continuation, the procedure for determination of simulated data is presented as well as the selected evolutionary optimization methods.

2.1. DC Motor Model and Simulations

The DC motor is presented schematically in Figure 1, which also shows a working machine that may or may not be connected to the motor. It is dependent on the requirement of which parameters of the drive should be determined. The presented approach can be used for a DC motor with separate excitation, with parallel excitation, and for a DC motor with permanent magnets.

The values on Figure 1 are: R_a ohmic resistance of the DC motor, L_a inductance of the DC motor, u_a voltage at the DC motor, i_a current of the DC motor, e induced voltage at the DC motor armature coil, J_m inertia of the motor, ω angular speed at the axis of the motor, and J_wm inertia of the working machine.

Step responses of the motor can be simulated using two Differential Equations (1)–(7). The first Differential Equation describes the electrical subsystem and is written in (1), the second describes the mechanical subsystem and is written in (2).

$$u_a=i_a·R_a+L_a·\frac{d{i}_a}{dt}+e$$

(1)

$$T_m-T_{load}=J\frac{d\omega }{dt}$$

(2)

J is the equivalent inertia of all parts in the drive reduced to the motor speed, T_m is the torque of the motor and T_load is the torque of the load, which can only be friction, or friction and load of the working machine together.

Equations (1) and (2) are coupled, because induced voltage e depends on the speed of the motor, and the torque of the motor T_m depends on the current of the motor. e and T_m are written in (3) and (4).

$$e=c_m·\omega$$

(3)

$$T_m=c_m·i_a$$

(4)

c_m is assumed to be a constant (in the continuation it is called the motor constant), because c_m=k_m·Φ, where k_m is a constant of the DC motor and magnetic flux Φ is assumed to be a constant in the case of the considered motors. T_load written in (2) can be divided into more parts, presented in (5).

$$T_{load}=T_{\mathrm{la}}+T_{\mathrm{lb}}·\omega +T_{\mathrm{lc}}·{\omega }^2$$

(5)

In case no load from the working machine is present, factors T_la and T_lb represent Coulomb and viscous friction, and T_lc represents air resistance at the fan on the motor axis, if present. In case the working machine produces load, factors T_la, T_lb, and T_lc represent friction, air resistance at the fan on the motor axis and load of the working machine together. By insertion of (3), (4) and (5) into (1) and (2), (6) and (7) are obtained.

$$u_a=i_a·R_a+L_a·\frac{d{i}_a}{dt}+c_m·\omega$$

(6)

$$c_m·i_a-\left(T_{\mathrm{la}}+T_{\mathrm{lb}}·\omega +T_{\mathrm{lc}}·{\omega }^2\right)=J\frac{d\omega }{dt}$$

(7)

Based on Equations (6) and (7), we are searching for seven parameters. The meaning of the searched parameters depends on the three different test cases: (a) Only the DC motor is tested, (b) The DC motor and a working machine which does not produce a load are tested and (c) The DC motor and a working machine which produces load are tested. The meaning of the seven searched parameters, dependent on the test case, are presented in

Table 1. Meaning of the searched parameters dependent on the test case.

Parameter	Test Case (a) Only Motor	Test Case (b) Motor and Working Machine Which Does Not Produce Load	Test Case (c) Motor and Working Machine Which Produces Load
Ra (Ω)	motor resistance
La (H)	motor inductance
cm (Vs)	motor constant
J (kgm2)	motor inertia	drive inertia	drive inertia
Tla (Nm)	} friction of the motor	friction of the drive (motor + working machine)	common load (friction of the drive + load of the working machine)
Tlb (Nm·s)
Tlc (Nm·s2)

.

Simulation of the motor step responses is made with numerical solving of the differential equations written in (6) and (7). A simple first-order Euler method [46] and more complex fourth-order Runge–Kutta method [46] are chosen for numerical solving. Derivatives must be expressed for both methods, the derivative of current from (6) and derivative of speed from (7), and they are written in (8) and (9).

$$\frac{d{i}_a}{dt}=\frac{1}{L_a}·\left(u_a-i_a·R_a-c_m·\omega \right)=f\left(t,i_a,\omega \right)$$

(8)

$$\frac{d\omega }{dt}=\frac{1}{J}·\left[c_m·i_a-\left(T_{la}+T_{lb}·\omega +T_{lc}·{\omega }^2\right)\right]=g\left(t,i_a,\omega \right)$$

(9)

2.1.1. Simulation with the Use of Euler’s Method

Euler’s method is a first-order method for numerical solving of differential equations [46]. The basic equation of Euler’s method is written in (10).

$$y^{(k+1)}=y^{(k)}+h·f\left(x^{(k)},y^{(k)}\right)$$

(10)

y and z are searching functions, k is the value counter, h is a calculation step and f and g are expressed derivatives of functions y and z. Considering (8) and (9), which describe the problem, and considering Euler’s method (10), expressions (11) and (12) are written, which describe the solving of the given system.

$$i_a^{(k+1)}=i_a^{(k)}+\Delta t·f\left(t^{(k)},i_a^{(k)},{\omega }^{(k)}\right)=i_a^{(k)}+\Delta t·\left[\frac{1}{L_a}·\left(u_a-i_a^{(k)}·R_a-c_m·{\omega }^{(k)}\right)\right]$$

(11)

$$\begin{array}{ll}{\omega }^{(k+1)}& ={\omega }^{(k)}+\Delta t·g\left(t^{(k)},i_a^{(k)},{\omega }^{(k)}\right)\\ & ={\omega }^{(k)}+\Delta t·\left[\frac{1}{J}·\left[c_m·i_a^{(k)}-\left(T_{la}+T_{lb}·{\omega }^{(k)}+T_{lc}·{\left({\omega }^{(k)}\right)}^2\right)\right]\right]\end{array}$$

(12)

In (11) and (12) the values of i_a and ω in time instant k+1 are calculated, based on the previous time instant k. Δt is a calculation time step.

2.1.2. Simulation with the Use of Fourth-Order Runge–Kutta Method

Runge–Kutta’s method is a fourth-order method for numerical solving of differential equations [46]. The basic equations of the fourth-order Runge–Kutta method are written in (13).

$$\begin{array}{l}{K}_1=f\left(x^{(k)},y^{(k)}\right);K_2=f\left(x^{(k)}+\frac{1}{2}h,y^{(k)}+\frac{1}{2}h·K_1\right);\\ K_3=f\left(x^{(k)}+\frac{1}{2}h,y^{(k)}+\frac{1}{2}h·K_2\right);K_4=f\left(x^{(k)}+h,y^{(k)}+h·K_3\right);\\ y^{(k+1)}=y^{(k)}+h·\left(\frac{1}{6}{K}_1+\frac{1}{3}{K}_2+\frac{1}{3}{K}_3+\frac{1}{6}{K}_4\right);\end{array}$$

(13)

y and z are searching functions, k is the value counter, h is a calculation step and f and g are expressed derivatives of functions y and z. K₁ to K₄ and L₁ to L₄ are derivatives of y and z. Considering (8) and (9), which describe the problem, and considering the fourth-order Runge–Kutta method (13), expressions (14)-(23) are written, which describe the solving of the given system.

$$K_1=f\left(t^{(k)},i_a^{(k)},{\omega }^{(k)}\right)=\frac{1}{L_a}·\left(u_a-i_a^{(k)}·R_a-c_m·{\omega }^{(k)}\right)$$

(14)

$$L_1=g\left(t^{(k)},i_a^{(k)},{\omega }^{(k)}\right)=\frac{1}{J}·\left[c_m·i_a^{(k)}-\left(T_{la}+T_{lb}·{\omega }^{(k)}+T_{lc}·{\left({\omega }^{(k)}\right)}^2\right)\right]$$

(15)

$$\begin{array}{ll}{K}_2& =f\left(t^{(k)}+\frac{1}{2}\Delta t,i_a^{(k)}+\frac{1}{2}\Delta t·K_1,{\omega }^{(k)}+\frac{1}{2}\Delta t·L_1\right)\\ & =\frac{1}{L_a}·\left(u_a-\left(i_a^{(k)}+\frac{1}{2}\Delta t·K_1\right)·R_a-c_m·\left({\omega }^{(k)}+\frac{1}{2}\Delta t·L_1\right)\right)\end{array}$$

(16)

$$\begin{array}{ll}{L}_2& =g\left(t^{(k)}+\frac{1}{2}\Delta t,i_a^{(k)}+\frac{1}{2}\Delta t·K_1,{\omega }^{(k)}+\frac{1}{2}\Delta t·L_1\right)\\ & =\frac{1}{J}·\left[c_m·\left(i_a^{(k)}+\frac{1}{2}\Delta t·K_1\right)-\left(T_{la}+T_{lb}·\left({\omega }^{(k)}+\frac{1}{2}\Delta t·L_1\right)+T_{lc}·{\left({\omega }^{(k)}+\frac{1}{2}\Delta t·L_1\right)}^2\right)\right]\end{array}$$

(17)

$$\begin{array}{ll}{K}_3& =f\left(t^{(k)}+\frac{1}{2}\Delta t,i_a^{(k)}+\frac{1}{2}\Delta t·K_2,{\omega }^{(k)}+\frac{1}{2}\Delta t·L_2\right)\\ & =\frac{1}{L_a}·\left(u_a-\left(i_a^{(k)}+\frac{1}{2}\Delta t·K_2\right)·R_a-c_m·\left({\omega }^{(k)}+\frac{1}{2}\Delta t·L_2\right)\right)\end{array}$$

(18)

$$\begin{array}{ll}{L}_3& =g\left(t^{(k)}+\frac{1}{2}\Delta t,i_a^{(k)}+\frac{1}{2}\Delta t·K_2,{\omega }^{(k)}+\frac{1}{2}\Delta t·L_2\right)\\ & =\frac{1}{J}·\left[c_m·\left(i_a^{(k)}+\frac{1}{2}\Delta t·K_2\right)-\left(T_{la}+T_{lb}·\left({\omega }^{(k)}+\frac{1}{2}\Delta t·L_2\right)+T_{lc}·{\left({\omega }^{(k)}+\frac{1}{2}\Delta t·L_2\right)}^2\right)\right]\end{array}$$

(19)

$$\begin{array}{ll}{K}_4& =f\left(t^{(k)}+\Delta t,i_a^{(k)}+\Delta t·K_3,{\omega }^{(k)}+\Delta t·L_3\right)\\ & =\frac{1}{L_a}·\left(u_a-\left(i_a^{(k)}+\Delta t·K_3\right)·R_a-c_m·\left({\omega }^{(k)}+\Delta t·L_3\right)\right)\end{array}$$

(20)

$$\begin{array}{ll}{L}_4& =g\left(t^{(k)}+\Delta t,i_a^{(k)}+\Delta t·K_3,{\omega }^{(k)}+\Delta t·L_3\right)\\ & =\frac{1}{J}·\left[c_m·\left(i_a^{(k)}+\Delta t·K_3\right)-\left(T_{la}+T_{lb}·\left({\omega }^{(k)}+\Delta t·L_3\right)+T_{lc}·{\left({\omega }^{(k)}+\Delta t·L_3\right)}^2\right)\right]\end{array}$$

(21)

$$i_a^{(k+1)}=i_a^{(k)}+\Delta t·\left(\frac{1}{6}{K}_1+\frac{1}{3}{K}_2+\frac{1}{3}{K}_3+\frac{1}{6}{K}_4\right)$$

(22)

$${\omega }_a^{(k+1)}={\omega }_a^{(k)}+\Delta t·\left(\frac{1}{6}{L}_1+\frac{1}{3}{L}_2+\frac{1}{3}{L}_3+\frac{1}{6}{L}_4\right)$$

(23)

In (22) and (23), the values of i_a and ω in time instant k+1 are calculated based on the time instant k. Δt is a calculation time step. In the case of the fourth-order Runge–Kutta method, many more calculations are needed to get i_a and ω at a current time instant, as in the case of Euler’s method.

The discussed problem is an optimization problem. The simulated data must be as similar as possible to the input data (measured or simulated). The Objective Function (OF) is defined as the square of differences between the input and simulated data of current and speed, and it is written in (24).

$$\mathrm{OF}=\frac{1}{N}{\displaystyle \sum _{i=1}^N\left({\left(\frac{i_{a\_simulated\_i}-i_{a\_input\_i}}{i_{a\_input\_\max }}\right)}^2+{\left(\frac{{\omega }_{simulated\_i}-{\omega }_{input\_i}}{{\omega }_{input\_\max }}\right)}^2\right)}$$

(24)

N is the number of points in the step responses. OF is divided by N to have the possibility to compare the same response measured with different numbers of points. Current and speed are normalized, divided with the maximum input value, to eliminate the different influence of current and speed on the OF due to the different size of current and speed data.

2.2. Methods for Parameter Determination

Evolutionary optimization methods [11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45] are used for parameter determination. They are population-based search algorithms. The space of all possible solutions is searched through many solutions, which have various suitability (fitness). Better solutions have higher probability to be transferred into the next generation, where they might be changed by crossover and/or mutation. Evolutionary optimization methods have the very exceptional ability to balance exploitation and exploration [15], and, with that, the ability to avoid local optima. They are also suitable for solving technical problems, such as to determine the parameters of different models describing real material properties [47,48,49], or also determine the parameters of models describing real devices. In the presented work, we used and compared GA [11,12,13,14,15,16], DE [17,18,19,20,21,22,23,24,25,26], TLBO [27,28,29,30,31,32,33,34,35] and ABC [36,37,38,39,40,41,42,43,44,45]. In GA, the mutation rate used was 0.2, and the used fraction of population kept was 0.5. For both DE strategies, which are DE/rand/1/exp and DE/best/1/bin, the used amplification of the differential variation was 0.6, and the used crossover probability was 0.8. Two different strategies were used, because a strategy that is the best for a given problem, may not be the best when applied for a different problem. The duplicate elimination phase was omitted in the TLBO used in our work. Therefore, the number of Fitness Evaluations (FEs) consumed was determined statically as: FEs = 2 × population size × iterations. The limit value, which is a control parameter for the bee population in ABC, was set at 100. Since in ABC a scout bee might not be employed for every iteration, the number of FEs cannot be determined statically [39].

In this work, we are dealing with seven parameters (D = 7). The population size is ten times the number of parameters, which is 70 for seven parameters (NP = 70). For stopping criterion for all tests, we used 140,000 FEs, which is 2000 iterations in the case of GA and DE, 1000 iterations in the case of TLBO, and less or equal to 2000 iterations in the case of ABC. The OF used for the quality estimation of simulated response, compared to input data/response, was written previously in (24).

3. Simulated Data Used as Input

To compare the efficiency and efficacy of the selected methods: GA, DE/rand/1/exp, DE/best/1/bin, TLBO, and ABC, simulated data were used as input data. Simulated input data were obtained using the motor model presented in Section 2.1, described with coupled Differential Equations (8) and (9). To perform the simulation the fourth-order Runge–Kutta method was used, described with expressions (14)–(23), presented in Section 2.1.2. For the first tests, simulated input data were chosen, due to the following advantages compared to the real measured data:

• The weaknesses of the model are avoided. For example, L_a is not a completely constant value, but it is considered to be a constant value.
• The phenomena, which are not covered by the model, such as armature reaction, is eliminated from the input data.
• The exact values of all seven parameters, which should be obtained as results, are known, because they are used as input for preparation of the simulated current and speed step responses.
• The same method (the fourth-order Runge–Kutta method) is used to prepare simulated input data and to evaluate the OF (24) during the parameter determination process. With that, the influence of the selected method on the results is eliminated.

The data used for the simulated input are presented in

Table 2. Data used for the simulated input.

Voltage ua	Start Time tstart	End Time tend	Time Step Δt	Initial Current ia (t = tstart)	Initial Speed ω (t = tstart)
220 V	0	0.05 s	1 × 10−4 s	0 A	0 s−1

.

Based on t_start, t_end and Δt, presented in Table 2, it can be seen that i_a and ω were calculated for 500 time instants. The parameters used to get simulated data are presented in

Table 3. Parameters used to get simulated data.

	Simulated Input Data
Parameter	SD1	SD2	SD3	SD4
Ra (Ω)	42.5	42.5	42.5	42.5
La (H)	0.08	0.08	0.008	0.08
cm (Vs)	0.4781	0.4781	0.4781	0.4781
J (kgm2)	2·10−5	6 × 10−5	2 × 10−5	2 × 10−6
Tla (Nm)	0.01	0.01	0.01	0.01
Tlb (Nm·s)	3.27 × 10−5	3.27 × 10−5	3.27 × 10−5	3.27 × 10−5
Tlc (Nm·s2)	8.55 × 10−8	8.55 × 10−8	8.55 × 10−8	8.55 × 10−8

. Four different sets of input data, current, and speed step responses were made. They are marked as simulated data sets 1 to 4 (SD1–SD4).

For SD2, bigger inertia was used, for SD3, lower L_a, and with that, lower electrical time constant was used and for SD3 smaller inertia.

The used simulated data responses are presented in Figure 2.

The limits of the parameters were set due to the physical meaning of the parameters. The used limits are presented in

Table 4. Limits of the parameters.

Parameter	Lower Limit	Upper Limit
Ra (Ω)	0	100
La (H)	0	1
cm (Vs)	0	5
J (kgm2)	0	1
Tla (Nm)	0	1
Tlb (Nm·s)	0	1 × 10−3
Tlc (Nm·s2)	0	1 × 10−6

.

3.1. Results Obtained Using Simulated Input Data

Calculations of the DC motor and drive parameters were made for all four simulated input data current and speed responses, which are SD1-SD4, presented in Figure 2. For parameter determination, all five selected methods were used, which are GA, DE/rand/1/exp, DE/best/1/bin, TLBO, and ABC. Due to the stochastic behavior of the evolutionary methods, 50 independent runs are made for each combination of input data/used method. The Best (B), Worst (W), Mean (M) OF, and Standard Deviation (SD) are written in

Table 5. OF B, W, M value, and SD for 50 independent runs of different methods.

Expression				Method
		GA	DE/Rand/	DE/Best/	TLBO	ABC
			1/exp	1/bin
	B	3.2054 × 10−4	4.8980 × 10−19	4.8980 × 10−19	3.9564 × 10−15	6.8520 × 10−10
SD1	W	4.3669 × 10−1	4.8980 × 10−19	2.5072 × 10−2	2.7757 × 10−10	1.6968 × 10−7
	M	1.0765 × 10−1	4.8980 × 10−19	6.6203 × 10−4	1.6827 × 10−11	3.4562 × 10−8
	SD	1.1031 × 10−1	2.0882 × 10−27	3.6637 × 10−3	4.7340 × 10−11	4.3330 × 10−8
	B	7.8123 × 10−5	6.2556 × 10−19	6.2556 × 10−19	7.8460 × 10−13	1.0932 × 10−9
SD2	W	1.8151 × 10−1	6.2556 × 10−19	4.6771 × 10−2	7.3488 × 10−9	2.1087 × 10−7
	M	5.2484 × 10−2	6.2556 × 10−19	9.3542 × 10−4	5.1639 × 10−10	4.1001 × 10−8
	SD	5.0508 × 10−2	3.5547 × 10−27	6.5479 × 10−3	1.0976 × 10−9	4.5495 × 10−8
	B	6.9449 × 10−5	4.6666 × 10−19	4.6666 × 10−19	1.9151 × 10−12	6.1774 × 10−11
SD3	W	2.8716 × 10−1	4.6666 × 10−19	4.6005 × 10−2	2.6491 × 10−9	2.1340 × 10−9
	M	7.4399 × 10−2	4.6666 × 10−19	2.3516 × 10−3	3.3533 × 10−10	5.2114 × 10−10
	SD	6.716 × 10−2	2.2635 × 10−27	9.1419 × 10−3	5.0565 × 10−10	4.4003 × 10−10
	B	3.6061 × 10−3	3.0248 × 10−19	3.0248 × 10−19	4.523 × 10−16	9.7925 × 10−9
SD4	W	2.8045 × 10−1	3.0248 × 10−19	8.7846 × 10−2	4.2365 × 10−9	4.1963 × 10−7
	M	1.3085 × 10−1	3.0248 × 10−19	4.4826 × 10−3	1.5931 × 10−10	1.3155 × 10−7
	SD	8.1103 × 10−2	1.1859 × 10−27	1.3509 × 10−2	6.8381 × 10−10	9.5220 × 10−8

. In Table 5 the lowest mean OF for each simulated input data is marked with bold, and the second lowest mean OF for each simulated input data is marked with bold italics.

Based on the results presented in Table 5, it can be seen that for all four tested simulated input data, DE/rand/1/exp was much better than the other methods. The second best, also for all tested simulated input data, was TLBO. The mean values of calculated parameters using all five calculation methods are presented in

Table 6. Mean values of calculated parameters for SD1.

Parameter				Method
		Used Value	GA	DE/Rand/	DE/Best/	TLBO	ABC
		(Simulation)		1/exp	1/bin
Ra (Ω)	M	42.5	37.34	42.5	43.58	42.5	42.48
La (H)	M	8 × 10−2	5.692 × 10−2	8 × 10−2	8.116 × 10−2	8 × 10−2	8.002 × 10−2
ce (Vs)	M	0.4781	0.6782	0.4781	0.4773	0.4781	0.4781
J (kgm2)	M	2 × 10−5	8.302 × 10−5	2 × 10−5	1.978 × 10−5	2 × 10−5	2.004 × 10−5
Tla (Nm)	M	1 × 10−2	1.358 × 10−2	1 × 10−2	0.981 × 10−2	0.998 × 10−2	0.736 × 10−2
Tlb (Nms)	M	3.27 × 10−5	2.033 × 10−5	3.270 × 10−5	3.026 × 10−5	3.260 × 10−5	2.577 × 10−5
Tlc (Nms2)	M	8.55 × 10−8	4.778 × 10−8	8.55 × 10−8	10.98 × 10−8	8.850 × 10−8	11.34 × 10−8

for SD1.

Based on Table 6, following conclusions can be derived:

• Mean values of parameters calculated using DE/rand/1/exp are the same as the values used to generate input data. These results are the best.
• Mean values of the parameters calculated using TLBO are almost the same as the values used to generate input data. The small difference is only for the parameters describing load, which are T_la, T_lb, and Tlc.
• Mean values of the parameters calculated using GA, DE/best/1/bin and ABC are more or less different than the values used to generate input data.

Similar results considering the mean values of calculated parameters were also obtained in the case of SD2, SD3, and SD4. They are not shown, due to the close similarity to the case of SD1 and to make this paper short and concise.

Results of simulations confirmed that evolutionary methods are appropriate for solving of the presented problem.

4. Measured Data Used as an Input and Improvement of the Method for Controlled Drive and Current Limitation

The use of current and speed step responses is limited to small motors supplied by powerful supply units, due to the large start-up currents of the DC motors. Usually supply units are not powerful enough for direct start of the motor, and controlled drives are often used, in which the motor speed is controlled. The combined drive of the AC/DC supply unit, DC motor, pulse encoder, and working machine is presented in Figure 3.

To consider a controlled drive we did not add a controller into the model of the motor presented in Section 2.1 as for example in [7]. The controller influences voltage at the motor u_a, which is no longer a constant. To consider that in the case of Euler’s method (EM) presented in Section 2.1.1, the u_a used in (11) is replaced with measured u_a(t). In the case of the Runge–Kutta method (RK) presented in Section 2.1.2, the u_a used in (14), (16), (18) and (20) is replaced with measured u_a (t). Also, after each calculated i_a, it should be checked, and, if needed, corrected to a limit value of i_{a_limit}. Algorithms for parameter determination using EM and RK are presented in Figure 4.

4.1. Measured Data Used as a Test Example

For measured data, a laboratory drive was used, which is presented in Figure 5.

The drive presented in Figure 5 is combined of the following parts:

• Supply unit: SIEMENS SIMOREG DC-Master 6RA7013-6DV62-0-Z.
• DC motor: SIEMENS 1GG5104-0ED40-6VV1.
• Pulse encoder: HUBNER Berlin, P0G 9D 1024.
• Load simulation motor: SIEMENS 1LA7139-4AA10-Z FDB0
• Measurements were made with the use of a “Trace” function, which is a part of the SIEMENS “Drive Monitor” software used for supply unit support.

The measurement data must be prepared properly. The excitation should be switched on before the armature circuit. In this way, it is not necessary to consider the transient phenomenon when switching on the excitation. Only the part of the measured data after switching the armature on is used for the parameter determination, as presented in Figure 6.

Four different measured input data sets were made, marked as MD1, MD2, MD3, and MD4. Since it was a controlled drive, the following values had to be defined: Speed up time t_{speed_up}, final speed ω_final and current limit i_{a_limit}. The values used for measured input data sets are presented in

Table 7. Values used to get measured data.

	Measured Input Data
Value	MD1	MD2	MD3	MD4
tspeed_up (s)	0	0	0	2
ωfinal (s−1)	126	182	182	182
ia_limit (A)	11.44 (110% of Ia_rated)	11.44 (110% of Ia_rated)	14.56 (140% of Ia_rated)	14.56 (140% of Ia_rated)
Load (Nm)	no load	no load	no load	≈ 7.53 × 10−2·ω
Number of measured points	87	81	65	106

.

The measurement time step was 0.132 s. The measured voltage u_a(t), current i_a(t) and speed ω(t) are presented in Figure 7 for all four measured data sets.

The limits of the parameters were set due to the physical meaning of the parameters. We set limits appropriate for different motor sizes and loads, but it is always good to check if the limits are appropriate for the size of the drive. The used limits are presented in

Table 8. Limits of the parameters.

Parameter	Lower Limit	Upper Limit
Ra (Ω)	0	100
La (H)	0	100
cm (Vs)	0	5
J (kg·m2)	0	1
Tla (Nm)	0	20
Tlb (Nm·s)	0	9.55 × 10−2 (torque of 20 Nm at speed 209 s−1—motor rated speed is 182 s−1)
Tlc (Nm·s2)	0	4.56 × 10−6 (torque of 20 Nm at speed 209 s−1—motor rated speed is 182 s−1)

.

4.2. Results Obtained Using the Measured Input Data

Calculation results are presented for each measured input data set MD1-MD4 in separate Section 4.2.1, Section 4.2.2, Section 4.2.3 and Section 4.2.4, and results using input data sets are compared in Section 4.2.5. Comparison between different data makes sense, because they are obtained on the same drive, and it is expected to get the same results for both motor and drive parameters. Cross-validation is presented in Section 4.2.6, testing of result correctness at stationary points in Section 4.2.7, calculation times are presented in Section 4.2.8, and the influence of weights added into OF is analyzed in Section 4.2.9.

4.2.1. Results Obtained Using Input Data MD1

In the case of MD1, parameters were calculated using EM and RK. Results using different evolutionary methods, presented in Section 2.2, were compared, to find out which of the methods was the best for the presented problem. OF and mean values of parameters calculated using EM are presented in

Table 9. OF and Mean value of the calculated parameters for 50 independent runs using EM for MD1.

OF and				Method
Parameters		Known Value	GA	DE/rand/1/exp	DE/best/1/bin	TLBO	ABC
	B	-	3.0377 × 10−3	2.8906 × 10−3	2.8906 × 10−3	2.8906 × 10−3	2.8906 × 10−3
OF	W	-	3.8698 × 10−2	2.8906 × 10−3	3.0919 × 10−2	2.6806 × 10−2	2.8907 × 10−3
	M	-	1.0227 × 10−2	2.8906 × 10−3	3.4514 × 10−3	3.3694 × 10−3	2.8906 × 10−3
	SD	-	8.3375 × 10−3	0.0	3.9239 × 10−3	3.3480 × 10−3	1.5865 × 10−8
Ra (Ω)	M	5.66	18.54	11.56	11.33	12.92	11.56
La (H)	M	not known	3.53	9.72 × 10−1	1.95	1.12	9.72 × 10−1
cm (Vs)	M	not known	1.288	1.337	1.338	1.324	1.337
J (kgm2)	M	≈ 4 × 10−2	2.83 × 10−2	4.64 × 10−2	4.55 × 10−2	4.55 × 10−2	4.64 × 10−2
Tla (Nm)	M	≈ 0.9	1.79 × 10−1	3.13 × 10−16	0.0	3.29 × 10−5	7.06 × 10−18
Tlb (Nms)	M	≈ 0	2.83 × 10−3	8.95 × 10−17	4.24 × 10−4	4.95 × 10−4	6.27 × 10−6
Tlc (Nms2)	M	≈ 0	2.39 × 10−5	4.99 × 10−5	4.71 × 10−5	4.60 × 10−5	4.99 × 10−5

, and using RK in

Table 10. OF and Mean value of the calculated parameters for 50 independent runs using RK for MD1.

OF and				Method
Parameters		Known	GA	DE/rand/1/exp	DE/best/1/bin	TLBO	ABC
	B	-	2.6083 × 10−3	2.2531 × 10−3	2.2531 × 10−3	2.2531 × 10−3	2.2531 × 10−3
OF	W	-	1.7874 × 10−2	2.2531 × 10−3	3.3070 × 10−3	2.2531 × 10−3	2.2555 × 10−3
	M	-	9.7051 × 10−3	2.2531 × 10−3	2.3423 × 10−3	2.2531 × 10−3	2.2534 × 10−3
	SD	-	4.2762 × 10−3	0.0	28456 × 10−4	0.0	4.2741 × 10−7
Ra (Ω)	M	5.66	7.12	5.06	5.41	5.06	5.07
La (H)	M	not known	9.98 × 10−1	2.44 × 10−1	2.50 × 10−1	2.44 × 10−1	2.44 × 10−1
cm (Vs)	M	not known	1.364	1.369	1.367	1.369	1.369
J (kgm2)	M	≈ 4 × 10−2	3.00 × 10−2	4.68 × 10−2	4.63 × 10−2	4.68 × 10−2	4.68 × 10−2
Tla (Nm)	M	≈ 0.9	3.64 × 10−1	7.99 × 10−1	5.11 × 10−1	7.99 × 10−1	7.86 × 10−1
Tlb (Nms)	M	≈ 0	1.93 × 10−3	7.70 × 10−18	2.87 × 10−21	2.71 × 10−17	9.06 × 10−7
Tlc (Nms2)	M	≈ 0	1.76 × 10−3	8.00 × 10−19	1.84 × 10−5	1.43 × 10−18	7.81 × 10−7

.

Considering Table 9 and Table 10, the best results were obtained using DE/rand/1/exp in the case of EM, and using DE/rand/1/exp and TLBO in the case of RK, where the calculated values for all 50 starts were the same (SD = 0). Results obtained using EM were not correct, considering known values, because EM was not precise enough. In the case of DE/rand/1/exp deviation of calculated R_a was 104%, of J was 16%, and calculated T_la was 0 (it should be around 0.9). Calculated values using RK were correct. In the case of DE/rand/1/exp deviation of calculated R_a was 11%, of J was 17% and of T_la was 11%. The measured and calculated responses (using parameters obtained with DE/rand/1/exp written in Table 9 and Table 10) are presented in Figure 8. To make the Figure more transparent, only the first two seconds of the responses are displayed, and visible differences are marked between calculated responses using EM and RK.

Also, from Figure 8, it can be seen that differences between measured and calculated responses are bigger in the case of EM than in the case of RK.

4.2.2. Results Obtained Using Input Data MD2

OF and mean values of calculated parameters for MD2 using EM are presented in

Table 11. OF and Mean value of the calculated parameters for 50 independent runs using EM for MD2.

OF and				Method
Parameters		Known Values	GA	DE/rand/1/exp	DE/best/1/bin	TLBO	ABC
	B	-	3.6370 × 10−3	3.4006 × 10−3	3.4006 × 10−3	3.4006 × 10−3	3.4006 × 10−3
OF	W	-	1.8426 × 10−2	3.4006 × 10−3	3.4056 × 10−3	3.4006 × 10−3	3.4007 × 10−3
	M	-	7.9539 × 10−3	3.4006 × 10−3	3.4008 × 10−3	3.4006 × 10−3	3.4006 × 10−3
	SD	-	3.6910 × 10−3	1.7347 × 10−18	9.7750 × 10−7	1.7347 × 10−18	5.9380 × 10−9
Ra (Ω)	M	5.66	14.96	11.45	11.45	11.45	11.45
La (H)	M	not known	1.69	9.76 × 10−1	9.76 × 10−1	9.76 × 10−1	9.76 × 10−1
cm (Vs)	M	not known	1.331	1.350	1.350	1.350	1.350
J (kgm2)	M	≈ 4 × 10−2	3.87 × 10−2	4.89 × 10−2	4.89 × 10−2	4.89 × 10−2	4.89 × 10−2
Tla (Nm)	M	≈ 0.9	2.66 × 10−1	3.58 × 10−16	0.0	6.75 × 10−16	0.0
Tlb (Nms)	M	≈ 0	2.00 × 10−3	9.61 × 10−17	2.02 × 10−4	3.44 × 10−16	7.78 × 10−7
Tlc (Nms2)	M	≈ 0	1.13 × 10−5	2.77 × 10−5	2.66 × 10−5	2.77 × 10−5	2.77 × 10−5

, and using RK in

Table 12. OF and Mean value of the calculated parameters for 50 independent runs using RK for MD2.

OF and				Method
Parameters		Known Value	GA	DE/rand/1/exp	DE/best/1/bin	TLBO	ABC
	B	-	2.9697 × 10−3	2.6384 × 10−3	2.6384 × 10−3	2.6384 × 10−3	2.6384 × 10−3
OF	W	-	3.6245 × 10−3	2.6384 × 10−3	3.9111 × 10−3	2.6384 × 10−3	2.6433 × 10−3
	M	-	1.0895 × 10−2	2.6384 × 10−3	2.7298 × 10−3	2.6384 × 10−3	2.6389 × 10−3
	SD	-	6.1378 × 10−3	1.3010 × 10−18	3.0860 × 10−4	1.3010 × 10−18	1.1636 × 10−6
Ra (Ω)	M	5.66	7.04	4.97	5.20	4.97	4.97
La (H)	M	not known	9.96 × 10−1	2.46 × 10−1	2.49 × 10−1	2.46 × 10−1	2.46 × 10−1
cm (Vs)	M	not known	1.375	1.376	1.374	1.376	1.376
J (kgm2)	M	≈ 4 × 10−2	3.70 × 10−2	5.00 × 10−2	4.91 × 10−2	5.00 × 10−2	5.00 × 10−2
Tla (Nm)	M	≈ 0.9	3.94 × 10−1	9.34 × 10−1	7.47 × 10−1	9.34 × 10−1	9.18 × 10−1
Tlb (Nms)	M	≈ 0	1.44 × 10−3	1.10 × 10−17	9.60 × 10−5	2.53 × 10−17	1.35 × 10−5
Tlc (Nms2)	M	≈ 0	1.0926 × 10−5	1.31 × 10−19	5.40 × 10−6	3.34 × 10−19	3.65 × 10−7

.

Considering Table 11 and Table 12, the best results were obtained using DE/rand/1/exp and TLBO in the case of EM and RK, where the calculated values for all 50 starts were the same. Results obtained using EM were not correct, considering the known values. In the case of DE/rand/1/exp deviation of calculated R_a was 102%, of J was 22% and calculated T_la was 0 (it should be around 0.9). Calculated values using RK were correct. In the case of DE/rand/1/exp deviation of calculated R_a was 12%, of J was 25% and of T_la was 3.8%. The first two seconds of measured and calculated responses (using the parameters obtained with DE/rand/1/exp written in Table 11 and Table 12) are presented in Figure 9.

Also, from Figure 9, it can be seen that differences between measured and calculated responses are bigger in the case of EM than in the case of RK.

4.2.3. Results Obtained Using Input Data MD3

OF and mean values of calculated parameters in case of MD3 are presented in

Table 13. OF and Mean value of the calculated parameters for 50 independent runs using EM for MD3.

OF and				Method
Parameters		Known Value	GA	DE/rand/1/exp	DE/best/1/bin	TLBO	ABC
	B	-	5.8927 × 10−3	5.7731 × 10−3	5.7731 × 10−3	5.7731 × 10−3	5.7731 × 10−3
OF	W	-	4.9909 × 10−2	5.7731 × 10−3	4.2488 × 10−2	5.7731 × 10−3	5.7763 × 10−3
	M	-	1.2440 × 10−2	5.7731 × 10−3	6.5079 × 10−3	5.7731 × 10−3	5.7732 × 10−3
	SD	-	8.0686 × 10−3	8.6736 × 10−19	5.1400 × 10−3	8.6736 × 10−19	4.4973 × 10−7
Ra (Ω)	M	5.66	16.42	11.92	11.68	11.92	11.92
La (H)	M	not known	2.95 × 10−1	1.28	1.84	1.28	1.28
cem(Vs)	M	not known	1.313	1.342	1.344	1.342	1.342
J (kgm2)	M	≈ 4 × 10−2	2.98 × 10−2	4.42 × 10−2	4.33 × 10−2	4.42 × 10−2	4.42 × 10−2
Tla (Nm)	M	≈ 0.9	3.12 × 10−1	6.45 × 10−16	0.0	2.03 × 10−15	3.88 × 10−4
Tlb (Nms)	M	≈ 0	2.26 × 10−3	9.91 × 10−18	1.13 × 10−4	6.82 × 10−17	5.60 × 10−6
Tlc (Nms2)	M	≈ 0	1.52 × 10−5	3.11 × 10−5	3.11 × 10−5	3.11 × 10−5	3.10 × 10−5

for EM and

Table 14. OF and Mean value of the calculated parameters for 50 independent runs using RK for MD3.

OF and				Method
Parameters		Known Value	GA	DE/rand/1/exp	DE/best/1/bin	TLBO	ABC
	B	-	4.0166 × 10−3	3.7766 × 10−3	3.7766 × 10−3	3.7766 × 10−3	3.7766 × 10−3
OF	W	-	2.1531 × 10−2	3.7766 × 10−3	3.8196 × 10−3	3.7766 × 10−3	3.7778 × 10−3
	M	-	7.7614 × 10−3	3.7766 × 10−3	3.7817 × 10−3	3.7766 × 10−3	3.7766 × 10−3
	SD	-	3.9337 × 10−3	0.0	1.3977 × 10−5	0.0	1.8574 × 10−7
Ra (Ω)	M	5.66	5.79	4.87	4.85	4.87	4.87
La (H)	M	not known	9.43 × 10−1	4.51 × 10−1	4.49 × 10−1	4.51 × 10−1	4.51 × 10−1
cm (Vs)	M	not known	1.371	1.374	1.374	1.374	1.374
J (kgm2)	M	≈ 4 × 10−2	3.66 × 10−2	4.44 × 10−2	4.47 × 10−2	4.44 × 10−2	4.44 × 10−2
Tla (Nm)	M	≈ 0.9	5.38 × 10−1	1.01	8.90 × 10−1	1.01 × 10−1	1.01 × 10−1
Tlb (Nms)	M	≈ 0	1.59 × 10−3	5.38 × 10−18	0.0	2.65 × 10−17	0.0
Tlc (Nms2)	M	≈ 0	1.13 × 10−5	6.59 × 10−20	3.62 × 10−6	3.60 × 10−19	3.66 × 10−8

for RK.

Considering Table 13 and Table 14, the best results were obtained using DE/rand/1/exp and TLBO in the case of EM and RK, where the calculated values for all 50 starts were the same. Results obtained using EM were not correct, considering the known values. In the case of DE/rand/1/exp deviation of calculated R_a was 111%, of J was 11% and calculated T_la was 0 (it should be around 0.9). Calculated values using RK were correct. In the case of DE/rand/1/exp deviation of calculated R_a was 14%, of J was 11% and of T_la was 12%. The first two seconds of measured and calculated responses (using the parameters obtained with DE/rand/1/exp written in Table 13 and Table 14) are presented in Figure 10.

Also, from Figure 10, it can be seen that the difference between measured and calculated responses is bigger in the case of EM than in the case of RK.

4.2.4. Results Obtained Using Input Data MD4

OF and mean values of the calculated parameters in the case of MD4 are presented in

Table 15. OF and Mean value of the calculated parameters for 50 independent runs using EM for MD4.

OF and				Method
Parameters		Known Value	GA	DE/rand/1/exp	DE/best/1/bin	TLBO	ABC
	B	-	3.1654 × 10−3	2.9627 × 10−3	2.9627 × 10−3	2.9627 × 10−3	2.9674 × 10−3
OF	W	-	1.0235 × 10−2	2.9627 × 10−3	2.9627 × 10−3	2.9627 × 10−3	3.2046 × 10−3
	M	-	4.9229 × 10−3	2.9627 × 10−3	2.9627 × 10−3	2.9627 × 10−3	3.0880 × 10−3
	SD	-	1.7022 × 10−3	4.3368 × 10−19	4.3368 × 10−19	4.3368 × 10−19	5.9401 × 10−5
Ra (Ω)	M	5.66	5.20	6.97	6.97	6.97	6.56
La (H)	M	not known	4.16	2.82	2.82	2.82	3.08
cm (Vs)	M	not known	1.351	1.251	1.251	1.251	1.271
J (kgm2)	M	≈ 4 × 10−2	9.80 × 10−2	1.19 × 10−1	1.19 × 10−1	1.19 × 10−1	1.13 × 10−1
Tla (Nm)	M	≈ 0.9	6.93 × 10−1	1.27 × 10−15	0.0	2.69 × 10−14	5.27 × 10−2
Tlb (Nms)	M	≈ 7.5 × 10−2	3.06 × 10−2	2.66 × 10−17	8.14 × 10−26	1.77 × 10−16	1.14 × 10−2
Tlc (Nms2)	M	≈ 0	1.36 × 10−4	2.91 × 10−4	2.91 × 10−4	2.91 × 10−4	2.33 × 10−4

for EM and

Table 16. OF and Mean value of the calculated parameters for 50 independent runs using RK for MD4.

OF and				Method
Parameters		Known Value	GA	DE/rand/1/exp	DE/best/1/bin	TLBO	ABC
	B	-	3.7368 × 10−3	3.2949 × 10−3	3.2949 × 10−3	3.2949 × 10−3	2.1111 × 10−3
OF	W	-	7.3916 × 10−3	3.2949 × 10−3	5.7046 × 10−3	3.2949 × 10−3	3.8064 × 10−3
	M	-	4.9182 × 10−3	3.2949 × 10−3	3.3464 × 10−3	3.2949 × 10−3	3.3975 × 10−3
	SD	-	9.3375 × 10−4	8.6736 × 10−19	3.3730 × 10−4	1.0294 × 10−9	3.1891 × 10−4
Ra (Ω)	M	5.66	3.04	5.38	5.29	5.38	4.98
La (H)	M	not known	3.73	3.46 × 10−1	4.26 × 10−1	3.46 × 10−1	4.00 × 10−1
cm (Vs)	M	not known	1.442	1.314	1.319	1.314	1.336
J (kgm2)	M	≈ 4 × 10−2	9.74 × 10−2	1.02 × 10−2	1.02 × 10−2	1.02 × 10−2	9.95 × 10−2
Tla (Nm)	M	≈ 0.9	1.67	7.85 × 10−1	7.22 × 10−1	7.86 × 10−1	1.369
Tlb (Nms)	M	≈ 7.5 × 10−2	2.87 × 10−2	5.13 × 10−2	5.19 × 10−2	5.13 × 10−2	4.12 × 10−2
Tlc (Nms2)	M	≈ 0	1.38 × 10−4	7.37 × 10−20	5.54 × 10−22	7.46 × 10−12	4.46 × 10−5

for RK.

Considering Table 15 and Table 16, the best results were obtained using DE/rand/1/exp, DE/best/1/bin and TLBO in the case of EM, and DE/rand/1/exp and TLBO in the case of RK, where the calculated values for all 50 starts were the same. Results obtained using EM were not correct, considering known values. In the case of DE/rand/1/exp deviation of calculated R_a was 23%, of J was 198% and calculated T_la was 0 (it should be around 0.9). Calculated values using RK were better. In the case of DE/rand/1/exp deviation of calculated R_a was 4.9%, of J was 155% and of T_la was 13%.

It is interesting to note that OF in the case of EM was lower than in the case of RK, although the parameters calculated using RK were better than the parameters calculated using EM. In the cases of MD1, MD2, and MD3 EM was not appropriate. In this case, the measurement deviations and errors were bigger for the RK method, which simulates a real motor better, cannot get as close to the measured signal as the EM method, which simulates a real motor worse. Also, 155% deviation of calculated inertia in the case of RK can be explained by measurement deviations.

The first four seconds of measured and calculated responses (using the parameters obtained with DE/rand/1/exp written in Table 15 and Table 16) are presented in Figure 11.

Also, from Figure 11 it can be seen that the differences between measured and calculated responses are different in cases of EM and RK.

4.2.5. Comparison of Results Obtained Using All Four Input Data Sets

As previously mentioned, comparison between different data makes sense, because they are obtained on the same drive, and it is expected to get the same results for both motor and drive parameters. With that, also the correctness can be checked of the unknown parameters. Mean values of the parameters calculated using DE/rand/1/exp and RK for MD1, MD2, MD3 are presented in

Table 17. Mean values of calculated parameters for 50 independent runs using DE/rand/1/exp and RK for all input data.

				Input Data
Parameters		Known Value	MD1	MD2	MD3	Known Value	MD4
Ra (Ω)	M	5.66	5.06	4.97	4.87	5.66	5.38
La (H)	M	not known	2.44 × 10−1	2.46 × 10−1	4.51 × 10−1	not known	3.46 × 10−1
ce (Vs)	M	not known	1.369	1.376	1.374	not known	1.314
J (kgm2)	M	≈ 4 × 10−2	4.68 × 10−2	5.00 × 10−2	4.44 × 10−2	≈ 4 × 10−2	1.02 × 10−2
Tla (Nm)	M	≈ 0.9	7.99 × 10−1	9.34 × 10−1	1.01	≈ 0.9	7.85 × 10−1
Tlb (Nms)	M	≈ 0	7.70 × 10−18	1.10 × 10−17	5.38 × 10−18	≈ 7.5 × 10−2	5.13 × 10−2
Tlc (Nms2)	M	≈ 0	8.00 × 10−19	1.31 × 10−19	6.59 × 10−20	≈ 0	7.37 × 10−20

.

Considering Table 17, the first five parameters should be the same for all input data. The last two parameters T_lb and T_lc are different, because in the cases of MD1, MD2, and MD3, no load was present, and in case of MD4 a load of approximately 14.6 Nm at angular speed 182 s⁻¹ was present (the load was a linear function of speed).

Conclusions considering the calculated values:

• R_a: The maximum deviation of calculated R_a was less than 14%, which is a good result.
• L_a: Calculated L_a was between 244 mH and 451 mH, the difference was due to the deviations in measurements. As can be seen from the result, a small deviation in the measurement causes a significant change of L_a.
• c_m: Calculated c_m was between 1.314 and 1.376, which is a good result, because the difference between calculated values was less than 5%.
• J: The real inertia of the drive was approximately 0.04 kgm² (inertia of some small parts is not known). The calculated inertia was in the range of 20% in the cases of MD1, MD2, and MD3. Only in the case of MD4, was the calculated inertia too small, due to the deviations of the measurements (see a more detailed explanation in Section 4.2.4).
• T_la: Parameter T_la was a friction of the drive, which was approximately 0.9 Nm. The calculated values for all four input data were in the scope of 15%.
• T_lb, T_lc: Parameters T_lb and T_lc were negligibly small in cases MD1, MD2, and MD3. In the case of MD4 load was present, which was linearly dependent on speed. The total load at angular speed 182 s⁻¹ was approximately 14.6 Nm. The calculated load considering parameters T_la, T_lb, and T_lc was 10.1 Nm, which was a deviation of 31%.

The obtained results were very satisfactory, considering the facts that the input data were measured data with measurement deviations and errors, that the motor model did not consider the reactance of the armature, and that real L_a was not completely a constant value. Although the results were very satisfactory, we tried to improve them, which is presented in the continuation of the paper in Section 4.3.

4.2.6. Cross-Validation of the Obtained Results

From Table 17 it can be seen that for MD1, MD2, and MD3, all seven parameters should be the same, because the same drive without load is used. However, different conditions were used for the measurement (different final speed and current limit, as can be seen from Table 7), which offers a test of the correctness of the results using cross-validation.

The parameters obtained using MD1 are input data for the motor simulation, using the values for MD2 and MD3 presented in Table 7 and OF calculation. In the same way, parameters obtained using MD2, are input data for the motor simulation using values for MD1 and MD3, and the parameters obtained using MD3 are input data for the motor simulation using values for MD1 and MD2. The cross-validation procedure and calculated OFs are presented in Figure 12.

In Figure 12, R_MD1, R_MD2, and R_MD3 are calculated sets of the parameters using MD1, MD2, and MD3. Analyzing MD1 test data, based on Figure 12, it can be seen that OF = 2.2531 × 10⁻³ for MD1. In the case of the motor simulation using values for MD1 and R_MD2, OF = 2.4928 × 10⁻³, which is only a slightly higher value. In the case of the motor simulation using values for MD1 and R_MD3, OF = 5.3863 × 10⁻³, which is a less than two times higher value, and it can be assumed that cross-validation confirms the correctness of the results. Based on Figure 12, similar can be concluded analyzing MD2 and MD3 test data.

The simulated motor speed and current using R_MD1, R_MD2, and R_MD3 for MD1 values are presented in Figure 13, for MD2 values are presented in Figure 14 and for MD3 values are presented in Figure 15.

In Figure 13, Figure 14 and Figure 15, it can be seen that the simulated motor speed and current are similar, although different parameters (R_MD1, R_MD2 and R_MD3) are used. For all three Figures, bigger differences appear at current responses and smaller at speed responses.

4.2.7. Testing of Result Correctness at Stationary Points

Testing of result correctness at stationary points was made for results obtained using RK, presented in Table 17. The whole measured data for MD1 are presented in Figure 16.

From Figure 16 it can be seen that after the start of the motor, stationary operation occurs. Although small changes of measured voltage are present, due to the measurement deviations and influence of the controller, operating points after the start of the motor can be treated as stationary. In the case of stationary operation derivatives di_a/dt=0 and dω/dt=0, and, with that Equations (6) and (7) are changed into (25) and (26).

$$u_a=i_a·R_a+L_a·\frac{d{i}_a}{dt}+c_m·\omega \Rightarrow u_a=i_a·R_a+c_m·\omega \Rightarrow u_{a\_\mathrm{calc}}=i_{a\_\mathrm{meas}}·R_a+c_m·{\omega }_{\mathrm{meas}}$$

(25)

$$\begin{array}{rr}{c}_m·i_a-\left(T_{\mathrm{la}}+T_{\mathrm{lb}}·\omega +T_{\mathrm{lc}}·{\omega }^2\right)=J\frac{d\omega }{dt}& \Rightarrow c_m·i_a-\left(T_{\mathrm{la}}+T_{\mathrm{lb}}·\omega +T_{\mathrm{lc}}·{\omega }^2\right)=0\\ & \Rightarrow i_{a\_\mathrm{calc}}=\frac{T_{\mathrm{la}}+T_{\mathrm{lb}}·{\omega }_{\mathrm{meas}}+T_{\mathrm{lc}}·{\omega }_{\mathrm{meas}}^2}{c_m}\end{array}$$

(26)

In the case of the correct calculated parameters presented in Table 17, calculated u_{a_calc} using (25) and calculated i_{a_calc} using (26) should be as similar as possible to measured values. Parameters L_a and J are not included in (25) and (26), and their values are not the subject of correctness tests using stationary points. Calculated current and speed and their deviation at selected stationary points for MD1, MD2, MD3, and MD4 are presented in

Table 18. Calculated current and speed using (25) and (26) and their deviation at selected points for MD1, MD2, MD3, and MD4 for calculated parameters obtained with RK.

Stationary Point—Measured and Calculated Values
Data	t (s)	ia_meas (A)	ωmeas (1/s)	ua_meas (V)	ua_calc (V)	Deviation of ua_calc (%)	ia_calc (A)	Deviation of ia_calc (%)
MD1	2	0.61	125.8	177.0	175.3	1.0	0.58	4.3
	4	0.54	124.8	172.8	173.6	0.5	0.58	8.1
	6	0.60	125.5	171.6	174.8	1.9	0.58	2.7
	8	0.56	124.5	171.4	173.3	1.1	0.58	4.2
	10	0.56	125.0	172.6	173.9	0.8	0.58	4.22
MD2	2	0.78	181.9	252.8	254.1	0.5	0.68	13.0
	4	0.74	181.2	251.0	253.1	0.8	0.68	8.3
	6	0.72	182.9	258.4	255.2	1.2	0.68	5.7
	8	0.78	183.0	258.0	255.6	0.9	0.68	13.0
	10	0.69	182.8	262.0	255.0	2.7	0.68	1.6
MD3	2	0.72	181.0	253.0	252.2	0.3	0.74	2.1
	4	0.76	183.2	255.3	255.4	0.0	0.74	3.3
	6	0.65	183.1	255.9	254.8	0.4	0.74	13.1
	7	0.73	183.1	257.5	255.1	0.9	0.74	0.7
	8	0.67	182.7	257.1	254.3	1.1	0.74	9.7
MD4	4	7.87	182.5	281.5	282.1	0.2	7.72	1.9
	6	7.68	182.3	283.9	280.8	1.1	7.71	0.4
	8	7.79	182.3	280.9	281.4	0.2	7.71	1.0
	10	7.69	182.2	283.6	280.8	1.0	7.71	0.3
	12	7.78	182.4	282.3	281.5	0.3	7.72	0.8

.

In Table 18, times 2, 4, 6, 8, and 10 s were selected for MD1 and MD2. In the case of MD3, times 2, 4, 6, 7, and 8 s were selected, because the signal was measured for a shorter final time. In the case of MD4, times 4, 6, 8, 10, and 12 were selected, because the start of the motor was longer than 2 s.

From Table 18 it can be seen that deviations of calculated voltages are up to 3%. Deviations of calculated speeds are, at 13 stationary points, less than 5%, and in the remaining seven, up to 13%. Based on the presented results it can be concluded that calculated parameters are correct for stationary points.

4.2.8. Calculation Times

Mean values of calculation times for all four measured input data, for EM and RK, and for all five tested evolutionary methods are presented in

Table 19. Mean values of calculation times for 50 independent runs using EM and RK for all measured input data.

Measured Data—Method			Method
		GA	DE/rand/1/exp	DE/best/1/bin	TLBO	ABC
MD1—E	M t(s)	18.8	25.5	25.1	52.0	62.0
MD1—RK	M t(s)	22.2	30.2	28.9	70.3	62.5
MD2—E	M t(s)	16.8	23.7	22.7	56.5	54.9
MD2—RK	M t(s)	20.5	27.7	26.6	66.3	58.0
MD3—E	M t(s)	16.2	21.7	21.2	58.8	46.1
MD3—RK	M t(s)	19.3	25.1	26.0	58.8	49.3
MD4—E	M t(s)	16.4	24.6	19.9	40.9	70.2
MD4—RK	M t(s)	21.4	30.3	28.8	80.0	71.5

.

Based on Table 19, the following conclusions can be written:

• Calculation times were mostly longer for a higher number of measured points (MD1 87, MD2 81, MD3 65, MD4 106 measured points), which was the expected result.
• Calculation times using RK were mostly approximately not more than 25% longer than calculation times using EM, although in the case of RK, ten expressions (14)–(23) were calculated for each time instant, and in the case of EM, only two expressions (11) and (12) were calculated for each time instant.
• Calculation times using different evolutionary methods were different. The fastest was GA. The calculation times of DE were 20 to 40% longer than the times of GA. Calculation times of TLBO and ABC were 100 to 150% longer than the calculation times of DE.

Considering the results presented in Section 3.1 (Simulated Input Data), and Section 4.2.1, Section 4.2.2, Section 4.2.3 and Section 4.2.4 (measured input data), DE/rand/1/exp can be suggested for the presented problem as the best method among those used in this study. TLBO was only slightly worse, but calculation times were much longer than in the case of DE. DE/best/1/bin and ABC were also acceptable, but not as good as DE/rand/1/exp. GA cannot be suggested as a method appropriate for the presented problem.

Improvements of the calculations using only DE/rand/1/exp are presented in the continuation of the paper. Improvements were needed because EM was not suitable, and even in the case of RK, we do not know if we got the most that can be obtained from the presented measurements.

4.2.9. Influence of Weights Used in OF on Results

The OF written in (24) is composed of two parts: Current and speed. Each of them is normalized. By adding weights, we can increase the impact of an individual part, and thus improve the results. Weights w₁ and w₂, to which w₁ + w₂ = 1 applies, are added to (24), and (27) is obtained.

$$\mathrm{OF}=\frac{1}{N}{\displaystyle \sum _{i=1}^N\left(w_1·{\left(\frac{i_{a\_simulated\_i}-i_{a\_input\_i}}{i_{a\_input\_\max }}\right)}^2+w_2·{\left(\frac{{\omega }_{simulated\_i}-{\omega }_{input\_i}}{{\omega }_{input\_\max }}\right)}^2\right)}$$

(27)

Nine combinations of weights w₁, w₂ are used, which are 0.1, 0.9; 0.2, 0.8; … 0.9, 0.1, and tested using MD1, MD2, MD3, and MD4. For each test, 50 independent runs were made using RK. In

Table 20. Calculated parameters with different weights using RK for MD1 (mean values for 50 independent runs).

	Parameters and OF
Weights w1, w2	Ra (Ω)	La (H)	La (H)	J (kgm2)	Tla (Nm)	Tlb (Nm·s)	Tlc (Nm·s2)	OF
0.1, 0.9	4.56	2.32 × 10−1	1.371	4.85 × 10−2	5.00 × 10−16	7.72 × 10−18	5.19 × 10−5	4.9217 × 10−4
0.2, 0.8	4.70	2.32 × 10−1	1.371	4.87 × 10−2	5.95 × 10−16	1.76 × 10−17	5.12 × 10−5	6.6101 × 10−4
0.3, 0.7	4.80	2.34 × 10−1	1.370	4.89 × 10−2	1.08 × 10−15	6.01 × 10−18	5.09 × 10−5	8.2256 × 10−4
0.4, 0.6	4.88	2.37 × 10−1	1.370	4.90 × 10−2	1.60 × 10−13	7.47 × 10−18	5.06 × 10−5	9.8028 × 10−4
0.5, 0.5	5.06	2.44 × 10−1	1.369	4.68 × 10−2	7.99 × 10−1	7.70 × 10−18	8.00 × 10−19	1.1266 × 10−3
0.6, 0.4	5.14	2.46 × 10−1	1.368	4.68 × 10−2	8.00 × 10−1	4.91 × 10−18	1.25 × 10−19	1.2647 × 10−3
0.7, 0.3	5.21	2.48 × 10−1	1.368	4.69 × 10−2	8.00 × 10−1	4.38 × 10−18	5.54 × 10−20	1.4006 × 10−3
0.8, 0.2	5.27	2.51 × 10−1	1.367	4.69 × 10−2	8.00 × 10−1	2.34 × 10−18	3.52 × 10−20	1.5345 × 10−3
0.9, 0.1	5.34	2.53 × 10−1	1.367	4.70 × 10−2	8.00 × 10−1	1.47 × 10−18	3.37 × 10−20	1.6667 × 10−3
Known value	5.66	not known	not known	≈ 4 × 10−2	≈ 0.9	≈ 0	≈ 0	/

,

Table 21. Calculated parameters with different weights using RK for MD2 (mean values for 50 independent runs).

	Parameters and OF
Weights w1, w2	Ra (Ω)	La (H)	La (H)	J (kgm2)	Tla (Nm)	Tlb (Nm·s)	Tlc (Nm·s2)	OF
0.1, 0.9	4.37	2.27 × 10−1	1.379	5.12 × 10−2	6.50 × 10−16	1.03 × 10−17	2.88 × 10−5	4.4913 × 10−4
0.2, 0.8	4.57	2.31 × 10−1	1.378	5.16 × 10−2	1.02 × 10−15	7.00 × 10−18	2.84 × 10−5	6.8282 × 10−4
0.3, 0.7	4.69	2.35 × 10−1	1.377	5.19 × 10−2	7.86 × 10−15	2.57 × 10−17	2.82 × 10−5	9.0588 × 10−4
0.4, 0.6	4.87	2.42 × 10−1	1.376	4.99 × 10−2	8.99 × 10−1	2.96 × 10−17	1.07 × 10−6	1.1189 × 10−3
0.5, 0.5	4.97	2.46 × 10−1	1.376	5.00 × 10−2	9.34 × 10−1	1.10 × 10−17	1.31 × 10−19	1.3192 × 10−3
0.6, 0.4	5.05	2.48 × 10−1	1.375	5.01 × 10−2	9.32 × 10−1	7.22 × 10−18	9.37 × 10−20	1.5156 × 10−3
0.7, 0.3	5.12	2.51 × 10−1	1.375	5.02 × 10−2	9.32 × 10−1	5.88 × 10−18	3.42 × 10−20	1.7089 × 10−3
0.8, 0.2	5.18	2.53 × 10−1	1.374	5.02 × 10−2	9.30 × 10−1	4.85 × 10−18	4.11 × 10−20	1.8999 × 10−3
0.9, 0.1	5.24	2.55 × 10−1	1.374	5.03 × 10−2	9.11 × 10−1	3.43 × 10−18	5.51 × 10−7	2.0902 × 10−3
Known value	5.66	not known	not known	≈ 4 × 10−2	≈ 0.9	≈ 0	≈ 0	/

,

Table 22. Calculated parameters with different weights using RK for MD3 (mean values for 50 independent runs).

	Parameters and OF
Weights w1, w2	Ra (Ω)	La (H)	La (H)	J (kgm2)	Tla (Nm)	Tlb (Nm·s)	Tlc (Nm·s2)	OF
0.1, 0.9	4.63	3.70 × 10−1	1.375	4.71 × 10−2	5.82 × 10−2	9.40 × 10−18	2.88 × 10−5	5.4941 × 10−4
0.2, 0.8	4.69	4.03 × 10−1	1.375	4.63 × 10−2	3.19 × 10−1	2.93 × 10−17	2.08 × 10−5	8.9597 × 10−4
0.3, 0.7	4.78	4.25 × 10−1	1.375	4.52 × 10−2	7.21 × 10−1	1.73 × 10−17	8.69 × 10−6	1.2326 × 10−1
0.4, 0.6	4.86	4.42 × 10−1	1.374	4.45 × 10−2	1.01	6.96 × 10−18	4.48 × 10−19	1.5621 × 10−3
0.5, 0.5	4.87	4.51 × 10−1	1.374	4.44 × 10−2	1.01	5.38 × 10−18	6.59 × 10−20	1.8883 × 10−3
0.6, 0.4	4.89	4.57 × 10−1	1.374	4.44 × 10−2	1.01	7.68 × 10−18	5.46 × 10−20	2.2128 × 10−3
0.7, 0.3	4.90	4.63 × 10−1	1.374	4.43 × 10−2	1.01	4.67 × 10−18	4.66 × 10−20	2.5363 × 10−3
0.8, 0.2	4.91	4.67 × 10−1	1.374	4.43 × 10−2	1.02	5.85 × 10−18	2.51 × 10−20	2.8590 × 10−3
0.9, 0.1	4.91	4.70 × 10−1	1.374	4.43 × 10−2	1.02	3.67 × 10−18	3.32 × 10−20	3.1812 × 10−3
Known value	5.66	not known	not known	≈ 4 × 10−2	≈ 0.9	≈ 0	≈ 0	/

and

Table 23. Calculated parameters with different weights using RK for MD4 (mean values for 50 independent runs).

	Parameters and OF
Weights w1, w2	Ra (Ω)	La (H)	La (H)	J (kgm2)	Tla (Nm)	Tlb (Nm·s)	Tlc (Nm·s2)	OF
0.1, 0.9	4.44	2.64 × 10−1	1.357	1.07 × 10−1	1.21 × 10−15	5.77 × 10−2	1.74 × 10−20	4.5890 × 10−4
0.2, 0.8	4.64	2.82 × 10−1	1.348	1.08 × 10−1	2.27 × 10−15	5.72 × 10−2	3.12 × 10−20	7.7123 × 10−4
0.3, 0.7	4.82	2.95 × 10−1	1.340	1.08 × 10−1	4.05 × 10−14	5.68 × 10−2	6.08 × 10−20	1.0743 × 10−3
0.4, 0.6	5.08	3.17 × 10−1	1.328	1.06 × 10−1	3.53 × 10−1	5.43 × 10−2	6.69 × 10−20	1.3678 × 10−3
0.5, 0.5	5.38	3.46 × 10−1	1.314	1.02 × 10−2	7.85 × 10−1	5.13 × 10−2	7.37 × 10−20	1.6474 × 10−3
0.6, 0.4	5.76	3.89 × 10−1	1.297	9.76 × 10−2	1.30	4.77 × 10−2	6.33 × 10−20	1.9089 × 10−3
0.7, 0.3	5.38	3.15	1.319	8.55 × 10−2	2.91	2.53 × 10−2	8.15 × 10−5	2.0945 × 10−3
0.8, 0.2	5.34	3.99	1.323	7.12 × 10−2	3.92	3.50 × 10−2	4.79 × 10−19	2.0480 × 10−3
0.9, 0.1	6.28	4.60	1.282	6.06 × 10−2	4.94	2.77 × 10−2	1.53 × 10−19	1.8430 × 10−3
Known value	5.66	not known	not known	≈ 4 × 10−2	≈ 0.9	≈ 7.5 × 10−2	≈ 0	/

the mean values of calculated parameters are presented for MD1, MD2, MD3, and MD4 separately.

In Table 20, Table 21, Table 22 and Table 23 OFs are also presented, although a lower value does not mean a better result and the values are not comparable.

Based on Table 20, Table 21, Table 22 and Table 23 the following conclusions can be written:

• In the case of MD1, correct results (due to the known values written at the bottom of Table 20) are obtained for weights w₁ ≥ 0.5. For w₁ ≤ 0.4, parameter T_la has an incorrect value. For higher w₁, slightly better results were obtained, for example R_a at w₁ = 0.5 is 5.06 and at w₁ = 0.9 it is 5.34.
• In the case of MD2, correct results were obtained for weights w₁ ≥ 0.4. For w₁ ≤ 0.3 parameter, T_la has an incorrect value. For higher w₁, slightly better results were obtained, for example, R_a at w₁ = 0.4 is 4.87 and at w₁ = 0.9 it is 5.24.
• In the case of MD3 correct results were obtained for weights w₁ ≥ 0.4. For w₁ ≤ 0.3 parameter, T_lc has an incorrect value. For higher w₁, almost no difference in results was obtained, for example, R_a at w₁ = 0.4 is 4.86 and at w₁ = 0.9 it is 4.91.
• In the case of MD4, only results for w₁ = 0.4 and 0.5 are acceptable. For w₁ ≤ 0.3, T_la is too small, and for w₁ ≥ 0.6, T_1a is too big.

Based on the written findings, only weight w₁ = 0.5 was appropriate for all the presented test cases. However, general conclusions cannot be made only based on this one example. The basic OF written in (24) is used in the continuation of the paper.

4.3. Motor Model Simulations with Time Interval Division

The measurement interval, the time step, could not be reduced, because it depends on the hardware. However, it is well known that the solutions of differential equations can be better if a smaller time step is used [49]. To improve results, we reduced the time step in the case of the motor model simulations (OF calculation), and only calculated current and speed at time instants in which measured values were known were used for OF calculation. For the tests, the number of time step divisions (ND) was set to 5 and 10. The used methods are marked with EM5 for Euler’s method using ND = 5, EM10 for Euler’s method using ND = 10, RK5 for the Runge–Kutta method using ND = 5 and RK10 for Runge–Kutta using ND = 10. The times between two sequential measurement time instants, marked with t_i, were calculated using (28).

$$t_i=t_{k-1}+\frac{t_k-t_{k-1}}{ND}·i;i=1,\dots ,ND-1$$

(28)

t_k−1 and t_k are sequential measurement time instants. Also, voltage is an input value, which is measured, and it changes with the time. The voltages between two sequential measured voltages, marked with u_{a_i}, were determined using linear interpolation, defined with (29).

$$u_{a\_i}=u_{a\_k-1}+\frac{u_{a\_k}-u_{a\_k-1}}{ND}·i;i=1,\dots ,ND-1$$

(29)

u_{a_k−1} and u_{a_k} are sequential measured voltages.

The algorithm showing the calculation procedure is presented in Figure 17. Parts which were added to the algorithms of EM and RK without a time step division (presented in Figure 4) are marked with red.

OF and mean values of the calculated parameters using EM, EM5, EM10, RK, RK5, and RK10 for MD1 are presented in

Table 24. OF and Mean value of the calculated parameters for 50 independent runs using EM, EM5, EM10, RK, RK5, and RK10 for MD1.

OF and				Method
Parameters		EM	EM5	EM10	RK	RK5	RK10
	B	2.8906 × 10−3	2.1738 × 10−3	2.1541 × 10−3	2.2531 × 10−3	2.1453 × 10−1	2.1464 × 10−1
OF	W	2.8906 × 10−3	2.1738 × 10−3	2.1541 × 10−3	2.2531 × 10−3	2.1453 × 10−1	2.1464 × 10−1
	M	2.8906 × 10−3	2.1738 × 10−3	2.1541 × 10−3	2.2531 × 10−3	2.1453 × 10−1	2.1464 × 10−1
	SD	0.0	1.7347 × 10−18	1.3010 × 10−18	0.0	8.6736 × 10−19	2.1684 × 10−18
Ra (Ω)	M	11.56	5.63	4.94	5.06	4.25	4.25
La (H)	M	9.72 × 10−1	2.64 × 10−1	2.07 × 10−1	2.44 × 10−1	1.57 × 10−1	1.57 × 10−1
ce (Vs)	M	1.337	1.366	1.370	1.369	1.373	1.373
J (kgm2)	M	4.64 × 10−2	4.68 × 10−1	4.71 × 10−2	4.68 × 10−2	4.77 × 10−2	4.76 × 10−2
Tla (Nm)	M	3.13 × 10−16	7.96 × 10−1	8.00 × 10−1	7.99 × 10−1	8.03 × 10−1	8.03 × 10−1
Tlb (Nms)	M	8.95 × 10−17	7.84 × 10−18	6.53 × 10−18	7.70 × 10−18	6.27 × 10−17	3.22 × 10−18
Tlc (Nms2)	M	4.99 × 10−5	1.22 × 10−19	1.35 × 10−19	8.00 × 10−19	6.04 × 10−20	5.09 × 10−20

.

Base on Table 24 it can be seen that results using EM were improved. In the cases of EM5 and EM10 calculated R_a and T_la were correct, although they were not correct using EM. Considering OFs, EM5 was an improvement on EM, and EM10 was only a small improvement on EM5. RK5 was a small improvement on RK, and RK10 was not an improvement on RK5. It must be noted that lower OF values does not automatically mean more precisely calculated parameters, because in the measured signal, measurement deviations and errors are present.

Measured and calculated responses using EM5 and EM10 are shown in Figure 18, and using RK5 and RK10 in Figure 19. Only the first second of the responses is displayed to make Figures more transparent.

Comparing Figure 8 and Figure 18, the improvement achieved with EM5 and EM10 can be seen. The simulated responses were closer to the measured responses.

OF and mean values of the calculated parameters using EM, EM5, EM10, RK, RK5, and RK10 for MD2 are shown in

Table 25. OF and Mean value of the calculated parameters for 50 independent runs using EM, EM5, EM10, RK, RK5, and RK10 for MD2.

OF and				Method
Parameters		EM	EM5	EM10	RK	RK5	RK10
	B	3.4006 × 10−3	2.4870 × 10−1	2.4192 × 10−1	2.6384 × 10−3	2.3706 × 10−1	2.3676 × 10−1
OF	W	3.4006 × 10−3	2.4870 × 10−1	2.4664 × 10−1	2.6384 × 10−3	2.3706 × 10−1	2.4366 × 10−1
	M	3.4006 × 10−3	2.4870 × 10−1	2.4306 × 10−1	2.6384 × 10−3	2.3706 × 10−1	2.3690 × 10−1
	SD	1.7347 × 10−18	2.6021 × 10−18	2.0131 × 10−5	1.3010 × 10−18	2.3706 × 10−1	9.6658 × 10−6
Ra (Ω)	M	11.45	5.72	5.00	4.97	4.38	4.37
La (H)	M	9.76 × 10−1	2.69 × 10−1	2.11 × 10−1	2.46 × 10−1	1.61 × 10−1	1.61 × 10−1
ce (Vs)	M	1.350	1.373	1.376	1.376	1.378	1.378
J (kgm2)	M	4.89 × 10−2	4.86 × 10−2	4.95 × 10−2	5.00 × 10−2	4.98 × 10−2	4.96 × 10−2
Tla (Nm)	M	3.58 × 10−16	9.19 × 10−1	7.01 × 10−1	9.34 × 10−1	9.23 × 10−1	9.05 × 10−1
Tlb (Nms)	M	9.61 × 10−17	9.43 × 10−18	6.73 × 10−20	1.10 × 10−17	7.19 × 10−18	8.90 × 10−18
Tlc (Nms2)	M	2.77 × 10−5	7.32 × 10−20	6.59 × 10−6	1.31 × 10−19	4.75 × 10−20	5.48 × 10−7

.

Base on Table 25, it can be seen that results using EM were improved. In the cases of EM5 and EM10 calculated R_a and T_la were correct, although they were not correct using EM. Considering OFs, EM5 was an improvement on EM, and EM10 was only a small improvement on EM5. RK5 was a small improvement on RK, and RK10 was not an improvement on RK5.

The first second of the measured and calculated responses using EM10 and RK10 is shown in Figure 20.

Comparing Figure 9 and Figure 20, especially the improvement achieved with EM10 can be seen.

OF and mean values of the parameters calculated using EM, EM5, EM10, RK, RK5, and RK10 for MD3 are presented in

Table 26. OF and Mean values of the calculated parameters for 50 independent runs using EM, EM5, EM10, RK, RK5, and RK10 for MD3.

OF and				Method
Parameters		EM	EM5	EM10	RK	RK5	RK10
	B	5.7731 × 10−3	4.0635 × 10−3	3.9179 × 10−1	3.7766 × 10−3	3.7887 × 10−1	3.7887 × 10−1
OF	W	5.7731 × 10−3	4.0635 × 10−3	3.9179 × 10−1	3.7766 × 10−3	3.7887 × 10−1	3.7887 × 10−1
	M	5.7731 × 10−3	4.0635 × 10−3	3.9179 × 10−1	3.7766 × 10−3	3.7887 × 10−1	3.7887 × 10−1
	SD	8.6736 × 10−19	2.6021 × 10−18	0.0	0.0	3.7887 × 10−1	3.0358 × 10−18
Ra (Ω)	M	11.92	6.05	5.42	4.87	4.74	4.74
La (H)	M	1.28	5.96 × 10−1	5.38 × 10−1	4.51 × 10−1	4.79 × 10−1	4.79 × 10−1
ce (Vs)	M	1.342	1.369	1.372	1.374	1.375	1.375
J (kgm2)	M	4.42 × 10−2	4.45 × 10−1	4.41 × 10−1	4.44 × 10−2	4.45 × 10−2	4.45 × 10−1
Tla (Nm)	M	6.45 × 10−16	7.33 × 10−1	1.01	1.01	1.01	1.01
Tlb (Nms)	M	9.91 × 10−18	1.08 × 10−17	1.62 × 10−17	5.38 × 10−18	9.13 × 10−18	1.81 × 10−17
Tlc (Nms2)	M	3.11 × 10−5	8.71 × 10−6	1.09 × 10−7	6.59 × 10−20	2.46 × 10-19	2.51 × 10−19

.

Based on Table 26, it can be seen that results were improved using EM. In the cases of EM5 and EM10 calculated R_a and T_la were correct, although they were not correct using EM. Considering OFs, EM5 was an improvement on EM, and EM10 was only a small improvement on EM5. In the cases of RK5 and RK10 no improvements were made.

The first second of the measured and calculated responses using EM10 and RK10 is shown in Figure 21.

Comparing Figure 10 and Figure 21, especially the improvement achieved with EM10 can be seen.

OF and mean values of the calculated parameters using EM, EM5, EM10, RK, RK5, and RK10 for MD4 are shown in

Table 27. OF and Mean value of the calculated parameters for 50 independent runs using EM, EM5, EM10, RK, RK5, and RK10 for MD4.

OF and				Method
Parameters		EM	EM5	EM10	RK	RK5	RK10
	B	2.9627 × 10−3	3.3114 × 10−3	3.2788 × 10−3	3.2949 × 10−3	3.2719 × 10−1	3.2719 × 10−3
OF	W	2.9627 × 10−3	3.3114 × 10−3	3.2911 × 10−3	3.2949 × 10−3	3.2719 × 10−1	3.2719 × 10−3
	M	2.9627 × 10−3	3.3114 × 10−3	3.2906 × 10−3	3.2949 × 10−3	3.2719 × 10−1	3.2719 × 10−3
	SD	4.3368 × 10−19	1.7347 × 10−18	2.4102 × 10−6	8.6736 × 10−19	2.1684 × 10−18	8.6736 × 10−19
Ra (Ω)	M	6.97	5.84	5.60	5.38	5.36	5.36
La (H)	M	2.82	3.48 × 10−1	2.63 × 10−1	3.46 × 10−1	1.99 × 10−1	1.98 × 10−1
ce (Vs)	M	1.251	1.295	1.305	1.314	1.315	1.315
J (kgm2)	M	1.19 × 10−1	1.02 × 10−1	1.02 × 10−1	1.02 × 10−1	1.02 × 10−1	1.01 × 10−1
Tla (Nm)	M	1.27 × 10−15	6.35 × 10−1	7.34 × 10−1	7.85 × 10−1	8.29 × 10−1	8.29 × 10−1
Tlb (Nms)	M	2.66 × 10−17	5.13 × 10−2	5.12 × 10−2	5.13 × 10−2	5.11 × 10−2	5.11 × 10−2
Tlc (Nms2)	M	2.91 × 10−4	1.44 × 10−19	2.63 × 10−19	7.37 × 10−20	8.78 × 10−20	1.49 × 10−19

.

Based on Table 27 it can be seen that results using EM were improved, although the OF value was higher. In this case, the measurement deviations and errors were bigger than the EM5 and EM10 methods, which simulated a real motor better, could not get as close to the measured signal as the EM method, which simulated a real motor worse. In the cases of EM5 and EM10, calculated R_a and T_la were correct, although they were not correct using EM. Considering OFs, EM10 was only a small improvement on EM5. In the cases of RK5 and RK10 only very small improvements were made.

The first four seconds of the measured and calculated responses using EM10 and RK10 are shown in Figure 22.

Division of the time interval in the case of the motor simulations results in a higher number of calculations, and, with that, with a longer calculation time. Calculation times using E, EM5, EM10, RK, RK5, and RK10 for all measured input data are presented in

Table 28. Mean values of the calculation times for 50 independent runs using EM, EM5, EM10, RK, RK5, and RK10 for all measured input data.

Measured Data—Method				Method
		EM	EM5	EM10	RK	RK5	RK10
MD1	M t(s)	25.5	39.31	48.9	30.2	54.9	79.9
MD2	M t(s)	23.7	38.4	45.7	27.7	52.3	74.8
MD3	M t(s)	21.7	33.4	39.7	25.1	45.5	63.3
MD4	M t(s)	24.6	42.4	50..2	30.3	62.7	91.0

.

Based on Table 28, it can be seen that calculation times using EM5 were, on average, 1.5 times longer than the times of EM, and using EM10 were, on average, 2 times longer than the times of EM. Calculation times using RK5 were, on average, 2 times longer than the times of RK, and using RK10 were, on average, 3 times longer than the times of RK.

4.4. Analysis of Memory Assistance

Memory Assistance offers the possibility to shorten calculation times [50], which is always welcome. The presented problem is interesting for Memory Assistance analysis, because using different approaches different calculation times are achieved (without division of the time interval and with 5 or 10 divisions of the time interval).

We used three different Memory Assistances, which are the following:

• STMA: In each iteration, population members are compared with members from only one previous iteration. After each iteration, the whole population is saved into memory. Only additional memory for one population set is used using the presented approach. In the presented case, this is ((7 parameters + OF) × 70 population members) memory locations.
• Long-Term Memory Assistance—Strategy 1 (LTMA-S1): Each population member is compared with all members written into memory, obtained from all previous iterations. Each member which is not found in the memory is added to the memory. Strategy 1 means that the search in the memory starts from first added to the last added (when the same population member is found, the search is finished). Theoretically ((7 parameters + OF) × 70 population members × 2000 iterations) locations of additional memory can be used.
• Long-Term Memory Assistance—Strategy 2 (LTMA-S2): It is the same as LTMA-S1, only the strategy of the search is changed. It starts from the last added to the first added.

LTMA-S2 is better than LTMA-S1, probably because the latter solutions are less different and duplicates occur sooner. Using both strategies it will be seen if memory processing time is negligible or large compared to the calculation time of the Evolutionary Method.

The STMA algorithm is presented in Figure 23, and the LTMA algorithm for both strategies is presented in Figure 24.

The expression “record” used in Figure 23 and Figure 24 means one set of parameters + OF value (R_a, L_a, c_m, J, T_la, T_lb, T_lc, OF). It is important to set appropriate precision, differences of the parameters for which duplication may be identified. These values determine the decimal accuracy of the calculated parameters, and they depend on the physical background of the problem. Precisions used for duplication identification are presented in

Table 29. Precisions which determine duplication identification.

Parameter	Maximum Difference
Ra (Ω)	<1 mΩ
La (H)	<1 mH
cm (Vs)	<1 mVs
J (kgm2)	<10−4 kgm2
Tla (Nm)	<1 mNm
Tlb (Nm·s)	<5 × 10−6 Nm·s	(torque of 1 mNm at speed 209 s−1—motor rated speed is 182 s−1)
Tlc (Nm·s2)	<2 × 10−8 Nm·s2	(torque of 1 mNm at speed 209 s−1—motor rated speed is 182 s−1)

.

Results presenting calculation times and repetitions obtained with EM, EM5, EM10, RK, RK5, and RK10 using STMA, LTMA-S1, and LTMA-S2 are presented in

Table 30. Mean value of the calculation Times (Time) and Duplications (Dupl.) in respect to no additional memory use for 50 independent runs using EM, EM5, EM10, RK, RK5, and RK10 for MD1.

MD1	Method
	EM		EM5		EM10		RK		RK5		RK10
Memory Use	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)
STMA	104.8	51.5	84.0	54.2	67.0	54.2	96.1	50.1	73.6	54.9	62.7	56.5
LTMA-S1	334.3	81.5	210.8	84.9	159.5	84.4	318.1	80.9	148.6	85.4	104.6	86.1
LTMA-S2	188.2	79.3	104.2	84.0	86.2	84.2	189.6	75.8	83.3	84.1	56.8	85.4

for MD1,

Table 31. Mean value of the calculation Times (Time) and Duplications (Dupl.) in respect to no additional memory use for 50 independent runs using EM, EM5, EM10, RK, RK5, and RK10 for MD2.

MD2	Method
	EM		EM5		EM10		RK		RK5		RK10
Memory Use	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)
STMA	102.9	51.4	78.1	55.1	69.7	57.5	94.3	53.1	70.1	56.2	61.1	57.6
LTMA-S1	354.7	84.0	220.6	84.0	173.5	85.8	297.1	84.0	149.1	86.1	107.7	86.4
LTMA-S2	176.0	82.1	106.9	82.7	84.7	84.3	145.3	82.6	78.2	85.1	58.4	85.3

for MD2,

Table 32. Mean value of the calculation Times (Time) and Duplications (Dupl.) in respect to no additional memory use for 50 independent runs using EM, EM5, EM10, RK, RK5, and RK10 for MD3.

MD3	Method
	EM		EM5		EM10		RK		RK5		RK10
Memory Use	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)
STMA	109.0	55.2	96.8	36.3	88.4	38.9	100.8	56.3	74.3	56.5	64.7	57.7
LTMA-S1	361.9	85.4	290.4	81.3	231.1	82.6	307.0	85.8	180.0	85.7	133.4	86.0
LTMA-S2	162.8	83.7	159.2	77.9	134.4	78.6	138.1	85.2	85.0	84.9	64.3	84.9

for MD3, and

Table 33. Mean value of the calculation Times (Time) and Duplications (Dupl.) in respect to no additional memory use for 50 independent runs using EM, EM5, EM10, RK, RK5, and RK10 for MD4.

MD4	Method
	EM		EM5		EM10		RK		RK5		RK10
Memory Use	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)	Time (%)	Dupl. (%)
STMA	102.5	40.2	88.7	33.8	51.1	34.6	99.9	32.6	85.3	33.1	68.5	33.9
LTMA-S1	431.4	77.8	210.6	72.5	175.2	74.0	434.3	71.8	203.2	72.8	104.6	73.6
LTMA-S2	215.8	74.3	164.7	68.8	132.4	70.6	240.1	67.1	127.0	68.7	91.3	69.3

for MD4. Time (%) is the relative value of used time according to the calculation time without Memory Assistance. Duplications (%) is the relative number of duplications according to all possible duplications, which are, in our problem, (70 population members × 2000 iterations = 140,000) duplications.

Based on Table 30, Table 31, Table 32 and Table 33, the following remarks can be written:

• Duplications for STMA were between 50–60% for MD1 and MD2, between 35–60% for MD3 and between 30–40% for MD4.
• Duplications for LTMA (S1 and S2) were between 75–86% for MD1, MD2, and MD3 and between 67–78% for MD4.
• Based on previous remarks, it can be seen that in the case of LTMA, the number of duplications was much higher than in the case of STMA, which was expected. It is interesting that duplications depend not only on precision, stopping criteria, type of problem, etc., but also on input data.
• Calculation times using STMA were, in the cases of EM and RK, in the range from 94.3% up to 109% of times without Memory Assistance, so the time was only a little shorter or even longer. In the cases of EM5, EM10, RK5, and RK10, times were shortened up to 51.1% of times without Memory Assistance.
• Calculation times using LTMA-S1 were longer than times without Memory Assistance. They were in the range from 104.6% in the case of RK10 (MD1 and MD4), up to 431.4% in the case of EM (MD4). The presented problem is not time-consuming enough to be suitable for LTMA-S1.
• Calculation times using LTMA-S2 were longer than times without Memory Assistance in the cases of EM, EM5 and RK, in the range from 104.2% in the case of EM5 (MD1) up to 240.1% in the case of RK (MD4). In the cases of EM10 and RK5, times were in some cases shorter and in some longer. They were in the range from 78.2 up to 132.4%. The shortest times were obtained in the case of RK10, which was the most time-consuming method. They were in the range from 56.8 up to 64.3% for MD1, MD2, and MD3, which were even shorter than the times in the case of STMA. Only in the case of MD4, where the number of duplications is smaller, was the calculation time 91.3%, which was a longer time than in the case of STMA.

The use of STMA can be suggested for the presented problem. To get a feeling about the relation between OF calculation time and memory handling time, Figure 25 shows the mean OF calculation times and the mean memory handling times (entry and search) for one fitness evaluation.

5. Conclusions

The extension of the step responses approach into the approach also using measured voltage and considering limitation of the current, which was one of the main contributions of this work, was successful. Such an approach offers the possibility to determine the parameters of big controlled drives without modeling of the controller.

In the case of simulated input data, the computed parameters (in the case of DE/rand/1/exp) were exactly the values used for preparation of the input data, and, thus, confirmed the correctness of the approach. In the case of measured input data, the calculated parameters (for known values) were mainly in the range of 15%, in some cases up to 25%, and only in the case of calculated torque (MD4) was the deviation 30.1%. Due to the fast time responses and measurement deviations and errors, it can be assumed that deviations were within normal limits.

The tested evolutionary methods (GA, DE/rand/1/exp, DE/best/1/bin, TLBO, and ABC) showed the same properties in the case of simulated input data (ideal input without armature reaction and measuring errors) as in the case of measured input data (real input). The DE/rand/1/exp was the most robust. For each run in the scope of individual test data, results were the same. TLBO was only slightly worse, but calculation times were 2 times longer. Results obtained using DE/best/1/bin and ABC were also acceptable. Only results obtained using GA were too poor. Based on this, DE/rand/1/exp can be suggested as the best method for the presented problem.

Comparison between EM and RK showed that EM was not good enough, while RK was appropriate. The suitability of the method depends on the measurement time interval, electrical time constant of the motor and mechanical time constant of the drive, and it is strongly problem dependent. It was expected that the fourth-order RK method was more appropriate than first-order EM, but, on the other hand, calculation times using RK were longer and programming of RK is more demanding.

Results were improved by division of the time interval in the case of the motor model simulations. Also, results using EM5 and EM10 were correct, although results using EM were incorrect. Based on OF analysis, it was seen that in some cases, no improvement was achieved, improvement was achieved using ND = 5 and no more using ND = 10, or improvement was achieved using ND = 5 and only a small, irrelevant additional improvement was achieved using ND = 10. At the start, with a new problem with unknown time constants, it is difficult to determine the appropriate number of divisions. It is suggested to make calculations by increasing the number of divisions until no correction is present.

Analyses of Memory Assistance showed that STMA was appropriate for the presented problem. A good feature of STMA is that only a small amount of additional memory was used (one set of population members). The use of LTMA-S2 was successful only in some cases of EM10 and RK5, and in all cases of RK10. The use of LTMA-S1 was not successful. In the case of a similar problem (similar precision, number of measured points, stopping criteria, used method, ND) the use of STMA can be suggested.

It must be mentioned that in the case of the presented or similar problems, three types of errors are present: (I) The error of the mathematical model, (II) The error of the calculation method and (III) The computer round-off errors [51,52,53]. We tried to reduce error (I) by properly prepared measured data, as presented in Figure 6. In the presented work, the main emphasis was on the reduction of error (II) by decreasing the time step of numerical integration to increase the degree of matching between the parameters of the physical and mathematical models.

We did not consider torsional dynamic properties, masses with torsional stiffness and damping effect, and this can be a topic for further research. However, we considered constant, linear, and square parts of the load, but this does not cover all the physical phenomena during the start of the motor.

In the future, findings about the used evolutionary methods and the use of Memory Assistance, presented in the paper, will be used to solve optimization problems in Electromagnetics and to shorten calculation times. Also, the possibilities to adopt the presented method to parameter determination of other types of motor will be investigated.