Convex Relaxation of Optimal Power Flow, Part I: Formulations and Equivalence

I Introduction

For our purposes an optimal power flow (OPF) problem is a mathematical program that seeks to minimize a certain function, such as total power loss, generation cost or user disutility, subject to the Kirchhoff’s laws as well as capacity, stability and security constraints. OPF is fundamental in power system operations as it underlies many applications such as economic dispatch, unit commitment, state estimation, stability and reliability assessment, volt/var control, demand response, etc. There has been a great deal of research on OPF since Carpentier’s first formulation in 1962 . An early solution appears in and extensive surveys can be found in e.g. .

Power flow equations are quadratic and hence OPF can be formulated as a quadratically constrained quadratic program (QCQP). It is generally nonconvex and hence NP-hard. A large number of optimization algorithms and relaxations have been proposed. A popular approximation is a linear program, called DC OPF, obtained through the linearization of the power flow equations e.g. . See also for a more accurate linear approximation. To the best of our knowledge solving OPF through semidefinite relaxation is first proposed in as a second-order cone program (SOCP) for radial (tree) networks and in as a semidefinite program (SDP) for general networks in a bus injection model. It is first proposed in as an SOCP for radial networks in the branch flow model of . See Remark 6 below for more details. While these convex relaxations have been illustrated numerically in and , whether or when they will turn out to be exact is first studied in . Exploiting graph sparsity to simplify the SDP relaxation of OPF is first proposed in and analyzed in .

Convex relaxation of quadratic programs has been applied to many engineering problems; see e.g. . There is a rich theory and extensive empirical experiences. Compared with other approaches, solving OPF through convex relaxation offers several advantages. First, while DC OPF is useful in a wide variety of applications, it is not applicable in other applications; see Remark 10. Second a solution of DC OPF may not be feasible (may not satisfy the nonlinear power flow equations). In this case an operator may tighten some constraints in DC OPF and solve again. This may not only reduce efficiency but also relies on heuristics that are hard to scale to larger systems or faster control in the future. Third, when they converge, most nonlinear algorithms compute a local optimal usually without assurance on the quality of the solution. In contrast a convex relaxation provides for the first time the ability to check if a solution is globally optimal. If it is not, the solution provides a lower bound on the minimum cost and hence a bound on how far any feasible solution is from optimality. Unlike approximations, if a relaxed problem is infeasible, it is a certificate that the original OPF is infeasible.

This two-part tutorial explains the main theoretical results on semidefinite relaxations of OPF developed in the last few years. Part I presents two power flow models that are useful in different situations, formulates OPF and its convex relaxations in each model, and clarifies their relationship. Part II presents sufficient conditions that guarantee the relaxations are exact, i.e. when one can recover a globally optimal solution of OPF from an optimal solution of its relaxations. We focus on basic results using the simplest OPF formulation and does not cover many relevant works in the literature, such as stochastic OPF e.g. , distributed OPF e.g. , new applications e.g. , or what to do when relaxation fails e.g. , to name just a few.

Many mathematical models have been used to model power networks. In Part I of this two-part paper we present two such models, we call the bus injection model (BIM) and the branch flow model (BFM). Each model consists of a set of power flow equations. Each models a power network in that the solutions of each set of equations, called the power flow solutions, describe the steady state of the network. We prove that these two models are equivalent in the sense that there is a bijection between their solution sets (Section II). We formulate OPF within each model where the power flow solutions define the feasible set of OPF (Section III). Even though BIM and BFM are equivalent some results are much easier to formulate or prove in one model than the other; see Remark 2 in Section II.

The complexity of OPF formulated here lies in the nonconvexity of power flow equations that gives rise to a nonconvex feasible set of OPF. We develop various characterizations of the feasible set and design convex supersets based on these characterizations. Different designs lead to different convex relaxations and we prove their relationship (Sections IV and V). When a relaxation is exact an optimal solution of the original nonconvex OPF can be recovered from any optimal solution of the relaxation. In Part II we present sufficient conditions that guarantee the exactness of convex relaxations.

Branch flow models are originally proposed for networks with a tree topology, called radial networks, e.g. . They take a recursive structure that simplifies the computation of power flow solutions, e.g. . The model of also has a linearization that offers several advantages over DC OPF in BIM; see Remark 10. The linear approximation provides simple bounds on the branch powers and voltage magnitudes in the nonlinear BFM (Section VI). These bounds are used in to prove a sufficient condition for exact relaxation.

We make algorithmic recommendations in Section VII based on the results presented here.

This extended version differs from the journal version only in the addition of two Appendices. Appendix VIII provides some mathematical preliminaries and Appendix IX proofs of all main results. Even though all proofs can be found in their original papers, we provide proofs here because (i) it is convenient to have all proofs in one place and in a uniform notation, and (ii) some of the formulations and presentations here are slightly different from those in the original papers.

I-B Notations

II Power flow models

In this section we describe two mathematical models of power networks and prove their equivalence. By a “mathematical model” we mean a set of variables and a set of equations relating these variables. These equations are motivated by the physical system, but mathematically, they are the starting point from which all claims are derived.

The bus injection model (BIM) is defined by the following power flow equations that describe the Kirchhoff’s laws:

Let the set of power flow solutions $V$ for each $s$ be:

For convenience we include $V_{0}$ in the vector variable $V:=(V_{j},j\in N^{+})$ with the understanding that $V_{0}:=1\angle 0^{\circ}$ is fixed.

Bus types. Each bus $j$ is characterized by two complex variables $V_{j}$ and $s_{j}$, or equivalently, four real variables. The buses are usually classified into three types, depending on which two of the four real variables are specified. For the slack bus 0, $V_{0}$ is given and $s_{0}$ is variable. For a generator bus (also called $PV$-bus), Re$(s_{j})=p_{j}$ and $|V_{j}|$ are specified and Im$(s_{j})=q_{j}$ and $\angle V_{j}$ are variable. For a load bus (also called $PQ$-bus), $s_{j}$ is specified and $V_{j}$ is variable. The power flow or load flow problem is: given two of the four real variables specified for each bus, solve the $n+1$ complex equations in (1) for the remaining $2(n+1)$ real variables. For instance when all $n$ buses $j\neq 0$ are all load buses, the power flow problem solves (1) for the $n$ complex voltages $V_{j},j\neq 0$, and the power injection $s_{0}$ at the slack bus 0. This can model a distribution system with a substation at bus 0 and $n$ constant-power loads at the other buses. For optimal power flow problems $p_{j}$ and $|V_{j}|$ on generator buses or $s_{j}$ on load buses can be variables as well. For instance economic dispatch optimizes real power generations $p_{j}$ at generator buses; demand response optimizes demands $s_{j}$ at load buses; and volt/var control optimizes reactive powers $q_{j}$ at capacitor banks, tap changers, or inverters. These remarks also apply to the branch flow model presented next.

II-B Branch flow model

The branch flow model (BFM) in is defined by the following set of power flow equations:

where (2b) is the Ohm’s law, (2c) defines branch power, and (2a) imposes power balance at each bus. The quantity $z_{ij}|I_{ij}|^{2}$ represents line loss so that $S_{ij}-z_{ij}|I_{ij}|^{2}$ is the receiving-end complex power at bus $j$ from bus $i$.

For convenience we include $V_{0}$ in the vector variable $V:=(V_{j},j\in N^{+})$ with the understanding that $V_{0}:=1\angle 0^{\circ}$ is fixed.

II-C Equivalence

Even though the bus injection model (1) and the branch flow model (2) are defined by different sets of equations in terms of their own variables, both are models of the Kirchhoff’s laws and therefore must be related. We now clarify the precise sense in which these two mathematical models are equivalent. We say two sets $A$ and $B$ are equivalent, denoted by $A\equiv B$, if there is a bijection between them .

III Optimal power flow

As mentioned in Remark 1 an optimal power flow problem optimizes both variables $V$ and $s$ over the solution set of the BIM (1). In addition all voltage magnitudes must satisfy:

where $\underline{v}_{j}$ and $\overline{v}_{j}$ are given lower and upper bounds on voltage magnitudes. Throughout this paper we assume $\underline{v}_{j}>0$ to avoid triviality. The power injections are also constrained:

where $\underline{s}_{j}$ and $\overline{s}_{j}$ are given bounds on the injections at buses $j$.

OPF constraints. If there is no bound on the load or on the generation at bus $j$ then $\underline{s}_{j}=-\infty-\textbf{i}\infty$ or $\overline{s}_{j}=\infty+\textbf{i}\infty$ respectively. On the other hand (4) also allows the case where $s_{j}$ is fixed (e.g. a constant-power load), by setting $\underline{s}_{j}=\overline{s}_{j}$ to the specified value. For the slack bus 0, unless otherwise specified, we always assume $\underline{v}_{0}=\overline{v}_{0}=1$ and $\underline{s}_{0}=-\infty-\textbf{i}\infty$, $\overline{s}_{0}=\infty+\textbf{i}\infty$. Therefore we sometimes replace $j\in N^{+}$ in (3) and (4) by $j\in N$.

We can eliminate the variables $s_{j}$ from the OPF formulation by combining (1) and (4) into

Then OPF in the bus injection model can be defined just in terms of the complex voltage vector $V$. Define

Let the cost function be $C(V)$. Typical costs include the cost of generating real power at each generator bus or line loss over the network. All these costs can be expressed as functions of $V$. Then the problem of interest is: OPF:

III-B Branch flow model

OPF variants. OPF as defined in (7) and (9) is a simplified version that ignores other important constraints such as line limits, security constraints, stability constraints, and chance constraints; see extensive surveys in and a recent discussion in on real-life OPF problems. Some of these can be incorporated without any change to the results in this paper (e.g. see for models that include shunt elements and line limits). Indeed a shunt element $y_{j}$ at bus $j$ can be easily included in BIM by modifying (1) into:

or included in BFM by modifying (2a) into:

III-C OPF as QCQP

Before we describe convex relaxations of OPF we first show that, when $C(V):=V^{H}CV$ is quadratic in $V$ for some Hermitian matrix $C$, OPF is indeed a quadratically constrained quadratic program (QCQP) by converting it into the standard form. We will use the derivation in for OPF (7) in BIM. OPF (9) in BFM can similarly be converted into a standard form QCQP.

Define the $(n+1)\times(n+1)$ admittance matrix $Y$ by

$Y$ is symmetric but not necessarily Hermitian. Let $I_{j}$ be the net injection current from bus $j$ to the rest of the network. Then the current vector $I$ and the voltage vector $V$ are related by the Ohm’s law $I=YV$. BIM (1) is equivalent to:

where $e_{j}$ is the $(n+1)$-dimensional vector with 1 in the $j$th entry and 0 elsewhere. Hence, since $I=YV$, we have

where $Y_{j}:=e_{j}e_{j}^{H}Y$ is an $(n+1)\times(n+1)$ matrix with its $j$th row equal to the $j$th row of the admittance matrix $Y$ and all other rows equal to the zero vector. $Y_{j}$ is in general not Hermitian so that $V^{H}Y_{j}^{H}V$ is in general a complex number. Its real and imaginary parts can be expressed in terms of the Hermitian and skew Hermitian components of $Y_{j}^{H}$ defined as:

Let their upper and lower bounds be denoted by

Let $J_{j}:=e_{j}e_{j}^{H}$ denote the Hermitian matrix with a single 1 in the $(j,j)$th entry and 0 everywhere else. Then OPF (7) can be written as a standard form QCQP:

IV Feasible sets and relaxations: BIM

Since OPF is a nonconvex QCQP there is a standard semidefinite relaxation through the equivalence relation: for any Hermitian matrix $M$, $V^{H}MV=$ tr $MVV^{H}=$ tr $MW$ for a psd rank-1 matrix $W$. Applying this transformation to the QCQP formulation (10) leads to an equivalent problem of the form:

for appropriate Hermitian matrices $C_{l}$ and real numbers $b_{l}$. This problem is equivalent to (10) because given a psd rank-1 solution $W$, a unique solution $V$ of (10) can be recovered through rank-1 factorization $W=VV^{H}$. Unlike (10) which is quadratic in $V$ this problem is convex in $W$ except the nonconvex rank-1 constraint. Removing the rank-1 constraint yields the standard SDP relaxation.

We start with some basic definitions on partial matrices and their completions; see e.g. for more details. Fix any connected undirected graph $F$ with $n$ vertices and $m$ edges connecting distinct vertices.In this subsection we abuse notation and use $n,m$ to denote general integers unrelated to the number of buses or lines in a power network. A partial matrix $W_{F}$ is a set of $2m+n$ complex numbers defined on $F$:

$W_{F}$ can be interpreted as a matrix with entries partially specified by these complex numbers. If $F$ is a complete graph (in which there is an edge between every pair of vertices) then $W_{F}$ is a fully specified $n\times n$ matrix. A completion $W$ of $W_{F}$ is any fully specified $n\times n$ matrix that agrees with $W_{F}$ on graph $F$, i.e.,

Given an $n\times n$ matrix $W$ we use $W_{F}$ to denote the submatrix of $W$ on $F$, i.e., the partial matrix consisting of the entries of $W$ defined on graph $F$. If $q$ is a clique (a fully connected subgraph) of $F$ then let $W_{F}(q)$ denote the fully-specified principal submatrix of $W_{F}$ defined on $q$. We extend the definitions of Hermitian, psd, and rank-1 for matrices to partial matrices, as follows. A partial matrix $W_{F}$ is Hermitian, denoted by $W_{F}=W_{F}^{H}$, if $[W_{F}]_{jk}=[W_{F}]_{kj}^{H}$ for all $(j,k)\in F$; it is psd, denoted by $W_{F}\succeq 0$, if $W_{F}$ is Hermitian and the principal submatrices $W_{F}(q)$ are psd for all cliques $q$ of $F$; it is rank-1, denoted by rank $W_{F}=1$, if the principal submatrices $W_{F}(q)$ are rank-1 for all cliques $q$ of $F$. We say $W_{F}$ is $2\times 2$ psd (rank-1) if, for all edges $(j,k)\in F$, the $2\times 2$ principal submatrices

are psd (rank-1), denoted by $W_{F}(j,k)\succeq 0$ (rank $W_{F}(j,k)=1)$. $F$ is a chordal graph if either $F$ has no cycle or all its minimal cycles (ones without chords) are of length three. A chordal extension $c(F)$ of $F$ is a chordal graph that contains $F$, i.e., $c(F)$ has the same vertex set as $F$ but an edge set that is a superset of $F$’s edge set. In that case we call the partial matrix $W_{c(F)}$ a chordal extension of the partial matrix $W_{F}$. Every graph $F$ has a chordal extension, generally nonunique. In particular a complete supergraph of $F$ is a trivial chordal extension of $F$.

For our purposes chordal graphs are important because of the result [62, Theorem 7] that every psd partial matrix has a psd completion if and only if the underlying graph is chordal. When a positive definite completion exists, there is a unique positive definite completion, in the class of all positive definite completions, whose determinant is maximal. Theorem 2 below extends this to rank-1 partial matrices.

IV-B Feasible sets

Then the constraints (5) and (3) imply that the partial matrix $W_{G}$ satisfies The constraint (12a) can also be written compactly in terms of the admittance matrix $Y$ as in : $\displaystyle\underline{s}\ \leq\ \text{diag }\left(WY^{H}\right)\ \leq\ \overline{s}$

Following Section III-C these constraints can also be written in a (partial) matrix form as:

We say that a partial matrix $W_{G}$ satisfies the cycle condition if for every cycle $c$ in $G$

When $\angle[W_{G}]_{jk}$ represent voltage phase differences across each line then the cycle condition imposes that they sum to zero (mod $2\pi$) around any cycle. The next theorem, proved in [57, Theorem 3] and , implies that $W_{G}$ has a psd rank-1 completion $W$ if and only if $W_{G}$ is $2\times 2$ psd rank-1 on $G$ and satisfies the cycle condition (13), if and only if it has a chordal extension $W_{c(G)}$ that is psd rank-1. The theorem also holds with psd replaced by negative semidefinite.

Consider the following conditions on $(n+1)\times(n+1)$ matrices $W$ and partial matrices $W_{c(G)}$ and $W_{G}$:

Fix a graph $G$ on $n+1$ nodes and any chordal extension $c(G)$ of $G$. Assuming $W_{jj}>0$, $\left[W_{c(G)}\right]_{jj}>0$ and $\left[W_{G}\right]_{jj}>0$, $j\in N^{+}$, we have:

Given an $(n+1)\times(n+1)$ matrix $W$ that satisfies (14), its submatrix $W_{c(G)}$ satisfies (15).

Given a partial matrix $W_{c(G)}$ that satisfies (15), its submatrix $W_{G}$ satisfies (16) and the cycle condition (13).

Given a partial matrix $W_{G}$ that satisfies (16) and the cycle condition (13), there is a completion $W$ of $W_{G}$ that satisfies (14).

Informally Theorem 2 says that (14) is equivalent to (15) is equivalent to (16)$+$(13). It characterizes a property of the full matrix $W$ (rank $W=1$) in terms of its submatrices $W_{c(G)}$ and $W_{G}$. This is important because the submatrices are typically much smaller than $W$ for large sparse networks and much easier to compute. The theorem thus allows us to solve simpler problems in terms of partial matrices as we now explain.

Fix any chordal extension $c(G)$ of $G$ and define the set of Hermitian partial matrices $W_{c(G)}$:

Finally define the set of Hermitian partial matrices $W_{G}$:

Note that the definition of psd for partial matrices implies that $W_{c(G)}$ and $W_{G}$ are Hermitian. The assumption $\underline{v}_{j}>0,j\in N^{+}$ implies that all matrices or partial matrices have strictly positive diagonal entries.

Theorem 4 suggests three equivalent problems to OPF. We assume the cost function $C(V)$ in OPF depends on $V$ only through the partial matrix $W_{G}$ defined in (11). For example if the cost is total real line loss in the network then $C(V)=\sum_{j}\text{Re }s_{j}=\sum_{j}\sum_{k:(j,k)\in E}\text{Re}\left([W_{G}]_{jj}-[W_{G}]_{jk}\right)y_{jk}^{H}$. If the cost is a weighted sum of real generation power then $C(V)=\sum_{j}\left(c_{j}\,\text{Re }s_{j}+p_{j}^{d}\right)$ where $p_{j}^{d}$ are the given real power demands at buses $j$; again $C(V)$ is a function of the partial matrix $W_{G}$. Then Theorem 4 implies that OPF (7) is equivalent to

IV-C Semidefinite relaxations

This is a second-order cone and hence OPF-socp is indeed an SOCP in the rotated form.

Literature. SOCP relaxation for OPF seems to be first proposed in for the bus injection model (1), and in for the branch flow model (2) as explained in the next section. By defining a new set of variables $v_{j}:=|V_{j}|^{2}$, $R_{jk}:=|V_{j}||V_{k}|\cos(\theta_{j}-\theta_{k})$, and $I_{jk}:=|V_{j}||V_{k}|\sin(\theta_{j}-\theta_{k})$ where $\theta_{j}:=\angle V_{j}$, rewrites the bus injection model (1) in the complex domain as a set of linear equations in these new variables in the real domain and the following quadratic equations:

IV-D Solution recovery

Then it can be checked that $V$ is in (6) and feasible for OPF.

IV-E Tightness of relaxations

Let $C^{\text{opt}},C^{\text{sdp}},C^{\text{ch}},C^{\text{socp}}$ be the optimal values of OPF (7), OPF-sdp (21), OPF-ch (22), OPF-socp (23) respectively. Theorem 4 and Theorem 5 directly imply

$C^{\text{opt}}\geq C^{\text{sdp}}=C^{\text{ch}}\geq C^{\text{socp}}$. If $G$ is a tree then $C^{\text{opt}}\geq C^{\text{sdp}}=C^{\text{ch}}=C^{\text{socp}}$.

Tightness. Theorem 5 and Corollary 6 imply that for radial networks one should always solve OPF-socp since it is the tightest and the simplest relaxation of the three. For mesh networks there is a tradeoff between OPF-socp and OPF-ch/OPF-sdp: the latter is tighter but requires heavier computation. Between OPF-ch and OPF-sdp, OPF-ch is usually preferable as they are equally tight but OPF-ch is usually much faster to solve for large sparse networks. See for numerical studies that compare these relaxations.

IV-F Chordal relaxation

Theorem 2 through Corollary 6 apply to any chordal extension $c(G)$ of $G$. The choice of $c(G)$ does not affect the optimal value of the chordal relaxation but determines its complexity. Unfortunately the optimal choice that minimizes the complexity of OPF-ch is NP-hard to compute.

V Feasible sets and relaxations: BFM

We now present an SOCP relaxation of OPF in BFM proposed in in two steps. We first relax the phase angles of $V$ and $I$ in (2) and then we relax a set of quadratic equalities to inequalities. This derivation pinpoints the difference between radial and mesh topologies. It motivates a recursive version of BFM for radial networks (Section VI) and the use of phase shifters for convexification of mesh networks (Part II ).

i.e., $\beta(x)$ is in the range space of $B$ (mod $2\pi$). A solution $\theta(x)$, if exists, is unique in $(-\pi,\pi]^{n}$. Define the set

where the $n\times n$ submatrix $B_{T}$ corresponds to links in $T$ and the $(m-n)\times n$ submatrix $B_{\perp}$ corresponds to links in $T^{\perp}:=G\setminus T$. Similarly partition $\beta(x)$ into

In that case $\theta(x)=\mathcal{P}\left(B_{T}^{-1}\beta_{T}(x)\right)$ is the unique solution of (26) in $(-\pi,\pi]^{n}$, where $\mathcal{P}(\phi)$ projects $\phi$ to $(-\pi,\pi]^{n}$.

V-B SOCP relaxation

Let $C^{\text{opt}}$ be the optimal cost of OPF (9) in the branch flow model. Let $C^{\text{opf}}$, $C^{\text{nc}}$, $C^{\text{socp}}$ be the optimal costs of OPF (30), OPF-nc (31), OPF-socp (32) respectively defined above. Theorem 9 implies

V-C Equivalence

VI BFM for radial networks

Case I: Links point away from bus 0. Model (24) reduces to:

Use the boundary condition (34d), $S_{n}=f_{n}(s_{0})=0$, to solve for the scalar variable $s_{0}$. The other variables $x_{j}$ can then be computed from (35). This method can be extended to a general radial network with laterals . See also for techniques for solving the nonlinear equations (35), and for a different recursive approach called the forward/backward sweep for radial networks.

Case II: Links point towards bus 0. Model (24) reduces to:

VI-B Linear approximation and bounds

For $j\in N^{+}$, $v_{j}\leq v_{j}^{\text{lin}}$.

For $j\in N^{+}$, $\hat{v}_{j}\leq\hat{v}_{j}^{\text{lin}}$.

Linear approximations. For radial networks the linear approximations (37) and (38) of BFM have two advantages over the (linear) DC approximation of BIM. First they have a simple recursive structure that leads to simple bounds on power flow quantities. Second DC approximation assumes $r_{jk}=0$, fixes voltage magnitudes, and ignores reactive power, whereas (37) and (38) do not. This is important for distribution systems where $r_{jk}$ are not negligible, voltages can fluctuate significantly and reactive powers are used to regulate them. On the other hand (37) and (38) are applicable only for radial networks whereas DC approximation applies to mesh networks as well. See also for a more accurate linearization of BIM that addresses the shortcomings of DC OPF.

VII Conclusion

We have presented a bus injection model and a branch flow model, formulated several relaxations of OPF, and proved their relations. These results suggest a new approach to solving OPF summarized in Figure 2.

For radial networks we recommend solving OPF-socp in either BIM or BFM though there is preliminary evidence that BFM can be more stable numerically. For mesh networks we recommend solving OPF-ch for small networks and OPF-socp followed by a heuristic search for a feasible point for large networks. Also see Remarks 7 and 8.

The key for this solution strategy is that the relaxations are exact so that an optimal solution of the original OPF can be recovered. In Part II of this paper we summarize sufficient conditions that guarantee exact relaxation.

In this appendix we summarize some basic concepts in optimization, matrix completion and chordal relaxation that we use in this two-part tutorial. For notations see Section I. More details can be found in, e.g., .

Quadratic constrained quadratic program (QCQP) is the following problem:

Any psd rank-1 matrix $X$ has a unique spectral decomposition $X=xx^{H}$. Using $x^{H}C_{l}x=\text{tr }C_{l}xx^{H}=:\text{tr }C_{l}X$ we can rewrite a QCQP as the following equivalent problem where the optimization is over Hermitian matrices:

SDP is a convex program and can be efficiently computed. We call (41) an SDP relaxation of QCQP (39) because the feasible set of (40) is a subset of the feasible set of SDP (41). A strategy for solving QCQP (39) is to solve SDP (41) for an optimal $X^{\text{opt}}$ and check its rank. If rank $X^{\text{opt}}=1$ then $X^{\text{opt}}$ is optimal for (40) as well and an optimal solution $x^{\text{opt}}$ of QCQP (39) can be recovered from $X^{\text{opt}}$ through spectral decomposition $X^{\text{opt}}=x^{\text{opt}}(x^{\text{opt}})^{H}$. If rank $X^{\text{opt}}>1$ then, in general, no feasible solution of QCQP can be directly obtained from $X^{\text{opt}}$ but the optimal objective value of SDP provides a lower bound on that of QCQP.

To derive the Lagrangian dual of SDP (41), form the Lagrangian, for $y:=(y_{l},l=1,\dots,L)\geq 0$,

Then the primal problem (41) is equivalent to $\min_{X\succeq 0}\,\max_{y\geq 0}\ L(X;y)$ and its dual is $\max_{y\geq 0}\,\min_{X\succeq 0}\ L(X;y)$ (if we allow their objective values to be $\pm\infty$). Hence the dual objective function is

A pair $(X^{\text{opt}},y^{\text{opt}})$ is a primal-dual optimal if and only if

Primal feasibility: $X^{\text{opt}}\succeq 0$ and tr$\,C_{l}X^{\text{opt}}\,\leq\,b_{l}$, $l=1,\dots,L$.

Dual feasibility: $y^{\text{opt}}\geq 0$ and $C_{0}\,+\,\sum_{l=1}^{L}\,y_{l}^{\text{opt}}\,C_{l}\ \preceq\ 0$.

Complementary slackness: tr$\,\left(C_{0}\,+\,\sum_{l}\,y_{l}^{\text{opt}}\,C_{l}\right)X^{\text{opt}}\ =\ 0$.

A special case of SDP is a second-order cone program (SOCP):

For optimal power flow problems, we use SOCP in the following rotated form:

In this paper we formulate optimal power flow (OPF) problems as QCQPs and describe SDP and SOCP relaxations of OPF. The third relaxation we will discuss is chordal relaxation based on the notion of chordal extension of a network graph. We now review some basic concepts in graph theory, partial matrices and completions, and show that a chordal relaxation is indeed a semidefinite program.

-B Graph, partial matrix and completion

Consider a graph $G=(N,E)$ with $N:=\{1,\dots,n\}$. $G$ can either be undirected or directed with an arbitrary orientation. Two nodes $j$ and $k$ are adjacent if $j\sim k\in E$. A complete graph is one where every pair of nodes is adjacent. A subgraph of $G$ is a graph $F=(N^{\prime},E^{\prime})$ with $N^{\prime}\subseteq N$ and $E^{\prime}\subseteq E$. A clique of $G$ is a complete subgraph of $G$. A maximal clique of $G$ is a clique that is not a subgraph of another clique of $G$.

By a path connecting nodes $j$ and $k$ we mean either a set of distinct nodes $(j,n_{1},\dots,n_{i},k)$ such that $(j\sim n_{1}),(n_{1}\sim n_{2}),\dots,(n_{i}\sim k)$ are edges in $E$ or this set of edges, depending on the context. A cycle $(n_{1},\dots,n_{i})$ is a path such that $(n_{1}\sim n_{2}),\dots,(n_{i}\sim n_{1})$ are edges in $E$. By convention we exclude a pair of adjacent nodes $(j,k)$ as a cycle. We will only consider connected graphs in which there is a path between every pair of nodes.

A cycle in $G$ that has no chord (an edge connecting two nodes that are non-adjacent in the cycle) is called a minimal cycle. $G$ is chordal if all its minimal cycles are of length 3 (recall that an edge $(j,k)$ is not considered a cycle). A chordal extension of $G$ is a chordal graph on the same set of nodes as $G$ that contains $G$ as a subgraph. Every graph has a chordal extension; e.g. the complete graph on the same set of nodes is a trivial chordal extension.

Fix a graph $G=(N,E)$ with $N:=\{1,\dots,n\}$ and $E\subseteq N\times N$. For our purposes here we assume $G$ is undirected so that $(j,k)\in E$ if and only if $(k,j)\in E$. A $G$-partial matrix (or simply a partial matrix if $G$ is clear from the context) is a set of complex numbers:

One can treat a partial matrix $X_{G}$ as entries of an $n\times n$ matrix $X$ whose entries $X_{jk}$ are unspecified if $(j,k)\not\in E$. See Figure 3(a) below for an example. Given a partial matrix $X_{G}$ we call an $n\times n$ matrix $X$ a completion of $X_{G}$ if $X_{jj}=[X_{G}]_{jj},j\in N$, and $X_{jk}=[X_{G}]_{jk},(j,k)\in E$, i.e., $X$ agrees with $X_{G}$ on $G$.We abuse the $X_{G}$ notation: given $G$, $X_{G}$ is a partial matrix defined on $G$, and given an $n\times n$ matrix $X$, $X_{G}$ is the submatrix $(X_{jj},j\in N,X_{jk},(j,k)\in E)$ of $X$ defined by $G$. The meaning should be clear from the context.

Consider any $n\times n$ matrix $X$. Given any $k\leq n$ nodes $(n_{1},n_{2},\ldots,n_{k})$ let $X(n_{1},\ldots,n_{k})$ denote the $k\times k$ principal submatrix of $X$ defined by:

Any maximal clique $q:=(n_{1},n_{2},\ldots,n_{k})$ of $G$ with $k$ nodes defines a (fully specified) $k\times k$ principal submatrix denoted by $X(q):=X(n_{1},\ldots,n_{k})$. In particular each edge $(i,j)\in E$ is a clique and defines a $2\times 2$ principal submatrix $X(i,j)$, which we use heavily in discussing optimal power flow problems. These notions are extended to partial matrices with $X$ replaced by $X_{G}$.

We extend the notions of Hermitian, psd, rank-1, and trace to partial matrices as follows. We say that a partial matrix $X_{G}$ is Hermitian, denoted by $X_{G}=X_{G}^{H}$, if $[X_{G}]_{kj}=\left([X_{G}]_{jk}\right)^{H}$. An $n\times n$ matrix $X$ is psd if and only if all its principal submatrices (including $X$ itself) is psd. We extend the notion of psd to partial matrices using this property, by saying that a partial matrix $X_{G}$ is psd if all its “principal submatrices” that are fully specified are psd. Formally $X_{G}$ is psd, denoted by $X_{G}\succeq 0$, if $X_{G}(q)\succeq 0$ for all maximal cliques $q$ of $G$. Note that if $X_{G}$ is psd then it is Hermitian by definition. Similarly we say that a partial matrix $X_{G}$ is rank-1, denoted by rank $X_{G}=1$, if $X_{G}(q)$ is rank-1 for all maximal cliques $q$ of $G$. We say $W_{G}$ is $2\times 2$ psd on $G$ if, for all $(j,k)\in E$, the $2\times 2$ matrices $W_{G}(j,k)$ are psd, i.e.,

We say $W_{G}$ is $2\times 2$ rank-1 on $G$ if, for all $(j,k)\in E$, $W_{G}(j,k)$ are $2\times 2$ rank-1 matrices, i.e., they are not the zero matrices and

Finally we say that an $n\times n$ matrix $C$ is defined on graph $G$ if $C_{jk}=0$ if $(j,k)\not\in E$. We extend the operation tr to partial matrices $X_{G}$: if $C$ and $X_{G}$ are defined on the same graph $G$ then

Suppose the matrices $C_{l}$ in (41), $l=0,\dots,L$, are all defined on $G$, i.e., for all $l$, $[C_{l}]_{jk}=0$ if $(j,k)\not\in E$. Then given any $n\times n$ matrix $X$, tr $C_{l}X=$ tr $C_{l}X_{G}$ where $X_{G}$ is the submatrix of $X$ defined by $G$. Conversely, given a partial matrix $X_{G}$ that satisfies (41c), any completion $X$ of $X_{G}$ satisfies (41c). Even though both the objective function (41a) and the constraints (41c) depend only on the partial matrix $X_{G}$, the constraint $X\succeq 0$ in (41c) depends also on entries not in $X_{G}$. Indeed the number of complex variables in $X$ is $n^{2}$ while the number of complex variables in $X_{G}$ is only $n+2|E|$, which is much smaller than $n^{2}$ if $G$ is large but sparse. Hence instead of solving for a full psd matrix $X$ directly as in SDP (41) we would like to compute a partial matrix $X_{G}$ that has a psd completion $X$ that satisfies (41c)–(41c). If the completion $X$ is rank-1 then it also solves the problem (40) and hence yields a solution to the original QCQP (39) through spectral decomposition of $X$. Theorem 2 provides an exact characterization of when this is possible.

Two questions naturally arise in this approach: (i) How to formulate a semidefinite relaxation based on a given a chordal extension $F$ of $G$? (ii) How to choose a good chordal extension $F$ of $G$ so that the resulting relaxation can be solved efficiently? We next illustrate the issues involved in these two questions through an example. See for more details.

-C Chordal relaxation

We call this problem a chordal relaxation of QCQP (39). Recall that we assume $C_{l}$, $l=0,\dots,L$, are all defined on $G$, i.e., $[C_{l}]_{jk}=0$ if $(j,k)\not\in E$. This implies that tr\,C_{l}X=\tr\,C_{l}X_{F}=\tr$\,C_{l}X_{G}$. Then chordal relaxation (44) is equivalent to SDP (41) in the sense that given any feasible solution $X_{F}$ of (44), there is a psd completion $X$ that is feasible for (41) and has the same cost, and vice versa. This is a consequence of [62, Theorem 7] that says every psd partial matrix has a psd completion if and only if the underlying graph is chordal. See also Theorem 5 and Corollary 6.

The first step in constructing the chordal relaxation (44) is to list all the maximal cliques $q_{k}$. Even though listing all maximal cliques of a general graph is NP-hard it can be done efficiently for a chordal graph. This is because a graph is chordal if and only if it has a perfect elimination ordering and computing this ordering takes linear time in the number of nodes and edges . Given a perfect elimination ordering all maximal cliques $q_{k}$ can be enumerated and $X_{F}(q_{k})$ constructed efficiently . For optimal power flow problems the computation depends only on the topology of the power network, not on operational data, and therefore can be done offline.

We now show that (44) is indeed an SDP by converting it into the standard form (41) with the introduction of auxiliary variables, following the procedure described in . This conversion also illustrates the difficulty in choosing a good chordal extension $F$ (see Remark 11 below).

The (fully specified) matrices $X_{F}(q_{k})$ in (44c) can be treated as principal submatrices of an $n\times n$ matrix $X$. They may not however be integrated directly into a common $n\times n$ matrix variable $X$ because different $X_{F}(q_{k})$ may share entries. We now explain the issue and its resolution using the example in Figure 1. They are the same in the general case with more cumbersome notations; see .

Suppose we have chosen the chordal extension $F$ in Figure 3(b) with two overlapping cliques $q_{1}$ and $q_{2}$ as explained in the caption of the figure. To decouple the two matrices $X_{F}(q_{1})$ and $X_{F}(q_{2})$, define the $3\times 3$ matrix

where the decoupling variables $u_{jk}$ are constrained to be:

Define the $7\times 7$ block-diagonal matrix

Then the chordal relaxation (44) can be written in the standard form (41) in terms of these $7\times 7$ block-diagonal Hermitian matrices:

for appropriate choices of $C_{l}^{\prime}$, $l=0,\dots,L$. The constraint $X^{\prime}\succeq 0$ in (47d) is equivalent to the requirement (46) on its submatrices and $C_{r}^{\prime}$ in (47d) is chosen to enforce the requirement (45). Hence the chordal relaxation (44) is indeed an SDP.

There are two conflicting factors in choosing a good chordal extension $F$. First an $F$ that contains fewer number of maximal cliques $q$ generally involves larger cliques, leading to larger submatrices $X_{F}(q)$; for example the complete graph $F$ has a single maximal clique but the corresponding $X_{F}(q)=X$ has $n^{2}$ entries and the chordal relaxation (44) offers no computational advantage over solving (in fact it is exactly) the original SDP (41). This argues for a chordal extension $F$ with smaller, possibly more, maximal cliques $q$. Second, however, having more maximal cliques $q$ tends to require more decoupling variables $u_{jk}$. Every decoupling variable $u_{jk}$ introduces an extra equality constraint in (47d), thus increasing the required computational effort. For instance the transformed problem based on the chordal extension in Figure 3(b) involves 2 maximal cliques of sizes 3 and 4, and 4 additional equality constraints in (47d). The transformed problem based on the chordal extension in Figure 3(c), on the other hand, requires 3 maximal cliques each of size 3, and 8 additional equality constraints.

In summary even though the ambient dimension of the new variable $X^{\prime}$ is generally larger than that of the original $n\times n$ matrix variable $X$ ($7\times 7$ as opposed to $5\times 5$ for the example in Figure 3(b)), the chordal relaxation (44) can typically be solved much more efficiently than SDP (41) if $G$ is large and sparse; for OPF examples, see . Choosing a good chordal extension $F$ of $G$ is important but nontrivial. See for methods to compute efficient chordal extensions and sparse SDP solutions.

-D Proof of Theorem 1: equivalence

-E Proof of Theorem 2: rank-1 characterization

We will prove (1) $\Rightarrow$ (2) $\Rightarrow$ (3) $\Rightarrow$ (1). If $W$ is psd rank-1 then all its principle submatrices are psd and of rank 1 (the submatrix cannot be of rank 0 because, by assumption, $W_{jj}>0$ for all $j\in N^{+}$). This implies that its submatrix $W_{c(G)}$ is psd and rank-1. Hence (1) $\Rightarrow$ (2).

Fix a partial matrix $W_{c(G)}$ that is psd and rank-1 and consider its submatrix $W_{G}$. Since each link $(j,k)\in E$ is a clique of $c(G)$ the $2\times 2$ principle submatrix $W_{G}(j,k)$ is psd and rank-1. Therefore to prove that (2) $\Rightarrow$ (3), it suffices to show that $W_{G}$ satisfies the cycle condition (13). We now prove the following statement by induction on $3\leq k\leq n+1$: for all cycles $(j_{1},\dots,j_{k})$ of length $k$ in $c(G)$,

where $j_{k+1}:=j_{1}$. For $k=3$, a cycle $(n_{1},n_{2},n_{3})$ is a clique of $c(G)$ and therefore the following principle submatrix of $W_{c(G)}$:

Suppose (48) holds for all cycles in $c(G)$ of length up to $k>3$. Consider now a cycle $(j_{1},\dots,j_{k+1})$ of length $k+1$ in $c(G)$. Since $c(G)$ is chordal there is a chord, say, $(j_{1},j_{l})\in E$ for some $1<l<k+1$. Since both cycles $(j_{1},\dots,j_{l})$ and $(j_{1},j_{l},\dots,j_{k+1})$ satisfy (48) we have

where $j_{k+2}:=j_{1}$. Since $W_{G}$ is Hermitian, adding the above equations yields

proving (48) for $k+1$. This completes the proof of (2) $\Rightarrow$ (3).

Without loss of generality let $\angle V_{0}=0^{\circ}$; for $j=1,\dots,n$, set

-F Proof of Corollary 3: uniqueness of completion

-G Proof of Theorem 5: BIM feasible sets

-H Proof of Theorem 7: BFM feasible sets

Using the connectedness of $G$ and the definition of $B$, one can argue that $\alpha$ must be an integer vector for $k+B\alpha$ to be integral. We say $\sigma(\theta,k)$ is a solution of (49) if every vector in $\sigma(\theta,k)$ is a solution of (49), and $\sigma(\theta,k)$ is the unique solution of (49) if it is the only equivalence class of solutions. The lemma below implies that if a solution of (49) exists then it is unique.

there is at most one $\sigma(\theta,k)$ with $\theta\in(-\pi,\pi]^{n}$, that is the unique solution of (49) when it exists.

To prove the first claim, suppose $(\theta,k)$ is a solution of (26) for some $k=k(\theta)$. We need to show that (24), (49) together with (50) imply (2). Now (24a) is equivalent to (2a). Moreover (24c) and (50) imply (2c). To prove (2b), substitute (2c) into (49) to get

-I Proof of Theorem 8: BFM cycle condition

with $\theta\in(-\pi,\pi]^{n}$. Since $T$ is a spanning tree, the $n\times n$ submatrix $B_{T}$ is invertible. Moreover (52) has a unique solution if and only if $B_{\perp}B_{T}^{-1}(\beta_{T}+2\pi k_{T})=\beta_{\perp}+2\pi k_{\perp}$, or if and only if

for some integer vector $\hat{k}_{\perp}:=k_{\perp}-B_{\perp}B_{T}^{-1}k_{T}$. ((54) below implies that $\hat{k}_{\perp}$ is indeed an integer vector.)

Finally consider the unique solution $(\theta_{*},k(\theta_{*}))$ of (49) with $\theta_{*}\in(-\pi,\pi]^{n}$. By (52) we have $\theta_{*}=B_{T}^{-1}\beta_{T}+2\pi B_{T}^{-1}k_{T}(\theta_{*})$. The definition of $k(\theta_{*})$ and the fact $\theta_{*}\in(-\pi,\pi]^{n}$ imply that $\theta_{*}=\mathcal{P}\left(B_{T}^{-1}\beta_{T}\right)$; see the discussion preceding Lemma 14. This completes the proof of Theorem 8. ∎

-J Proof of Theorem 9: radial networks

-K Proof of Theorem 11: equivalence

where the second equality follows from (57) and the last equality from (24b). But

This and (4) imply (12a) as desired. To prove $W_{G}(j,k)\succeq 0$ for each $(j,k)\in E$, we have

Theorem 8 therefore implies that the cycle conditions (13) and (26) are equivalent under $g$ and $g^{-1}$.

This completes the proof of Theorem 11. ∎

-L Proof of Lemma 12: voltage bound

The proofs of 12(1)–(3) and Lemma 13 are obvious and omitted. The proof of Lemma 12(4) makes use of Lemma 13 and is provided here.

It is easy to see from (37) and (38) that $S_{jk}^{\text{lin}}=-\hat{S}_{kj}^{\text{lin}}$ and $v^{\text{lin}}=\hat{v}^{\text{lin}}$.

To show $v=\hat{v}$, define the following functions:

To obtain (61c) from (LABEL:eq:bdf.app1), note that the right-hand side of (LABEL:eq:bdf.app1) is

where the second equality follows from (61b). This together with (LABEL:eq:bdf.app1) imply (61c).

Hence Lemma 13(4) implies that $v=\hat{v}\leq\hat{v}^{\text{lin}}=v^{\text{lin}}$ as desired. ∎

References