Traction Load Modeling and Parameter Identification Based on Improved Sparrow Search Algorithm

: In this paper, a traction load model parameter identification method based on the improved sparrow search algorithm (ISSA) is proposed. According to the load characteristics of the AC traction power supply system under transient disturbance, the model structure of the traction load is equated to the composite load model structure of the static load shunt induction motor’s dynamic load. The traditional sparrow search algorithm is improved to enhance its accuracy and convergence. The generalization ability of the model was tested, and the accuracy of the proposed model was verified. Using the ISSA to determine the load model from the measured data, the results can verify the effectiveness of the ISSA for comprehensive load model parameter identification. Comparing the ISSA with the traditional SSA and PSO algorithms, it shows that the ISSA has better accuracy and convergence.


Introduction
Power systems consist of three main components: the generating unit, the transmission network, and the electrical load.With the continuous updating of computer technology, digital simulation of power systems [1] has become an important tool for power system planning, design, and computational analysis at present.The accuracy of the power system model required for simulation has attracted increasing attention from scholars.So far, the modeling studies of generators and transmission networks are quite mature with the continuous efforts of scholars.Compared with other components, the modeling of electric loads is still very difficult because of the variety of loads.The complex composition, complexity of dynamic and static characteristics, and other factors all lead to randomness, dispersion, and time-varying nature of loads.A large number of experiments and actual operating conditions of power systems already show that the load model has an important influence on the tide calculation, voltage stability analysis, and transient stability analysis of power systems.The establishment of a suitable and accurate load model is of great significance for the digital simulation of power systems.
Traction loads are characterized by randomness, asymmetry, and low power factor.The problems of power quality, such as harmonics and negative sequence currents of traction loads, are not the main focus of transient stability studies; the study of their transient characteristics mainly depends on the performance of the system after its dynamic characteristics are changed by transient disturbances [2,3].Thus, an accurate traction load equivalent composite load model should be established, and its transient stability could be studied by means of professional power system transient simulation analysis.
Currently, load modeling methods include statistical synthesis [4], overall measurement and identification [5], and fault simulation [6].The basic idea of the overall measurement and identification method is to treat the load as a whole system, first determining the load model structure and then using the data collected in the field and the effective identification algorithm to identify the relevant parameters of the load model structure, check the accuracy of the model, and then output the load model after the error meets the accuracy requirements [7].The overall measurement and identification method is an effective way to solve complex integrated load modeling without the need to understand the complex composition within the load.Parameter identification is the core of this method, and the evolving system identification theory provides a strong theoretical basis for this method.Therefore, choosing a suitable parameter identification method is the key to the accuracy of load modeling.
Load model parameter identification can be divided into linear and nonlinear methods.Linear methods mainly include the least squares method, Kalman filter method [8], etc.The nonlinear parameter identification method is mainly based on the optimization algorithm, and the main idea is to find a set of optimal parameter solutions to minimize the preset objective function.At present, the main nonlinear methods proposed are the gradient method, random search method, and simulated evolution method.In [9], a genetic algorithm is first introduced into the parameter identification of the induction motor load model.In [10], a genetic algorithm (GA) and the Levenberg-Marquardt (L-M) algorithm are combined to identify parameters of the comprehensive load model, which improves the accuracy and efficiency of identification.In [11], a combination of chaos and quantum particle swarm algorithms was used to carry out the parameter identification of the integrated load model, and it indicated that the CQDPSO algorithm has better performance in terms of convergence speed and accuracy.
In recent years, many scholars have proposed many new intelligent optimization algorithms [12][13][14][15] by analyzing different biological populations and physical phenomena.The sparrow search algorithm is a new group intelligence algorithm proposed by scholar Xue [16] in 2020, and it has been widely used by scholars because of its excellent convergence and accuracy.Scholar Li [17] tested the typical swarm intelligence algorithms proposed in recent years and compared the experimental performance of these algorithms in terms of convergence speed, accuracy, and stability through 22 standard CEC test functions; the results show that the sparrow search algorithm has better convergence, accuracy, and variance.However, it is still easy to converge prematurely and fall into a local optimum when solving complex nonlinear problems.Therefore, it is necessary to optimize its existing algorithm process to improve the algorithm accuracy and convergence performance.In [18], the authors introduced the random walk strategy into the SSA, proposed an adaptive sparrow search algorithm, and applied it to the model parameter optimization identification of proton exchange membrane fuel cell stack.The results show that the proposed algorithm has the best efficiency compared with the comparison algorithm.In [19], the improved sparrow search algorithm is applied to the IEEE33 nodes of the distributed-generation optimization configuration model, and the DG configuration scheme can reduce the active power loss and voltage deviation of the distribution network to a greater extent.In [20], the population is initialized by using the barycenter reversal mechanism.Aiming at the shortage of the global search ability of the discoverer, the learning coefficient is introduced into the position update part, and the improved algorithm is applied to the photovoltaic microgrid system.The results show that compared with the comparison algorithm, the algorithm can track the maximum power point more accurately with good robustness.However, the paper only used the particle swarm optimization algorithm for comparative experiments, and there is no benchmark function experiment based on it, which cannot fully explain the universality of the improvement strategy.In [21], an improved sparrow search algorithm is proposed by combining the sparrow search algorithm with the UAV planning problem.The convergence speed and detection ability of the adaptive inertia weight balance algorithm are adopted, and the Cauchy-Gaussian mutation strategy is introduced to improve the anti-stagnation ability of the algorithm.The experimental results show that the algorithm can more effectively solve the UAV route planning problem.In [22], the Kent chaos mapping, Student-t distribution, and Lévy flight strategy are combined with the basic sparrow search algorithm and applied to the unknown load robotic arm parameter identification.
In order to verify the application effectiveness of the power system load parameter identification method combined with the SSA, we proposed an improved sparrow search algorithm by adding tent chaos mapping, Levy flight strategy, Gaussian perturbation after population update, and tent chaos perturbation on the basis of the conventional sparrow search algorithm.In order to verify the performance of the algorithm, the paper performs parameter identification of the proposed traction load model.The method is compared with the conventional sparrow search algorithm and particle swarm optimization (PSO) algorithm.After satisfying the error accuracy, we check the generalization ability of the identified load model.The article verifies the practicality of the ISSA to solve the optimization problem of load modeling parameter identification.

Modeling the Load
The main power component of the traction power supply system of the electrified railroad is the high-power traction motor on the traction locomotive, which is a threephase asynchronous motor by nature and has obvious dynamic load characteristics.In addition, there are many small power components including air conditioning, lighting, and other static loads.Therefore, based on the load characteristics, the traction side can be equated to a composite load model with a parallel static load of motors when viewed from the grid side to the traction side, as in Figure 1.The load model uses a third-order induction motor dynamic load model in parallel with an exponential static load model [23].According to the research of previous scholars, induction motor transient models [24,25] are usually divided into electromechanical transient models, mechanical transient models, and voltage transient models according to the different amounts of parameters, such as dynamic response accuracy of active power and voltage stability index accuracy.A comparison of these three motor models is presented in Table 1.Among the many models proposed above, the third electromechanical transient model is generally suitable for the induction motor transient analysis problem, considering the balance of calculation volume and calculation accuracy.Since the work of this paper mainly focuses on the active power response and reactive power response of the load after the transient voltage disturbance, the dynamic load part of this paper adopts the third induction motor electromechanical transient model, and the static load is described by the exponential static model with few parameters and acceptable accuracy.It is shown schematically in Figure 2.

Dynamic Load Section
In order to increase the calculation accuracy and the calculation volume, and also reduce the number of parameters to be identified, the third model of induction motor adopts the polar coordinate form [26], and its state equation and power output equation are as shown in Equations ( 1) and (2).

•
Equation of state.
, where E′ is the induction electric potential amplitude, δ is the angle be- hind transient reactance, S ω is the angular velocity of stator, r ω is the angular velocity of rotor, m X is the magnetizing reactance, s X is the rotor reactance, r X is the rotor reactance,

/ ( )
T is the load torque constant, and d P and d Q are the dynamic load active and reactive power, respectively.
In Equation ( 1), E′ , δ , ω are the state variables.Since the differential of the state quantity is 0 in the steady-state process, it can be determined that the initial value of the state variable is solved as shown in Equation ( 3) as follows: In the iterative process of solving the optimal parameters of the load model, each iteration involves the calculation of the dynamic process of dynamic load state quantities.In this paper, the state differential equation is solved using the fourth-order Runge-Kutta methods in the iterative process.The principle of the formulation is described in Equation (4) as follows: , ) After solving the equation of state by the fourth-order Runge-Kutta methods, the state variable iteration equation can be obtained as: Therefore, the above Equation ( 5) and the state variable initial value of Equation ( 3) can be combined to iteratively calculate the induction motor in the time domain.

Static Load Section
The exponential function load model is chosen for the static part, ignoring the variation of frequency, and its expression is described as follows: where s P , s Q are the active and reactive power of the static load, respectively, 0 s P , 0 s Q are the active power and reactive power of the steady state before the disturbance, respectively, 0 U is the steady-state voltage, and v p , v q are the active power index and reactive power index, respectively.

Parameter Identification
From the above two subsections, the total power absorbed by the composite load is shown in Equation ( 7) as follows: In summary, the parameters to be determined for this integrated load model consist of five dynamic load parameters The objective function of parameter identification is generally taken as a non-negative monotonic increasing function of the output error, the squared error sum of the measured response of the system, and the calculated response of the model is used as the objective function.
Here, k P and k Q are the measured active power and reactive power, ˆk P and ˆk Q are the model calculates the output active and reactive power, and N is the number of sample data sets.

Traditional Sparrow Search Algorithm
The sparrow search algorithm (SSA) is a new group intelligent optimization algorithms proposed by Xue [16] in 2020 to simulate sparrow population searching.The background and rationale of this algorithm are mainly based on the foraging and anti-predatory behaviors of sparrow populations.This algorithm enables sparrows within a population to update the positions of discoverers, joiners, and perceivers by continuously comparing fitness values until the optimal solution of the objective function is found.
The sparrow search algorithm is divided into three sparrow populations based on their characteristics and adaptability (energy contained in themselves), as described below:

•
Producer: Producers are characterized by high fitness values and a wide search range, and they are responsible for finding food for the entire population and providing foraging directions; • Scrounger: Scroungers have low fitness values, but they always watch the Producers and leave their current location to compete for food if they sense that the Producers have found better food; • Perceiver: Perceivers originate from the Producers and the Scroungers.They can realize the update of the position by perceiving the danger.
The general content of the sparrow search algorithm is as follows.In the case of an optimization problem, the population composed of sparrows is as follows: ... ... ...
Where d denotes the dimension of the parameter to be optimized for this problem, and n indicates the number of sparrows.It can be determined that the fitness value of the sparrow is expressed as follows: 1,1 ) In the SSA, the location of the Producer is updated as described below: exp( ) where t is the number of current iterations, j is the dimension of the parameters, Max item is the maximum number of iterations, , i j X is the position coordinate of the i sparrow in the j dimension, represent the alarm value and the safety threshold, respectively, Q repre- sents a random number that obeys normal distribution, and L represents a matrix with one row and d columns, of which all its elements are 1.There are two situations about the location change of the Producers: when 2

R ST <
, there are no predators around in this environment, and the Producers can carry out extensive search.When 2

R ST ≥
, a sparrow in the population has spotted predators and sends a warning signal to other sparrows in the entire population, at which point all sparrows should immediately fly to other safe locations.
In the SSA, the location of the Scrounger is updated as described below: where P X represents the best location discovered by the Producers, and worst X represents the current global worst position.

( ) T T
A A AA + − = , in which A represents a ma- trix with one row and d columns and where each element in the matrix is randomly assigned 1 or −1.When

/ 2 i n >
, it means that the number i Scrounger with low fitness did not obtain food, so at this time it should change the position to obtain more energy.
In the SSA, for Perceivers who are aware of danger, we usually give 10-20% of the sparrows in the whole population.Its location update is described below.
Here, best X is the current global best position, and ~(0,1) N β is a random number that represents the step control parameter.
[ 1,1] K ∈ − represents the direction of sparrow movement, which is also a step size control parameter and a random number.i f represents the fitness value of the current individual sparrow.ε represents an infin- itesimal constant, which is used to avoid the occurrence of 0 in the denominator of the above formula.There are two situations for the position change of the perceiver.When i g f f > , the sparrow is at the edge of the population and easy to be attacked by predators.
When i g f f = , it means that the sparrow in the middle of the population is aware of the danger and needs to move to the position of other sparrows to avoid the danger.

Tent Chaotic Mapping
This paper uses a tent chaotic map [27] to initialize the sparrow population.A chaotic system is a kind of seemingly disordered but complex and orderly non-random motion, which has the characteristics of randomness, ergodicity, and regularity.Common chaotic maps include a circle map and logistic map, etc. Tent mapping and logistic mapping are topologically conjugate.A related study [28] shows that a tent map has better ergodic uniformity than a logistic map through comparative experiments.Therefore, a tent map is selected as the chaotic map in this study.The iteration formula of tent mapping is as follows: In this study, by setting the parameters of three dimensions and performing tent chaotic mapping between (0,1), it is found that when 0.3 is taken, the variables have relatively good randomness, ergodicity, and regularity, and the iteration diagram of the tent mapping is as shown in Figure 3.

Lévy Flight Strategy
When the Scroungers notice that the Producers have searched for a large amount of food, they will flock to the new location together.This situation will easily cause the joiners to gather to a certain search area at the same time, which will cause the global search ability to decrease and even make the algorithm fall into a local optimum in serious cases.To solve the above problem, we introduce the random step size of Lévy flight [29] when the Scroungers update its position and use its step size s and the uncertainty of its di- rection to enhance the search ability of the algorithm, avoiding premature convergence of the search and falling into a local optimum.The random step s size of the Lévy flight strategy is described as follows.
The location update following the introduction of the Lévy flight strategy by the Scroungers is described as follows.
We selected three sets of parameters for their Lévy flight iterations, and it can be seen from Figure 4 that each parameter has good randomness and uncertainty in its Lévy flight step and direction during the iterative process, realizing the diversity of population position changes.

Gaussian Variation and Tent Chaos Perturbation
Gaussian variation refers to a perturbation of the current global optimal solution searched for by a random number satisfying a Gaussian distribution during each iteration of the search for the optimal solution, which is described by the perturbation formula as follows: where best X represents the global optimal solution of the search at the current iteration number, and best X ′ represents the new value of the global optimum after Gaussian per- turbation.From the Gaussian distribution properties, it is known that the introduction of Gaussian variation enhances the local search ability of the sparrow, which is conducive to inducing the population to jump out of the local optimal solution and find the global optimal solution.
The tent chaos perturbation is derived from the tent chaos mapping, which uses chaotic variables d Z acting on parameter variables in d dimensions and is formulated as follows: min (max min ) where max d and min d are the maximum and minimum values of the d dimension parameters, respectively.The tent chaotic perturbation iterative update equation is described as follows: where X ′ represents the position of the individual sparrow to be disturbed, X is the generated perturbation variable, and new X ′ represents the position of an individual sparrow after chaotic perturbation.
According to the characteristics of Gaussian distribution, the main disturbance area of Gaussian disturbance is the local area near the original individual, which is conducive to the algorithm efficiently finding the minimum point of the current search area.Meanwhile, the tent chaotic disturbance can make the original individual have better randomness and ergodicity after disturbance, which can help the algorithm jump out of local optimization and prevent "premature".Therefore, this paper uses tent chaos perturbation in the early stage of algorithm optimization to search in a large range, which is as much as possible to avoid "premature", and it uses Gaussian perturbation to search in a small range more carefully in the late stage of optimization to speed up the convergence.

Improved Sparrow Search Algorithm Solving Process
According to the improvement of the algorithm in Section 3.2, the parameter identification process of the improved sparrow-based search algorithm in this study is implemented as follows:

•
Step 1: set the number of groups of parameters to be identified (population size N ), number of dimensions of parameters to be identified (dimensions d ), number of Producers PD , number of Perceiver SD , safety threshold ST , maximum number of iterations max Iter , and the objective function .( ) Obj f x ;

•
Step 2: apply tent chaos mapping to initialize the population location is the set of optimal identification parameters sought.
If not, go to step 4 and continue iteratively.
The flowchart of the ISSA is shown in Figure 5:

Data Acquisition
These research data come from the recording and broadcasting data of the PMU measurement device based on the Internet of things at the power grid side obtained by a traction substation in Tibet due to fault disturbance.The traction power supply method in this area adopts the direct power supply method with a return current line; the traction transformer type adopts three-phase V/v [30] wiring, and the locomotive is an D HX 1 [31] AC electric locomotive.The recording device installed in the substation records the instantaneous values of bus voltage and current before and after the substation disturbance to realize the input and output of parameter identification.
In the acquired recording data, the bus voltage perturbation causes the bus voltage to plunge by about 50%, and for the voltage perturbation curve for input data ( ) U k and output data ( ) P k and ( ) Q k , the data after normalization is shown in Figure 6 below:

Parameter Identification
After obtaining the processed voltage input data, the ISSA is used to identify the parameters of the traction load.The flow chart of parameter identification is shown in Figure 7: , the objective function (also the fitness value) is Equation (8): After the algorithm is set, the parameters of the model are identified.The comparison of the active power and reactive power responses with the original data is shown in Figure 8: Figure 8 shows that the output of the identified load model can accurately fit the measured power response curve, especially the reactive power.
According to the research in [7], the load model evaluation index can be evaluated by using the relative error coefficient: where ( ) y k is the measured power, and ˆ( ) y k is the calculated power (it can be active power or reactive power), and if the relative error is less than 5% [7], the load model is acceptable.
The results and errors of the parameters obtained after applying the ISSA to the identification of the model are shown in Table 2.  From Table 2, it can be seen that the identified parameters are all within reasonable ranges, the active and reactive errors are less than 5%, and the model is acceptable.The above results verify the effectiveness of the ISSA for the parameter identification of the composite load model.

Algorithm Comparison
In order to verify that the improved ISSA has the advantages of high accuracy and good convergence, the ISSA is compared with the SSA and PSO algorithms for experiments to jointly identify the parameters of this load model.The ISSA, SSA, and PSO algorithms are set as follows:  3: The global optimal fitness value reflects the change law of the global optimal solution in the iterative process of the algorithm and can reflect the convergence characteristics of the algorithm.Figure 8 takes the iteration times as the x-axis and the global optimal fitness value as the y-axis.The iterative change curves of the optimal fitness values of the three algorithms are shown in Figure 9: From Table 3 and Figure 9, it can be seen that: • PSO: Although the accuracy and convergence speed of PSO is acceptable, it is easy to fall into local optimization.The values of the load model identified by PSO are acceptable, but the relative error of its active power is slightly higher; • SSA: Although the search accuracy and convergence speed of the SSA are better than that of the PSO algorithm, it cannot jump out of the local optimal solution.The values of the load model identified by the SSA are acceptable, but the relative error of its active power is slightly higher; • ISSA: After improvement, ISSA has the ability to search quickly and jump out of the local optimal solution, so its search accuracy and convergence speed are better than the SSA and PSO.Moreover, its active power response, especially reactive power response, is described more accurately.
Therefore, the ISSA has better accuracy and convergence than the SSA and PSO.

Model Generalization Capability Study
The generalization capability of a load model [32] refers to the ability of the model built using known data to interpret the load location change data, also known as interpolation extrapolation capability.A load model has real practical value only when it has good extrapolation and interpolation capabilities.This study evaluates the generalization ability of the load model from two aspects: load response fitting curve and response residual The response residual [33] can reflect the description ability of interpolation and extrapolation of the model and is one of the important indicators to evaluate the applicability of the model.The average value of the response residuals is as follows: where X is the calculated power data of the original response of the model itself, and Ŷ is the actual power data of the interpolation or extrapolation response.Through the recording and broadcasting device, the disturbance data with a voltage disturbance of less than 50% is obtained, and the load model in Table 2 identified by ISSA is tested to verify the interpolation ability of the model.The results are shown in Figure 10: It can be seen from Figure 10 that although there is some error between the interpolation response value calculated by the model and the change curve of the measured value, the response tracking fitting curve is basically consistent, which indicates that the model has good nonlinear interpolation ability.
This study compares the interpolation ability of the models identified by the SSA and PSO algorithms by calculating the values of the response residuals.The mean values of residuals for each algorithm are shown in Table 4. overall model identified by the ISSA has a better interpolation capability, especially the reactive power interpolation response.

Extrapolation Capability
Through the recording and broadcasting device, the disturbance data with a voltage disturbance of less than 50% is obtained, and the load model in Table 2 identified by the ISSA is tested to verify the extrapolation ability of the model.The results are shown in Figure 11: It can be seen from Figure 11 that although there is some error between the extrapolation response value calculated by the model and the change curve of the measured value, the response tracking fitting curve is basically consistent, which indicates that the model has good nonlinear extrapolation ability.
From Table 5, it can be seen that the models identified by different algorithms have different extrapolation capabilities for active and reactive power descriptions, and the overall model identified by the ISSA has a better interpolation capability.Reactive Power(p.u.)

Extrapolated: Reactive power response
In general, the models identified by the ISSA have relatively better generalization capabilities, including interpolation and extrapolation capabilities, than those identified by the SSA and PSO.

Conclusions
In this paper, an improved sparrow search algorithm, ISSA, is proposed and applied to the parameter identification of the traction-integrated load model.The following conclusions are obtained through the application of algorithm identification analysis: 1.The load disturbance data are obtained from the recording and broadcasting data, a suitable traction load model is established, and its parameters are identified by the ISSA.The identification results meet the relative error of the dynamic load model, which verifies the correctness and effectiveness of the ISSA applied to the parameter identification of the comprehensive load model; 2. The results show that compared with the other two algorithms, the ISSA has stronger ergodicity in searching individuals and better performance in convergence speed and accuracy, as well as being able to constantly jump out of local optima; 3. The generalization ability of the identified load model was studied.The results show that the load model identified by the ISSA has good interpolation and extrapolation ability.Therefore, the ISSA can improve the accuracy of load modeling and show practical value.
Since the balance degree of negative sequence voltage on the 220 kV power grid side is considerable, this paper does not consider too much the negative sequence current and harmonic in the research, and future work should analyze and model these two aspects.

Figure 2 .
Figure 2. A schematic of the exponential-induction motor model.

iX• Step 4 : 2 R• Step 7 :
and gen- erate N d -dimensional sparrow individuals;•Step 3: After setting the objective function .( ) Obj f x , the current fitness value i f is calculated for each individual sparrow (the objective function value is taken as the fitness value in this study), and then the current optimal fitness value g f and the current worst fitness value w f are determined in order.Record the positions g X and w X corresponding to the fitness values; Compare the magnitude of the random value and the security threshold ST to update the position of the Producers according to Equation(11).Update the Scroungers position according to the Equation (16) after introducing Lévy's flight strategy.Update the Perceivers position according to Equation (13);•Step 5: According to the size of the current individual fitness value, perform Gaussian variation and tent chaos perturbation on the sparrow position after each iteration position update, and then one iteration is completed; • Step 6: Updating individual fitness values i f of sparrow populations, reorder the new population fitness to determine the current global optimal fitness value g f and the global worst fitness value w f and corresponding positions g X and w X ; Determine if the algorithm has reached the maximum number of iterations max Iter , and if the iterative maximum is reached, the optimal fitness value best f of the sparrow population and its corresponding sparrow position best X are output, where the optimal fitness value best f is the objective solution of the requested opti- mal objective function, best X

Figure 5 .
Figure 5.The flowchart of the ISSA.

Figure 7 .
Figure 7. Schematic diagram of parameter identification process.

Figure 8 .
Figure 8. Fitting curve between measured data and calculated data.
is set up, the parameters of the load model are identified.The identification results and errors are shown in Table

Figure 10 .
Figure 10.Interpolation fitting curve between measured data and calculated data.

Figure 11 .
Figure 11.Extrapolation fitting curve of measured data and calculated data.

Table 1 .
Comparison of three induction motor models.

Table 2 .
Parameter identification values and errors.

Table 3 .
Parameter identification values and errors for the three algorithms.

Table 4 .
Interpolation response residual.From Table4, it can be seen that the models identified by different algorithms have different interpolation capabilities for active and reactive power descriptions, and the