Development of Accurate Lithium-Ion Battery Model Based on Adaptive Random Disturbance PSO Algorithm

The performance behavior of the lithium-ion battery can be simulated by the batterymodel and thus applied to a variety of practical situations. Although the particle swarm optimization (PSO) algorithm has been used for the battery model development, it is usually unable to find an optimal solution during the iteration process. To resolve this problem, an adaptive random disturbance PSO algorithm is proposed. The optimal solution can be updated continuously by obtaining a new random location around the particle’s historical optimal location.There are two conditions considered to perform the model process. Initially, the test operating condition is used to validate the model effectiveness. Secondly, the verification operating condition is used to validate the model generality.The performance results show that the proposedmodel can achieve higher precision in the lithium-ion battery behavior, and it is feasible for wide applications in industry.


Introduction
A battery management system (BMS) is an electronic system used to manage a rechargeable battery or battery pack, and it is widely applied to many applications that use a battery or batteries, such as portable electronic devices, electric vehicles, and power grids [1].To ensure safe and efficient operation, it is essential for a BMS to be able to predict the static and dynamic behavior of the batteries.
Parameter identification is an essential step in battery modeling, and its results directly affect the accuracy and reliability of the model.Identification methods are usually divided into two types [12,13]: online and offline, corresponding to online and offline modeling, respectively.
Online identification methods adjust the parameters of the model in real time based on the condition and current state of the battery.Furthermore, the BMS makes use of such parameters and other information such as current, voltage, and temperature to evaluate the state of charge (SOC), state of health (SOH), and so on [14,15], which are required for real time control [16,17].However, the processor should require higher calculation speed.Offline identifications can adopt mass experimental data that reflects the characteristics of the batteries.As a result, the identified parameters have higher precision and adaptability, and this makes such offline identifications more suitable for battery or battery pack modeling.
At present, there are two major types of offline identifications.The first type of methods is generally called traditional identification methods, like fitting method based on Least Squares [8], subspace identification [16], multiple linear regression method [18], and so on.This type of methods is simple and intuitive; however, the identified parameters have larger errors, and hence it is usually used in applications with lower accuracy demand [19,20].The second type of methods is generally called bionic intelligent optimization algorithms, like PSO [21], genetic algorithm (GA), [22] and so on.Compared with the first type of methods, the second one has obvious advantages in accuracy and reliability, and it has become a popular method for parameter identification.However, when GA is used in parameter identification [22][23][24], there are certain issues that cannot be avoided, such as the high computation time and easily falling into local optimum, i.e., local extremum.The particle swarm optimization (PSO) algorithm may have low precision and often fall into local optimum [24][25][26][27][28][29].According to [25][26][27], the convergence speed and the inertia weight  control can be improved.In [28], both GA and PSO were combined, and in [29], the double quantum PSO was adopted to identify the system parameters, which can increase the traverse ability of particles and simplify the evolution equation without using a velocity vector.However, the abovementioned method may not solve the problem in easily falling into the local optimum.In [30][31][32], the chaotic optimization algorithm (COA) and the improved PSO algorithm were combined to identify the parameters of the battery, load, and solar cell models.In abovementioned papers [30][31][32], upon accepting that the PSO algorithm has an issue of falling into local optimum, the COA is adopted to find the new searching swarm and continue the search, in order to increase global convergence and calculation precision.In [30], the difference between the covariance matrix of the particle location and the predetermined threshold has been used as the basis, while in [32] the difference between the variance of the population's fitness and the predetermined threshold is used as the basis to judge whether the solution falls into local optimum.In addition, the predetermined threshold, which has certain influence on calculation precision, is set according to the experience.
For the development of the battery model using PSO algorithm, the particles may hover around local optimal solution during the iteration process without reaching the real optimal location.For this reason, an adaptive random disturbance PSO (ARDPSO) is proposed and its performances are validated by using a classical single model and a multiplex model as target optimization functions.Also, this algorithm can be used for identifying the parameters of the battery model thus achieving higher calculation precision.
The remaining contents of this paper are arranged as follows: Section 2 introduces the standard PSO algorithm.Section 3 describes ARDPSO algorithm and validates its performance.ECM parameters identification using ARDPSO is illustrated in Section 4. The experiment results and discussion are presented in Section 5. Section 6 concludes this paper.

Standard PSO Algorithm
Particle swarm optimization (PSO) is an evolutionary computation technique, especially searching for optimization for continuous nonlinear, constrained and unconstrained, and nondifferentiable multimodal functions [21].It can optimize a problem by iteratively moving particles around in the search space over the particle's position and velocity.Each particle's movement is updated as better positions in the search space [33].Assume  particles appear in a D-dimensional search space, and the sets can be expressed as where  denoted by   represents the group containing  particles with -dimensional vector in the search space.The term V represents the velocity of the particles.The term   represents the local optimal solution, that is, the individual best solution of particle .The term Gbest represents the global optimal solution.
In each iteration, every particle updates its position and velocity in search space according to its individual best solution   and the global optimal solution of the group Gbest, as follows: where  is the current number of iterations;  1 and  2 are two positive constants called acceleration factors. 1 and  2 are random numbers in the range 0-1.For a particle swarm optimization, a better global search is needed from a starting phase to help the algorithm converge to a target area quickly, and then a stronger local search is used to get a high precision value.Therefore, the modification of improved standard PSO introduces the inertia weight  in formula (2), and its velocity is represented as follows: For getting a high precision solution,  needs to be kept as a variable value, generally a decreasing value.In the improvements of the standard PSO, the linear PSO (LPSO) [25] is very representative, which uses a linearly decreasing inertia weight, given by  = ( max −  min ) (  −   ) +  min (5) in which  is current inertia weight,  min is the minimum value (that is, the final value) of inertia weight,  max is the maximum value (that is, the initial value) of inertia weight,  is the current iteration number, and  is the maximum number of allowable iterations.

ARDPSO Algorithm
3.1.The Principle of ARDPSO.Although some related PSO algorithms in convergence speed and inertia weight control have been reported in the literature, the particles may still encounter problems such as hovering around the real optimal location but being unable to locate it.It means that the local optimal solution of the particle is not updated, then resulting in the global optimal solution not to be updated.Consequently, the distance between the searched solution motion curve of the global optimal location real optimal location initial global optimal location final global optimal location and the real optimal solution will not become closer.For illustrating such a problem, a typical experiment uses the standard PSO and LPSO algorithms to find a solution that minimizes the single mode target function (, ) =  2 + 2 and obviously the real optimal solution is (0, 0).Without loss of generality, the number of particles is set as 50, and the maximum number of iterations is set as 100 in the tests.Note that the performances with different initial settings such as the number of the particles, iterative times, the initial range of the particles, and velocity may produce different results.
Figures 1(a) and 1(c) show the motion curves of the global optimal solution of the standard PSO and LPSO, respectively.On the other hand, Figures 1(b) and 1(d) show the corresponding fitness for standard PSO and LPSO, respectively.
From Figure 1, the following can be seen.
(1) The global optimal solution (Gbest) may stop updating during a certain period of iterations.For example, as shown in Figure 1(d), from the 10 th to 15 th iteration, the global optimal solution stops updating, while the global optimal solution continues to update after the 16 th iteration.
(2) During the iteration process, the global optimal solution of the standard PSO and LPSO gradually tends to the real optimal position (0, 0).However, it is noted that the standard PSO algorithm stops updating the global optimal solution at the 24 th iteration, and its fitness is 1.702324344665708e-05.On the other hand, the LPSO algorithm stops updating the global optimal solution at the 78 th iteration, and its fitness is 1.6916e-12.
(3) The inertia weight  shown in formula (4) in LPSO is the key factor to achieve a better fitness value than standard PSO by increasing the number of iterations to be closer to a real optimal solution.
However, both standard PSO and LPSO would encounter such problems: (1) During the iterations, the global optimal solution may not be updated for each iteration.
(2) The global optimal solution may stop updating before the maximum number of allowable iterations is achieved.
For this reason, the ARDPSO algorithm is proposed to give the particle more opportunities to continue to update its local optimal solution and thus to find the global optimal more accurately.The global optimal solution updating from the ARDPSO algorithm is shown in Figure 2(a).During the iteration process, the moving of the particle will update the global optimal location if a new location has a smaller distance from the real optimal location than the global optimal location.Figure 2(b) illustrates the procedure of the ARDPSO algorithm performance as follows.
(1) In each iteration process, for each particle, if the fitness of its current solution (that is,  +1  ) is inferior to the fitness of its local optimal solution (that is,    ),    will not be updated in this iteration.It is noted that the current solution ; the most likely reason for this is that (a)    is the real optimal solution and (b) the moving speed is too rapid so that  +1  exceeds the real optimal solution.If it is due to the first reason, then the algorithm has already achieved the optimal solution.In the second reason, a relatively slow moving speed needs to be set.
If    in the () ℎ iteration is not updated and if ( +1  ) in the ( + 1) ℎ iteration is inferior to (   ), it means that the current location of the particle  +1  comparing with    has a relatively long distance from the real optimal solution.In this situation, a new solution should be got around    .During the evolution of the global optimal solution, the distance between the global optimal solution and the real optimal solution decreases with increasing iteration time, and the coverage of the random disturbance thus becomes less.In Figure 2(b), it can be seen that the coverage of the ( + 1) ℎ random disturbance is less than that of the () ℎ random disturbance.
A new random location generation function can be defined as formula If the fitness of its current solution (that is,  +1  ) is inferior to the fitness of its local optimal solution, then it will use formula (6) to get a new location.In formula (6), exp(− * /) is used to scale the coverage of the random disturbance, commonly, a is a constant, n presents the current iteration number, and N presents the total iteration numbers.The function rand() is a random value between 0 and 1.
(2) By comparing the fitness of the randomly obtained new solution to that of the abovementioned current solution, the better solution is selected as the new current location of the particle.
(3) If the new current location of the particle is superior to its local optimal solution, thus the local optimal solution is updated.
As above, the process of the ARDPSO is concluded as follows.
Step 1. Get the moving speed: Step 2. Get the new location: And then, for all of the particles, we have the following.

Testing and Analyzing the ARDPSO.
In order to verify the availability and wide applicability of the proposed ARDPSO algorithm, nine benchmarks from the BBO repository are adopted to test its performance, shown in Table 1, in which  1 −  5 are single mode functions,  6 is a step function, and  7 −  9 are multiplex mode functions [34][35][36].
Without loss of generality, the parameters setting for each benchmark in Table 1 are as follows: (i) The dimensions of these functions are all set as 2.
(i) The solution range of the benchmarks is based on Table 1.(ii) Generate the initial position of the particles randomly within its solution range.
In addition, we have the following.
(i) For LPSO and ARDPSO, the maximum inertia weight is set as 0.9, while the minimum one is 0.5.ii.For ARDPSO, the parameter  in the formula ( 5) is set as 1.
Following the above rules, in the MATLAB 2015 (b) environment, standard PSO, LPSO, and ARDPSO are used to test 9 benchmark functions in Table 1, respectively, where each benchmark function was tested 30 times.The number of local optimal solution (that is,   ) updating is shown in Figure 3.The abscissa values 1-9 correspond to benchmarks 1-9, and the ordinate value indicates the number of updating process.Table 2 shows the average updating times for the local optimal solution of the particle.The number of global optimal solution (that is, Gbest) updating is shown in Figure 4. Similarly, in Figure 4, the abscissa values 1-9 correspond to benchmarks 1-9, and the ordinate value indicates the number of updating process.Table 4 shows the minimum error, maximum error, and mean error between the fitness of optimal solution and the real optimal solution.
The performance results as above indicate the following: (1) For local optimal solution (Pbest) from Figure 3 and Table 2, all benchmarks except the step function  6 show that the average updating times by standard PSO and LPSO are  about 5 to 18, respectively.However, the average updating times by ARDPSO are about 64 for  1 −  5 and  8 , while they are about 36 for  7 and  9 .It is found that the benchmark  6 may achieve relatively small updating times in standard PSO, LPSO, and ARDPSO.
(2) For global optimal solution (Gbest) from Figure 4 and Table 3, for all the benchmarks except the step function  6 , the average updating times by standard PSO are about 6 to 12, and they are about 28 to 32 by LPSO.The average updating times by ARDPSO are more than 60 for  1 −  5 and  8 , while they are about 33 for  7 and  9 .Similarly, the benchmark  6 has relatively small average updating times in standard PSO, LPSO, and ARDPSO.
(3) From the error statistics in Table 4, both LPSO and ARDPSO algorithms are superior to the standard PSO algorithm in all performance evaluations.For step function  6 , both LPSO and ARDPSO can find the real optimal result, reaching an average error=0.For  1 − 5 and  8 , APDPSO has much lower errors than LPSO.For  7 and  9 , the relatively small updating times for local and global optimal solution imply that the global optimal solution will not be updated if the real optimal solution is found.
The following can be concluded: (1) For step function  6 , many solutions may get the optimal fitness, and it means that there are more chances to find the real solution.In this situation, the Pbset and Gbest are relatively small for standard PSO, LPSO, and ARDPSO.
(2) From formula ( 2) and ( 4), it can be seen that LPSO adds an inertia weight , compared to the standard PSO. decreases with the increasing iterations, and this gives the particles more chances to update its Pbest as well as the Gbest for better solution.
(3) When the particle could not update its   , the current solution ( +1  ) is not a better solution.A too rapid moving speed may make it happen and thus lose the real optimal solution.In the proposed ARDPSO model, when the particle could not update its   , it can obtain a new solution  +1  rand around the   .If the   is not the real optimal solution, choosing the better one from  +1  and  +1  rand as the new position of the particle can encourage   to be updated.From the results shown in Tables 2-4, it can be seen that if the real optimal position is not found, the Pbest and Gbest updating times by ARDPSO may increase far higher than that of LPSO.

ECM Parameters Identification Based on the Proposed ARDPSO
In this section, the ARDPSO algorithm is used to identify the parameters of the ECM.

Preparation for the Experiment.
The following experiments use NEWARE BTS-4008 as the power battery test system and 18650 NMC batteries (the positive material is LiCoxNiyMnzO 2 ) as the experiment object.The selected 18650 NMC battery was manufactured in Tianjin with a rated capacity of 2000 mAh and a rated voltage of 3.7 V.The rated charging and discharging cut-off voltages are 4.2 V and 2.5 V, respectively.
The open circuit voltage-(OCV-) state of charge (SOC) curve of the above tested battery is obtained according to the method described in [37].
The error criterion is a crucial factor to determine the precision of parameters identification so that the root square error is commonly used for battery models parameters identification [9].However, many factors may influence the precision of the battery model, like the OCV-SOC relationship, the model itself, SOC initial value, and so on.For this reason, a weighted fitness function shown in formula ( 7) is defined as where  represents the total number of measurement intervals,   is the terminal voltage of the battery system that can be obtained by measurement,   is the voltage that is calculated from the battery model, and  represents the  ℎ sampling.In the following experiments,  =  = 0.5.

ECM Model Used by the Experiments.
At present, the simple model, the first-order RC model, and the secondorder RC model are the most applied battery models [9,38] in industry.The simple model shown in (8) has a simple structure and operating principle.The Thévenin ECM shown in ( 9) is known as the first-order RC model.The second-order RC model shown in ( 10) is an extension model from the firstorder RC model.
In above formulae,   represents the terminal voltage of the battery,   () represents the value of the OCV corresponding to the value of the SOC,   represents the load current, and  represents the internal resistance of the battery.In (9),   represents the polarization voltage and its discrete form can be shown using the following formula [39].
where  =     , Δ represents the sampling interval of the system, and   and   represent the polarization resistance and the polarization capacitor, respectively.In formula (10),  1 ,  2 , respectively, represent the electrochemical polarization voltage and the polarization voltage, and the discrete form of  1 ,  2 can be expressed as where  1 =  1  1 ,  2 =  2  2 , Δ represents the sampling interval of the system,  1 and  1 represent the electrochemical polarization resistance and the electrochemical polarization capacitor, respectively, and  2 and  2 represent the concentration polarization resistance and the concentration polarization capacitor, respectively.In formulae (10) and (11),  represents the  ℎ sampling.
For the three ECMs, the dynamic behavior of the battery cannot be effectively modeled by the simple model due to the lack of consideration in the battery polarization effect.For example, when the battery is at rest after it is charged/discharged for a certain time, its SOC value is a constant with 0 current.Then, the simple model will output a constant terminal voltage, but in fact its terminal voltage is variable.The first-order RC model is based on the connection with one RC network in series with the simple model, and the delay characteristic of the first-order RC network is used to simulate the polarization effect of the battery.Therefore, during the battery rest, the first-order RC model can output a variable.However, it is unsuitable for the voltage transient process of the battery.The second-order RC model that connects two RC networks in series with the simple model can simulate the electrochemical polarization and concentration polarization separately.Accordingly, it could better model the dynamic behavior of the battery.In this study, the influence from different parameter identification methods is mainly considered.

Parameter Identification Based on RDPSO Algorithm.
The RDPSO algorithm evaluated with the fitness function shown in ( 7) is used to identify the parameters of the battery model, and the process is indicated in Figure 5. Preparing battery date: current, voltage, SoC, OCV randomly is more superior?
The new x i got of random disturbance Calculating fitness of x i

Calculating Gbest
Initialization:  (iii) Figure 6(c) shows the value of the SOC from each measurement interval.In the test operating condition, the battery rests for two hours so that the terminal voltage of the battery can be used as initial OCV.Consequently, the initial SOC value can be obtained via the OCV-SOC relationship.Other SOC values can be obtained through Amper-Hour integral [8] using the known initial SOC value and the current shown in Figure 6(b).(iv) Figure 6(d) shows the value of OCV corresponding to the SOC value for each measurement interval.

Effectiveness Verification.
The simple, first-order RC, and second-order RC models were used to test the least square method [8], LPSO, and ARDPSO in the parameters identification under the condition shown in Figure 6.In Note that each variable in formula (13) has the same definition as that of formula (7).The above three models based on the same parameters have similar output voltages and error distributions, and only the results of the second-order RC model based on LPSO and ARDPSO are shown in Figure 8.The measured voltage of the battery and the simulation result are shown in Figure 8(a).The simulation errors of the battery model are shown in Figure 8(b).
From the experiment results, the following can be concluded: (1) When the three abovementioned models are identified under the same operating conditions, the global optimal location is further optimized by the ARDPSO with the added random disturbance function, as compared to the LPSO.As a result, the output voltage obtained from the ARDPSO is much closer to the true battery voltage than the LPSO.
(2) Under the same operating conditions, the same identification method is used to obtain the parameters of the three models.The second-order RC model that uses two RC networks to simulate the electrochemical polarization and concentration polarization has more precision than the firstorder RC model with one RC network.However, the simple model does not consider the polarization characteristics; comparing with the simple model, the first-order RC model is more precise.
(3) From Figure 8(a), it is shown that, during the battery discharging process, when the value of SOC is small, the voltage drops relatively fast for the same current due to the battery inherent characteristic.It would result in a large error in this case.Consequently, most of the literatures use only the 20%-80% of SOC information [9,40].It can be seen from Figure 8(b) that LPSO and ARDPSO have relatively small errors as the value of SOC is in between 20% and 80%.
(4) Table 5 reveals that the RMS error, the max error, and the min error of ARDPSO are smaller than those of the LPSO; therefore, it can be seen that the novel method combination with the adaptive random disturbance can produce a better performance.
(5) From Table 5, comparing with the least square algorithm, LPSO and ARDPSO could get lower Acc Error.However, the result of least square algorithm is very sensitive to the initial value setting, and an improper initial setting may cause a large error.

Generality Verification.
To validate the generality of the ARDPSO for each of ECM models, the parameters identified under the test operating condition are used as the input to obtain the output voltage of ECM model.Figure 9 shows the output voltage error distribution in different ECM models.The following can be seen from Figures 9 and 10: (1) For the three abovementioned models, the results of the verification operating condition (as shown in Figure 7) using the parameters identified under test operating condition (as shown in Figure 6) indicate that the output voltage of the ECMs almost coincides with the real battery voltage, shown in Figure 7(b).As a result, throughout the whole simulation process, the maximum error is low enough within a small range, as shown in Figure 9.
(2) Among the three models, the second-order model can better simulate the static and dynamic behavior of the battery.Thus the second-order RC model has smaller maximum error, RMS error, and average error than simple model and the first-order RC model.
(3) The parameters identified by the proposed ARDPSO under one operating condition are suitable for other operating conditions for ECM models, verifying the effectiveness and generality of the ARDPSO.That is, the battery model can be used by BMS to predict the states of the battery.

Conclusions
The PSO has been widely used in many applications like identification of ECM model.It and its extended algorithm such as LPSO could update both local and global optimal solutions by moving particles to achieve the target.However, it is found that the solution either local or global optimum may not keep updating for a period of time during the particles movement.The ARDPSO algorithm is proposed to continue to update the optimal solutions.Test results from multiple benchmark functions have verified that the ARDPSO can improve the updating process for both local and global optimal solutions.Accordingly, the ARDPSO can reach higher solution precision than the standard PSO and LPSO.
As ECM model parameters can affect the static and dynamic behaviors of the battery model, the ARDPSO therefore introduces a new weighted fitness function to identify ECM parameters.Based on the evaluation tool using the maximum error, RMS error, and the average error, it is obvious that the parameters of ECM model have been identified accurately under the test operating condition.Besides, it indicates the ARDPSO promises a better performance than the LPSO.For future work, the black box algorithms such as neural network and support vector machine will be used to model the battery and further compared with the ARDPSO algorithm in the state of charge (SOC) and state of health (SOH).

Figure 1 :
Figure 1: Optimization procedure: (a) the motion curve of the global optimal solution of the standard PSO; (b) the fitness curve of the global solution of the standard PSO; (c) the motion curve of the global solution of the LPSO; (d) the fitness curve of the standard LPSO.

Figure 2 :
Figure 2: The principle of the ARDPSO: (a) the global optimal solution updating; (b) the convergent procedure of ARDPSO.

Figure 3 :
Figure 3: Statistical results diagram of historical local optimal updating (in which the red point presents the standard PSO, the green point presents the LPSO, and the blue point presents ARDPSO; the abscissa values 1-9 correspond to benchmarks 1-9, and the ordinate value indicates the number of updating times).

Figure 4 :
Figure 4: Statistical results diagram of global optimal updating (the abscissa values 1-9 correspond to benchmarks 1-9, and the ordinate value indicates the number of updating times).
x i is set to Pbest i Calculating fitness of new x i got randomly Calculating fitness of x i Calculating Pbest i for x i

Figure 5 :
Figure 5: Flow diagram of system identification based on ARDPSO.

5. 1 .
The Working Conditions.Based on the known capacity and the OCV-SOC curve of the battery, the LPSO and ARDPSO are used to identify the parameters of the abovementioned three ECM models.In this section, two different operating conditions are selected for experiments.The first one is the test operating condition shown in Figures6(a)-6(d), and it is used to validate the effectiveness of the ARDPSO.The second one is the verification operating condition shown in Figures7(a)-7(d), and it is used to validate the generality of the ARDPSO.In other words, the battery model parameters identified by one operating condition can suit other operating conditions.In both operating conditions, four types of battery information are provided, shown in Figures6 and 7.In Figures6(a)-6(d) and 7(a)-7(d), the abscissa is time and the measurement interval is 1 s; moreover, Figures 7(a)-7(d) have the same meaning with corresponding (a)-(d) shown in Figure 6.In these figures, we have the following.(i)Figure6(a) shows the current for each measurement interval.(ii)Figure6(b) shows the battery voltage obtained from each measurement interval.

Figure 9 :
Figure 9: Error distribution of the battery models with the parameters identified by ARDPSO: (a) simple model; (b) first-order RC model; (c) second-order RC model.

Figure 10 :
Figure 10: Error results of the battery model parameters by ARDPSO.

Figure 10
Figure 10 concludes the statistical results of the maximum error, RMS error, and the average error.The following can be seen from Figures9 and 10:(1) For the three abovementioned models, the results of the verification operating condition (as shown in Figure7) using the parameters identified under test operating condition (as shown in Figure6) indicate that the output voltage of the ECMs almost coincides with the real battery voltage, shown in Figure7(b).As a result, throughout the whole simulation process, the maximum error is low enough within a small range, as shown in Figure9.

Table 2 :
Average updating times of the historical local optimal solution.

Table 3 :
Average updating times of the global optimal solution.

Table 4 :
Comparison of performance results.

Table 5 :
Error result of least square, LPSO, and ARDPSO in test operating conditions.