A Scatter Search Heuristic for the Optimal Location , Sizing and Contract Pricing of Distributed Generation in Electric Distribution Systems

In this paper we present a scatter search (SS) heuristic for the optimal location, sizing and contract pricing of distributed generation (DG) in electric distribution systems. The proposed optimization approach considers the interaction of two agents: (i) the potential investor and owner of the DG, and (ii) the Distribution Company (DisCo) in charge of the operation of the network. The DG owner seeks to maximize his profits from selling energy to the DisCo, while the DisCo aims at minimizing the cost of serving the network demand, while meeting network constraints. To serve the expected demand the DisCo is able to purchase energy, through long-term bilateral contracts, from the wholesale electricity market and from the DG units within the network. The interaction of both agents leads to a bilevel programming problem that we solve through a SS heuristic. Computational experiments show that SS outperforms a genetic algorithm hybridized with local search both in terms of solution quality and computational time.


Introduction
Distributed generation (DG), known as the small-scale production of electricity near consumers, has become increasingly common in distribution networks in the last decade.This trend has been driven by different factors.Such as environmental concerns, fossil fuel resource depletion, and new advances in small-scale generation technologies [1].Distribution systems were not originally designed to operate with DG; therefore, its insertion in the network must be carefully planned out.In this regard, the integration of DG in distribution networks has been the focus of several studies [2][3][4][5][6][7][8][9][10][11][12].Most of these studies aim at harvesting the potential benefits of DG such as: reduction of power losses [2][3][4], improvement of voltage profile [5], improvement of stability [6], and deferral of investments [7].In recent years there has been an important increase in the participation of DG in distribution systems.This has been motivated by several factors that include new advances in microgeneration and storage technologies, new information and communications technologies, a rapid growth of demand, and the process of liberalization of energy markets.Despite of these facts, there are still technical, economic, and regulatory barriers for the widespread participation of DG.A review of factors that can contribute to the definite boost of DG is presented in [13], and a survey regarding critical issues in smart grid technologies and their integration is presented in [14].Also, several studies have been conducted based on sensitivity analysis to face the challenge of improving voltage profile and active/reactive power flows for low voltage distribution systems with DG.Examples of such studies can be consulted in [15][16][17].
Potential benefits of DG depend on several factors such as its size, location, and network parameters.The studies conducted for the optimal location and/or sizing of DG can be broadly classified into three groups: (i) analytical expressions; (ii) exact algorithms; and (iii) metaheuristic algorithms.In [3,4], the authors present analytical methods for the optimal location of DG in distribution systems aiming to reduce power losses.The main disadvantage of such methods lies in the fact that they often resort to oversimplifications of the network.Exact methods for the optimal location and sizing of DG are presented in [18,19].These methods are based on classical mathematical programming and within certain conditions of convexity and differentiability can find optimal solutions for the location and/or sizing of DG.As with analytical expressions, for the sake of their success, mathematical-programming based methods must resort to some simplification in the modeling of the network.These simplifications are often applied to the power balance and power flow equations.If no simplifications of the network are considered, the optimal sizing and location of DG becomes a non-linear non-convex optimization problem.These problems are better handled with metaheuristic techniques rather than with classic mathematical approaches.Some metaheuristic techniques applied to the optimal location and sizing of DG include the use of particle swarm optimization [20], genetic algorithms [21], bee colony optimization [22], and differential evolution [23].Hybrid metaheuristics have also been proposed in [24,25].A detailed description of the aforementioned methods is out of the scope of this paper; the interested reader is referred to [26,27] for literature reviews on the optimal placement and sizing of DG.Particularly, in [27] the authors review the applications of different heuristics and metaheuristics for the optimal placement and sizing of DG.Despite the common usage of population heuristics (genetic algorithms, particle swarm optimization, differential evolution, etc.), to the best of our knowledge scatter search (SS) has not been applied to this problem before.
Consequently, the main contribution of this paper is the proposal of a SS heuristic for the optimal location, sizing, and contract pricing of DG in power systems.The rest of the paper is organized as follows.Section 2 presents the mathematical model used for the optimal location, sizing, and contract price of DG.Then, Section 3 summarizes the main elements of the proposed heuristic.Section 4 presents the computational results that compare SS against a memetic algorithm.Finally, Section 5 provides the conclusions and future research possibilities.

Mathematical Model
Most of the approaches regarding optimal location and/or sizing of DG tackle this decision from the standpoint of the Distribution Company (DisCo).By contrast, in this paper we consider the interaction between the DisCo and the DG owner.Consequently, we take into account a market structure in which the DisCo is free to purchase energy either from the wholesale energy market, from DG units within its network, or both.In this case, the DisCo aims at meeting the network demand at the lowest possible cost while keeping network constraints.On the other hand, the DG owner aims at maximizing his profits from the energy sales to the DisCo.This interaction leads to a Stackelberg leader-follower game that we model as a bilevel programming problem [19].
The methodology proposed in this paper is suitable for an energy market in which the DisCo and DG investor are willing to engage in business.Also, we assume that power losses (or a percentage of them) are paid by the DisCo.This is the case of most deregulated electricity markets in which the regulator promotes mechanisms to force DisCos to reduce power losses.Furthermore, the action-reaction nature of the proposed approach implies dispatchability of DG resources.That is, the DisCo, based on an optimal dispatch calculation, might be willing to purchase certain amount of power from the DG units, at a given period, for which the DG must be available.

Decision-Making Problem of the Distribution Company
To meet the expected demand, the DisCo must decide on how much energy to purchase from the wholesale energy market and how much from the DG units within its network.At first glance, the easiest solution would be to buy energy from the cheapest source.However, the energy price is not the only factor that the DisCo must take into account when purchasing energy, since a power injection in an unsuitable bus might lead to problems such as overvoltage or an increase in power losses.Therefore, in order to account for network constraints, the DisCo must solve an optimal dispatch considering the presence of DG.For a given location and contract price, the solution of the optimal dispatch would yield the right amount of power to be purchased from the DG units and the wholesale energy market in order to minimize costs while meeting network constraints.

Decision-Making Problem of the Distributed Generation Owner
The objective function of the DG owner, considered in the proposed model, is the maximization of profits from the energy sales to the DisCo.In this case the DG owner must decide over the location, size, and contract price of the DG units.Also, he must take into account the reaction of the DisCo.Regarding contract price, negotiating a low energy price might result in more energy sales but not necessarily in more profits.As regards location, the owner of the DG must be aware of the fact that the impact of DG is highly dependent on the bus where it is located.Regardless of the energy price, the DisCo might not be willing to buy energy from the DG units if it has a negative impact on the network.Finally, the sizing of the DG is important, since it determines the investment costs and has a direct impact on the profits of the DG owner.

Bilevel Modeling Framework
The decision-making problems of both the DisCo and DG owner can be combined into a bilevel programming framework as shown in Figure 1.In this case, the upper level optimization problem corresponds to the DG owner; he must decide over the size, location, and contract price of the DG units (action).Once this information is known by the DisCo (located in the lower level optimization problem), he uses it as parameters for an optimal dispatch and determines the amount of energy he would be willing to purchase from the DG units (reaction).The DG owner now knows how much energy would be bought from his units, and is able to compute his expected profits.With this information, he might want to change his strategy testing new locations, contract prices, and sizes of DG.The process is repeated until an equilibrium is reached.That is, until there are not incentives for any of the agents to change their decisions.The idea behind this leader-follower game lies in the fact that the leader must consider the reaction of the follower before setting his strategy.That is, the DG owner takes into account the impact of the DG units in the network before deciding on their parameters.Equations ( 1)-( 20) describe the mathematical model of the proposed methodology.

Decision-Making Problem of the Distribution Company
To meet the expected demand, the DisCo must decide on how much energy to purchase from the wholesale energy market and how much from the DG units within its network.At first glance, the easiest solution would be to buy energy from the cheapest source.However, the energy price is not the only factor that the DisCo must take into account when purchasing energy, since a power injection in an unsuitable bus might lead to problems such as overvoltage or an increase in power losses.Therefore, in order to account for network constraints, the DisCo must solve an optimal dispatch considering the presence of DG.For a given location and contract price, the solution of the optimal dispatch would yield the right amount of power to be purchased from the DG units and the wholesale energy market in order to minimize costs while meeting network constraints.

Decision-Making Problem of the Distributed Generation Owner
The objective function of the DG owner, considered in the proposed model, is the maximization of profits from the energy sales to the DisCo.In this case the DG owner must decide over the location, size, and contract price of the DG units.Also, he must take into account the reaction of the DisCo.Regarding contract price, negotiating a low energy price might result in more energy sales but not necessarily in more profits.As regards location, the owner of the DG must be aware of the fact that the impact of DG is highly dependent on the bus where it is located.Regardless of the energy price, the DisCo might not be willing to buy energy from the DG units if it has a negative impact on the network.Finally, the sizing of the DG is important, since it determines the investment costs and has a direct impact on the profits of the DG owner.

Bilevel Modeling Framework
The decision-making problems of both the DisCo and DG owner can be combined into a bilevel programming framework as shown in Figure 1.In this case, the upper level optimization problem corresponds to the DG owner; he must decide over the size, location, and contract price of the DG units (action).Once this information is known by the DisCo (located in the lower level optimization problem), he uses it as parameters for an optimal dispatch and determines the amount of energy he would be willing to purchase from the DG units (reaction).The DG owner now knows how much energy would be bought from his units, and is able to compute his expected profits.With this information, he might want to change his strategy testing new locations, contract prices, and sizes of DG.The process is repeated until an equilibrium is reached.That is, until there are not incentives for any of the agents to change their decisions.The idea behind this leader-follower game lies in the fact that the leader must consider the reaction of the follower before setting his strategy.That is, the DG owner takes into account the impact of the DG units in the network before deciding on their parameters.Equations ( 1)-( 20) describe the mathematical model of the proposed methodology.

Upper Level Optimization Problem
The objective function of the upper level optimization problem consists on the maximization of profits of the DG investor given by (1).In this expression, J represents the candidate locations for the DG units.The first term of the objective function corresponds to the income due to energy sales, where CP DGj is the contract price in $/MWh, associated to the j th DG unit; P DGj (t) is the power dispatched of the j th DG unit in time interval t in MW; and C DGj is the production cost of the j th DG unit in $/MWh.Whereas the second term corresponds to the investment cost, where IC DGj (t) is the portion of the installation cost of the j th DG unit corresponding to time interval t; P Max DGj is the size of the j th DG unit in MW and ∆t is the length of the time interval.The objective function given by ( 1) is subject to the minimum (P Min ) and maximum (P Max ) size limits of DG units, and the actual size selected (P Max DGj ) (Equation ( 2)); and, the fact that the DG owner is willing to operate k DG units (Equation ( 3)).Equations ( 4)-( 6) define the decision variables for the location (y i ), size (P Max DGj ) and contract price (CP DGj ) of each DG unit j ∈ J. Constraint (6) avoids selling energy from DG units that are not under operation.In this constraint, CP max represents the maximum admissible contract price.Finally, the reaction of the DisCo expressed in (7) represents the lower level optimization problem: Subject to: Reaction of the Distribution System Operator (7)

Lower Level Optimization Problem
The lower level optimization problem (reaction of the DisCo) is given by Equations ( 8)- (20).For a given candidate solution of the DG investor regarding location, sizing, and contract pricing of the DG units, the DisCo reacts by purchasing a certain amount of energy from the DG units.The DisCo aims to minimize the cost of meeting the expected demand (8).The first term of Equation ( 8) corresponds to the cost of energy provided by the wholesale energy market; ρ SE k (t) is the energy price at substation k in time interval t in $/MWh and P SEk (t) corresponds to the power provided to the network through substation k in time interval t.The second term of Equation ( 8) corresponds to the payments associated with the DG units within the DisCo's network.Equations ( 9) and ( 10) correspond to active and reactive power balance constraints.P G n (t) and P D n (t) are the active power generated and demanded in bus n, in time interval t, respectively, while P n (t) corresponds to net active power injected in bus n in time interval t.Q G n (t) and Q D n (t) are the reactive power generated and demanded in bus n, in time interval t, respectively, while Q n (t) corresponds to net reactive power injected in bus n in time interval t.Equations ( 11) and ( 12) correspond to the expressions of active and reactive power injected in bus n, respectively.G nm and B nm are the real and imaginary entries in position (m, n) of the system admittance matrix; θ nm (t) corresponds to the angle between buses n,m in time interval t and V n (t) corresponds to the voltage magnitude in bus n in time interval t.Equations ( 13) and ( 14) represent active power limits supplied by substations and DG units, upper scripts Max and Min stand for maximum and minimum limits, respectively.Equations ( 15) and ( 16) represent limits on reactive power and voltage magnitudes, respectively.Equation ( 17) represents apparent power limits in every line lnm.Equation ( 18) stands for the active and reactive components of apparent power.Equations ( 19) and ( 20) represent the mathematical expressions of active and reactive power flow on lines, respectively: For a given location, sizing and contract pricing of the DG units (upper level solution), the optimization problem given by ( 8)-( 20) is basically an optimal dispatch.

Illustrative Example
From an economic point of view, the straightforward solution for the DisCo would be to purchase energy from the cheapest source.However, the solution to the problem given by ( 1)-( 20) is not that simple, since the DisCo must take into account not only the price of energy but also its impact in the network (power losses, voltage profile, congestion, etc.).An effective way to consider these aspects is by means of an optimal dispatch based on an AC optimal power flow.Figure 2 illustrates such a situation.Consider the 5-bus distribution system depicted in Figure 2a.A demand of 10 MW and 5 Mvar is considered for all buses, and an impedance of 0.01 + 0.1 j Ω is considered for all lines.Suppose that there is only one single DG unit to be allocated to either in bus 2 or bus 5 that provides the total demand of the bus.Both DG locations are advantageous to the DisCo in terms of improvement of voltage profile (see Figure 2b) and the reduction of power losses (see Figure 2c).However, it is evident that the DisCo would prefer the DG to be located in bus 5.If the DisCo is willing to pay a higher price for the energy supplied in bus 5 than for the one supplied in bus 2, there is an incentive for the DG investor to raise the price of the contract and propose the location of the DG unit in bus 5.
Within the bilevel programming scheme given by ( 1)-( 20), the decision of the upper level problem would be to install a DG unit in bus 5 at a given energy price.With this information the DisCo evaluates the impact of such a DG unit solving the lower level problem given by ( 8)- (20).Note that in this case the DisCo considers technical (network limits) and economic issues (price offers).The solution of the lower level optimization problem determines how much energy the DisCo would be willing to buy from the DG owner.With this information, the DG owner determines his profits.Regarding DG size, there is a limit from which the DG creates counter flows and causes technical problems in the network.Such a situation, as well as other technical limitations can be accounted for as additional constraints in the lower level optimization problem.

Solution Approach
While genetic algorithms rely mostly on random decisions in their selection, crossover, and mutation operators, SS is a population-based heuristic that relies on deterministic procedures [28].In power systems, SS applications include the planning of transmission network expansion [29], electric power load forecasting [30], and bi-objective environmental/economic dispatch [31].
Developed in the 1970s, recent applications of SS show its effectiveness to solve nonlinear optimization problems [32], and in optimization problems where the evaluation of the objective function implies costly computations (e.g., black-box optimization [33] and simulation optimization [34]).This advantage makes SS a good candidate for the solution of the bilevel optimization problem in models (1)-( 20) since every feasible solution of (1)-( 7) requires the solution of an optimal dispatch to find the reaction of the DisCo (problems ( 8)-( 20)).

Scatter Search General Structure
SS combines and improves systematically the elements stored in a small reference set (RefSet) of good solutions found during the search process.Algorithm 1 depicts the general structure of the proposed SS for DG optimization.

Solution Approach
While genetic algorithms rely mostly on random decisions in their selection, crossover, and mutation operators, SS is a population-based heuristic that relies on deterministic procedures [28].In power systems, SS applications include the planning of transmission network expansion [29], electric power load forecasting [30], and bi-objective environmental/economic dispatch [31].
Developed in the 1970s, recent applications of SS show its effectiveness to solve nonlinear optimization problems [32], and in optimization problems where the evaluation of the objective function implies costly computations (e.g., black-box optimization [33] and simulation optimization [34]).This advantage makes SS a good candidate for the solution of the bilevel optimization problem in models (1)-( 20) since every feasible solution of (1)-( 7) requires the solution of an optimal dispatch to find the reaction of the DisCo (problems ( 8)-( 20)).

Scatter Search General Structure
SS combines and improves systematically the elements stored in a small reference set (RefSet) of good solutions found during the search process.Algorithm 1 depicts the general structure of the proposed SS for DG optimization.As Algorithm 1 shows, SS initially builds a set P of diverse solution using the diversification generation method (Lines 1-7).Once the number of desired diverse solutions (n) is reached, SS builds the initial Re f Set (Lines 8-11) with the best b/2 solutions in P.Moreover, to favor diversity, the remaining b/2 solutions added to Re f Set are those in P more distant to solutions already in Re f Set.Then, the main loop (Lines 13-22) of SS iterates until no new solutions are inserted in Re f Set.This main loop follows three steps for all pairs of solutions in Re f Set (that includes at least one new solution): (i) combines the two solutions to obtain a new solution S (Line 17); (ii) applies the solution improvement method to S (Line 18); and (iii) tries to update the Re f Set by inserting S on it (Line 19).SS finishes when no new solutions are inserted in the Re f Set.SS returns the best solution of the DG optimization problem in Re f Set (the solution with the largest profit for the DG owner) (Line 23).

Scatter Search Components
As it can be seen in Algorithm 1, the proposed SS uses four methods, namely: (i) diversification generator method (DiverseGeneration); (ii) solution combination method (Combine); (iii) improvement method (Improve); and (iv) reference set update method (U pdate).This section describes the details of each one of these methods.Likewise, two other important components for the SS implementation are discussed: the solution representation and the distance metric needed to maintain the diversity of the solutions in Re f Set.

Solution Representation and Objective Function Evaluation
To represent the solutions of the DG optimization problem we use the two-row array depicted in Figure 3a.The position in the array correspond to the elements in J, the first row corresponds to the size of the DG located in the j-th element of J (P Max DGj ).To reduce the search space, the possible sizes of the DGs are constrained to a discrete set of options T = {P S1 = P Min , P 2 , . . ., P |T| = P Max }.Therefore, as can be seen in Figure 3, the first row indicates the index in T and not the value of P Max DGj .On the other hand, the second row represents directly the contract price (CP DGj ) in monetary units.Figure 3b presents an example of a solution with 10 candidate buses in J where the DG owner seeks to locate 3 DGs (k = 3) with three possible sizes (|T| = 3).In this solution DGs are located in the 3rd, 4th and 10th elements of J.The DG located at the 3rd candidate bus is of maximum size, whereas the two others are of the smallest possible size.Therefore, as can be seen in Figure 3, the first row indicates the index in and not the value of .On the other hand, the second row represents directly the contract price ( ) in monetary units.Figure 3b presents an example of a solution with 10 candidate buses in where the DG owner seeks to locate 3 DGs ( = 3) with three possible sizes (| | = 3).In this solution DGs are located in the 3rd, 4th and 10th elements of .The DG located at the 3rd candidate bus is of maximum size, whereas the two others are of the smallest possible size.Figure 4 illustrates the evaluation process of the example solution of Figure 3.The leftmost part of Figure 4 presents the DG owner decisions for problems (1)-( 7) (i.e., location, size, and contract price of their DG units).With this information, the DisCo (rightmost part of Figure 4) solves an optimal power dispatch (Equations ( 8)-( 20)) to decide how much energy to buy from the wholesale market and the DG units to serve the expected demand at minimum cost.The DisCo returns this information to the DG owner so that he can calculate his profits (Equation ( 1)).This process repeats every time a new solution is found by the SS.
Energies 2017, 10, 1449 8 of 15 Figure 4 illustrates the evaluation process of the example solution of Figure 3.The leftmost part of Figure 4 presents the DG owner decisions for problems (1)-( 7) (i.e., location, size, and contract price of their DG units).With this information, the DisCo (rightmost part of Figure 4) solves an optimal power dispatch (Equations ( 8)-( 20)) to decide how much energy to buy from the wholesale market and the DG units to serve the expected demand at minimum cost.The DisCo returns this information to the DG owner so that he can calculate his profits (Equation ( 1)).This process repeats every time a new solution is found by the SS.

Distance Measure
SS seeks to balance the quality and the diversity of the initial set of solutions (Lines 9 and 10 of Algorithm 1).This creates the need to have a distance measure to compare the two solutions.To accomplish this task we use a simple distance measure that only compares the location of the DGs, ignoring the sizes and contract prices.The main motivation to ignore these two variables in the distance measure is the fact that the improvement operator modifies these variables.Therefore, after improvement, two solutions with the same locations for the DGs may result in the same local optima (with identical sizes and prices).For a given pair of solutions and , the distance measure ( , ) operates in two steps.First, according to Equation ( 21), it uses the values of the first row of each solution to derive the binary variables ( ), indicating whether or not a DG is located in candidate bus .Then, in the second step it uses the Hamming distance of the location variables to calculate ( , ) using Equation (22).Where ( ) indicates the value of the variable for the candidate bus is solution ( ), respectively.Using ( , ) in Line 10 of Algorithm 1, we add sequentially the solution ∈ that is more distant to all solutions already in until completing solution (i.e.,

Distance Measure
SS seeks to balance the quality and the diversity of the initial set of solutions P (Lines 9 and 10 of Algorithm 1).This creates the need to have a distance measure to compare the two solutions.
To accomplish this task we use a simple distance measure that only compares the location of the DGs, ignoring the sizes and contract prices.The main motivation to ignore these two variables in the distance measure is the fact that the improvement operator modifies these variables.Therefore, after improvement, two solutions with the same locations for the DGs may result in the same local optima (with identical sizes and prices).For a given pair of solutions S and S , the distance measure dist(S, S ) operates in two steps.First, according to Equation ( 21), it uses the values of the first row of each solution to derive the binary variables (y j ), indicating whether or not a DG is located in candidate bus j.Then, in the second step it uses the Hamming distance of the location variables to calculate dist(S, S ) using Equation (22).Where y S j (y S j ) indicates the value of the y variable for the candidate bus j is solution S (S ), respectively.Using dist(S, S ) in Line 10 of Algorithm 1, we add sequentially

Diversification Generator Method
The first step of the SS generates a large set of diverse solutions P of size n, from which it selects the initial solutions in Re f Set (Lines 2-7 of Algorithm 1).We consider three different diversification generator methods (DiverseGeneration).The first one is a completely random (Rand) method that selects randomly the location, size and price of the solutions in P. By contrast, the second one (Sist) generates diverse solutions in a deterministic way as follows.For the location, we use the procedure described in [35] that generates (systematically) binary vectors with k values different from zero.Then, for the size we use the maximum size (P |S| ) and for the contract price also the maximum admissible price (CP max ).The third method (SistRand) combines the two previous ones by generating the location systematically according to [35], but generating the size and the contract price randomly.

Solution Combination Method
To combine pairs of solutions in Re f Set (Line 17 of Algorithm 1), we use only the location decisions (i.e., the y variables).Therefore, we use again Equation (21) to derive the location decisions of a given solution S.Then, with the binary vectors of y variables we apply a deterministic variant of the combination method for the knapsack problem described in [36] (p.63).This combination method uses the objective function of the solutions ( f (S)) to bias the selection towards the candidate locations of high-quality solutions.
Method (Combine(S, S )) generates a new solution following four steps: (i) it calculates a score (score(j)) for each candidate location in J according to Equation ( 23); (ii) then, it sorts the locations in non-decreasing order of the score; (iii) it picks the first k to locate there the DGs; and finally; and (iv) it sets the size and contract price to random values as in the randomized diversification generation method:

Improvement Method
Since the distance operator, diversification generation and combination methods focus on the location variables, we decided to apply the improvement method (Line 18 of Algorithm 1) to the other decisions of the DG optimization problem.Therefore, the improvement method (Improve(S)) applies a local search to the DG size (P Max DGj ) or contract price (CP DGj ) of solution S. Since the evaluation of a candidate solution implies the call of an (time-consuming) optimal dispatch, each time the improvement method is called, it operates with the same probability in only one of these decisions for a randomly selected DG of the solution.If the improvement is applied to the size (contract price), the method evaluates the increase/decrease of ±0.5 MW (±0.5 $/MWh) systematically until no more improvements can be achieved.For the initial solution, both increase and decrease of the size (contract price) are considered, and the direction (increase/decrease) with the larger improvement (if any) is selected for the next iterations.The local search stops as soon as no improvement is obtained or the minimum/maximum values for the variables are reached.

Reference Set Update Method
When trying to insert a new solution in the Re f Set (Line 19 of Algorithm 1), procedure U pdate(Re f Set, S) returns a Boolean value (new), that is true only if S is inserted in Re f Set.The following condition must be valid to insert S in Re f Set : ∃S ∈ Re f Set : f (S) > f (S ).That is, solution S must be better than at least one solution already contained in Re f Set.As pointed out before, SS tries to balance the quality and the diversity of the solutions in Re f Set.Therefore, once the update is successful (i.e., new = true) the distance measure is used to select which solution is replaced.In this case, the solution removed from Re f Set (S out ) is the closest to S with a worse objective function (i.e., S out = argmin S ∈Re f Set: f (S)> f (S ) {dist(S, S )}).The update of the reference set is then performed with the following operation: Re f Set = Re f Set ∪ {S}\{S out }.

Tests and Results
The SS heuristic was implemented in Matlab.We used the Matpower software [37] to solve the optimal dispatch (Equations ( 8)-( 20)) required to evaluate the reaction of the DisCo for each solution found during the search.After some preliminary experiments, we set the size of the initial population of diverse solutions (P) to n = 20 and the size of the Re f Set to b = 6.These small values provide short running times without compromising the solution quality.All of the experiments reported in this section were performed in a computer with an Intel Pentium processor with 2 GB of RAM running under Windows 7.
To evaluate the performance of SS we use a distribution system with 34 buses.Figure 5 depicts the structure of this system.Appendix A provides the detailed information of the buses, lines and demands of this system.To evaluate the performance of SS we use a distribution system with 34 buses.Figure 5 depicts the structure of this system.Appendix A provides the detailed information of the buses, lines and demands of this system.

Diversification Generator Method Comparison
This first experiment compares the performance of the three diversification generation methods ( , ), as described in Section 3.2.3.We performed three runs of each variant of the SS, with different seeds for the random number generator.Tables 1 and 2 present the results of this experiment.As it can be seen in Table 1, the randomized diversification generation method achieves the best results.Even its worst results outperform the other variants.As expected, the deterministic variant (SS-Sist) obtains the same results in all of the runs.The detailed information of Table 2 shows a convergence towards the same contract price and size in all of the best solutions, being the main difference the location of the DG units.This result can be explained by the effect of the improvement method, which searches for the best size and contract price, but does not affect the location of the DG units.Moreover, this result also supports the SS design that only takes into account the location of the DG units while combining and comparing the solutions in the RefSet.

Diversification Generator Method Comparison
This first experiment compares the performance of the three diversification generation methods (Sist, Rand and SistRand), as described in Section 3.2.3.We performed three runs of each variant of the SS, with different seeds for the random number generator.Tables 1 and 2 present the results of this experiment.As it can be seen in Table 1, the randomized diversification generation method achieves the best results.Even its worst results outperform the other variants.As expected, the deterministic variant (SS-Sist) obtains the same results in all of the runs.The detailed information of Table 2 shows a convergence towards the same contract price and size in all of the best solutions, being the main difference the location of the DG units.This result can be explained by the effect of the improvement method, which searches for the best size and contract price, but does not affect the location of the DG units.Moreover, this result also supports the SS design that only takes into account the location of the DG units while combining and comparing the solutions in the RefSet.Figure 6 depicts the results of the best runs of the three variants.As it can be seen, SS-Sist obtains very good results in the initial iterations, but is unable to improve the quality of the solutions in Re f Set during many iterations.On the other hand, the randomized variants begin with worst solutions but steadily improve their quality.

Comparison against Other Metaheuristics
In a second experiment, we compare the results of SS against those obtained with a memetic algorithm (i.e., a genetic algorithm hybridized with local search [38]).As a benchmark, we use the memetic algorithm (MA) developed in [39].This MA extends the work of [21] to include the size of the DG as a decision variable.The general control of this MA is based on the work of Beasley and Chu [40]; as the SS, this MA uses a steady-state population replacement strategy, where only one chromosome is replaced in each iteration; furthermore, both heuristics share the solution representation and the improvement procedure (local search in the MA).On the other hand, there are several differences between the MA and the SS.First, while the MA uses a classical one-point crossover, the SS uses a deterministic voting procedure.Second, the MA relies on a (randomized) mutation operator to preserve the diversity of its population, whereas SS uses the distance measure and the reference set updating procedure as a mechanism to obtain a diverse .Third, while SS has a small reference set ( = 6), the MA uses a large population with 100 chromosomes.Furthermore, to evaluate the impact of the solution improvement component (local search), we used as additional benchmark a simple genetic algorithm (GA).This simple method is obtained by removing the local search from the MA.Table 2 compares the results of the MA, the GA, and the proposed SS.For comparison we select the best and worst variants of the SS (SS-Rand and SS-Sist).

Comparison against Other Metaheuristics
In a second experiment, we compare the results of SS against those obtained with a memetic algorithm (i.e., a genetic algorithm hybridized with local search [38]).As a benchmark, we use the memetic algorithm (MA) developed in [39].This MA extends the work of [21] to include the size of the DG as a decision variable.The general control of this MA is based on the work of Beasley and Chu [40]; as the SS, this MA uses a steady-state population replacement strategy, where only one chromosome is replaced in each iteration; furthermore, both heuristics share the solution representation and the improvement procedure (local search in the MA).On the other hand, there are several differences between the MA and the SS.First, while the MA uses a classical one-point crossover, the SS uses a deterministic voting procedure.Second, the MA relies on a (randomized) mutation operator to preserve the diversity of its population, whereas SS uses the distance measure and the reference set updating procedure as a mechanism to obtain a diverse Re f Set.Third, while SS has a small reference set (b = 6), the MA uses a large population with 100 chromosomes.Furthermore, to evaluate the impact of the solution improvement component (local search), we used as additional benchmark a simple genetic algorithm (GA).This simple method is obtained by removing the local search from the MA.Table 2 compares the results of the MA, the GA, and the proposed SS.For comparison we select the best and worst variants of the SS (SS-Rand and SS-Sist).Three runs of both methods were performed.Results for the MA and GA are taken from [39].
As it can be seen in Table 3, the proposed SS outperforms both the MA and the GA, both in solution quality and running time.The best SS solution offers 5.8% more profit to the DG owner than the best MA solution.Even the worst variant of SS (Sist) obtained solutions with 4.7% more profit.Remarkably, SS took a small fraction of the running time of the MA to reach better results (only 0.13% of the time of MA).According to [39], the MA took 3.5 h of running time on average to solve the DG optimization problem.Clearly, removing the local search from the MA produces a much faster method.The GA only took 20 min of running time on average.However, the quality of the solutions deteriorates.The best GA solution offers only 21% of the profit of the worst SS solution.Moreover, both methods were implemented in Matlab, but the MA experiments were run in a much faster computer with an Intel Core i5 processor (with 4 GB of RAM).These results show that SS offers a good option to solve DG optimization problems where the objective function evaluation consumes long running times.

Conclusions
This paper tackles the optimal location, sizing, and contract pricing of DG in electric distribution systems.To model this problem we follow a bilevel programming approach that considers the interaction between the DisCo and the DG owner.
These joint decisions (location, sizing, and contract price for the DGs) represent a challenge since the complexity of electric distribution systems hinders the application of exact optimization methods.Therefore, population metaheuristics like genetic algorithms and particle swarm optimization, among others, are one of the preferred options to solve this problem.As an alternative, this paper proposes a SS that relies (mostly) on deterministic mechanisms to preserve the quality and diversity of a small set of reference solutions.
Computational experiments on a test system with 34 buses reveal the effectiveness of the proposed SS.The comparison against a memetic algorithm and a genetic algorithm shows that the proposed heuristic is able to find better solutions in much shorter running times.The use of a small set of solutions and the deterministic diversity preserving mechanism of SS represent an advantage in problems where the objective function evaluation is a time-consuming task (as the DG optimization problem tackled in this paper).
Future research directions include the use of SS for multiobjective variants of the problem considering other performance measures (power losses, reliability, voltage stability, etc.), in addition to the DG owner profit.Likewise, the use of SS for the solution of other optimization problems Energies 2017, 10, 1449 13 of 16 that arise in the planning and operation of electric distribution systems (network configuration and reconfiguration, reliability optimization, etc.) appears as a promising alternative.

Figure 2 .
Figure 2. (a) 5-bus distribution system; (b) impact of distributed generation (DG) location in voltage profile; and (c) impact of DG location in power losses.

Figure 2 .
Figure 2. (a) 5-bus distribution system; (b) impact of distributed generation (DG) location in voltage profile; and (c) impact of DG location in power losses.
solution S d ∈ P that is more distant to all solutions already in Re f Set until completing b solution (i.e., S d = argmax S ∈P\Re f Set {min S∈Re f Set dist(S, S )}): Energies 2017, 10, 1449 10 of 15

15 Figure 6 .
Figure 6.Progress of the best solution in RefSet for different scatter search variants.

Figure 6 .
Figure 6.Progress of the best solution in RefSet for different scatter search variants.

Figure A1 .
Figure A1.Aggregated demand curve and price in the wholesale energy market.Figure A1.Aggregated demand curve and price in the wholesale energy market.
Set with the best b/2 solutions of P, Re f Set = {S 1 , . . .S b/2 } 10 Add to Re f Set the most diverse b/2 solutions in P with respect to those already in Re f Set * = argmax S∈Re f Set { f (S)}

Table 1 .
Results for scatter search (SS) with different diversification generation methods.

Table 1 .
Results for scatter search (SS) with different diversification generation methods.

Table 2 .
Best solutions of each SS variant for the 34-bus distribution system.

Table 3 .
Results for SS, memetic algorithm (MA) and genetic algorithm (GA) on the 34-bus distribution system.

Table A2 .
Line data for the 34-bus distribution system.Aggregated demand curve and price in the wholesale energy market.