Research on Three-phase Optimal Power Flow for Distribution Networks Based on Constant Hessian Matrix

: The optimal power flow (OPF) problem for active distribution networks with distributed generation (DG) and a variety of discretely adjustable devices (e.g., on-load tap-changers, OLTCs) is essentially a non-convex, nonlinear, mixed-integer optimization problem. In this paper, the quadratic model of three-phase OLTCs is proposed by adding branch currents as unknown variables, which guarantee a constant Hessian matrix throughout iterations. This paper proposes a three-phase OPF model for active distribution networks, considering a three-phase DG model. The OPF model is solved by an interior point method incorporating a quadratic penalty function as opposed to a Gaussian penalty function. Furthermore, a voltage regulator is also incorporated into the OPF model to form an integrated regulation strategy. The methodology is tested and validated on the IEEE 13-bus three-phase unbalanced test system.


Introduction
The optimal power flow (OPF) for active distribution networks with distributed generation (DG) and a variety of discretely adjustable devices is a non-convex, non-linear, mixed-integer optimization problem.It involves both discrete and continuous variables.By optimizing the operation strategy of the adjustable devices, including shunt capacitors, adjustable distributed energy resources, on-load tap-changers (OLTCs), etc., the OPF aims to minimize the operation cost which takes into account network losses.DG and energy storage units bring new challenges to traditional OPF problems [1][2][3][4][5][6][7].Reference [6] studied how the uncertainty in wind outputs and the correlation among multiple wind farms would affect the OPF.Reference [6] did not consider discrete adjustable devices, which are prevalent in active distribution networks.Reference [7] proposed an OPF model for active distribution networks with OLTCs, where a piecewise linear model is adopted to model OLTCs.The OPF is solved by a second-order cone relaxation method [7].Reference [8] optimizes the operation of an active distribution network where a non-coupled model is adopted for DG, yielding approximate results.Reference [9] proposed a three-phase steady-state model for DG -the model can be applied to OPF for active distribution networks.
The primal-dual interior point method has been widely used to solve the OPF problem for traditional power systems because of its advantages such as efficiency and fast convergence [10][11][12][13][14][15][16][17][18][19][20].Reference [13] proposed a discrete variable processing method based on the Gaussian penalty function, but it does not compare the computation efficiency of the function with that of other penalty functions, e.g., the quadratic penalty function as introduced in [20].By introducing a virtual node into the OLTC model, reference [14] transformed the OPF model into a quadratic optimization in the Cartesian coordinate system, thus improving computation efficiency.This paper proposes a new asymmetrical injection model for DG with three phase coupling and a new quadratic model for three-phase OLTCs.The three-phase OPF for active distribution networks is solved by a predictor-corrector primal-dual interior point method (PCPDIPM) incorporating a penalty function.The network voltages are fine-tuned on each phase by a voltage regulator, and an integrated regulation strategy is proposed and incorporated into the OPF problem so that the voltages are fine-tuned on each phase to further optimize the network loss.Case studies demonstrate that the quadratic models for three-phase OLTCs and DG improve both the computation efficiency and the accuracy of the OPF.Case studies also compared the advantages and disadvantages of 1) the Gaussian penalty function; and 2) the quadratic penalty function with continuous variable discretization process.The impacts of the voltage regulator on nodal voltages and network losses are also demonstrated.

Quadratic model of three-phase OLTCs
In the Cartesian coordinate system, the OPF model is a nonlinear optimization problem with an order higher than quadratic if considering the turns ratio of the OLTC as a control variable.When solving the OPF problem by the interior point method, the Hessian matrix is updated in each iteration, bringing about significant computation burden [14] [21].In order to solve this problem, a new three-phase OLTC model is proposed in this paper.By adding the branch currents as state variables, the Hessian matrix becomes constant throughout the iterations of the OPF.
Fig. 1 shows the OLTC in a Wye-delta (Yd) configuration, a virtual node m is added to branch ij.This transforms the OLTC into an ideal transformer im (with an adjustable turns ratio k) connected in series with an equivalent impedance mj (with impedance R + jX).The nodal voltages and branch currents are also shown in Fig. 1.

Fig.1. Three-Phase Model of OLTCs
With the ideal transformer, there is Based on the conservation of energy, there is Substituting (1) in (2) yields, In the three-phase system, the line voltage is set as the base voltage.So, for the Yd transformer, the standard transformation ratio is 1: 3 under the per-unit system.
Considering the equivalent impedance branch mj, there is

Quadratic model of three-phase DGs
The model in Fig. 2 represents different types of DG units, including those with the interface of voltage-source converters, wind generation based on doubly fed induction generators, synchronous generators, etc.

Fig.2. A generalized sequence component model of DG
Where I 0 d +jI 0 q , I + d +jI + q , and I - d +jI - q represent zero sequence current, positive sequence current, and negative sequence current from DG to the network, respectively.U 0 d +jU 0 q , U + d +jU + q , and U - d +jU - q represent zero sequence voltage, positive sequence voltage, and negative sequence voltage at the coupling point, respectively.R 0 +jX 0 and R - +jX -represent zero sequence impedance and negative sequence impedance, respectively.U 0 rd +jU 0 rq and U - rd +jU - rq represent the zero sequence and negative sequence of the excitation voltage, respectively.The DG control system consists of three parts: active power control, reactive or voltage control, and the control of the unbalanced components.The active power control keeps the positive sequence active power unchanged.The reactive control maintains a constant reactive power injection to the grid.DG excitation voltage contains only positive sequence component, with negative and zero sequence components being zero.Consequently, the equality constraints are established, as shown in equations (5, 6, 7).

Re (
)( ) Im ( )( ) Where, Re[] and Im[] correspond to the real part and the imaginary part of the expression, respectively.P sp and Q sp correspond to the positive sequence values of the target active power and reactive power, respectively.
The grid interface of a three-phase coupled DG is modelled as both a voltage controlled voltage source and a current controlled current source, as illustrated in Fig. 3.The sequence values of the DG injection currents and the terminal voltages are unknown variables in the quadratic model of DG.The variables include 0

Fig.3. Sequence-phase coupled interface
The transformations from phase values to sequence values are given by 3 0 0 Where, 0 abc T → + − is phase to sequence transformation matrix.'•' corresponds to complex quantities.
The constraints (5) -( 9) consider the coupling characteristics of three-phase powers.They are different from the steady-state model in [22] which does not consider sequence control under unbalanced condition.The equations involving DG are either linear or quadratic -this guarantees a constant Hessian matrix in the OPF.
A popular model of DG is the three-phase noncoupling model [8], described in (10).

(
)( ) Where, P P and Q P correspond to the active and reactive power outputs from DG on phase P, respectively.U P re and U P im correspond to the real part and the imaginary part of the nodal voltage of DG on phase P, respectively.I P re and I P im correspond to the real part and the imaginary part of the injection current from DG on phase P, respectively.
The three-phase coupled model of DG is accurate model, and the three-phase non-coupled model is an approximate model.So, the optimal power flow results of the three-phase coupled model are more accurate than those of the approximate three-phase non-coupled model.
The above formulas constitute the DG constraints in the three-phase OPF model.Since the highest order of these formulas are quadratic, the Hessian matrix is constant throughout iterations.

Three-phase model of voltage regulator
A voltage regulator is connected to the grid for the fine-tuning of nodal voltages.A three-phase Y-connected voltage regulator consists of three single-phase voltage regulators, each with a tap changer that changes the tap position.The configuration of the three-phase voltage regulator is illustrated in Fig. 4.

Fig.4. Three-phase star-shaped voltage regulator
The following constraints apply: In formula (13), a plus ('+') corresponds to 'step up' and a minus (' -') corresponds to 'step down'.'Tap' denotes the position of voltage regulator tap (e.g., 32-level adjustments within a range of ±10%).The tap is controlled to regulate the voltage of load center with a line voltage drop compensator.The line voltage drop compensator is coupled with the distribution line through a voltage transformer (turns ratio N PT : 1) and a current transformer (turns ratio C TP : C TS ), as shown in Fig. 5.The impedance of the compensator represents the equivalent impedance from the regulator to the load center: The voltage at the load center is If the voltage level of the load center is 120V and the bandwidth is 2V, the voltage will change by 0. The objective function of the OPF model is the minimum loss of the network, as shown in equation (18).
Where, bn is the number of branches in the network.U br,P and I br,P correspond to the branch voltage and branch current on phase P of the br th branch, respectively.And the branch voltage is the difference between the nodal voltages at both ends of the branch.
The equations include the KCL and KVL constraints (as shown in equation ( 19)), transformer branch constraints (as shown in equation (1) (3) (4)), DG branch constraints (such as (5) -( 9)), the regulator branch constraints (as shown in equation (11) Where, U br,P,re and U br,P,im correspond to the real part and the imaginary part of the branch voltage on phase P of the br th branch, respectively.U n,P,re and U n,P,im correspond to the real part and the imaginary part of the nodal voltage on phase P of the n th node, respectively.I br,P,re and I br,P,im correspond to the real part and the imaginary part of the branch current on phase P of the br th branch, respectively.A is the node-branch incidence matrix The inequality constraints include state variable constraints and control variable constraints.The state variable constraints include the generator active and reactive power constraints, the node voltage amplitude constraint, and the line transmission power P ij constraint (as shown in equation ( 20)).The control variable constraints including the constraints of the OLTC turns ratio K T , the voltage regulator tap position K V , reactive power capacitor compensation capacity Q C , as shown in equation (21) Thus, the equations in the OPF model are either linear or quadratic which guarantees a constant Hessian matrix when the interior point method is used to solve the model.In this paper, we use the predictor-corrector primaldual interior point method [23] to solve the optimization model with the advantages of efficiency and fast convergence.

Continuous variable discretization process
The OPF for active distribution networks involves discrete control variables (such as OLTCs turns ratio, capacitor compensation capacity, etc.).The objective function incorporates a penalty function to cope with the discrete variables.This introduces a virtual loss to the objective function which will reduce the error brought by rounding-off.At present, the quadratic penalty function [20] and the Gaussian penalty function [13] are two widely used penalty functions, as given by ( 22) and ( 23), respectively.
Where v is penalty factor.The value of the function drops to zero when the distance from the vector When the penalty factor is large enough, the discrete control variable will be optimized to the corresponding discrete value, so that the penalty function value becomes zero and the objective function to reach the minimum value.
Compared with the quadratic penalty function, the Gaussian penalty function decreases faster than the quadratic penalty function with the decrease of the distance between the vectors x and b; the former is more sensitive to the variation of discrete variables.However, the Gaussian penalty function is a non-quadratic function, which makes the Hessian matrix change throughout iterations -this compromises computation efficiency.

Integrated regulatory strategy
In order to minimize network losses, this paper proposes an OPF integrated regulatory strategy with voltage regulators.The integrated regulatory strategy is defined as an OPF strategy for a variety of adjustable devices to coordinate and optimize in this paper -the voltage regulator fine-tunes nodal voltages on each phase, and constitutes integrated regulatory strategy with the generator, OLTCs, shunt capacitors and DGs.
The voltage regulator is a special transformer with a zero equivalent impedance.The voltage regulator is modelled by a quadratic model where the tap position is chosen as a control variable.The Hessian matrix remains constant under the integrated regulatory strategy.This is the reason why the integrated regulatory strategy is suitable for the OPF model proposed in this paper.

Example
According to the three-phase OPF model for distribution networks based on constant Hessian matrix, this paper uses Symbolic Math Toolbox to verify the model on MATLAB.And the OPF analysis is carried out on the modified IEEE13 test system as shown in Fig. 6.The mutual-impedance and the mutual-admittance are omitted in the figure.The control variables in the example are generator reactive power, DG output, the turns ratio of the OLTC, voltage regulator tap position and capacitors adjustment, in which the turns ratio of the OLTC and the capacitors are discrete control variables.The range of the turns ratio is set at 0.90~1.10,divided into 8 taps, its step Page 4 of 7

IET Generation, Transmission & Distribution
This article has been accepted for publication in a future issue of this journal, but has not been fully edited.Content may change prior to final publication in an issue of the journal.To cite the paper please use the doi provided on the Digital Library page.size is 0.025.The maximum reactive power of the capacitors is 0.02 p.u. and the step size is 0.01 p.u..The range of all nodal voltages is 0.9 ~1.1 p.u..The voltage regulator tap has a 32-stage adjustment with ± 10% regulation range, the voltage level of load center is 120V, and the bandwidth is 2V.

Analyze the quadratic model of three-phase
OLTCs In order to verify the rationality of the proposed three-phase OLTCs model in the OPF problem for active distribution networks, the OPF calculation based on the predictor-corrector primal-dual interior point method embedding quadratic penalty function (PCPDIPM-QPF) is used to compare the traditional non-quadratic model for three-phase OLTCs [21] and the quadratic model for threephase OLTCs proposed in this paper.The results are shown in Table 1.  1 shows that the results (network loss) of the two methods using different models are the same, but the quadratic model has a faster calculation speed.The reason is that, the Hessian matrix of the non-quadratic model of OLTCs is not constant in the OPF calculation process, and is updated in each iteration, which makes the calculation speed very slow.The quadratic model of OLTCs guarantee a constant Hessian matrix throughout iterations, thus reducing the calculation time greatly.
In addition, the computation time of the two models is relatively long.This is because the Hessian matrix is generated using Matlab's automatic differentiation function, which is convenient but not efficient.This, however, does not affect the conclusion of this paper.The case study involves a total of 17 iterations.The one-time calculation of the Hessian matrix takes 6.67s on the computer with Intel core i3-3240M CPU 3.40GHz, 4GB memory.Therefore, a conservative estimate of the computation time of the nonquadratic model is 17 times of the time for calculating the Hessian matrix, i.e., 113.84s.Therefore, the computation time in Table 1 is justified.

Analyze the continuous variable discretization process
In order to verify the convergence properties of the quadratic penalty function and the Gaussian penalty function in the OPF calculation, the PCPDIPM-QPF and predictorcorrector primal-dual interior point method embedding Gaussian penalty function (PCPDIPM-GPF) are designed.The OPF for the modified IEEE13 three-phase system is calculated by the two algorithms.This paper proves the advantages and disadvantages of each method by comparing the optimization results (network loss), iteration number and calculating time (the accuracy is e -10 ).
Table 2 contains the results of the OPF calculation.It shows that, the PCPDIPM-GPF is better than the PCPDIPM-QPF in terms of the optimization results (the network loss is reduced by 2.4%).But the iteration times of the PCPDIPM-GPF is slightly higher than that of the PCPDIPM-QPF, and the PCPDIPM-GPF takes longer time.The reason is that the Gaussian penalty function is a higher order function and the Hessian matrix is a non-constant matrix throughout iterations.So the Hessian matrix is updated as the iterations increasing, which increases the computation time.Fig. 7 shows the comparison of the PCPDIPM-QPF and PCPDIPM-GPF.It can be seen that the dual gap decreases as the iterations increasing, and the time of the PCPDIPM-QPF is about 1/5 of the PCPDIPM-GPF when the convergence effect of the two penalty functions reaches the desired value.

Analyze the integrated regulatory strategy
The integrated regulatory strategy of OPF is implemented by the modified IEEE13 three-phase system with the voltage regulator.It verifies the effect of the voltage regulator on each phase voltages fine-tuning in active distribution networks, and further reduce the network loss in the OPF calculation.The results are shown in Table 3.As shown in Table III, the network loss of the network with voltage regulator is reduced by 0.32% compared with that without the voltage regulator.The threephase voltage profile of all nodes, before and after integrated regulatory, is shown in Fig8, Fig9 and Fig10.Thus, the voltage regulator has an effect on fine-tuning the nodal voltages.It can be seen that, the voltage regulator can regulate the nodal voltages at each phase and further reduce the network loss in the OPF calculation.

Conclusion
1) This paper proposes a quadratic model for threephase OLTCs.The model is incorporated into the OPF model, where the Hessian matrix becomes constant throughout iterations.This results in a considerable improvement in computation efficiency.
2) This paper solves the three-phase OPF problem for active distribution networks by using the interior point method.The OPF model incorporates both the quadratic model for the three-phase OLTC (as explained above) and the asymmetric model for the three-phase DG.The methodology is validated by case studies.
3) The case studies demonstrate the advantages and disadvantages of the interior point method incorporating the quadratic penalty function and the Gaussian penalty function.
4) An integrated regulation strategy supported by a voltage regulator is incorporated into the OPF.Case studies demonstrate that the voltage regulator fine-tuned the network voltages, thus further reducing network losses.

Acknowledgments
This work was supported by National Key Research and Development Program (2016YFB0900101) and Chinese Scientific Fund (2017QC145).

Fig. 5 .
Fig.5.Line voltage drop compensator3.Optimal Power Flow Model and Algorithm3.1.OPF modelThis paper establishes three-phase OPF model for distribution networks by taking the nodal voltages and branch currents as the state variables.The objective function of the OPF model is the minimum loss of the network, as shown in equation(18).
This article has been accepted for publication in a future issue of this journal, but has not been fully edited.Content may change prior to final publication in an issue of the journal.To cite the paper please use the doi provided on the Digital Library page.

Table 2
The optimization results