Optimal Operation of Energy Storage Units in Distributed System Using Social Spider Optimization Algorithm

Optimal placement and performance of energy storage units are two important issues for power systems. Energy storage units help in energy supply and elevated reliability in the network. Energy storage system (ESS) cause many of the variables of planning, designing, and exploiting the network to alter, and thus planning and operational methods for the network will change as well. These units also assist in reducing the peak load with the aim of profit maximization. In this paper, a method is presented to optimize the performance of energy storage units in the power market. To solve the problem, a new optimization algorithm called social spider optimization (SSO) algorithm is employed. Cost reduction, losses reduction and voltage profile improvement are three main scenarios which are considered. The proposed method is applied to a 33-bus standard distribution system and results show that this novel method is efficient to use for ESSs optimal operation.

Many researchers have tried to optimize the performance of the ESSs for different applications. The important thing about the ESS is its planning to meet the load demands. In [6] a novel approach is introduced for optimal operation of distribution networks at the presence of DG and ESS. Reliability and costs are objective of this study and a hybrid evolutionary algorithm called Hybrid Grey Wolf Optimizer-Particle Swarm Optimization (HGWO-PSO) algorithm, is used to optimize. An analysis of the operation mode of energy storage is presented in [7]. The analysis is performed through a mixed integer linear programming (MILP) model programmed via the General Algebraic Modelling System (GAMS) to reduce operation costs.
A multi-objective methodology for locating, sizing and operation of ESSs in distribution network for a day is presented in [8]. The objective functions of the problem are minimization of power losses and the ESS cost that uses non dominated sorting genetic algorithm II and tabu search optimization algorithm. In [9] a method to find the optimal ESS capacity for a peakshaving application is proposed. This method is based on Dynamic programming to improve voltage profile and energy cost. Optimal storage operation problem as a finite-horizon dynamic program is formulated in [10] for a grid with renewable DGs and AC/DC load to minimize expected energy cost. Ref. [11] presented mixed-integer second-order cone programing and mixed-integer linear programing models for solving the problem of optimal operation of distribution networks with ESSs. Optimal operation of the grid with DG and ESS for multiple objective functions such as minimization of losses, imbalanced power at the substation, total energy costs, and peak demand is presented in [12].
Study on optimal operation of hybrid storage systems include battery and super-capacitor is presented in [13,14] for distribution system. Optimal operation and energy management of a hybrid power system in smart grid framework is presented in [15]. In this research wind turbine, photovoltaic and ESS are integrated and model predictive control (MPC) is used for optimal management.
A distribution network with ESS is optimized with epsilon-constrained differential evolution, which is categorized as a metaheuristic optimization method, in [16] to minimize operating costs under nonlinear conditions and various energy system connections. A novel stochastic planning framework is proposed in [17] to determine the optimal ESS capacity in an isolated micro-grid. In this paper optimal operation of the ESSs in a distribution system is presented. In this study SSO algorithm is used as a novel method to minimize cost, losses and improve voltage profile. The rest of this paper is organized as follows: problem formulation includes decision variables, objective functions and constraints are illustrated in section II. Section III is dedicated to backward/forward power flow. The SSO algorithm is explained in section IV. Simulation results and conclusion are presented in section V and VI, respectively.

A. Decision Variables
In this study, contribution of storage units is considered in electricity markets to achieve maximum financial benefit within a specified period. Determining the optimal time of storage unit charge and discharge is the main decision variable. This time should be determined in such a way that the problem objective function is optimized.

B. Objective Functions
The ESSs are organized and contributed in the market to achieve more profit. The storage units' profit can be expressed as follows: (1) Where, Pbuy t is the power purchased from the grid, Psale t represents the power sold to the grid, and λgrid t shows power cost at t.

C. Constraints C.1. Buses voltage limits
The permitted range of the grid voltage level changes is very limited, and standards usually allow for small variations around the nominal value. According to the limitation defined for voltage of buses, the following equation can be expressed: Where, vi is the bus voltage, vi min and vi max are minimum and maximum voltages of i th bus, respectively and Nbus is the total bus number.

C.2. Permitted line current flow
The thermal capacity of the lines is simply determined based on current flow; thus, the following equation can be used for the current flow limitation: Where, Ii is the current of i th line, Ii max is the maximum allowed current for i th line, and NLine is the total number of lines.

C.3. Constraint on ESS
State of charge (SOC) is defined as the ratio of battery remaining energy to total energy saved. Each EES has its specified SOC which must always be within the allowed range, as:

SOC
SOC SOC  min max (4) Where, SOCmin and SOCmax are minimum and maximum values of the storage SOC, respectively.

C.4. Load modeling
With respect to the load static features, three general models are proposed in power system studies, including constant power (powers are independent of voltage), constant current (powers are proportional to bus voltage) and constant impedance (powers are proportional to square bus voltage). In this paper, constant power model is considered.

Power Flow
Traditional power flow methods in transmission networks are less effective for distribution systems due to radial structure and high R/X ratio. Also, by increasing the DG unit's penetration, the distribution network changes from a passive to an active grid [18]. Therefore, instead of common power flow methods, backward/forward sweep method based on current summation is used.

A. Backward/forward Sweep Method
The backward/forward sweep method is appropriate for radial construction, whose effectiveness has already been proven. Simple structure, fast convergence and suitability for online and offline problems are the main advantages of this method [19]. This method has three step which are explained below:

A.1. Backward sweep
In the first iteration, all bus voltages are equal to source voltage (1 p.u.). In the next iteration, bus voltages are calculated from the previous step. With respect to bus voltages, the load current can be calculated as follows: Where, ILi, Pi, Qi and Vi are the current, active power, reactive power and voltage of i th load, respectively. Thereafter, it is necessary to calculate the current through lines starting from the farthest line relative to the reference line. For example, line j can be expressed as: Optimal Operation of Energy Storage Units in Distributed System Using = ∑ ∈ i=1… Nbus (6) Where, Nbus is the bus number, ILi is the current of j th line, and D represent all the lines connected to bus i. Consequently, backward sweep is finished and the current of lines is updated.

A.2. Forward sweep
In the forward sweep starting with source bus, whose voltage is known, and considering impedance and current of lines, the voltage of i th bus is calculated as follows: Where, Vi is the voltage of i th bus, Vi-1 is the first bus voltage at i th line, ILi shows the current of i th line and Nbus is the bus number. Upon completion of this process, the forward sweep is finished and all bus voltages are updated.

A.3. Convergence criterion survey
After backward/forward sweep, the convergence criterion needs to be calculated to assess the need for further iterations.
Vi old is the voltage of i th node in the previous iteration, Vi new denotes the voltage of i th node in the last iteration and ε is the permissible voltage deviation value. If the above equation persists, the power flow will end; otherwise, iterations will continue.

B. Power flow in distributed network with DG sources
Distributed systems with DGs are similar to multi source networks, so buses with DG are modelled as PQ or PV bus [18]. In this paper, DG units are controlled like a PQ bus and hence, they are considered as a negative load. The PQ buses are modeled as follows:

Social Spider Optimization (SSO) Algorithm
Optimization is a challenging area that has attracted increasing attention in recent decades. Optimization finds the best solution depending on the problem, limits and constraints [20]. Recently, metaheuristic algorithms have gained more popularity as they are superior to traditional gradient-based approaches.
Although SSO may also be considered a metaheuristic algorithm, it has fundamental differences with them. SSO is different in terms of agents' trainability, the manner of information dissemination, and the exploration frameworks used in the algorithm. SSO indicates lower convergence potential in some problems, but it is also able to solve multi-objective optimization problems with numerous local optimization points [21]. In many algorithms based on collective intelligence including artificial bee colony (ABC) and group search optimization (GSO), populations are categorized into different types. Different independent agents perform different tasks and the entire population cooperate with each other to search and explore the search space. However, in SSO, all of the independent agents (spiders) are the same and equal. Each agent performs all tasks which should be carried out by several types of species of population in other algorithms [21]. If we want to integrate SSO into the traditional and conventional framework, a new characteristic will be added to the algorithm, whereby different factors can be converted to other factors simply and without user control, which can relatively improve the algorithm efficiency.

A. Spider
Spiders are optimization agents in SSO. At the beginning of running the algorithm, a number of predetermined spiders lie on the web. Each spider S has a memory, and stores independent information including the location of the spider s on the web, the extent of qualification of the current position of the spider S, vibration of the spider S target in the previous repetition, the Mehdi Tabasi, et al.
number of repetitions when spider S altered its target's vibration for the last time, displacement of spider S in the previous repetition, dimensional cover used by spider S for running its displacement in the previous repetition [21].
The position of the spider as well as the degree of eligibility is used for stating the special status of spider S, and other information is employed for navigation of spider S towards new locations.
Based on observations, spiders have a very accurate sense towards vibrations. In addition, they are able to distinguish between different vibrations propagated on the same web, and perceive the intensity associated with each of them. In this way, the spiders lying on a web share their personal information with others to develop a public social awareness.

B. Vibration
In SSO, two characteristics including the location of the vibration source and its intensity are used to define a vibration. The location of the source is defined based on the search space of the optimization problem, while the intensity of each vibration is defined within the range of [0, +∞). By transferring each spider to a new position, a vibration is developed in its current position. The position of spider A at the time of t is defined as PA(t), or more simply as PA. Further, G(PA, PB, t) is used to represent the intensity of vibration sensed by a spider at the time of t and position of PB, whose source is located at PA. Therefore, G(PS, PS, t) means intensity of vibration produced by spider S at the source position [21]. The intensity of vibration at the source position is dependent on the eligibility of its position F(PS): Where, C is a small constant value. According to (10), it is concluded that all possible values for the intensity of vibration and optimization problem are positive and positions with greater eligibility have greater intensity compared to positions with unsuitable eligibility values. As a type of energy, vibration is damped in a distance. This phenomenon is taken into account in SSO design. The distance between spiders A and B is defined as D(PA, PB), which is obtained by the following Relation: The Where, (0, ) a r  is adjusted and controlled by the user and controls the damping rate of the vibration intensity along the distance. For a weaker damping, ra should be considered large.

C. Search Pattern
In SSO, there are three consecutive stages of initialization, iteration, and finalization. At the stage of initialization, the objective function and search space as well as the numerical value used for the algorithm are defined. After adjusting these values, the algorithm starts generating an initial population of spiders. As the number of spiders along the SSO simulation remains unchanged, a memory with a constant size is considered for storing their information. The positions of spiders are generated randomly and then eligibility values of each of them are calculated and stored. The first target vibration related to each spider present in the population is adjusted and determined at the current position of that spider, and the intensity of the vibration is considered zero. All of the initial data stored for spiders are also valued as zero. In this way, the stage of initialization is terminated, and the algorithm begins the iteration stage, and deals with searching for the problem solution by the generated spiders. At this stage, a specific number Global optimal value. The eligibility values for each spider along every iteration are assessed. Once all vibrations were generated, the algorithm simulates propagation of these vibrations. In this process, every spider S will receive the number |population| different vibrations generated by other spiders. The information received from these vibrations includes the source of vibration and the extent of its damped intensity. V is used to represent the number |population| of vibration [21].
As soon as the vibration's V is received, the spider S selects the most powerful vibration vs best among the vibration V and compares its intensity with the intensity of the target vibration vs tar stored in its memory. If the intensity of vs best is higher, spider S stores vs best as vs tar and Cs (number of iterations since s has last changed its target vibration) is reset to zero. Otherwise, the main vtar is kept, and Cs increase one unit. PS i and PS tar represent the positions of the source for generation of V and vtar, in which i={1, 2, …, |population|}.
Thereafter, the algorithm takes some measure causing spider S to move randomly towards vs tar . Every spider has a m dimensional cover, which is a binary vector consisting of 0-1, with a length of D, with D showing the number of dimensions of the optimization problem. First, all cover values are zero. In each iteration, spider S alters its cover with the probability of 1-; Here ∈ (0˒ 1) and is defined by the user, which represents probability of cover variation. If it is decided to change the cover, each bit of this vector is replaced by a probability of pm equal to 1, and remains with a probability of 1-pm unchanged and equal to zero. The pm is also defined by the user within the range of (0, 1). Every bit of the cover alters independently, with absolutely no relationship with the previous cover. If all bits are zero, a random value of the cover changes into value of one. Similarly, if all bits are equal to one, a random bit will be considered as zero. Once the dimensional cover was determined, the following new position called it generated based on the cover of spider S. The value of ith dimension of the following position, ˒ , is generated as follows: Where, r is a random number within the range of [1, |population|], which is produced separately, while ms,I represents the ith dimension of the dimension m cover associated with the spider S.
After generating PS fo , spider S performs a random movement towards this position. This random movement is conducted and controlled by the following relation: Where, r is the element-wise multiplication and R is a vector of decimal random numbers between zero and one. Before movement towards PS fo , spider S moves towards its previous movement direction (the direction related to the previous iteration). The length of traversal in this direction is equal to a random percentage of the previous movement. Hence, spider S moves along any dimension and with random coefficients within (0, 1) range towards PS fo . These random coefficients are generated independently and individually for each dimension. After running this random movement, spider S stores its manner of displacement in the current iteration to be used in the next iteration. In this way, the sub stage of random movement is finished. The final sub stage of the stage of iterations is controlling establishment of constraints. During the stage of random movement, spiders may leave the web's zone, causing transgression of the optimization problem constraints. Reflective method can be used to develop the constraints, and one can generate a free position from bounded constraints, i.e. PS (t+1) using the following relation: Where, ̅ is the upper bound of the search space in the ith dimension, while represents the low bound of its corresponding dimension. The r is a random decimal number within the range of (0,1). The stage of iteration continues iteratively until the final criterion is fulfilled. The final Mehdi Tabasi, et al.
criterion can be defined as achieving the maximum number of iterations, achieving a specific error, or any other suitable criterion. By the completion of the iteration stage, the best solution found with the best eligibility is displayed as the output of the algorithm. The above three stages constitute a complete the SSO algorithm.  By introducing problem and the SSO algorithm, proposed optimization method is simulated by MATLAB software under different scenarios.

Results and Discussion
The proposed method is applied to a standard 33-bus radial distribution network with 4 laterals and 33 nodes, as shown in Figure 1. The main data of this network are presented in Appendix. This test network has 32 individual branches and its total active and reactive power  Table 1 presents energy price and hourly load during 24 hours.

A. Scenario I
Maximizing the operational profit from storage units during 24 hours is the main purpose of the first scenario. The SSO algorithm population and maximum iteration are considered equal to 100 and 200, respectively.
The convergence curve of the first scenario is illustrated in Figure 2. The SSO algorithm has the ability to converge rapidly to the optimal global answer. The ESSs cost is calculated as about -1016 $ which means ESSs have earned 1016$ during a day through energy sales and purchases.       Optimal Operation of Energy Storage Units in Distributed System Using Figure 6 shows the charge and discharge rate of ESSs during 24 hours. The blue sections (positive) represent charge state and red sections (negative) indicate discharge state.
When the energy price is cheap or the load is not heavy, ESSs store energy and vice versa, leading to decreased ESSs operational costs and more profits. For example at hour 24, the price of energy is too high, so three ESSs are produced energy (red section).
According to scenario I, ESSs satisfies SOC limits and they produce in expensive time. This figures demonstrate behavior of ESSs for 24 hours so as to obtain high profit.

B. Scenario II
Reducing power losses is the main purpose of the second scenario. The SSO algorithm population and maximum iteration are similar to the previous scenario. This scenario is investigated in two modes, which are with and without storage units. Figure 7 shows the convergence curve of the second scenario. According to this figure, the SSO has been converged within a short time.
Total losses for 24 hours a day are 3109 kW, which is equal to 3101 kW before presence of storage units. After optimal charging and discharging, losses grow by about 0.25% which can be ignored. This scenario indicates that ESSs do not have a significant effect on power losses. It should be noted that situation depends on the number of storage units, their capacity as well as charge and discharge rates. Similar to scenario I, the SOC of installed storage units at buses 10, 18 and 33 is shown in Figure 8, Figure 9 and Figure 10, respectively. These figures confirm that the SOC satisfies constraints, with its change rates being equal to 8 kW.   The ESSs charging and discharging rate for one day is shown in Figure 11.  Optimal Operation of Energy Storage Units in Distributed System Using According to proposed information in scenario II, the ESSs cannot reduce significantly power losses. Furthermore, since the aim of this scenario is loss reduction, results may different with scenario I. similar to previous scenario, all constraints are satisfied.

C. Scenario III
In this scenario, the problem is solved with the aim of improving the voltage profile. The SSO algorithm population and maximum iteration are considered equal to 100 and 200, respectively. This scenario is investigated in two modes; with and without storage units. Figure 12 displays convergence curve of the third scenario. The algorithm is converged after 42 iterations.   The SOC of connected ESSs at buses 10, 18 and 33 is shown in Figure 14, Figure 15 and Figure 16, respectively. Figure 17 shows the charge (blue) and discharge (red) rate of ESSs during 24-hours.     Improvement of voltage profile with the presence of ESSs is marginally. Moreover, the cost of energy and profit is not the matter of third scenario, so it can be seen ESSs operates in accordance with voltage level as well as constraints. For example, the installed ESSs at bus 10 and 33 don't operate at 24th hour but the ESS at bus 18 produce energy at same time. The behavior of ESSs in this hour is opposite with scenario I.

Conclusion
In this paper, contribution of storage units was examined in the form of an optimization problem with the aim of acquiring the maximum profit for a 24-hour period. The technical limitations including grid constraints, charging and discharging rate of the storage sources, and SOC of these forces were also taken into account in definition of the problem. Furthermore, to solve the problem, a new smart optimization algorithm called social spider optimization (SSO) algorithm was used.
The results obtained from SSO algorithm indicated that it is efficient and can be converged to the global optimal solution with a high speed. Optimal charging and discharging of storage units cause profitability. This is because when electricity price is cheaper, units usually start charging and vice versa. This causes reduction of operational costs of units, bringing about profitability. Although presence of storage units causes increased grid loss, this increase is usually trivial, and is dependent on the number, placement of units as well as their charging and discharging rates. In addition, the effect of storage units on improving voltage profile is trivial, though through optimal placement of ESS, one can improve voltage profile and reduce losses.