{'fraction_non_alphanumeric': 0.04210991508570346, 'fraction_numerical': 0.0156886265349285, 'mean_word_length': 4.200262639527248, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Distributed Algorithm for Optimal Power Flow on Unbalanced Multiphase Distribution Networks

The optimal power flow (OPF) problem is funda- mental in power distribution networks control and operation that underlies many important applications such as volt/var control and demand response, etc.. Large-scale highly volatile renewable penetration in the distribution networks calls for real-time feed- back control, and hence the need for distributed solutions for the OPF problem. Distribution networks are inherently unbalanced and most of the existing distributed solutions for balanced networks do not apply. In this paper we propose a solution for unbalanced distribution networks. Our distributed algorithm is based on alternating direction method of multiplier (ADMM). Unlike existing approaches that require to solve semidefinite programming problems in each ADMM macro-iteration, we exploit the problem structures and decompose the OPF problem in such a way that the subproblems in each ADMM macro- iteration reduce to either a closed form solution or eigen-decomposition of a 6x6 hermitian matrix, which significantly reduce the convergence time. We present simulations on IEEE 13, 34, 37 and 123 bus unbalanced distribution networks to illustrate the scalability and optimality of the proposed algorithm.


I. INTRODUCTION
The optimal power flow (OPF) problem seeks to minimize a certain objective, such as power loss and generation cost subject to power flow physical laws and operational constraints. It is a fundamental problem that underpins many distribution system operations and planning problems such as economic dispatch, unit commitment, state estimation, volt/var control and demand response, etc.. Most algorithms proposed in the literature are centralized and meant for applications in today's energy management systems that centrally schedule a relatively small number of generators. The increasing penetrations of highly volatile renewable energy sources in distribution systems requires simultaneously optimizing (possibly in realtime) the operation of a large number of intelligent endpoints. A centralized approach will not scale because of its computation and communication overhead and we need to rely on distributed solutions.
Various distributed algorithms for OPF problem have been proposed in the literature. Some early distributed algorithms, including [1], [2], do not deal with the non-convexity issue of OPF and convergence is not guaranteed for those algorithms. Recently, convex relaxation has been applied to convexify the OPF problem, e.g. semi-definite programming (SDP) relaxation [3]- [6] and second order cone programming (SOCP) relaxation [7]- [9]. When an optimal solution of the original OPF problem can be recovered from any optimal solution of the SOCP/SDP relaxation, we say the relaxation is exact. It is shown that both SOCP and SDP relaxations are exact for radial networks using standard IEEE test networks and many practical networks [4], [6], [8], [10]. This is important because almost all distribution systems are radial. Thus, optimization decompositions can be applied to the relaxed OPF problem with guaranteed convergence, e.g. dual decomposition method [11], [12], methods of multiplier [13], [14], and alternating direction method of multiplier (ADMM) [10], [15], [16].
There are at least two challenges in designing distributed algorithm that solves the OPF problem on distribution systems. First, distribution systems are inherently unbalanced because of the unequal loads on each phase [17]. Most of the existing approaches [11]- [16], [18] are designed for balanced networks and do not apply to unbalanced networks.
Second, the convexified OPF problem on unbalanced networks consists of semi-definite constraints. To our best knowledge, all the existing distributed solutions [1], [2], [10] require solving SDPs within each macro-iteration. The SDPs are computationally intensive to solve, and those existing algorithms take significant long time to converge even for moderate size networks.
In this paper, we address those two challenges through developing an efficient distributed algorithm for the OPF problems on unbalanced networks based on alternating direction methods of multiplier(ADMM). The advantages of the proposed algorithm are twofold: 1) instead of relying on SDP optimization solver to solve the optimization subproblems in each iteration as existing approaches, we exploit the problem structures and decompose the problem in such a way that the subproblems in each ADMM macro-iteration reduce to either a closed form or a eigen-decomposition of a 6 × 6 hermitian matrix, which greatly speed up the convergence time. 2) Communication is only required between adjacent buses.
We demonstrate the scalability of the proposed algorithms using standard IEEE test networks [19]. The proposed algorithm converges within 3 seconds on the IEEE-13, 34, 37, 123 bus systems. To show the superiority of using the proposed procedure to solve each subproblem, we also compare the computation time for solving a subproblem by our algorithm and an off-the-shelf SDP optimization solver (CVX, [20]). Our solver requires on average 3.8 × 10 −3 s while CVX requires on average 0.58s.
A preliminary version has appeared in [21]. In this paper, we improve the algorithm in [21] in the following aspects: 1) We consider more general forms of objective function and power injection region such that the algorithm can be used in more applications. In particular, we provide a sufficient condition, which holds in practice, for the existence of efficient solutions to the optimization subproblems. 2) Voltage magnitude constraints, which are crucial to distribution system operations, are considered in the new algorithm. 3) We study the impact of network topologies on the rate of convergence in the simulations.
The rest of the paper is structured as follows. The OPF problem on an unbalanced network is defined in section II. In section III, we develop our distributed algorithm based on ADMM. In section IV, we test its scalability on standard IEEE distribution systems and study the impact of network topologies on the rate of convergence. We conclude this paper in section V.

II. PROBLEM FORMULATION
In this section, we define the optimal power flow (OPF) problem on unbalanced radial distribution networks and review how to solve it through SDP relaxation.
We denote the set of complex numbers with C, the set of ndimensional complex numbers with C n and the set of m × n complex matrix with C m×n . The set of hermitian (positive semidefinite) matrix is denoted by S (S + ). The hermitian transpose of a vector (matrix) x is denoted by x H .
The trace of a square matrix x ∈ C n×n is denoted by tr(x) := n i=1 x ii . The inner product of two matrices (vectors) x, y ∈ C m×n is denoted by x, y := Re(tr(x H y)). The Frobenius (Euclidean) norm of a matrix (vector) x ∈ C m×n is defined as x 2 := x, x . Given x ∈ C n×n , let diag(x) ∈ C n×1 denote the vector composed of x's diagonal elements.

A. Branch flow model
We model a distribution network by a directed tree graph T := (N , E) where N := {0, . . . , n} represents the set of buses and E represents the set of distribution lines connecting the buses in N . Index the root of the tree by 0 and let N + := N \ {0} denote the other buses. For each bus i, it has a unique ancestor A i and a set of children buses, denoted by C i . We adopt the graph orientation where every line points towards the root. Each directed line connects a bus i and its unique ancestor A i . We hence label the lines by E := {1, . . . , n} where each i ∈ E denotes a line from i to A i . Note that E = N + and we will use N + to represent the lines set for convenience.
Let a, b, c denote the three phases of the network. For each bus i ∈ N , let Φ i ⊆ {a, b, c} denote the set of phases. In a typical distribution network, the set of phases for bus i is a subset of the phases of its parent and superset of the phases of its children, i.e.
For each line i ∈ N + connecting bus i and its ancestor A i , the set of phases is Some notations are summarized in Fig. 1. A variable without a subscript denotes the set of variables with appropriate components, as summarized below.
Branch flow model is first proposed in [22], [23] for balanced radial networks. It has better numerical stability than bus injection model and has been advocated for the design and operation for radial distribution network, [8], [14], [18], [24]. In [6], it is first generalized to unbalanced radial networks and uses a set of variables (v, s, , S). Given a radial network T , the branch flow model for unbalanced network is defined by: where P i (v Ai ) denote projecting v Ai to the set of phases on bus i and P i (S j − z j j ) denote lifting the result of S j − z j j to the set of phases Φ i and filling the missing phase with 0, Given a vector (v, s, , S) that satisfies (1), it is proved in [6] that the bus voltages V i and branch currents I i can be uniquely determined if the network is a tree. Hence this model (1) is equivalent to a full unbalanced AC power flow model. See [6, Section IV] for details.

B. OPF and SDP relaxation
The OPF problem seeks to optimize certain objective, e.g. total power loss or generation cost, subject to unbalanced power flow equations (1) and various operational constraints. We consider an objective function of the following form: For instance, • to minimize total line loss, we can set for each φ ∈ Φ i , i ∈ N , • to minimize generation cost, we can set for each i ∈ N , where α φ i , β φ i > 0 depend on the load type on bus i, e.g. α φ i = 0 and β φ i = 0 for bus i where there is no generator and for generator bus i, the corresponding α φ i , β φ i depends on the characteristic of the generator.
For each bus i ∈ N , there are two operational constraints on each phase φ ∈ Φ i . First, the power injection s φ i is constrained to be in a injection region I φ i , i.e.
The feasible power injection region I φ i is determined by the controllable loads attached to phase φ on bus i. Some common controllable loads are: • For controllable load, whose real power can vary within [p i , p i ] and reactive power can vary within [q i , q i ], the injection region I i is For instance, the power injection of each phase φ on substation bus 0 is unconstrained, thus p i , q i = −∞ and p i , q i = ∞. • For solar panel connecting the grid through a inverter with nameplate s φ i , the injection region I i is Second, the voltage magnitude needs to be maintained within a prescribed region. Note that the diagonal element of v i describes the voltage magnitude square on each phase φ ∈ Φ i . Thus the constraints can be written as where v φφ i denotes the φ th diagonal element of v i . Typically the voltage magnitude at substation buses is assumed to be fixed at a prescribed value, i.e. v φ 0 = v φ 0 for φ ∈ Φ 0 . At other load buses i ∈ N + , the voltage magnitude is typically allowed to deviate by 5% from its nominal value, i.e. v φ To summarize, the OPF problem for unbalanced radial distribution networks is: (1) and (5) − (7) The OPF problem (8) is nonconvex due to the rank constraint (1d). In [6], an SDP relaxation for (8) is obtained by removing the rank constraint (1d), resulting in a semidefinite program (SDP): Clearly the relaxation ROPF (9) provides a lower bound for the original OPF problem (8) since the original feasible set is enlarged. The relaxation is called exact if every optimal solution of ROPF satisfies the rank constraint (1d) and hence is also optimal for the original OPF problem. It is shown empirically in [6] that the relaxation is exact for all the tested distribution networks, including IEEE test networks [19] and some real distribution feeders.

III. DISTRIBUTED ALGORITHM
We assume SDP relaxation is exact and develop in this section a distributed algorithm that solves the ROPF problem. We first design a distributed algorithm for a broad class of optimization problem through alternating direction method of multipliers (ADMM). We then apply the proposed algorithm on the ROPF problem, and show that the optimization subproblems can be solved efficiently either through closed form solutions or eigen-decomposition of a 6 × 6 matrix.

A. Preliminary: ADMM
ADMM blends the decomposability of dual decomposition with the superior convergence properties of the method of multipliers [25]. It solves optimization problem of the form 1 : where f (x), g(y) are convex functions and K x , K y are convex sets. Let λ denote the Lagrange multiplier for the constraint x = y. Then the augmented Lagrangian is defined as where ρ ≥ 0 is a constant. When ρ = 0, the augmented Lagrangian degenerates to the standard Lagrangian. At each iteration k, ADMM consists of the iterations: Specifically, at each iteration, ADMM first updates x based on (12a), then updates y based on (12b), and after that updates the multiplier λ based on (12c). Compared to dual decomposition, ADMM is guaranteed to converge to an optimal solution under less restrictive conditions. Let which can be viewed as the residuals for primal and dual feasibility, respectively. They converge to 0 at optimality and are usually used as metrics of convergence in the experiment.
Interested readers may refer to [25,Chapter 3] for details.
In this paper, we generalize the above standard ADMM [25] such that the optimization subproblems can be solved efficiently for our ROPF problem. Instead of using the quadratic penalty term ρ 2 x − y 2 2 in (11), we will use a more general quadratic penalty term: and Λ is a positive diagonal matrix. Then the augmented Lagrangian becomes The convergence result in [25,Chapter 3] carries over directly to this general case.

B. ADMM based Distributed Algorithm
In this section, we will design an ADMM based distributed algorithm for a broad class of optimization problem, of which the ROPF problem is a special case. Consider the following optimization problem: are matrices with appropriate dimensions. A broad class of graphical optimization problems (including ROPF) can be formulated as (15). Specifically, each node i ∈ N is associated with some local variables stacked as x i , which belongs to an intersection of R i + 1 local feasible sets K ir and has a cost objective function f i (x i ). Variables in node i are coupled with variables from their neighbor nodes in N i through linear constraints (15c). The objective then is to solve a minimal total cost across all the nodes. The goal is to develop a distributed algorithm that solves (15) such that each node i solve its own subproblem and only exchange information with its neighbor nodes N i . In order to transform (15) into the form of standard ADMM (10), we need to have two sets of variables x and y. We introduce two sets of slack variables as below: 1) x ir . It represents a copy of the original variable x i for 1 ≤ r ≤ R i . For convenience, denote the original x i by x i0 . 2) y ij . It represents the variables in node i observed at node j, for j ∈ N i .
Then (15) can be reformulated as where x and y represent the two groups of variables in standard ADMM. Note that the consensus constraints (16d) and (16e) force all the duplicates x ir and y ij are the same. Thus its solution x i0 is also optimal to the original problem (15). (16) falls into the general ADMM form (10), where (16b) corresponds to K y , (16c) corresponds to K x , and (16d) and (16e) are the consensus constraints that relates x and y.
Following the ADMM procedure, we relax the consensus constraints (16d) and (16e), whose Lagrangian multipliers are denoted by λ ir and µ ij , respectively. The generalized augmented Lagrangian then can be written as where the parameter Λ ir and M ij depend on the problem we will show how to design them in section III-C. Next, we show that both the x-update (12a) and y-update (12b) can be solved in a distributed manner, i.e. both of them can be decomposed into local subproblems that can be solved in parallel by each node i with only neighborhood communications.
First, we define the set of local variables for each node i, denoted by A i , which includes its own duplicates x ir and the associated multiplier λ ir for 0 ≤ r ≤ R i , and the "observations" y ji of variables from its neighbor N i and the associated multiplier µ ji , i.e.
In the x-update at each iteration k, the optimization subproblem that updates x k+1 is where the constraint K x is the Cartesian product of K ir , i.e.
The objective can be written as a sum of local objectives as shown below where the last term is independent of x and Then the problem (19) in the x-update can be written explicitly as where both the objective and constraint are separable for 0 ≤ r ≤ R i and i ∈ N . Thus it can be decomposed into i∈N (R i + 1) independent problems that can be solved in parallel. There are R i + 1 problems associated with each node i and the r th (0 ≤ r ≤ R i ) one can be simply written as whose solution is the new update of variables x ir for node i.
In the above problem, the constants y k ij , µ k ij ∈ A j are not local to i and stored in i's neighbors j ∈ N i . Therefore, each node i needs to collect (y ij , µ ij ) from all of its neighbors prior to solving (22). The message exchanges is illustrated in Figure  2a.
In the y-update, the optimization problem that updates y k+1 is where the constraint set K y can be represented as a Cartesian product of |N | disjoint sets, i.e. The objective can be written as a sum of local objectives as below.
where the last term is independent of y and Then the problem (23) in the y-update can be written explicitly as which can be decomposed into |N | subproblems and the subproblem for node i is s.t.

j∈Ni
A ij y ji = 0, whose solution is the new update of {y ji | j ∈ N i } ∈ A i . In (24), the constants x j0 ∈ A j are stored in i's neighbor j ∈ N i . Hence, each node i needs to collect x j0 from all of its neighbor prior to solving (24). The message exchanges in the y-update is illustrated in Figure 2b. The problem (24) can be solved with closed form solution. we stack the real and imaginary part of the variables {y ji | j ∈ N i } in a vector with appropriate dimensions and denote it asỹ. Then (24) takes the following form: where M is a positive diagonal matrix,Ã is a full row rank real matrix, and c is a real vector. M, c, A are derived from (24). There exists a closed form expression for (25) given bỹ In summary, the original problem (15) is decomposed into local subproblems that can be solved in a distributed manner using ADMM. At each iteration, each node i solves (22) in the x-update and (25) in the y-update. There exists a closed form solution to the subproblem (25) in the y-update as shown in (26), and hence whether the original problem (15) can be solved efficiently in a distributed manner depends on the existence of efficient solutions to the subproblems (22) in the x-update, which depends on the realization of both the objectives f i (x i ) and the constraint sets K ir .
Next, we show the ROPF problem (9) is a special case of (15), hence can be solved in a distributed manner using the above method. In particular, we show the corresponding subproblems in the x-update can be solved efficiently.

C. Application on OPF problem
We assume the SDP relaxation is exact and now derive a distributed algorithm for solving ROPF (9). Using the ADMM based algorithm developed in Section III-B, the global ROPF problem is decomposed into local subproblems that can be solved in a distributed manner with only neighborhood communication. Note that the subproblems in the y-update for each node i can always been solved with closed form solution, we only need to develop an efficient solution for the subproblems (22) in the x-update for the ROPF problem. In particular, we provide a sufficient condition, which holds in practice, for the existence of efficient solutions to all the optimization subproblems. Compared with existing methods, e.g. [10], [13]- [16], that use generic iterative optimization solver to solve each subproblem, the computation time is improved by more than 100 times.
The ROPF problem defined in (9) can be written explicitly as over v, s, S, Then (27) takes the form of (15) with R i = 1 for all i ∈ N , where (27b)-(27c) correspond to (15c) and (27d)-(27f) correspond to (15d). Then we have the following theorem, which provides a sufficient condition for the existence of an efficient solution to (22). Theorem 3.1: Suppose there exists a closed form solution to the following optimization problem for all i ∈ N and φ ∈ Φ i given any constantŝ φ and ρ, then the subproblems for ROPF in the x-update (22) can be solved via either closed form solutions or eigen-decomposition of a 6 × 6 hermitian matrix. Proof: We will prove Theorem 3.1 through elaborating the procedure to solve (22).
Recall that there is always a closed form solution to the optimization subproblem (24) in the y-update, if the objective function f φ i s φ and injection region I φ i satisfy the sufficient condition in Theorem 3.1, all the subproblems can be solved efficiently.
Remark 3.1: In practice, the objective function f φ i (s), usually takes the form of f φ i (s) := αi 2 p 2 + β i p, which models both line loss and generation cost minimization as discussed in Section II-B. For the injection region I φ i , it usually takes either (3) or (4). It is shown in Appendix B that there exist closed form solution for all of those cases. Thus (31) can be solved efficiently for practical applications.
Following the procedure in Section III-B, we introduce two set of slack variables: x ir and y ij . Then the counterpart of (16) is where we put superscript (·) (x) and (·) (y) on each variable to denote whether the variable is updated in the x-update or y-update step. Note that each node i does not need full information of its neighbor. Specifically, for each node i, only voltage information v (y) Aii is needed from its parent A i and branch power S ii : children j ∈ C i based on (32). Thus, y ij contains only partial information about x i0 , i.e.
On the other hand, only x i0 needs to hold all the variables and it suffices for x i1 to only have a duplicate of v i , i.e.
As a result, x ir , y ii in (32g) and x i0 , y ij in (32h) do not consist of the same components. Here, we abuse notations in both (32g) and (32h), which are composed of components that appear in both items, i.e.
x i0 − y ij ii ) r = 1 . Let λ denote the Lagrangian multiplier for (32g) and µ the Lagrangian multiplier for (32g). The detailed mapping between constraints and those multipliers are illustrated in Table I. Next, we will derive the efficient solution for the subproblems in the x-update. For notational convenience, we will skip the iteration number k on the variables. In the x-update, there are 2 subproblems (22) associated with each bus i. The first problem, which updates x i0 , can be written explicitly as: where H i0 (x i0 ) is defined in (20) and for our application, Mij is chosen to be i0 ), which is defined below, can be written as the Euclidean distance of two Hermitian matrix, which is one of the key reasons that lead to our efficient solution. Therefore, H i0 (x i0 ) can be further decomposed as The last step in (35) is obtained using square completion and the variables labeled with hat are some constants. Hence, the objective (33a) in (33) can be decomposed into two parts, where the first part H (1) i0 . Note that the constraint (33c)-(33d) can also be separated into two independent constraints. Variables i0 ) only depend on (33c) and s (x) i0 depends on (33d). Then (33) can be decomposed into two subproblems, where the first one (36) solves the optimal (S ) and the second one (37) solves the optimal s (x) i0 . The first subproblem can be written explicitly as which can be solved using eigen-decomposition of a 6 × 6 matrix via the following theorem. Theorem 3.2: Suppose W ∈ S n and denote X(W ) := arg min X∈S+ X − W 2 2 . Then X(W ) = i:λi>0 λ i u i u H i , where λ i , u i are the i th eigenvalue and orthonormal eigenvector of matrix W , respectively.
Proof: The proof is in Appendix A. Denote .
Then (36) can be written abbreviately as which can be solved efficiently using eigen-decomposition based on Theorem 3.2. The second problem is Recall that if f i (s , then both the objective and constraint are separable for each phase φ ∈ Φ i . Therefore, (37) can be further decomposed into |Φ i | number of subproblems as below.
which takes the same form as of (31) in Theorem 3.1 and thus can be solved with closed form solution based on the assumptions. For the problem (39) that updates x i1 , which consists of only one component v i1 , it can be written explicitly as min (20) and for our application, Λir is chosen to be Then the closed form solution is given as: To summarize, the subproblems in the x-update for each bus i can be solved either through a closed form solution or a eigendecomposition of a 6 × 6 matrix, which proves Theorem 3.1.
In the y-update, the subproblem solved by each node i takes the form of (24) and can be written explicitly as which has a closed form solution given in (26) and we do not reiterate here. Finally, we specify the initialization and stopping criteria for the algorithm. Similar to the algorithm for balanced networks, a good initialization usually reduces the number of iterations for convergence. We use the following initialization suggested by our empirical results. We first initialize the auxiliary variables {V i | i ∈ N } and {I i | i ∈ E}, which represent the complex nodal voltage and branch current, respectively. Then we use these auxiliary variables to initialize the variables in (32). Intuitively, the above initialization procedure can be interpreted as finding a solution assuming zero impedance on all the lines. The procedure is formally stated in Algorithm 1.
For the stopping criteria, there is no general rule for ADMM based algorithm and it usually hinges on the particular problem. In [25], it is argued that a reasonable stopping criteria is that both the primal residual r k defined in (13a) and the dual residual s k defined in (13b) are below 10 −4 |N |. We adopt this criteria and the empirical results show that the solution is accurate enough. The pseudo code for the complete algorithm is summarized in Algorithm 2.
Algorithm 1 Initialization of the Algorithm 3 π for i ∈ N 2: Initialize s φ i using any point in the injection region i0 = s i for i ∈ N 5: y ij = x i0 for j ∈ N i and i ∈ N 6: x i1 = x i0 for i ∈ N 7: function DFS(i,φ) 8: for j ∈ C i do 10: In the x-update, each agent i solves both (33) and (39) to update x i0 and x i1 . 6: In the y-update, each agent i solves (40) to update y ji for j ∈ N i . 7: In the multiplier update, update λ, µ by (12c). 8: end while

IV. CASE STUDY
In this section, we first demonstrate the scalability of the distributed algorithm proposed in section III-C by testing it on the standard IEEE test feeders [19]. To show the efficiency of the proposed algorithm, we also compare the computation time of solving the subproblems in both the x and y-update with off-the-shelf solver (CVX). Second, we run the proposed algorithm on networks of different topology to understand the factors that affect the convergence rate. The algorithm is implemented in Python and run on a Macbook pro 2014 with i5 dual core processor.

A. Simulations on IEEE test feeders
We test the proposed algorithm on the IEEE 13, 34, 37, 123 bus distribution systems. All the networks have unbalanced three phase. The substation is modeled as a fixed voltage bus (1 p.u.) with infinite power injection capability. The other buses are modeled as load buses whose voltage magnitude at each phase can vary within [0.95, 1.05] p.u. and power injections are specified in the test feeder. There is no controllable device in the original IEEE test feeders, and hence the OPF problem degenerates to a power flow problem, which is easy solve. To demonstrate the effectiveness of the algorithm, we replace all the capacitors with inverters, whose reactive power injection ranges from 0 to the maximum ratings specified by the original capacitors. The objective is to minimize power loss across the network, namely f φ i (s φ i ) = p φ i for φ ∈ Φ i and i ∈ N . We mainly focus on the time of convergence (ToC) for the proposed distributed algorithm. The algorithm is run on a single machine. To roughly estimate the ToC (excluding communication overhead) if the algorithm is run on multiple machines in a distributed manner, we divide the total time by the number of buses.
In Table II, we record the number of iterations to converge, total computation time required to run on a single machine and the average time required for each node if the algorithm is run on multiple machines excluding communication overhead. From the simulation results, the proposed algorithm converges within 2.5 second for all the standard IEEE test networks if the algorithm is run in a distributed manner.
Moreover, we show the advantage of using the proposed algorithm by comparing the computation time of solving the subproblems between off-the-shelf solver (CVX [20]) and our algorithm. In particular, we compare the average computation time of solving the subproblem in both the x and y update.
In the x-update, the average time required to solve the subproblem (40) is 9.8 × 10 −5 s for our algorithm but 0.13s for CVX. In the y-update, the average time required to solve the subproblems (33)-(39) are 3.7 × 10 −3 s for our algorithm but 0.45s for CVX. Thus, each ADMM iteration takes about 3.8×10 −3 s for our algorithm but 5.8×10 −1 s for using iterative algorithm, a more than 100x speedup.

B. Impact of Network Topology
In section IV-A, we demonstrate that the proposed distributed algorithm can dramatically reduce the computation time within each iteration. The time of convergence (ToC) is determined by both the computation time required within each iteration and the number of iterations. In this subsection, we study the number of iterations, namely rate of convergence.
Rate of convergence is determined by many different factors. Here, we only consider the rate of convergence from two factors, network size N , and diameter D, i.e. given the termination criteria in Algorithm 2, the impact of network size and diameter on the number of iterations. The impact from other factors, e.g. form of objective function and constraints, is beyond the scope of this paper.
To illustrate the impact of network size N and diameter D on the rate of convergence, we simulate the algorithm on two extreme cases: 1) Line network in Fig. 3a, whose diameter is  the largest given the network size, and 2) Fat tree network in Fig. 3b, whose diameter is the smallest given the network size.
In Table III, we record the number of iterations for both line and fat tree network of different sizes. For the line network, the number of iterations increases notably as the size increases. For the fat tree network, the trend is less obvious compared to line network. It means that the network diameter has a stronger impact than the network size on the rate of convergence.

V. CONCLUSION
In this paper, we have developed a distributed algorithm for optimal power flow problem on unbalanced distribution system based on alternating direction method of multiplier. We have derived an efficient solution for the subproblem solved by each agent thus significantly reducing the computation time. Preliminary simulation shows that the algorithm is scalable to all IEEE test distribution systems.