An Iterative Approach to Finding Global Solutions of AC Optimal Power Flow Problems
ii An Iterative Approach to Finding Global Solutionsof AC Optimal Power Flow Problems
Ling Zhang,
Student Member, IEEE, and Baosen Zhang,
Member, IEEE,
Abstract —The existence of multiple solutions to AC optimalpower flow (ACOPF) problems has been noted for decades.Existing solvers are generally successful in finding local solutions,which are stationary points but may not be globally optimal.In this paper, we propose a simple iterative approach to findglobally optimal solutions to ACOPF problems. First, we callan existing solver for the ACOPF problem. From the solutionand the associated dual variables, we form a partial Lagrangian.Then we optimize this partial Lagrangian and use its solution as awarm start to call the solver again for the ACOPF problem. Byrepeating this process, we can iteratively improve the solutionquality, moving from local solutions to global ones. We showthe effectiveness our algorithm on standard IEEE networks. Thesimulation results show that our algorithm can escape from localsolutions to achieve global optimums within a few iterations.
I. I
NTRODUCTION
Optimal power flow (OPF) is a fundamental resource al-location problem in power system operations that minimizesthe cost of power generation while satisfying demand. TheACOPF formulation of the problem uses nonlinear power flowequations, resulting in nonlinear and nonconvex optimizationproblems [1]–[3]. The consequence of the nonconvexity ofACOPF we study in this paper is the presence of multiplesolutions. The existence of multiple local solutions of theOPF problem has been well-known for several decades [4].A common assumption is that OPF problems tend to have asingle “practical” solution, and therefore the fact that multiplesolutions can exist do not impact day-to-day operations [5],[6]. However, an increasingly large body of work have pointedto that multiple solutions do occur under reasonable conditionsand cannot easily ruled out [2], [7], [8]. For example, [7] showshow modifications of the standard IEEE benchmarks can leadto each having more than one local solutions.Most ACOPF problems are solved via variations of non-linear optimization algorithms, including Newton-Raphson,sequential programming, interior points and others (see [2],[9], [10] and the references within). These algorithms arein general only able to certify a solution is locally optimal.These solutions satisfy some necessary optimality conditions(e.g., the KKT conditions) but may not globally optimal, thatis, the least cost solution. An open question in the field isto develop algorithms that can find global optimal solutions.In addition, understanding and distinguishing the local andthe global optimal solutions can lead to important theoreticaldiscoveries about the ACOPF problem.
L. Zhang and B. Zhang are with the Department of Electrical andComputer Engineering, University of Washington, WA, 98195 USA e-mail: { lzhang18,zhangbao } @uw.edu The authors are partially supported by NSFgrants ECCS-1942326 and ECCS-2023531 In this paper, we propose a simple algorithm that caneffectively escape from local solutions to seek global optimalsolutions. The process is outlined in Fig. 1. First, we solve theACOPF problem using an existing solver (e.g., IPOPT [11]).From the solution and its associated dual variables, we form apartial Lagrangian. This partial Lagrangian serves to “flatten”the geometry of the optimization problem. We then optimizethis partial Lagrangian, which can lead to a different solution.Using this second solution as a warm start, we again callthe solver for the ACOPF problem. Repeating this iterativeprocess, we can successively improve the solution quality,moving from local solutions to globally optimal ones.Fig. 1: Outline of the solution process.We show the performance of our algorithm on standard 9-bus, 22-bus and 39-bus networks. We show that our algorithmcan escape from local solutions [7] and find the global optimalsolution a single iteration for most instances. We also sketchthe intuition behind the algorithm using a two bus system.Our approach can be seen as a way to provide good warmstarts to nonlinear optimization solvers. Existing approachesalong this line either randomizes [12] or uses a previoussolution as the starting point [13]. The former tends to betime consuming, while the latter does not guarantee theresulting solution would be globally optimal [7]. The work in[14] suggests that solutions can escape local minima if loadundergoes random fluctuations. Our approach can be seen asproviding an explicit and deterministic algorithm to search forglobal solutions by using dual variables at local solutions.II. M
ODEL AND P ROBLEM F ORMULATION
Consider a power system network where n buses are con-nected by m edges. For bus i , let P Gi and Q Gi denote theactive and reactive output of the generated and P Di and Q Di denote the active and reactive load. We use P fij and Q fij todenote the active and reactive power flowing from bus i tobus j . Note that if buses i and j are not connected, then theflows are zero. We use θ ij as a shorthand for θ i − θ j . a r X i v : . [ ee ss . S Y ] F e b i The OPF problem is to minimize the total of active powergeneration costs while satisfying a number of constraints: min V , θ (cid:80) i c i ( P Gi ) (1a)s.t. P Gi = P Di + (cid:80) Nj =1 P fij (1b) Q Gi = Q Di + (cid:80) Nj =1 Q fij (1c) P fij = V i g ij − V i V j ( g ij cos( θ ij ) − b ij sin( θ ij )) (1d) Q fij = V i ˆ b ij − V i V j ( b ij cos( θ ij ) + g ij sin( θ ij )) (1e)V i ≤ V i ≤ ¯ V i (1f)P Gi ≤ P Gi ≤ ¯ P Gi (1g)Q Gi ≤ Q Gi ≤ ¯ Q Gi (1h) ( P fij ) + ( Q fij ) ≤ ( S max ij ) (1i)where ˆ b ij = b ij − . b Cij and b Cij is the line charging suscep-tance. The constraints (1b) and (1c) enforce power balance,(1d) and (1e) are the AC power flow equations, (1f) limits thebus voltage magnitudes, (1g) and (1h) represent the active andreactive limits and (1i) are the line flow limits. The cost at bus i is c i ( · ) and can be linear, quadratic or other functions.A pair of ( V , θ ) that minimizes the objective cost among allfeasible solutions to (1) is called a global solution. In practice,a nonlinear programming (NLP) solver may only return alocal solution, which is a feasible point that satisfy some localoptimality conditions (e.g., KKT). Since a local solution is notnecessarily global, in the following subsection, we proposean iterative approach to find global solutions by alternativelysolving (1) and its partial Lagrangian. Any existing NLP canbe used, and and we use IPOPT in this paper.III. A LGORITHM
Our algorithm starts with a call to a NLP solver withan initial guess, denoted by θ init , V init . For example, thiscan be the standard flat start. Then we assume the solverreturns a feasible solution. At this solution, we record the dualvariables associated with the power balance equations (1b) and(1c), denoted as ¯ µ P and ¯ µ Q . Using these dual variables, weform the following partial Lagrangian by dualizing the powerbalance equations: L ( V , θ , µ P , µ Q ) = (cid:80) i c i P Gi + (cid:80) i µ Pi ( P Di + (cid:80) Nj =1 P fij − P Gi )+ (cid:80) i µ Qi ( Q Di + (cid:80) Nj =1 Q fij − Q Gi ) . We then minimize the partial Lagrangian by solving min V , θ L ( V , θ , µ P , µ Q ) (3)s.t. (1 d ) − (1 i ) . We solve the problem in (3) at the same initial point ( V init , θ init ) that was used to solve the original primal problemin (1). Denote this solution to (3) by ( ¯ V , ¯ θ ) . Then we startthe solver again to solve (1) but with the initial point ( ¯ V , ¯ θ ) .This process can be repeated until the solutions stop changingor up to a predefined number of iterations.It tures out that the initial point ( ¯ V , ¯ θ ) that is given bysolving the partial Lagrangian is often a much better startingpoint than the original choice of ( V init , θ init ) . Therefore, by repeating these steps, we can iteratively improve the quality ofthe starting point and get the global solution. The algorithm issummarized below as Algorithm 1. We illustrate the intuitionbehind this algorithm in the next section using a two busnetwork and present the numerical results in Section V. Algorithm 1: Solving ACOPF iterativelyInputs: θ (0) init , V (0) init
1: At i -th iteration: Initialized at θ ( i ) init , V ( i ) init :2: Call IPOPT solver to solve problem (1), record ( ¯ µ P ( i ) , ¯ µ Q ( i ) )
3: Given ( ¯ µ P ( i ) , ¯ µ Q ( i ) ) , call IPOPT for the partial Lagrangian in (3),record the solutions as (¯ θ ( i ) , ¯ V ( i ) )
4: Call IPOPT for (1) initialized at (¯ θ ( i ) , ¯ V ( i ) ) record solutions (ˆ θ ( i ) , ˆ V ( i ) )
5: If the solution from line 4 corresponds to lowerobjective function value, then update initial points: θ ( i +1) init = ˆ θ ( i ) , V ( i +1) init = ˆ V ( i )
6: Otherwise, stay with the original initial points:7: Repeat the above process until solutions stop changingor reach the maximum number of iterations.IV. G
EOMETRY
Algorithm 1 is successful because of the geometry of thepartial Lagrangian is “better” than the geometry of the originalproblem. We use a two-bus network as an illustrative example.For simplicity, we set both voltage magnitudes to 1 p.u. andoptimize over the angle. Suppose bus 1 is a generator and thereference (slack) bus with linear cost at $1/MW, and bus 2 isthe load bus with angle − θ . The line admittance is g − jb .Given a load of l at bus 2, the ACOPF is: min θ g − g cos( θ ) + b sin( θ ) (4a)s.t. l + g − g cos( θ ) − b sin( θ ) = 0 . (4b)For a feasible load l , there are two solutions to (4b). To makevisualize how an nonlinear solver would approach (4), weconvert it to a penalized unconstrained problem [15]: L ρ = g − g cos( θ ) + b sin( θ ) (5) + ρ/ l + g − g cos( θ ) − b sin( θ )) , where ρ is a (large) penalty parameter. The function L ρ isplotted in Fig. 2 and there are two local minima, with the leftone being global. If a NLP is initialized in the right valley,the suboptimal solution would be found.Now consider the partial Lagrangian for (4), given by L µ = g − g cos( θ ) + b sin( θ ) (6) + µ ( l + g − g cos( θ ) − b sin( θ )) , where µ is the multiplier corresponding to the equality con-straint (4b) at a local solution. The blue line in Fig. 2 plots L µ when dual variable is computed at the right (suboptimal) localsolution. In contrast to L ρ , L µ has a unique minimum thatfalls into the left valley, and can be obtained by minimizing L µ from any initialization. Therefore, if we use the minimumof L µ as an initialization point for the nonlinear solver, wewould reach the global optimal solution. ii Fig. 2: Geometry of the penalized objective functions L ρ and thepartial Lagrangian L µ . The above example illustrates that the Lagrangian has amore advantageous geometry for optimization than the orig-inal primal problems. This is not surprising, since equalityconstraints tend to complicate the optimization landscape [2].The key insight is to observe that it suffices to use the dualvariables at suboptimal local solutions to form the Lagrangian.V. S
IMULATION R ESULTS
In this section we report the simulation results to validatethe effectiveness of our algorithm. The NLP solver used isIPOPT [11] and the convergence tolerance is set to . .We assume IPOPT returns a feasible solution, which may notbe a global optimum. We test our algorithm on different sizesof IEEE networks. For the 9-bus and 22-bus networks, weutilize the local solutions found in [7] as starting points tolaunch IPOPT. Specifically, there are 4 local solutions in the9-bus network, and 2 local solutions in the 22-bus network.For each network, starting from different local solutions, ouralgorithm can achieve global solutions within one iteration.For the 39-bus network, we take a set of initial pointsfor θ and V . To find the global solution, we do an exhaustivesearch within the bounds of each variable. In Fig. 3, weillustrate the improvement of solution quality for the 39-busnetwork. The x-axis represents the number of iterations thatAlgorithm 1 is ran, and y-axis represents the percentage ofglobally optimal solutions after each iteration. Without Algo-rithm 1, less than half of the solutions are globally optimal.One application of Algorithm 1 increases the percentage ofglobally optimal solutions to . After two iterations, onlyfour cases do not reach global optimas. All solutions areglobally optimal after three iterations.VI. C ONCLUSION
In this paper, we propose a simple algorithm to find globallyoptimal solutions to ACOPF problems iteratively. First, wesolve the ACOPF problem using an existing nonlinear solver.From the solution and its associated dual variables, we con-struct a partial Lagrangian. Optimizing this partial Lagrangianleads to a new solution. With this solution as an initial point,we again call the solver for the ACOPF problem. By repeatingthese steps, we can iteratively improve the solution quality,escaping from local solutions to achieve global optimums.We illustrate the intuition behind our algorithm using a two-bus network, which shows that the partial Lagrangian has a
Fig. 3: Using a set of random starting points, 47% of them leads tothe global optimal after a direct call to IPOPT. The fraction of globaloptimal solutions increases to 98%, 99.93% and 100% after one, twoand three iterations of algorithm 1, respectively. flatter optimization landscape compared to the original primalproblem. We validate the effectiveness of our algorithm onstandard 9-bus, 22-bus and 39-bus networks. Regardless of theinitial points, our algorithm always finds the global optimumwithin at most three iterations.R
EFERENCES[1] M. B. Cain, R. P. O’neill, A. Castillo et al. , “History of optimal powerflow and formulations,”
FERC , vol. 1, pp. 1–36, 2012.[2] D. K. Molzahn and I. A. Hiskens, “A survey of relaxations andapproximations of the power flow equations,”
Foundations and Trendsin Electric Energy Systems , vol. 4, no. 1-2, pp. 1–221, 2019.[3] I. A. Hiskens and R. J. Davy, “Exploring the power flow solution spaceboundary,”
IEEE transactions on power systems , vol. 16, no. 3, pp. 389–395, 2001.[4] J. A. Momoh, R. Adapa, and M. El-Hawary, “A review of selectedoptimal power flow literature to 1993. i. nonlinear and quadratic pro-gramming approaches,”
IEEE transactions on power systems , vol. 14,no. 1, pp. 96–104, 1999.[5] J. Momoh, R. Koessler, M. Bond, B. Stott, D. Sun, A. Papalexopoulos,and P. Ristanovic, “Challenges to optimal power flow,”
IEEE Transac-tions on Power systems , vol. 12, no. 1, pp. 444–455, 1997.[6] H. Wei, H. Sasaki, J. Kubokawa, and R. Yokoyama, “An interior pointnonlinear programming for optimal power flow problems with a noveldata structure,”
IEEE Transactions on Power Systems , vol. 13, no. 3, pp.870–877, 1998.[7] W. A. Bukhsh, A. Grothey, K. I. McKinnon, and P. A. Trodden, “Localsolutions of the optimal power flow problem,”
IEEE Transactions onPower Systems , vol. 28, no. 4, pp. 4780–4788, 2013.[8] D. Wu, D. K. Molzahn, B. C. Lesieutre, and K. Dvijotham, “A deter-ministic method to identify multiple local extrema for the ac optimalpower flow problem,”
IEEE Transactions on Power Systems , vol. 33,no. 1, pp. 654–668, 2017.[9] Z. Qiu, G. Deconinck, and R. Belmans, “A literature survey of optimalpower flow problems in the electricity market context,” in
IEEE/PESPower Systems Conference and Exposition , 2009, pp. 1–6.[10] F. Capitanescu, “Critical review of recent advances and further devel-opments needed in ac optimal power flow,”
Electric Power SystemsResearch , vol. 136, pp. 57–68, 2016.[11] A. W¨achter and L. T. Biegler, “On the implementation of an interior-point filter line-search algorithm for large-scale nonlinear programming,”
Mathematical programming , vol. 106, no. 1, pp. 25–57, 2006.[12] N. Costilla-Enriquez, Y. Weng, and B. Zhang, “Combining newton-raphson and stochastic gradient descent for power flow analysis,”
IEEETransactions on Power Systems , vol. 36, no. 1, pp. 514–517, 2021.[13] Y. Tang and S. Low, “Distributed algorithm for time-varying optimalpower flow,” in . IEEE, 2017, pp. 3264–3270.[14] J. Mulvaney-Kemp, S. Fattahi, and J. Lavaei, “Load variation enablesescaping poor solutions of time-varying optimal power flow,” in
PESGM ,2020.[15] D. P. Bertsekas, “Nonlinear programming,”