The Optimal Allocation and Operation of an Energy Storage System with High Penetration Grid-Connected Photovoltaic Systems

High-penetration grid-connected photovoltaic (PV) systems can lead to reverse power flow, which can cause adverse effects, such as voltage over-limits and increased power loss, and affect the safety, reliability and economic operations of the distribution network. Reasonable energy storage optimization allocation and operation can effectively mitigate these disadvantages. In this paper, the optimal location, capacity and charge/discharge strategy of the energy storage system were simultaneously performed based on two objective functions that include voltage deviations and active power loss. The membership function and weighting method were used to combine the two objectives into a single objective. An energy storage optimization model for a distribution network considering PV and load power temporal changes was thus established, and the improved particle swarm optimization algorithm was utilized to solve the problem. Taking the Institute of Electrical and Electronic Engineers (IEEE)-33 bus system as an example, the optimal allocation and operation of the energy storage system was realized for the access of high penetration single-point and multi-point PV systems in the distribution network. The results of the power flow optimization in different scenarios were compared. The results show that using the proposed approach can improve the voltage quality, reduce the power loss, and reduce and smooth the transmission power of the upper-level grid.


Introduction
In the context of the rapid development of the global economy, the continuous growth of energy demands, the gradual depletion of fossil fuels and the increasing environmental pressure, many countries and regions are vigorously advocating for the reform and transformation of current energy structures, and have introduced various policies to promote the use and development of clean energy. Photovoltaic (PV) power generation technology has been rapidly developed and applied during this period [1][2][3]. With the increasing penetration rate and number of distributed PVs in the distribution network, adverse effects such as power backflows, the bus voltage being out of limit, power loss increases, harmonic increases and power supply reliability reductions may occur [4,5]. At the same time, the PV power output varies over time due to intermittency, randomness and volatility. The load changes also have time sequence characteristics, and the peak and valley of the PV power output and those of the load in the distribution network are not completely matched [6]. The double In the literature [24,25], a charging and discharging operational strategy optimization for energy storage was realized, but this study was carried out under the premise that the location and capacity of the energy storage are pre-established. Mahani et al. [26] optimized the planning and operation of the energy storage location and capacity in the distribution network with a high penetration rate of renewable energy to realize the peak adjustment and reduce the reverse power flow and price arbitrage. Grisales-Noreña et al. [27] selected the minimum energy loss as the optimization objective, and used a master-slave strategy to solve the problem to obtain the optimal position and the charge/discharge operational control strategies for energy storage batteries and capacitors banks. In the literature [26,27], the voltage regulation of a distribution network before and after the configuration of the energy storage system was not considered.
To our knowledge, there are few studies on energy storage optimization allocation and operation that simultaneously consider the joint optimization of the installation location, capacity, charge/discharge power of an energy storage system, and power flow. Usually, the studies are only part of them, while other factors are given in advance. In addition, in the present study, the location of the PV system in the power grid is composed of one point or several points, and multi-point cases are not considered.
The objective of the paper is to propose a new optimization approach that used an improved particle swarm optimization algorithm to determine the energy storage systems' locations, capacities and charge/discharge strategy simultaneously to minimize the voltage deviations and total active power loss for the access of high penetration single-point and multi-point PV systems in the distribution network. The membership function and weighting method are used to combine the two objectives into a single objective. The differential inertial weights and global guided cross search mechanism are utilized in the improved algorithm. The Institute of Electrical and Electronic Engineers (IEEE)-33 bus feeder system is examined. The power flow optimization results under different scenarios in the distribution network are achieved to verify its efficiency.

Mathematical Model of the Optimal Allocation and Operation for an Energy Storage System
In order to optimize the power flow of the distribution network with PV systems, the optimization objectives are to minimize the bus voltage deviation and active power loss; the power balance, bus voltage, energy storage power, charge/discharge balance and remaining capacity of the energy storage system are the constraints. The optimal location and capacity of the energy storage system are determined by the improved particle swarm optimization algorithm.

Objective 1: The Minimization of Bus Voltage Deviations
Voltage stability is one of the main purposes of power flow optimization. It is of great significance for the safe and stable operation of the distribution network to maintain the voltage of each bus at a certain level in the network. PV power generation has randomness and volatility. Although PV power under certain penetration rates improves the voltage quality to some extent, it also aggravates the deviation of the bus voltage. Therefore, the sum of the voltage deviations of all buses in the distribution network for one day was taken as the objective function of the optimal energy storage allocation, as shown in (1).
where N is the total number of buses in the distribution network, and T is the total number of selected moments in this research, T = 24. U i,t is the voltage amplitude of bus i at t, and U N is the rated voltage of the distribution network. The average value of active power loss during the study period was selected as the objective function of optimization, as shown in (2).
where P t,L is the active power loss at time t. The calculation method is the same as that in [28].

The Double Objectives Become a Single Objective
In this paper, our objective is that the voltage deviation of each bus and the power loss of the grid be as small as possible, so it is necessary to consider objective 1 and objective 2 at the same time. Since the dimensions of objective 1 and objective 2 are different, they cannot be added directly. Therefore, this paper will use the membership function of fuzzy mathematics to calculate the membership degree of each objective [29], then the weighted method will be used to transform objective 1 and objective 2 into a single objective. The general idea is: first, the optimal solution of each single objective under all constraints is obtained, and then the membership degree of each objective is calculated according to the membership function, and then the memberships are added linearly. Finally, the multi-objective optimization problem is transformed into a single objective optimization problem. According to the above method, the problems of the different dimensions for objective 1 and objective 2 are solved at the same time. The steps are as follows: In the first step, the optimal solution of objective 1 under the constraint condition is first obtained; that is, the minimum value of the voltage deviation is obtained as the lower limit of voltage deviation f 1L . At the same time, the value of objective 2, active power loss, is calculated as the upper limit of power loss f 2U .
In the second step, the optimal solution of objective 2 under the constraint conditions, the minimum value of active power loss, is obtained as the lower limit of active power loss f 2L . At the same time, the value of objective 1, the value of voltage deviation, is calculated as the upper limit of voltage deviation f 1U .
In the third step, fuzzy processing is carried out for the objective function, and an ascending half trapezoid distribution is selected for the membership function, as shown in (3). The distribution is shown in Figure 1.
where µ( f j ) is the membership of the objective function j.  The membership of objective function (1) and (2) is then linearly weighted to obtain the single objective function as follows:

Power Balance Constraints
where P sum is the active power injected into the distribution network from the upper-level grid; m is the number of energy storage systems; P c,j is the power of the j-th energy storage system, where the specified positive value is the discharge; n is the number of grid-connected PV systems; P g,k is the power of the k-th PV systems; N is the number of buses in the distribution network; and P l,i is the load power of bus i.

Bus Voltage Constraints
where U min is the minimum value of the specified bus voltage, which is 0.95 U N ; U max is the maximum bus voltage, which is 1.05 U N ; and U i,j is the voltage value of bus i at time j.

Power Constraints of the Energy Storage System
P c,min < P c,j < P c,max where P c,min is the minimum power of the energy storage system; P c,max is the maximum power of the energy storage system; and P c,j is the power of the energy storage system at j moment.

Charge/Discharge Balance Constraints of the Energy Storage Systems
T t=1 P c,t * ∆t = 0 where P c,t is the power of the energy storage system at moment t, and ∆t is the time interval between two moments.

The Remaining Power Constraints of the Energy Storage System
where E ct is the remaining power of the energy storage unit at moment t, and E cmax and E cmin are the upper and lower limits of the remaining capacity of the energy storage system, which are 0.9 and 0.2 times the rated capacity of the energy storage system, respectively. That is, the state of charge (SOC) values of the energy storage system range from 20% to 90%.

Solving Algorithm
In the particle swarm optimization (PSO) algorithm, the particles update their speed and position by tracking their own historical optimal solution and the global optimal solution of the population. The update method is shown in (10) and (11) [30].
where w is the inertial weight, c 1 and c 2 are the acceleration factors, and r 1 and r 2 are random numbers between 0 and 1; p id is the d-dimensional component of the i-th particle in the optimal position vector of the k-th iteration, and p In order to solve the problem of the fixed inertial weights and make it easy to fall into a local optimal solution in the iteration process of a conventional PSO algorithm, an improved PSO algorithm is proposed in this paper.

Differential Inertial Weight
The balance between the exploration ability and development ability of the intelligent optimization method is the key to its successful operation, and a reasonable value of inertial weight ω is an important way to achieve balance between these two capabilities in the PSO algorithm. In the PSO algorithm, the value of inertial weight ω has different forms, such as fixed weight, which is a better value obtained from experience; time-varying weight, which decreases as the number of iterations increases; fuzzy weight; which is a kind of dynamic adjustment based on fuzzy systems; and random weight, which is a random value in a certain range [31]. However, none of these inertial weights take into account the characteristics of the particles in the iteration process.
The difference between the position of a particle and the corresponding position of the global optimum particle of the population reflects the difference between the current solution and the optimum solution. Thus, increasing the value of ω when the difference is large can increase the flight speed and improve the global search ability of the particle. When the difference is small, reducing the value of ω can reduce the flight speed and improve the particle's local search ability. In the improved particle swarm optimization (IPSO) of this paper, the size of ω is linearly adjusted according to the difference between the position of the particle and the corresponding variable of the optimal particle position of the population.
The difference between the position of the i-th particle in the k-th iteration and the global optimal position of the population is shown in (12). The inertial weight of the i-th particle in the update speed of its position in the k-th iteration is defined as the differential inertial weight, as shown in (13).
where x max and x min are the maximum and minimum values of the particle position variables, respectively and ω max and ω min are the maximum and minimum inertial weights, respectively. The differential inertial weight curve when the maximum inertial weight ω max is set to 1 and ω min is set to 0.5 is shown in Figure 2. respectively and ω max and ω min are the maximum and minimum inertial weights, respectively.
The differential inertial weight curve when the maximum inertial weight max ω is set to 1 and min ω is set to 0.5 is shown in Figure 2.

Global Guided Cross-Search Mechanism
Because PSO is prone to falling into a local optimal solution during its iterations, a global guided cross-search mechanism to update the location of the particles is introduced [32] using (14). Guiding the new solution to a global optimal solution improves the convergence speed of the algorithm: where rand is the random value of the uniform distribution between 0 and 1, and b is a random value between −1 and 1. Increasing the value of the coefficient cr is conducive to enhancing the development ability of the algorithm, while reducing the value of the coefficient cr is beneficial to the search ability of the algorithm. In this paper, cr is 0.6 [28].

Coding Method
Because the charging and discharging power of the energy storage system at each moment is affected by the load change, the PV power output and the charge/discharge of the energy storage system have constraints. In the optimization process, the installation location of the energy storage system and the charging or discharging power at each moment are taken as the optimization variables, and the coding method is as follows.
The position x of each particle in the solution process is the variable to be optimized, as shown in (15), including the installation location and the charging or discharging power at each time of the energy storage system in the distribution network; each particle is a (m + 1) *T-dimensional vector.
where lj x is the position variable of the j-th energy storage system, which is an integer; j t + p x is the charging or discharging power of the j-th energy storage system at moment t; m is the number of energy storage systems in the distribution network; and T is the total number of moments.

Global Guided Cross-Search Mechanism
Because PSO is prone to falling into a local optimal solution during its iterations, a global guided cross-search mechanism to update the location of the particles is introduced [32] using (14). Guiding the new solution to a global optimal solution improves the convergence speed of the algorithm: where rand is the random value of the uniform distribution between 0 and 1, and b is a random value between −1 and 1. Increasing the value of the coefficient cr is conducive to enhancing the development ability of the algorithm, while reducing the value of the coefficient cr is beneficial to the search ability of the algorithm. In this paper, cr is 0.6 [28].

Coding Method
Because the charging and discharging power of the energy storage system at each moment is affected by the load change, the PV power output and the charge/discharge of the energy storage system have constraints. In the optimization process, the installation location of the energy storage system and the charging or discharging power at each moment are taken as the optimization variables, and the coding method is as follows.
The position x of each particle in the solution process is the variable to be optimized, as shown in (15), including the installation location and the charging or discharging power at each time of the energy storage system in the distribution network; each particle is a (m + 1) *T-dimensional vector.
where x lj is the position variable of the j-th energy storage system, which is an integer; x pj+t is the charging or discharging power of the j-th energy storage system at moment t; m is the number of energy storage systems in the distribution network; and T is the total number of moments.
In this case, the differential inertial weight in the update speed of the different location variables for each particle is calculated as follows: Sustainability 2020, 12, 6154 where ∆x li is the difference between the position variable of the energy storage system and the global optimal position of the population for the i-th particle in k iterations; ω (k) li is the differential inertial weight of the i-th particle in the update speed of the energy storage system position variable in k pi is the difference between the power variable of the energy storage system and the global optimal position of the population for the i-th particle in k iterations; and ω (k) pi is the differential inertial weight in the update speed of the power variable of the storage system for the i-th particle in k iterations.
x lmax and x lmin are the maximum and minimum values of the location variables of the energy storage system, respectively.
x pmax and x pmin are the maximum and minimum power variables of the energy storage system, respectively.

Solution Steps
Step 1: read the network data. Read the parameters of the distribution network to be optimized, the locations of the grid-connected PV systems, the load demand and the PV power output data at each moment, etc.
Step 2: initialize the parameters of the IPSO. Create the initialization conditions, determine the total number of particles N, the particle dimensions D, the number of iterations k = 1, and the maximum number of iterations MaxDT, as well as the size of parameters ω max , ω min , c 1 , c 2 , r 1 and r 2 .
Step 3: initialize the particles. According to the coding method shown in (15)-(17), the energy storage position variable uses integer encoding, the energy storage power variable uses real encoding, and the particle position x i and velocity v i are randomly generated within the constraint range of the network parameters. The particle position vector includes the installation position variable of the energy storage system and the charging or discharging power variable of the energy storage at each moment. Correspondingly, the particle velocity vector includes the velocity of the energy storage system position variable (which should be an integer) and the velocity of the energy storage system power.
Step 4: call the power flow program to calculate the value of objective function (1) and objective function (2), thereby satisfying the constraint conditions. According to the fuzzification method mentioned in Section 2.1, the double objectives are changed into a single objective. The return value is the initial fitness of each particle, which is the current optimal value of the particle. The current global optimal value is obtained by sorting and comparing the current optimal value of each particle. The particle whose objective function value is equal to the global optimal value is the global optimal particle. Set iteration number k = 1.
Step 5: enter the main cycle. According to (18)-(21), calculate the differential inertial weight of the energy storage position variable and the power variable in each iteration; then, calculate the new position and velocity variables of the particles in each iteration according to (10) and (11). The position and velocity variables of the particles are updated according to the global guided cross mutation operation, as shown in (14). Next, call the power flow calculation program, calculate and compare, update the optimal fitness and global optimal fitness of each particle, and record the corresponding positions of the particles. Step 6: the convergence judgment. If the iteration number k reaches the maximum iteration number MaxDT, the iteration ends, and the optimal solution is output; otherwise, k = k + 1, and the loop command continues to be executed in step 4.
The calculation flow chart is shown in Figure 3.
position and velocity variables of the particles in each iteration according to (10) and (11). The position and velocity variables of the particles are updated according to the global guided cross mutation operation, as shown in (14). Next, call the power flow calculation program, calculate and compare, update the optimal fitness and global optimal fitness of each particle, and record the corresponding positions of the particles.
Step 6: the convergence judgment. If the iteration number k reaches the maximum iteration number MaxDT, the iteration ends, and the optimal solution is output; otherwise, k = k + 1, and the loop command continues to be executed in step 4.
The calculation flow chart is shown in Figure 3.

Output results
According to (15), (16), and (17), encode each variable to initialize the particles. The particle position xi and velocity vi are randomly generated, and the particle position vector includes the installation position and power variable of the energy storage system.
N Call the power flow program for calculation. The obtained fuzzy objective function value is then used as the initial fitness of each particle. By sorting and comparison, the current global optimal value and optimal particle are acquired.
According to (18), (19), (20), and (21), the differential inertial weights of the energy storage position variable and power variable are calculated, and then the new particle position and velocity are obtained according to (10) and (11). Then update the position and velocity variables of the particles according to the global guided cross mutation operation in (14).
Call the power flow calculation program, calculate and compare the fitness of each moment, update the optimal fitness value of each particle and the global optimal fitness value, and record the corresponding position of the particles.

k=k+1
The optimal particle is then obtained, (that is, the optimal operational planning).

Start
Read the parameters of the distribution network to be optimized, the locations of the grid-connected PV system, the PV power output and the load demand at each moment.
Initialize the parameters of the IPSO: N, D, MaxDT, w max , w max , c 1 , c 2 , r 1 and r 2.

Results and Discussion
The IEEE-33 bus radial distribution system was examined to test the validity of the optimization ideas proposed in this paper. The distribution network diagram is shown in Figure 4. Bus 1 is the balanced bus, and is also the bus connected to the upper-level power grid. The voltage amplitude of bus 1 is 12.66 kV, the nominal value is 1.0 p.u., and the phase angle is 0. The allowable voltage

Results and Discussion
The IEEE-33 bus radial distribution system was examined to test the validity of the optimization ideas proposed in this paper. The distribution network diagram is shown in Figure 4. Bus 1 is the balanced bus, and is also the bus connected to the upper-level power grid. The voltage amplitude of bus 1 is 12.66 kV, the nominal value is 1.0 p.u., and the phase angle is 0. The allowable voltage deviation range of each bus in the distribution network in this area is ±5% U N , the total active power load is 3715 kW, and the total reactive power load is 2300 kvar. The total active power loss is 202.6878 kW, and the total reactive power loss is 135.6112 kvar without grid-connected PV systems. The typical daily load curve of the distribution network is shown in Figure 5. The data in Figure 5 is from reference [33], as shown in Table A1 of Appendix A, which shows the total active power load of all of the buses (except bus 1) in the IEEE-33 bus system at each time, and the maximum active power load is 3715kW. The line and load data of the IEEE-33 bus system are shown in Table A2 of Appendix A. In the IEEE-33 bus system, the active power load of each bus is considered to be the maximum active power of the corresponding bus. The active power changes of each bus at different times are considered to be consistent with the time series changing trend of a typical daily load curve ( Figure 5). The active power load of each bus at each time in the distribution network is shown in Figure 6. The data of Figure 6 is shown in Table A3 of Appendix A. The typical characteristic curve of a PV system is shown in Figure 7. The data of Figure 7 is from the PV power output characteristics of an actual PV power station at each time, as shown in Table A4 of Appendix A. The maximum PV power output is 3715 kW. The initial conditions of the IPSO algorithm are shown in Table 1. deviation range of each bus in the distribution network in this area is ±5% UN, the total active power load is 3715 kW, and the total reactive power load is 2300 kvar. The total active power loss is 202.6878 kW, and the total reactive power loss is 135.6112 kvar without grid-connected PV systems. The typical daily load curve of the distribution network is shown in Figure 5. The data in Figure 5 is from reference [33], as shown in Table A1 of Appendix A, which shows the total active power load of all of the buses (except bus 1) in the IEEE-33 bus system at each time, and the maximum active power load is 3715kW. The line and load data of the IEEE-33 bus system are shown in Table A2 of Appendix A. In the IEEE-33 bus system, the active power load of each bus is considered to be the maximum active power of the corresponding bus. The active power changes of each bus at different times are considered to be consistent with the time series changing trend of a typical daily load curve ( Figure  5). The active power load of each bus at each time in the distribution network is shown in Figure 6. The data of Figure 6 is shown in Table A3 of Appendix A. The typical characteristic curve of a PV system is shown in Figure 7. The data of Figure 7 is from the PV power output characteristics of an actual PV power station at each time, as shown in Table A4 of Appendix A. The maximum PV power output is 3715 kW. The initial conditions of the IPSO algorithm are shown in Table 1.     deviation range of each bus in the distribution network in this area is ±5% UN, the total active power load is 3715 kW, and the total reactive power load is 2300 kvar. The total active power loss is 202.6878 kW, and the total reactive power loss is 135.6112 kvar without grid-connected PV systems. The typical daily load curve of the distribution network is shown in Figure 5. The data in Figure 5 is from reference [33], as shown in Table A1 of Appendix A, which shows the total active power load of all of the buses (except bus 1) in the IEEE-33 bus system at each time, and the maximum active power load is 3715kW. The line and load data of the IEEE-33 bus system are shown in Table A2 of Appendix A. In the IEEE-33 bus system, the active power load of each bus is considered to be the maximum active power of the corresponding bus. The active power changes of each bus at different times are considered to be consistent with the time series changing trend of a typical daily load curve ( Figure  5). The active power load of each bus at each time in the distribution network is shown in Figure 6. The data of Figure 6 is shown in Table A3 of Appendix A. The typical characteristic curve of a PV system is shown in Figure 7. The data of Figure 7 is from the PV power output characteristics of an actual PV power station at each time, as shown in Table A4 of Appendix A. The maximum PV power output is 3715 kW. The initial conditions of the IPSO algorithm are shown in Table 1.    deviation range of each bus in the distribution network in this area is ±5% UN, the total active power load is 3715 kW, and the total reactive power load is 2300 kvar. The total active power loss is 202.6878 kW, and the total reactive power loss is 135.6112 kvar without grid-connected PV systems. The typical daily load curve of the distribution network is shown in Figure 5. The data in Figure 5 is from reference [33], as shown in Table A1 of Appendix A, which shows the total active power load of all of the buses (except bus 1) in the IEEE-33 bus system at each time, and the maximum active power load is 3715kW. The line and load data of the IEEE-33 bus system are shown in Table A2 of Appendix A. In the IEEE-33 bus system, the active power load of each bus is considered to be the maximum active power of the corresponding bus. The active power changes of each bus at different times are considered to be consistent with the time series changing trend of a typical daily load curve ( Figure  5). The active power load of each bus at each time in the distribution network is shown in Figure 6. The data of Figure 6 is shown in Table A3 of Appendix A. The typical characteristic curve of a PV system is shown in Figure 7. The data of Figure 7 is from the PV power output characteristics of an actual PV power station at each time, as shown in Table A4 of Appendix A. The maximum PV power output is 3715 kW. The initial conditions of the IPSO algorithm are shown in Table 1.

Initial Condition Numerical
Learning factor c1, c2 1.4962 Maximum inertial weight max The optimal configuration and operational strategy of an energy storage system, considering the time-series changes of the 'source and load', with a high penetration single-point and multi-point grid-connected PV system access to the distribution network, are analyzed in this paper. According to the study in [28], in an IEEE-33 bus system, when a single point PV system is connected to the distribution network, the installation location located at bus 7 can minimize the total active power loss. When multi-point PV systems are connected, the grid-connected capacity can be improved to satisfy the voltage quality. Therefore, under a high penetration single-point PV grid connection, the PV system is set to be connected at bus7. In this paper, the percentage of the maximum total PV power output to the maximum total active power load of the distribution network is defined as the penetration rate of grid-connected PV systems. The maximum total active power load in the IEEE-33 bus system is 3715 kW, and the maximum total PV power output in a typical PV daily characteristic curve is also 3715 kW, so 3715 kw is regarded as the 100% PV penetration rate. The power output characteristic curve of a grid-connected PV system at the single point of bus 7 is shown in Figure 7. For high penetration and multi-point PV system grid connections, the total PV penetration is set at 200%. In this case, the total PV power output at each time is twice that of the 100% penetration rate. Then, it is evenly distributed to every bus except for bus 1, so that each bus is connected to the PV system with the same capacity and has the time sequence changes in proportion to the typical daily characteristic curve of the PV system; the PV power output of each bus is shown in Figure 8. The data of Figure 8 is shown in Table A5 of Appendix A. The load of each bus changes according to the time sequence of the load curve in Figure 6. The analysis will be divided into five scenarios:

Initial Condition Numerical
Learning factor c 1 , c 2 1.4962 Maximum inertial weight ω max 1 Minimum inertial weight ω min 0 Maximum number of iterations MaxDT The optimal configuration and operational strategy of an energy storage system, considering the time-series changes of the 'source and load', with a high penetration single-point and multi-point grid-connected PV system access to the distribution network, are analyzed in this paper. According to the study in [28], in an IEEE-33 bus system, when a single point PV system is connected to the distribution network, the installation location located at bus 7 can minimize the total active power loss. When multi-point PV systems are connected, the grid-connected capacity can be improved to satisfy the voltage quality. Therefore, under a high penetration single-point PV grid connection, the PV system is set to be connected at bus7. In this paper, the percentage of the maximum total PV power output to the maximum total active power load of the distribution network is defined as the penetration rate of grid-connected PV systems. The maximum total active power load in the IEEE-33 bus system is 3715 kW, and the maximum total PV power output in a typical PV daily characteristic curve is also 3715 kW, so 3715 kW is regarded as the 100% PV penetration rate. The power output characteristic curve of a grid-connected PV system at the single point of bus 7 is shown in Figure 7. For high penetration and multi-point PV system grid connections, the total PV penetration is set at 200%. In this case, the total PV power output at each time is twice that of the 100% penetration rate. Then, it is evenly distributed to every bus except for bus 1, so that each bus is connected to the PV system with the same capacity and has the time sequence changes in proportion to the typical daily characteristic curve of the PV system; the PV power output of each bus is shown in Figure 8. The data of Figure 8 is shown in Table A5 of Appendix A. The load of each bus changes according to the time sequence of the load curve in Figure 6. The analysis will be divided into five scenarios: Sustainability 2020, 13, x FOR PEER REVIEW 12 of 22 Scenario 1: no PV system is connected; only load timing changes are considered. Scenario 2: the PV system is connected to bus 7 with 100% penetration, and no energy storage system is connected. At the same time, the PV system and load power changes are considered.
Scenario 3: the PV system is connected to bus 7 with 100% penetration, and is connected to the energy storage system, taking into account the time sequence changes of the PV system and load power.
Scenario 4: the multiple-point PV systems are connected to the distribution network with 200% penetration, and no energy storage system is connected. At the same time, the PV system and load power timing changes are considered.
Scenario 5: multiple-point PV systems are connected to the distribution network with 200% penetration, and the energy storage system is connected, taking into account the time sequence changes of the PV system and load power.

Bus Voltage Deviations
Using power flow calculations for different scenarios, the total voltage deviation values of the distribution network in scenarios 1-5 are 54.448 kV, 104.81 kV, 49.278 kV, 252.06 kV and 0 kV, respectively. No matter how great the PV penetration is, the voltage deviation will increase when connected to the PV system and decrease when connected to the energy storage system simultaneously. When the energy storage system is not connected, the higher the PV grid-connected penetration is, the larger the voltage deviation value will be. After the energy storage system is connected, the voltage deviation can be strongly reduced. The bus voltage curves of different scenarios are shown in Figure 9. As seen in Figure 9, in traditional distribution networks without PV system access, under a change of load, the voltage of each bus changes smoothly at each moment. The voltage level at the end of the line is low, and the bus number exceeds the lower limit rate by 37.5%. When the PV system is connected but the energy storage system is not connected, the voltage deviation is obvious and is caused by the randomness and deviation of the PV output. After accessing the energy storage system, the increase in the bus voltage near the PV access point is significantly reduced, and the voltage of each bus is more average and largely within a permissible range. Therefore, with a high PV penetration, the reasonable allocation of energy storage can effectively reduce voltage deviations, increase the voltage level, and ensure that the voltage is near the rated voltage and does not exceed the upper limit. This process can effectively improve the voltage quality. Scenario 1: no PV system is connected; only load timing changes are considered. Scenario 2: the PV system is connected to bus 7 with 100% penetration, and no energy storage system is connected. At the same time, the PV system and load power changes are considered.
Scenario 3: the PV system is connected to bus 7 with 100% penetration, and is connected to the energy storage system, taking into account the time sequence changes of the PV system and load power. Scenario 4: the multiple-point PV systems are connected to the distribution network with 200% penetration, and no energy storage system is connected. At the same time, the PV system and load power timing changes are considered.
Scenario 5: multiple-point PV systems are connected to the distribution network with 200% penetration, and the energy storage system is connected, taking into account the time sequence changes of the PV system and load power.

Bus Voltage Deviations
Using power flow calculations for different scenarios, the total voltage deviation values of the distribution network in scenarios 1-5 are 54.448 kV, 104.81 kV, 49.278 kV, 252.06 kV and 0 kV, respectively. No matter how great the PV penetration is, the voltage deviation will increase when connected to the PV system and decrease when connected to the energy storage system simultaneously. When the energy storage system is not connected, the higher the PV grid-connected penetration is, the larger the voltage deviation value will be. After the energy storage system is connected, the voltage deviation can be strongly reduced. The bus voltage curves of different scenarios are shown in Figure 9. As seen in Figure 9, in traditional distribution networks without PV system access, under a change of load, the voltage of each bus changes smoothly at each moment. The voltage level at the end of the line is low, and the bus number exceeds the lower limit rate by 37.5%. When the PV system is connected but the energy storage system is not connected, the voltage deviation is obvious and is caused by the randomness and deviation of the PV output. After accessing the energy storage system, the increase in the bus voltage near the PV access point is significantly reduced, and the voltage of each bus is more average and largely within a permissible range. Therefore, with a high PV penetration, the reasonable allocation of energy storage can effectively reduce voltage deviations, increase the voltage level, and ensure that the voltage is near the rated voltage and does not exceed the upper limit. This process can effectively improve the voltage quality.

Total Active Power Loss Changes at Different Times in the Distribution Network
The average active power loss values for all times of the day in Scenario 1-5 are 139.81 kW, 110.65 kW, 94.35 kW, 131.03 kW and 73.66 kW, respectively. When the PV and energy storage systems are not connected, the active power loss is the largest. The average active power loss will be reduced after the PV system is connected, and the average active power loss will be further reduced after the energy storage system is connected. When not connected to the energy storage system, the higher the PV penetration rate is, the greater the average active power loss will be. After the energy storage system is connected to the distribution network, the multi-point PV grid connection with a high penetration rate will reduce the active power loss more than the single point PV grid connection. The change curve of the total active power loss of the distribution network at different times under different scenarios is shown in Figure 10. It can be seen that the active power loss at each time of the traditional distribution network is large, and the change trend is similar to the load change trend of the distribution network. Under a PV system with a 100% penetration rate and single point access in a traditional distribution network, the total active power loss of the distribution network can be significantly reduced, and the total active power loss of the distribution network can be further reduced after the storage system is configured, where the maximum decrease rate of the active power loss is 47.05%. Using a PV system with a 200% penetration rate in the traditional distribution network, the total active power loss fluctuates greatly at different times in the distribution network. Sometimes the active power loss increases, and the maximum increase rate is 53.19%. In some cases, the active power loss decreases, with the maximum reduction rate reaching 65.22% due to the simultaneous consideration of the time series changes of load and PV power and the high penetration rate of the PV system in the distribution network. The total active power loss of the distribution network at each time after the allocation of the energy storage systems is significantly reduced, with a maximum reduction rate of 66.09%. In most periods, the power loss is lower than the power loss of a single-point PV grid connection after the energy storage is configured.

Total Active Power Loss Changes at Different Times in the Distribution Network
The average active power loss values for all times of the day in Scenario 1-5 are 139.81 kW, 110.65 kW, 94.35 kW, 131.03 kW and 73.66 kW, respectively. When the PV and energy storage systems are not connected, the active power loss is the largest. The average active power loss will be reduced after the PV system is connected, and the average active power loss will be further reduced after the energy storage system is connected. When not connected to the energy storage system, the higher the PV penetration rate is, the greater the average active power loss will be. After the energy storage system is connected to the distribution network, the multi-point PV grid connection with a high penetration rate will reduce the active power loss more than the single point PV grid connection. The change curve of the total active power loss of the distribution network at different times under different scenarios is shown in Figure 10. It can be seen that the active power loss at each time of the traditional distribution network is large, and the change trend is similar to the load change trend of the distribution network. Under a PV system with a 100% penetration rate and single point access in a traditional distribution network, the total active power loss of the distribution network can be significantly reduced, and the total active power loss of the distribution network can be further reduced after the storage system is configured, where the maximum decrease rate of the active power loss is 47.05%. Using a PV system with a 200% penetration rate in the traditional distribution network, the total active power loss fluctuates greatly at different times in the distribution network. Sometimes the active power loss increases, and the maximum increase rate is 53.19%. In some cases, the active power loss decreases, with the maximum reduction rate reaching 65.22% due to the simultaneous consideration of the time series changes of load and PV power and the high penetration rate of the PV system in the distribution network. The total active power loss of the distribution network at each time after the allocation of the energy storage systems is significantly reduced, with a maximum reduction rate of 66.09%. In most periods, the power loss is lower than the power loss of a singlepoint PV grid connection after the energy storage is configured.

The Active Power Changes of the Upper-Level Grid Injected into the Distribution Network at Different Times
The active power changes of the upper-level grid injected into the distribution network under different scenarios are shown in Figure 11. As shown in Figure 11, the power curve of the traditional distribution network is consistent with the typical daily load curve, i.e., the load power and power loss are supplied by the upper-level grid. After accessing the PV system, the injected power of the upper-level grid curve is the opposite of the power curve of the PV power output. At noon, when the PV power output is at its maximum, the distribution network no longer needs the power supply of the upper-level grid, and an excessive PV penetration rate will reverse the power to the upper-level grid. After accessing the energy storage system, the power injection of the upper-level grid will be less than that of the traditional distribution network and smoother than that without the energy storage system. The configured energy storage system can effectively reduce and smooth the transmission power of the upper-level grid.
Sustainability 2020, 13, x FOR PEER REVIEW 15 of 22 The active power changes of the upper-level grid injected into the distribution network under different scenarios are shown in Figure 11. As shown in Figure 11, the power curve of the traditional distribution network is consistent with the typical daily load curve, i.e., the load power and power loss are supplied by the upper-level grid. After accessing the PV system, the injected power of the upper-level grid curve is the opposite of the power curve of the PV power output. At noon, when the PV power output is at its maximum, the distribution network no longer needs the power supply of the upper-level grid, and an excessive PV penetration rate will reverse the power to the upper-level grid. After accessing the energy storage system, the power injection of the upper-level grid will be less than that of the traditional distribution network and smoother than that without the energy storage system. The configured energy storage system can effectively reduce and smooth the transmission power of the upper-level grid.

Energy Storage Optimized Allocation and Operating Results
The location variable of the energy storage can be obtained by the optimization result of the IPSO algorithm. The optimal access position of the energy storage system in scenario 3 is bus 7, and the optimal access position of the energy storage system in scenario 5 is bus 11. Figure 12 illustrates the charging and discharging of the power of the energy storage systems at different times in scenarios 3 and 5. According to the optimization results, the optimal capacity of the energy storage system can be obtained by the power and time period of the charge/discharge of the energy storage system. Combined with the charge/discharge power and the duration at each time, the charge or discharge power energy of the energy storage system in each time period can be obtained; then, the cumulative charge/discharge power energy of the continuous charge/discharge period of the energy storage system can be calculated, taking the maximum value of the absolute value as the optimal capacity of the energy storage system. According to the calculation results in scenario 3, the optimal configurable capacity of the energy storage system is 0.81 MWh. According to the calculation results of Scenario 5, the optimal configuration capacity of the energy storage system is 22 MWh, and the energy storage capacity that needs to be configured is larger due to the high penetration of PV systems and the low transmission power of the upper-level grid. The energy storage SOC curves for scenarios 3 and 5 are shown in Figure 13.

Energy Storage Optimized Allocation and Operating Results
The location variable of the energy storage can be obtained by the optimization result of the IPSO algorithm. The optimal access position of the energy storage system in scenario 3 is bus 7, and the optimal access position of the energy storage system in scenario 5 is bus 11. Figure 12 illustrates the charging and discharging of the power of the energy storage systems at different times in scenarios 3 and 5. According to the optimization results, the optimal capacity of the energy storage system can be obtained by the power and time period of the charge/discharge of the energy storage system. Combined with the charge/discharge power and the duration at each time, the charge or discharge power energy of the energy storage system in each time period can be obtained; then, the cumulative charge/discharge power energy of the continuous charge/discharge period of the energy storage system can be calculated, taking the maximum value of the absolute value as the optimal capacity of the energy storage system. According to the calculation results in scenario 3, the optimal configurable capacity of the energy storage system is 0.81 MWh. According to the calculation results of Scenario 5, the optimal configuration capacity of the energy storage system is 22 MWh, and the energy storage capacity that needs to be configured is larger due to the high penetration of PV systems and the low transmission power of the upper-level grid. The energy storage SOC curves for scenarios 3 and 5 are shown in Figure 13. Sustainability 2020, 13, x FOR PEER REVIEW 16 of 22

Conclusions
In this paper, an optimal approach for the simultaneous location, capacity and charge/discharge of an energy storage system under high penetration single point and multi-point grid-connected PV conditions was proposed. The two objectives were to minimize the active power loss and voltage deviations, which were combined into a single objective using the membership function and weighting method, and the improved particle swarm algorithm was used to solve the problem. In the improved particle swarm optimization algorithm, the differential inertial weight was used to balance the exploration ability and developmental ability of the algorithm, and a global guided cross search mechanism was used to improve the convergence speed. For the IEEE-33 bus system, taking into account the timing changes of the load and PV power output, the optimal allocation and operation results of the energy storage under high penetration single point and multi-point gridconnected PV conditions were achieved, which were compared with the optimization results of different scenarios. The results show that the proposed method can determine the optimal configuration and operation strategy for an energy storage system with high penetration gridconnected PV systems, thereby improving the voltage level, reducing the power loss, reducing and

Conclusions
In this paper, an optimal approach for the simultaneous location, capacity and charge/discharge of an energy storage system under high penetration single point and multi-point grid-connected PV conditions was proposed. The two objectives were to minimize the active power loss and voltage deviations, which were combined into a single objective using the membership function and weighting method, and the improved particle swarm algorithm was used to solve the problem. In the improved particle swarm optimization algorithm, the differential inertial weight was used to balance the exploration ability and developmental ability of the algorithm, and a global guided cross search mechanism was used to improve the convergence speed. For the IEEE-33 bus system, taking into account the timing changes of the load and PV power output, the optimal allocation and operation results of the energy storage under high penetration single point and multi-point gridconnected PV conditions were achieved, which were compared with the optimization results of different scenarios. The results show that the proposed method can determine the optimal configuration and operation strategy for an energy storage system with high penetration gridconnected PV systems, thereby improving the voltage level, reducing the power loss, reducing and

Conclusions
In this paper, an optimal approach for the simultaneous location, capacity and charge/discharge of an energy storage system under high penetration single point and multi-point grid-connected PV conditions was proposed. The two objectives were to minimize the active power loss and voltage deviations, which were combined into a single objective using the membership function and weighting method, and the improved particle swarm algorithm was used to solve the problem. In the improved particle swarm optimization algorithm, the differential inertial weight was used to balance the exploration ability and developmental ability of the algorithm, and a global guided cross search mechanism was used to improve the convergence speed. For the IEEE-33 bus system, taking into account the timing changes of the load and PV power output, the optimal allocation and operation results of the energy storage under high penetration single point and multi-point grid-connected PV conditions were achieved, which were compared with the optimization results of different scenarios. The results show that the proposed method can determine the optimal configuration and operation strategy for an energy storage system with high penetration grid-connected PV systems, thereby improving the voltage level, reducing the power loss, reducing and smoothing the transmission power of the upper-level grid, and improving the storage capacity of grid-connected PV systems.
In the process of the optimal allocation and operation of the energy storage system, the optimization objectives were only to improve the voltage quality and reduce the network loss from the technical point of view; the economic issues, such as equipment investment cost, operation cost and profit, will be considered in the future by authors.     Table A3. The data of Figure 6 (kW).