Methodology for Multi-stage, Operations- and Uncertainty-Aware Placement and Sizing of FACTS Devices in a Large Power Transmission System

07/07/2017 ∙ by Vladimir Frolov, et al. ∙ Los Alamos National Laboratory 0

We develop new optimization methodology for planning installation of Flexible Alternating Current Transmission System (FACTS) devices of the parallel and shunt types into large power transmission systems, which allows to delay or avoid installations of generally much more expensive power lines. Methodology takes as an input projected economic development, expressed through a paced growth of the system loads, as well as uncertainties, expressed through multiple scenarios of the growth. We price new devices according to their capacities. Installation cost contributes to the optimization objective in combination with the cost of operations integrated over time and averaged over the scenarios. The multi-stage (-time-frame) optimization aims to achieve a gradual distribution of new resources in space and time. Constraints on the investment budget, or equivalently constraint on building capacity, is introduced at each time frame. Our approach adjusts operationally not only newly installed FACTS devices but also other already existing flexible degrees of freedom. This complex optimization problem is stated using the most general AC Power Flows. Non-linear, non-convex, multiple-scenario and multi-time-frame optimization is resolved via efficient heuristics, consisting of a sequence of alternating Linear Programmings or Quadratic Programmings (depending on the generation cost) and AC-PF solution steps designed to maintain operational feasibility for all scenarios. Computational scalability and application of the newly developed approach is illustrated on the example of the 2736-nodes large Polish system. One most important advantage of the framework is that the optimal capacity of FACTS is build up gradually at each time frame in a limited number of locations, thus allowing to prepare the system better for possible congestion due to future economic and other uncertainties.



There are no comments yet.


page 1

This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.




Number of power lines in operation

Number of buses in the system

Number of loading scenarios representing given time frame

Number of scenarios representing planning horizon

Number of time frames representing horizon

Index of a decision point

Index of a scenario at time frame t

Occurrence probability of a scenario

at time frame

Vector of initial line inductances


Vector of maximum (minimum) active power generator outputs


Vector of maximum (minimum) reactive power generator outputs


Vector of active (reactive) power demands

Vector of line apparent power limits


Vector of maximum (minimum) allowed voltages

Cost per Ohm of a series FACTS device

Cost per MVAr of a shunt FACTS device

Planning horizon

Optimization variables (operational, scenario dependent):

Vector of bus voltage magnitudes

Vector of bus voltage angles

Vector of generator active power injections

Vector of generator reactive power injections

Vector of line inductances modified by SC devices

Vector of series FACTS settings

Vector of shunt FACTS settings

Optimization variables (investment, scenario independent):

Vector of series FACTS capacities built at decision point t

Vector of shunt FACTS capacities built at decision point t

I Introduction

Energy market deregulation and massive installation of renewables are among significant stress factors forcing transmission systems across the world to operate closer to their limits. When a transmission system is constrained building new lines seems a natural remedy [1, 2, 3, 4]. However, this option is costly and severely limited in many countries due to social and environmental concerns, hence rising the question: if the system can be upgraded creatively and gracefully by installing Flexible Alternating Current Transmission System (FACTS) devices which are both less expensive and whose spatial footprint is much smaller?

This question was addressed in the literature [5, 6, 7, 8, 9] suggesting that FACTS devices may indeed be effective in increasing transmission capacity and keeping the system safe and operational as the demand growth. However, it also became clear from the studies that this sustainable upgrade needs to be creative – types, locations and capacities of the newly introduced FACTS devices must be chosen carefully in order to exploit their benefits. Several objectives, including decreasing operational [7, 10, 11] and investment costs [5, 9], reducing transmission losses [6, 10], increasing power system loadability and managing system congestion [12, 10, 13, 14], reducing load curtailment [15] and improving voltage profile [9] and voltage stability index [8, 14], have been considered. Formally, these problems were stated as mixed Integer Nonlinear Programming (MINLP) optimizations which exact solutions are limited to only very small (typically not exceeding tens of nodes) power system models. Sensitivity analysis [15, 16], Mixed Integer Linear Programming (MILP) [17]

and genetic algorithms

[18, 19, 20, 21] were suggested to resolve the MINLP formulations approximately. The sensitivity based methods consist in computing a few indicators to identify the most critical lines which are thus suggested as good option for FACTS placement. However, this ad-hoc methodology is obviously not optimal and it also does not offer any suggestion on how to size the devices. Genetic algorithms aim at finding optimal solutions but come with unacceptably high computational cost. Relaxation and approximation techniques were suggested to convert MINLP to MILP [22]. Approximation techniques are computationally less expensive as they use the simplified line flow based model [9] or DC power flow model [17]. However, such approximations have serious limitations especially for planning installation of shunt FACTS devices.

One significant complication in resolving MINLP is related to combinatorial explosion in the number of choices for possible FACTS placements. It is thus desirable to enforce sparcity of the installation. A principal way to achieve sparse placement in a computationally feasible way was suggested in [23, 24]. An alternative way to achieve sparse placement was also discussed in [22].

Uncertainty of the economic growth is another notable obstacle for efficient and practical implementation of this and other related planning problems. A solution for resolving this problem was suggested [25] where the uncertainty was modeled via exogenously described multiple operational samples of loads and related nonlinear AC power flow solutions was embedded in the large-scale optimization explicitly.

In this manuscript we extend the approach of [25] to resolve the last remaining issue in the optimal multiple-scenario aware, AC-based and sparse placement and sizing of the FACTS devices in a large transmission system. Here, we choose to represent future not in one time step, as was done in [25], but in multiple time steps. In other words here we complete general description of the framework started in [23, 24, 25] and propose a comprehensive resolution for finding optimal locations of FACTS devices in a large transmission system by preparing the system for future loading gradually through multi-stage, properly paced investments. Main highlights of our comprehensive approach are as follows:

  1. Planning horizon is represented by multiple decision points (multiple time frames). At each new time frame a set of new FACTS devices can be installed, and they are assumed available for operations immediately such that respective operational values do not exceed the installed capacities. Therefore, installation of FACTS devices can be paced. In this manuscript we work with a finite number, , of the time-horizon sub-intervals. Notice, however, that extension of the approach to the case of a receding horizon, thus accounting at each step for updated forecast, is straightforward. (We plan to conduct extensive discussion of the receding horizon experiments in a future journal version of the manuscript.)

  2. Future operational conditions are represented through multiple loading scenarios and associated probabilities broken into time frames. Our framework is set in the way that the scenarios are stated as an exogeneous input, which allows us to separate the problem of scenario generation from the intrinsic optimization details. Given the exogeneously prescribed scenarios, optimal installation of FACTS is resolved within the optimization framework by accounting for both investment variables and operational variables, characterizing installation decisions and operational implementations (per scenario) respectively. It is important to stress that an optimal constraints (i.e. feasible) all the scenarios. (This is in a contrast with the worst case planning approach.)

  3. Both capital and operational expenditures are optimized simultaneously. To the best of our knowledge no prior works have considered optimizing them at the same time. But for practical planning horizons operational cost is much bigger than the cost of FACTS installation. Thus relatively small additional investment allows to save a significant amount of money by reducing congestion additionally to resolving infeasibility of particular load scenarios.

  4. A novel optimization iterative heuristics is developed, which is a combination of analytic linearization of non-linear constraints, a solution of Quadratic Programming (QP) or Linear Programming (LP) (depending on generation cost) problem for finding of investment variables and operational settings for all scenarios and Alternating Current (AC) Power Flow (PF) resolution (for each scenario) to update previously found states.

  5. The developed heuristics for finding optimal locations of FACTS devices considering multiple loading scenarios and multiple time frames can be applied to large power systems. In other words, the developed approach is scalable. To the best of our knowledge the system considered in the literature addressing optimal locations of FACTS placement consists of a maximum of 1228 buses [12]. The results of the proposed methodology is demonstrated on 2736-bus Polish system, thereby proving its scalability. Moreover, algorithm provides upper bound solutions of the objective function with the gap less than 0.1

The material in the manuscript is organized as follows. Our basic optimization model is introduced in Section II. The solution algorithm is described in Section III. Section IV is at the core of the manuscript - it describes in details our experiments. We conclude and discuss path forward in Section V. Appendixes describe constructions and steps needed to run our algorithm. (This material is largely a repetition from [25] reproduced here for completeness.)

Ii Optimization model

This section describes our optimization framework for operations-aware installation of FACTS devices taking into account multiple future decision points (or multiple time intervals).

Assume that the planning time horizon is , is the number of time intervals, is number of given loading configurations (scenarios) per each time frame (the number of scenarios per time frame may also vary with the time frame). In this setting we aim to place and size the Series Compensation (SC) and Static Var Compensation (SVC) devices, where an SC device, installed at a line, modifies inductance of the line (thus allowing to reroute apparent power), while an SVC device, installed at a node, injects or consumers reactive power at the node thus helping to balance the voltage locally.

Fig. 1 illustrates the setting. Since scenarios are generated within each time interval independently, the total number of paths accounted for within our optimization formulation is , where a path is a sequence of scenarios (each per time interval). Notice that even though the number of paths is exponential in , the total number of the operational constraints in the optimization formulation, , scales linearly in .

The overall problem is to minimize a combination of the sum (over the time intervals) of the investment cost and the sum of operational costs over all the scenarios taking into account (a) operational constraints for every scenarios (per time interval) and (b) investment constraints requiring that the operational variables (for every scenario per tiem interval) do not exceed the respective installed capacities. Operational settings can be different for different scenarios but installed capacities of the devices are the same for all the scenarios representing given time interval.

Fig. 1: Illustration of the relation between the number of scenarios, (each defined in the time interval, ) and the number of time intervals, . is the total number of constraints imposed (per line or per node) within our optimization formulation. (See text for additional explanations.)

Mathematically the optimization problem is stated as follows:

subject to:


where labels the scenarios; upper index labels the time intervals, ; and denotes the set of nodes and the set of (undirected) edges, of the grid-graph, where a node can be of the load type, , or of the generator type, .

The objective function in (II) consists of three terms. The first two, sparsity promoting terms [23, 24], express the capital investment costs of the installation of the two types of FACTS devices (investment can be performed at each decision point ). The third term stands for the operational cost in which the summation is over all the scenarios for each time frame and over the time frames ( is a shortcut for, ) accounting for respective occurrence probability multiplied by the number of years (service period). Therefore, the optimization (II) is nothing but an operational aware planning.

Each scenario in (II) is stated in terms of the set of operational variables. Description of the optimization constraints in (II) is as follows. (3) and (4

) represent total available capacity at the decision moment

. (5) and (6) ensures that the already installed capacities are inherited in the future time frames. (7) bounds actual line inductances, which are adjusted according to the operational value of the installed series compensation for each scenario, within their respective installed capacities (represented by (14)). (8) and (9) represent active and reactive power balances at each bus of the network. Components of the vectors () and () are assumed equal to zero at the buses containing, respectively, no generators or loads. (10) and (11) represent the net active () and reactive () power injections at the system buses. The term expresses shunt compensation by SVC adjusted to a scenario bounded by the respective installed capacities (represented by (15)). Active and reactive power generation limits are set by (12) and (13). Voltage and thermal line flow constraints are represented by (16) and (17). and stand, respectively, for the apparent power flows from to and to from along the line . One also accounts for the budget constraint per time step, (18), and/or for the maximum built capacity constraint per time step, (19).

The main challenge in resolving the optimization is related to nonlinearity of the Power Flow relations (10), (11) and also to nonlinearity of the line thermal limits (17). Available non-linear solvers, such as IPOPT, are not effective in resolving the nonlinearities for large systems efficiently. This has motivated us to develop an heuristic algorithm consisting in sequential linearization of the nonlinear constrains discussed in the following Section (and, specifically, in Subsection III-A.

Iii The Algorithm

This Section describes the algorithm which allows us to resolve efficiently (and in spite of its complexity) the optimization problem just stated. Our algorithm consists of the following steps:

  1. Scenarios are generated for each time frame according to the methodology suggested and described in details in [25]

    . Briefly, one picks the base case, re-scale it for different time frame (taking into account the economic growth), and then introduce fluctuations around the re-scaled solutions to represent the forecasted load uncertainty. The fluctuations are chosen to be Gaussian with the standard deviation proportional to the mean.

  2. Generation is initialized (for each load scenario) according to scheme explained in Appendix C.

  3. If some of the constraints (7)-(17) are violated the initial state of the system is outside of the feasible domain defined by them. The non-linear constraints (10), (11) and (17) are linearized around the current state. This allows to construct current linearized version of the non-linear optimization problem (II)-(17).

  4. The resulting linearized problem is solved by QP (or LP, depends on generation cost functions) using one of the available algorithms of the CPLEX solver [26].

  5. AC power flow (AC-PF) is solved to update the state obtained at the previous step. This step is needed to prepare a feasible solution for the next iteration.

  6. Steps 2-5 are repeated till either no constraints remain violated or the target precision is reached or the maximum allowed number of iterations is reached.

It is important to emphasize that, by construction, the algorithm maintains a feasible physical states at each iteration of its main loop including linearization, solution of the current QP optimization and back projection to the non-linear PF equations (achieved through the AC PF step).

Below we present details of the main steps of the algorithm.

Iii-a Linearization

The Power Flow (PF) equality constraints (10), (11) and the thermal limit constraints (17) are nonlinear. We choose to add to the optimization formulation auxiliary variables – active and reactive power flows expressed via voltages, phases and system parameters (see Appendix A for modeling details). This allows to localize non-linearities in local relations between the line power flows and respective voltages and phases at the two ends of the lines. Details of this technical trick, aimed at improving performance of the CPLEX solver, are as follows.

Introduce the following operational state.


and substitute the left hand side of Eq. (17) by its Taylor expansion around the current state, ,


thus arriving at the two equations representing a line

where the newly introduced auxiliary variables should also be substituted by the respective linearized expressions

Power balance constraints (8), (9), (10), (11) will be exact linear in terms of auxiliary variables.

Three comments/clarifications are in order. First, notice that the operational variables are adjusted independently for each scenario, thus enabling devices’ efficient utilization. Second, to manage possible degeneracy of the resulting system of linear constraints and limit the change of reactive flows at every we add to Eqs. (II)-(19) the following soft constraints for reactive power dispatch at each QP/LP step of the procedure

(In the case of LP, when the degeneracy is stronger, we also add similar constraints imposed on the active power.) Finally, third, to speed up the QP/LP computations one uses a cutting plane (constraint management) procedure. We split the whole set of constraints (17) into “active” and “inactive” sets including the constraints which were overloaded and, respectively, not overloaded, at the current state (of the previous iteration) or at any of the preceding steps. Only active constraints are explicitly accounted for in the optimization, while the validity of the inactive constraints is verified post-factum and the active/passive split is updated at every LP/QP step.

Iii-B QP/LP implementation

Standard CPLEX solver is called at each QP/LP step which outputs operational variables for each scenario along with investment variables for each time frame, and .

Iii-C AC-PF feasibility

The QP/LP step is followed by the AC-PF step, which is needed to maintain the AC PF feasibility destroyed by the linearization. Overall, combination of the QP/LP and AC-PF steps allow to maintain solution and resolve contingencies of the system simultaneously and gracefully.

Iv Case studies

The AC PF and optimization algorithms are implemented in Julia/JuMP. (See [27] and references therein.) QP/LP optimizations (called at internal steps of our algorithm) are resolved by CPLEX [26]. When possible we utilize IPOPT [28], called from JuMP, to solve the optimization problem (in its original, nonlinear formulation). The (brute-force) IPOPT solution is computationally expansive and it is used as a ground truth (to validate our heuristic algorithm). Computational performance of the algorithms is analyzed on a Macbook Pro laptop (Core i7 3.3 GHz (2 Cores), 16 Gb of RAM).

Iv-a Algorithm validation

The algorithm is validated by comparison with the IPOPT solution. Since IPOPT is not able to resolve the 2376 bus-large Polish model (even with a single scenario) we perform the initial validation study on the 30 bus IEEE model. (Both the Polish model and the IEEE 30 bus model are available within the MathPower package [29].) Table I presents results of the (IPOPT vs our heuristics) comparison, where the planning horizon is taken to be year, , and for each time interval. Budget constraint per time interval is set to .

to — X[c] — X[c] — X[c] — X[c] — X[c] — X[c] — t. int.: & 1 & 2 & 3 & 4 & 5
sc.: & 10 & 10 & 10 & 10 & 10
alpha: & 1.02 & 1.04 & 1.06 & 1.08 & 1.15
dev.: & 0.001 & 0.003 & 0.005 & 0.01 & 0.02
IPOPT & & & & &
SVC (bus 8), MVAR: & 4.00006 & 7.99995 & 10.5666 & 10.5666 & 10.5666
SC (line 10), : & 0 & 0 & 1.661575 & 1.661618 & 1.661625
main & & & & &
SVC (bus 8), MVAR: & 4.0 & 8.0 & 10.5241 & 10.5241 & 10.5241
SC (line 10), : & 0 & 0 & 2,321978 & 2,321978 & 2,321978

TABLE I: Comparison of our heuristic algorithm against the (brute-force) IPOPT algorithm for the IEEE 30-bus model.

IPOPT objective is . Main algorithm objective is . Our main algorithm gives upper bounds solutions of the objective function with the gap less than 0.1.

Iv-B Scalability analysis

Our next step (after completion of the aforementioned validation study) was to perform a comparative analysis of the algorithms’ (computational) scalability.

First, we fix the number of scenarios ( per time interval) and study dependence on the number of time intervals. The results are shown in Fig.  2. Then, we consider one time interval and study dependence on the number of scenarios. These results are shown in Fig.  3. Both tests are done still on the 30 bus model with the quadratic generation cost (de-fault in the Mathpower package). In both cases we compare performance of the brute-force IPOPT solver with performance of our main algorithm and of the main algorithm reinforced by the cutting plane.

Comparing performance of the IPOPT in the two setting one observes a drastic difference. In the case of Fig. 2 IPOPT shows a surprisingly good performance, outperforming in speed both of our algorithms. (We relate this good performance of the IPOPT to using the primal simplex option within the IPOPT.) However the situation is reversed in the case of Fig. 3 where, moreover, performance of the IPOPT degrades exponentially with the number of scenarios. Our main algorithm and the reinforced (by cutting plane) algorithm show similar scaling performance in both cases (with the reinforced algorithm performing slightly better). Notice that juxtaposing here our algorithms against each other has a sense because of the comparable number of constraints contributing the two optimization settings.

Fig. 2: Computational time of IPOPT (orange), of our main algorithm (red) and of our algorithm sped up with the cutting plane (dark blue) are shown as functions of the number of time intervals, , for the 30 bus model. In all the tests shown the number of scenarios (per time interval) was . Loading level is unity initially and it increases (imitating economic growth) by the factor per time step. Scenarios are generated with the deviation factor . (See Appendix B for details.) Budget constraint (of per time step) is applied. Our algorithms (with and without cutting plane sped up) are limited to 20 iterations.
Fig. 3: Computational time of IPOPT (orange), of our main algorithm (red) and of our algorithm sped up with the cutting plane (dark blue) are shown as functions of the number of samples in the case of a single time interval for the 30 bus model. Loading level is set to . Scenarios are generated with the deviation factor . (See Appendix B for details.) No budget constraints are applied. Our algorithms (with and without cutting plane sped up) are limited to 20 iterations.
Fig. 4: Computational time of our main algorithm reinforced by the cutting plane shown vs the number of time intervals for the 2736 bus-large Polish model in the case of a single scenario (per time interval). The optimization cost is linear (according to the base case documented in Mathpower) and thus LP is used at each iteration step of our reinforced algorithm solved in iterations. In this case the loading level is set to unity in each time frame. Budget constraint of per time step is applied.

Moving to the scaling analysis of the Polish model, one first of all note that in this case the IPOPT fails to converge. To illustrate performance of our reinforced algorithm we focus on analyzing dependence on the number of time intervals. The results of our scaling experiments with the Polish model are shown in Fig. (4), where the dashed line show a (rather satisfactory) quadratic match (for dependence of the overall computational time on the number of intervals.) In this case we use the Primal Simplex CPLEX solver at each LP step of our algorithm. (This is LP and not QP, as in the 30 bus model, because the de-fault generation cost is linear in the Polish model of Mathpower.)

We conclude this Section with a number of preliminary, and not yet fully conclusive but calling for further investigation, remarks. First of all, we have observed that developing efficient computational strategy for our linearization algorithm/heuristics became the task which is rather sensitive to the functional form of the generation cost and the choice of variables. If the cost is linear (in generated power) introducing auxiliary (line flow) variables and using the Primal Simplex solver is the winning strategy. However, the same approach in the case of quadratic cost (QP step replacing LP step) leads to slower convergence for the Primal Simplex algorithm while the barrier algorithm fails to converge at all.

Iv-C Gradual investments to resolve congestion

We apply our newly developed algorithm to study effect of the gradual investment, available only within the multi-time period framework, on the overall cost. We study the Polish model in the case of a single scenario with the optimization horizon of one year broken in periods. The (single) loading scenario, chosen to be stressed but still feasible (it is only away from the boundary of the AC OPF infeasibility - see Appendix B for details), stays the same over time. The congestion cost of the initial loading scenario (yet no investments in FACTS) is 17000 . The investment budget is limited to per (one month) time interval. The results of optimal investment generated by our reinforced algorithm are illustrated in Figs. 5,6,7. We observe that only SC devices were installed at lines, of which only two would be overloaded (if the line limits are, first, ignored while solving AC OPF and then checked for the overload). (Lines which are both overloaded, and thus contained in the active set of our cutting plance algorithm, and which are also selected for optimal SC installation are shown blue in Fig. 5. Lines which are shown green were not overloaded but chosen for SC installation. Lines which are shown red were overloaded but were not chosen for SC installation.)

We observe that the optimal installation is gradual. Moreover, all (constrained) available money are spent at each (time interval) decision. Fig. 6 shows how the congestion reduces with time, thus leading to reduction of the operational cost (blue bars) as time progresses. Orange line marks result of the AC OPF before investments start. Red line marks result for the (initial, i.e. before investments) AC OPF with the thermal limits ignored.

Fig. 7 shows that distribution of investments over lines and time period is nontrivial, therefore utilizing the newly available SC-capacities with other operational degrees of freedom.

We also show in Fig. 6 and Fig. 7 the result of optimal investment when the entire year (time horizon) is considered as one time interval. We observe that in the latter case the entire installation budget is used immediately such that the final solutions in the case one and intervals are the same and equal to the available budget. Notice, however, that making one investment upfront for the entire year, as opposed to breaking it into , periods is preferable because the overall (integrated over the year) cost of generation is reduced.

Fig. 5: Snapshot of the final solution (after all the investments are made). Red marks lines from the active set of cutting plane. Thin green marks lines with installed SCs which are not overloaded (initially). Dashed blue marks lines which are both in the active set (overloaded) and chosen for SC installation.
Fig. 6: Dependence of the optimal cost on time (blue bars) provided by our reinforced algorithm for the Polish model in the case of a single scenario (the same for different time intervals) and the year-long horizon split in periods (months). Orange line shows initial ACOPF cost (without investments). Red line shows ACOPF cost with thermal limits ignored. Yellow line corresponds to the optimal solution found in the case when the entire horizon is treated as a single time-interval.
Fig. 7: Installed capacity of the SC devices at corresponding lines shown in percentage of the initial line inductances. This is the case of a single scenario therefore operational values coincide with respective capacities. For example, inductance of line is reduced to .

V Conclusion

In this manuscript, we propose a new optimization framework for algorithmic resolution of placement and sizing of FACTS devices in a large transmission grids. Our algorithmic solution of the problem takes into account non-linear power flow equations. The most important features of the newly developed algorithms include, scalability, allowing to resolve congestion over practical (thousands of buses) size transmission systems, and also the ability to resolve multiple scenarios and over multiple time intervals simultaneously. The optimization can be considered as generalizing standard AC-OPF over multiple scenarios and multiple time intervals with an added cost of installation. This optimization accounts for both installation and operations, thus allowing the installed devices to adjust operationally to a particular scenario within the bounds set by the installed capacity.

Many interesting cases were analyzed experimentally. In all the cases considered the output (optimal placement) is spatially sparse also showing strong non-locality (in the sense that placement of a device may resolve congestion in a distant region of the grid). It is also observed that under highly loaded conditions FACTS devices are beneficial in reducing the total cost of generation. Optimal installation of the devices helps to resolve infeasibilities that are projected to become even more severe in the future.

Main technical achievement reported in this paper is the development of an efficient heuristics for solving the non-linear, non-convex and multi-time-interval optimization. The developed algorithm builds a convergent sequence of convex optimizations with linear constraints. Each constraint is represented explicitly through exact analytical linearization of the original nonlinear constraints (e.g. representing power flows and apparent power line limits) over all the degrees of freedom (including FACTS corrections) around the current operational point for particular loading scenario over particular time interval. In order to represent uncertainty in the projected growth of the system (loads) a custom scenario sampling split over multiple time intervals is introduced. Practicality of our approach for resolving the problem of investment (new installation) planning is illustrated on the IEEE 30-bus model and 2736-bus Polish model. It is evident from the experimental results that the approach is capable of both improving the system’s economy (reduce congestion price and generation cost) and also of resolving feasibility issues by introducing additional degrees of freedom (associated with the newly installed FACTS devices).


The work at LANL was supported by funding from the U.S. Department of Energy’s Office of Electricity as part of the DOE Grid Modernization Initiative.

Appendix A Transmission system modeling

Figure 8 shows a pictorial representation of model of a transmission line.

Power lines:

Fig. 8: model of a transmission line [29].

The model is parametrized by series impedance , total charging susceptance , transformation ratio and shift angle . A transmission line has two ends: sending (often called “from” end) and receiving (often called “to” end).

Explicit expressions for the apparent powers injected at “from” and “to” ends of a line in terms of voltages and phases are:


where and . Note that not only , but also expressions for and are non-symmetric with respect to the change of and indices. The symmetry is broken as the transformer is positioned at “from” side of the line.

Generators and loads:

Generators and loads are modeled as apparent power injections or consumptions:


Generation cost is a quadratic function of active power generation.

FACTS devices:

FACTS devices are described by two variables - first is capacity, second is setting.

Series Compensation (SC) device capacity and setting for the line :


Static Var Compensation (SVC) device capacity and setting for the load node :


Expressions (22) and (23) are used for analytic linearization around current states.

Appendix B Generation of scenarios

The method of scenario generation/sampling is used to include the uncertainty related to system load for the planning period. Standard power system load growth over the time horizon is modeled via the LD curve. The base LD curve is illustrated in Figure 9. In our simulations we usually have just a single base case loading profile for a power case. In order to define where it is situated on the corresponding LD curve we define loading level to represent condition which is from ACOPF infeasibility by uniform load rescaling (we call it top point conditions).

We use the base LD curve, first, to generate LD curvess for consecutive years, re-scaling the base LD curve by the load growth factor of a year. Second, each early LD curve is split into piece-wise-constant parts. (We choose in our experiments.) Finally, each piece of an LD curve is used to generate a scenarios according to a random (thus called sampling) procedure described below. This scheme of scenario generation/sampling models variations in the distribution of loads thus simulating power system behavior during an extended period of time in the future.

Fig. 9: Piece-wise-constant approximation of the LD curve.

We assume (and this assumption is maintained in all of our experimental tests) that each of the generated (sampled) load scenario is ACOPF feasible when the line constraints are ignored. (In other words, we consider the setting when there is enough of generation capacity even for the stressed cases.) Depending on the sampled scenario, 3 situations may arise.

  1. ACOPF is feasible and congestion price is zero (low loading level).

  2. ACOPF is feasible and congestion price is positive (higher loading level representing peak conditions).

  3. ACOPF is infeasible due to either congestion of lines and/or voltage constraints but the system has enough generation capacity. ACOPF without apparent power limits on lines (and without voltage constraints if infeasible) is feasible (overloaded conditions which are possible in the future).

The aim of planning installation of FACTS devices at the right locations with their corresponding capacities is to reduce generation cost for point 2, and to improve or extend feasibility domain of the system for the point 3. Extra years of service can hence be added to the existing grid by making it more flexible, thereby allowing to delay investments into new lines and generators.

B-a Scenarios sampling for each segment

The loading level for a segment is represented by:


Future loading configurations are obtained from the base case by re-scaling all active and reactive loads by uniformly. The resulting vector of loads for a segment is thus given by:


Loading configurations are generated, for each segment and each , through modification of initial . It is done by adding Gaussian correction to each load with zero expected value and a respective standard deviation:


where, is given by:


The choice of parameters used in our experimental test to sample the scenarios is described in Table II.

to — X[c] — X[c] — X[c] — X[c] — & & &
1 & 5,50 & 0,940 & 0,064
2 & 19,50 & 0,845 & 0,041
3 & 25,00 & 0,775 & 0,045
4 & 25,00 & 0,685 & 0,080
5 & 18,80 & 0,590 & 0,068
6 & 6,20 & 0,51 & 0,078

TABLE II: Implementation of the LD curve scheme

Congestion Analysis Correction

If we study a case where for given load configuration standard AC-OPF outputs solution which is not congested, i.e. solution for which each constraint (on line flows or voltages) is satisfied with a margin, then this scenario does not require any FACTS installations. If the whole segment (from the procedure described in the preceding Subsection) is of this “zero-congestion” type, then obviously do not need to generate many samples representing the segment. Instead, we pick one re-scaled base scenario to represent the whole segment.

Appendix C Initializing the generation profile

Generation capacity is assumed to be large enough, i.e. not limiting, for the loading level considered. The initial profile of the generation for each load scenario have to be determined to run the algorithm. The initial generation profile is derived following the following two steps. 1) Re-solving ACOPF with the thermal limits ignored. 2) Setting up proportional response of the generators. The second step includes a) search for the smallest load re-scaling factor lowering the load and thus making the resulting case feasible; b) resolving ACOPF with this new re-scaled loading; c) increasing the generation and load proportionally to match the initial system loading (voltages are kept equal to these provided by the ACOPF solution); finally, d) solving ACPF to obtain generation matching the initial load.


  • [1] N. Alguacil, A. L. Motto, and A. J. Conejo, “Transmission expansion planning: a mixed-integer lp approach,” IEEE Transactions on Power Systems, vol. 18, no. 3, pp. 1070–1077, Aug 2003.
  • [2] R. Bent, G. Toole, and A. Berscheid, “Transmission network expansion planning with complex power flow models,” IEEE Transactions on Power Systems, vol. 27, no. 2, pp. 904–912, 2012.
  • [3] A. H. Fuldner, “Upgrading transmission capacity for wholesale electric power trade,” :, accessed: 2016-11-16.
  • [4] H. Nagarajan, R. Bent, P. Van Hentenryck, S. Backhaus, and E. Yamangil, “Resilient Transmission Grid Design: AC Relaxation vs. DC approximation,” arxiv:1703.05893, 2017.
  • [5] K. Aoki, M. Fan, and A. Nishikori, “Optimal VAR planning by approximation method for recursive mixed-integer linear programming,” IEEE Trans. Power Syst., vol. 3, no. 4, pp. 1741–1747, Nov. 1988.
  • [6] Y. T. Hsiao, C. C. Liu, H. D. Chiang, and Y. L. Chen, “A new approach for optimal VAR sources planning in large scale electric power systems,” IEEE Trans. Power Syst., vol. 8, no. 3, pp. 988–996, Aug. 1993.
  • [7] D. Chattopadhyay, K. Bhattacharya, and J. Parikh, “Optimal reactive power planning and its spot-pricing: an integrated approach,” IEEE Trans. Power Syst., vol. 10, no. 4, pp. 2014–2020, Nov. 1995.
  • [8] F. Dong, B. H. Chowdhury, M. L. Crow, and L. Acar, “Improving voltage stability by reactive power reserve management,” IEEE Trans. Power Syst., vol. 20, no. 1, pp. 338–345, Feb. 2005.
  • [9] G. Y. Yang, G. Hovland, R. Majumder, and Z. Y. Dong, “TCSC allocation based on line flow based equations via mixed-integer programming,” IEEE Trans. Power Syst., vol. 22, no. 4, pp. 2262–2269, Nov. 2007.
  • [10] A. Lashkar Ara, A. Kazemi, and S. A. Nabavi Niaki, “Multiobjective optimal location of FACTS shunt-series controllers for power system operation planning,” IEEE Trans. Power Del., vol. 27, no. 2, pp. 481–490, Dec. 2011.
  • [11] S. J. Galloway, I. M. Elders, G. M. Burt, and B. Sookananta, “Optimal flexible alternative current transmission system device allocation under system fluctuations due to demand and renewable generation,” IET Gen. Trans. Distrib., vol. 4, no. 6, pp. 725–735, June 2010.
  • [12] R. Minguez, F. Milano, R. Zarate-Miano, and A. J. Conejo, “Optimal network placement of SVC devices,” IEEE Trans. Power Syst., vol. 22, no. 4, pp. 1851–1860, Nov. 2007.
  • [13] E. Ghahremani and I. Kamwa, “Optimal placement of multiple-type FACTS devices to maximize power system loadability using a generic graphical user interface,” IEEE Trans. Power Syst., vol. 28, no. 2, pp. 764–778, Aug. 2013.
  • [14] R. S. Wibowo, N. Yorino, M. Eghbal, Y. Zoka, and Y. Sasaki, “FACTS devices allocation with control coordination considering congestion relief and voltage stability,” IEEE Trans. Power Syst., vol. 26, no. 4, pp. 2302–2310, Nov 2011.
  • [15] J. G. Singh, P. Guha Thakurta, and L. Soder, “Load curtailment minimization by optimal placements of SVC/STATCOM,” Int. Trans. Elect. Energy Syst., vol. 25, no. 11, pp. 2769–2780, Sept. 2014.
  • [16]

    E. A. Leonidaki, D. P. Georgiadis, and N. D. Hatziargyriou, “Decision trees for determination of optimal location and rate of series compensation to increase power system loading margin,”

    IEEE Trans. Power Syst., vol. 21, no. 3, pp. 1303–1310, Jul. 2006.
  • [17] F. G. M. Lima, F. D. Galiana, I. Kockar, and J. Munoz, “Phase shifter placement in large-scale systems via mixed integer linear programming,” IEEE Trans. Power Syst., vol. 18, no. 3, pp. 1029–1034, Aug. 2003.
  • [18] S. Rahimzadeh, M. Tavakoli, and B. A. H. Viki, “Simultaneous application of multi-type FACTS devices to the restructured environment: achieving both optimal number and location,” IET Gen. Trans. Distrib., vol. 4, no. 3, pp. 349–362, Mar. 2010.
  • [19] S. Gerbex, R. Cherkaoui, and A. J. Germond, “Optimal location of multi-type FACTS devices in a power system by means of genetic algorithms,” IEEE Trans. Power Syst., vol. 16, no. 3, pp. 537–544, Aug. 2001.
  • [20] P. Paterni, S. Vitet, M. Bena, and A. Yokoyama, “Optimal location of phase shifters in the French network by genetic algorithm,” IEEE Trans. Power Syst., vol. 14, no. 1, pp. 37–42, Feb. 1999.
  • [21] L. Ippolito and P. Siano, “Selection of optimal number and location of thyristor-controlled phase shifters using genetic based algorithms,” IEE Proc. Gen. Trans. Distrib., vol. 151, no. 5, pp. 630–637, Sept. 2004.
  • [22] C. Duan, W. Fang, L. Jiang, and S. Niu, “FACTS devices allocation via sparse optimization,” IEEE Trans. Power Syst., vol. 31, no. 2, pp. 1308–1319, Mar. 2016.
  • [23] V. Frolov, S. Backhaus, and M. Chertkov, “Efficient algorithm for locating and sizing series compensation devices in large power transmission grids: I. model implementation,” New Journal of Physics, vol. 16, no. 10, 2014.
  • [24] ——, “Efficient algorithm for locating and sizing series compensation devices in large power transmission grids: II. solutions and applications,” New Journal of Physics, vol. 16, no. 10, 2014.
  • [25] V. Frolov, P. Guha Thakurta, S. Backhaus, J. Bialek, and M. Chertkov, “Optimal Placement and Sizing of FACTS Devices to Delay Transmission Expansion,” arxiv:, 2016.
  • [26] “ILOG CPLEX optimization studio,”
  • [27] “[online] julia/jump package for steady-state power network optimization,”, accessed: 2017-06-17.
  • [28] “[online] interior point optimizer,”, accessed: 2017-06-18.
  • [29] R. D. Zimmerman, C. E. Murillo-Sanchez, and R. J. Thomas, “MATPOWER: Steady-State Operations, Planning and Analysis Tools for Power Systems Research and Education,” IEEE Trans. Power Syst., vol. 26, no. 1, pp. 12–19, Feb 2011.