Two-Layer Optimization Model for the Siting and Sizing of Energy Storage Systems in Distribution Networks

: One of the most important issues that must be taken into consideration during the planning of energy storage systems (ESSs) is improving distribution network economy, reliability, and stability. This paper presents a two-layer optimization model to determine the optimal siting and sizing of ESSs in the distribution network and their best compromise between the real power loss, voltage stability margin, and the application cost of ESSs. Thereinto, an improved bat algorithm based on non-dominated sorting (NSIBA), as an outer layer optimization model, is employed to obtain the Pareto optimal solution set to offer a group of feasible plans for an internal optimization model. According to these feasible plans, the method of fuzzy entropy weight of vague set, as an internal optimization model, is applied to obtain the synthetic priority of Pareto solutions for planning the optimal siting and sizing of ESSs. By this means, the adopted fuzzy entropy weight method is used to obtain the objective function’s weights and vague set method to choose the solution of planning ESSs’ optimal siting and sizing. The proposed method is tested on a real 26-bus distribution system, and the results prove that the proposed method exhibits higher capability and efficiency in finding


Introduction
In recent years, distribution networks using renewable energy such as wind power and photovoltaic power have become one of the research hotspots at home and abroad. After the increasing distributed generations (DGs) are connected to the distribution network, the randomness and intermittency of renewable energy access will have a great impact on the safe and stable operation of distribution networks. Therefore, incorporating ESSs within the distribution network has a positive impact on improving system reliability, such as improved overall energy efficiency, maximized voltage stability, and reduced real power losses [1,2]. Moreover, ESSs are of promising technology with peak shaving, and they facilitate the integration of high penetration levels of renewable energy sources. Despite the above advantages, if the site and size of ESSs are improperly selected, it may cause insufficient voltage stability and economic benefits will decrease [3][4][5]. Therefore, one of the main problems associate with the use of ESSs in distribution networks is to find their best site and size in order to maximize their impact on the grid.
Considering the multidimensional, constrained, and nonlinear characteristic of optimal siting and sizing, the solution techniques for EESs are attained via optimization methods. Much research information on decision-making. Based on the above, NSIBA and the method of fuzzy entropy weight of vague sets are integrated into a two-layer optimization model that considers several constraints in the process of ascertaining ESSs' most cost-effective scheme and gives the planner the capability of making the final decision. The effectiveness of the proposed algorithm is validated using the real 26bus system. To test the effectiveness and reliability of the proposed method, the results are compared with NSGA II and NSPSO; moreover, it demonstrates and verifies that the proposed algorithm can improve a distribution network's economy and power quality. The rest of this paper is structured as follows. The mathematic formulation is presented in Section 2. The outer layer optimization model of the NSIBA optimization algorithm is constructed in Section 3. Section 4 presents the fuzzy entropy weight of vague set as the internal optimization model. The simulation results are discussed in Section 5, and the conclusion is given in Section 6.

Mathematical Problem Formulation
The main objective of current work is to find out the optimal siting and sizing of ESSs in distribution networks, together with minimum network power loss, voltage stability margin, and the application cost of ESSs. Each of the above factors can be considered as an objective function (OF) that is subject to equality and inequality constraints as well as to boundary restrictions imposed by the planner. In mathematical terms, a multi-objective optimization problem can be formulated as follows.
Here x = (x1, x2, …, xn) represents the solution to the n-dimensional problem; OF is the target space to the m-dimensional problem; h(x) and g(x) are the equality and inequality constraints, respectively. It is impossible to find a solution to minimize all objectives at the same time when the OFs are in conflict; therefore, the dominance relationship and Pareto optimal solution set are introduced as follows.

Objective Function
Considering the evaluation indexes of network power loss, voltage stability margin, and the application cost of ESSs, the mathematical model of multi-objective optimization of ESSs is established as follows: It is impossible to guarantee that all objective functions can achieve a relative minimum at the same time. Therefore, the Pareto optimal solution set is introduced to search the frontier solution in the target region.
(1) Network power loss. In the process of ESSs' siting and sizing, minimize the power loss as much as possible is the first and foremost objective. The power losses depend on the impedances of line and transformer. Incorporating ESSs within the distribution network has a positive impact on reducing the network power loss and improving the electric efficiency of the distribution network. Therefore, under ESSs' optimal configuration, the scheme with low total network power loss should be considered. The equivalent models of transformer and line are shown in Figure 1. The objective function of network power loss can be expressed as Here, Ploss is the total losses of the distribution network. Ptran is the rated active power of the transformer. Pline is the active loss of line. Rij is the resistance of the distribution line connecting the ith and jth buses. Ui is the ith bus voltage. N is the total number of buses in the distribution network. θi is the ith phase angle. Pi and Qi are the active and reactive powers of the ith bus, respectively.
(2) Voltage stability margin. The bus voltage of the network often experiences fluctuations with the increase of load and DGs. ESSs connected to the network properly are conducive to improving the overall voltage profile. The objective function of the voltage stability margin can be stated as follows: Here, Ue is the rated voltage. Up is the rated voltage allowable deviation, namely, Up = 0.05. (3) Application cost. The ESSs' application cost includes investment and operation costs. In order to pursue a higher economy, smaller operating and investment costs must be considered. Therefore, the general formula of the objective function can be described as follows:

Equality Constraints
Due to the advantages of fewer iterations, simplicity, and flexibility, the back/forward sweep method is widely adopted in distribution networks. Therefore, under operating frequency conditions, the conserved formula of the active and reactive powers at a certain bus can be expressed as Here, ESS i P is ESSs' output active and reactive powers of the ith bus, respectively. PLi and QLi represent the active and reactive powers losses of the load, respectively. Gij and Bij are conductance and susceptance, respectively. θij is the angle of the distribution line connecting the ith and jth buses.

Inequality Constraints
The node voltage, ESSs installed capacity, branch current, and SOC are bounded between two extreme levels forced by physical limitations.
Here, ESS i S is the total capacity of the ith ESSs. Iij is the current of the distribution line connecting the ith and jth buses. In order to protect the service life of the ESSs, the SOC must be limited to prevent the overcharge and overdischarge of ESSs.
Here, Pmax and Pmin are the maximum and minimum power of ESSs. In addition, in order to extend the useful life of the ESSs, we set the usage range of the ESSs' SOC as 20% to 90%.

Outer Layer Optimization Model
ILSSIWBA is a swarm intelligence optimization algorithm based on local iterative search and random inertia weights proposed by Chao Gan [21]. The basic framework of ILSSIWBA is similar to the previous bat algorithm (BA), which is popular in solving optimization problems. However, differing from the BA, ILSSIWBA introduces the local optimal solution and the global inertia weight, which can make the optimal solution more stable. This paper combines ILS strategy, SIW strategy, and balance strategy with a non-dominated sorting strategy. NSIBA is proposed to obtain the Pareto optimal front, which can solve the multi-objective optimization problem with faster convergence speed. This algorithm is described in Figure 2.

ILS Strategy
ILS strategy is introduced to make the target jump out of the local optimal solution [22], which is defined as follows.
Here * x is the current optimal solution; * * x is the perturbation solution. Because multiobjective functions cannot be compared directly, the judgment condition is improved: it is equivalent to single-objective comparison when there are at least two OFs values that are less than the comparison value.

SIW Strategy
Stable solutions are obtained by introducing the SIW strategy, which updates the bat pulse frequency, position, and velocity as follows [23].
Here, f is the frequency of bat; fmax and fmin are the maximum and minimum value limits of the frequency; t i X and t i V are the position and velocity of bat in the tth iterations, respectively; ω is the random inertia weight; μmax and μmin are the maximum and minimum influencing factors of the inertia weight, respectively; σ is the deviation coefficient; rand is a uniformly distributed random number. This strategy introduces a weight coefficient for the bat speed update and adjusts the weight by using random variables, which is beneficial to solve the problem of optimal solution instability in traditional BA.

Balance Strategy
In order to balance the local and global solutions, new emissivity and volume update formulas are adopted, in which the pulse emissivity controls the search of the bat in local and global, and the volume controls the acceptance of the new solution. The volume and pulse emission rate are updated by Here, A t is the bat volume in t iterations; r0 and A0 are the initial pulse emissivity and volume, respectively; r ꝏ and A ꝏ are the maximum values of initial pulse emissivity and volume, respectively.

Non-Dominant Sorting and Elite Preservation Strategy
Fast non-dominant is a key feature in stratifying the population according to the level of the Pareto optimal solution set and making the target closer to the Pareto optimal front. The procedures of non-dominant sorting and elite strategy are given as follows [24]: Crowding distance is introduced to characterize the distance between two bats, to make the distribution of bats more uniform in space, which can be formulated as follows: (5) Elite strategy. Elite preservation strategy combines the parents and individuals generated by improved bat algorithm and genetic algorithm to form a set with the size of 3N and selects N individuals with better performance to form the offspring generation. In the selection operation, the solution with the smaller non-dominated Rank and the bigger crowding distance have priority to be chosen than others when the non-dominated Rank is the same. The main idea of elite strategy is given in Figure 3.

Internal Layer Optimization Model
After the Pareto optimal front is obtained through the outer optimization model, the optimal solution can be selected by setting different weights. Due to the difference and ambiguity of weight caused by different design schemes, the fuzzy entropy weight of the vague set is used to reduce the errors caused by decision making.

Analytic Hierarchy Process (AHP)
AHP is a subjective weight method, which determines the weight through the historical experience or subjective preference of the decision made. The calculation process is given below.
(1) Construct a judge matrix. In order to make the scheme have a unified standard, the numbers 1-9 and their inverses are used as scales to define the judge. The comparative important scale of criteria is given in Table 1. Two factors are equally important 3 One factor is weakly important than another 5 One factor is strongly important than another 7 One factor is demonstrably important than another 9 One factor is absolutely important than another 2, 4, 6, 8 The median of two adjacent judgments (2) Arithmetic mean estimated weight vector.
Here, ωAHP is the weight value corresponding to each OF; the elements in the judgment matrix are normalized and averaged to obtain the weight vector.

Fuzzy Entropy Weight
The entropy weight method uses information entropy to describe the objectivity of information, and objectively determines the weight of attributes based on the degree of information difference in the decision matrix. The calculation process is given below: (1) Normalization is represented as follows: (2) The entropy is represented as follows:

Fuzzy Entropy Weight of Vague Set
Due to the ambiguity and uncertainty of the information, the linear weighting method often loses some useful information in multi-objective decision making. Therefore, Gau proposed the concept of the vague set based on the fuzzy set [25], which divides the membership function of each element into support and opposition. The vague set of A can be formulated as follows: Here, U is the universe of discourse; ui denotes a generic element of U; tA(ui) is lower bound on the grade of membership on the evidence for ui and fA(ui) is a lower bound on the negation of ui derived from the evidence against ui. πA is the degree of hesitation of the vague set A. The larger the value of πA is, the more unknown information it may contain. The vague set can be interpreted as tA(ui) and fA(ui) are the number of votes in favor and against; πA is the abstention.
As a generalization of fuzzy sets, vague fuzzy sets consider both membership and nonmembership information, so they can more fully express fuzzy information in multi-objective decisions. The decision process is as follows: Here, P + and P − are the best and worst solution for each attribute to reach candidate solution, respectively; the true and false membership of Pij relative to the positive ideal index P + and the negative ideal index P − are shown in (35). The comprehensive vague membership is shown in (36).
(2) Determine the comprehensive vague value of each ideal solution in the Pareto optimal front, which can be described as follows: (3) Sort the solutions by a new score function to choose the best solution.
The new score function not only emphasizes the difference between true and false membership degrees but also strengthens the influence of unknown information on decision-making. Obviously, the greater Si is, the more credible of the scheme is. The flowchart based on the two-layer optimization model is shown in Figure 4.  When the optimal configuration analysis of the obtained scheme is carried out, the comprehensive evaluation mainly based on the system's operating economy has attracted much attention. Therefore, the three OFs constructed in the study are converted into economic indicators, and the economic optimization configuration of ESSs are analyzed by calculating the following costs [26][27][28] Here, G1 is the cost of network power loss; Cep is the cost of unit network power loss; Ploss,t is the network power loss at time t; Δt is equal to 1 h; G2 is the cost of the voltage stability margin; Ddev,t is the voltage stability margin at time t; ω is the conversion coefficient between the total line loss and voltage stability margin.
ESSs in distribution networks not only have fast corresponding speed but also have the advantages of "spread arbitrage", which can alleviate the power supply tension during the peak of power consumption. Therefore, in conjunction with peak and valley power prices, it can reduce economic costs.
Here G4 is the benefits of economical recycling; Cebuy is the peak-valley unit price spread; the loss of the ESSs in its life cycle is expressed as nloss, which is equal to 0.95.

Simulation Results and Discussion
The proposed algorithm is implemented on a 26-node system of a 10-KV substation that is depicted in Figure 5. The system contains 12 transformers, where the transformer nodes are set as candidate nodes for ESSs.

Algorithm Performance Analysis
NSGA II, NSPSO, and NSIBA are applied in this paper to compare the performance for solving multi-objective optimization problems. Algorithm parameters are set as follows: The maximum iteration is 300; the bat population is 50; the initial and maximum value of the pulse emissivity are 0.1 and 0.7, respectively; the initial and minimum values of the volume are 0.9 and 0.6, respectively; maximum and minimum values of inertia weight are 0.9 and 0.4, respectively; the coefficient of deviation is 0.2; lifetime a = 10 year; discount rate r = 10%; investment cost Ins DG C = 172 USD/kW and operating cost Oper DG C = 257 USD/kW; the cost of unit network power loss Cep= 0.07 USD/(kW h); the peak-valley unit price spread Cebuy = 0.083 USD/(kW h).
The convergence curves of OFs are depicted in Figure 6. As shown in Figure 6, OF1 and OF2 of NSIBA have a sharp decline in the first 50 iterations. It can be affirmed that NSIBA has quicker convergence speed and better search accuracy than other algorithms.

Generations
Generations Generations In order to explain the performance of NSIBA further, Figure 7 shows the Pareto optimal front obtained by different approaches. Comparing the OF values, the solutions by NSIBA has a significant advantage than other algorithms, which not only distributes more uniformly in the Pareto optimal front but also has a wider range of solution set distribution. Therefore, the NSIBA algorithm can search for more possible solutions and avoid falling into a local optimum.

Scheme Economic Analysis
In the actual operation of a distribution network, affordability is gradually increasing as the load increases. Therefore, this section compares the economics of ESSs with capacity expansion equipment such as transformers and analyzes the optimal configuration of the ESSs, which is listed in Table 2. In order to realize the comparison in the same conditions, it supposes that the capacity of the transformer is increased from 315 to 500 kVA on the bus where the ESSs are installed. The expansion of the transformer is generally increased by 50% of the original equipment capacity. The transformer upgrade cost is about USD 6857.14, which is represented by G3 in Table 2. From the comparison results, we can see that G1 and G2 within ten years has a small difference in the two scenarios. Although the investment cost G3 of ESSs is higher than that of traditional transformer expansion, the cost-saving through "ESSs high storage and low generation" is enough to make up for investment cost and achieve capital recovery. Therefore, in comparison to the expansion of traditional transformers, the configuration of ESSs can bring more economic benefits to the network.

Scheme Effectiveness Analysis
In order to enrich the background of the example, 200 kW photovoltaic systems are connected at bus 6 and 7, and a 200-kW wind power system is connected at bus 24 and 25. The typical daily characteristic curves of load, photovoltaic, and wind power are shown in Figure 8. The typical daily system voltage at each bus is shown in Figure 9. Comparing Scenario 1 with Scenario 2, due to the randomness of power output by DGs, the system voltage has a large deviation, and the total line loss of the system has increased greatly. Comparing Scenario 2 with Scenario 3, the adjustability of ESSs can greatly suppress the voltage fluctuation, which is of great significance to the safe and stable operation of the system. Thus, the ESSs performs very well in terms of suppressing bus voltage fluctuation and load fluctuation. The voltage and current of each bus at 20 h are shown in Figure 10. The voltage amplitude increases toward the value of the rated voltage, and the current amplitude is smaller than the corresponding value in Scenario 1. The voltage and current quality of the network with ESSs have been significantly improved, which not only help improve the accommodation ability of the network, but also verify the practicability of the proposed two-layer optimization model.

Conclusions
According to the network power loss, voltage stability margins, and application costs of ESSs, this paper proposes a two-layer optimization model for the optimal siting and sizing of ESSs. Based on a theoretical analysis and simulation verification, the following conclusions can be drawn: (1) This paper establishes a two-layer optimization model by integrating NSIBA and fuzzy entropy weight of the vague set. The proposed NSIBA can keep the population evolving during the search and jump out of the local optimal solution due to the integration of ILS strategy, SIW strategy, and balance strategy within it. The proposed fuzzy entropy weight of the vague set can strengthen the influence of unknown information on decision-making due to the integration of the new score function. The two-layer optimization model can fairly and reasonably give decision-makers a comprehensive plan. (2) The proposed method is tested on a real 26-bus distribution system. The convergence curve and Pareto optimal front show that NSIBA has the quicker convergence speed and better search accuracy than NSGA II and NSPSO. The simulation in different scenarios demonstrates that the proposed scheme not only can achieve arbitrage of USD 2612.22 by "ESSs high storage and low generation", but also improve the power quality when accessing the DGs. Hence, it may be concluded that the proposed two-layer optimization model is a good choice over the other algorithm for determining the optimal siting and sizing of ESSs in distribution networks. It will be very useful for decision-making in the optimal configuration of ESSs using simulation tools.
Author Contributions: F.Z. and T.S. conceived and designed the study. L.Z. performed the simulation. P.Z. and X.X. provided the simulation case. Y.C. and T.S. wrote the paper. F.Z. and T.S. reviewed and edited the manuscript. All authors have read and agreed to the published version of the manuscript.
Funding: This research was funded by the science and technology project of State Grid Hubei Electric Power Co. LTD through grant number 5215C018002H and the Natural Science Foundation of Fujian Province, China, through grant number 2019J01249.

Conflicts of Interest:
The authors declare no conflict of interest.