Optimal Dispatch of a Virtual Storage Plant Using Inexact ADMM

A virtual storage plant (VSP) is identified as an effective approach to aggregating distributed storage devices participating in power network supports with similar capability as bulk storage systems. In this study, we develop a distributed control framework for cost-effective storage coordination in the distribution networks, in which the energy storage units are coordinated to contribute to a given power reference at the aggregated level while regulating the local network voltages in the presence of renewable generations. The salient features of the proposed VSP control roots from the successful employment of an inexact alternating direction method of multiplier (ADMM) algorithm, in which the primal updates have analytical solutions in closed form using proximal operators, which significantly reduces the computation efforts of individual storage agents, and renders fast storage dispatch. The proposed control is favorable for near real-time storage dispatch in an optimal manner, and its effectiveness is demonstrated using realistic distributed networks in the simulations.


INTRODUCTION
Distributed energy resources (DERs) transform the conventional role of passive distribution networks in modern power grids. These DERs include dispatchable units, such as diesel generators, energy storage systems, and nondispatchable units, for example, rooftop PV and small wind turbine systems. Although the deployment of renewable energy significantly helps with reduced carbon emissions, with their increased penetration, new network regulation issues appear. Distributed energy storage devices, located closed to the end-users and covering wide geographical areas, are identified as an effective measure to accommodate the intermittency of renewables without compromising the quality of power delivery (Morstyn et al., 2016).
Energy storage systems are widely available at different levels of a power system. As utility-scale storage systems, pumped-hydro storage, compressed-air energy storage, etc., have a long history of participating in energy markets and providing ancillary services to the transmission levels. In the near decades, distributed storage devices are integrated with smaller capacities but larger unit numbers. By the end of 2019, the residential energy storage systems coupled with rooftop PVs reached 2GW, representing a 57% annual increase in Europe; electrical vehicles, whose batteries can potentially contribute to grid supports, attained a 184% annual increase by September 2020 in the United Kingdom (Energy Storage News, 2019). The aggregation of distributed storage devices has a competitive capacity and presents more benefits over bulk storage. Such aggregated storage system is designed as a virtual storage plant (VSP), which seeks to use available storage resources providing flexibility to the large network interconnections.
The coordination of storage devices in a VSP relies on the cyber network for storage communications. This coincides with the concept of virtual power plants (VPPs) for the integration of general DERs (Kim et al., 2019;Cheng et al., 2017). The typical centralized control, in which a control center has access to every storage unit, will soon become inadequate regarding the expansions of the scales of storage in future distribution levels. On the other hand, distributed control collaborates individual agents to settle a given target based on neighboring communication, not only reducing the required bandwidth of the cyber layers but also enhancing robustness against failures and noises of the integrated systems (Olfati-Saber et al., 2007). Existing literature studies on distributed controls of power networks can be categorized into consensus-based controls (Wang et al., 2018Li et al., 2017) and optimization-based distributed controls (Yang et al., 2013a;Dall'Anese et al., 2018;Sulc et al., 2014;Tang et al., 2019;Li et al., 2019). Compared to the consensus controls, distributed optimization can generally achieve more sophisticated objectives, although deriving distributed solutions are not straightforward as that in a centralized setup. Yang et al. (2013b) and Xu et al. (2015) exploited the summation invariant of agent states under special cyber topologies, such that all agents are aggregated to provide a fixed amount of power and meanwhile minimize the overall cost; a similar method is employed for optimal VPP dispatch considering cyber attack in Li et al. (2018). Although these consensus controls are simple and effective in synchronizing agent states, system constraints cannot be easily handled by purely consensus-based methods. The power limits of DERs are enforced in Yang et al. (2013b), and then the algorithm stability is established by generalizing the unconstrained consensus cases. In addition, the optimal VSP operation should account for local voltage variations due to storage dispatch and load variations in the distribution level, but most of the studies (Yang et al., 2013a;Wang et al., 2019;Li et al., 2018;Zheng et al., 2018) were only concerned with the cyber layer and ignored the physical grid that holds the controllable units.
Distributed optimization, on the counterpart, can realize optimal storage dispatch and respect device and network constraints in a systematic manner Zhao and Ding, (2021). Primal-dual algorithms are usually used to derive distributed solutions for a centralized optimization problem by exploring the sparsity in the problem setup. For example, Mallada et al. (2017) developed the optimal frequency control in multi-area power networks in which the proposed control keeps the distributed structure of the typical automatic generation control (AGC) and incorporates frequency dynamics as part of the primal updates. For voltage control in distribution networks, Zhu and Liu (2016) modified the primal updates by a positive definite matrix, which maintains the direction of primal descent and leads to decentralized computation in radial networks. Tang et al. (2019) further employed the dual ascent by minimizing primal variables in each iteration for accelerated convergence but assumptions are made for homogeneous line parameters and 2-hop neighboring communications. A more relevant study by Dall'Anese et al. (2018) investigated optimal operations of a distribution network feeder as a VPP; the local voltage variations are considered, and the primaldual algorithm is derived for the regularized Lagrangian showing better convergence. Dall'Anese et al. (2018), Zhu and Liu (2016), and Tang et al. (2019) used gradient-based methods for variable updates, whose slow converging rates are undesirable in timecritical and real-time applications. Combining the features of dual decomposition and the method of the multiplier, the alternating direction method of multiplier (ADMM) improves algorithm robustness and significantly reduces the needed number of iterations for convergence (Boyd et al., 2010;Feijer and Paganini, 1974). Distributed voltage regulation using ADMM is investigated in Sulc et al. (2014), in which the optimization is formulated in a consensus form with copied variables for fully distributed implementation. Similar concepts are employed for area voltage regulations of a partitioned distribution network (Xu and Wu, 2020), in which the dual update is accelerated to outperform the convergence of the original ADMM.
In Sulc et al. (2014), Xu and Wu (2020), and Zheng et al. (2018), the primal updates of ADMM were to find the minimizer of a constrained optimization problem. This process can be computational-intensive for high-dimensional problems. This disadvantage is first tackled for the linear regression issues in Mateos et al. (2010) using a consensus ADMM, which is further extended in Chang et al. (2015) and Chang (2016) considering coupled constraints by optimizing the dual problems based on inexact ADMM. In this investigation, a VSP is controlled to deliver the requested amount of power at the aggregated level by costeffective storage coordination. This aggregated power reference can either follow an AGC signal or sustain for a specified amount and duration, for example, in the fast frequency reserve (FFR) defined by the National Grid in the United Kingdom (Zhao et al., 2020). In this study, the concerned VSP dispatch is formulated as an optimization problem considering both storage charging and discharging. We use voltage feedback and represent the formulation in an incremental form that accounts for voltage variations against storage dispatch and renewable generations. A proximal dual-consensus ADMM algorithm (PDC-ADMM) (Chang et al., 2015;Chang, 2016) is employed to solve the problem in a distributed manner, whose benefit is that the   minimizer of the primal variables can be represented in a closedform under coupled affine constraints. This results in a lowcomplexity algorithm that reduces the overall computation cost of ADMM. A similar problem setup appears in Chen and Li (2018), in which the primal update still depends on effective solvers to minimize a constrained optimization problem, and the impact of undispatchable units on network operation is ignored. The contribution of this study is that we develop a distributed control framework for optimal storage cooperation through VSP, considering local voltage regulations during storage dispatch based on the inexact PDC-ADMM, in which the variable update is accomplished by per-agent estimates and enjoy analytical forms that reduce the overall computation efforts. The rest of the study is organized as follows: the System Modeling and Problem Formulation section describes the modeling of the cyber-physical systems involving the communication and distribution networks and formulates the optimization problem for VSP dispatch; the Inexact ADMM section illustrates the classical ADMM and the inexact variant of the ADMM algorithms using the proximal operators; Case Study discusses simulation results showing promising performance of the proposed approach in realistic distribution networks with real PV generation data; Conclusions are given in Data Availability Statement at the end.

SYSTEM MODELING AND PROBLEM FORMULATION
In this section, models of the considered cyber-physical system are presented. We assume that the distribution network has a radial topology as for most real-world distribution networks. The consequent DistFlow model is linear with small approximation errors from the nonlinear power flow models (Zhu and Liu, 2016). If more general networks with meshed topologies are

DistFlow Model for Radial Distribution Networks
A radial distribution network is represented by a directed graph G p . Except for the substation node having a constant voltage, buses in the distribution networks are denoted by the set N , among which the buses with storage are collected in the set S, while all the other buses are collected in the set L N \S. Note that the buses in L involve loads and renewables. The set of line segments connecting buses in N are denoted by the edge set Ε p {(i, j) ⊂ N × N }. For every line segment (i, j) ∈ Ε p , the DistFlow model is established as follows: where P ij , Q ij denote the active and reactive power flows on line (i, j); the net injections of the bus j are P n j P g j − P l j , Q n j Q g j − Q l j , where P g j , Q g j denote active/reactive powers from DERs and P l j , Q l j denote the powers of load consumptions. The power flows from the substation to the end buses lead to voltage drops as V i − V j described in (1c), where r l ij , x l ij are resistance and reactance of the distribution lines, N j denotes the set of neighbors of note j. The DistFlow model in (1) can be represented in a more compact form as (Zhu and Liu, 2016) follows: In Equation 2, V is the vector of bus voltages and denoted as V (V j ) j∈N ; similarly, P n (P j ) j∈N , Q n (Q j ) j∈N are vectors of net bus injections, and V 0 v 0 1 |N | where v 0 is the voltage magnitude of node 0 and 1 |N | denotes a vector of all one with dimension |N |, where |N | is the cardinality of N . The matrix M is the incidence matrix of the graph (N , E p ), and D r diag(r l ij ) ij∈E p , D x diag(x l ij ) ij∈E p are diagonal matrices of the line parameters.

Cyber Network
The cyber network communicating distributed storage devices are represented as a graph G c {S, E c } in this study, where the graph edges in E c represent the communication lines among storage agents in S. The connectivity of the graph is described by the adjacency matrix A c (a ij ) i,j∈S where a ij 1 for all j ∈ S i and 0 otherwise, where S i is the neighboring set of node i. We assume that the cyber graph is undirected and connected.

Problem Formulation
The optimization problem is formulated for the optimal storage coordination of a VSP located in a radial distribution network. The VSP is controlled to deliver a given amount of active power while regulating the local network voltage using its reactive power control capability. Considering every time slot t in a time horizon T and ∀i ∈ S, we have min .
In the optimization problem (3), P b i,t , Q b i,t are active and reactive power outputs of the ith storage unit at the current time step. We assign a quadratic cost for every storage unit that is described by the parameters α P i,t , γ P i,t , α Q i,t , γ Q i,t as in (Eq. 3a), and the absolute values indicate that both charging and discharging yield costs are related to storage usage. In (Eq. 3b), R l i (R l ij ) j∈S , X l i (X l ij ) j∈S ∈ R |S| are column vectors of R l , X l corresponding to the ith storage response on other storage buses, which can be understood as the bus voltage sensitivities with respect to the power injections at the ith buses. Here, we separately represent the power injections on the storage buses and load buses in L, whose net powers P n j,t , Q n j,t combine the effects of renewables and loads as in (Eq. 3b); V min , V max are the bus voltage limits. The VSP is optimized to provide the aggregated power P A t as in (Eq. 3c). We mainly focus on inverter-based storage technologies, and the local constraints (3d)-(3e) characterize the limits of storage power capacity where S b i is the apparent power limit of the power converter. Eq. 3f imposes the limits on storage energy capacity based on a general storage model, in which η is the coefficient for charging/discharging, and CA b i is the storage energy capacity. Note that (Eq. 3f) can be incorporated into (Eq. 3e) with time-varying power limits depending on the storage SoC (Dall'Anese et al., 2018); for example, storage with a full charge has P b,max i,t 0 and P b,min i,t −S b i . To facilitate the algorithm implementation, the quadratic constraint (3d) is further linearized as follows: where the accuracy loss is 1.5% when κ 8 (Jabr, 2017). Moreover, considering that only the power injections on storage buses are known by agents, the network voltage constraint (3b) is rewritten in an increment form of storage power, as follows: It is assumed that the algorithm is fast enough such that the powers of loads and renewables are nearly constant within each time slot. From (Eq. 5) we denotẽ which is the vector of voltage magnitudes of storage buses before the agent actions and can be measured at the beginning of the current time step t.
We denote and then the optimization problem (3) is arranged in the canonical form with the incremental variables x i,t , considering storage outputs at t − 1 are known by each agent. where The vector 1 i has the ith entity equal to 1 and 0 for the others. The 1-norm ■ 1 in the objective function (7a) accounts for storage costs in both charging and discharging modes; the coupled equality constraints 3(c) are represented in inequality, which with (Eq. 3b) is incorporated into (7b); (7c) are local constraints for storage capacity limits while Δr i,t ≥ 0 is a slake variable.

INEXACT ADMM
In this section, we introduce the basic principle of ADMM, based on which a proximal dual consensus ADMM is presented to solve convex optimization problem with coupled constraints. The inexact ADMM has the closed form for primal updates and thus reduces computational overheads of individual agents (Chang et al., 2015).

Alternating Direction Method of Multiplier
Consider the following optimization problem: min .
x L ρ x, r k , y k , r k+1 arg min . r L ρ x k+1 , r, y k , The updates of the primal variables are carried out jointly with respect to the augmented Lagrangian in a sequential fashion, while the dual update takes on a gradient step with coefficient ρ. ADMM brings robustness and yields convergence without the restrictive assumptions for primal-dual and dual ascent algorithms (Boyd et al., 2010).

Inexact ADMM
In the classical ADMM, a central coordinator is usually needed to collect the primal variables in (Eq. 9c), and then the updated dual variable is broadcasted to every agent. To make the ADMM completely distributed, a consensus ADMM is proposed by Mateos et al. (2010), which is further developed to solve the problem with coupled constraints via its Lagrange dual problem (Chang et al., 2015;Chang, 2016). For an optimization problem in the form of (7), in which we make the generalization that Δx i,t ∈ R n , Δr i,t ∈ R m and the other matrices have compatible sizes. The objection function of (7) can be separated for individual agents and denoted as where g i,t (Δx i,t ) and h i,t (Δx i,t ) are the smooth and non-smooth parts, respectively. The Lagrange dual problem of (7) is considered and arranged in a consensus form, as follows: min .
y i,t t ij , y j,t t ij , ∀j ∈ S i , where y i,t , z i,t are dual variables associated with the coupled constraints and local constraints. Since the dual variables of the coupled constraints are kept by individual agents, constraint (10b) imposes their consensus, in which t ij is an auxiliary variable. Constraint (10c) is a dummy constraint to make a strongly convex subproblem facilitating the following derivation (Chang, 2016). The augmented Lagrangian of (10) is as follows: where we drop the time stamp t to simplify notions, and variables associated with the consensus and dummy constraints, τ, σ are positive coefficients of the Lagrange multipliers. Applying the ADMM procedure for (11), we have Deriving solutions of (12a)-(12c) and combining with (12d)-(12e) gives t k ij t k ji (y k i + y k j )/2, z k i l k i , ] k i 0. So, (12) can be simplified as follows: (13b) Using the minmax theorem and the strong convexity of (13a) (Chang, 2016), y i , z i can be conveniently obtained by expanding φ i and combining the linear and quadratic terms, which gives where p k i j∈S i (λ k +,ij + λ k −,ji ) and [■] + max(0, ■). The iterative minimization of primal variables in (14a) does not need to be very accurate since it is an intermediate step in ADMM procedures (Mateos et al., 2010). In this regard, the inexact ADMM approximates (14a) by its first-order Taylor expansion at (Δx k−1 i,t , Δr k−1 i,t ) (Chang, 2016). Denote the smooth part of (14a) as follows: To derive an analytical solution, (14a) is simplified as follows: where β i > 0 is some penalty coefficient. Expanding the gradients ∇ xg To solve (17a) with non-smooth 1-norm terms, we employ the proximal operator for each coordinate of Δx i,t . Denoting the jth where Γ i (j) is the jth diagonal element of Γ i and∇ xg is separatable, and notice that (Eq. 18) is a proximal operator for a scaler function h j i,t (Δx i,t (j)), which is precomposed by an affined form. Thus, we rearrange (18) as follows: where the last equation uses the scaling and translation property of the proximal operator (denoted as prox.) for scalar functions (Beck, 2017). The proximal operator of the absolute value is called soft thresholding (Boyd et al., 2010), and thus, (19) can be computed analytically as follows: From (Eq. 20), we recover the optimal VSP control as defined in the System Modeling and Problem Formulation section. This gives the control commands for storage agents: Combining (21) with the primal and dual updates in (17b), (14b)-(14d) makes up the control law of storage devices in a VSP. Assuming a zero-duality gap, and the primal and dual optimal is attainable, the only requirement for the algorithm convergence is that the penalty coefficient β i should be larger than some constant determined by the modulus and Lipschitz constant of the smooth part of the objection function g i,t (Chang, 2016).

CASE STUDY
The performance of the proposed method in optimal storage dispatch is demonstrated in this section. A VSP has an aggregator which receives high-level power reference for participating in transmission-level operations. Note that individual storage agents need to know this VSP power reference P A t , which can be either estimated by a first-order consensus algorithm or broadcasted by the VSP aggregator with low communication costs. We consider three study cases in the simulations, i.e., the IEEE 33-bus system, IEEE-68 bus system (Schneider et al., 2017), and a 199-bus system originally presented in Zhang et al. (2007). All the three distribution networks have radial topology, but the proposed algorithm can be employed for general network scenarios using linear approximation or voltage sensitivity matrix derived in fastdecoupled power flow. The diagram of the IEEE 33-bus system is shown in Figure 1 with an illustration of the VSP concept.

Simulation Setup
The simulations use real PV generation data sampled at 5 min, which are collected from a 75 MW solar power plant in Colorado, United States (Solar Power Data for Integration). The PV generation data from different days are used (Solar Power Data for Integration Studies, 2006) for 25 PV systems covering 24 h, and their profiles are shown in Figure 2. The PV powers are scaled-down by 150 times to suit the capacities of the distribution networks. The buses where the storage and PVs are located are randomly selected, while the loads are assumed to be constant in the simulations. The VSP dispatch commands are each randomly generated in 10 min. The voltage limits for cases 1 and 2 are 0.95-1.05 p.u., and 0.99-1.01 p.u. for case 3. Figure 3 shows the voltage profile of the IEEE 33 and IEEE 69 bus systems with integrated PV systems. This indicates that without storage response, a significant voltage increase is observed during peak PV generation hours. Table 1 gives the storage and control parameters in the simulations. The storage agents are labeled by the same bus number where they are located. We assume a ring topology for storage communication in the VSP cyber network, that id each storage agent communicates to its 2 most adjacent neighbors in the undirect graph. The sampling time of the algorithm is set to 1 min, which gives enough time for algorithm convergence at each time step. The storage power capacities and costs are randomly selected from the ranges provided in Table 1. We assume that all the storages have enough energy in the considered duration. The condition for algorithm termination is the accuracy of the objection function (obj k − obj p )/obj p ≤ 1e −5 , where obj p is the cost of a centralized optimization for the same problem.

Simulation Results
Case I: The IEEE 33-bus system is rated at 12.7 kV with 3.75 MW/2.3MVar loads in total. The PV generations are selected from the data profile and integrated into the network at random buses. Figures 4A, B compares the VSP aggregated powers and costs under the proposed distributed control and a centralized optimization using YALMIP. This illustrates that the proposed VSP control can track the given power reference in a timely fashion and yield almost identical storage costs with respect to a centralized optimization result. The voltage profiles of the system are given in Figures 4C,D. Again, the proposed distributed control, i.e., Figure 4D generates similar results as the centralized scenario in Figure 4C, both of which successfully restrict the voltage within 0.95-1.01p.u., compared to the significant voltage variation observed in Figure 2A. Generally, the computational time of the distributed approach scales linearly with respect to the dimension of the decision variables, while the centralized approach has a high-order polynomial relation between the computational time and the dimension.
Case II: The IEEE 69 bus system is used in this study case. The radial distribution network is rated at 12.7 kV. The total loads are 3.8 MW and 2.7 MVar. Again, the simulation results in Figure 5 indicate that the proposed distributed control accomplishes nearly identical results as a centralized optimization in VSP dispatch and voltage regulations, while the distributed control has fewer requirements on the communication capability of the cyber network.
Case III: The data of the 119-bus system is obtained from Zhang et al. (2007) and the diagram is shown in Figure 6. The nominal voltage of the system is 11 kV, and all the initial tie switches are open to making a radial distribution network. The total power loads are 22.7 MW and 17.0 MVar. In Figure 7, we present the active powers and reactive powers of the total 40 storage units in the VSP. The first column gives the results from a centralized optimization using CPLEX, while the second column is the result using the proposed PDC-ADMM method. It can be observed that the centralized control maintains the storage power references in each dispatch duration, while the proposed distributed control can gradually approach the optimal dispatch using feedback from the previous time step. Figure 8 validates the performance of the proposed VSP control in terms of the power reference tracking and voltage regulations similar to the results in IEEE 33-bus and IEEE 69-bus systems. Note that the proposed control framework generalizes the functionalities of storage solely by providing voltage supports in distribution networks.
The centralized optimization can provide the global optimal regarding the storage dispatch in the convexified optimization problem. It provides a lower bound on the cost of storage coordination. The proposed distributed control framework uses an inexact approach to reduce the computational efforts of each storage agent. So, the accuracy of the optimization is compromised in exchange for the speed of storage dispatch. In addition, for the sake of practical implementation, we set a fixed value for the maximum number of iterations. So, the distributed control approach might output the calculated storage setpoints even when a suboptimal solution is obtained. These are the reasons that cause slight differences in the cost and storage dispatch between the centralized and distributed approach, but the control performance regarding the reference tracking and voltage regulation is almost in the two control scenarios.

CONCLUSION
This study presents a distributed control framework for the optimal coordination of distributed energy storage devices providing active power at the aggregated level in response to the request from system operators. The proposed control also covers the functionality of local voltage regulation that is widely defined by grid codes in the distribution networks. An improved ADMM algorithm is employed to solve the formulated optimization problem considering storage charging/discharging and its active/reactive power control capabilities. This leads to a fully distributed storage dispatch with analytical control laws based on the inexactness of intermediate ADMM calculations. The closed-form solution significantly reduces the computation overhead in algorithm iterations. The simulation results using MATLAB validate that the proposed VSP control can track the given power reference in a cost-effective manner while maintaining a flat voltage profile against renewable variations. In future work, the proposed control framework will be extended to a robust optimization framework that can handle uncertainties from the renewables and loads.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusion of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
DC: algorithm, mathematical formulation, draft preparation; ZW: algorithm and modeling; XW: idea and supervision; YW: manuscript writing and proofreading; WW: review and editing; KZ: data provision and figure plot; DY: proofreading.