Cheap Conic OPF Models for Low-Voltage Active Distribution Networks

The increasing focus on the active participation of low-voltage (LV) active distribution networks (DNs) in electricity markets requires the real-time optimal control of these DNs. To achieve this goal, a cheap semi-definite programming (SDP)-based optimal power flow (OPF) model for active neutral-equipped DNs, hosting both wye- and delta-connected loads, is proposed in this paper, aiming at overcoming the high computational requirement of the primal SDP-based OPF model. The coupled power injections between conductors are explicitly represented for each conductor by utilizing the network admittance matrix-based approach. Furthermore, three novel propositions (P1, P2 and P3) are proposed for the modelling of the constant current component of ZIP end-users in the context of the proposed OPF model. Moreover, the impact of the voltage-angle deviation on the exactness of the P1- and P2-based models is discussed. Simulations are carried out on several LV active DNs for various parameters of ZIP end-users, and the quality of the proposed OPF model is verified through the %optimality gap, power mismatch, voltage violation and root-mean-square error criteria. It is successfully shown that the proposed OPF model provides an optimal and feasible solution for all load types (wye, delta, mixed wye-delta) under a large range of ZIP load parameters. Furthermore, among the three propositions, the P3-based OPF model appears to be the most accurate in terms of determining an optimal and feasible solution. Finally, the reduced computational time of the cheap conic model allows its real-time implementation for medium- and large-sized DNs for which the primal multi-phase SDP-based model is practically difficulty to realize.

SOCP Second-Order-Cone-Programming. ZIP Constant Impedance, constant Current and constant Power.

I. INTRODUCTION
The optimal power flow (OPF) problem formulation based on the convex relaxation and approximation techniques has gained significant attention in the past few years due to the fact that these techniques ensure global optimality and problem infeasibility through several sufficient conditions and criteria. In this context, various relaxations [1]- [8] based on the concepts of lift-and-project, convex envelopes and hierarchies of moment have been proposed mainly for the single-phase equivalent representation of distribution networks (DNs). The global optimality of semi-definite programming (SDP)-based relaxation for radial and mesh DNs is proved in [9]- [11], whereas the exactness of SDP-based and second-order-cone-programming (SOCP)-based OPF models are discussed in [12], [13]. A detailed survey on these relaxations can be found in [14]. Due to the large computational requirements of SDP-based relaxation, the sparsity of power networks is exploited in [15]- [18], which significantly reduces the computational time (CT) of SDP-based relaxation and, therefore, allows its practical realization for large networks.
To determine the optimal active and reactive power operational point of single-and three-phase distributed energy resources (DERs) connected to medium-voltage (MV) and low-voltage (LV) DNs, it is necessary to solve a multi-phase OPF problem. Since these DNs exhibit significant unbalance loading, solving a single-phase OPF model for them provides an incomplete picture of the state variables of these DNs. However, due to the increased complexity and non-convexity of a multi-phase OPF problem in comparison to a single-phase OPF model, very few researchers have studied it for three-phase networks. Among the relaxation techniques, the single-phase SDP-based relaxation in its primal 1 form is extended in [19] and solved through an alternative direction method of multipliers (ADMM)based distributed algorithm. The single-phase bus injection model (BIM)-based and branch flow model (BFM)-based SOCP relaxations are extended in [20], [21] to three-phase unbalanced DNs, which leads to cheap 2 BIM-SDP-based and BFM-SDP-based OPF relaxations. It must be noted that the resultant relaxations involve positive semi-definite (PSD) constraints in their formulation, even for a radial system, due to the existence of mutual coupling among network phases. Unlike a single-phase equivalent DN, for which the SDP and SOCP models have been proved to be equivalent due to the transformation of SDP constraints into SOCP constraints [22], [23], the same has not yet been proved for three-phase or multi-wire models. Consequently, the primal multi-phase SDP model is considered a generalization of cheap multi-phase SDP-based OPF models. However, it must be noted that the PSD constraints are defined over smaller matrices, which involve either node voltages (BIM) or a node voltage along with the branch flow quantities (BFM), in cheap OPF relaxations. Consequently, the computational efficiency of these cheap SDP-based OPF models does not degrade by much relative to the computational requirements of single-phase SOCP-based relaxations. In fact, the CT reduces significantly when the PSD constraints are defined only on the edge minor, three-cycle minor or four-cycle minor of a full matrix variable.
Quite recently, a chordal-relaxation-based SDP model was proposed in [24], which exploits the chordal sparsity of a three-phase SDP-based relaxation to speed up the computational process. Furthermore, in [25], the chordal relaxation is embedded inside a convex iterative algorithm, which, on the one hand, increases the computational efficiency for large DNs and, on the other hand, guarantees the global optimality of the recovered solution when the trace of the regularization terms becomes zero.
In all the above-mentioned proposed single-and three-phase OPF relaxations, except [21], phase-ground constant power injections are taken into consideration by assuming only the wye-connected end-users. In [21], a BFM-SDP-based relaxation is solved for the mixed wyeand delta-connected loads. However, the optimization variable (OV) related to the latter load type has a large eigen-value ratio (EVR), which indicates that the proposed approach fails to provide a feasible solution for delta end-users. In [26], an attempt is made to incorporate a constant impedance, constant current and constant power (ZIP) load into the SDP-based OPF model by introducing additional voltagemagnitude-based OVs at load buses. However, the proposed scheme is developed only for the single-phase equivalent representation of a DN, and it has been shown in our latest work [27] that the developed approach cannot be directly generalized for multi-phase DNs.
In the case of distribution networks equipped with neutral conductor(s) (DNNs) and hosting shunt elements either between the phase-neutral (wye) or phase-phase (delta) conductors, the assumption of a phase-ground power injection is no longer valid. Consequently, the existing three-phase OPF models cannot be extended to 4 or more wire networks by simply specifying the voltage drop and power balance constraints corresponding to the additional conductors. This restriction is due to the fact that in these networks, power injection is coupled between any two conductors, which is not the case for a single-or three-phase network configuration. Consequently, for neutral-equipped DNs, a detailed primal SDP-based OPF model was recently introduced by the authors in [27], [28]; this model successfully incorporates both wye-and delta-connected loads and distributed generators (DGs) without reducing the network configuration from four conductors to three phases. However, as shown in [27], [28], the primal SDP-based OPF model suffers from a high CT and cannot be realized for medium and large real LV active DNs. Since LV DNs are rapidly transforming from passive to active networks and, therefore, can play a significant role in the near future in the day-ahead and ancillary services markets, it is therefore necessary to optimally control and manage these networks on a real-time basis.
To realize the real-time implementation of an OPF model for LV active DNs, a cheap SDP-based model provides an alternative to the primal SDP-based relaxation because of its low computational requirements. Therefore, the following goals are set and successfully achieved in this work.
A. CONTRIBUTIONS 1) Extension of a centralized three-phase cheap SDPbased OPF model, which utilizes the BFM-based power system representation, for DNNs containing both wye-and delta-connected ZIP loads 2) Proposal of three novel techniques for the approximate modelling of the constant current component of a ZIP end-user to incorporate it into the proposed OPF model 3) Development of tighter limits for the complex phase-ground voltage variables involved in propositions 1 and 2, and proposal of a novel approach for determining the bounds on neutral conductor voltage variables 4) Demonstration of the impact of the voltage-angledeviation bounds on the exactness of the proposed OPF models The rest of the paper is organized as follows. Section II introduces the nomenclature and provides a quick recap of the correction current injection (CCI)-methodology-based [29], [30] power decoupling approach presented in [28]. In section III, a non-linear OPF model, based on the BFM representation, is initially presented, followed by the three novel modelling techniques related to the constant current component. Moreover, a complete multi-wire cheap SDP-based OPF model is presented. In section IV, a quality assessment of the developed OPF model is thoroughly carried out, whereas section V concludes this paper.

II. MODELLING OF MULTI-PHASE RADIAL DISTRIBUTION NETWORK AND DECOUPLED POWER INJECTION REPRESENTATION
A. NOMENCLATURE Consider a radial active DN comprising a set N of n + 1 nodes as N = {0, 1, 2, · · · , n} and a set E of branches as E = {(i, j) ⊂ N × N } having γ conductors between any two nodes. Let (j : i ∼ j) ∈ E denote a line connecting node i to node j. Let 0 represent the substation node, and define N + as N \{0}. Let φ ij ⊂ {a ij , b ij , c ij , n 1 ij , · · · , n w ij } and φ i ⊂ {a i , b i , c i , n 1 i , · · · , n w i } be sets containing the phase and neutral conductor of a line (i, j) ∈ E and the phase and neutral conductor of a node i ∈ N , respectively.
Let ω i ⊂ {n 1 i , · · · , n w i } and ω ij ⊂ {n 1 ij , · · · , n w ij } be sets containing the neutral conductor of a node i ∈ N and line (i ∼ j) ∈ E, respectively, and define η i = φ i \ω i and η ij = φ ij \ω ij as sets containing only phase conductors of a node i ∈ N and a line (i ∼ j) ∈ E, respectively. Consider as sets containing the phase-neutral and phase-phase connection of a node i ∈ N , respectively. Lines are modelled as π-equivalent representations, and complex series impedance of a line (i, j) ∈ E is represented by Z ij ∈ C |φ ij |×|φ ij | , where | X | defines the cardinality of a general set X . Consider V ϕ i as the complex phase-ground voltage of a phase ϕ ∈ φ i of node i, and let 0] T as the slack node voltage vector, defined as a positive sequence voltage triplet for phases {a, b, c} and as 0 for neutrals {n 1 , · · · , n w }. Let I ϕ i be the complex current injected into the phase ϕ ∈ φ i of a node i, and let I ϕ ij be the complex current flowing on the phase ϕ ∈ φ ij of a line (j : i ∼ j). In vector form, let I ij = [I a ij I b ij I c ij I n 1 ij , · · · , I n w ij ] T be the complex current flowing on a line between nodes i and j. Let λ ij = I ij I H ij be a hermitian current matrix for each line (j : i ∼ j). Consider as the complex phase-ground (for both DGs and loads) power injection vector at node i. Let S ij = V i I H ij be a complex power flow matrix for each line (j : i ∼ j). Let x and x denote the lower and upper limits of a scalar/vector variable x, respectively, and finally, let x * denote the conjugate of a complex variable x.

B. DECOUPLED POWER INJECTIONS BASED ON THE CORRECTION CURRENT INJECTION APPROACH
Consider a typical branch of a 4-wire DN incorporating loads, DGs and grounding impedance, as shown in Fig. 1. Based on the CCI approach [29], [30], any shunt element can be modelled by a combination of constant admittance and a suitable current injector, if required. Consequently, the apparent power of a single-phase shunt element, referred to as the r-th FIGURE 1. DN representation with constant wye and delta load admittances and correction current injections. VOLUME 8, 2020 iteration of a load flow, can be written as: where i ∈ N , α ∈ {υ, σ }, υ ∈ ψ i , σ ∈ χ i , I is the CCI term and Y i(0) is the constant admittance, which refers to the rated apparent power and voltage values. For wye-and delta-connected end-users, it is defined as follows: where subscript (0) refers to the rated value. On the other hand, the voltage dependency of each component of a ZIP load can be defined as follows: By substituting (3a)-(3c) in (1), the apparent power related to each component of a ZIP end-user can be expressed as follows: represents the constant admittance term related to each component of a ZIP load at node i. Note that the apparent powers expressed in (4a)-(4c) do not explicitly represent the power injection in each conductor due to the presence of a coupled admittance between the conductors. This indicates that the explicit power injection corresponding to each conductor can be obtained by splitting the admittance between conductors. To achieve this, the following load matrices Y i and i are introduced.
whereĨ Y andĨ are the incidence matrices and are defined as follows: T are 3 × 1 vectors containing constant admittances of wye-and delta-connected shunt elements, respectively, as shown in (2). Through Y i and i , it eventually becomes possible to obtain the explicit power injections as expressed below for a ZIP end-user connected at node i.
where I c , P c are load matrices that are independent of the ZIP parameters, whereas load matrices Z , I and P are formed by multiplying the k Z , k I and k P coefficients by the admittances of the respective load components. Once explicit injections are obtained, they can be represented in terms of the OVs, which are based on the chosen modelling approach. In [27], [28], the primal SDP-based approach is used to develop an OPF model. In this paper, the cheap SDP-based OPF model is adopted for this purpose, which is comprehensively reported in the next section.

III. NON-LINEAR AND CHEAP SDP-BASED OPF MODELS
The BFM-based power system representation was introduced in [31] and later adopted by [32] and [20] in the context of the OPF problem formulation for the single-phase equivalent and three-phase representation of electrical networks, respectively.

M1: Non-Linear BFM-Based OPF Model for DNNs
Var : subject to : In this paper, the non-linear (NL) and cheap SDP-based OPF models are proposed for a 4-wire DN. The resultant models are based on the BFM-based power system representation and are shown in models M1 and M2, respectively. Notably, the non-convexity in M1 appears due to the lower bound in (10d) and the bilinear product in (10f), whereas in M2, the non-convexity is removed by specifying the PSD constraints in (11g)-(11h). Furthermore, note that all terms, except the load injection S i term, involved in the constraints of both models either are based on the variables of these models or come from the network data.

M2: Cheap BFM-SDP Based OPF Relaxation for DNNs
subject to : In the case of an NL model, the load injection term is equal to the summation of the power injection terms (7)- (9). However, in the case of the cheap SDP-based OPF model, the load power injections (7)-(9) have to be transformed into M2 variables. Since in M2, the OV W ii = V i V H i (∀i = 1, · · · , n) is defined for each node, the modelling of constant impedance, constant power and the first term of constant injections can easily be done through this variable. However, for the latter constant current term, which involves the absolute voltage term { I i (|V i | ⊗ |V i (0) |)}, W ii cannot be used directly. Consequently, to model this component, three novel propositions are proposed in section III-B, which ultimately complete the modelling of the constant current component of a ZIP load.
The representation of each ZIP component in terms of the OVs of M2 is presented in the subsequent sections.

A. REPRESENTATION OF CONSTANT IMPEDANCE AND CONSTANT POWER LOADS
The presence of the product of same-phase and differentphase voltage terms in W ii allows for the quick representation of constant impedance and constant power components without requiring any new OV. To this extent, the decoupled power injections related to these components can subsequently be defined as shown in (12)-(13): where W 0 = V int V H int and V int is defined as the vector of the initial voltage at each node and the Z /P/Pc Y / i matrices can be determined by multiplying k Z and k P by the Y / i matrices as defined in (5)-(6).

B. REPRESENTATION OF CONSTANT CURRENT LOADS
The constant current injection (8) contains both the product and absolute value of the voltage variable terms; therefore, additional mathematical modelling is required to represent this component mainly due to the presence of the latter term. The representation of the former term is pretty straight forward and follows the steps described in the previous section for the constant impedance and constant power components. To this end, it can be modelled as: On the other hand, to incorporate the former term, three novel independent mathematical models are developed due to the fact that the approach presented in [28] for the modelling of the absolute part cannot be directly extended to the cheap SDP-based modelling framework.
To begin, consider the last term in (8) for a load connected between the phase a and neutral n of a node i.
where R and IM represent the real and imaginary component of a complex variable, respectively. Note that (16) is the firstorder Taylor series approximation of (15), evaluated at V a 0 . The decoupled injections at phase a and neutral n can now be expressed as follows: Since (17)-(18) contain linear complex-voltage variables that are absent in W ii , additional variables are needed to model this term, which leads to the following three propositions. Please note that the terms involving the subtraction of the phase/neutral voltage in the imaginary part are dropped from the subsequent derivation. This is carried out due to the fact that the computational software does not produce a hard zero and that by introducing extremely small coefficients in the optimization problem (OP), infeasible or suboptimal solutions are obtained, as comprehensively demonstrated in [28].

1) PROPOSITION 1 (P1)
In this proposition, an additional matrix variable ii (19) is introduced at each node i containing constant current loads.
Based on (19), the injection related to (17)- (18) can be expressed as: whereL ϕ,I i,l is calculated as follows: where e ϕ i /c ϕ i is a standard canonical basis vector of R |φ|×1 corresponding to phase ϕ and I D Y / i is a diagonal matrix formed by the diagonal elements of I Y / i . N IL represents the set of nodes containing constant current loads, and Tr represents the trace of a matrix. It should be kept in mind that after obtaining the variableL ϕ,I i,l , it has to be appended with an additional row and column of zeros to have the same dimensions as ii .
Since ii is not used in any equality/inequality constraint of M2, it has to be coupled with the other OVs of this OPF model. Consequently, to enforce the coupling between ii and W ii , the following additional constraints are introduced.
However, constraints (22a)-(22c) do not enforce any coupling between the first-row, first-column variables (FR-FC-Vs) of ii and W ii , which is extremely important due to the involvement of these variables in representing the power injection (20). Consequently, additional bounds (23)- (24) are enforced on the terms of ii in this regard, which restrict the real and imaginary part of each phase voltage within the specified limits. These bounds, in turn, are based on the chosen value of the voltage-angle deviation δ defined for each phase-ground voltage vector. The bound calculation approach is comprehensively explained in the next section.
where i represents the sub-vector ii(2:|η|,1) and V R/IM i is a 3 × 1 vector containing maximum and minimum bounds on real and imaginary components of the phase voltages. Based on (22a)- (24), the diagonal and FR-FC-Vs of ii are bounded.
Please note that in [26], a similar approach to P1 is proposed to model the constant current component in the framework of a single-phase equivalent representation of a power system. However, the OVs introduced in [26] consist of voltage-magnitude terms, and it was successfully demonstrated in [27] that this scheme cannot be generalized for a multi-wire DN configuration. Furthermore, the incorporation of additional constraints (23)-(24) is a novelty of this work since it has been observed that the P1-based OPF model provides a completely meaningless result in the absence of these bounds. Moreover, in [33], a γ -order complex moment relaxation is presented (25), which defines the voltage vector variable associated with an n-bus system as follows: For the first-and second-order complex moment of relaxation, the voltage vector for a 2-bus system can be defined as: Eqs. (26)- (27) indicate that the voltage vector introduced in the P1-based model is completely different from the voltage vector introduced in the complex moment relaxation; therefore, it is not justifiable to compare both methodologies.

2) BOUNDS ON LINEAR COMPLEX VOLTAGE VARIABLES
In practice, bounds on the magnitude of a phase voltage are readily available and are specified by the technical standards. However, as can be seen in (23)-(24), the real and imaginary bounds are required for a complex-voltage variable but are not usually available. Resultantly, a methodology is presented in [34] that derives the bounds on a complex-voltage variable V i and, subsequently, limits its real and imaginary components between −|V i | and |V i |. However, as shown in Fig. 2, this approach encloses a large feasible region (FR) (green region) for V i and allows the solver to calculate the value of V i anywhere between 0 • and 360 • . However, this is not the case in real DNs, where each phase-voltage vector remains within a small region around the phase-voltage vector position in a balanced network scenario (taken as a reference position in this study), i.e., (0 • , -120 • , 120 • ). This region, in turn, can be defined by a voltage-angle deviation δ, which represents the maximum deviation of the phase-voltage vector from its reference position. Since phase voltages are defined with respect to the ground, it is expected that the maximum value of δ at each phase of a node will remain in a small symmetrical range, i.e., ±10 • , and as a result, the bounds presented in [34] can be tightened by making use of the δ parameter. It is extremely important to emphasize that in W ii and ii , phase-ground voltage variables are involved, and these variables must not be confused with the phase-neutral voltage variables. In an unbalanced DN, the phase-neutral voltage angle can assume any value depending upon the degree of unbalance in a network. However, in the presented OP, phase-ground voltage variables are involved; therefore, the assumed range of δ can be used with a higher level of confidence.
Based on δ, the following lower and upper bounds on the real and imaginary voltage components are derived for each phase-voltage variable of ii .
• Phase A: • Phase B: • Phase C: Based on (28)- (33), it can be noted that the FR for each phase voltage is tightened, which leads to a more accurate solution than the solution obtained by setting the bounds described in [34].
The last point to be discussed in this section is the determination of bounds that can be put on the terms involving neutral-voltage variables in ii . The technical standards usually provide the maximum limit on voltage unbalance but do not provide any guidelines on the maximum amount of voltage appearing on a neutral conductor. However, the information on this voltage term may be of significant importance, for example, in the assessment and allocation of losses [35], [36]. Consequently, in this work, the following procedure is adopted to set the minimum and maximum bounds on the real and imaginary components of a neutral voltage. To determine the bounds, a CCI-based load flow model is initially solved, and the real and imaginary components of the neutral voltage are obtained at each node. Based on these components, the parameter σ R/IM is set to the max(|max(V nt R/IM )|, |min(V nt R/IM )|) value, which leads to the bounding of each component of the neutral voltage within the limit (−σ R/IM , σ R/IM ). This is equivalent to imposing the condition that the neutral voltage after solving the OPF problem is not higher than the initial condition.

3) PROPOSITION 2 (P2)
Since P1 involves an additional matrix ii variable at each node containing a constant current load, setting additional McCormick constraints on the bilinear product of off-diagonal terms ( ) of ii leads to a quadratic matrix inequality. This inequality cannot yet be solved due to the lack of a solver for handling this type of constraint. Consequently, to enforce the coupling between the voltage variables corresponding to the different phases of the same node, scalar complex-voltage variables V ϕ i {∀ϕ ∈ φ i , i ∈ N IL } are introduced in P2 since the bilinear product of two scalar variables can, subsequently, be replaced by the McCormick envelopes [37], which ultimately leads to the formation of a convex region. With the introduction of new scalar variables, the following constraints are enforced in P2 to link V ϕ i with the terms of W ii , as shown below: where M stands for the McCormick envelope. Eq. (34) sets the limit on the voltage magnitude of each complex variable V ϕ i by linking it to the diagonal terms of W ii . Please note that the equality constraint cannot be set in (34) since it leads to a quadratic constraint that is non-convex. Eq. (35) links the off-diagonal terms of W ii with the bilinear product of scalar variables. Interestingly, it has been observed that if the OP is solved by taking into consideration only (34), then the real value for each complex-voltage variable is obtained. However, the collective application of (34)- (35) leads to the correct complex value of the phase and neutral voltages. It is also important to emphasize that the bilinear product of phase and neutral voltages is also replaced by the convex region imposed by the McCormick envelopes, where the upper and lower bounds on the neutral voltage are obtained as defined in III-B2. Based on this proposition, the power injection corresponding to the absolute term in (17) can be written as: where represents only the first element of the vector on the right-hand side of (36) and the variableL ϕ,I i,l is calculated as shown in (21).

4) PROPOSITION 3 (P3)
In propositions P1 and P2, additional variables and constraints are introduced, which increase the complexity of the VOLUME 8, 2020 P1-and P2-based OPF models. Consequently, to keep the proposed OPF model as simple as possible, a further approximation is introduced in (16), which models the constant current component using only W ii without introducing any new OV.
In the following discussion, the derivation is shown only for phase a; however, it is equally applicable to the injections corresponding to other phase and neutral conductors. For phase a, (16) can be re-written as: Note that an approximation is introduced in (39) under the assumption that for lightly or moderately unbalanced DNs, the imaginary part of the phase a voltage will be small compared to its real part, i.e., V i R , and therefore, these terms V 2ag i IM and 2V ag i R V ag i IM can be dropped. Nevertheless, the former term is kept in (39) to obtain the product of variables V ag i V * ag i , which in turn is readily available in W ii . Furthermore, it can also be noted that a further approximation is introduced in (41) by applying the firstorder Taylor series to (40), evaluated at V a 0 , which leads to a linear equation that can be easily realized through W ii . In (41), µ is a constant that depends on the point at which the Taylor series is evaluated. In this work, its value is chosen to be 0.5. Based on (41), the approximate injection related to the absolute term in (17) can be written as: Notably, to express (20), (36) and (42) for other phases b and c, their voltage vectors must be rotated towards the reference position of phase a to have the equivalent power injection in these phases correspond to the same load attached to them as connected to phase a. It can also be noted that no additional constraints and variables are needed in this proposition. However, the introduction of two approximations, i.e., (39) and (41), can lead to a suboptimal solution in the case of highly unbalanced DNs serving a large constant current load since the imaginary component of a phase voltage cannot be ignored under these circumstances. The complete injection related to a constant current component can now be expressed as a summation of (14) and S ϕ,Î i : where S ϕ,Î i is a C |φ i |×1 vector representing the injection obtained from the P1- (20), P2-(36) or P3-based (42) approaches.

C. MODELLING OF THE GROUNDING IMPEDANCE
The grounding impedance can be incorporated in the M2 by formulating a ground matrix gnd ii = Y gnd * Y H gnd , having dimensions equal to γ ×γ , at each node i, which contains only a neutral-ground admittance value at the diagonal position corresponding to a neutral conductor. Based on this matrix, the power injections related to the impedance connected between a neutral conductor and ground can be defined as:

D. CONVEXIFIED CHEAP OPF MODEL FOR NEUTRAL-EQUIPPED DISTRIBUTION NETWORKS
The cheap SDP-based OPF model for multi-phase active DNNs containing both wye-and delta-connected loads can now be developed since explicit (decoupled) power injections for each phase and neutral conductor have been derived in the previous sections for a full ZIP load, which ultimately lead to the model M3.

E. OBJECTIVE FUNCTIONS
For the evaluation of the proposed OPF models, minimization of the slack-bus power injection and power system losses by objective functions (OFs) is carried out in this work. In terms of the OVs of M3, these OFs can be expressed as follows:

IV. NUMERICAL TESTS
The quality assessment of the proposed OPF model M3 in terms of the recovered solution, voltage-angle deviation limit and CT is reported in this section for several combinations of load parameters [38], which are shown in Table 1. The proposed approach M3 is applied to the LV CIGRE [39] and Italian-37 (IT-37) [35] active DNs, with their degrees of unbalance reported in Table 3. For all simulation scenarios, V and V are set to 0.90 and 1.05 pu, respectively. The simulations are carried out using the MATLAB-based tool box YALMIP [40] along with the MOSEK solver on a DELL 64-bit OS, core i7 with a processor speed of 2.80 GHz and 16 GB of RAM.

A. CRITERIA FOR THE NUMERICAL EXACTNESS
Two metrics, namely, the optimality gap (OG) and the power mismatch at load (PQ) buses & voltage violation at all nodes, are used to check the optimality and feasibility of the recovered solution, respectively.

1) OPTIMALITY CHECK
The %OG criterion is defined as: The %OG parameter (46) measures the difference between the objective values of the proposed OPF model M3 and its non-convex model M1, which is solved by the non-linear solver IPOPT [41] in this work. The solution from model M3 with an %OG less than 1%, as selected in [5], is considered an optimal solution in this study.

2) FEASIBILITY CHECK
To check the feasibility of the obtained solution, the power mismatch criterion checks the power balance constraint at each load bus once the resultant approach is solved. In the case of a feasible recovered solution, the active and reactive power mismatch at each load bus must be zero. Similarly, the voltage violation criterion checks the recovered voltages at all nodes to ensure that they are within the specified bounds; consequently, the phase voltage at any node outside the specified bounds is considered an infeasible solution.
In this work, the threshold value for the power mismatch criterion is chosen to be 10 −2 [W/VAr], which is lower than the usually specified tolerance value of 10 −1 [W/VAr] in power flow solution packages such as PSS/E.

3) QUALITY OF THE RECOVERED SOLUTION
The quality of the recovered solution from model M3 is cross-checked through the root-mean-square error (RMSE) (47), which compares its value with the results obtained from model M1. The value of this parameter is calculated for phase voltages and branch currents in SI units, and the threshold for characterizing the solution as an accurate one is set as 1 V and 1 A for these state variables, respectively.
where x = {|V|, I br }, M1 and M3 represent models 1 and 3, and P is the length of the vector containing all node voltages/branch currents. It must be noted that the solution-recovery algorithm presented in [20] for three-phase BFM-SDP relaxation is adopted in this work for the 4-wire network OPF model M3.
Finally, please note that the multiple-point grounded configuration is considered for the test DNs by setting R gnd = 1 at each node. Table 3 shows the %OG of each proposition (P1-P3)-based OPF model M3 with respect to the NL model M1. To facilitate the reading process and better clarify the results, Table 2 lists the font styles used in Table 3, as well as their interpretation. Furthermore, please be aware that the results reported in Table 3 are based on a voltage-angle deviation δ of 5 • . The following observations can be made from the mentioned results.  3) The dominance of the proposed P1-, P2-and P3-based OPF approaches over each other strongly depends on the sign of the k I coefficient, i.e., whether the constant current load injects (−k I ) or absorbs (+k I ) power into/from the system. With respect to the −k I coefficient, the P3-based model provides an optimal solution in all these LMs and becomes inexact only when the NL model fails to provide a meaningful solution. Furthermore, it dominates the P1-and P2-based OPF models, as can be noted in the bold-underline cases (OF2), where these OPF models exhibit a large %OG, whereas the P3-based OPF model provides a more accurate (tighter) solution than that of the P1-and P2-based OPF models. However, few cases, as identified by the italic style, are observed where the P2-based model provides a more accurate solution than that of the P3-based model. The slight increase in the %OG for the P3-based model compared to that for the P2-based model is due to the fact that two more approximations are introduced in the former model and, therefore, a less accurate solution is obtained. Nevertheless, the exactness of the P3-based model even in these cases makes it a stronger model among all propositions due to the fact that the P1-and P2-based models can fail for some LMs having −k I coefficients. 4) With respect to the same coefficient, i.e., −k I , the P1-and P2-based OPF models provide almost the same solution in the majority of test cases. However, few cases, as indicated by the underline style, show the superiority of the P1-based model over the P2-based model, especially in the case of OF2, where the P2-based OPF model, although providing an exact solution, shows a larger %OG than that of the P1-based OPF model. This is an unexpected observation since the P2-based OPF approach is expected to be tighter (more accurate) than the P1-based model due to the introduction of additional McCormick envelope terms. The reason for this unanticipated behaviour lies in the chosen value of δ, which appears to be too tight for the scalar complex-voltage variable of proposition 2. However, it has been observed that by slightly relaxing the angle limit from 5 • to 7 • , the solutions, although not shown here, obtained from both models become equal again. 5) With respect to the +k I coefficient, the P1-based OPF model fails for all these LMs, as indicated by the extremely large %OG. The P2-based OPF approach reduces the %OG and provides an optimal solution in the majority of cases, as indicated by the bold-italic style, and therefore exhibits a superior performance to that of the P1-based model. 6) With respect to the +k I coefficient, the P3-based approach further tightens model M3 and provides a more accurate solution than that of the P2-based model. Furthermore, it remains exact and provides a meaningful solution even in the case of large-power-absorbing constant current LMs (such as LM 14) for which a large %OG is observed in the case of the P2-based OPF model, as indicated by the bold-underline style. Nevertheless, few cases are also observed for this LM, where P3 fails, as well in the case of OF1. However, even in these cases, it still provides a more accurate solution than those of the other two propositions. 7) To conclude the discussion with respect to the propositions for constant current component modelling, it is justified to claim that the P3-based OPF model appears to be the strongest among all the presented OPF models. It tends to remain exact under all LMs, having power absorbing as well as power injecting constant current loads. However, its general dominance cannot be established due to the presence of a few cases where the P2-based model dominates it. On the other hand, the P2-based OPF model dominates the P1-based OPF approach in the case of positivecoefficient-based constant current loads, whereas for the negative-coefficient-based loads, the P1 approach can be adopted due to its easy and simple implementation in comparison to the P2-based model.

2) FEASIBILITY AND QUALITY ANALYSIS
In this section, the feasibility and quality of the recovered solution from the P3-based OPF model are discussed since this proposition provides the most optimal solution among   Table 3. 1) With respect to OF1, Fig. 3 shows that the power mismatch for all LMs and load types remains well below the specified threshold level in the case of the  IT-37 network. Similarly, the phase voltage values at all nodes remain within the specified bounds, and no phase-voltage violation is observed in any case. Furthermore, the RMSE values indicate that a highly accurate solution is obtained in this case since the error remains below the 1 V/1 A criterion. Although in the case of a wye-connected load the P3-based model shows a large %OG of 2.96% (above the specified threshold value) for LM 14, the feasibility analysis shows that the obtained solution is still feasible. Nevertheless, as expected, it is less accurate, as evident from the higher RMSE values for phase voltages and branch currents. In this case, the obtained values can be used as the initial values of a multi-wire load flow solver to obtain a more accurate solution. 2) With respect to the same OF for the CIGRE network, the power mismatch criterion is satisfied for all LMs and load types. However, the phase voltage values violate the lower bound for all load types, as shown in Fig. 4. In the case of the Y and Y − load types, voltage violation is observed for LM 4, whereas in the case of the load type, LMs 4 and 6 provide an infeasible solution. In Table 3, for these LMs, either the problem is infeasible or shows a large %OG; consequently, the feasibility analysis results further confirm the recovery of an infeasible solution.  type. Since LM 05 gives an optimal solution, it can be observed in Fig. 7 that the branch currents recovered from both the NL and P3-based models are almost equal to each other. However, this is not the case for LM 14, which has a large %OG; consequently, a significant difference can be observed in Fig. 8 between the recovered current values.

C. IMPACT OF THE VOLTAGE-ANGLE-DEVIATION LIMIT ON THE EXACTNESS OF THE P1-AND P2-BASED OPF MODELS
The voltage-angle-deviation limit δ plays a significant role in the exactness, and lack thereof, of the P1-and P2-based models when applied to LMs that consist of a positive +k I value, as can be seen from the results reported in Table 4. It can be noted that for these LMs, tightening the angle limit from 15 • to 5 • reduces the %OG for both propositions; in particular, the P2-based OPF technique becomes exact for a δ value of 5 • , as indicated by the results reported in the bold font style. In the case of the P1-based OPF model, the small angle limit does reduce the %OG. However, the obtained solution is still far from the solution obtained from the NL VOLUME 8, 2020 OPF model, and consequently, it cannot be considered as an optimal solution. On the other hand, it has also been noted that in the case of LMs comprising constant current components injecting active power into the DNs (−k I values), the impact of the angle limit on the exactness of the P1-and P2-based OPF models is minimal. For these LMs, the difference between the %OG values obtained for all angle limit values is almost negligible; therefore, a moderate value of the δ variable can be selected if an OP involves only these LMs. However, it must be kept in mind that setting an extremely tight angle limit can make the P2-based model less accurate than the P1-based model, as observed in a few cases reported using the underline style format in Table 3.

D. INEXACTNESS OF P1-BASED OPF MODEL UNDER POWER ABSORBING CONSTANT CURRENT LOAD
The sub-optimality of the P1-based OPF approach in the case of +k I components largely depends on the choice of OF. To understand this, it is very important to clarify first what the sign of the load coefficient means. Normally, loads are considered passive and absorb power from a network. For all these loads, the active and reactive power coefficients have positive values, which make the power injection equations (7)- (9) positive; per the standards, these injections are considered as a withdrawal of power from the system. However, it was reported in [38] that the majority of constant current household appliances have a negative coefficient; therefore, these loads inject power into a network. Furthermore, the power injection/absorption of a constant current load is voltage dependent and is therefore highly influenced by the node voltage to which it is connected. Now, under OF2 (power loss minimization) for +k I loads, the solver tends to go towards the lower voltage limits placed on the additional variables ( ii ) of the P1-based model for constant current load modelling to reduce the power absorption of these components. The reduction in the power drawn by these loads subsequently reduces the power flow in the lines/branches of a network, and as a result, the total system losses are reduced, which is the objective of the OP. Similarly, in the case of OF1, which addresses the minimization of slack bus power injection, lowering the voltage on the load buses leads to a reduced power absorption by the constant current loads, which ultimately reduces the power injection into the system from the slack bus.
To better understand this point, consider the results presented in Table 5, which reports the OF2 values obtained in the case of power absorption and injection from a purely constant current load, i.e., LM 2, along with the voltages recovered from ii and the V ϕ i variables used in the P1-and P2-based OPF models, respectively. Furthermore, the minimum and maximum values of real and imaginary bounds corresponding to δ = 5 • are also presented in Table 6, which serve as reference values to cross-check the recovered voltages. Please note that the voltages are presented only for node 8 of the CIGRE network due to the fact that the same trend is observed for other node voltages.
In the case of the +k I coefficient, the real and imaginary components of complex voltages appearing in the first   column of ii are almost identical to the maximum and minimum bounds presented in Table 6, whereas the voltages recovered from W ii of the P1-based model are quite far away from the voltages of ii and do not lie at the extreme bounds of the FR. In terms of polar coordinates, both recovered voltages share the same angle position, but their magnitudes are quite different. The voltage magnitudes in the case of ii are almost equal to the minimum bound value, whereas in the case of W ii of the P1-based model, they are almost identical to the values recovered from the NL OPF model, as shown in Table 7. Since the constant current loads are modelled through the ii variable, the voltage values calculated by the solver lead to a reduced power absorption by these loads, which in turn reduces the power system losses to an unrealistic value. On the other hand, the voltages obtained from V In the case of the −k I coefficient, a constant current load injects power into a system and therefore can be considered as a local power source. In this case, the solver tends to increase the power injections from these loads to reduce the power drawn from the slack bus to bring down the overall system losses. As a result, the optimal voltage values calculated for ii and V ϕ i in the P1-and P2-based OPF models, respectively, are almost identical to each other and equal to the voltage values obtained from the W ii matrices of these models and the NL OPF model variables. These results clearly show why, in the case of negative coefficients, the P1-based OPF approach works well and gives the same optimal solution as obtained by the other OPF models.

E. COMPUTATIONAL PERFORMANCE OF THE PROPOSED OPF APPROACH
The real strength of the cheap SDP-based OPF model in comparison to the primal SDP-based model lies in its superior computational performance, as can be clearly noticed from the CT reported in Table 8 and Fig. 9. For a medium DN (IT-37 bus), the multi-phase primal SDP-based OPF model takes an exceptionally large CT; consequently, the real-time implementation of this OPF model in the context of a distribution management system (DMS) seems difficult to realize. On the other hand, a remarkable reduction in the CT can be observed when the cheap SDP-based OPF model is solved for this DN. It can also be noted that the CT is invariably independent of the load types and LMs, and therefore, no additional techniques such as sparsity exploitation [15], [33] or distributed algorithms [19] are required for the real-time implementation  of this model. Furthermore, for large real DNs [42], the primal SDP-based OPF model is still difficult to solve due to the immature SDP solving technology, as the large dimensions of the OV involved in this OPF model cause the solver to hit its memory limit, as discussed in [28]. Consequently, in these cases, the cheap SDP-based OPF model provides an alternative for the real-time control of these networks.

V. CONCLUSIONS
The cheap SDP-based OPF models for DNNs hosting both wye-and delta-connected ZIP loads provide an alternative solution to the primal SDP-based OPF methodology for realizing the real-time optimal control and management of these networks. The coupled power injections, which are considered the main issue in the development of an OPF model for DNNs, are decoupled between the conductors through the network admittance matrix-based approach. Furthermore, the incorporation of the constant current component of a ZIP load in the OPF model is facilitated by the three novel propositions.
The proposed OPF technique provides an exact solution for all LMs that do not include a constant current component. For LMs consisting of a constant current component, the P3-based OPF model shows promising results compared with the P1-and P2-based OPF models for both +k I and −k I . However, its general dominance cannot be established due to the fact that the P2-based model dominates it, as observed in a few cases with −k I components. On the other hand, the P2-based approach dominates the P1-based model in the case of +k I components since the latter approach fails to completely provide a meaningful solution. For these LMs, the P2-based approach significantly reduces the optimality gap; however, it is even surpassed by the P3-based approach in terms of its ability to provide an optimal solution.
The voltage-angle-deviation limit strongly affects the quality of the P1-and P2-based models since setting a large value for the angle deviation parameter leads to a weaker conic model and, subsequently, to suboptimal results.
With respect to the CT, the developed OPF approach provides an optimal solution in only a fraction of a second even for medium-sized DNs (the primal multi-phase SDP-based technique can take up to three orders of magnitude more time for these networks). Although the primal multi-phase SDPbased relaxation dominates the cheap SDP-based OPF model due to the involvement of a more stringent PSD constraint, the provision of an exact solution in a significantly reduced amount of time by the cheap SDP-based OPF technique allows its practical realization for active DNs to ensure their optimal management. MASSIMILANO COPPO (Member, IEEE) received the Ph.D. degree in electrical energy engineering from the University of Padova, Italy, in 2016. He is currently a Research Associate with the University of Padova. His main research interests include the modeling and simulation of power systems for smart grid management and energy market participation, network stability, and power quality analysis related to the integration of distributed resources in electrical networks.
FABIO BIGNUCOLO (Member, IEEE) is a Research Fellow with the Department of Industrial Engineering, University of Padova. He has authored or coauthored more than 50 articles presented at national and international conferences or published in esteemed international journals. His research interests especially concern computer applications in electrical power engineering, the regulation of distribution networks hosting dispersed generators, innovative control architectures, and the modeling of components and plants. Recently, he has also worked on network applications for electrochemical storage units, aiming at providing ancillary services to transmission and distribution grids.
ROBERTO TURRI (Senior Member, IEEE) was born in Padova, Italy, in 1958. He received the Dr.Ing. degree in electrical engineering from the University of Padova, Italy, in 1984, and the Ph.D. degree from the University of Wales, in 1987. He was with the Physics Department, University College of Swansea, Wales, U.K. In 1990, he joined the Electrical Engineering Department, University of Padova, where he is currently working as an Associate Professor in power systems. His main research interests include power system analysis and simulation, smart grids, and the assessment and mitigation of human exposure to low-frequency electromagnetic fields.