A Benders’ Decomposition Approach for Renewable Generation Investment in Distribution Systems

: A model suitable to obtain where and when renewable energy sources (RES) should be allocated as part of generation planning in distribution systems is formulated. The proposed model starts from an existing two ‐ stage stochastic mixed ‐ integer linear programming (MILP) problem including investment and scenario ‐ dependent operation decisions. The aim is to minimize photovoltaic and wind investment costs, operation costs, as well as total substation costs including the cost of the energy bought from substations and energy losses. A new Benders’ decomposition framework is used to decouple the problem between investment and operation decisions, where the latter can be further decomposed into a set of smaller problems per scenario and planning period. The model is applied to a 34 ‐ bus system and a comparison with a MILP model is presented to show the advantages of the model proposed.


Introduction
Distributed generation (DG) has been used to produce energy in remote and isolated places, where the distance between the demand and the producer is short [1]. At present, this trend is changing as it has been proven that DG provides technical, economic, and environmental improvements [2]. The benefits of using DG are based on the reduction of network losses, the voltage level improvement, or the dependence reduction of energy, fuel prices, and traditional generation. All these advantages cause CO2 emission reduction [3,4]. On the other hand, DG has also disadvantages, mainly related to technical aspects due to the fact that the existing networks have not been designed to incorporate this type of generation [2]. Some drawbacks are reverse current flows, the need for network redesign, or frequency instability [5].
In addition, planning of distribution networks must account for renewable generation uncertainty to meet the future demand in any possible future scenario. In that case, the decisions to be made could be to invest in the network, in substations, in DG generation, or in any combination of them [6]. When the chosen option is DG, it is necessary to determine the type of location, by solving the Distribution Generation Planning (DGP) problem. An example of this sort of problems can be seen in [7], where a linear model that optimizes size and location of DG is used. The objective function maximizes DG real power. In [6], a model is proposed for minimizing the investment and operation DG costs, the cost of the electricity bought from the substation, and the cost of network losses of a distribution company. The model is formulated as a mixed-integer-nonlinear one. Other authors have worked with particle swarm optimization methods [8], where a multi-objective model considers investment and operating costs of new generation, the cost of the energy purchased, and CO2 emission cost.
1. The stochastic multi-stage DGP problem is formulated using Bendersʹ decomposition that decomposes the problem per both scenario and planning period and 2. A comparison with a standard MILP model is provided and results with a 34-bus test case are shown, where the significant computational advantage of using Bendersʹ decomposition with respect to the MILP model is also shown.

Notation
Subscripts t, k, and ω below refer to the values in year t, time block k, and scenario ω, respectively. Superindex υ indicates the values in the υ-th Benders' iteration. All indices, sets, constants, parameters and variables used in the document are shown in Table 1. Value of the r-th block associated with the reactive power through branch , [MVAr ].

Objective Function
The goal of the model is to minimize the total system costs (TSC) considering DG. The model utilizes a two-stage stochastic mixed-integer linear programming model. The investment variables that do not depend on the scenarios are established in the first stage. In the second stage, dependent or stochastic operation variables that depend on the scenarios are determined.
TSC are composed of two terms (Equation (1)). The first term corresponds to the first-stage variables that determine the number of new units such as wind turbines, PV modules, and transformers to install. The second term corresponds to the second-stage variables whose values are obtained after the outcome of scenario ω is known. See [18] for details. (1) Total costs are updated using the present worth factor, 1/ 1 . The interest rate, i, is used to calculate the annualized investment cost payment rates. The payment contributions for the three technologies, namely, transformers, wind turbines, and PV modules, are obtained in Equations (2), (3), and (4), respectively.
After having determined the payment rates, the investment costs are calculated as shown in Equations (5) and (6).
The operation and maintenance total costs (Equation (7)) take into account the cost of losses (Equation (8)), the penalty for the energy not supplied (Equation (9)), the cost of purchase of energy from substations (Equation (10)), as well as renewable DG candidates' operation and maintenance costs (Equation (11)).
The maximum value of the power allowed to install at each bus n is limited by Equation (15).
, , The investments are limited by Equations (16) and (17). Equation (16) refers to the annual investment payment limit that describes the budget available for each investment period, while the portfolio investment (Equation (17)) is related to the amount of money that is available for investment in the long term.
, ∈ , , The power flow constraints (Equations (18)-(41)) refer to both active and reactive power equations. The distribution system is composed of two types of buses (load and substation buses) where each kind of bus has a different load flow expression.
The voltage is related to the different electrical values (Equation (22)) and is limited by the maximum and minimum values (Equation (23)).
The substation available power is subject to yearly updates (44). New investments are made considering the already existing power and the new possible power of the candidate transformers (45), (46). ; ; ; ; The renewable power (wind turbines and PV modules) available is updated every year. Note that the vector of the candidate nodes appears in Equations (47) to (50). , , , , The maximum amount of generation that is available is limited by the generation levels of each of the technologies available (Equations (51), (52)).
The active power injected into the distribution system must keep a constant power factor (Equation (53)). In this way the reactive power is also limited (Equation (54)). , , , , , , , , The reactive renewable power that is injected in the distribution system must have a constant power factor (Equations (55) and (56)).

Solution Procedure: Benders' Decomposition
Model (1)-(56) is set as a MILP problem whose size depends on the number of time blocks, scenarios, and years considered in the study. In order to obtain informed expansion decisions, a large number of time blocks and scenarios are needed to obtain a good representation of the loads and renewable power outputs, as well as a large number of years. However, in such a case the MILP problem (1)-(56) may become computationally intractable. Nevertheless, note that if variables ; , ; , and ; , i.e., the expansion decision variables, are fixed to given values, then, problem (1)-(56) can be decomposed into a set of smaller problems, one for each time block, scenario, and year. This allows us to use Benders' decomposition to solve the problem.
Model (1)-(56) uses binary variables for the linearization of power flow equations. Since Benders' decomposition cannot have binary variables in the sub-problem, a procedure to fix the values of these variables must be used. Further details are shown in step 1.
In the following sections superindex (υ) indicates the values in the υ-th Benders' iteration. The proposed Benders' decomposition algorithm comprises the steps below: 1.
Step 1. Sub-problem solution: The variables ; , ; , and ; , ∀ , are set, to their optimal values from the previous step and solve the following sub-problems, one for each year , time block , and scenario .
Next, the lower bound of the objective function is updated (1): (75) Then, the algorithm continues at Step 1. For clarity purposes, the flowchart of the algorithm is shown in Figure 1. Note that the Benders' decomposition algorithm described in this section can solve each sub-problem separately, and hence, the computational time of the proposed algorithm grows linearly with the number of scenarios, time blocks, and years, being its scalability excellent.

Case Study Data
This section describes the data used for the simulations of the problem explained in Section 3. The Benders' decomposition method is applied to a 34-bus three-phase radial feeder [21], which has 1 substation, 29 buses with load, and 5 buses with no load. The topology of the system is shown in Figure 2. The data are based on those provided in [18] and [21], but adding new scenarios. In [18]  ; ; Step 0 Step 1 Step 2 Step 3 there are 3 levels of demand, wind production, and solar production per each of the 8 time blocks. Now, there are 4 levels and 12 time blocks. The wind turbines chosen have a capacity of 100 kW and the PV modules have power of 2.5 kW. The candidate buses for each technology are shown in Figure 2 and the maximum power (wind turbine and PV modules) that may be installed at each bus is 250 kW. The maximum number of units is limited to 2 for wind turbines and to 85 for PV modules, and the substation can be expanded adding transformers of sizes ranging from 1 MVA up to a maximum of 5 MVA. Demand increases 2% each year along the planning horizon (20 years). The voltage data of the distribution system is 1 pu in the substation node and the minimum and maximum allowable voltage values are 0.95 pu and 1.05 pu, respectively. The values defined for the interest and discount rates are 8% and 12.5% [6], respectively. Investment data of new devices are shown in Table 2. The operation and maintenance costs of the new renewable production units are €0.007/kWh [8]. The annual budget is €150,000 and the maximum portfolio investment for the life time of the devices is €5,500,000. Data on demand, wind speed, and PV factors used per time block are shown in Table 3, where wind and PV production levels are also displayed. The different levels are combined with each other in every time block. There are four levels of demand, wind, and PV production, hence, the total number combinations (scenarios) per time block is 64. Hence, there are 12 time blocks, 6 per season (winter (October-March) and summer (April-September)), accounting for 768 scenarios. Note that energy prices are not used in scenario generation and each price corresponds with a defined demand factor. These prices are shown in Table 3 and increase 1% each year with respect to the base year. All the levels in each time block have a probability of 1/4. Therefore, the weight of each of the scenarios within each time block is 1/64 [22]. Two blocks are used for piecewise linearization, where the cost of energy not supplied is €15,000/MWh. The tolerance ε of the Benders' algorithm is specified as €1.

Results Discussion
Two case studies have been simulated to test the model, where the constraints related to installed power and limits of investment are different. The results of Benders' algorithm for each case study case are compared with the MILP model given in Equations (1)-(56). These results are obtained using CPLEX 11 under GAMS [23] on an Intel Xeon E7-4820 computer with 4 processors at 2 GHz and 128 GB of RAM. Table 4 presents the numerical results of each simulation in each study case for 768 scenarios. The relative gap is set to 0.01 in all simulations for Benders' simulation.  Table 5). The expansion of the substation is carried out in the ninth year, with a new transformer. The new PV devices are installed at the end of the branches because it reduces the costs associated with energies losses. The location of the PV modules within the network is displayed in Table 6. The CPU time decreases by 97.9% when Benders' decomposition is the method chosen for the simulations.  (16) and (17)) are not considered. This new constraint scenario allows the investment in, not only new transformers and PV modules, but also in wind units (see Table 7). Two expansions of the substation are made, in years 9 and 16, PV modules are installed in all candidate nodes, and six wind units, in total, are also installed in the last year. In this case, the new installed capacity is 4,725 kW. The location of the PV modules within the network is displayed in Table 8. The CPU time decreases by 94.4% using Benders' algorithm.

--
The results are better in Case b than Case a because the problem is less constrained. The introduction of renewable energy in distribution systems reduces the total operation and maintenance (O&M) system costs (see Table 9). The reduction of O&M costs is due a reduction of losses costs and purchase energy by substation. These results highlight the advantages of investing in renewable generation in the long term. Finally, Table 10 shows the CPU times required to solve the DGP problem using both the MILP model (1)-(56) and the Benders' algorithm for different number of scenarios. The CPU time required to solve the MILP model (1)-(56) increases drastically with the number of scenarios (Figures 3 and 4). However, its increment is approximately linear if Benders' algorithm is considered. This causes that, up to 64 scenarios, MILP model solves the problem faster than Benders' algorithm but, from 216 scenarios on, Benders' algorithm becomes a much more efficient way to solve it.

Conclusions
This paper has considered the DGP problem in a stochastic environment, where both PV and wind technologies as well as demand are subject to random changes. The use of Benders' decomposition to solve the two-stage stochastic investment problem has allowed us to further decompose the problem by scenarios and planning periods, making it a fully decomposable one. The model has been tested for a 34-bus example with excellent results. In terms of computing time, the increase in the number of scenarios makes the differences between MILP and Benders' models evident. Up to 64 scenarios, the MILP model is much faster. Benders' model needs to perform iterative processes (loops) that increase the computing time. This trend changes as the number of scenarios increases and Benders' method becomes faster than MILP. In general, MILP's computing time behaves in a non-linear way, whereas Benders' model is more linear in terms of computing time. This proves the significant computational advantage of Benders' with respect to a conventional MILP model.
Note that applying Benders' decomposition may also allow extending the investment problem to address other relevant issues in distribution systems, such as switching or network reconfiguration. This would be non-viable using the original mixed-integer linear programming problem due to its high computational burden.
The work developed in this article can help investors to decide the kind and size of renewable technologies and the place where install the new devices in the distribution system. As improvements to the proposed problem, other kinds of producers (biomass, hydraulic) can be introduced, incorporate reliability in generation and electric vehicle. In addition, making a comparative with other uncertainty management methods such as K-means.