History of Optimal Power Flow and Formulations • December 2012

Page 1
History of Optimal Power Flow and Formulations • December 2012

Page 2

History of Optimal Power Flow and Formulations 

Optimal Power Flow Paper 1

Mary B. Cain, Richard P. O’Neill, Anya Castillo

mary.cain@ferc.gov; richard.oneill@ferc.gov; anya.castillo@ferc.gov

December, 2012 

Revised August 2013

Abstract: 

The purpose of this paper is to present a literature review of the AC Optimal Power 

Flow (ACOPF) problem and propose areas where the ACOPF could be improved. 

The ACOPF is at the heart of Independent System Operator (ISO) power markets, 

and is solved in some form every year for system planning, every day for day-ahead 

markets, every hour, and even every 5 minutes. It was first formulated in 1962, and 

formulations have changed little over the years. With advances in computing power 

and solution algorithms, we can model more of the constraints and remove 

unnecessary limits and approximations that were previously required to find a 

solution in reasonable time. One example is nonlinear voltage magnitude 

constraints that are modeled as linear thermal proxy constraints. In this paper, we 

refer to the full ACOPF as an ACOPF that simultaneously optimizes real and reactive 

power. Today, 50 years after the problem was formulated, we still do not have a fast, 

robust solution technique for the full ACOPF. Finding a good solution technique for 

the full ACOPF could potentially save tens of billions of dollars annually. Based on 

our literature review, we find that the ACOPF research community lacks a common 

understanding of the problem, its formulation, and objective functions. However, we 

do not claim that this literature review is a complete review—our intent was simply 

to capture the major formulations of the ACOPF. Instead, in this paper, we seek to 

clearly present the ACOPF problem through clear formulations of the problem and 

its parameters. This paper defines and discusses the polar power-voltage, 

rectangular power-voltage, and rectangular current-voltage formulations of the 

ACOPF. Additionally, it discusses the different types of constraints and objective 

functions. This paper lays the groundwork for further research on the convex

approximation of the ACOPF solution space, a survey of solution techniques, and 

computational performance of different formulations. 

Disclaimer: The views presented are the personal views of the authors and not the Federal Energy 

Regulatory Commission or any of its Commissioners.
History of Optimal Power Flow and Formulations • December 2012

Page 3

Table of Contents

1. Introduction.......................................................................................4

2. History of Power System Optimization ....................................7

3. Conventions, Parameters, Sets and Variables........................13

4. Admittance Matrix and AC Power Flow Equations...............16

5. ACOPF Formulations.......................................................................22

6. Literature Review of Formulations ...........................................28

7. Conclusions ........................................................................................32

References 
History of Optimal Power Flow and Formulations • December 2012

Page 4

1. Introduction

The heart of economically efficient and reliable Independent System 

Operator (ISO) power markets is the alternating current optimal power flow 

(ACOPF) problem. This problem is complex economically, electrically and 

computationally. Economically, an efficient market equilibrium requires multi-part 

nonlinear pricing. Electrically, the power flow is alternating current (AC), which 

introduces additional nonlinearities. Computationally, the optimization has 

nonconvexities, including both binary variables and continuous functions, which 

makes the problem difficult to solve. The power system must be able to withstand 

the loss of any generator or transmission element, and the system operator must 

make binary decisions to start up and shut down generation and transmission 

assets in response to system events. For investment planning purposes, the problem 

needs binary investment variables and a multiple year horizon. 

Even 50 years after the problem was first formulated, we still lack a fast and 

robust solution technique for the full ACOPF. We use approximations, 

decompositions and engineering judgment to obtain reasonably acceptable 

solutions to this problem. While superior to their predecessors, today’s 

approximate-solution techniques may unnecessarily cost tens of billions of dollars 

per year. They may also result in environmental harm from unnecessary emissions 

and wasted energy. Using EIA data on wholesale electricity prices and U.S. and 

World energy production, Table 1 gives a range of potential cost savings from a 5% 

increase in market efficiency due to improvements to the ACOPF.(EIA 2012). Small 

increases in efficiency of dispatch are measured in billions of dollars per year. Since 

the usual cost of purchasing and installing new software for an existing ISO market 

is less than $10 million dollars (O’Neill et. al. 2011), the potential benefit/cost ratios 

of better software are in the range of 10 to 1000.
History of Optimal Power Flow and Formulations • December 2012

Page 5

TABLE 1: POTENTIAL COST SAVINGS OF INCREASED EFFICIENCY OF DISPATCH (EIA 2012)

2009 gross 

electricity 

production 

(MWh)

Production cost 

($billion/year) 

assuming 

$30/MWh energy 

price

Savings 

($billion/year) 

assuming 5% 

increase in 

efficiency

Production cost 

($billion/year) 

assuming 

$100/MWh 

energy price

Savings 

($billion/year) 

assuming 5% 

increase in 

efficiency

U.S. 3,724,000 112 6 372 19

World 17,314,000 519 26 1731 87

An ultimate goal of ISO market software, and a topic of future research, is the 

security-constrained, self-healing (corrective switching) AC optimal power flow 

with unit commitment over the optimal network. The optimal network is flexible, 

with assets that have time-varying dynamic ratings reflecting the asset capability

under varying operating conditions. The optimal network is also optimally 

configured – opening or closing transmission lines becomes a decision variable, or 

control action, rather than an input to the problem, or state. When possible, the 

security constraints are corrective rather than preventive. With preventive security 

constraints, the system is operated conservatively to survive loss of any 

transmission element or generator. In contrast, corrective constraints reconfigure 

the system with fast-acting equipment such as special protection systems or 

remedial action schemes immediately following loss of a generator or transmission 

element, allowing the system to be reliably used closer to its limits. This problem 

must be solved weekly in 8 hours, daily in 2 hours, hourly in 15 minutes, each five 

minutes in 1 minute and for self-healing post-contingency in 30 seconds. Currently, 

the problem is solved through varying levels of approximation, depending on 

application and time scale, but with increases in computing power it may be 

possible to reduce the number of approximations and take advantage of parallel 

computing.

Today, the computational challenge is to consistently find a global optimal 

solution with speeds up to three to five orders of magnitude faster than existing 

solvers. There is some promising recent evidence that this could be a reality in five 

to ten years. For example, in the last two decades mixed-integer programming (MIP) 

has achieved speed improvements of 107; that is, problems that would have taken 
History of Optimal Power Flow and Formulations • December 2012

Page 6

10 years in 1990 can be solved in one minute today. As a consequence, MIP is 

replacing other approaches in ISO markets. Implementation of MIP into the dayahead and real-time markets, with the Commission’s encouragement, has saved 

American electricity market participants over one-half billion dollars per year 

(FERC 2011). More will be saved as all ISOs implement MIP and the new 

formulations it permits in the next several years. 

Due to idiosyncrasies in design, current software oversimplifies the problem 

in different ways, and requires operator intervention to address real-time problems 

that do not show up in models. This operator intervention unnecessarily alters 

settlement prices and produces suboptimal solutions. The Joint Board on Economic 

Dispatch for the Northeast Region stated in 2006 that improved modeling of system 

constraints such as voltage and stability constraints would result in more precise 

dispatches and better market signals, but that the switch to AC-based software 

would increase the time to run a single scenario from minutes to over an hour, 

making use of ACOPF impractical, even for the day-ahead market (FERC 2006). One 

example is the Midwest Independent System Operator (MISO), where operators 

have to commit resources before the unit commitment and economic dispatch 

software models are run to address local voltage issues that MISO has had difficulty 

modeling in its market software (FERC 2012). PJM Interconnection (PJM) employs 

an approach, called Perfect Dispatch, that ex-post solves the real-time market 

problem with perfect information (PJM 2012). The Perfect Dispatch solution is used 

to train operators, where they can compare the “perfect dispatch,” which is based on 

“perfect” after-the-fact information to the actual dispatch, which is based on the 

information available at the time. ISO models solve proxies or estimates for reactive 

power and voltage constraints, where they calculate linear thermal constraints to 

approximate quadratic voltage magnitude constraints. The details of transmission 

constraint modeling and transmission pricing have been neglected, but need to be 

considered to improve the accuracy of ACOPF calculations. Transmission 

constraints can be modeled in terms of current, real power, apparent power, voltage 

magnitude differences, or angle differences. The choice of constraint depends on the 

type of model, data availability, and physical limit (voltage, stability, or thermal 
History of Optimal Power Flow and Formulations • December 2012

Page 7

limit). Surrogate constraints can be calculated based on the line flow equations, but 

these calculations have inherent assumptions. One example is the Arizona-Southern 

California outage in 2011, where some line limits were modeled and monitored as 

real power transfer limits while others were modeled as current transfer limits 

(FERC/NERC 2012). This paper seeks to better understand the ACOPF problem 

through clear formulations of the problem, theoretical properties of the problem 

and its parameters, approximations to the nonlinear functions that are necessary to 

make the problem solvable, and to produce computational results from large and 

small test problems using various solvers and starting points. Discrete variables 

such as equipment states, generator commitments, and transmission switching 

further complicate the ACOPF, but we do not discuss these in this paper. With the 

increased measurements and controls inherent in smart grid upgrades, the potential 

savings are greater, although the problem may become more complex with more 

discrete devices to model. 

In the rest of the paper, we provide a brief history of power system 

optimization, present notation and nomenclature, formulate the admittance matrix 

and power flow equations, formulate constraints, present different formulations of 

the ACOPF, and present a literature review of ACOPF formulations.

2. History of Power System Optimization

Power system optimization has evolved with developments in computing 

and optimization theory. In the first half of the 20th century, the optimal power flow 

problem was “solved” by experienced engineers and operators using judgment, 

rules of thumb, and primitive tools, including analog network analyzers and 

specialized slide rules. Gradually, computational aids were introduced to assist the 

intuition of operator experience. The optimal power flow problem was first 

formulated in the 1960’s (Carpentier 1962), but has proven to be a very difficult 

problem to solve. Linear solvers are widely available for linearized versions of the 

optimal power flow problem, but nonlinear solvers cannot guarantee a global 

optimum, are not robust, and do not solve fast enough. In each electricity control 

room, the optimal power flow problem or an approximation must be solved many 

times a day, as often as every 5 minutes. 
History of Optimal Power Flow and Formulations • December 2012

Page 8

There are three types of problems commonly referred to in power system 

literature: power flow (load flow), economic dispatch, and optimal power flow. 

Three other classes of power system optimization, specifically unit commitment, 

optimal topology, and long-term planning, involve binary and integer variables, and 

are not discussed in this paper; but combined with the insights on formulations in 

this paper, could be promising areas for future research. 

Table 2 compares the major characteristics of the power flow, economic 

dispatch, and optimal power flow problems. The power flow or load flow refers to 

the generation, load, and transmission network equations. Power flow methods find 

a mathematically but not necessarily physically feasible or optimal solution. The 

power flow equations themselves do not take account of limitations on generator 

reactive power limits or transmission line limits, but these constraints can be 

programmed into many power flow solvers.

A second type of problem, economic dispatch, describes a variety of 

formulations to determine the least-cost generation dispatch to serve a given load

with a reserve margin, but these formulations simplify or sometimes altogether 

ignore power flow constraints. 

A third type of problem, the optimal power flow, finds the optimal solution to 

an objective function subject to the power flow constraints and other operational 

constraints, such as generator minimum output constraints, transmission stability 

and voltage constraints, and limits on switching mechanical equipment. Optimal 

power flow is sometimes referred to as security-constrained economic dispatch 

(SCED); most implementations of SCED include only thermal limits, and proxies for 

voltage limits. There are a variety of formulations with different constraints, 

different objective functions, and different solution methods that have been labeled 

optimal power flow; these are discussed in the formulations section later in this 

paper. Formulations that use the exact AC power flow equations are known as 

“ACOPF.” Simpler versions, known as DCOPF, assume all voltage magnitudes are 

fixed and all voltage angles are close to zero. DC stands for direct current, but is a bit 

of a misnomer; a DCOPF is a linearized form of a full alternating current network 

(ACOPF) and not a power flow solution for a direct current network. We use the 
History of Optimal Power Flow and Formulations • December 2012

Page 9

general term OPF to include both ACOPF and DCOPF. The ACOPF is often solved 

through decoupling, which takes advantage of the structure of the problem, where 

real power (P) and voltage angle (θ) are tightly coupled and voltage magnitude (V) 

and reactive power (Q) are tightly coupled, but the P-θ and V-Q problems are 

weakly coupled due to the assumptions that the phase angle differences between 

adjacent buses are rather small, and high-voltage transmission networks have much 

higher reactance compared to resistance. The decoupled OPF divides the ACOPF 

into two linear subproblems, one with power and voltage angle and another with 

voltage magnitude and reactive power. In this paper, we use the term ACOPF to 

refer to the full ACOPF that simultaneously optimizes real and reactive power, and 

decoupled OPF to refer to the decoupled problems that separately optimize real and 

reactive power and iterate between the two to reach an optimal solution.
History of Optimal Power Flow and Formulations • December 2012

Page 

10

TABLE 2: MAJOR TYPES OF POWER SYSTEM PROBLEMS

General 

problem 

type

Problem name Includes 

voltage angle 

constraints?

Includes 

bus voltage 

magnitude 

constraints?

Includes 

transmission 

constraints?

Includes 

losses?

Assumptions Includes 

generator 

costs?

Includes 

contingency 

constraints?

OPF ACOPF, or Full 

ACOPF

Yes Yes Yes Yes Yes No

OPF DCOPF No No; all 

voltage 

magnitudes 

fixed

Yes Maybe Voltage magnitudes are 

constant

Yes No

OPF Decoupled OPF Yes Yes Yes Yes Power-voltage angle are 

independent of voltage 

magnitude-reactive 

power

Yes No

OPF SecurityConstrained 

Economic 

Dispatch (SCED)

Yes No Yes Yes Voltage magnitudes are 

constant

Yes Yes

Power 

flow

Power Flow, or 

Load Flow

No, but can be 

added

Yes No, but can 

be added

Yes No No

Economic 

dispatch

Economic 

Dispatch

No No No Depends No transmission 

constraints

Yes No

OPF Security 

Constrained OPF 

(SCOPF)

Yes Depends Yes Yes Depends Yes Yes
History of Optimal Power Flow and Formulations • December 2012

Page 11

We now discuss early research of the three types of problems in power 

system optimization: economic dispatch, power flow, and optimal power flow. 

As early as the 1930’s, the economic dispatch problem was solved by hand or 

specially-developed slide rule using the principle of equal incremental loading, 

taking as long as 8 hours to complete (Happ 1977). Early computations of economic 

dispatch were slow. Kirchmayer estimated that it would take 10 minutes of 

computational time to produce the schedules for a 10 generator system at a given 

system price (Kirchmayer 1958). In contrast, RTOs today solve systems of hundreds 

of generators in a matter of seconds. In the survey of economic dispatch methods up 

through the 1970’s, Happ provides an overview of the evolution of economic 

dispatch formulations and different ways to account for losses. 

Prior to digital computers, as early as 1929, the power flow problem was 

solved with analog network analyzers that simulated power systems (Sasson 1967). 

Ward and Hale published the first automated digital solution to the power flow 

problem in 1956 (Ward 1956). Sasson and Jaimes provide a survey and comparison 

of early load flow solution methods, which are various iterative methods based on 

the nodal admittance matrix (Y matrix) or its inverse, the nodal impedance matrix 

(Z matrix) (Sasson 1967). Early researchers, including Carpentier, used the GaussSeidel method. The Newton-Raphson method became the commonly used solution 

method during the 1960’s (Peschon et. al. 1968), after Tinney and others developed 

sparsity techniques to take advantage of the structure of the admittance matrix in 

the OPF problem. The admittance matrix is sparse, meaning it has many zero 

elements; this is because power system networks are not densely connected. 

Sparsity techniques can be used to reduce data storage and increase computation 

speed (Stott 1974). 

Early research on OPF used classical Lagrangian techniques for the 

optimality conditions, but neglected bounds on variables (Squires 1961). In 1962, 

Carpentier published the optimality conditions for an OPF, including variable 

bounds, based on the Kuhn-Tucker conditions; this is generally considered the first 

publication of a fully formulated OPF (Carpentier 1962). Carpentier assumes that 

the applicable functions display “suitable convexity” for the Kuhn-Tucker (now 
History of Optimal Power Flow and Formulations • December 2012

Page 12

referred to as the Karush-Kuh-Tucker or KKT) conditions to apply (Carpentier 

1962). Given the structure of the power flow equations, this may be a big 

assumption (Hiskens 2001 and Schecter 2012). Carpentier includes the full AC 

power flow equations, generator real and reactive power constraints, bus voltage 

magnitude constraints, and bus voltage angle difference constraints for buses 

connected by transmission elements. 

Huneault and Galliana provide an extensive survey of optimal power flow 

literature up to 1991, surveying over 300 articles and citing 163 (Huneault 1991). 

They conclude, “The history of optimal power flow (OPF) research can be 

characterized as the application of increasingly powerful optimization tools to a 

problem which basically has been well-defined since the early 1960’s.” The paper 

outlines the evolution of OPF literature, grouped by solution method. The solution 

methods include various forms of gradient methods, linear programming, quadratic 

programming, and penalty methods. The authors conclude that “commercially 

available OPF algorithms all satisfy the full nonlinear load flow model and a full set 

of bounds on variables.” The authors further conclude that the OPF remains a 

difficult mathematical problem. The present algorithms cannot compute quickly 

enough, and are prone to serious ill-conditioning and convergence problems.

Another area of research, security-constrained OPF, accounts for

transmission contingency constraints and poses additional computational 

challenges (Carpentier 1979, Stott 1987). Our discussion in this paper focuses on 

ACOPF. Future research could extend the formulations to include contingency 

constraints that are required to maintain the system after an outage. This 

formulation increases the size of the problem formulation by a factor equal to the 

number of contingencies studied.

Researchers have identified challenges to solving the OPF, including 

modeling discrete variables, local minima, lack of uniform problem definition, 

solution reliability and computing time. Some of these have been solved: both 

Tinney et al. and Momoh et al. discussed the challenges in modeling discrete 

variables in OPF solutions (Tinney 1988), (Momoh 1997). Today, with advances in 

mixed integer programming (MIP), discrete variables can be modeled and included 
History of Optimal Power Flow and Formulations • December 2012

Page 13

in ACOPF solutions. Other challenges persist today: Koessler states that the “lack of 

uniformity in usage and definition” has been a challenge to users and developers in 

OPF, and specifically discusses local minima, which indicate that the problem is 

nonconvex (Momoh 1997). Huneault and Galliana found that algorithms available in 

1991 could not compute OPF solutions quickly and reliably enough, and that the 

OPF, like many nonlinear problems, is prone to ill-conditioning and difficult 

convergence (Huneault 1991).

3. Conventions, Parameters, Sets and Variables

Notation and Nomenclature 

When n and m are subscripts, they index buses; k indexes the transmission 

elements. When j is not a superscript, j = (-1)1/2; i is the complex current. When j is a 

superscript, it is the ‘imaginary’ part of a complex number. Matrices and vectors are 

represented with upper case letters. Scalars and complex numbers are in lower case

letters. For column vectors A and B of length n, where ak and bk are the kth

components of A and B respectively, the Hadamard product ‘∙’ is defined so that A∙B 

= C, where C is a column vector also of length n, with kth component ck = akbk.

The complex conjugate operator is * (superscript) and * (no superscript) is an 

optimal solution.

We assume balanced, three-phase, steady-state conditions. All variables are 

associated with a single-line model of a balanced, three-phase system. A common 

practice in power system modeling is the per-unit (p.u.) representation, where base 

quantities for voltage magnitude, current, power, and impedance (or admittance) 

are used to normalize quantities in a network with multiple voltage levels. Such 

normalization is a convenience. We use the convention that an injection occurs 

when the real part of the complex number is positive and a withdrawal occurs when 

the real part of the complex number is negative.

The topology of the network consists of locations known as buses or nodes 

and transmission elements connecting paired buses. The network is an undirected 

planar graph. 
History of Optimal Power Flow and Formulations • December 2012

Page 14

Indices and Sets 

n, m are bus (node) indices; n, m ε {1, ..., N} where N is the number of buses. (m is an 

alias for n)

k is a three-phase transmission element with terminal buses n and m. 

k ε {1, ..., K} where K is the number of transmission elements; k counts from 1 to the 

total number of transmission elements, and does not start over for each bus pair nm.

K’ is the set of connected bus pairs nm (|K’| ≤| K|).

Unless otherwise stated, summations (∑) are over the full set of indices. 

Variables

pn is the real power injection (positive) or withdrawal (negative) at bus n

qn is the reactive power injection or withdrawal at bus n

sn = pn + jqn is the net complex power injection or withdrawal at bus n

We distinguish between the real, reactive, or complex power injected at a specific 

bus (pn, qn, and sn) and the real, reactive, or complex power flowing in a transmission 

element between two buses:

pnmk is the real power flow from bus n to bus m on transmission element k

qnmk is the reactive power flow from bus n to bus m on transmission element k

snmk is the apparent complex power flow from bus n on transmission element k. snmk

= sr

nmk + jsjnmk = pnmk + j qnmk

θn is the voltage angle at bus n 

θnm = θn - θm is the voltage angle difference from bus n to bus m

θ – δ is the power angle.

i is the current (complex phasor); we distinguish between current injected at a 

specific bus and current flowing in a transmission element between two buses:

in is the current (complex phasor) injection (positive) or withdrawal (negative) at 

bus n where in = ir

n + jijn

inmk is the current (complex phasor) flow in transmission element k at bus n (to bus 

m). inmk = ir

nmk + jijnmk

vn is the complex voltage at bus n. vn = vr

n + jvjn
History of Optimal Power Flow and Formulations • December 2012

Page 15

ynmk is the complex admittance on transmission element k connecting bus n and bus 

m (If buses n and m are not connected directly, ynmk= 0.); yn0 is the self-admittance 

(to ground) at bus n. 

V = (v1, ..., vN)T is the complex vector of bus voltages; V = Vr + jVj

I = (i1, ..., iN)T is the complex vector of bus current injections; I = Ir + jIj

P = (p1, ..., pN)T is the vector of real power injections

Q = (q1, ..., qN)T is the vector of reactive power injections

G is the N-by-N conductance matrix

B is the N-by-N susceptance matrix

Note that elements of G and B will be constant for passive transmission elements 

such as transmission lines, but can be variable when active transmission elements 

such as phase shifting transformers, switched capacitors/reactors, or power 

electronic flexible AC transmission system (FACTS) devices are included.

Y = G + jB is the N-by-N complex admittance matrix

gnm, bnm, and ynm represent elements of the G, B, and Y matrices respectively.

Functions and Transformations

Re( ) is the real part of a complex number, for example, Re(ir

n + jijn) = irn

Im( ) is the real part of a complex number, for example, Im(ir

n + jijn) = ijn

| | is the magnitude of a complex number, for example, |vn| = [(vr

n)2+(vjn)2]1/2

abs( ) is the absolute value function.

The transformation from rectangular to polar coordinates for complex voltage is:

vr

n = |vn|cos(θn)

vj

n = |vn|sin(θn)

(vr

n)2 + (vjn)2 = [|vn|sin(θn)]2 + [|vn|cos(θn)]2 = |vn|2[sin(θn)2 + cos(θn)2] = |vn|2

We drop the bus index n and let θ be the voltage angle and δ be the current angle. 

For real power,

p = vrir + vj

ij = |v|cosθ|i|cosδ + |v|sinθ|i|sinδ = |v||i|[cosθcosδ + sinθsinδ]

 = |v||i|(0.5[cos(θ-δ)+ cos(θ+δ)] + 0.5[cos(θ-δ)- cos(θ+δ)] )

 = |v||i|cos(θ-δ)

For reactive power, 

q = vjir - vrij = |v|sinθ|i|cosδ - |v|cosθ|i|sinδ = |v||i|[sinθcosδ - cosθsinδ]
History of Optimal Power Flow and Formulations • December 2012

Page 16

 = |v||i|.5[sin(θ + δ) + sin(θ - δ)] - |v||i|.5[sin(θ + δ) - sin(θ - δ)]

q = |v||i|sin(θ - δ)

θ – δ is the power angle.

Parameters

rnmk or rk is the resistance of transmission element k. 

xnmk or xk is the reactance of transmission element k. 

smaxk is the thermal limit on apparent power over transmission element k at both 

terminal buses. 

θminnm, θmaxnm are the maximum and minimum voltage angle differences between n

and m

pminn, pmaxn are the maximum and minimum real power for generator n

qminn, qmaxn are the maximum and minimum reactive power for generator n

C1 = (c11, ..., c1N)T and C2 = (c21, ..., c2N)T are vectors of linear and quadratic objective 

function cost coefficients respectively.

4. Admittance Matrix and AC Power Flow Equations

In this section, we develop the admittance matrix and the current-voltage

flow equations (IV equations), which are a different formulation of the commonly

used power flow equations. In the following sections, we develop the additional 

constraints that bound the solutions. 

We define the conductance (G), susceptance (B) and admittance (Y) matrices, 

with elements gnm, bnm, and ynm respectively, and Y = G + jB. We start with a simple 

admittance matrix defined by resistance, r, and reactance, x. We assume shunt 

susceptance is negligible. The elements of G, B and Y matrices are derived as follows: 

gnmk = rnmk/(rnmk2 + xnmk2) for n ≠ m

bnmk =-xnmk/(rnmk2 + xnmk2) for n ≠ m

ynmk = gnmk +jbnmk for n ≠ m

ynmk = 0 for n = m

ynm = ∑k ynmk for n ≠ m

ynn = yn0 -∑n ≠ m ynm

Transformers. The admittance matrix above does not include transformer 

parameters. For an ideal in-phase transformer (assuming zero resistance in 
History of Optimal Power Flow and Formulations • December 2012

Page 17

transformer windings, no leakage flux, and no hysteresis loss), the ideal voltage 

magnitude (turns ratio) is anmk =|vm|/|vn| and θn = θm, where n is the primary side 

and m is the secondary side of the transformer. Since θn = θm, 

anmk =|vm|/|vn| = vm/vn = -inm/imn

The current-voltage (IV) equations for ideal transformer k between buses n and m

are: 

inmk = anmk2ynmkvn - anmkynmkvm

imnk = -anmkynmkvn + ynmkvm

For the phase shifting transformer (PAR) with a phase angle shift of φ,

vm/vn = tnmk = anmkejφ

inm/imn = tnmk* = -anmke-jφ

The current-voltage (IV) equations for the phase shifting transformer k between 

buses n and m are: 

inmk = anmk2ynmkvn - tnmk*ynmkvm

imnk = -tnmkynmkvn + ynmkvm

Admittance Matrix. If there are no transformers or FACTS devices, G is positive 

semidefinite and B is negative semidefinite. A matrix where ynn ≥ abs(∑m ynm) is 

called diagonally dominant and strictly diagonally dominant if ynn > abs(∑m ynm). 

If there are no transformers and yn0 = 0, G and B are weighted Laplacian 

matrices of the undirected weighted graph that describes the transmission network. 

Much is known about the weighted Laplacian matrices. Y is a complex weighted 

Laplacian matrix. The admittance matrix is Y = G+jB, where G and B are real 

symmetric diagonally dominant matrices. A symmetric diagonally dominant matrix 

has a symmetric factorization, for example, B = UUT where each column of U has at 

most two non-zeros and the non-zeroes have the same absolute value. 

For large problems, the admittance matrix, Y= G+jB, is usually sparse. The 

density of both G and B is (N+2K’)/N2 where K’ is the number of off-diagonal nonzero entries (the aggregate of multiple transmission elements between adjacent 

buses) and N is the number of buses. For example, in a topology with 1000 buses 

and 1500 transmission elements, G and B would have a density of 

(1000+3000)/10002 = .004. The lowest density for a connected network is the 
History of Optimal Power Flow and Formulations • December 2012

Page 18

spanning tree. It has N-1 transmission elements and the density is (N+2(N-1))/N2. 

For large sparse systems, (N+2(N-1))/N2 ≈ 3/N. 

Transformers and FACTS devices change the structure of the Y matrix. If 

there are transformers and FACTS devices, let 

ynmk if no transformer

ynmk = { anmk2ynmk if an ideal transformer

tnmk*ynmk, or -tnmkynmkvn if a phase shifting transformer

as appropriate off-diagonal element, then ynn = yn0 + ∑k,m ynmk, ynm = ∑k ynmk, and Y is 

the matrix [ynm]. If there are only ideal in-phase transformers, the Y matrix is 

symmetric. If there are phase shifting transformers, the symmetry of the Y matrix is 

lost. 

AC Power Flow Equations

Kirchhoff’s Current Law. Kirchhoff’s current law requires that the sum of the 

currents injected and withdrawn at bus n equal zero: 

in = ∑k inmk (1)

If we define current to ground to be yn0(vn – v0) and v0 = 0, we have: 

in = ∑k ynmk(vn - vm) + yn0vn (2)

inmk = ynmk(vn - vm) = gnmk(vr

n - vrm)-bnmk(vjn - vjm) + j(bnmk(vrn - vrm)+gnmk(vjn - vjm))

ir

nmk = gnmk(vrn - vrm) - bnmk(vjn - vjm)

ij

nmk = bnmk(vrn - vrm) + gnmk(vjn - vjm)

Current is a linear function of voltage. Rearranging,

in = vn(yn0 + ∑k ynmk) - ∑k ynmkvm (3)

In matrix notation, the IV flow equations in terms of current (I) and voltage (V) in 

(3) are

I = YV = (G + jB)(Vr + jVj) = GVr - BVj + j(BVr + GVj) (4)

where Ir = GVr - BVj and Ij = BVr + GVj

In another matrix format, (4) is

I = (Ir, Ij ) = Y(Vr, Vj)T or

I = (Ir, Ij ) = G -B Vr

 where Y = 

G -B

B G Vj B G
History of Optimal Power Flow and Formulations • December 2012

Page 19

If a and φ are constant, the I =YV equations are linear. If not, the linearity is lost 

since some elements of the Y matrix will be functions of V.

Power Flow Equations. The traditional power-voltage power flow equations defined 

in terms of real power (P), reactive power (Q) and voltage (V) are

S = P + jQ = diag(V)I* = diag(V)[YV]* = diag(V)Y*V* (5)

The power injections are 

S = V•I* = (Vr + jVj )•(Ir - jIj ) = (Vr•Ir + Vj•Ij ) + j(Vj•Ir - Vr•Ij) (6)

where 

P = Vr•Ir + Vj•Ij (7)

Q = Vj•Ir - Vr•Ij (8)

The power-voltage power flow equations (5) and (6) are quadratic. The IV flow 

equations (4) are linear.

Constraints. First, we introduce the physical constraints of generators, load, and 

transmission. 

Generator and Load Constraints. The lower and upper bound constraints for 

generation (injection) and load (withdrawal) are:

Pmin ≤ P ≤ Pmax (9)

Qmin ≤ Q ≤ Qmax (10)

In terms of V and I, the injection constraints are:

Vr•Ir + Vj

•Ij ≤ Pmax (11)

Pmin ≤ Vr•Ir + Vj

•Ij (12)

Vj

•Ir - Vr•Ij ≤ Qmax (13)

Qmin ≤ Vj•Ir - Vr•Ij (14)

Inequalities (11)-(14) along with other thermal constraints on equipment 

enforced at each generator bus constitute a four-dimensional reactive capability 

curve, also known as a “D-curve’ since it is shaped like the capital letter D, in the PQ

space. Additional D-curves defining the tradeoff between real and reactive power 

constitute a convex set and can be easily linearized (FERC 2005). Equations (11)-

(14) are nonconvex quadratic constraints. Since here we model a single period, 

ramp rates are unnecessary.
History of Optimal Power Flow and Formulations • December 2012

Page 20

Voltage Magnitude Constraints. The two constraints that limit the voltage 

magnitude in rectangular coordinates at each bus m are 

(vr

m)2 +(vjm)2 ≤ (vmaxm)2 (15)

(vminm)2 ≤ (vr

m)2 +(vjm)2 (16) 

Again, each nonlinear inequality involves only the voltage magnitudes at bus m. In 

matrix terms, the voltage magnitude constraints are: 

Vr•Vr + Vj•Vj ≤ (Vmax)2 (17)

(Vmin)2 ≤ Vr •Vr + Vj•Vj (18) 

Vmin and Vmax are determined by system studies. The voltage magnitude bounds are 

generally in the range, [.95, 1.05] per unit. High voltages are often constrained by 

the capabilities of the circuit breakers. Low voltage magnitude constraints can be 

due to operating requirements of motors or generators. 

Line Flow Thermal Constraints. Smaxk is a thermal transmission limit on k based on 

the temperature sensitivity of the conductor and supporting material in the 

transmission line and transmission elements. Transmission assets generally have 

three thermal ratings: steady-state, 4-hour and 30-minute. These ratings vary with 

ambient weather. The apparent power at bus n on transmission element k to bus m

is:

snmk = vninmk* = vny*nmk(vn - vm)*. = vny*nmkv*n - vny*nmkv*m)

The thermal limit on snmk is

(sr

nmk)2+ (sjnmk)2 = |snmk|2 ≤ (smaxk)2 (19)

These constraints are quadratic in sr

nmk and sjnmk and quartic in vrn, vjn, vrm, vjm. Since 

vn = vr

n + jvjn and ynmk = gnmk + jbnmk,

vny*nmkv*m = (vr

n + jvjn )(gnmk + jbnmk )(vrm + jvjm)

Expanding, we obtain 

vny*nmkv*m =(vr

n + jvjn )[gnmkvrm - bnmkvjm - j(gnmkvjm + bnmkvrm)]

Expanding again, we obtain 

vny*nmkv*n = gnmk(vr

nvrm+vjnvjm)+bnmk(vjnvrm-vrnvjm)

 +j[gnmk(vj

nvrm-vrnvjm) - bnmk(vrnvrm+vjnvjm)] (20)
History of Optimal Power Flow and Formulations • December 2012

Page 21

In matrix notation,

Re(vny*nmkv*m) = [vr

n, vjn] gnmk -bnmk vr

m

bnmk gnmk vj

m

Im(vny*nmkv*m) = [vr

n, vjn] -bnmk -gnmk vr

m

gnmk -bnmk vj

m

Similarly, vny*nmkv*n = (vr

n + jvjn )(gnn - jbnn)(vrn - jvjn)

Expanding, we obtain

= (vr

n + jvjn)[gnnvrn - bnnvjn - j(gnnvjn + bnnvrn)] 

Expanding and collecting terms,

= gnn(vr

nvrn+vjnvjn)-jbnn(vrnvrn+vjnvjn)] (21)

In matrix notation,

Re(vny*nmkv*n) = [vr

n, vjn] gnn 0 vr

n

0 gnn vj

n

Im(vny*nmkv*n) = [vr

n, vjn] -bnn 0 vr

n

0 -bnn vj

n

Inequality (19) becomes a quadratic constraint.

Line Flow Constraints as Current Limitations. As current increases, lines sag and 

equipment may be damaged by overheating. The constraints that limit the current 

magnitude in rectangular coordinates at each bus n on k are 

(ir

nmk)2 +(ijnmk)2 ≤ (imaxnmk)2 (23)

Again, the nonlinearities are convex quadratic and isolated to the complex current at 

the bus. Generally, the maximum currents, imaxnmk, are determined by material 

science properties of conductors and transmission equipment, or as a result of 

system stability studies.

Line Flow Constraints as Voltage Angle Constraints. The power flowing over an 

AC line is approximately proportional to the sine of the voltage phase angle 

difference at the receiving and transmitting ends. For stability reasons, the voltage 

angle difference for terminal buses n and m connected by transmission element k 

can be constrained as follows: 

θminnm ≤ θn - θm ≤ θmaxnm (24)
History of Optimal Power Flow and Formulations • December 2012

Page 22

In the rectangular formulation, the arctan function appears in some constraints. 

θminnm ≤ arctan(vj

n/vrn) - arctan(vjm/vrm) ≤ θmaxnm (25)

The theoretical steady-state stability limit for power transfer between two 

buses across a lossless line is 90 degrees. If this limit were exceeded, synchronous 

machines at one end of the line would lose synchronism with the other end of the 

line. In addition, transient stability and relay limits on reclosing lines constrain 

voltage angle differences. The angle constraints used in the ACOPF should be the 

smallest of these angle constraints, which depend on the equipment installed and 

configuration. However, many test cases do not include any voltage angle or line 

flow constraints. In general, system engineers design and operate the system 

comfortably below the voltage angle limit to allow time to respond if the voltage 

angle difference across any line approaches its limit. 

5. ACOPF Formulations

We begin with a discussion of objective functions, then a note on bus types, 

and finally discuss different formulations of the ACOPF. The formulations of the 

ACOPF presented here include all the constraints, but may take different approaches 

to modeling the constraints. As discussed above, current, voltage magnitude, and 

voltage angle constraints can be calculated that are surrogates for each other. We 

discuss constraints further in (O’Neill 2012). 

Objective Function. Various authors formulate the ACOPF with different objective 

functions. They include minimizing generation costs, maximizing market surplus, 

minimizing losses, minimizing generation (equivalent to minimizing losses), and 

maximizing transfers. Without demand functions, minimizing generation costs and 

maximizing market surplus are equivalent. 

A full ACOPF that accurately models all constraints and controls with an 

objective function of minimizing cost would inherently meet the objectives of 

minimizing generator fuel costs, minimizing generation output, minimizing losses,

minimizing load shedding, and minimizing control actions. 

When it is not feasible to run a full ACOPF due to time constraints, computing 

power, or lack of a robust solution algorithm, a common substitute is to decouple 

the problem and iterate between a DCOPF that minimizes costs by varying real 
History of Optimal Power Flow and Formulations • December 2012

Page 23

power, then fix the generator outputs from the DCOPF and run an ACOPF that 

minimizes losses by varying reactive power of generators, capacitors, etc. For 

economically dispatching resources in an ACOPF that fully models voltage and 

stability constraints, minimizing cost is the correct objective function; objective

functions of minimizing losses, minimizing generation, and maximizing transfers for 

an ACOPF are inconsistent with economic principles, and result in sub-optimal 

dispatch. We do not discuss the details of decoupled OPF here, but save it for a 

future review of solution algorithms.

Stott et al. discuss badly-posed problems when an OPF formulation does not 

adhere to the normal engineering principles of power system operation (Stott 

1987). They mention a few examples in decoupled formulations: minimizing losses 

with generator real power output as variables would move away from a minimumcost solution; imposing limits on MW reserves with only generator voltage controls 

and transformer voltage tap controls, but no real power control to meet the reserve 

limit. They state that it is helpful to associate each control, constraint, and objective 

in a decoupled OPF with either or both the active and reactive power subproblems. 

They further note that some objective functions and constraints are not algebraic or 

differentiable, and that multiple solutions are likely to exist, in particular when 

there are many reactive power controls (such as switched capacitors, FACTS 

devices, or generators) in network loops. 

It is possible to formulate an objective function that includes the cost of 

reactive power. For a generator the cost of generation is a function of the apparent 

power generated, c(S) = cP(P) + cQ(Q), where S = (P2 + Q2)1/2. If we assume that the 

cost of reactive power is small compared to the cost of real power and if the cost 

function, c(S), is linear in S, an approximation of c(S) is 

c(S) ≈ cP(P) + cQ(|Q|).

Bus-type. In P, Q, |V|, θ space, there are four quantities at each bus: voltage 

magnitude (V), voltage angle (θ), real power (P), and reactive power (Q). In a power 

flow solution without optimization, buses are classified into three bus types: PQ, PV 

and slack. PQ buses generally correspond to loads and PV buses to generators. 

Generator buses are called PV buses because power and voltage magnitude are 
History of Optimal Power Flow and Formulations • December 2012

Page 24

fixed; load buses are known as PQ buses because real and reactive power are fixed, 

that is, Pmin = Pmax and Qmin = Qmax; slack or reference buses have a fixed voltage 

magnitude and voltage angle. For a power flow to solve, the slack bus needs to have 

sufficient real and reactive power to make up for system losses and hold the slack 

bus voltage magnitude at 1; for this reason, a bus with a large generator is

commonly chosen as a slack bus. Table 3 compares the different types of buses.

Table 3: Bus classification used in power flow problems

Bus Type Fixed quantities Variable quantities Physical 

interpretation

PV real power, voltage 

magnitude

reactive power, voltage 

angle

generator

PQ real power, reactive 

power

voltage magnitude, 

voltage angle

load, or generator with 

fixed output

Slack voltage magnitude, 

voltage angle

real power, reactive 

power

an arbitrarily chosen 

generator

In a power flow, the slack bus serves partly to ensure an equal number of 

variables and constraints; without a designated slack bus, the system would be 

over-determined, with more equations than unknowns. Stott states that the need for 

a slack bus also arises because the system I2R losses are not precisely known in 

advance of the load-flow calculation for linear DC models and therefore cannot be 

assigned to a particular generator dispatch (Stott 1974). Some models use a 

distributed slack bus where generators at several different buses provide system 

slack.

We note that an ACOPF that iterates between a simplified OPF and an AC 

power flow may need a slack bus for the power flow iterations, but even then the 

voltage magnitude at the slack bus does not have to be fixed.

When using an iterative method such as Newton or Gauss-Seidel to solve the 

power flow equations, the convergence tolerance is generally set based on the 

“mismatch” terms. Mismatch refers to the difference between known values at each 

bus, such as P and Q at load buses, and the values P(x) and Q(x) computed with the 

power flow equations at each iteration. 
History of Optimal Power Flow and Formulations • December 2012

Page 25

Since the ACOPF is an optimization problem, where the number of variables 

does not have to equal the number of constraints, specifying a slack or reference bus 

is unnecessary. In fact, Carpentier noted this as early as 1962 (Carpentier 1962).1 In 

all optimization formulations herein, we forgo the bus type designation. In an 

optimization context, these categorizations seem overly prescriptive, and could 

unnecessarily over-constrain the problem. For example, fixing the reference voltage 

magnitude at 1.0 per unit when in normal operations generators vary voltage 

magnitude between 0.95 and 1.05 per unit could result in a sub-optimal solution. 

Most modern solvers pre-process the problem, removing variables that have equal 

lower and upper bounds and replacing them with a constant. 

ACOPF Power-Voltage (PQV) Formulation. Most of the ACOPF literature uses the 

polar power-voltage formulations based on the early work of Carpentier during the 

1960’s (Carpentier 1962). 

Polar Power-Voltage Formulation. The polar power-voltage (polar PQV) ACOPF 

(polar ACOPF-PQV) replaces quadratic equality constraints in (32) with the polar 

formulation of (27)-(28):

 Network-wide objective function: Min c(S) (26)

 Network-wide constraints: 

Pn = ∑mk VnVm(Gnmkcosθnm + Bnmksinθnm) (27)

Qn = ∑mk VnVm(Gnmksinθnm - Bnmkcosθnm) (28)

Vmin ≤ V ≤ Vmax (29)

 θminnm ≤ θn - θm ≤ θmaxnm. (30)

 

1 Rough translation of (Carpentier 1962): If voltage and angle are taken as variables in place of P and Q, the restriction of 

fixing the reference voltage can be lifted; voltage and angle are in effect independent variables that fix the state of the network, 

and it suffices to write an objective function that is minimized with respect to these variables. The arbitrarily chosen reference 

bus disappears and the problem is the most general that one can pose.
History of Optimal Power Flow and Formulations • December 2012

Page 26

In this formulation, (27) and (28) represent 2N nonlinear equality constraints with 

quadratic terms and sine and cosine functions that apply throughout the network. 

In this formulation, voltage angle difference constraints are linear. In the 

rectangular formulation discussed below, arctan functions appear in the angle 

difference constraints.

Rectangular Power Voltage Formulation. The rectangular power-voltage

formulation, shown below, is less common in the literature. The rectangular powervoltage (rectangular PQV) ACOPF (rectangular ACOPF-PQV) formulation is shown 

below.

Network-wide objective function: Min c(S) (31)

Network-wide constraint: P + jQ = S = V•I* = V•Y*V* (32)

Bus-specific constraints

Pmin ≤ P≤ Pmax (33)

Qmin ≤ Q ≤ Qmax (34)

(|snmk|)2 ≤ (smaxk)2 for all k (35)

(29) is replaced by:

Vr•Vr + Vj•Vj ≤ (Vmax)2 (36)

(Vmin)2 ≤ Vr •Vr + Vj•Vj (37)

(30) is replaced by:

θminnm ≤ arctan(vj

n/vrn) - arctan(vjm/vrm) ≤ θmaxnm (38)

In this formulation, (32) represents 2N quadratic equalities that apply throughout 

the network; (33)-(34) are simple variable bounds at each bus; (35) and (37) 

represents convex quadratic inequalities at each bus; (37) represents a nonconvex 

quadratic inequalitiy at each bus; and (37) and (38) represents nonconvex 

inequalities between each set of connected buses. 

ACOPF Current Injection (IV) Formulation. Current injection formulations use 

power flow equations based on current and voltage rather than power flow 

equations based on power and voltage discussed above. We only consider the 

rectangular current-voltage (rectangular IV) ACOPF (rectangular ACOPF-IV)

formulation due to the advantages in expressing the current injections as linear 
History of Optimal Power Flow and Formulations • December 2012

Page 27

equality constraints; however, the polar current-voltage formulation could be easily 

derived.

The IV formulation has 6N variables (P, Q, Vr, Vj, Ir, Ij) and the VΘ has 4N variables 

(P, Q, |V|, Θ).

Rectangular ACOPF-IV formulation. The rectangular ACOPF-IV formulation is

shown below.

Network-wide objective function: Min c(S) (40)

Network-wide constraint: I = YV (41)

Bus-specific constraints:

P = Vr•Ir + Vj•Ij ≤ Pmax (42)

Pmin≤ P = Vr•Ir + Vj•Ij (43)

Q = Vj•Ir - Vr•Ij ≤ Qmax (44)

Qmin ≤ Q = Vj

•Ir - Vr•Ij (45)

Vr•Vr + Vj

•Vj ≤ (Vmax)2 (46)

(Vmin)2 ≤ Vr •Vr + Vj•Vj (47)

(inmk)2 ≤ (imaxk)2 for all k (48)

θminnm ≤ arctan(vj

n/vrn) - arctan(vjm/vrm) ≤ θmaxnm (49)

In this formulation, (41) represents 2N linear equality constraints that apply 

throughout the network. This is in contrast to the PQV formulations where 

quadratic and trigonometric constraints apply throughout the network and linear 

constraints are isolated at each bus. Equations (42) to (45) are local quadratic 

nonconvex constraints. Equations (46) and (48) are local convex quadratic 

inequality constraints, but (47) are non-convex local quadratic inequality 

constraints. Overall, the constraint set is still nonconvex, but we hypothesize that 

this formulation may be easier to solve than the power-voltage formulations, since 

the nonlinearities are isolated to each bus and each transmission element, while the 

constraints that apply throughout the network are linear. In general, linear solvers 

solve problems faster than nonlinear solvers. As discussed previously, the voltage 

angle limit (49) could be replaced with an analogous current limit and the problem 

becomes locally quadratic with linear network equations, and (48) and (49) are 

essentially redundant constraints. 
History of Optimal Power Flow and Formulations • December 2012

Page 28

Polar PQV Rectangular PQV Rectangular IV

Network 

constraints

2N nonlinear equality 

constraints with 

quadratic terms and 

sine and cosine 

functions

2N quadratic equalities 2N linear equality 

constraints

Voltage angle 

difference 

constraints

Linear Nonconvex (arctan) Linear if replaced 

with current or 

apparent power 

constraint

Bus constraints Linear Noncovex quadratic 

inequalities

Locally quadratic, 

some nonconvex, 

some convex

6. Literature Review of Formulations

Most literature uses the polar power-voltage formulation, while a smaller 

group of papers use the rectangular power-voltage formulation. Some have also 

proposed hybrid and alternative formulations. So, rather than attempt to review the 

vast literature on the traditional formulation based on power and reactive power 

equations, we focus on alternative formulations in this section.

Stott et al. criticize that much OPF research since the classical formulations of 

Carpentier, Dommel and Tinney have addressed similar formulations without 

considering the additional requirements needed for practical real-time applications, 

partly because OPF problems are still stretching the limits of applied optimization 

technology, and also that utilities have been slow to adopt software to calculate OPF 

“on-line,” or in near-real-time (Stott 1987). They further note that it is a mistake to 

analytically formulate OPF problems by regarding them as simple extensions of 

conventional power flow; once the power flow problem is formulated as an 

optimization problem with degrees of freedom, problems that appear easy to solve 

can turn out to be badly posed, for example with conflicting objective function, 

controls, and constraints. For OPF, they note that researchers have not agreed on

“rules of solvability,” which are the engineering criteria needed for an OPF solution 
History of Optimal Power Flow and Formulations • December 2012

Page 29

to be physically valid, especially for voltage and reactive power, and that these 

“rules of solvability” have hardly if ever been mentioned in the vast literature on 

OPF. They also identify several common problems with the OPF formulation. Most of 

these relate to modeling voltage characteristics of generation, load, and 

transformers, but also include problems with incompatibility of objective, controls, 

and constraints. For example, one incompatibility problem uses an objective of 

minimizing losses with generator real power outputs as variables, rather than fixing 

generator real power outputs at the minimum cost dispatch and adjusting reactive 

power settings to minimize losses (Stott 1987). 

A few researchers have developed a current injection formulation for the 

power flow or optimal power flow equations. Current injection and reactive current 

are terms used in the literature for a formulation similar to the IV formulation 

discussed earlier in this paper. Additionally, some literature uses the term “in 

phase” for the real component of current (Ir) and “quadrature” for the imaginary 

component of current (Ij); in this context, quadrature refers to being 90 degrees out 

of phase. Most of these papers identify challenges modeling generator, or PV buses, 

where the real power injection and voltage magnitude are known but the reactive 

power injection is not. Several authors have identified ways to model PV buses. We 

discuss these formulations here.

Dommel et al. present a power flow formulation using current injections and 

a mix of polar and rectangular coordinates, where each PQ bus is represented by 

two equations for the real and imaginary components of current mismatches in 

terms of complex voltage in rectangular coordinates, while PV buses are 

represented by a single active power mismatch equation and associated voltage 

angle deviation (Dommel 1970). Tinney later mentions that a current injection 

algorithm with a constant nodal admittance matrix could not be used for general 

power flow applications because a satisfactory method of modeling PV buses had 

not been developed (Tinney 1991). Other authors allude to difficulties modeling PV 

buses using current injections, and much of literature using current injection 

formulations applies to radial distribution networks where PV buses are less 

common. For some solution techniques, modeling PV buses with current injection 
History of Optimal Power Flow and Formulations • December 2012

Page 30

equations introduces singularities into some matrices in the solution technique. 

Substitutions introduce dependencies in the Jacobian, meaning that the entire 

Jacobian would have to be recalculated at each step (Gómez Romero 2002). Various 

authors have proposed substitutions and approximations to model PV buses in a 

current injection formulation.

Stadlin and Fletcher discuss a “voltage versus reactive current” model for 

voltage and reactive control that is well suited for use with a linear programming 

algorithm (Stadlin 1982). This paper does not directly discuss an OPF, but provides 

a model that could be used in a linear programming optimization for reactive 

dispatch and voltage control. The model would be used after a real power dispatch 

model, such as a decoupled power flow, was run, and would assume fixed real 

power generation, except at the swing bus. This model uses real and reactive 

current (computed as P/V and Q/V, respectively). The authors use an incremental 

current model rather than an incremental power model because the Jacobian 

matrices of a current model are less sensitive to bus voltage variations. In addition, 

the sensitivity coefficient of voltage to reactive current is much less sensitive than 

the sensitivity coefficient of voltage to reactive power. The authors fix the swing bus 

voltage angle at zero, but allow the voltage magnitude to float. The authors note 

their assumptions result in a more accurate “decoupled” relationship between 

incremental reactive current and voltage than is given by the B matrix used in B-θ 

decoupled OPF, and that this more accurate and more linear model reduces the 

iterations in an optimization algorithm. The sensitivity coefficients in the B matrix 

are accurate only in a small range of voltage, requiring recalculation of the B matrix 

for large changes in voltage; Stadlin and Fletcher’s model is accurate and linear over 

a larger voltage operating range than a B-θ model. Stadlin and Fletcher wanted to 

define a model which remains nearly linear for changes in voltage and reactive 

variables so that efficient linear programming techniques could be applied. 

Da Costa and Rosa note that for the rectangular formulation, generation or 

PV buses have different equations than load or PQ buses. At load buses, active and 

reactive power mismatches are known. At generation buses, reactive power 

mismatches are not known but voltage magnitude constraints are known, because 
History of Optimal Power Flow and Formulations • December 2012

Page 31

in a traditionally formulated power flow, generator reactive power output is 

variable (Da Costa 2008). Therefore, a voltage magnitude constraint is added to 

each load bus, resulting in a different Jacobian matrix. 

Da Costa et al. present a rectangular formulation of a Newton-Raphson 

power flow based on current injections, for both PQ and PV buses (Da Costa 1999, 

Lin 2008). In this formulation, the Jacobian matrix has the same structure as the 

nodal admittance matrix, except for PV buses. For PV buses, the authors introduce a 

new dependent variable, ΔQ, and an additional constraint on voltage magnitude 

deviation. The voltage magnitude constraint is linearized:

ΔVn = 0 ≈ (Vr

n/Vn)ΔVrn + (Vjn/Vn)ΔVjn, where Vn is the voltage magnitude at bus 

n, Vr

n is the real component of voltage at bus n, and Vjn is the imaginary component 

of voltage at bus n. 

Da Costa and Rosa note that the current injection equations are linear for 

electrical networks with only PQ buses and a constant impedance load model (Da 

Costa 2008). 

Jiang et al. published a power-current hybrid rectangular OPF formulation. 

They divide buses into two types, those with non-zero injections, and those with 

zero injections (Jiang 2009). For buses with non-zero injections, the power 

mismatch formulation is used, while the current mismatch formulation is used for 

buses with zero injections. The authors note that in the current mismatch 

formulation, which is similar to the IV formulation presented above, the first-order 

derivatives of the equations are constants and the second-order derivatives are 

zeros. By dividing the buses into two groups, the hybrid method saves computation 

time for the Jacobian and Hessian matrices. 

Meliopoulos and Tao use a formulation referred to as “Quadratic Power 

Flow,” with current conservation equations from Kirchhoff’s current law in 

rectangular coordinates instead of power flow equations, and add operational 

constraints to the model only when they are violated in the previous iteration

(Meliopoulos 2011). The equations modeling generators, constant power loads, and 

transformers are quadratic equations separated into real and imaginary parts. The 

objective function is to minimize the sum of a penalty factor times the sum of 
History of Optimal Power Flow and Formulations • December 2012

Page 32

current mismatches and the total generator costs. The model includes a slack bus as 

the “mismatch current source” where the voltage magnitude is a state variable and 

the real and imaginary components of complex voltage are control variables, while a 

PV bus has the voltage magnitude as a control variable and real and imaginary 

components of complex voltage as state variables. The authors linearize to eliminate 

integer state variables. The quadratic constraints are linearized when they are 

added to the model. 

7. Conclusions

This paper has presented a literature review of different formulations of the 

ACOPF and discussed areas for future research where the ACOPF could be 

improved. The ACOPF problem is inherently difficult due to nonconvexities, 

multipart nonlinear pricing, and alternating current. We do not yet have practical 

approaches to solving nonconvex problems. The ACOPF is a well-structured 

problem, and has developed during 50 years of research. Academia and industry 

have developed various approaches to solving the ACOPF, with different 

formulations, algorithms, and assumptions. The traditional approach has been to 

linearize the full ACOPF problem and decompose it into subproblems. The ACOPF is 

not a hypothetical problem – it is solved every 5 minutes through approximations 

and judgment. After 50 years, there is not yet a commercially viable full ACOPF. 

Many possibilities and ways to examine the ACOPF remain. Today’s solvers do not 

return the gap between the given and globally optimal solution; if we make a rough 

estimate that today’s solvers are on average off by 10%, and world energy costs are 

$400 billion, closing the gap by 10% is a huge financial impact.
History of Optimal Power Flow and Formulations • December 2012

Page 33

References

J. Carpentier, “Contribution á l’étude du dispatching économique,” Bulletin de la 

Société Française des Électriciens, ser. 8, vol. 3, pp. 431-447, 1962.

J. Carpentier, “Optimal power flows,” International Journal of Electrical Power and 

Energy Systems, Vol. 1, Issue 1, pp. 3-15, Apr. 1979.

V.M. da Costa, N. Martins, J.L.R. Pereira, “Developments in the Newton Raphson 

Power Flow Formulation Based on Current Injections,” IEEE Transactions on Power 

Systems, Vol. 14, No. 4, pp. 1320-1326, Nov. 1999.

V.M. da Costa and A.L.S. Rosa, “A Comparative Analysis of Different Power Flow 

Methodologies,” IEEE Transmission & Distribution Conference and Exposition: Latin 

America,” Bogota, Aug. 13-15, 2008.

H.W. Dommel, W.F. Tinney, and W.L. Powell, “Further Developments in Newton’s 

method for power system applications,” IEEE Winter Power Meeting, Conference 

Paper No. 70 CP 161-PWR, Jan. 1970, cited in V.M. da Costa, N. Martins, J.L.R. Pereira, 

“Developments in the Newton Raphson Power Flow Formulation Based on Current 

Injections,” IEEE Transactions on Power Systems, Vol. 14, No. 4, p. 1320, Nov. 1999. 

U.S. Department of Energy, Energy Information Administration (EIA), International 

Energy Statistics, 2012, available at 

http://www.eia.gov/cfapps/ipdbproject/IEDIndex3.cfm?tid=2&pid=2&aid=2; also 

Wholesale Market Data, 2012, available at 

http://www.eia.gov/electricity/wholesale/.

FERC Staff, “Principles for Efficient and Reliable Reactive Power Supply and 

Consumption,” available at 

http://www.ferc.gov/eventcalendar/files/20050310144430-02-04-05-reactivepower.pdf, Feb. 2005.

FERC Joint Boards on Security Constrained Economic Dispatch, “Study and 

Recommendations Regarding Security Constrained Economic Dispatch (SCED) in 

the Northeast by the Joint Board on Economic Dispatch for the Northeast Region,”

FERC Docket AD05-13, page 15, available at 

http://www.ferc.gov/industries/electric/indus-act/joint-boards/northeastrecom.pdf, 2006.
History of Optimal Power Flow and Formulations • December 2012

Page 34

FERC Staff, “Recent ISO Software Enhancements and Future Modeling Plans,” 

available at http://www.ferc.gov/industries/electric/indus-act/rto/rto-iso-soft2011.pdf, Nov. 2011.

FERC and NERC staff, “Arizona-Southern California Outages on September 8, 2011: 

Causes and Recommendations,” pp. 23-59, Apr. 2012.

FERC Order, Midwest Independent System Operator, 138 FERC ¶ 61,235, FERC 

Docket ER12-678 (2012).

A. Gómez Expósito, E. Romero Ramos, “Augmented Rectangular Load Flow Model,” 

IEEE Transactions on Power Systems, Vol. 17, No. 2, pp. 271-274, May 2002.

H.H. Happ, “Optimal Power Dispatch – A Comprehensive Survey,” IEEE Transactions 

on Power Apparatus and Systems, Vol. PAS-96, no. 3, pp. 841-844, May/Jun. 1977.

I.A. Hiskens and R.J. Davy, “Exploring the Power Flow Solution Space Boundary,” 

IEEE Transactions on Power Systems, Vol. 16, No. 3, pp. 389-395, Aug. 2001.

M. Huneault and F.D. Galiana, “A Survey of the Optimal Power Flow Literature,” IEEE 

Transactions on Power Systems, Vol. 6 No. 2, May 1991, pp. 762-770.

Q.Y. Jiang, H.D. Chiang, C.X. Guo, Y.J. Cao, “Power-current hybrid rectangular 

formulation for interior-point optimal power flow,” Institute of Engineering and 

Technology (IET) Generation, Transmission and Distribution, Vol. 3, Iss. 8, pp. 748-

756, 2009.

L.K. Kirchmayer, Economic Operation of Power Systems, Wiley and Sons 1958. 

W.M. Lin, C.H. Huang, T.S. Zhan, “A Hybrid Current-Power Optimal Power Flow 

Technique,” IEEE Transactions on Power Systems, Vol. 23, No. 1, pp. 177-185, Feb. 

2008. 

J.A. Momoh, R.J. Koessler, M.S. Bond, B. Stott, D. Sun, A. Papalexopoulos, P. 

Ristanovic, “Challenges to Optimal Power Flow,” IEEE Transactions on Power 

Systems, Vol. 12, No. 1, pp. 444-447, Feb. 1997.

R.P. O’Neill, “The IV Formulation of the ACOPF and its Linearizations,” FERC 

Technical conference to discuss opportunities for increasing real-time and dayahead market efficiency through improved software , Docket No. AD10-12-003, 
History of Optimal Power Flow and Formulations • December 2012

Page 35

available at http://www.ferc.gov/EventCalendar/Files/20120626084555-

Tuesday_SessionTB_ONeill.ppt.pdf, June 26, 2012.

R. P. O’Neill, T. Dautel and E. Krall, Recent ISO Software Enhancements and Future 

Software and Modeling Plans, Staff Report, Federal Energy Regulatory Commission, 

available at http://www.ferc.gov/industries/electric/indus-act/rto/rto-iso-soft2011.pdf, November 2011, and Private communication with ISO personnel.

J. Peschon, D.S. Piercy, W.F. Tinney, O.J. Tveit, M. Cuénod, “Optimum Control of 

Reactive Power Flow,” IEEE Transactions on Power Apparatus and Systems, Vol. 

PAS-87, No. 1, pp. 40-48, Jan. 1968.

PJM, “PJM’s Striving for Perfect Dispatch Nets Nearly $200 Million in Savings,” News 

Release, Jan. 19, 2012, available at https://www.pjm.com/~/media/aboutpjm/newsroom/2012-releases/20120119-perfect-dispatch-news-release2012.ashx.

J. Riquelme Santos, A. Gómez Expósito, J.L. Martínez Ramos. “Slack Bus Selection to 

Minimize the System Power Imbalance in Load-Flow Studies.” IEEE Transactions on 

Power Systems, 19(2), May 2004.

A.M. Sasson and F.J. Jaimes, “Digital Methods Applied to Power Flow Studies,” IEEE 

Trans. on Power Apparatus and Systems, vol. 86, No. 7, pp. 860-867, Jul. 1967.

A. Schecter, “Exploration of the ACOPF Feasible Region for the Standard IEEE Test 

Set,” FERC Technical conference to discuss opportunities for increasing real-time 

and day-ahead market efficiency through improved software , Docket No. AD10-12-

003, available at http://www.ferc.gov/EventCalendar/Files/20120626084625-

Tuesday,%20Session%20TC-1,%20Schecter%20.pdf, June 26, 2012.

R.B. Squires, “Economic Dispatch of Generation Directly from Power System 

Voltages and Admittances,” AIEE Trans. vol. 79, pt. III, pp. 1235-1244, 1961.

W.O. Stadlin and D.L. Fletcher, “Voltage Versus Reactive Current Model for Dispatch 

and Control,” IEEE Transactions on Power Apparatus and Systems, Vol. PAS-101, 

NO. 10, pp. 3751-3760, Oct. 1982.

B. Stott, “Review of Load-Flow Calculation Methods,” Proceedings of the IEEE, vol. 

62 no. 7, page 916, Jul. 1974.
History of Optimal Power Flow and Formulations • December 2012

Page 36

B. Stott, O. Alsac, and A.J. Monticelli, “Security Analysis and Optimization,” 

Proceedings of the IEEE, vol. 75, no. 12, pp. 1623-1644, 1987.

Y. Tao and A.P.S. Meliopoulous , “Optimal Power Flow via Quadratic Power Flow,” 

IEEE Power Systems Conference and Exposition, Phoenix, May 2011. 

W.F. Tinney, J.M. Bright, K.D. Demaree, B.A. Hughes, “Some Deficiencies in Optimal 

Power Flow,” IEEE Transactions on Power Systems, Vol. 3 No. 2, pp. 676-683, May 

1988.

W.F. Tinney, “A presentation to the workshop in engineering mathematics and 

computer sciences,” EPRI publication EAR/EL-7/07, pages 1-15, August 1991, cited 

in W.-M. Lin, C.-H. Huang, T.-S. Zhan, “A Hubrid Current-Power Optimal Power Flow 

Technique,” IEEE Transactions on Power Systems, Vol. 23, No. 1, p. 177, Feb. 2008. 

J.B. Ward and H.W. Hale, “Digital Computer Solution of Power Flow Problems,” 

Trans. AIEE (Power Apparatus and Systems), vol. 75, pp. 398-404, Jun. 1956.