Fair Coordination of Distributed Energy Resources with Volt-Var Control and PV Curtailment

1]Daniel Gebbrancor1
[1]School of Electrical and Information Engineering, The University of Sydney, Sydney, Australia
[cor1]Corresponding author

2]Sleiman Mhanna

[2]Department of Electrical and Electronic Engineering, The University of Melbourne, Melbourne, Australia

1]Yiju Ma

3]Archie C. Chapman

1]Gregor Verbič

[3]School of Information Technology and Electrical Engineering, The University of Queensland, Brisbane, Australia

This paper presents a novel distributed optimal power flow (DOPF) method for fair distributed energy resource (DER) coordination in the context of mandated rooftop PV inverter control modes.
In practice, inverter reactive power control is increasingly required by grid connection codes, which often unfairly curtail PV generation of prosumers towards the end of low-voltage feeders.
Similarly, optimization-based DER coordination methods that aim solely for technically-efficient DER coordination do not consider the distribution of PV curtailment across customers.
To address these concerns, we develop a tractable multi-objective DOPF method for optimal DER coordination that (i) curtails PV generation fairly across prosumers, and (ii) incorporates a standard piecewise-linear volt-var control reactive power control function without using integer variables.
Three equity principles representing different interpretations of fairness are implemented in our coordination method; namely, egalitarian, proportional and uniform dynamic PV curtailment redistribution.
The performance of our approach is demonstrated on low-voltage distribution feeders of different sizes (5, 10, 25, 50 and 100 prosumers) using two network topologies: line topology without lateral spurs and tree topology with lateral spurs. Each network considers three levels of PV penetration, giving 30 test systems in total.

The results demonstrate the effectiveness of the proposed DOPF method for fair DER coordination: PV curtailment is equitably distributed among prosumers with a computational burden on par with conventional DOPF approaches. Moreover, different fairness methods result in different patterns of curtailment, which a regulator may choose between.

distributed optimal power flow (DOPF) distributed energy resources (DER) ADMM prosumers demand response Volt-Var Control fair PV curtailment.

[A,8]$\mathcal{X}$Feasible set of network variables for OPF.
[A,a1]$\mathcal{Z}_h$Feasible set of variables for prosumer $h$ for OPF.
[A,9]$\hat{\mathcal{X}}$Feasible set of network variables for DOPF.
[A,a2]$\hat{\mathcal{Z}}_h$Feasible set of variables for prosumer $h$ for DOPF.
[A,6]$\mathcal{H}$Set of prosumers in the power network.
[A,7]$\mathcal{T}$Time horizon of the problem.
[A,5]$\mathcal{B}$Set of buses in the power network.

[A,1]$\boldsymbol{x}$Set of network variables for OPF.
[A,3]$\boldsymbol{z}_h$Set of prosumer variables for OPF.

[A,2]$\hat{\boldsymbol{x}}$Set of network variables for DOPF.
[A,4]$\hat{\boldsymbol{z}}_h$Set of prosumer variables for DOPF.

[B,1]${q}_{\text{g},0,t}$Active power ($\SI{}{\kilo\watt}$) from reference bus at time $t$.
[B,2]$q_{\text{g},0,t}$Reactive power ($\SI{}{\kilo\VAr}$) from reference bus at $t$.
[B,3]$v_{i,t}$Voltage (p.u.) at bus $i$ at $t$.
[B,4]$\theta_{ij,t}$Voltage phase angle between buses $i$ and $j$ at $t$.

[B,5]$p_{h,t}/q_{h,t}$Active/reactive power demand ($\SI{}{\kilo\watt}/\SI{}{\kilo\VAr}$) of prosumer $h$ at $t$.
[B,6]$v_{h,t}$Voltage (p.u.) at prosumer $h$ at $t$.
[B,7]$p^{\text{bat}}_{h,t}$Battery active power ($\SI{}{\kilo\watt}$) of prosumer $h$ at $t$.
[B,8]$p^{\text{PV}}_{h,t}$PV active power generation ($\SI{}{\kilo\watt}$) of prosumer $h$ at $t$.
[B,9]$y_{h,t}$PV curtailment ($\SI{}{\kilo\watt}$) of prosumer $h$ at $t$.
[B,a1]$p_{h,t}^{\text{ch}}/p_{h,t}^{\text{ch}}$Charging/discharging of prosumer's $h$ battery at $t$.
[B,a2]$SoC_{h,t}$State-of-charge ($\SI{}{\kilo\watt\hour}$) of prosumer's $h$ battery, at $t$.

[B,a3]$\overline{y}_{t}$Maximum amount of PVC ($\SI{}{\kilo\watt}$) across all prosumers, at $t$.
[B,a4]$\underline{y}_{t}$Minimum amount of exported power ($\SI{}{\kilo\watt}$) across all participating prosumers, at $t$.
[B,b5]$\lambda^{\bullet}_{h,t}$Dual variable of prosumer $h$ at $t$, associated with given $\bullet$ (voltage, active or reactive power).

[C,1]$b_{ij}$Susceptance (p.u.) of branch $ij$.
[C,2]$g_{ij}$Conductance (p.u.) of branch $ij$.
[C,3]$c_{0}/c_{1}/c_{2}$Constant ($\SI{}{\$}$), linear ($\SI{}{\$\per\kilo\watt}$) and quadratic ($\SI{}{\$\per\kilo\watt\squared}$) coefficient terms of the network cost function.

[C,6]$v_{\text{r},t}/\theta_{\text{r},t}$Reference voltage (magnitude and angle), in p.u.

[C,7]$c^{\text{ToU}}_{t}$Time-of-use tariff ($\SI{}{\$\per\kilo\watt}$) at time interval $t$ of prosumer cost functions.
[C,8]$c^{\text{FiT}}$Feed-in-tariff ($\SI{}{\$\per\kilo\watt}$) of prosumer cost functions.
[C,9]$\overline{s}_{h}$Apparent power rating ($\SI{}{\kilo\VA}$) of prosumer's $h$ inverter.
[C,a1]$\overline{q}_{h}$Reactive power set rating ($\SI{}{\kilo\VAr}$) of prosumer's $h$ inverter.
[C,a2]${V}_{\bullet}$Voltage setting ($\SI{}{\volt}$) of Volt-Var Control equation.
[C,a3]$p_{h,t}^{\text{d}}$Active power demand ($\SI{}{\kilo\watt}$) of prosumer $h$ at $t$.
[C,a4]$\tilde{p}_{h,t}^{\text{PV}}$PV power generation ($\SI{}{\kilo\watt}$) of prosumer $h$ at $t$.
[C,a5]$\eta_{h}^{\bullet}$Battery charging/discharging efficiency of prosumer $h$.
[C,a6]$\Delta t$Time interval ($\SI{}{\hour}$) within $\mathcal{T}$.

[C,a7]$k$Iteration number.
[C,a8]$\rho_{\bullet}$ADMM penalty parameter.
[C,a9]$\alpha$Weighting coefficient for PVC penalty function.
[C,b1]$\beta$Weighting coefficient for reactive power VVC relaxation.
[C,b2]$\gamma$Positive small number to avoid division by 0.
[C,b3]$\epsilon^{\text{abs}}/\epsilon^{\text{rel}}$Absolute and relative feasibility tolerances for ADMM termination criteria.
[C,b4]$\mu/\tau$ADMM residual balancing parameters.

[C,b6]$\underline{\bullet}/\overline{\bullet}$Minimum/maximum magnitude operator.
[C,b7]$\bullet^{+}/\bullet^{-}$Non-negative terms of given $\bullet = \bullet^{+} - \bullet^{-}$.
[C,b8]${\lVert \bullet \rVert}_2$2-norm (Euclidean distance) of a set.

[C,b9]$\bullet^{k}$Previous variables values, used in current ADMM iteration $k$.
[C,c1]$\bullet^{k+1}$Calculated variables values obtained in current ADMM iteration $k$.

[E]ACAlternating current.
[E]ADMMAlternating direction method of multipliers.

[E]DERDistributed energy resources.
[E]DNSPDistribution network service provider.
[E]DOPFDistributed optimal power flow.

[E]MINLPMixed-integer nonlinear programming.

[E]OPFOptimal power flow.

[E]PVCPhotovoltaic curtailment.
[E]SOCSecond-order cone.
[E]VVCVolt-Var Control.

§ INTRODUCTION

The centralized power system operation paradigm is rapidly shifting towards a more distributed one, driven by an increasing uptake of prosumer-owned distributed energy resources (DER), which necessitates new tools for DER management [AEMO, Energy Networks Australia, 2018]. In this context, distributed AC optimal power flow (DOPF) approaches [Scott et al., 2019, Andrianesis and Caramanis, 2019, Attarha et al., 2020, Gebbran et al., 2020] have emerged as an effective tool for DER coordination, for three main reasons: (i) they explicitly consider network constraints, (ii) they permit a prosumer-based decomposition, which preserves prosumer prerogative and privacy, and (iii) they are computationally scalable. However, the current DOPF approaches fail to consider reactive power compensation, which is increasingly required by grid connection codes, and tend to curtail more PV generation from prosumers at the end of low-voltage feeders.

§.§ Background

The AC optimal power flow problem is typically solved using interior-point methods (IPM). Although they only guarantee local optimality, IPM are desirable as they can converge to feasible solutions. In contrast, convex relaxations often converge to infeasible solutions [Molzahn et al., 2017], but are also desirable as they provide a lower bound on the optimal solution, which can be used to assess the quality of the solution obtained from IPM. However, the AC OPF problem quickly becomes intractable when considering a large number of DER and the longer optimization horizons required to accommodate inter-temporal couplings (e.g. battery storage). This therefore encouraged distributed solution approaches.

Conceptually, distributed approaches decompose the optimization problem into several smaller subproblems, where each subproblem can be solved by an individual agent. The solution is reached by iteratively exchanging messages between the agents until the problem converges. Distributed optimization approaches can be hierarchical or fully distributed. The former include a higher level central coordinator, whereas the latter have no central entity.

Distributed optimization techniques can be broadly categorized into: (i) methods based on Lagrangian decomposition, and (ii) techniques for solving the Karush-Kuhn-Tucker conditions in a distributed manner. We provide a brief overview of each category next, but interested readers are referred to [Molzahn et al., 2017] for a comprehensive survey of distributed optimization and control algorithms applied to power systems, including mathematical formulations and convergence properties.

§.§.§ Lagrangian-based methods

Lagrangian-based methods use the Lagrangian function of an optimization problem with a separable structure.
A decomposition is obtained by duplicating (consensus) variables shared between subproblems, e.g., voltages and/or power flows.
Consistency between the copies is enforced by additional equality constraints.

The simplest Lagrangian-based method is dual decomposition [Everett, 1963]. The basic idea is to relax the coupling constraints by multiplying them by Lagrange multipliers (dual variables) and adding them to the objective. The solution procedure then iterates between minimization of the Lagrangian over the primal variables and maximization over the dual variables in a procedure called dual ascent. Because of the decomposable structure, the minimization of the Lagrangian can be parallelized; however the dual update needs to be computed centrally. Convergence of dual ascent is often poor because the dual can be non-smooth.

The poor convergence of the dual decomposition can be improved by augmenting the Lagrangian with a quadratic penalty on the consensus constraints, which smooths the Lagrangian but destroys the separability of the problem. The most popular Augmented Lagrangian-based method is the alternating direction method of multipliers (ADMM) [Boyd et al., 2011], which uses alternate minimizations over the duplicated variables, which restores the separability.

More recently, the Augmented Lagrangian alternating direction inexact Newton (ALADIN) method was proposed to solve distributed optimization problems [Engelmann et al., 2019]. Similar to ADMM, ALADIN combines a coordination step between agents after solving local optimization problems based on the Augmented Lagrangian. The coordination step, however, entails solving an equality-constrained quadratic program, obtained from sensitivity evaluations of all the subproblems (gradients of local objective functions, Hessian approximations of local objective functions and constraints, and Jacobians of constraints). In part due to this improved coordination step, ALADIN offers locally quadratic convergence properties and may converge faster than ADMM, even for non-convex OPF [Engelmann et al., 2019].

In a different approach, the analytical target cascading (ATC) method [Tosserams et al., 2006] makes use of a hierarchical separation in a tree structure, where parent and child subproblems share coupling variables. The coupling constraints are enforced by a quadratic penalty relaxed function. Furthermore, it is possible to use the Augmented Lagrangian modeled within the ATC to be solved at each subproblem, effectively solving the ATC hierarchical subproblems using ADMM.

The auxiliary problem principle (APP) [Cohen, 1980] technique also decomposes an optimization problem into smaller parts, solving the subproblems based on the Augmented Lagrangian and sharing the common variables, and eventually ensuring consensus between them. However, the APP uses linearized cross-terms when modeling the Augmented Lagrangian, which simplifies it and
decouples the subproblems in a way that does not require a central coordinator.

§.§.§ KKT conditions-based methods

An alternative approach to Lagrangian relaxation is the distributed solution of the Karush-Kuhn-Tucker (KKT) optimality conditions.
In the optimality condition decomposition (OCD) [Conejo et al., 2002], each agent performs one Newton step to the KKT conditions for its subproblem in each iteration and shares the resulting primal and dual variables with the neighboring agents. Each primal and dual variable is assigned to a specific subproblem, where all other variables are fixed at their last update. The coupling for these fixed variables is modeled using linear penalties, added to the objective function, which does not require a central coordinator [Molzahn et al., 2017].

The consensus+innovation (C+I) algorithm [Kar et al., 2014] also solves the KKT conditions in a distributed manner without a central coordinator, but allows all variables in each subproblem to vary. Each update takes into account the coupling between the neighboring Lagrange multipliers, in tandem with a local innovation term, enforcing the demand and supply power balance. The limit point of the algorithm satisfies the KKT conditions.

§.§ Decomposition of the Optimal Power Flow

Different decomposition approaches to solve the distributed OPF problem have been studied for more than 20 years [Kim and Baldick, 1997].
They range from region-wise decompositions [Erseghe, 2014, Guo et al., 2017] to more recent component-wise decomposition [Peng and Low, 2018, Mhanna et al., 2019, Mhanna et al., 2018].

Approaches for DER coordination typically use prosumer-based decomposition, decomposing the problem at the connection point between prosumers and the network [Scott et al., 2019, Andrianesis and Caramanis, 2019, Attarha et al., 2020, Gebbran et al., 2020].

In prosumer-based decomposition, the home energy management problem [Azuatalam et al., 2019] used to schedule prosumer-owned DER is solved in one piece. Because the only information that is exchanged with the aggregator is the grid demand, but not the individual device schedules, this preserves prosumer prerogative and privacy.

Furthermore, because the prosumer subproblems can be solved on single-board computers at prosumer premises [Gebbran et al., 2020], this approach can be easily implemented in practice, as demonstrated in a recent Australian trial [Scott et al., 2019].
Existing approaches that use prosumer-based decomposition also solve the network problem in one piece; this results in a significantly lower number of iterations compared to component-based decomposition used, for example, in [Peng and Low, 2018, Mhanna et al., 2019, Mhanna et al., 2018]. The total solution time is well below the DER dispatch intervals, which are typically five minutes or more [Scott et al., 2019, Gebbran et al., 2020].
Finally, solving the network problem in one piece by the aggregator aligns well with the hierarchical structure of DER aggregation where the aggregator assumes the role of a central coordinator.

Comparison of the existing approaches for DER coordination and voltage control with the proposed method.

[Scott et al., 2019, Andrianesis and Caramanis, 2019, Attarha et al., 2020, Gebbran et al., 2020]     [Farivar et al., 2012]     [Tonkoski et al., 2011, Domínguez-García et al., 2012, Haque et al., 2019, Alyami et al., 2014, Liu et al., 2020]     [Zhou et al., 2020]     [Lusis et al., 2019]      [Dall'Anese et al., 2014, Su et al., 2014, Almasalma et al., 2019, Peng and Low, 2018, Engelmann et al., 2019]     This paper

DOPF with Volt-Var control and fair PV curtailment                                   $\times$

DOPF with prosumer-based decomposition     $\times$                              $\times$

Active or reactive inverter dispatch          $\times$     $\times$

Active and reactive inverter dispatch     $\times$                 $\times$     $\times$     $\times$     $\times$

Local voltage control with Volt-Var Control                        $\times$                $\times$

Fair PV Curtailment              $\times$          $\times$           $\times$

§.§ Literature gap

Theoretically, ALADIN allows for a prosumer-based decomposition. However, the computation of the sensitivities for the consensus step makes ALADIN impractical for formulations containing nonsmooth functions, such as piecewise-linear functions and problems with integer variables in general. As shown next, these are required in our problem formulation.
APP, OCD and C+I, on the other hand, do not require a central coordinator, which does not align with the hierarchical prosumer aggregation structure. Further, they require each agent to estimate voltages and power flows in each iteration. This is redundant for prosumer-based approaches, where the aggregator computes the network OPF, obtaining these variables' optima at each iteration, and not just estimates.
Hence, there are no gains from using OCD or C+I methods in prosumer-based DOPF formulations. Moreover, C+I algorithms have not been applied to non-convex AC OPF [Molzahn et al., 2017].
Similarly, ATC offers no benefits because the prosumer-based decomposition has only two layers of subproblems, and does not benefit from further partitioning. In case the prosumer-based DOPF was combined with a higher layer OPF (e.g., a medium- and low-voltage coordination problem partitioned in three levels: MV, LV and prosumers), it is possible that ATC could offer benefits, especially if using ADMM to solve parent-child subproblems.

Given the drawbacks of the existing approaches mentioned above, all the existing prosumer-based DOPF formulations for DER coordination [Scott et al., 2019, Andrianesis and Caramanis, 2019, Gebbran et al., 2020, Attarha et al., 2020] use ADMM. Nevertheless, in the context of Volt-Var Control (VVC) and PV curtailment, ADMM-based approaches suffer from two main pitfalls. First, they do not consider standard inverter reactive power compensation, which is crucial for voltage regulation in medium- and low-voltage networks [Turitsyn et al., 2011], and is increasingly required by grid connection codes [IEE, 2018, AS_, 2015, VDE, 2018].
The standardized VVC curves defined in [IEE, 2018, AS_, 2015, VDE, 2018] are piecewise-linear functions, which adds to the complexity of the already hard, non-convex optimization problem. Therefore, most works that consider reactive power dispatch (in prosumer-based decomposition as well as in other methods) do not account for these standardized functions [Scott et al., 2019, Andrianesis and Caramanis, 2019, Attarha et al., 2020, Gebbran et al., 2020, Lusis et al., 2019, Dall'Anese et al., 2014, Su et al., 2014, Almasalma et al., 2019, Peng and Low, 2018]. The only exception is [Zhou et al., 2020], which uses a control strategy to account for the piecewise-linear VVC in local voltage control instead of using OPF.

Second, a high PV penetration can result in infeasibility of the OPF problem unless PV generation is curtailed. Still, simply allowing the OPF to curtail PV power, as in [Attarha et al., 2020, Gebbran et al., 2020], puts prosumers at the end of the feeder at a disadvantage [Liu et al., 2020]. In countries like Australia, PV penetration has reached levels that require active network management [AEMO, Energy Networks Australia, 2018, Liu et al., 2020], which calls for a solution that distributes curtailment in an equitable manner. A fair redistribution of PV curtailment has been addressed both from a control perspective [Tonkoski et al., 2011, Domínguez-García et al., 2012, Haque et al., 2019] and from an optimization perspective [Liu et al., 2020, Lusis et al., 2019]. However, both approaches [Liu et al., 2020, Lusis et al., 2019] disregard DER coordination. Moreover, the existing methods use proportional PV curtailment only. We propose two additional PV curtailment methods, using different notions of fairness.

§.§ Contributions

Table <ref> summarizes the existing approaches for voltage control in low-voltage networks in the context of DER coordination. They treat different aspects of the problem in isolation, but a unified method is still missing.
Against this background, the major contribution of this paper is a novel OPF formulation for coordination of prosumer-owned DER with a standardized (piecewise-linear) VVC function and a fair PVC, and its implementation in a computationally efficient DOPF. To the best of our knowledge, our formulation is the first to include both a standardized VVC function within an OPF problem, and a fair PV curtailment within a DOPF. In a nutshell, this work is a unified distributed framework that incorporates advantages of previous works and further pushes the state-of-the-art for DER coordination algorithms, while avoiding known disadvantages, as shown in Table I.

In more detail, the proposed method features:

* Standard inverter Volt-Var control (VVC) piecewise-linear function for reactive power compensation embedded in the prosumer subproblem;

* Fair PV curtailment (PVC) using three different notions of fairness resulting in different curtailment patterns, implemented in the DOPF problem;

* Tractable[By “tractable” we mean solvable in practice, not polynomial-time computable as in computational complexity theory.] computational performance that scales well owing to its distributed nature  despite the mixed-integer nonlinear nature of the problem.[In more detail, the problem belongs to the class of mixed-integer nonlinear programming (MINLP) problems that have a non-convex continuous relaxation, which, when compared to MINLP problems with a convex continuous relaxation, are more challenging to solve using current state-of-the-art MINLP technology.]

§.§ Paper structure

This paper is structured as follows: Section <ref> discusses inverter operation modes for voltage control in distribution networks, including the VVC formulation used in the proposed DOPF approach. Section <ref> formulates the DER coordination problem. Section <ref> discusses fair PVC strategies using the three proposed equity principles. Section <ref> details the problem decomposition with PVC redistribution and the resulting formulation of the ADMM-inspired DOPF problem. Section <ref> discusses the algorithm and implementation details. Section <ref> presents detailed results for the 50-bus test system with tree topology, gives a comparison across all 50 simulation cases focusing on the algorithmic performance, and discusses practical implementation considerations.

Finally, Section <ref> gives concluding remarks.

§ PV INVERTERS FOR VOLTAGE CONTROL

In traditional distribution networks, voltage regulation only catered for diurnal energy demand fluctuations, which was done at the substation using transformer tap-changers and capacitor banks [Turitsyn et al., 2011]. With the emergence of distributed PV generation, the slowly responding utility equipment is no longer fit for purpose. Not only does PV generation significantly reshape net demand patterns, it also causes power fluctuations due to rapid changes in solar insolation, which requires voltage control along the feeder.

PV inverters can address these challenges, and are increasingly being used for feeder voltage control. Various standards [IEE, 2018, AS_, 2015, VDE, 2018] offer different modes of operation:

* Fixed Power Factor mode (no voltage regulation);

* Volt-Var Control (VVC) regulates inverter reactive power based on the voltage at the point of network coupling;

* Volt-Watt Control regulates the active power using $P$-$V$ droop control above a certain voltage;

* Power Rate Limit mode limits the rate of change of inverter active power output;

* Voltage Balance mode balances voltages between phases by injecting unbalanced three-phase current.

With the increasing PV penetration, curtailment is required to keep voltages within limits. A network company can impose a limit on the amount of active power that can be injected into the grid, but this is clearly suboptimal. A better solution is to use dynamic curtailment using $P$-$V$ droop control. This can be unfair because the over-voltage at which prosumers are required to curtail PV is not evenly distributed across the grid [Liu et al., 2020].
To achieve a fair PV curtailment,  different approaches have been proposed, including using voltage sensitivity information computed centrally [Tonkoski et al., 2011] or in a decentralized fashion [Alyami et al., 2014] to set adaptive droop control for inverters, and using a consensus algorithm [Haque et al., 2019] to coordinate the curtailment control across prosumers. More recently, the authors in [Liu et al., 2020, Lusis et al., 2019] used a centralized DC OPF formulation to compute fair PV curtailment. In [Lusis et al., 2019], they analyzed the dispatch of optimal active and reactive inverter power considering a fairer curtailment, but without actually distributing all curtailment equally and fairly. This is addressed in [Liu et al., 2020], which also presents a thorough analysis of proportional and economically fair PVC using different metrics.

An alternative approach is to cast the voltage control problem as an OPF problem with inverter reactive power compensation and active power curtailment (APC) as decision variables [Farivar et al., 2012, Dall'Anese et al., 2014, Su et al., 2014, Almasalma et al., 2019]. Authors in [Farivar et al., 2012] and [Ma et al., 2020] formulate a relaxed second-order cone (SOC) OPF, where [Farivar et al., 2012] calculates optimal inverter reactive power compensation to minimize line losses, whereas [Ma et al., 2020] includes inverter VVC within the OPF to account for reactive compensation when formulating an investment game to assess network PV hosting capacity.
Other approaches control both reactive power compensation and active power curtailment: [Su et al., 2014] formulates a SOC OPF for an unbalanced four-wire distribution network;
a Jacobi-Proximal ADMM is employed by [Almasalma et al., 2019] on a decentralized peer-to-peer network, based on voltage sensitivity coefficients instead of solving the actual OPF problem; and [Dall'Anese et al., 2014] formulates an AC OPF, solving the resulting problem using ADMM to minimize power losses and voltage deviations.

Yet another approach is proposed in [Zhou et al., 2020], where the authors study voltage control in distribution networks as a distributed control problem. They show that careful design of local Volt-Var controller is required to ensure convergence. They also demonstrate that the power system dynamics with existing controls can be interpreted as distributed algorithms for solving the optimization problems for reactive power dispatch. The paper, however, does not consider prosumer-owned DER, nor APC.

Our work is based upon VVC and APC.
Unlike [Dall'Anese et al., 2014, Farivar et al., 2012, Su et al., 2014] and similar to [Ma et al., 2020, Zhou et al., 2020], our formulation uses a standard VVC function, using a local voltage-based response included in various international standards [IEE, 2018, AS_, 2015, VDE, 2018] and adopted by all manufacturers in Australia [Australia's Clean Energy Council, 2020].
VVC is typically represented as a deadband graph, with a region close to the nominal voltage under which no response is enacted.
The amount of VAr injection or absorption based on voltage and maximum available reactive power at each connection is given by:

\begin{equation}
q_{t} = \phi \left( v_t \right) = \begin{cases}
\overline{q}, & \text{if} \quad v_{t} \leq {V}_{1}\\
\overline{q} \frac{ V_{2}-v_{t} }{ V_{2}-V_{1} }, & \text{if} \quad {V}_{1} < v_{t} < {V}_{2}\\
0, & \text{if} \quad {V}_{2} \leq v_{t} \leq {V}_{3}\\
-\overline{q} \frac{ v_{t}-V_{3} }{ V_{4}-V_{3} }, & \text{if} \quad {V}_{3} < v_{t} < {V}_{4}\\
-\overline{q}, & \text{if} \quad v_{t} \geq {V}_{4},
\end{cases} \label{eq:HEM_Reactive}
\vspace{-0mm}
\end{equation}

where $\overline{q}$ is the maximum reactive power output (normally defined as a percentage of the inverter's maximum capacity $\overline{s}$), $v_t$ is the voltage amplitude at the point of common coupling, and ${q}_{t}$ is the reactive power output, at any given time $t$. A graphical representation of (<ref>) is shown in Fig. <ref>.

Generic Volt-Var Control (VVC) function.

Different standards suggest or enforce different values for $V_1$ to $V_4$, and implementation follows the requirements of the local distribution network service provider (DNSP) [Australia's Clean Energy Council, 2020].

§ DISTRIBUTED ENERGY RESOURCES COORDINATION

We formulate our proposed approach for DER coordination as a multiperiod OPF problem. Conceptually, it consists of two sources of costs. On the one hand, prosumers schedule their DER to minimize energy expenditure[When ${c}^{\text{ToU}} \! > \! {c}^{\text{FiT}}$, as is the case in Australia, this corresponds to PV self-consumption.]. On the other hand, the distribution system operator coordinates prosumer actions to minimize energy-based operating costs while meeting network and operational constraints. In our formulation, the network cost component is modeled as the cost of power generation, but it could also be modeled as the marginal cost of losses, or the cost of drawing power from the grid.

This is conceptually different from bi-level optimization, where the prosumer problem is modeled explicitly, and its solution is included as a constraint in the aggregator’s problem.
Solving bi-level optimization problems requires reformulating the lower-level problem (the prosumer problem in our case) with its optimality KKT conditions, which are then included in the upper level (aggregator) problem, resulting in a MILP problem; however, this can only be done for convex problems [Conejo and Ruiz, 2020]. Our problem is inherently non-convex due to the AC network constraints and also requires the use of binary variables for the VVC functions, resulting in an MINLP with a non-convex relaxation. Moreover, existing bi-level models (e.g. [Riaz et al., 2019, Zugno et al., 2013]) do not consider network constraints. Solving bi-level models with non-convex mixed-integer problems is an open research problem [Conejo and Ruiz, 2020].

The objective function of our problem formulation is

\begin{equation} \label{eq:Central_OPF_1st}
\begin{split}
&\underset{\boldsymbol{x}, \: \boldsymbol{z}} {\mbox{minimize}} \quad F\left( \boldsymbol{x},\boldsymbol{z}\right)  \coloneqq f \big(  \boldsymbol{x} \big)  + \sum_{h \in \mathcal{H}} g_h\left( \boldsymbol{z}_h\right)  \\
& = \underset{t \in \mathcal{T}}{{\sum}} \left( {c}_{2}{ \big(  {p}_{\text{g},0,t}^{+}\big)  }^2 + {c}_{1} {p}_{\text{g},0,t}^{+} + {c}_{0} + \underset{h \in \mathcal{H}}{\sum} \big(  {c}^{\text{ToU}}_{t}{ {p}_{h,t}^{+} } - {c}^{\text{FiT}}{ {p}_{h,t}^{-} \big)  } \right) ,
\end{split}
\end{equation}

where $f \left( \boldsymbol{x}_{} \right) $ represents the network OPF objective function,
$g_{h} \left( \boldsymbol{z}_{h} \right) $ are prosumer objective functions for each household $h$;
$\mathcal{H}$ is the set of prosumers,
$\boldsymbol{x}$ is the set of network variables (active and reactive power flows, and voltages, for each $t \in \mathcal{T}$),
$\boldsymbol{z}_{h}$ is the set of internal variables of prosumer $h$ for each $t \in \mathcal{T}$ (e.g., battery charging and discharging),
and $\boldsymbol{z} \coloneqq { \left\{ \boldsymbol{z}_h \right\} }_{h\in \mathcal{H}}$.
The network objective $f \left( \boldsymbol{x}_{} \right)$ can include, for example, loss minimization, peak load reduction, or minimizing the use of backup diesel as in [Scott et al., 2019].
The prosumer objective $g_{h} \left( \boldsymbol{z}_{h} \right) $ accounts for a fixed time-of-use (ToU) tariff for purchasing energy, and a feed-in-tariff (FiT) for selling excess PV generation to the grid.

For the sake of simplicity, we assume a balanced three-phase network, which can be modeled by its single-phase equivalent. The model can be readily extended to include unbalanced networks with a combination of single- and three-phase connections as in [Scott et al., 2019].
The network constraints for the OPF are given for each bus $i \in \mathcal{B}$, and for each time interval $t \in \mathcal{T}$ (index 0 denotes the reference bus):

\begin{align}
& p_{\text{g},i,t} - {p}_{h,t} = v_{i,t} \underset{j \: \in \: \mathcal{B}}{\sum} v_{j,t} \left (g_{ij}\cos\theta_{ij,t}+b_{ij}\sin\theta_{ij,t} \right), \label{eq:OPF_P-flow}\\
& q_{\text{g},i,t} - {q}_{h,t} = v_{i,t} \underset{j \: \in \: \mathcal{B}}{\sum} v_{j,t} \left (g_{ij}\sin\theta_{ij,t}-b_{ij}\cos\theta_{ij,t} \right), \label{eq:OPF_Q-flow}\\
% & {p}_{g,t} - {p}_{h,t} = g_{ij}v_{i,t}^2 + v_{i,t} \underset{j \: \in \: \mathcal{B}}{\sum} v_{j,t}(g_{ij}\cos\theta_{ij,t}+b_{ij}\sin\theta_{ij,t}), \label{eq:OPF_P-flow}\\
% %
% & {q}_{g,t} - {q}_{h,t} = b_{ij}v_{i,t}^2 + v_{i,t} \underset{j \: \in \: \mathcal{B}}{\sum} v_{j,t}(g_{ij}\sin\theta_{ij,t}-b_{ij}\cos\theta_{ij,t}), \label{eq:OPF_Q-flow}\\
& v_{0,t} = 1, \quad \theta_{0,t} = 0, \label{eq:OPF_V-slack}\\
& \underline{v}_{i} \leq v_{i,t} \leq \overline{v}_{i}, \label{eq:OPF_V-lim}\\
%& \underline{\alpha}^{\text{g}}{p}_{g,t-1} \leq {p}_{g,t} \leq \overline{\alpha}^{\text{g}}{p}_{g,t-1}, \label{eq:OPF_G-ramp}\\
& \underline{p}_{\text{g},0,t} \leq {p}_{\text{g},0,t} \leq \overline{p}_{\text{g},0,t}, \quad %
\underline{q}_{\text{g},0,t} \leq {q}_{\text{g},0,t} \leq \overline{q}_{\text{g},0,t}, \label{eq:OPF_G-limits}
\end{align}

where $p_{\text{g},i,t}, q_{\text{g},i,t} = 0 \; \forall \; i \in \mathcal{B} \setminus 0$,
${p}_{h,t}, {q}_{h,t}$ are the total net active and reactive power of prosumer $h$ connected to bus $i$,[Note that prosumers have their own dedicated buses. In other words, at most one prosumer is connected to at a certain bus $i$, hence $\mathcal{H} \subset \mathcal{B}$.] and $\theta_{ij,t} = \theta_{i,t} - \theta_{j,t}$ is the angle difference between $i$ and its neighboring bus $j$. Constraints (<ref>) and (<ref>) model the power flow equations, (<ref>) models the reference bus, and (<ref>) and (<ref>) capture voltage and generator limit constraints.

On the other hand, constraints of prosumers consist of the power balance equation

\begin{equation} \label{eq:HEM_Balance}
{p}_{h,t} = {p}_{h,t}^{\text{bat}} + {p}_{h,t}^{\text{d}} - {p}_{h,t}^{\text{PV}}, \quad \forall \; t \in \mathcal{T}, h \in \mathcal{H},
\end{equation}

where ${p}_{h,t}$ is the total power of household $h$ exchanged with the grid, with $\underline{p}_{h,t} \leq {p}_{h,t} \leq \overline{p}_{h,t}$, ${p}_{h,t}^{\text{bat}}$ is the scheduled battery charging power such that $\underline{p}_{h,t}^{\text{bat}} \leq {p}_{h,t}^{\text{bat}} \leq \overline{p}_{h,t}^{\text{bat}}$, ${p}_{h,t}^{\text{d}}$ is the household fixed demand (non-controllable), and ${p}_{h,t}^{\text{PV}}$ is the PV generation power output. Additionally, let ${p}_{h,t} = {p}_{h,t}^{+} - {p}_{h,t}^{-}$, where ${p}_{h,t}^{+}$ and ${p}_{h,t}^{-}$ are non-negative terms that represent the imported and exported active power.[Note that because the second term in (<ref>) is a convex piecewise linear function, at least one of the variables ${p}_{h,t}^{+}$ and ${p}_{h,t}^{-}$ can be zero at time slot $t$. This therefore obviates the need to use binary variables.

Now, we enable a reduction in PV generation by introducing a schedulable curtailment variable $y_{h,t}$.
The inverter's capacity is also considered, which may require additional active power curtailment according to the reactive power response.
The following equations are introduced $ \forall \: t \in \mathcal{T}, h \in \mathcal{H}$:

\begin{align}
& 0 \leq {p}_{h,t}^{\text{PV}} \leq \tilde{p}_{h,t}^{\text{PV}}, \label{eq:HEM_PV_lim} \\
& {y}_{h,t} = \tilde{p}_{h,t}^{\text{PV}} - {p}_{h,t}^{\text{PV}}, \label{eq:HEM_PVC_Gen} \\
& \overline{s}_{h}^{2} \geq { \left( {p}_{h,t}^{\text{PV}} \right) }^2 + {q}^{2}_{h,t}, \label{eq:HEM_Inv_Cap}%%
\end{align}

where $\tilde{p}_{h,t}^{\text{PV}}$ is a parameter denoting the maximum PV generation for prosumer $h$ at $t$, and ${y}_{h,t} \geq 0$ denotes the amount of curtailed PV power.[In our model, we assume batteries and PV have independent inverters. A hybrid inverter would require slight modifications in (<ref>).] For batteries, the constraints are, $ \forall \: t \in \mathcal{T}, h \in \mathcal{H}$:

\begin{align}
& {p}_{h,t}^{\text{bat}} = {p}_{h,t}^{\text{ch}} - {p}_{h,t}^{\text{dis}}, \label{eq:HEM_bat_Dec}\\
& {SoC}_{h,0} \le {SoC}_{h,T}, \label{eq:HEM_SoC0}\\
& {SoC}_{h,t} = {SoC}_{h, t - \Delta t } + \left( \eta_{h}^{\text{ch}} {p}_{h,t}^{\text{ch}} - {{p}_{h,t}^{\text{dis}}}/{\eta_{h}^{\text{dis}}} \right)\Delta t, \label{eq:HEM_bat_SoC}%%
\end{align}

where ${p}_{h,t}^{\text{ch}}, {p}_{h,t}^{\text{dis}} \geq 0$ are the battery charging and discharging powers respectively, ${SoC}_{h,t}$ is the scheduled battery state of charge, with $\underline{SoC}_{h,t} \leq {SoC}_{h,t} \leq \overline{SoC}_{h,t},$[The inclusion of (<ref>) avoids the full depletion of the battery, without considering the subsequent time horizon. Replacing this constraint would be recommended when implementing the algorithm in a rolling horizon basis.] $\eta_{h}$ is the battery charge or discharge efficiency, and $\Delta t$ is the time interval.

Let OPF constraints (<ref>) define a feasible set $\mathcal{X}$ for network variables $\boldsymbol{x}$. Likewise, let prosumer constraints (<ref>), (<ref>), and (<ref>) define a feasible set $\mathcal{Z}_{h}$ for the variables $\boldsymbol{z}_h$ (prosumer's active and reactive power demand, battery charging, state of charge, and PV curtailment). Henceforth, $\boldsymbol{x} \in \mathcal{X}$ and $\boldsymbol{z}_h \in \mathcal{Z}_{h}$, and $\mathcal{Z}$ is the feasible set for all prosumer variables $\boldsymbol{z}$.

We can rewrite the above problem in its compact form, accounting for the VVC coupling between prosumers and the network, which yield the following mixed-integer nonlinear programming (MINLP) problem:[As no two prosumers share the same bus and $\mathcal{H} \subset \mathcal{B}$, as discussed previously, we have $v_{h,t} = v_{i,t}$ for prosumer $h$ connected at bus $i$.]

\begin{align}
& \underset{\boldsymbol{x} \in \mathcal{X}, \: \boldsymbol{z} \in \mathcal{Z}} {\mbox{minimize}} \quad F \left( \boldsymbol{x},\boldsymbol{z} \right) \label{eq:Central_OPF_2nd_Obj} \\
%\intertext{subject to}
%& q_{h,t} = \phi(v_{h,t}), \quad \forall \; h \in \mathcal{H}, \; t \in \mathcal{T}, \label{eq:Central_OPF_2nd_VVC}
& \hspace{0.01cm} \text{subject to} \hspace{0.27cm} q_{h,t} = \phi \left( v_{h,t} \right) , \quad \forall \; h \in \mathcal{H}, \; t \in \mathcal{T}. \label{eq:Central_OPF_2nd_VVC}%
\end{align}

This DER coordination problem is an all-encompassing model that can be used to compute a DER-based OPF with VVC, either centrally or in a distributed fashion. Previous works have utilized similar problems, but [Attarha et al., 2020] and [Gebbran et al., 2020] do not consider the VVC given by (<ref>), whereas [Scott et al., 2019] and [Andrianesis and Caramanis, 2019] do not account for the VVC nor the PVC given by (<ref>), and [Ma et al., 2020] accounts for the VVC in a SOC OPF but does not consider schedulable DER.

However, this problem still lacks consideration of the fair distribution of PV curtailment across prosumers, which we discuss next.

§ FAIR PV CURTAILMENT REDISTRIBUTION

To control the amount of power curtailed at prosumer PV generation, a penalty function is added to (<ref>):

\begin{align}
&\underset{\boldsymbol{x} \in \mathcal{X}, \: \boldsymbol{z} \in \mathcal{Z}} {\mbox{minimize}} \quad  F \left( \boldsymbol{x},\boldsymbol{z} \right) + \alpha \Psi \text{(} . \text{)}, \label{eq:Central_OPF_PVC_Obj}\\
& \hspace{0.01cm} \text{subject to} \hspace{0.27cm} q_{h,t} = \phi \left( v_{h,t} \right) , \quad \forall \; h \in \mathcal{H}, \; t \in \mathcal{T} \label{eq:Central_OPF_PVC_VVC},%
\end{align}

where $\alpha > 0$ is a numerical weight that determines how much a reduction in $\Psi \text{(} . \text{)}$ is favored over the cost minimization $F \left( \boldsymbol{x},\boldsymbol{z} \right) $ or vice versa. For $\alpha=0$ we recover (<ref>).

Next, constraints are introduced to problem (<ref>) to model the type of curtailment, encoding different equity principles: egalitarian, proportional, and uniform dynamic.

§.§ Egalitarian PVC redistribution

In this case, PVC is redistributed equally among all prosumers. An example of this case is shown in Fig. <ref>b). It can be mathematically translated into an egalitarian objective, where the end-goal is to achieve the lowest PVC across all prosumers. The constraint to be included in (<ref>) is:

\begin{equation} \label{eq:PVC_case1}
\overline{y}_{t} \geq {y}_{h,t}, \qquad \forall \; h \in \mathcal{H}, \; t \in \mathcal{T},
\end{equation}

where $\overline{y}_t$ is the highest amount of curtailed PV ${y}_{h,t}$ at any given house, at time period $t$. Correspondingly, we use $\Psi \left( \boldsymbol{\overline{y}} \right) = \sum_{t} \overline{y}_{t}$ as the penalty function.

§.§ Proportional PVC redistribution

In this case, each prosumer is curtailed proportionally to how much they would be able to export on an uncongested operation (which is different from their actual PV generation).[Note that it is also possible to determine a proportional PVC based on installed PV capacities [Liu et al., 2020]. However, the self-interested allocation of PV has been shown to be inefficient [Ma] et al., 2020]. This is shown in Fig. <ref>c). Once again, we use $\Psi \left( \boldsymbol{\overline{y}} \right) = \sum_{t} \overline{y}_{t}$, but for the proportional PVC redistribution, the corresponding constraint to be considered in (<ref>) is:

\begin{equation} \label{eq:PVC_case2_const}
%\overline{y}_{t} \geq \frac{{p}_{h,t}^{-} + {y}_{h,t}}{\underset{h \in \mathcal{H}}{\sum} ({p}_{h,t}^{-} + {y}_{h,t})} {y}_{h,t}, \qquad \forall \; h \in \mathcal{H}, \; t \in \mathcal{T}, % DID NOT WORK WHEN IMPLEMENTED CENTRALLY!
\overline{y}_{t} \geq \frac{{y}_{h,t}}{{p}_{h,t}^{-} + {y}_{h,t}}, \qquad \forall \; h \in \mathcal{H}, \; t \in \mathcal{T},
\end{equation}

where ${p}_{h,t}^{-}$ is the exported power of house $h$ during period $t$. A proportional factor now divides ${y}_{h,t}$, corresponding to the maximum available export power of each prosumer.

§.§ Uniform dynamic PVC redistribution

This case employs a uniform dynamic PVC redistribution that reduces the highest power-exporting prosumers first. Observe in Fig. <ref>d) that the optimisation algorithms enforces a smallest possible uniform export limit. The limit is determined based on the current operating condition, hence the name dynamic. This stands in contrast to the current industry practice of imposing a fixed export limit irrespective of the operating condition. We use $\Psi \left(\underline {\boldsymbol{y}} \right) = \sum_{t} \left( - \underline{y}_{t} \right)$. Its associated constraint is:

\begin{equation} \label{eq:PVC_case3}
\underline{y}_{t} \leq {p}^{{-}}_{h,t}, \quad \forall \; { \left\{ h \in \mathcal{H} \mid \left( {p}^{{-}}_{h,t}+{y}_{h,t} \right) > \underline{y}_{t} \right\} }, \; t \in \mathcal{T},
\end{equation}

where $\underline{y}_t$ is the lowest amount of exported power ${p}^{{-}}_{h,t}$ at any given house subject to curtailment, at time period $t$. In other words, prosumers exporting less than $\underline{y}_t$ are not curtailed.

Examples of PVC redistribution for five households where black indicates power exported and red power curtailed: a) shows the unfair base case (i.e., curtailed as in [Attarha et al., 2020, Gebbran et al., 2020]) where three prosumers only export energy (black), and two prosumers are completely curtailed (red), b) shows an egalitarian redistribution, c) depicts the proportional redistribution, and d) demonstrates the uniform dynamic PVC redistribution. The gray bars show the base case. Note the solid colors retain the same total length in all cases.

§ DISTRIBUTED COORDINATION PROBLEM

The resulting DER coordination problem (<ref>) and the corresponding PVC constraint (<ref>)-(<ref>), account for a fair PVC sharing among prosumers on top of the reactive power compensation from the VVC.

Three main complications arise when solving this problem centrally. First, the privacy of all prosumers is violated because prosumers have to share their energy consumption details with the central coordinator. Second, this problem is difficult to solve, consisting of a large mixed-integer problem with a non-convex continuous relaxation that is computationally challenging in itself [Bienstock and Verma, 2019]. Third, the function introduced by the VVC function (<ref>) is piecewise linear (PWL) and non-convex. By including it in the OPF formulation (<ref>), the problem becomes an MINLP with a number of binary variables proportional to $HT$.

The problem in its central form (<ref>) cannot be solved in a distributed fashion because it is not directly decomposable. Two steps are necessary to decouple common variables between the aggregator and prosumers in order to bestow a decomposable structure, which gives rise to $1+H$ subproblems.

§.§ Problem decomposition

First, constraints associated with VVC (<ref>) and PVC (<ref>)-(<ref>) need to be included in either the network feasible set $\mathcal{X}$ or the prosumer feasible sets $\mathcal{Z}$. As the PVC constraint is associated with the PVC penalty function, it is included in the network set $\mathcal{X}$. Now, the new feasible set for the network is $\hat{\mathcal{X}}$. Similarly, as the VVC-related constraint is related to the inverter's reactive power output, we include this equation in the prosumer sets $\mathcal{Z}$, which then becomes $\hat{\mathcal{Z}}$. Another reason for this attribution of the VVC equation to prosumer sets is to avoid further complicating the already NP-hard network subproblem by adding the mixed-integer variables associated with the PWL function.

Second, since some variables appear in both $\hat{\mathcal{X}}$ and $\hat{\mathcal{Z}}$, the general approach in ADMM-based solutions is to duplicate these variables by assigning one copy for the network and one copy for the respective prosumer. Next, consensus between these copies is captured by, $\forall \ h \in \mathcal{H}, \ t \in \mathcal{T}$,

\begin{align}
\hat{p}_{h,t} &= {p}_{h,t}, \label{eq:Coupling_p}\\
\hat{y}_{h,t} &= {y}_{h,t}, \label{eq:Coupling_PVC}\\
\hat{q}_{h,t} &= {q}_{h,t}, \label{eq:Coupling_q}\\
{v}_{i,t} &= {v}_{h,t}, \: { \Big{\{} i \in \mathcal{B} \mid \exists \; h \; \text{connected to} \; i \Big{\}} }, \label{eq:Coupling_V}%
%{v}_{i,t} &= {v}_{h,t}, \quad \left\{ i \in \mathcal{B} \; | \; \exists \; h \; \text{connected to} \; i \right\}, \label{eq:Coupling_V}%
%i \in B so that agent exists h connected to B
\end{align}

where the left-hand terms are copies for the network problem, $\left\{ \hat{p}_{h,t}, \hat{y}_{h,t}, \hat{q}_{h,t}, {v}_{i,t} \right\} \subset {\boldsymbol{\hat{x}}} \in \hat{\mathcal{X}}$, and the right-hand terms are copies for the prosumer problem, $\left\{ {p}_{h,t}, {y}_{h,t}, {q}_{h,t}, {v}_{h,t} \right\} \subset \hat{\boldsymbol{z}} \in \hat{\mathcal{Z}}$. This confers a decomposable structure on the problem.
In more detail, duplicating the common variables allows us to rewrite problem (<ref>) as

\begin{align}
&\underset{\boldsymbol{\hat{x}} \: \in \: \mathcal{\hat{X}}, \: \boldsymbol{\hat{z}} \: \in \: \hat{\mathcal{Z}}} {\mbox{minimize}} %&
\quad F \left( \boldsymbol{\hat{x}},\boldsymbol{\hat{z}} \right) + \alpha \Psi \text{(} . \text{)}, \label{eq:Decomposed_OPF_Obj}\\
%\quad F(\boldsymbol{\hat{x}}) + \alpha \Psi(.), \label{eq:Decomposed_OPF_Obj}\\
& \hspace{0.01cm} \text{subject to} \hspace{1.30cm} \text{(\ref{eq:Coupling})} \label{eq:Decomposed_OPF_const},%
\end{align}

where $\boldsymbol{\hat{x}} \in \hat{\mathcal{X}}$ is the new set of network variables, and $\boldsymbol{\hat{z}} \in \hat{\mathcal{Z}}$ is the new set of prosumer variables.

In this form, the sets of variables $\mathcal{\hat{X}}$ and $\hat{\mathcal{Z}}$ are decoupled, and (<ref>) is separable if (<ref>) is relaxed (shown next).
We will exploit this decomposable structure to solve the resulting problem in a distributed fashion.

The ADMM [Boyd et al., 2011] capitalises on the decomposable structure in (<ref>) and performs alternating minimizations over sets $\mathcal{\hat{X}}$ and $\hat{\mathcal{Z}}_{h}$. At each iteration $k$, each subproblem calculates the next set of variables denoted by $k+1$. The prosumer subproblems utilize the next set of variables $k+1$ calculated by the network OPF subproblem.

To begin with, we write the Augmented (partial) Lagrange function of the problem:

\begin{equation} \label{eq:Lagrangian}
\small
\begin{split}
%L_{t} (\hat{p}_{h,t},{p}_{h,t},{\lambda}_{h,t}) := & \Big( \frac{\rho}{2} (\hat{p}_{h,t} - {p}_{h,t})^2 + \lambdabar_{h,t} (\hat{p}_{h,t} - {p}_{h,t}) \Big),\\
& L \coloneqq F \left( \boldsymbol{\hat{x}},\boldsymbol{\hat{z}} \right) + \alpha \Psi \text{(} . \text{)} + \underset{t \in \mathcal{T}}{\sum} \underset{h \in \mathcal{H}}{\sum} \Big( \frac{\rho^{\text{p}}}{2} \left( \hat{p}_{h,t} - {p}_{h,t} \right)^2 \\
& \hspace{3mm} + \lambda^{\text{p}}_{h,t} \left( \hat{p}_{h,t} - {p}_{h,t} \right) + \frac{\rho^{\text{y}}}{2} \left( \hat{y}_{h,t} - {y}_{h,t} \right)^2 + \lambda^{\text{y}}_{h,t} \left( \hat{y}_{h,t} - {y}_{h,t} \right) \\
& \hspace{3mm} + \frac{\rho^{\text{q}}}{2} \left( \hat{q}_{h,t} - {q}_{h,t} \right)^2 + \lambda^{\text{q}}_{h,t} \left( \hat{q}_{h,t} - {q}_{h,t} \right) \Big) \\
%& \qquad + \lambda_{q,t} (\hat{q}_{t} - {q}_{t}) \big) \\
& \hspace{3mm} = F \left( \boldsymbol{\hat{x}},\boldsymbol{\hat{z}} \right) +  \alpha \Psi \text{(} . \text{)} +\underset{h \in \mathcal{H}}{\sum} \big( L^{\text{p}}_{h} + L^{\text{y}}_{h} + L^{\text{q}}_{h} \big) \\
& \hspace{3mm} = F \left( \boldsymbol{\hat{x}},\boldsymbol{\hat{z}} \right) +  \alpha \Psi \text{(} . \text{)} + \underset{h \in \mathcal{H}}{\sum} L_{h},
%& \forall \; h \in \mathcal{H}, \; t \in \mathcal{T},%
\normalsize
\end{split}
\end{equation}

where $\rho > 0$ is a penalty parameter associated with each type of variable, and ${\lambda}$ is the dual variable (Lagrangian multiplier) associated with each coupling constraint.[The dual variables for active power are used as locational marginal prices (LMP) in other approaches to the DER coordination problem. See [Papavasiliou], 2018] for a discussion on LMP.

Before we reach the final ADMM-inspired formulation, we need to discuss three additional steps due to the nature of the problem.

§.§ VVC relaxation

The voltage at prosumer $v_{h,t}$ is not determined in any way by the household, but calculated in the network subproblem.

Therefore, instead of taking the Augmented Lagrangian of the voltage, we instead
(i) treat $v^{k+1}_{h,t}$ as an input to the prosumer subproblems, which yields $q^{k+1}_{h,t}$ with (<ref>),
and (ii) include a modified penalty term to the network Lagrangian.

From the first point, the PWL function for the VVC (<ref>) is no longer part of the prosumer optimization subproblem. Instead, it is calculated externally, as is shown in Fig. <ref>. As a result, $q^{k+1}_{h,t}$ becomes a parameter in the prosumer subproblem. We define for each household $\Tilde{L}_h = L^{\text{p}}_{h}+L^{\text{y}}_{h}$.

Algorithm flowchart for DOPF with prosumer-based VVC.

From the second point, we introduce a new penalty parameter that damps oscillations of ${v}_{i,t}$ in subsequent iterations:

\begin{equation}
\Phi \left( \boldsymbol{\hat{q}} \right) = \beta \underset{t \in \mathcal{T}}{\sum} \underset{h \in \mathcal{H}}{\sum} { \big( \hat{q}_{h,t} - \hat{q}_{h,t}^{k} \big) }^2, \label{eq:OPF_Obj_VVC}\\
\end{equation}

where $\beta$ is a weighting factor. This equation penalizes oscillations of reactive power $\hat{q}^{k+1}_{h,t}$, which occur in ill-conditioned problems. This can lead to a violation of constraint (<ref>) that does not appear in the Lagrangian, because it is not a variable in the prosumer subproblems.

This procedure greatly simplifies the network problem by removing the binary variables associated with the non-convex piece-wise linear VVC function.

§.§ PVC constraints

Since (<ref>) is not accounted for in the network subproblem, we need to limit the amount of curtailed power. We do so by adding the following constraint to $\hat{\mathcal{X}}$:

\begin{equation} \label{eq:PVC_ADMM_2}
{p}^{-}_{h,t} + \hat{y}_{h,t} \leq {p}^{-,k}_{h,t} + {y}_{h,t}^{k},  \quad \forall \; h \in \mathcal{H}, \; t \in \mathcal{T},
\end{equation}

which guarantees the sum of exported and curtailed power on the aggregator problem is never larger than the present sum of the exported and curtailed power. Therefore, if a prosumer is generating power but also consuming the same amount of power within its household, it will not be curtailed.

Moreover, a modification is required for the proportional PVC redistribution constraint; the prosumer constraints that model the power and curtailed power are no longer within the OPF (network) subproblem. Therefore, (<ref>) is rewritten using prosumer values from the previous iteration as parameters. We have $\forall \; h \in \mathcal{H}, \; t \in \mathcal{T}$:

\begin{equation} \label{eq:PVC_ADMM_3}
{\overline{y}_{t}} \geq \frac{\underset{h \in \mathcal{H}}{\sum} {y}_{h,t}^{k}} {\left( {p}_{h,t}^{-,k} + {y}_{h,t}^{k} + \gamma \right) \underset{h \in \mathcal{H}}{\sum} \left( {p}_{h,t}^{-,k} + {y}_{h,t}^{k} \right) } {{\hat{y}_{h,t}}},
\end{equation}

where $\gamma > 0$ is a very small value to prevent an undefined operation if no power is exported or curtailed.

§.§ Power flow equations

The curtailed power value is changed within the network subproblem, so a small modification of the power flow equation is needed to capture the resulting difference in each prosumer's power demand.

Hence, we update the left-hand side of the OPF constraint (<ref>), where $\forall \: i \in \mathcal{B}$:

\begin{equation} \label{eq:OPF_const_PF_updated}
\begin{split}
& {p}_{\text{g},i,t} - \left( \hat{p}_{h,t} + \big( \hat{y}_{h,t} - {y}_{h,t}^{k} \big) \right) \\
& = v_{i,t} \underset{j \: \in \: \mathcal{B}}{\sum} v_{j,t} \left( g_{ij}\cos\theta_{ij,t}+b_{ij}\sin \theta_{ij,t} \right),
\end{split}
\end{equation}

which reflects the new power demand from each household; that is, considering the changes within the PVC variable in addition to the updated variable for prosumer's demand $\hat{p}_{h,t}$.

§.§ Resulting ADMM-inspired algorithm

Finally, at a given iteration $k$, the next iterate $k+1$ is computed by solving the following subproblems for primal variables:

\begin{align}
& \boldsymbol{\hat{x}}^{k+1} \coloneqq \underset{\boldsymbol{\hat{x}} \: \in \: \mathcal{\hat{X}}}{\mbox{argmin}} \; \big[ F \left( \boldsymbol{\hat{x}},\boldsymbol{\hat{z}} \right) + \Phi \left( \boldsymbol{\hat{q}} \right) + \alpha \Psi \text{(} . \text{)} + \underset{h \in \mathcal{H}}{\sum} L_{h} \big], \label{eq:xk+1}\\
& {q}^{k+1}_{h,t} \coloneqq \phi \left( {v}^{k+1}_{h,t} \right), \hspace{1.2cm} \forall \; h \in \mathcal{H}, \; t \in \mathcal{T},\label{eq:qk+1} \\
& \boldsymbol{\hat{z}}^{k+1}_{h} \coloneqq \underset{\hat{\boldsymbol{z}}_{h} \: \in \: \hat{\mathcal{Z}}_{h}}{\mbox{argmin}} \; \big[ g \left( {\boldsymbol{\hat{z}}}_{h} \right) + \Tilde{L}_h \big] \hspace{6.2mm} \forall \; h \in \mathcal{H}, \label{eq:zk+1}\\
\intertext{ \normalsize {and update dual variables $\forall \; h \in \mathcal{H}, \; t \in \mathcal{T}$:}} %\small
& {\lambda}^{\text{p},k+1}_{h,t} \coloneqq {\lambda}^{\text{p},k}_{h,t} + \rho^{\text{p}} \left( {\hat{p}_{h,t}}^{k+1} - {p}_{h,t}^{k+1} \right) ,\label{eq:Lk+1} \\
& {\lambda}^{\text{y},k+1}_{h,t} \coloneqq {\lambda}^{\text{y},k}_{h,t} + \rho^{\text{y}} \left( {\hat{y}_{h,t}}^{k+1} - {y}_{h,t}^{k+1} \right) ,\label{eq:L1k+1} \\
& {\lambda}^{\text{q},k+1}_{h,t} \coloneqq {\lambda}^{\text{q},k}_{h,t} + \rho^{\text{q}} \left( {\hat{q}_{h,t}}^{k+1} - {q}_{h,t}^{k+1} \right) ,\label{eq:L2k+1}
%& \qquad {\Lambda}^{k+1}_{h,t} \coloneqq {\Lambda}^{k}_{h,t} + R({\hat{y}_{h,t}}^{k+1} - {y}_{h,t}^{k+1}), \quad \forall \; h \in \mathcal{H}, \; t \in \mathcal{T},\label{eq:L_PVC_k+1}
\end{align}

where (<ref>) is the subproblem solved by the aggregator (holding $\boldsymbol{p},\boldsymbol{y}, \boldsymbol{q}$ constant at $k$), (<ref>) is the VVC update, (<ref>) is the subproblem of each prosumer (holding $\hat{\boldsymbol{p}},\hat{\boldsymbol{y}}, \hat{\boldsymbol{q}}$ constant at $k+1$, results of the network subproblem and VVC update), and (<ref>)-(<ref>) are, respectively, the dual updates for $\boldsymbol{\lambda}^{\text{p}}, \boldsymbol{\lambda}^{\text{y}}$ and $\boldsymbol{\lambda}^{\text{q}}$.

Because problems (<ref>) and (<ref>) are decoupled, they can be solved in parallel, for example, by individual prosumers.

§ IMPLEMENTATION

The proposed coordination method is described in Algorithm 1. The problem (<ref>) is solved on a $\SI{32}{\giga\byte}$ RAM, Intel i7-7700, $\SI{3.60}{\giga\hertz}$ PC. All subproblems were implemented in Python using Pyomo [Hart et al., 2011] as a modeling interface, and solved using IPOPT
v3.12.11 [Wächter and Biegler, 2006], with linear solver MA27 [Smith, 2018]. All subproblems were solved sequentially.

DOPF with VVC and PVC Redistribution
[1]
Initialization: $k=1, \boldsymbol{\lambda}^{1} = 0, \boldsymbol{\rho}^{1} = 1$, $v_{h,t} = v_{i,t} = 1$  p.u., all other variables set to zero.
${\lVert\boldsymbol{r}^{k}\rVert}_2 > \epsilon^{\text{pri}}$ or ${\lVert\boldsymbol{s}^{k}\rVert}_2 > \epsilon^{\text{dual}}$
Compute $\boldsymbol{\hat{x}}^{k+1}$ using (<ref>) at aggregator.
Communicate $\boldsymbol{\hat{p}}_{h}^{k+1}, \boldsymbol{\hat{y}}_{h}^{k+1}, \boldsymbol{v}_{h}^{k+1}, \boldsymbol{\lambda}^{k}_{h}$ to all prosumers.
Calculate (<ref>), in parallel at each prosumer
Compute $\boldsymbol{\hat{z}}^{k+1}$ using (<ref>), in parallel at each prosumer
Communicate all $\boldsymbol{p}_{h}^{k+1}, \boldsymbol{y}_{h}^{k+1}, \boldsymbol{q}^{k+1}_{h}$ to aggregator
Update $\boldsymbol{\lambda}^{k+1}$ using (<ref>)-(<ref>)
Update $\boldsymbol{\rho}^{k+1}$ using  (<ref>)
$k = k + 1$

§.§ ADMM algorithmic details

We use primal and dual residuals to define the stopping criteria [Boyd et al., 2011]. They are, respectively:

\begin{align}
\begin{split}
\boldsymbol{r}^{k} \coloneqq & \big[ \hat{p}_{1,0}^{k+1} - {p}_{1,0}^{k+1}, ..., \hat{p}_{H,T}^{k+1} - {p}_{H,T}^{k+1}, \hat{y}_{1,0}^{k+1} - {y}_{1,0}^{k+1}, ...,\\
& \hat{y}_{H,T}^{k+1} - {y}_{H,T}^{k+1}, \hat{p}_{1,0}^{k+1} - {p}_{1,0}^{k+1}, ..., \hat{q}_{H,T}^{k+1} - {q}_{H,T}^{k+1} \big] \\
=& \boldsymbol{r}_{p}^{k} + \boldsymbol{r}_{y}^{k} + \boldsymbol{r}_{q}^{k}, \label{eq:Residuals_primal}
\end{split} \\
\begin{split}
\boldsymbol{s}^{k} \coloneqq & \big[ \rho_{p} \left( {p}_{1,0}^{k+1} - {p}_{1,0}^{k} \right), ..., \rho_{p} \left( {p}_{H,T}^{k+1} - {p}_{H,T}^{k} \right), \\
& \rho_{y} \left( {y}_{1,0}^{k+1} - {y}_{1,0}^{k} \right) , ..., \rho_{y} \left( {y}_{H,T}^{k+1} - {y}_{H,T}^{k} \right), \\
& \rho_{q} \left( {q}_{1,0}^{k+1} - {q}_{1,0}^{k} \right), ..., \rho_{q} \left( {q}_{H,T}^{k+1} - {q}_{H,T}^{k} \right) \big] \\
=& \boldsymbol{s}_{p}^{k} + \boldsymbol{s}_{y}^{k} + \boldsymbol{s}_{q}^{k}, \label{eq:Residuals_dual}%
\end{split}
%\boldsymbol{s}^{k} &= [{p}_{h,t}^{k+1} - {p}_{h,t}^{k}]
\end{align}

where (<ref>) represent the violations of the coupling <ref> associated with the Augmented Lagrangian (<ref>) at the current iteration, and (<ref>) represents the violation of the Karush-Kuhn-Tucker (KKT) stationarity constraints.

The termination criteria are then given by:

\begin{equation} \label{eq:termination}
{\lVert\boldsymbol{r}^{k}\rVert}_2 \leq \epsilon^{\text{pri}} \quad \text{and} \quad {\lVert\boldsymbol{s}^{k}\rVert}_2 \leq \epsilon^{\text{dual}},
\end{equation}

where $\epsilon^{\text{pri}}$ and $\epsilon^{\text{dual}}$ are feasibility tolerances determined by the following equations [Boyd et al., 2011]:

\begin{align}
\begin{split}
\epsilon^{\text{pri}} & = \sqrt{3HT} \epsilon^{\text{abs}} \\
& + \epsilon^{\text{rel}} \text{max} \big\{ {\lVert{\boldsymbol{\hat{p}}^{k}} + {\boldsymbol{\hat{y}}^{k}} + \boldsymbol{\hat{q}}^{k} \rVert}_2, {\lVert{\boldsymbol{p}^{k}} + {\boldsymbol{y}^{k}} + {\boldsymbol{q}^{k}}\rVert}_2 \big\}, \label{eq:Tolerances_primal}
\end{split}\\[10pt]
& \epsilon^{\text{dual}} = \sqrt{3HT} \epsilon^{\text{abs}} + \epsilon^{\text{rel}} {\lVert{\boldsymbol{\lambda}_{p}^{k}} + {\boldsymbol{\lambda}_{y}^{k}} + {\boldsymbol{\lambda}_{q}^{k}}\rVert}_2 \label{eq:Tolerances_dual}.%
\end{align}

In real-world applications, the standard ADMM may exhibit poor performance due to poor conditioning of the problem. Accelerated methods and adaptive penalty methods have been proposed to speed up the convergence of ADMM, but accelerated methods have been found to show no consistent or substantial improvement on non-convex AC OPF problems. Therefore, we adopt a residual balancing adaptive method [Mhanna et al., 2019], based on updating different $\rho$ independently, according to the relative magnitude of the primal and dual residuals [Boyd et al., 2011]:

\begin{equation} \label{eq:ADMM_mu}
\rho_{p}^{k+1} \coloneqq
\begin{dcases}
\rho_{p}^{k} \left(1+\tau^{\text{incr}} \right) & \text{if} \quad {\lVert\boldsymbol{r}^{k+1}_{p}\rVert}_2 > \mu^{\text{incr}}{\lVert\boldsymbol{s}^{k+1}_{p}\rVert}_2,\\
\rho_{p}^{k}{ \left (1+\tau^{\text{decr}} \right) }^{-1} \: & \text{if} \quad {\lVert\boldsymbol{s}^{k+1}_{p}\rVert}_2 > \mu^{\text{decr}}{\lVert\boldsymbol{r}^{k+1}_{p}\rVert}_2,\\
\rho_{p}^{k} & \text{otherwise},
\end{dcases}
\end{equation}

where $\tau^{\text{incr}}, \tau^{\text{decr}} > 0$, $\mu^{\text{incr}}, \mu^{\text{decr}} > 1$ are parameters that tune how the differences between primal and dual residuals affect $\rho^{k+1}_{p}$. The same formula applies for $\rho_{y}$ and $\rho_{q}$, accordingly.

For our simulations, $\epsilon^{\text{abs}} = 0.001, \epsilon^{\text{rel}} = 0.01$, $\tau^{\text{incr}} = 1.15, \tau^{\text{decr}} = 0.9$, $\mu^{\text{incr}} = \frac{1}{1.15}, \mu^{\text{decr}} = 1.15$ for proportional and uniform dynamic redistribution, $\mu^{\text{incr}} = \frac{1}{0.8}, \mu^{\text{decr}} = 0.8$ for egalitarian redistribution.

Implementations without any form of adaptive methods have shown a very poor performance for any of the coordinated cases, where the number iterations can be orders of magnitude above the adaptive method. This is common among poorly conditioned non-convex problems [Mhanna et al., 2019]. Furthermore, using separate weights $\rho_p, \rho_y$ and $\rho_q$ has also shown significant improvement over using a single weight, which is also a known characteristic of ADMM [Boyd et al., 2011].

The appropriate selection of the multi-objective weights $\alpha$ and $\beta$ may require close attention, and are proportional to the numerical conditioning of the problem. However, this is a general characteristic of multi-objective problems [Boyd et al., 2011]. For $\beta$, values outside the optimal range may lead to an infeasible or slowly converging problem. The optimal range for $\beta$ was found to be from $1$ to $10$.

Exceedingly high values of $\alpha$ can lead to slower convergence, and values too small fail to distribute the PVC across all prosumers. However, as noted in our results and mentioned in [Liu et al., 2020], there is a trade-off between fairly distributing PVC and maximizing energy export, which can be seen as a cost of fairness. Small values for $\alpha$ can be used to obtain different fairness trade-offs. The optimal range of $\alpha$ was found to be from $50$ to $100$. Lower values of $\alpha$ fail to properly distributed the PVC across all prosumers, while higher values increase the number of iterations $k$.

6-, 16-, 26-, 51- and 101-bus networks showing buses, lines, and the generator in red. The cyan, magenta, blue, black and red areas encompass, respectively, 5, 15, 25, 50 and 100 prosumers.

§.§ Test networks

The proposed method was tested on low-voltage distribution feeders of different sizes (6, 11, 26, 51 and 101 buses) using two network topologies: tree topology with lateral spurs and line topology without lateral spurs.
The tree topology with lateral spurs is shown in Fig. <ref>. The network is split into subareas of different sizes with, respectively, 5, 10, 25, 50 and 100 prosumers shown as black dots (the red dot is the connection to the upstream medium voltage network modeled as an equivalent generator).
The line topology is simply a straight line feeder, as in  [Peng and Low, 2018]. Unlike in [Peng and Low, 2018], we did not use the fat-tree topology because PV curtailment of a single prosumer does not affect other prosumers due to the connection to the medium voltage network in between.

A limit for the reactive power limit of $\SI{10}{\percent}$ was used in the simulations ($\overline{q}_h = 0.1 \overline{s}_h$). The inverter capacity $\overline{s}_h$ was set equal to the PV generation capacity, which required active power curtailment according to (<ref>) when the inverter provided reactive power.

The voltage reference is $v_{r} = \SI{230}{\volt}$, and the voltage limits $\underline{v_{i}} = \SI{216}{\volt}$ and $ \overline{{v}_i} = \SI{253}{\volt}$ according to the Australian standard ($\SI{230}{\volt} \pm \SI{6}{\percent}/\SI{10}{\percent} $).

Voltages for VVC are $V_1 = \SI{216}{\volt}, V_2 = \SI{225}{\volt}, V_3 = \SI{244}{\volt}$ and $V_4 = \SI{253}{\volt}$.

We used electricity demand and PV data with half-hourly resolution $\left( \Delta t = 0.5 \right)$ on a spring day (2011/11/07) from a recent Australian smart grid trial [Ratnam et al., 2015].

All prosumers have batteries with maximum capacity of $\SI{7.5}{\kilo\watt\hour}$, charging and discharging efficiency $\eta = 0.92$, and maximum charging or discharging rate of $\SI{3.75}{\kilo\watt}$. The time-of-use tariff ${c}^{\text{ToU}}_{t}$ was $\SI{0.12}{\$\per\kilo\watt\hour}$ off-peak (23:00 - 7:00), $\SI{0.22}{\$\per\kilo\watt\hour}$ shoulder (7:30 - 14:00, 20:00 - 22:30), and $\SI{0.52}{\$\per\kilo\watt\hour}$ peak (14:00 - 19:30), and the feed-in-tariff ${c}^{\text{FiT}}$ was $\SI{0.10}{\$\per\kilo\watt\hour}$. The battery state of charge ate the beginning of the horizon was set to $SoC_{h,0} = \SI{3}{\kilo\watt\hour}$. The APC of scenario B limits exporting up to $\SI{2}{\kilo\watt}$ for each prosumer.
The connection to the medium-voltage network is modeled as a generator with a quadratic cost function with ${c}_{2} =  \SI{0.025}{\$\per\kilo\watt\squared}$ for the three smaller networks and ${c}_{2} =  \SI{0.015}{\$\per\kilo\watt\squared}$ for the two bigger networks, and ${c}_{1} = {c}_{0} = 0$ for all networks.

§.§ Test cases

We applied the three fair PVC methods on each network with three levels of PV penetration: low, medium, and high.

The low and medium cases do not require curtailment, so the choice of the PVC method is not important.
This gives 50 test cases in total, as summarized in Table <ref>.
In summary, we considered the following seven scenarios:

* Uncoordinated: prosumers act in a selfish manner, seeking the best cost according to connection tariffs.

* Uncoordinated, VVC, fixed APC: inverters have VVC enabled. The active power export limit is set to $\SI{2}{\kilo\watt}$ for each prosumer.

* Coordinated: the standard DER coordination problem [Gebbran et al., 2020, Attarha et al., 2020] is utilized.

* Coordinated, VVC: (<ref>) implemented in a distributed manner.

* Coordinated, VVC, PVC #1: implementation of (<ref>) using an egalitarian PVC redistribution.

* Coordinated, VVC, PVC #2: implementation of (<ref>), using a proportional PVC redistribution.

* Coordinated, VVC, PVC #3: implementation of (<ref>), using an uniform dynamic PVC redistribution.

§ SIMULATION RESULTS

We present the results in two parts: we first discuss in detail the results for  the 51-bus network  with  tree topology and 50 prosumers (Fig. <ref> black), followed by a comparison across all 50 cases focusing on the algorithmic performance only.

§.§ 50-prosumer network with tree topology

Fig. <ref> shows a comparison of all seven scenarios detailed in Section <ref>, focusing on the occurrence of voltage violations (top graphs) and the curtailed PV generation (bottom).

Results for the 50-prosumer network with tree topology under the seven proposed scenarios. Figures a-g show the distribution of voltages $v_{h,t}$ over a day, and figures h-n show the amount of energy exported $\sum_{t} \boldsymbol{p}_h^-$, in black, and curtailed $\sum_{t} \boldsymbol{y}_h$, in red, for each prosumer over a day.

Scenarios A and B: As expected in networks with a high PV penetration, the voltage stress is very high. When uncoordinated, prosumers exceed the voltage limits: over $1.3$ p.u., and under $0.93$ p.u. as seen in Fig. <ref>a.
In scenario B, shown in Fig. <ref>b, we still have over-voltage because there is no coordination between agents. Due to the APC, no prosumer exceeds $\SI{25}{\kilo\watt\hour}$ in a day of exported energy, as seen in Fig. <ref>l.

Scenarios C and D: When coordinating prosumers, and curtailing PV generation, we can ensure proper operation within network limits. This is shown in  Fig. <ref>c, where voltages lie within the required limits $0.96$ - $1.1$ p.u. Note, however, that the PVC depends on the electrical distance; that is, prosumers at the end of the feeder are curtailed more. Fig. <ref>j shows that over half of the prosumers are curtailed for the most part of the day.
If we enable the VVC for all inverters, we have less voltage stress on the network, as shown in Fig. <ref>d, but we still have an unfair (albeit smaller) PVC, shown in  Fig. <ref>k.

Scenario E: The curtailment is spread across all prosumers roughly equally (Fig. <ref>l). There are minimal fluctuations of prosumers curtailment, since some prosumers may not export at any given time interval $t$. As proposed, prosumers who are not exporting energy at $t$ are not curtailed.
Observe in Fig. <ref>e how the voltage stress in the network is reduced.

Scenario F: The proportional redistribution spreads the PVC across all prosumers in proportion to the exported power at each interval $t$ (Fig. <ref>m). Once again, there are small variations for prosumers curtailment, because the contributions of individuals are different at each $t$.
Nevertheless, the proportion of energy curtailed versus the sum of curtailed and exported energy ranges between 15% and 30%, which can be attributed to different rooftop PV orientations and demand profiles.

Scenario G: The uniform dynamic redistribution curtails more aggressively prosumers who export more. Prosumers who export less are curtailed very little or not at all, as shown in Fig. <ref>n. The reduction in voltage stress on the network (Fig. <ref>g) is comparable to Scenarios E and F.

§.§.§ General observations

Observe in Table <ref> that the total amount of curtailed energy $\sum y_{h,t}$ is higher among cases with a fair PVC (E, F, and G) when compared to a similar case without a fair PVC (D). The fairness penalty drives the pattern of curtailment away from the most efficient configuration of Scenario D. However, the results show a smaller total curtailed energy when compared to the coordinated case without VVC (C). This is, in turn, proportional to the maximum reactive power $\overline{q}_h$.[Our selected value $\overline{q}_h = 0.1 \overline{s}_h$ is below all Australian DNSP's requirements (which range from $0.3$ to $0.6$) [Australia's Clean Energy Council], 2020], and is already enough to demonstrate the impact of VVC in networks with a high DER penetration.

The impact of fair curtailment and the reactive power compensation using VVC in the original cost function (<ref>) is shown in Table <ref> (columns $F_{\%C}$ and $F_{\%D}$). $F_{\%C}$ shows the difference of $F(\boldsymbol{x},\boldsymbol{z})$ between any given scenario and the coordination without neither VVC nor fair PVC (C), while $F_{\%D}$ shows the same difference with respect to Scenario D (with VVC but without fair PVC). The worst case scenario is Scenario G, where the increase in the cost function is $\SI{15.3}{\percent}$ compared to the optimal result of the DOPF without considering fair PVC (D). As discussed previously, all scenarios outperform the original case without VVC (C).

The difference of total curtailed energy $\sum y_{h,t}$ between scenarios E and F is negligible, but Scenario G has a larger amount of curtailment. This shows that the uniform dynamic PVC redistribution reduces the total exported energy the most.
However, the uniform dynamic PVC redistribution affects
the lowest-exporting prosumers the least. Conversely, the egalitarian redistribution spreads the PVC roughly equally among all energy-exporting prosumers, curtailing lowest-exporting prosumers the most. The proportional redistribution strikes a middle ground between these two.

The last column of Table <ref> shows the average coefficient of variation of the curtailment $\overline{\mathrm{CV}} y_{h,t}$ during the daylight hours of a day.
The variation in scenarios C, D and G are very high, which indicates the curtailment is spread unevenly. This is because prosumers are curtailed according to electrical distances for scenarios C and D, and higher exporting prosumers are curtailed more aggressively in scenario G. Scenario E shows a very low variation (evenly spread curtailment), while Scenario F lies somewhere between these two extremes.

Number of iterations, parallel computation time, sum of PVC, deviation from the optimal cost of the solution without VVC (scenario C) and with VVC (scenario D), and average coefficient of variation of PVC, for all distributed scenarios (C-G).

Scen.     $k$     $\text{t}_\text{Par} [\SI{}{\second}]$     $\sum y_{h,t}$  []     $F_{\%C}$     $F_{\%D}$     $\overline{\mathrm{CV}} y_{h,t}$

C     27     157.8     311.9     N/A     23.8     1.48

D     25     149.0     145.2     -31.2     N/A     1.87

E     19     156.4     201.3     -15.4     12.0     0.19

F     24     171.3     194.9     -18.5     9.7     0.92

G     19     165.6     227.0     -11.1     15.3     1.85

§.§.§ Computational performance

The four distributed scenarios were also simulated in their centralized form (<ref>). Bonmin [Bonami et al., 2008] is used to solve Problem (<ref>).

However, these simulations became intractable. The simulation did not converge for any of the VVC scenarios for the first part of the results (D, E, F, and G, 50-prosumer network) in more than 72 hours. The number of binary variables resulting from the VVC constraint was over a thousand for this network.

The number of iterations $k$ for the algorithm to converge
was actually smaller for scenarios E, F and G, which is shown in Table <ref>. The total parallel computation time $\text{t}_\text{Par}$ (the sum of the slowest parts when solving (<ref>)) is slightly larger, resulting from a longer computation time per iteration due to the larger number of variables. The number of iterations is comparable to other approaches that use prosumer-based decomposition [Scott et al., 2019, Andrianesis and Caramanis, 2019, Attarha et al., 2020].

§.§.§ Network and prosumer power flows

Cumulative PV curtailments $y_{h,t}$ (red) and cumulative power flows through the distribution transformer $p_{g,t}$ (black) for the 50-prosumer network with tree topology for all scenarios.

Fig. <ref> shows cumulative PV curtailments $y_{h,t}$ and cumulative power flows through the distribution transformer $p_{g,t}$.

Observe that the power flow through the distribution transformer varies little across scenarios, with slightly more PV curtailment in cases with PVC redistribution, as discussed in Section <ref>.

Fig. <ref> shows power flows for scenarios D, E, F and G for four selected prosumers: two at the feeder head and two at the end of the feeder, and two with a bigger and two with a smaller PV system.

Prosumers close to the feeder head (Figs. <ref>a and <ref>c) show no PVC in scenario D; that is, they are only curtailed when a fair PVC scheme is in effect.
Meanwhile, prosumers located at the end of the feeder (Figs. <ref>b and <ref>d) are curtailed more in scenario D (no PVC redistribution) when compared to scenarios E, F and G.
Prosumers with bigger PV systems (Fig. <ref>c and <ref>d) are curtailed more in scenarios F (proportional) and G (uniform dynamic) when compared to scenario E (egalitarian). Conversely, prosumers with smaller PV systems (Figs. <ref>a and <ref>b) are curtailed more in scenario E than in scenario F, and are not curtailed at all in scenario G.
These results agree with the results in [Liu et al., 2020], both in terms of the location of the prosumers and the size of their respective PV systems.

Power flows for four selected prosumers for the 50-prosumer network with tree topology for scenarios D, E, F and G. Prosumers are selected based on the size of their PV system (small in cases a and b, big in cases c and d) and their location on the feeder (close to the feeder head in cases a and c, at the end of the feeder in cases b and d).

Number of iterations for fifty simulations with different PV penetration, curtailment redistribution method, network size and topology.

PV Penetration
& Method    Network
Size     $k$
(Line topology)
(Tree topology with
radial spurs)

5*Low     5     13     15

15     12     13

25     26     17

50     27     16

100     26     18

5*Medium     5     23     8

15     7     6

25     29     7

50     30     8

100     28     12

5*High –
Egalitarian    5     28     32

15     31     22

25     35     19

50     36     19

100     45     25

5*High –
Proportional     5     22     18

15     33     18

25     29     16

50     27     24

100     27     19

5*High –
Dynamic     5     17     12

15     16     15

25     20     20

50     23     19

100     22     22

§.§ Comparison of algorithmic performance of all cases

The number of iterations $k$ for each of the 50 test cases is shown in Table <ref>.
Observe that the the number of iterations is not directly related to the size or the type of the network, which is also shown in [Gebbran et al., 2020].

In general, the line topology required a slightly larger number of iterations when compared to the tree topology with radial spurs, which agrees with the results in [Peng and Low, 2018]. Note, however, that [Peng and Low, 2018] uses a component-based decomposition so the total number of iterations is significantly higher. By contrast, we solve the network subproblem in one piece, which reduces the number of iterations significantly at the expense of the increase in the computational burden in each iteration.

The primal residuals ($\boldsymbol{r}_p, \boldsymbol{r}_y, \boldsymbol{r}_q$) across all fifty simulations were on average less than $\SI{0.0015}{\kilo\watt}$, and their maximum values were $\SI{0.050}{\kilo\watt}$, $\SI{0.048}{\kilo\watt}$ and $\SI{0.068}{\kilo\watt}$ for $\boldsymbol{r}_p, \boldsymbol{r}_y, \boldsymbol{r}_q$, respectively. These values are within the same range as in [Scott et al., 2019], which has been successfully deployed in a real-life trial.
Note that a rolling-horizon implementation used in [Scott et al., 2019] would likely result in a faster convergence due to the warm start of each subsequent horizon.

§.§ Practical implementation considerations

§.§.§ Convergence
It is known that distributed optimization approaches for the non-convex AC OPF problem do not theoretically guarantee a globally optimal solution [Molzahn et al., 2017]. However, as discussed in [Molzahn et al., 2017] and demonstrated extensively in [Mhanna et al., 2019, Mhanna et al., 2018], ADMM has been numerically shown to converge on a wide range of OPF cases. Our simulations featuring ten different network configurations under three different levels of PV penetration, for a total of fifty different simulation scenarios, all exhibit good convergence, as evidenced by the vanishing primal and dual residuals.
In other words, the KKT conditions are satisfied since the primal and dual residuals converge within the established tolerance [Boyd et al., 2011]. It is noteworthy that the KKT conditions are necessary to ascertain a (locally) optimal point in a non-convex problem, and sufficient to ascertain a globally optimal point in a convex problem. The robustness shown by the results of our ADMM-based algorithm for different test cases agree with previous results using ADMM for non-convex DOPF problems, whether in large systems decomposed by areas [Guo et al., 2017], in element-wise decomposition [Mhanna et al., 2018] or in prosumer-based decomposition [Scott et al., 2019, Andrianesis and Caramanis, 2019, Attarha et al., 2020, Gebbran et al., 2020].

§.§.§ Parameter Settings
Furthermore, parameters need to be reset between different network sizes/topologies and curtailment redistribution methods (beyond the discussion in subsection 6.1). This is a known requirement for non-convex distributed optimization problems in general, in particular for ADMM [Boyd et al., 2011, Mhanna et al., 2018, Mhanna et al., 2019]. However, once the parameter settings are established for a given system and redistribution method, any scenario can be simulated under the same algorithmic parameter settings. This means that for implementation purposes on a single network, the distribution system operator would only need to tune the parameters once for each desired curtailment redistribution method (or only once, if using a single curtailment method), and proceed to operating the algorithm with that set of parameters for all scenarios.

§.§.§ Hardware implementation
It has been shown that edge computing devices with low processing power are able to compute the prosumer subproblem quickly (in less than a second), which does not impede the deployment of this algorithm in actual distributed settings [Gebbran et al., 2020]. Alternatively, different distributed computation archetypes may be employed, such as grouping a small number of prosumers within an area and assigning one local computing agent to handle their subproblems, as has been implemented in a real-world trial in Australia [Scott et al., 2019].

§.§.§ Communication Requirements
The proposed algorithm requires communication between the aggregator and the prosumers.
The size of the messages exchanged between the agents is smaller than 2 KB, which can be accommodated using the existing last-mile network technologies, such as 4G, or networks tailored for the Internet of Things [Gebbran et al., 2018] such as LTE-M, NB-IoT and EC-GSM-IoT.
These technologies offer sufficiently small latencies relative to the computation time of the DOPF network subproblem so that the communication delay accounts for less than 5% of the total computation time [Gebbran et al., 2020].
Economical aspects and limiting factors such as low area coverage and poor internet connection would probably play a bigger tole in choosing the appropriate network technology.

§ CONCLUSION

This paper has proposed a novel approach for distributed energy resources (DER) coordination via an AC distributed optimal power flow inspired by the alternating direction method of multipliers (ADMM). Compared to existing methods, the proposed solution includes a standard inverter Volt-Var Control (VVC) function for reactive power compensation and a fair PV curtailment (PVC) distribution among prosumers, which have previously been only addressed in isolation. The resulting problem is a mixed-integer nonlinear problem due to the piecewise linear VVC function and the non-convex network problem, and is shown to be intractable if solved centrally. The proposed distributed approach, on the other hand, is tractable because it avoids using integer variables. Moreover, including fair PVC and inverter VVC in the DER coordination problem retains its privacy-preserving characteristic and the computational efficiency, both in terms of the computation time and the number of iterations.

We have demonstrated the effectiveness of the proposed method on 50 different test cases using low-voltage distribution networks of different sizes and topologies, and with different PV penetration levels. The simulation results on a large number of diverse cases demonstrate the robustness of the proposed approach, even though ADMM does not offer convergence guarantees for non-convex problems.

We have proposed three different approaches for PVC redistribution, each adopting different notions of fairness. The results demonstrate that the combination of inverter VVC and PVC effectively reduces voltage stress on the network without penalizing prosumers at the end of the feeder. This stands in contrast to conventional DER coordination approaches that curtail electrically-distant prosumers the most. The overall increase in the cost compared to an optimal solution with unfair allocation PVC (with VVC) is $\SI{15.3}{\percent}$ in the worst-case scenario.

The proposed method requires frequent message exchange between prosumers and the central coordinator, which requires a communication network with sufficiently low latency and high enough bandwidth.  Current network technologies have a sufficiently small latency, but the impact on the communication network traffic  for a very large number of agents needs further research.

[AEMO, Energy Networks Australia, 2018]
authorAEMO, Energy Networks Australia[].
titleOpen Energy Networks.
typeTech. Rep.; year2018.

[Scott et al., 2019]
authorScott[ P.], authorGordon[ D.],
authorFranklin[ E.], authorJones[ L.],
authorThiébaux[ S.].
titleNetwork-aware coordination of residential distributed
energy resources.
journalIEEE Transactions on Smart Grid

[Andrianesis and Caramanis, 2019]
authorAndrianesis[ P.], authorCaramanis[
titleOptimal grid-distributed energy resource
In: booktitle57th Annual Allerton Conference on
Communication, Control, and Computing (Allerton). year2019, p.

[Attarha et al., 2020]
authorAttarha[ A.], authorScott[ P.],
authorThiébaux[ S.].
titleAffinely adjustable robust ADMM for residential
DER coordination in distribution networks.
journalIEEE Transactions on Smart Grid

[Gebbran et al., 2020]
authorGebbran[ D.], authorMhanna[ S.],
authorChapman[ A.C.], authorHardjawana[
W.], authorVucetic[ B.], authorVerbič[
titlePractical considerations of DER coordination with
distributed optimal power flow.
In: booktitleProceedings of the 2020 International
Conference on Smart Grids and Energy Systems (SGES 2020).
year2020a, p. pages1–6.

[Molzahn et al., 2017]
authorMolzahn[ D.K.], authorDörfler[
F.], authorSandberg[ H.], authorLow[
S.H.], authorChakrabarti[ S.],
authorBaldick[ R.], et al.
titleA survey of distributed optimization and control
algorithms for electric power systems.
journalIEEE Transactions on Smart Grid

[Everett, 1963]
authorEverett[ H.].
titleGeneralized lagrange multiplier method for solving
problems of optimum allocation of resources.
journalOperations Research

[Boyd et al., 2011]
authorBoyd[ S.], authorParikh[ N.],
authorChu[ E.], authorPeleato[ B.],
authorEckstein[ J.].
titleDistributed optimization and statistical learning via
the alternating direction method of multipliers.
journalFoundations and Trends in Machine Learning

[Engelmann et al., 2019]
authorEngelmann[ A.], authorJiang[ Y.],
authorMühlpfordt[ T.], authorHouska[
B.], authorFaulwasser[ T.].
titleToward distributed OPF using ALADIN.
journalIEEE Transactions on Power Systems

[Tosserams et al., 2006]
authorTosserams[ S.], authorPapalambros[ P.],
authorEtman[ L.], authorRooda[ J.].
titleAn augmented Lagrangian relaxation for analytical
target cascading using the alternating direction method multipliers.
journalStructural and Multidisciplinary Optimization

[Cohen, 1980]
authorCohen[ G.].
titleAuxiliary problem principle and decomposition of
optimization problems.
journalJournal of Optimization Theory and Applications

[Conejo et al., 2002]
authorConejo[ A.J.], authorNogales[ F.J.],
authorPrieto[ F.J.].
titleA decomposition procedure based on approximate newton
journalMathematical Programming

[Kar et al., 2014]
authorKar[ S.], authorHug[ G.],
authorMohammadi[ J.], authorMoura[
titleDistributed state estimation and energy management in
smart grids: A consensus${+}$ innovations approach.
journalIEEE Journal of Selected Topics in Signal

[Kim and Baldick, 1997]
authorKim[ B.], authorBaldick[ R.].
titleCoarse-grained distributed optimal power flow.
journalIEEE Transactions on Power Systems

[Erseghe, 2014]
authorErseghe[ T.].
titleDistributed optimal power flow using ADMM.
journalIEEE Transactions on Power Systems

[Guo et al., 2017]
authorGuo[ J.], authorHug[ G.],
authorTonguz[ O.K.].
titleA case for nonconvex distributed optimization in
large-scale power systems.
journalIEEE Transactions on Power Systems

[Peng and Low, 2018]
authorPeng[ Q.], authorLow[ S.H.].
titleDistributed optimal power flow algorithm for radial
networks, I: Balanced single phase cas.
journalIEEE Transactions on Smart Grid

[Mhanna et al., 2019]
authorMhanna[ S.], authorVerbič[ G.],
authorChapman[ A.C.].
titleAdaptive ADMM for distributed AC optimal power
journalIEEE Transactions on Power Systems

[Mhanna et al., 2018]
authorMhanna[ S.], authorChapman[ A.C.],
authorVerbič[ G.].
titleComponent-based dual decomposition methods for the
OPF problem.
journalSustainable Energy, Grids and Networks

[Azuatalam et al., 2019]
authorAzuatalam[ D.], authorParidari[ K.],
authorMa[ Y.], authorFörstl[ M.],
authorChapman[ A.C.], authorVerbič[
titleEnergy management of small-scale PV-battery
journalRenewable and Sustainable Energy Reviews

[Gebbran et al., 2020]
authorGebbran[ D.], authorVerbič[ G.],
authorChapman[ A.C.], authorMhanna[ S.].
titleCoordination of prosumer agents via distributed
optimal power flow.
In: booktitleProceedings of the 19th International
Conference on Autonomous Agents and Multiagent Systems (AAMAS 2020).
year2020b, p. pages1–3.

[Farivar et al., 2012]
authorFarivar[ M.], authorNeal[ R.],
authorClarke[ C.], authorLow[ S.].
titleOptimal inverter VAR control in distribution
systems with high PV penetration.
In: booktitle2012 IEEE Power and Energy Society General
Meeting. year2012, p. pages1–7.

[Tonkoski et al., 2011]
authorTonkoski[ R.], authorLopes[
L.A.C.], authorEl-Fouly[ T.H.M.].
titleCoordinated active power curtailment of grid
connected PV inverters for overvoltage prevention.
journalIEEE Transactions on Sustainable Energy

[Domínguez-García et al., 2012]
authorDomínguez-García[ A.D.],
authorCady[ S.T.], authorHadjicostis[
titleDecentralized optimal dispatch of distributed energy
In: booktitle2012 IEEE 51st IEEE Conference on Decision
and Control (CDC). year2012, p. pages3688–3693.

[Haque et al., 2019]
authorHaque[ A.N.M.M.], authorXiong[ M.],
authorNguyen[ P.H.].
titleConsensus algorithm for fair power curtailment of
PV systems in LV networks.
In: booktitle2019 IEEE PES GTD Grand International
Conference and Exposition Asia (GTD Asia). year2019, p.

[Alyami et al., 2014]
authorAlyami[ S.], authorWang[ Y.],
authorWang[ C.], authorZhao[ J.],
authorZhao[ B.].
titleAdaptive real power capping method for fair
overvoltage regulation of distribution networks with high penetration of PV
journalIEEE Transactions on Smart Grid

[Liu et al., 2020]
authorLiu[ M.], authorProcopiou[ A.],
authorPetrou[ K.], authorOchoa[ L.],
authorLangstaff[ T.], authorHarding[ J.],
et al.
titleOn the fairness of PV curtailment schemes in
residential distribution networks.
journalIEEE Transactions on Smart Grid

[Zhou et al., 2020]
authorZhou[ X.], authorFarivar[ M.],
authorLiu[ Z.], authorChen[ L.],
authorLow[ S.].
titleReverse and forward engineering of local voltage
control in distribution networks.
journalIEEE Transactions on Automatic Control
year2020;:pages1–1noteEarly access.

[Lusis et al., 2019]
authorLusis[ P.], authorAndrew[ L.L.H.],
authorChakraborty[ S.], authorLiebman[
A.], authorTack[ G.].
titleReducing the unfairness of coordinated inverter
dispatch in PV-rich distribution networks.
In: booktitle2019 IEEE Milan PowerTech.
year2019, p. pages1–6.

[Dall'Anese et al., 2014]
authorDall'Anese[ E.], authorDhople[
S.V.], authorJohnson[ B.B.],
authorGiannakis[ G.B.].
titleDecentralized optimal dispatch of photovoltaic
inverters in residential distribution systems.
journalIEEE Transactions on Energy Conversion

[Su et al., 2014]
authorSu[ X.], authorMasoum[ M.A.],
authorWolfs[ P.J.].
titleOptimal PV inverter reactive power control and real
power curtailment to improve performance of unbalanced four-wire LV
distribution networks.
journalIEEE Transactions on Sustainable Energy

[Almasalma et al., 2019]
authorAlmasalma[ H.], authorClaeys[ S.],
authorDeconinck[ G.].
titlePeer-to-peer-based integrated grid voltage support
function for smart photovoltaic inverters.
journalApplied Energy
year2019;volume239:pages1037 – 1048.

[Turitsyn et al., 2011]
authorTuritsyn[ K.], authorŠulc[ P.],
authorBackhaus[ S.], authorChertkov[ M.].
titleOptions for control of reactive power by distributed
photovoltaic generators.
In: booktitleProceedings of the IEEE;
vol. volume99. year2011, p.

[IEE, 2018]
titleIEEE Standard for Interconnection and Interoperability of
Distributed Energy Resources with Associated Electric Power Systems
noteIEEE, Std. IEEE 1547-2018.

[AS_, 2015]
titleAustralia and New Zealand Standard for Grid connection of
energy systems via inverters.
noteStandards Australia and Standards New Zealand, Std.
AS/NZS 4777.2:2015.

[VDE, 2018]
titleVDE Application Rule for Generators connected to the
low-voltage distribution network.
noteVerband Deutscher Elektrotechniker, Std. VDE-AR-N

[Ma et al., 2020]
authorMa[ Y.], authorGebbran[ D.],
authorChapman[ A.C.], authorVerbič[
titleA photovoltaic system investment game for assessing
network hosting capacity allocations.
In: booktitleThe Eleventh ACM International Conference on
Future Energy Systems (e-Energy'20). year2020, p.

[Australia's Clean Energy Council, 2020]
authorAustralia's Clean Energy Council[].
titleClean Energy Council submission to the AS/NZS
4777.2:2015 Proposed Changes.
typeTech. Rep.; year2020.

[Conejo and Ruiz, 2020]
authorConejo[ A.J.], authorRuiz[ C.].
titleComplementarity, not Optimization, is the Language
of Markets.
journalIEEE Open Access Journal of Power and Energy

[Riaz et al., 2019]
authorRiaz[ S.], authorMarzooghi[ H.],
authorVerbič[ G.], authorChapman[
A.C.], authorHill[ D.J.].
titleGeneric demand model considering the impact of
prosumers for future grid scenario analysis.
journalIEEE Transactions on Smart Grid

[Zugno et al., 2013]
authorZugno[ M.], authorMorales[ J.M.],
authorPinson[ P.], authorMadsen[ H.].
titleA bilevel model for electricity retailers'
participation in a demand response market environment.
journalEnergy Economics
year2013;volume36:pages182 – 197.

[Bienstock and Verma, 2019]
authorBienstock[ D.], authorVerma[ A.].
titleStrong NP-hardness of AC power flows
journalOperations Research Letters

[Papavasiliou, 2018]
authorPapavasiliou[ A.].
titleAnalysis of distribution locational marginal prices.
journalIEEE Transactions on Smart Grid

[Hart et al., 2011]
authorHart[ W.E.], authorWatson[ J.P.],
authorWoodruff[ D.L.].
titlePyomo: modeling and solving mathematical programs in
journalMathematical Programming Computation

[Wächter and Biegler, 2006]
authorWächter[ A.], authorBiegler[ L.T.].
titleOn the implementation of a primal-dual interior point
filter line search algorithm for large-scale nonlinear programming.
journalMathematical Programming

[Smith, 2018]
authorSmith[ J.].
titleHSL archive: A collection of Fortran codes for
large scale scientific computation.

[Ratnam et al., 2015]
authorRatnam[ E.L.], authorWeller[ S.R.],
authorKellett[ C.M.], authorMurray[ A.T.].
titleResidential load and rooftop PV generation: an
Australian distribution network dataset.
journalInternational Journal of Sustainable Energy

[Bonami et al., 2008]
authorBonami[ P.], authorBiegler[ L.],
authorConn[ A.], authorCornuéjols[ G.],
authorGrossmann[ I.], authorLaird[ C.],
et al.
titleAn algorithmic framework for convex mixed integer
nonlinear programs.
journalDiscrete Optimization

[Gebbran et al., 2018]
authorGebbran[ D.], authorChapman[ A.C.],
authorVerbič[ G.].
titleThe Internet of Things as a facilitator of smart
building services.
In: booktitle2018 Australasian Universities Power
Engineering Conference (AUPEC). year2018, p.