Open Journal of Mathematical Optimization

We address the issue of computing a global minimizer of the AC Optimal Power Flow problem. We introduce valid inequalities to strengthen the Semideﬁnite Programming relaxation, yielding a novel Conic Programming relaxation. Leveraging these Conic Programming constraints, we dynamically generate Mixed-Integer Linear Programming (MILP) relaxations, whose solutions asymptotically converge to global minimizers of the AC Optimal Power Flow problem. We apply this iterative MILP scheme on the IEEE PES PGLib [2] benchmark and compare the results with two recent Global Optimization approaches.


Motivation and related works
The Alternating-Current Optimal Power Flow (ACOPF) is a seminal optimization problem related to the dispatching of electricity in a power network. The authorship of this problem is attributed to Carpentier [6], who introduced it in 1962 as "Economic Dispatch". Since then, this problem has not only interested the Power Systems research community, but also the community of Operations Research and Mathematical Programming [5,28]. Indeed, ACOPF was identified as a challenging and fruitful application of Nonlinear Programming (NLP) and Global Optimization methods. Thanks to Interior-Point algorithms, developed since the 1990s, the computation of ACOPF feasible solutions and local minimizers is accessible, even for instances of several thousand nodes [33].
To bound the optimality gap of feasible solutions found by such NLP algorithms, several convex relaxations have been introduced during the past decade. A review of the numerous relaxation techniques for the ACOPF problem is available in [28]. Leveraging NLP algorithms and convex relaxations techniques, several approaches emerged to solve the ACOPF problem to global optimality. We gather these works in four different categories.

Relaxation Strengthening and Bound Tightening.
Strengthening the classical convex relaxations [28] such as the rank relaxation helps improving the corresponding lower bounds. This strengthening is possible through additional valid inequalities coming from the polar formulation of the ACOPF problem [10,19,20] or derived from the Reformulation-Linearization Technique (RLT) [31]. Feasibility-Based Bound Tightening (FBBT) and Optimization-Based Bound Tightening (OBBT) techniques [3], the latter being based on the value of a known feasible solution, are also known to be particularly efficient for the ACOPF problem [9,32]. Even if these methods do not have a guarantee of convergence towards a global solution, the aforementioned articles report that they significantly reduce the optimality gap and even close the gap for some instances. Moment-Sum-of-Squares hierarchy. The celebrated Moment-Sum-of-Squares hierarchy of relaxations for polynomial optimization problems [23] has been applied to the ACOPF problem in several works [14,13,27,29]. The convergence of the relaxations' values towards the optimal value of the ACOPF problem is proven, at the price of the rapidly increasing size and computational cost of the resulting convex relaxations. In practice, only the first and second order relaxations are solvable for medium-scale ACOPF instances, using the sparse variant of the Moment-Sum-of-Squares hierarchy [24].
In Section 2, we present the ACOPF problem and an equivalent reformulation of it. Section 3 introduces our valid inequalities, the resulting Conic Programming relaxation, and the Bound Tightening procedure that we apply. The iterative MILP scheme is presented in Section 4 and the numerical experiments in Section 5.

Mathematical notation
For any complex number x ∈ C, x * = Re(x) − i Im(x) is its complex conjugate, |x| is its magnitude and ∠x its argument. We denote by C n×n the C-vector space of n × n matrices with complex entries. We denote by (E ab ) ab the canonical basis of this C-vector space. For any matrix M ∈ C n×n , its Hermitian transpose is M H , defined such as M H ab = M * ba for all b, a ∈ {1, . . . , n}. The R-vector space of Hermitian matrices, H n ⊂ C n×n , is the set of matrices M ∈ C n×n such that (s.t.) M = M H .

2
Mathematical Programming formulations for the ACOPF

Original formulation
A power grid is a network of buses interconnected by lines. We give an arbitrary orientation to each line, so as to distinguish its two extremities. Hence, the grid is modelled as a directed graph N = (B, L) with size n = |B|. The set L is s.t. L ∩ L R = ∅, where L R is the set of couples (b, a) s.t. (a, b) ∈ L. A line ∈ L is described by a couple (b, a) s.t. b ∈ B is the "from" bus (denoted by f), a ∈ B is the "to" bus (denoted by t). Electricity generating units are located at several buses in the network. We denote by G b the set of generators located at bus b ∈ B. The set of all generators is G = b∈B G b , whose cardinality is m = |G|. The parameters of the ACOPF problem are described in Table 1. Table 1 Parameters of the ACOPF problem Parameters Index set Description c1g ∈ R, c2g ∈ R+ g ∈ G Cost parameters sg, sg ∈ C g ∈ G Power injection bounds In the (S,V) formulation [5], the ACOPF problem has two types of decision variables: for any g ∈ G, S g ∈ C is the complex power injection of generator g, for any b ∈ B, V b ∈ C is the complex voltage at bus b.
It is traditionally assumed that the generation cost of each generator g ∈ G is a convex quadratic form in Re(S g ). The objective function is the sum of the generation costs g∈G c 1g Re(S g ) + c 2g Re(S g ) 2 , to be minimized. The decision variables are subject to different types of constraints: Injection Limits for Generators. For each g ∈ G, we have These inequalities between complex numbers designate the respective real inequalities for the real and for the imaginary parts. Voltage Magnitude Limits. For each b ∈ B, the voltage at b satisfies The ACOPF instances in the literature are scaled so that each bus has a nominal voltage normalized to 1. The lower and upper bounds respectively correspond to deviations of at most 50% around this nominal value.
With this notation, we write the Power Flow conservation at bus b ∈ B as Constraint (4) describes the equality between the net injection of power at b and the power transfer towards the adjacent buses. Thermal Limits for Lines. For each line (b, a) ∈ L, the operational limit in terms of apparent power is For (b, a) ∈ L R , this reads Line Phase Angle Difference Limits. For any (b, a) ∈ L ∪ L R , In summary, the ACOPF problem is the following Nonconvex Optimization problem

ACOPF reformulation
Definition 1. A tree decomposition T of the graph N = (B, L) is a tree where each node k ∈ T is associated with a set B k ⊂ B, and satisfying the following properties The union of the subsets B k equals the set B: k∈T B k = B, For every (b, a) ∈ L, there exists k ∈ T s.t {b, a} ⊂ B k , If b ∈ B k ∩ B for any k, ∈ T , then b ∈ B j for all nodes j of the tree T in the unique path between k and .
We consider a given tree decomposition T of the graph N , and we introduce the symmetric set E ⊂ B × B of arcs defined as E = k∈T B k × B k . As a matter of fact, the sets B k are cliques of the undirected graph induced by (B, E). In this respect, the sets B k are called cliques. We denote by H n (E) the set of partially defined matrices W , seen as vectors indexed by E, and s.t. W ba = W * ab for all (b, a) ∈ E. For any k ∈ T , we denote by W B k ,B k the matrix (W ba ) (b,a)∈B 2 k . With this notation, we reformulate (ACOPF) as While the clique-based SDP relaxation is well known, this clique-based reformulation of the ACOPF problem itself is not properly stated in the literature, as far as we know. Yet, we acknowledge that the proof of Theorem 2 is closely related to the developments presented in [7].
We prove the equivalence for the feasibility, which also proves the equivalence for the optimality since both problems share the same objective value. We take (S, V ) a feasible solution in (ACOPF) and we define W ∈ H n (E) as W ba = V b V * a for any (b, a) ∈ E. For any b ∈ B, we make the following observations: . To conclude about the feasibility of (S, W ) in (ACOPF W ), we state that Conversely, we consider any (S, W ) feasible in (ACOPF W ). Since W B k ,B k 0 and |W ba | 2 = W bb W aa for all (b, a) ∈ B 2 k , we can apply [7,Prop. 6] to state that rank W B k ,B k = 1 for all k ∈ T . By induction on the tree decomposition T , we prove that there exists V ∈ C n s.t.
The case |T | = 1 is trivial, since any rank-one positive semidefinite (PSD) matrix W can be written as W = V V H . We assume now that the induction hypothesis is true for any graph with a tree decomposition with size less or equal than p ∈ N * , and we consider a graph N with a tree decomposition T with size p + 1. We consider a leaf k of T , B k the corresponding clique, B = =k B the union of the other cliques, and C k = B k \ B. By property of a tree decomposition, since k is a leaf of T , T \ {k} is a tree decomposition of the graph ( B, E), where E denotes the edges in E that are not adjacent to C k . Applying the induction hypothesis, since T \ {k} has size p, there exists a complex vector V ∈ C |B| s.t.
By induction, this proves that there exists a vector V ∈ C n s.t. W ba = V b V * a for all (b, a) ∈ E. The feasibility of (S, V ) in (ACOPF) follows by substituting W ba by V b V * a in the constraints of (ACOPF W ).

Strengthening the SDP relaxation
In formulation (ACOPF W ), Constraints ( ) are the only nonconvex constraints. Removing them leads to the clique-based SDP relaxation [28,30]. Instead of merely deleting Constraints ( ), we add valid inequalities based on Voltage Magnitude and Phase Angle Difference bounds.

Conic Programming Outer-Approximation of Constraints ( )
that stands for |V b ||V a | and is subject to For all b ∈ B, we also define the following constraints Whereas Constraints (8)-(11) approximate the equality R 2 ba = W bb W aa , we also need to approximate |W ba | = R ba . For this purpose, we impose for all (b, a) ∈ E, For all (b, a) ∈ E \ (L ∪ L R ), we define θ ba = −2π and θ ba = 2π. In fact, we present in Section 3.3.3 how these Phase Angle Difference bounds may be tightened based on a Shortest Path algorithm. Then, we can define the angles φ ba = θ ba +θ ba 2 and δ ba = θ ba −θ ba 2 for any (b, a) ∈ E. With this notation, the following constraints are valid for any (b, a) ∈ E s.t. δ ba ≤ π 2 : cos(φ ba )Re(W ba ) + sin(φ ba )Im(W ba ) ≥ R ba cos(δ ba ).
Finally, for every k ∈ T , we require that where R B k B k denotes the matrix (R ba ) (b,a)∈B 2 k and L B k denotes the vector (L b ) b∈B k . Adding the decision vectors L ∈ R n and R ∈ R E to the optimization problem (ACOPF W ) and replacing Constraints ( ) by Constraints (8)- (14), we obtain a Conic Programming problem, that we denote (R).

Proposition 3. The Conic Programming problem (R) is a relaxation of (ACOPF W ).
Proof. We prove the validity of the Constraints (8)- (14). More specifically, we prove that for any couple (S, W ) ∈ C m × H n (E) feasible in (ACOPF W ), the quadruplet (S, W, L, R) is feasible in (R), where L and R are defined as L b = √ W bb and R ba = |W ba | for all (b, a) ∈ E. Since the objective function is the same in (R) and (ACOPF W ), this will prove that (R) is a relaxation of (ACOPF W ). Since , the triplet (R ba , L b , L a ) satisfies the Mc Cormick inequalities [26] with respect to (w.r.t.) these bounds, i.e., Constraints (8)- (9). Constraint (10) is satisfied since W bb ∈ R, as (S, W ) is feasible in (ACOPF W ), yielding R bb = |W bb | = W bb = L 2 b . Constraint (11) also being a Mc Cormick constraint (for b = a), it is satisfied by (R bb , L b ), as R bb = L 2 b . Constraint (12) just follows from the definition of R ba = |W ba |. For any (b, a) ∈ E, we define θ ba = ∠W ba ; considering the definition of φ ba and δ ba , we notice that |θ ba − φ ba | ≤ δ ba . For this reason, if δ ba ≤ π 2 , we obtain cos(|θ ba − φ ba |) ≥ cos(δ ba ), as cos is decreasing over [0, π 2 ]. Using the parity of cos, and multiplying by R ba ≥ 0, we obtain R ba cos(θ ba − φ ba ) ≥ R ba cos(δ ba ) Moreover, R ba cos(θ ba − φ ba ) = |W ba |(cos(φ ba ) cos(θ ba ) + sin(φ ba ) sin(θ ba )) = cos(φ ba )Re(W ba ) + sin(φ ba )Im(W ba ), explaining that (R ba , W ba ) satisfies Constraint (13), whenever δ. Finally, Constraint (14) just follows from the By construction, the relaxation (R) is tighter than the clique-based SDP relaxation, the value of which equals the value of the standard SDP relaxation [15], also known as rank relaxation. The following theorem shows how Constraints (8)-(13) help having |W ba | 2 ≈ W bb W aa when the Voltage Magnitude and Phase Angle Difference intervals are small. We recall the notation ∆ b = v b − v b and that we assume ∆ b ≤ 1 throughout the paper.

Connections to previous works
Previous works in the Power Systems community proposed valid inequalities to strengthen the SDP relaxation of the ACOPF problem [19,20,10]. In [10], the authors show that these valid inequalities are all dominated by the inequalities [10, (36a) and (36b)]. Using the parameter v σ b = v b + v b , the inequalities [10, (36a) and (36b)] read, with our notation, The following Proposition states that Constraints (8)- (13), that we introduce here to strengthen the SDP relaxation, dominate Equations ( †)-( ‡).

Proposition 5. For any
Proof. We take any (b, a) ∈ L ∪ L R and any quadruplet ( First, we combine Equations (19) as Multiplying Constraint (13) by (22), we deduce Equation ( †). We underline that Constraint (13) is indeed applicable since δ ba ≤ π 2 , as (b, a) ∈ L ∩ L R (see Table 1). Second, we combine the Equations (19) Multiplying Equation (23) Multiplying Constraint (13) by The advantage of Constraints (8)-(13) is to enforce a coupling between the convex envelopes of the quadruplets (Re(W ba ), Im(W ba ), W bb , W aa ) involving a same index b. This coupling is realized by the additional decision vectors L and R. In Appendix A, we present an illustrative example of two quadruplets (Re(W ba ), Im(W ba ), W bb , W aa ) and (Re(W bc ), Im(W bc ), W bb , W cc ) satisfying Equations ( †)-( ‡) introduced in [10], but for which there is no vector L and R s.t. Constraints (8)- (13) are satisfied. In this respect, we can state that Constraints (8)-(13) strictly dominate Equations ( †)-( ‡).

Bound Tightening procedures
We use Bound Tightening procedures to reduce the interval lengths ∆ b and δ ba and, thus, reduce the error bound in Theorem 4.

Feasibility-based Bound Tightening (FBBT)
The Power Flow limit for the line (b, a) ∈ L implicitly restricts the phase ∠V b V * a and, consequently, can help reduce the length of the interval [θ ba , θ ba ]. Dividing the inequality . Representing the ratio |V b | |Va| by a variable λ, we can formulate the following Convex Optimization problem Denoting by h its value, we deduce that arcsin(h) is an upper-bound on ∠V b − ∠V a . Hence, we can set θ ba ← min(θ ba , arcsin(h)) without changing the value of (ACOPF). If we minimize Im(u) under the same constraints to get a value h, we can set θ ba ← max(θ ba , arcsin(h)).
Va| to tighten θ ba and θ ba . This type of Bound Tightening is cheap, since it requires to solve a 2-variable optimization problem for each bound.

Optimization-Based Bound Tightening (OBBT)
We also apply a OBBT procedure to the Conic Programming relaxation (R), as performed in [32] with the QCP relaxation. We use any NLP algorithm to find an ACOPF feasible solution. With the corresponding upper-bound denoted obj, we add the constraint g∈G c 1g Re(S g ) + c 2g Re(S g ) 2 ≤ obj to Problem (R). We denote by F the resulting convex feasible set for the tuplet (S, W, L, R). Then, we update the following bounds: • For the Phase Angle Difference on line (b, a) ∈ L, we compute h ba = max (S,W,L,R)∈F Im(W ba ) and h ba = min (S,W,L,R)∈F Im(W ba ) and set

Shortest Path algorithm to tighten Phase Angle Difference bounds
Through FBBT and OBBT, we may individually improve the bounds θ ba and θ ba for any (b, a) ∈ E. To propagate the reduction of the Phase Angle Difference domains, we apply a Shortest Path algorithm. Indeed we notice that, for any s=0 θ bs+1bs . The Shortest Path between b 0 and b t in the directed weighted graph (B, E, θ) helps finding the lowest sum s=0 θ bs+1bs . The Shortest Path between b 0 and b t in the directed weighted graph (B, E, −θ) helps improving the lower-bound on ∠V bt − ∠V b0 to update θ btb0 . To compute Shortest Paths, we apply the Floyd-Warshall algorithm [11], which fits the context of a weighted directed graph, with weights of unspecified sign. May the Floyd-Warshall algorithm find a cycle with negative weight in the directed weighted graph (B, E, θ), it would certify the infeasibility of (ACOPF), since it would give a path b 0 , b 1 , . . . , b t Similarly, finding a cycle of negative weight in (B, E, −θ) certifies the infeasibility of (ACOPF).

A MILP-based Global Optimization algorithm
Leveraging the Conic Programming relaxation (R) and its solution, we generate a sequence of MILP problems whose values converge to the ACOPF value.

Linear Programming Outer-Approximations
The disadvantage of Problem (R) is its computational cost, that is higher, due to SDP constraints, than the cost of a Linear Programming (LP) or a convex QCP relaxation. Hence, it may not be computationally efficient to solve such a relaxation at every node of an exploration tree. The idea of our approach is to solve the relaxation (R) at the root node only, and use it to generate a LP relaxation with the same value. We denote by x ∈ R N the decision vector (Re(S), Im(S), Re(W ), Im(W ), L, R), and we notice that the Problem (R) may also be seen as with P ⊂ R N being a polytope and f 0 (x), f 1 (x), . . . , f M (x) continuous and convex functions. In Appendix B, we detail this polytope, the functions f j (x), and show that they share a common structure: for any j ∈ {0, . . . , M }, there exists an affine application π j : R N → R pj and a compact and convex set U j ⊂ R pj s.t. for all x ∈ P, f j (x) = max u∈Uj u π j (x). For any finite subset q U j ⊂ U j , we define the following polyhedral function q f j (x) = max u∈Ǔ j u π j (x). This function, called "cutting-plane model", is an underestimator of f j . If we relax the formulation (R) by replacing each function f j (x) by its polyhedral underestimator q f j (x) related to a given finite set q U j , we obtain the LP relaxation We show that based on a primal-dual solution of (R), we can compute finite sets q U 0 , . . . , q U M s.t. val(R) = val(R L ). For j ∈ {1, . . . , M } we define K j as the convex cone generated by U j . We define K 0 as the convex cone generated by {1} × U 0 . We also define λ and λ as a priori lower and upper bounds on the value of (R), that may be very rough estimates. We introduce a Lagrangian L function for the conic program (R), defined for any With this definition, we see that the conic program (R) is the min-max problem inf x∈P,λ∈ [λ,λ] Proposition 6. There is no duality gap between Problem (R) and Problem (30), i.e., they share the same value. Moreover, if Problem (30) has an optimal solution κ * ∈ K 0 ×K 1 ×· · ·×K M , written as κ Proof. The absence of duality gap follows from the Sion min-max theorem [22], since The primal optimization set P × [λ, λ] is convex and compact, The dual optimization set K 0 × K 1 × · · · × K M is convex, The Lagrangian L is continuous and convex w.r.t. (x, λ) for any κ ∈ K 0 × K 1 × · · · × K M , and, The Lagrangian L is continuous and concave w.r.t. κ for any (x, λ) ∈ P × [λ, λ].
At the price of finding an optimal primal-dual solution (x * , λ * , κ * ) of (R)-(30), we can build an LP relaxation with the same value. In practice, we obtain such a primal-dual solution for every tested instance.

Partitioning Voltage Magnitude intervals
For any b ∈ B, we may want to split the interval [v b , v b ] in subintervals. We introduce a tree J b and pairs and for any i Moreover we add the following constraint for every j ∈ J b , For every j ∈ J b and for all a ∈ B s.t. (b, a) ∈ E, we add the inequalities

Updating the partitions of the intervals
During the algorithm presented in Section 4.4, the partitions of the intervals [v b , v b ] and [θ ba , θ ba ] are not static but are made dynamically. The partition trees are all initialized as single-node graphs, and are then updated over the course of the algorithm. We present how these trees are updated, at any iteration t of the algorithm, where the current iterate is (S t , W t , L t , R t ). For a given bus b ∈ B, we update the partition tree J b by selecting the active leaf j, i.e. the only leaf j of We create three new leaves j 1 , j 2 , j 3 in the tree, which are attached to node j, and we partition the interval [v bj , v bj ] as follows: For a given pair (b, a) ∈ E, we update the partition tree J ba by selecting the active leaf j, i.e. the only leaf j of J ba s.t. ∠W t ba ∈ [θ baj , θ baj ]. We create three new leaves j 1 , j 2 , j 3 in the tree, which are attached to node j, and we partition the interval [θ baj , θ baj ] as follows: We define θ baj1 = θ baj and θ baj3 = θ baj ; and θ baj2 = θ baj3 = ∠W t ba . The construction procedure of the trees J b and J ba guarantees that (i) v bj − v bj , the length of the interval associated with a node j ∈ J b of depth D(j), is less than ∆ b 2 D(j) , (ii) the coefficient δ baj associated with a node j ∈ J ba is less than δ ba 2 D(j) . Proposition 7. We assume that the convex Constraints (8)- (13) and the MILP Constraints (31)-(38) are satisfied, but with a tolerance ρ ∈ [0, 1] for the nonlinear Constraints (10) and (12). Then, for any nodes j b ∈ J b , j a ∈ J a and j ba ∈ J ba that are active, i.e., s.t. α bj b = α aja = β baj ba = 1, we have We underline that the implication is still valid if δ ba = 0 and log 2 ( 2δ ba π ) = −∞.

The MILP-based iterative scheme
The following Global Optimization algorithm is executed based on (i) a local NLP solver (ii) a Conic Programming solver and (iii) a MILP solver. In this pseudo-code, we use the function (W ) = max (b,a)∈E |W ba | 2 − W bb W aa , which denotes the feasibility error in Constraints ( ). 0. Input: A target optimality gap targetOptGap ≥ 0, a tolerance ≥ 0, integers N 1 , N 2 ∈ N * and a sequence (ρ t ) t∈N with ρ t > 0. 1. Initialization: Compute an ACOPF feasible solution with a NLP solver and denote by obj its value (if the NLP solver fails, obj ← +∞). Solve the Conic Programming relaxation (R). If the gap is greater than targetOptGap, apply FBBT and OBBT to (R). Based on the optimal solution of Problem (R), generate the LP relaxation (R L ) with same value as (R) (see Subsect. 4.1). Set t ← 0 and LB t ← val(R).

3.2.
Solve the resulting MILP problem to global optimality to get x = (Re(S), Im(S), Re(W ), Im(W ), L, R). Theorem 9 states that, if the parameters targetOptGap and are set to zero and if (ρ t ) t∈N vanishes, the algorithm asymptotically recovers global minimizers of (ACOPF W ). Before stating this Theorem, we introduce a preliminary Proposition about the finite termination of the inner-loops.

Proposition 8. For any outer-iteration index t ∈ N * , for any tolerance ρ t > 0, the inner-loop has a finite number of iterations.
Proof. During outer-iteration t ∈ N * and the previous iterations, several auxiliary binary variables and associated linear constraints have been added to the relaxation (R). From the perspective of the decision vector x = (Re(S), Im(S), Re(W ), Im(W ), L, R), this yields a closed nonconvex set X . We also inherit finite sets ( q U j0 ) j∈{0,1,...,M } , the subscript 0 denoting the inner-iteration of index s = 0. The inner-iteration s ∈ N consists in solving to obtain a solution x s , and in defining U j(s+1) = {u js } ∪ U js for some u js ∈ argmax u∈Uj u π j (x s ) for all We reason by contrapositive and assume now that the inner-loop is not terminating in finite time, meaning that the generated MILP relaxation is feasible at each inner-iteration and the stopping condition of the inner-loop is never met. This second point means that ρ t ≤ max j∈{0,...,M } e js for all s ∈ N. We take any j ∈ {0, . . . , M } and show that e js → s 0. Since the sets P ∩ X and U j are compact and since the functions x → f j (x) and (x, u) → u π j (x) are continuous, we deduce that the sequence (e js ) s∈N is bounded. We take any limit point e * of this sequence, and we take a converging subsequence (e jψ(s) ) → s e * . Without loss of generality, since P ∩ X and U j are compact, we can choose the extraction ψ(s) so that x ψ(s) → s x * and u jψ(s) → s u * , for (x * , u * ) ∈ P ∩ X × U j . For s ∈ N * , we notice that ≤ u jψ(s) π j (x ψ(s) ) − u jψ(s−1) π j (x ψ(s) ).

Theorem 9.
If targetOptGap = = 0 and ρ t → t 0, then Either the algorithm stops due to the stopping criterion, and yields a global minimizer of (ACOPF W ), Or the algorithm stops due to the infeasibility of a relaxation, certifying the infeasibility of (ACOPF W ), Or the algorithm generates an infinite sequence of iterates (S t , W t , L t , R t ), and The sequence LB t monotonously converges to val(ACOPF W ) = val(ACOPF), The limit points of the sequence (S t , W t ) t∈N are global minimizers of (ACOPF W ).

Proof.
We consider the first case where the algorithm meets the stopping criterion at the beginning of a certain outer-iteration t. This means that either (a) obj = LB t , proving that the solution (S, V ) found by the NLP solver at step 1. is globally optimal in (ACOPF) and yields (S, V V H ) globally optimal in (ACOPF W ), or (b) the solution (S t , W t , L t , R t ) of the current MILP relaxation of (ACOPF W ) satisfies (W t ) = 0, i.e., (S t , W t ) is in fact feasible in (ACOPF W ) and thus optimal in (ACOPF W ) since it is the optimal solution of a relaxation. The second case is trivial: if the relaxation (R) or any MILP relaxation during the iterations is infeasible, this implies that (ACOPF W ) is also infeasible.
We consider now the third case, where the algorithm does not terminate. We invoke Proposition 8 to claim that for any outer-iteration t ∈ N, the inner-loop terminates in finite time. Hence, there is an infinite number of outer-iterations and we define the infinite sequence x t = (S t , W t , L t , R t ) t∈N , where x t is the solution of the MILP relaxation at the beginning of the outer-iteration t . For any (b, a) We let J t b (resp. J t ba ) denote the state of the tree J b (resp. J ba ) at the beginning of iteration t, and J b (resp. J ba ) the (potentially infinite) limit tree t J t b (resp. t J t ba ). We first show that and j a (t) ∈ J t at the active leaves to which three child nodes are attached during step 2.1 since (b t , a t ) presents the largest violation. For any j ∈ b J b , we recall that D(j) is the depth of j in the (unique) tree J b it belongs to. As x t is the output of the outer-iteration t − 1, it satisfies Constraint (10) and (12) with tolerance ρ t−1 , and we can apply Proposition 7 with ρ = ρ t−1 . This yields 9∆ bt ∆ at 2 max{D(j b (t)),D(ja(t))} + max 9ρ t−1 , (46) We notice that the sequence (j b (t)) t∈N is injective, since each j b (t) is a leaf in J t bt , but not in the trees J s bt for s ≥ t + 1. We deduce that D(j b (t)) → t ∞, otherwise by contrapositive, there would exist M ∈ N s.t. an infinite number of nodes of depth less or equal than M are created in the union of ternary trees b J b ; This is false since the number of nodes with depth less or equal than M is bounded by n M =0 3 . By the same argument, we have D(j a (t)) → t ∞. Combined with (46), we deduce that χ t btat → t 0 since ρ t → t 0 and because ∆ bt , ∆ at are bounded. For any t ∈ N and (b, a) ∈ E, we have 0 ≤ χ t ba ≤ χ t btat by definition of (b t , a t ), implying χ t ba → t 0. We apply the same approach to prove that ξ t ba → t 0 for any (b, a) ∈ E. For t ∈ N, we define ( b t , a t ) ∈ E s.t. ξ t btãt = max ba ξ t ba , and we define j(t) ∈ J t btãt the active leaf to which three child nodes are attached during step 2.2. We also define D( j(t)) as the depth of j(t) in J t btãt , which satisfies D( j(t)) → t ∞ by injectivity of ( j(t)) n∈N and since the number of nodes in (b,a)∈E J ba with depth less or equal than M is bounded by |E| M =0 3 . As D( j(t)) → t ∞, we know that it exists t 0 ∈ N s.t. D( j(t)) ≥ 2 for all t ≥ t 0 . Hence, for all t ≥ 0, D( j(t)) ≥ log 2 ( 4π π ) ≥ log 2 2δb tãt π , since δ btat ∈ [0, 2π]. Applying Proposition 7, we know that for any t ≥ t 0 , Combined with ρ t → t 0 and D( j(t)) → t ∞, we deduce that ξ t btãt → t 0. Additionally, since ξ t btãt = max ba ξ t ba , we have 0 ≤ ξ t ba ≤ ξ t btãt and, thus, ξ t ba → t 0 for any (b, a) ∈ E. We deduce that (W t ) → t 0, since (b,a)∈E χ t ba + ξ t ba , due to the triangle inequality. Hence, for any limit point (S, W ) of (S t , W t ), we thus have (W ) = 0. As ρ t → t 0, this also proves that (S, W ) satisfies the nonlinear convex constraints in (R). Hence, (S, W ) is feasible in (ACOPF W ). We denote by q f t 0 the cutting-plane model of the objective function at the beginning of iteration t; this function only depends on S, hence, we write q f t 0 (S) instead of q f t 0 (x). As the successive MILP relaxations over the iterations have nonincreasing feasible sets w.r.t. variables (S, W, L, R) and have nondecreasing sequence q f t 0 (S) as objective functions, the sequence q f t 0 (S t ) = LB t is nondecreasing. It is also bounded above by val(ACOPF W ) and, thus, converges to a value v * ≤ val(ACOPF W ). Since ρ t → t 0, q f t 0 (S t ) → t f 0 (S) = g∈G c 1g Re(S g ) + c 2g Re(S g ) 2 , for any limit point (S, W ). By uniqueness of the limit of q f t 0 (S t ) and since (S, W ) is feasible in (ACOPF W ), we deduce that v * = g∈G c 1g Re(S g ) + c 2g Re(S g ) 2 ≥ val(ACOPF W ). We conclude that v * = val(ACOPF W ) = val(ACOPF) and that (S, W ) is optimal in (ACOPF).

5
Experimental evaluation

Experimental setting
For all experiments, we use a 64-bit Ubuntu computer with 32 Intel(R) Xeon(R) CPU E5-2620 v4 @ 2.10GHz and 64 GB RAM. Along our algorithm, we use the commercial solvers MOSEK [1] and CPLEX [18] called through their Python APIs, as well as the academic solver IPOPT [33] called through the Pyomo interface [16]. We compute the tree decomposition with the approximate minimum degree (AMD) ordering routine of the chompack package. We consider a relative optimality gap of 0.01% for global optimality (GOPT) and use the parameters (N 1 , N 2 , ) = (4, 4, 10 −6 ). The FBBT and OBBT procedures are applied for all variables a maximum of 4 times, and with a time limit of 10 hours (TL 1 ). After each pass of FBBT and OBBT, we apply the Floyd-Warshall algorithm and we check whether the gap of the tightened conic relaxation reaches the target optimality gap.
If the maximum number of Bound Tightening iterations or time limit is reached, we enter the MILP iterative scheme with a time limit of 5 hours (TL 2 ). Our code is available at github.com/aoustry/SDP-MILP4OPF. This study focuses on the network instances from the IEEE PES PGLib AC-OPF v21.07 library [2] with less than 500 buses. As shown in Table 2, the instances of this benchmark are split in three categories depending on their characteristics: Typical Operating Condition (TYP) instances correspond to a reference scenario, Congested Operating Condition (API) correspond to situations with greater Power Demands, and Small Angle Condition (SAD) correspond to tighter constraints for the Phase Angle Difference.
We compare our approach with the standard SOCP and SDP relaxations [28], and with two other Global Optimization approaches [13,32]. We performed these comparative experiments on the same PGLib v21.07 instances, with the same hardware, the same time limit (TL 1 ) as our OBBT algorithm, and the same relative optimality gap tolerance. The concurrent approach from [32] is an OBBT algorithm based on a strengthened QC relaxation. We ran the Julia implementation of this algorithm provided in the PowerModels.jl package [8]. The QC relaxations are solved with IPOPT. The concurrent approach from [13] consists of an OBBT algorithm, based on a Determinant SDP relaxation strengthened with RLT constraints. We ran a C++ implementation of this algorithm, that is based on the Mathematical Modeling Language Gravity [17], and is available at the link indicated in [13]; The corresponding relaxations are also solved with IPOPT. We point out that these competing approaches [13,32] solely rely on open-source tools, whereas our implementation uses commercial solvers. Table 2 presents the optimality gap (in %) and the computational times obtained by the several approaches for the considered list of instances. For lack of space, we do not detail the computation time of the SOCP and SDP relaxations. To give an idea, the computation time of the SOCP relaxation is below 2s for all instances; for the SDP relaxation, the computation time goes from 0.2s for the smallest cases to about 40s for the largest cases. As regards the column "This work", in the optimality gap section of the table: the entry (GOPT) means that the Bound Tightening procedure based on the Conic Programming relaxation closed the gap; else, the entry a/b represents the gap after the Bound Tightening procedure (a) and the gap after the iterative MILP scheme (b).

Numerical results
This table shows that our algorithm reaches Global Optimality for 35 instances over 51. The optimality gap is below 0.5% for 43 instances over 51. For 4 instances only, the optimality gap at the end of the algorithm is above 2%. As regards the instances with less than 57 buses, they are all solved to Global Optimality in less than 220 seconds. For all these instances with less than 57 buses, except case5_pjm and case30_as_api, the Bound Tightening procedure based on the Conic Programming relaxation (R) manages to close the gap. For the instances case5_pjm and case30_as_api, the gap is closed by the iterative MILP scheme. For all the instances with more than 57 nodes where the MILP scheme is executed, the gap is admittedly not closed within the time limit, but it is reduced, except for case500_goc_api. For 44 over these 51 instances, our algorithm has the lowest gap; for 14 instances over 51, it has a strictly lower gap than the others approaches. For 5 instances (among those 14 instances), our approach is the only one to reach Global Optimality. Regarding the 7 instances where our approach has not the best gap: for 3 instances, only the QC relaxation-based Bound Tightening algorithm [32] yields a strictly lower gap than our approach; for the 4 other instances, only the Determinant-SDP relaxation-based Bound Tightening algorithm [13] yields a strictly lower gap than our approach.

Conclusion and perspectives
We introduce a Conic Programming relaxation for the AC Optimal Power Flow problem. This relaxation is a tightening of the classical Semidefinite Programming relaxation with additional variables and valid inequalities. These inequalities dominate previously introduced nonlinear cuts, used to strengthen convex relaxations. Our numerical experiments on a reference benchmark illustrate that this Conic Programming relaxation is particularly suitable for a Bound Tightening procedure: it closes the gap in many cases where a Bound Tightening based on a Quadratic Convex relaxation does not. We also introduce an iterative scheme based on Mixed-Integer Linear Programming, that converges asymptotically towards global minimizers of the AC Optimal Power Flow problem.
For the instances where the Bound Tightening procedure does not close the gap, this iterative scheme is able to reduce significantly the gap in most of the cases. A future line of research will consist in improving the scalability of the Optimization-Based Bound Tightening: parallelizing this procedure, or targeting the bounds to tighten, based on the graph structure. Another avenue to explore is the possibility of speeding-up the solution of the Mixed-Integer Linear Programming problem at a given iteration, by reusing the Branch-and-Bound trees of the problems solved during the previous iterations.
s.t. Constraints (8)-(13) are satisfied. Since W aa = v 2 a , we deduce from Constraint (11) that v 2 a + v a v a = R aa + v a v a ≤ (v a + v a )L a , i.e., that (v a + v a )v a ≤ (v a + v a )L a , and, thus, v a ≤ L a since (v a + v a ) > 0.
As v a ≥ L a by definition of L a , we observe that L a = v a . We use then R ba ≤ v a L b + v b L a − v a v b from Constraint (9) to deduce that R ba ≤ v a L b . Constraint (12) gives |W ba | ≤ R ba , meaning that L b ≥ |W ba | va = 1. As L 2 b ≤ W bb = 1 according to Constraint (10), L b = 1. As we did for a, we deduce from W aa = v 2 a and Constraint (11) Hence, R bc = v c L b = 1.1. As Constraint (13) imposes cos(φ bc )Re(W bc ) + sin(φ bc )Im(W bc ) ≥ R bc cos(δ bc ), we deduce that (φ bc )Re(W bc ) + sin(φ bc )Im(W bc ) ≥ 1.089. This is contradictory with the fact that (φ bc )Re(W bc )+sin(φ bc )Im(W bc ) = Re(W bc ) = 1.085. As a conclusion, there does (8)- (13) are satisfied simultaneously for the pairs (b, a) and (b, c).
This illustrates the interest of setting the trigonometric cuts (13) with a variable radius R ba , whereas previous works, to our knowledge, only use cuts associated to an extreme value of R ba .

B Nonlinear but convex objective and constraints in relaxation (R)
We recall that the decision vector in relaxation (R) is x = (Re(S), Im(S), Re(W ), Im(W ), L, R). First, we denote by P ⊂ R N the polytope defined by the following box constraints: For all g ∈ G, Re(S g ) ∈ [Re(s g ), Re(s g )] and Im(S g ) ∈ [Im(s g ), Im(s g )], For all (b, a) ∈ E, Now, we review the nonlinear terms in the objective and in the constraints of relaxation (R), as functions of x. We show that all these functions have the form f (x) = max u∈U u π(x) for all x ∈ P, with a given affine application π : R N → R p and a compact and convex set U ⊂ R p .