Distributed Reactive Power Regulation Considering Load Voltage and Power Loss

In this study, optimal reactive power regulation in distribution networks is achieved through the use of distributed reactive power regulators that can 1) perceive their own voltage magnitude and the P/Q flows in the connected branches, 2) communicate with nearby regulators, and 3) adjust the reactive power injections into the grid to minimize system power losses and maintain the bus voltages of nearby loads. Compared with many existing distributed reactive power regulation strategies, the proposed method can estimate and maintain the bus voltage of unmeasurable load buses within the limitations. Furthermore, this method releases the hardly achieved bus voltage angle requirement, which makes it practical for real-world.


Ⅰ. INTRODUCTION
Distributed renewable generators (DGs) are deployed worldwide in distribution networks to produce clean, inexpensive electrical power [1]. By using inverters, DGs can provide various ancillary services, such as harmonic compensation [2], voltage support [3], and reactive power (VAR) regulation [4]. In this study, we focus on reactive power regulation.
The main objective of reactive power regulation is to reduce the active power loss in the distribution network and to keep the node voltages within the security limits. Traditionally, these objectives are accomplished to some extent by local reactive regulation [5][6][7], but it has been suggested that these strategies may not be able to guarantee the desired regulation due to the lack of communication [8]. A centralized reactive power regulation strategy can overcome this barrier by solving a centralized optimal reactive power flow (OPF) problem with a central coordinator that receives all the required measurements of the grid [9][10][11]. Nevertheless, an online centralized reactive power regulation strategy requires detailed power flow along with feeders, which is difficult for the distribution system to achieve in real time. For offline strategies, it is also a great challenge to predict the power curve of stochastic DG outputs. These difficulties make centralized optimization impractical in the real world.
In terms of both effectiveness and practicability, the regulation strategies of distributed OPF methods are used to achieve optimal regulation with limited measurement and communication requirements.
The distributed OPF methods decompose the systematic optimization target into decoupled subproblems that can be solved by agents based on partial information. To achieve the decomposition, the OPF problem should be formulated for semidefinite programming (SDP) [12][13][14] or second-order cone programming (SOCP) [15], [16] so that the sparsity of the coefficient matrix can be utilized. A well-designed distributed OPF method exhibits performance close to that of the centralized OPF method with much fewer data requirements and is therefore preferred in practice.
The distributed OPF was achieved via linear approximation in many earlier works. In this case, the detailed power system model and power flow predictions are no longer needed, but measurements at all the buses in the distribution network are still required to maintain bus voltages [17], [18]. The latest reports further relax the above requirements, in which only DGs require phasor measurements at their own connection node and can communicate with adjacent ones; by approximating voltage phasors as a linear function of the injected reactive power, the function between reactive load demands and the total system losses can be formulated as a linear function of the imaginary parts of the DG node voltages, which can be calculated at each node [19][20][21][22][23]. However, synchronous phasor measurement devices in distribution systems are still costly, which prevents the above method from practical applications. Furthermore, these methods can maintain only the voltage of DG buses that can be measured.
To cope with the above problem, a distributed optimal reactive power regulation method that does not require the bus voltage phasor is proposed in this study, where DGs are able to perceive their own voltage magnitude and P/Q flows in the connected branches, communicate with nearby DG buses, and adjust the reactive power injections into the grid, while their common objective is to minimize the total power losses of the system with maintained bus voltages. The proposed method can work under practical conditions, and only RTU measurements such as branch flow and bus voltage magnitude are obtained. The load bus voltage, which cannot be obtained by DG buses, can also be maintained by the estimation method proposed in this paper.
The remainder of the paper is organized as follows. Section Ⅱ presents a generalized LinDistFlow model to formulate the distribution system reactive power regulation problem by considering the load voltage limitation and power loss. Section Ⅲ proposes a distributed method to solve the reactive power regulation problem based on branch flows. Section Ⅴ tests the performance of the strategy with a numerical simulation under different circumstances. Finally, Section Ⅵ presents the conclusions of this work and future research directions.

A. Reactive Power Regulation Model for Distribution Networks
To avoid an electromagnetic loop, distribution networks are radial in operation with one single system bus (SB) connected to the power transmission system. Hence, the topology of such networks can be described by a tree under graph theory, where each branch designates an electric line, while each node represents a bus bar.
Consider a radial distribution network that consists of = + + 1 buses corresponding to a tree with + + 1 nodes, in which the root node indexed by 0 denotes the SB, nodes indexed by = { 1 , 2 , … , } denote buses connected to DGs (DG buses), and nodes indexed by = { 1 , 2 , … , } denote buses connected to only loads (load buses). All the nodes are interconnected by − 1 branches, where the branch indexed by denotes the distribution line connecting node and its parent node. Then, the reactive power regulation problem considering the bus voltage limitations and reactive power constraints of DGs can be defined by the optimization problem in (1):

min
(1) subject to , < < , , ∀ ∈ < < , ∀ ∈ ∪ where denotes the total active power loss of the system, the decision variable ∈ denotes the reactive power generation of the DG connected to bus ∈ , and denotes the voltage magnitude of bus . and can be written as a quadratic function and a linear function of , respectively. Under the approximation model in the following parts, the problem can be solved much more easily.

B. LinDistFlow Approximation of the Reactive Power Regulation Model
The LinDistFlow model [24], which is the linearization of the DistFlow model, is a widely used approximation model of distribution systems with a mainline and no laterals, as shown in Fig. 1 [25]. Since the branch flow is much larger than the power loss, it can be approximated as the total power injections from the downstream nodes of the branch. Since the nodal voltages are close to the voltage of the system bus when the system is operating in a steady state, the voltage drop of a branch can be approximated as a linear function of the branch flow. Under such approximation, the LinDistFlow model can be formulated as (2): where the impedance of branch is defined by = + i . −1 / −1 denotes the branch flows of branch − 1, while / denotes the power injection on bus . denotes the voltage magnitude. Specifically, SB is indexed by 0, whose voltage is assumed to be 0 = 1.00 p.u.
In this study, we extend this model to radial systems with laterals, the topology of which is similar to the example shown in Fig. 2. We denote the path set PATH , = { , }, where is composed of all the branches (edges) in the path between buses < , >, and is composed of all the buses (vertexes) in the path set. For example, buses 7, 2 and 3 and branches 7 and 3 constitute PATH 7,3 in Fig. 2.
The matrix with a size of ( − 1) × ( − 1) is introduced to describe the topology of the system. The element at row column of is defined as: In this definition, , = 1 means node is a downstream node of branch , and branch is a part of the path from system bus to node . For example, only 6,6 equals 1 in the 6 th row of since only node 6 is the downstream node of branch 6; only 1,6 and 6,6 in the 6 th column equal 1 since only branch 1 and 6 are in the path from system bus to node 6. Thus, can describe the sensitivity of bus injections to branch flows, while can describe the sensitivity of branch voltage drops to the bus voltage: / denotes the active/reactive power vectors of all the branches, while / denotes the active/reactive power vectors of all the buses except for SB, denotes the voltage drop vector of all the branches, denotes the voltage magnitude vectors of all the buses except for SB, and is an all-one vector. Since the voltage drop of a branch is approximated as a linear function of the branch flow, the LinDistFlow model for radial systems with laterals can be given as: Using the approximation given in [22] and [23] that all the branches in the system are assumed to have the same impedance angle , the model can be further simplified: According to 3 and 4, the element at row column of can be formulated by: where 0 ∩ 0 is the common part of 0 and 0 . For example, 8,10 can be calculated by the total impedance magnitude of branches 1 and 2.
To separate the components related to in and , can be blocked as = [ T ], where is an × matrix, of which element is defined by = , , is an × matrix, of which element is defined by = , , and is an × matrix whose elements are defined as = , . Thus, the model can be reformulated as: In (8), is the contribution of the active power injection of the system to the active power loss, while is the contribution of the reactive power injection of the system to the active power loss. Since and T are irrelevant to the decision vector , the optimization problem (1) can be simplified to the following quadratic form: where , , , , , , , , . , and , are the vectors of all the limitations.
As a quadratic programming problem, Formulation (9) can be solved with a distributed dual ascent algorithm, as shown in Section III.

A. Centralized Dual Decomposition Solution
The Lagrangian of 9 can be formed as: ] T is the vector of the Lagrangian multipliers of the constraints. Thus, the partial derivatives of 10 can be calculated: Then, the optimization model can be solved with a dual ascent algorithm with two iterative steps: 1) Update the Lagrangian multipliers with dual gradient ascent. Here, can be calculated by (11), [·] + denotes a projection operator on the positive orthant, and denotes a suitable positive constant.
2) Minimize the Lagrangian to update . Since (10) has a quadratic form, the minimum value can be reached when ∂ ∂ = 0, and the optimal solution can be calculated: where ̃ is the contribution of ( + 1) to ( + 1). Alternately executing the given steps can drive the system toward the optimal configuration, but the calculation depends on the power injection and bus voltage measurements of all the buses in the system. To overcome this limitation, the sparsity of matrices −1 and −1 will be utilized, as discussed in the following subsections. The above sparsity allow us to define the neighborhood relationship between buses:

B. Distributed Solution Exploiting Matrix Sparsity
· Two buses are neighbors if and only if there is no other DG bus in the path between them.
· For each DG bus , define its neighborhood ℕ( ) as the collection of all its neighbors. In particular, itself is seen as a member of ℕ( ). Fig. 3 is an example of a neighborhood where the buses 1 , 2 , 3 , 5 , 6 , and 1 , 2 , 3 are neighbors of 2 , while 4 and 4 are not. Based on these definitions, we can transform (13) into a distributed form: where ̃ is the contribution of ( + 1) to ( + 1), and is the weighted sum of the load reactive injections in ℕ( ).
The equations in (14) indicate that the optimal reactive power at a single DG bus can be calculated from measurements of its neighborhood buses so that the communication cost can be effectively reduced. However, the power consumption and voltage of all the buses in the neighborhood are required, which could be relieved by the estimation methods in the following part.

C. Further Relief of the Data Requirements
To update the Lagrangian multipliers , and , in 11, the measurement of bus voltages is required. To calculate , ( ) in 14, the measurement of power injections of load buses is required. Since these data are hard to obtain for privacy and economic reasons, we will estimate them by the measurements in DG buses in the following part of this section.

1) Estimation of the weighted sum of the load reactive injections
, ( ) can be estimated by branch flows and the voltage drop between DGs in a neighborhood.
The neighbors of each DG bus and the branches interconnecting them constitute a subtree of the whole system, and the root node of the subtree is denoted as . For example, 3 = 2 is the root node of ℕ( 3 ) in Fig. 3. Thus, we can divide , in 14 into two parts: Since is the root node of ℕ( ), e 0, will be the subset of e 0, for every ∈ ℕ( ). Then, according to 7, we have: where is the branch flow from to its parent node, which can be seen as the reactive power injection of the equivalent node of and all its descendant nodes. According to 17 and 18, , ( ) can be accurately estimated.
According to 4 and 6, , ( ) can be estimated by the voltage drop between and , which can be calculated by: In most countries, retail customers may be required to correct their power factors to a standard (e.g., 0.95) or pay an additional electricity tariff for a low power factor. In that case, the power factors of load buses in a neighborhood will have similar power factors. Then, res( , ( )) can be ignored so that , ( ) can be approximated by ̂, ( ) , which can be calculated without the requirement of the measurement of load buses. Then, , ( ) can be estimated as , ( ) +̂, ( ).

2) Estimation of the load voltage
Since there is no DG connected to load buses, load buses will not undergo overvoltage. Thus, , can be fixed to all 0 vectors so that only the risk of low voltage needs to be considered for the load buses.
The load voltage will be estimated with a subneighborhood as a unit. If the bus with the lowest voltage can meet the voltage constraints, then the buses in the subneighborhood will not offend the voltage constraints; thus, we do not need to estimate the load voltage in a subneighborhood.
According to the network topology, subneighborhoods can be classified into three types, as shown in Fig. 4. Type I: As shown in Fig. 4-A, there is only one DG bus that injects power flow into a type I subneighborhood where all the buses at the edge of the subneighborhood are equipped with DG. In that case, there will be no reverse flows in the subneighborhood so that load buses will not be under voltage if all the DG buses are in the same subneighborhood. Hence, we do not need to estimate the load voltages in a type I subneighborhood so that the Lagrange multipliers corresponding to the lower bounds of the load buses in a type I subneighborhood can be set to 0.
Type II: As shown in Fig. 4-B, branch flows are not in the same direction in a type II subneighborhood, but every load bus is directly connected with at least one DG bus with a branch whose power flow can be measured directly by the DG bus. In that case, the load voltage in the subneighborhood can be calculated by the LinDistFlow model.
Type III: As shown in Fig. 4-C, the branch flows are not in the same direction, and not all load buses are connected with a DG bus directly in a type III subneighborhood. In this case, the load voltages of some buses are determined by the load distribution inside the subneighborhood.
To describe the load distribution, define d ℕ (  Then, d̂ℕ ( ) can be constrained by: In addition, all the load buses absorb power from the system, and the total demand in the subneighborhood can be determined by the branch flows: According to the LinDistFlow model, for all d̂ℕ ( ) that satisfy Eqs. (24)- (25), it is easy to prove that ̂= if ℕ( ) is a type II subneighborhood, and ̂ will higher than the lowest ℎ if ℕ( ) is a type I subneighborhood. That is, using (23) to estimate the lowest voltage in a type I or II subneighborhood will not introduce any error, regardless of the error of d̂ℕ ( ) .
For a type III subneighborhood, the error of d̂ℕ ( ) may introduce error to the estimation of ̂ if there is no branch directly connecting and a DG bus and the equations in (24) are not enough to determine d̂ℕ ( ) . The estimation of the lowest voltage in a type III neighborhood will be more conservative if the estimation of load distribution is more clustered at a location close to DG buses, and it will be more aggressive if the estimation of load distribution is more clustered at the location far from DG buses. As a compromise, we assumed that the power demand is distributed as evenly as possible in the unmeasurable loads under the constraints in Eqs. (24)- (25) min d̂ℕ ( ) 2 (d̂ℕ ( ) ) (26) subject to Eqs. (24)- (25) where 2 (d̂ℕ ( ) ) denotes the variance in all the elements in d̂ℕ ( ) .
Substituting the estimates of , ( ) and into 12 and 14, the reactive power regulation problem can be solved only with the measurement of DG buses.

IV. REAL-TIME CONTROL SYSTEM
According to the distributed solution method given in Section III, a real-time distributed reactive power regulation system can be designed for a distribution system.
In this system, all the DG buses are assumed to be equipped with reactive power regulators that can gather the local bus voltage, bus power injection, and branch flows measured by the CT on local buses. Regulators in a neighborhood can communicate with each other and regulate the reactive power injection of the local DG.
The control law for the regulator in each DG bus is given below:

Executed at time step 0：
Initialize the local Lagrangian multiplier vector (0) as an all-zero vector. Executed at each step t： 1. Gather the measurement of local power injection, local voltage and branch flows.
2. Gather the measurements of DG voltages, DG power injections and branch flows measured on DG buses as well as the Lagrangian multipliers of DG voltages and DG reactive power injections from the neighborhood.
3. Estimate the load distribution for all subneighborhoods in the neighborhood by 26

V. CASE STUDY
The proposed algorithm has been simulated on a testbed inspired by the IEEE 33-bus radial distribution system [26]. A schematic diagram of the test system is shown in Fig. 5. The default load injections of the system, as shown in Fig.  6, are not distributed uniformly in each neighborhood, and the power factors of the loads in each neighborhood are not the same. These characteristics reflect the actual situation of the real-world distribution system and will produce errors in the estimations given by 22 and 26. In addition to the proposed method, a centralized optimizer is used to give the optimum output of the case. The DORPF method given in [23], which requires the PMU data of all the DG buses without considering the load voltage constraints, is used as a positive control group, while a negative control group in which all the DG buses are used as pure active sources is also given. The simulation lasts 100 steps, while the load demands and the DG capacity are increased by 1.2 times the default value and considered static. The curves of the lowest bus voltage in the system and the active power losses are shown in Fig. 7 and 8, respectively. The simulation result shows that both the given method and the DORPF can maintain the bus voltages and reduce the power losses; the estimate used in part C of Section III does not introduce error to the regulation even though the power factor and load distribution in each neighborhood are not ideal.  By increasing the load demand in the system and the capacity of DG buses to 1.4 times the default value, the curves of the lowest bus voltage in the system and the active power losses are shown in Fig. 9 and 10, respectively. The simulation result shows that the DORPF method cannot maintain the under voltage of load buses, while the given method can maintain the load voltage within 80 steps. The active power loss after the voltage is maintained by the proposed method is slightly larger than the optimal solution given by the centralized optimizer, but the difference is acceptable relative to the total active loss. To determine whether the method works in more general scenarios with different load power factors and load distributions, we generate 1000 scenarios in which the active and reactive power injections of buses follow the normal distribution, where the mean value is the default value of a 33-bus system and the standard deviation is 4%, resulting in different power injections and power factors for each bus in each scenario. For each scenario, the DG reactive power will be regulated by 100 steps of the proposed method and the ORPF method, respectively, and we can validate the performance of the proposed method through comparison. When we set the load level of the system as 1.2x of the default value, both the given method and the DORPF method can maintain the load voltage and reduce the power losses. The lowest bus voltage and the active power losses for each sample are shown in Fig. 11 and 12, respectively. When we increase the load level to 1.4x, the proposed method can still maintain the load voltage, whereas DORPF cannot. Although the proposed method undergoes a higher power loss than DORPF to maintain the load voltage, the loss is still significantly reduced compared with that of the negative control group.
These results show that the proposed method performs similarly to the existing distributed method, which uses PMU under light-load conditions, in terms of its ability to reduce the active loss while maintaining the load voltage, and it can better maintain the load voltage than the existing method under heavy load conditions.

VI. CONCLUSIONS
In this paper, we propose a distributed reactive power regulation strategy that works in an environment in which no PMU data can be obtained. In the proposed strategy, distributed generators can perceive their own voltage magnitude and the P/Q flows in the connected branches, communicate with nearby DG buses, and then adjust the reactive power injections into the grid to minimize power losses and maintain bus voltages. Compared with the existing methods, the proposed method can be implemented more easily due to its low measurement requirements. Furthermore, it can keep not only the DG buses but also the unmeasurable load buses within the voltage constraints.