Asynchronous Distributed Power Control of Multi-Microgrid Systems Based on the Operator Splitting Approach
11 Asynchronous Distributed Power Control ofMulti-Microgrid Systems Based onthe Operator Splitting Approach
Zhaojian Wang, Shengwei Mei,
Fellow, IEEE , Feng Liu,
Senior Member, IEEE ,Peng Yi, Ming Cao,
Senior Member, IEEE
Abstract —Forming (hybrid) AC/DC microgrids (MGs) hasbecome a promising manner for the interconnection of variouskinds of distributed generators that are inherently AC or DCelectric sources. This paper addresses the distributed asyn-chronous power control problem of hybrid microgrids, consider-ing imperfect communication due to non-identical sampling ratesand communication delays. To this end, we first formulate theoptimal power control problem of MGs and devise a synchronousalgorithm. Then, we analyze the impact of asynchrony on optimalpower control and propose an asynchronous iteration algorithmbased on the synchronous version. By introducing a randomclock at each iteration, different types of asynchrony are fittedinto a unified framework, where the asynchronous algorithmis converted into a fixed-point problem based on the operatorsplitting method, leading to a convergence proof. We furtherprovide an upper bound estimation of the time delay in thecommunication. Moreover, the real-time implementation of theproposed algorithm in both AC and DC MGs is introduced. Bytaking the power system as a solver, the controller is simplifiedby reducing one order and the power loss can be considered.Finally, a benchmark MG is utilized to verify the effectivenessand advantages of the proposed algorithm.
Index Terms —Asynchronous control, distributed power con-trol; hybrid AC/DC microgrids, time delay.
I. I
NTRODUCTION
Multi-Microgrid systems or Microgrids (MGs) are clustersof distributed generators (DGs), energy storage systems andloads, which are generally categorized into three types: AC,DC and hybrid AC/DC MGs [1]–[3]. A hybrid AC/DC MG hasthe great advantage of reducing processes of multiple inverseconversions in the involved individual AC or DC grid [4]. Inthis paper, we address the distributed power control problemof hybrid AC/DC MGs considering asynchrony.Traditionally, a hierarchical control structure is utilizedin MGs for power control, which is composed of primarycontrol, secondary control and tertiary control [5], usuallyin a centralized way. Such a centralized control architecture,
This work was supported by the National Natural Science Foundation ofChina ( No. 51677100, U1766206, No. 51621065 ). (
Corresponding author:Feng Liu )Z. Wang, F. Liu and S. Mei are with the State Key Laboratory ofPower System and Department of Electrical Engineering, Tsinghua University,Beijing, 100084, China (e-mail: [email protected]).P. Yi is with the Department of Electrical and Systems Engineering,Washington University in St. Louis, 1 Brookings Drive, St. Louis, MO, 63130,USA.M. Cao is with the Faculty of Science and Engineering, University ofGroningen, Groningen 9747 AG, The Netherlands. however, may face great challenges raised by ever-increasinguncertain and volatile renewable generations that require fastresponse of controllers [6]. On the other hand, as MGs usuallybelong to different owners, privacy concerns may prevent thecontrol center acquiring information from individual MGs. Inthis context, breaking the hierarchy of MG control architecturebecomes an emerging research topic, supported by the newidea that real-time coordination could be embedded in the localsteady-state optimization of individual MGs by exchanginginformation only between neighboring MGs. This essentiallyadvocates a distributed control paradigm [7]–[10].Different distributed strategies have been developed in liter-atures for optimal real-time coordination, which can roughlybe divided into two categories in terms of methodology:consensus based methods [11]–[17], and (sub)gradient baseddecomposition methods [18]–[21]. In the consensus basedcontrol, the agents manage to estimate the global variableusing a consensus algorithm [10], [22]. Specifically, in powersystems, the global variable could be the generation ratioand the marginal cost. The former implies that all generatorshave the same generation ratio with respect to its maximalcapability [12]–[14]. The latter implies that all generatorsshare the same marginal cost, and hence the generation con-figuration is economically optimal [15]–[17]. Even thoughthe consensus methods are easy to be implemented, they aredifficult to address complicated (global) constraints. In thissituation, (sub)gradient based decomposition methods couldbe applied, where the optimization problem is solved by dual(sub)gradient ascent [18]–[21], [23].Although the aforementioned methods have achieved greatsuccess, there are still two issues to be hurdled before a practi-cal implementation. First, most of them have considered onlysynchronous distributed control. Therefore, all MGs must carryout computation simultaneously, implying that a global clockis necessary to ensure the instants for control actions gettingstrictly synchronized. This is computationally inefficient andimpractical since in each iteration all MGs have to wait forthe slowest one to finish before executing their local actionsin the next iteration. In fact, asynchrony widely exists inpower systems, such as time delays and non-identical samplingrates [24], [25]. Hence, the synchrony requirement limits theapplication of distributed control. Second, the load demand inexisting literature is usually assumed to be known, especiallyin those papers addressing optimal economic dispatch problem[17], [19], [21]. However, the load demand is very difficult to a r X i v : . [ c s . S Y ] M a y measure and is always time varying when demand responseand electric vehicles exist. In addition, the fast varying situa-tion requires the real-time implementation of the algorithm.This somewhat makes the existing method difficult to beapplied. Thus, it is significant to investigate the controller thatadapts to asynchronous and real-time implementation.The purpose of this paper is to propose an asynchronousalgorithm for the optimal power control of hybrid AC/DCMGs. In addition, we will also introduce the real-time imple-mentation of the algorithm, where the power system is takenas a solver to compute some variables automatically. In thisway, the algorithm can be greatly simplified and the powerloss can be considered. Main contributions of this paper areas follows. • We devise an asynchronous distributed algorithm to solvethe optimal power control problem of hybrid AC/DCMGs. Different from most existing works, in this paper,different kinds of asynchrony are fitted into one unifiedform i.e., the time interval between two consecutiveiterations, by introducing a virtual random clock; • We give a convergence proof of the distributed algorithmby converting it into a fixed-point iteration using theoperator splitting approach. The upper bound of commu-nication delay that guarantees the convergence is given,which is approximately proportional to the square root ofthe number of MGs; • We provide a real-time implementation of the proposedalgorithm, where the power system is taken as a solver tocompute some variables automatically. This simplifies thealgorithm by reducing one oder. Moreover, we providemethods to estimate the unknown load demand in ACand DC MGs with physical system variables respectively.In this way, the impact of power loss such as line andinverter loss can be considered.The rest of this paper is organized as follows. In SectionII, some notations and preliminaries are introduced. SectionIII formulates the power dispatch problem in hybrid MGs andproposes the synchronous algorithm. In Section IV, differenttypes of asynchrony are introduced and an asynchronousalgorithm is proposed. The optimality of its equilibrium pointand convergence of the asynchronous algorithm are provedin Section V. The implementation method in hybrid MGs isintroduced in Section VI. We confirm the performance of thecontroller via simulations on a benchmark low voltage MGsystem in Section VII. Section VIII concludes the paper.II. N
OTATIONS AND P RELIMINARIES
Notations : A hybrid MG system is composed of a clusterof AC and DC MGs connected by lines. Each MG is treatedas a bus with both generation and load. Denote AC MGs by N ac = { , , . . . , n ac } , and DC MGs by N dc = { n ac +1 , n ac + 2 , . . . , n ac + n dc } . Then the set of MG buses is N = N ac ∪ N dc . Let E ⊆ N × N be the set of lines, where ( i, k ) ∈ E if MGs i and k are connected directly. Then theoverall system is modeled as a connected graph G := ( N , E ) .Besides the physical connection among MGs, we also definea communication graph for MGs. Denote by N i the set ofinformational neighbors of MG i over the communication graph, implying MGs i, j can communicate if and only if j ∈ N i . Denote by N i the set of two-hop neighbors of MG i over the communication graph. The cardinality of N i isdenoted by | N i | . The communication graph is also assumedto be undirected and connected, which could be different fromthe physical graph. Denote by L the Laplacian matrix ofcommunication graph. Preliminaries : In this paper, R n ( R n + ) is the n -dimensional(nonnegative) Euclidean space. For a column vector x ∈ R n (matrix A ∈ R m × n ), x T ( A T ) denotes its transpose. Forvectors x, y ∈ R n , x T y = (cid:104) x, y (cid:105) denotes the inner productof x, y . (cid:107) x (cid:107) = √ x T x denotes the Euclidean norm of x .For a positive definite matrix G , denote the inner product (cid:104) x, y (cid:105) G = (cid:104) Gx, y (cid:105) . Similarly, the G -matrix induced norm (cid:107) x (cid:107) G = (cid:112) (cid:104) Gx, x (cid:105) . Use I to denote the identity matrix withproper dimensions. For a matrix A = [ a ij ] , a ij stands for theentry in the i -th row and j -th column of A . Use (cid:81) ni =1 Ω i to denote the Cartesian product of the sets Ω i , i = 1 , · · · , n .Given a collection of y i for i in a certain set Y , y denotes thecolumn vector y := ( y i , i ∈ Y ) with a proper dimension with y i as its components.Define the projection of x onto a set Ω as P Ω ( x ) = arg min y ∈ Ω (cid:107) x − y (cid:107) (1)Use Id to denote the identity operator, i.e., Id( x ) = x , ∀ x .Define N Ω ( x ) = { v | (cid:104) v, y − x (cid:105) ≤ , ∀ y ∈ Ω } . We have P Ω ( x ) = (Id + N Ω ) − ( x ) [26], [27, Chapter 23.1].For a single-valued operator T : Ω ⊂ R n → R n , a point x ∈ Ω is a fixed point of T if T ( x ) ≡ x . The set of fixedpoints of T is denoted by F ix ( T ) . T is nonexpansive if (cid:107)T ( x ) − T ( y ) (cid:107) ≤ (cid:107) x − y (cid:107) , ∀ x, y ∈ Ω . For α ∈ (0 , , T is called α -averaged if there exists a nonexpansive operator R such that T = (1 − α )Id + α R . We use A ( α ) to denotethe class of α -averaged operators. For β ∈ R , T is called β -cocoercive if β T ∈ A ( ) .III. S YNCHRONOUS D ISTRIBUTED A LGORITHM
In this section, we introduce the economic dispatch problemin MGs and propose a synchronous algorithm.
A. Economic dispatch model
The power dispatch is to achieve the power balance in MGswhile minimizing the generation cost, which can be formulatedas the following optimization problem min P gi (cid:88) i ∈N f i ( P gi ) (2a)s.t. (cid:88) i ∈N P gi = (cid:88) i ∈N P di (2b) P gi ≤ P gi ≤ P gi (2c)where f i ( P gi ) = a i ( P gi ) + b i P gi , with a i > , b i > . P gi , P di are the power generation and load demand of MG i respectively. P gi , P gi are the lower and upper bounds of P gi respectively. The objective function (2a) is to minimize thetotal generation cost of the MGs. Constraint (2b) is the powerbalance over MGs. And (2c) is the generation limit of eachMG. For the optimization problem (2), we make the followingassumption.
Assumption The Slater’s condition [28, Chapter 5.2.3] of(2) holds, i.e., problem (2) is feasible due to affine constraints.
B. Synchronous Algorithm
The Lagrangian function of (2) is L = (cid:88) i ∈N f i ( P gi ) + µ (cid:32)(cid:88) i ∈N P gi − (cid:88) i ∈N P di (cid:33) + (cid:88) i ∈N γ − i ( P gi − P gi ) + (cid:88) i ∈N γ + i ( P gi − P gi ) where µ, γ − i , γ + i are Lagrangian multipliers. Here µ is a globalvariable, but will be estimated by individual MGs locally.Define the sets Ω i := (cid:110) P gi | P gi ≤ P gi ≤ P gi (cid:111) , Ω = (cid:89) Ni =1 Ω i (3)Then, we give the synchronous distributed algorithm forpower dispatch (SDPD). In this case, the update of agent i atiteration k is given as, which takes the form of Krasnosel’skiˇi-Mann iteration [27, Chapter 5.2]. ˜ µ i,k = µ i,k + σ µ (cid:18) − (cid:88) j ∈ N i ( µ i,k − µ j,k )+ (cid:88) j ∈ N i ( z i,k − z j,k ) + P gi,k − P di (cid:19) (4a) ˜ z i k ,k = z i k ,k − σ z (cid:18) (cid:88) j ∈ N i (˜ µ i,k − ˜ µ j,k ) − (cid:88) j ∈ N i ( µ i,k − µ j,k ) (cid:19) (4b) ˜ P gi,k = P Ω i (cid:16) P gi,k − σ g (cid:0) f (cid:48) i ( P gi,k ) + 2˜ µ i,k − µ i,k (cid:1)(cid:17) (4c) µ i,k +1 = µ i,k + η k (˜ µ i,k − µ i,k ) (4d) z i,k +1 = z i,k + η k (˜ z i,k − z i,k ) (4e) P gi,k +1 = P gi,k + η k (cid:16) ˜ P gi,k − P gi,k (cid:17) (4f)where σ µ , σ z , σ g , η k are positive constants, and σ µ , σ z , σ g aresupposed to be chosen such that Φ in (10) (given in SectionV.B) is positive definite.In (4), the load demand P di is usually difficult to know.We will provide a practical method to estimate P gi − P di instead of directly measuring P di in the implementation, asexplained in Section VI. Later in Section IV, we will showthat the SDPD is simply a special case of the asynchronousalgorithm. Therefore, its properties, such as the optimalityof the equilibrium point and the convergence, are immediateconsequence of the results of asynchronous algorithm, whichare skipped here.IV. D ISTRIBUTED A SYNCHRONOUS A LGORITHM
In this section, we first introduce several typical types ofasynchrony existing in MGs. Then, we devise an asynchronousalgorithm by modifying Algorithm SDPD.
MG1MG2MG3MG4 MG1MG2MG3MG4 t t t t t t t t t t t t t t t idle idleidleidleidleidle (a) Synchronous distributed computing (b) Asynchronous distributed computing idle Fig. 1. Synchronous versus asynchronous computation
A. Asynchrony in Microgrids
In SDPD, each MG gathers information, computes locallyand conveys new information to its neighbors over the com-munication graph. In this process, asynchrony may arise ineach step. When gathering information, individual MGs mayhave different sampling rates, which results in non-identicalcomputation rates accordingly. In addition, other imperfectcommunication situations such as time delay caused by con-gestion or even failure are very common in power systems,which essentially result in asynchrony.In synchronous computation, an MG has to wait for theslowest neighbor to complete the computation by insertingcertain idle time. Communication delay, congestion or evenpackage loss can further lengthen the waiting time. Thisprocess is illustrated in Fig.1(a). Thus, the slowest MG andcommunication channel may cripple the system in the syn-chronous execution. In contrast, the MGs with asynchronouscomputation do not need to wait and computes continuouslywith little idling, as shown in Fig.1(b). Even if some of itsneighbors fail to update in time, the MG can use the previouslystored information. That means, the MG could execute aniteration without the latest information from its neighbors.
B. Asynchronous Algorithm
In this subsection, we propose an asynchronous distributedalgorithm for power dispatch (ASDPD) based on SDPD.Different from the iteration number k in (4), here each MGhas its own iteration number k i , implying that a local clock is used instead of the global clock. At each iteration k i , MG i computes in the following way. ˜ µ i,k i = µ i,k i − τ kii + σ µ (cid:18) − (cid:88) j ∈ N i (cid:0) µ i,k i − τ kii − µ j,k j − τ kjj (cid:1) + (cid:88) j ∈ N i (cid:0) z i,k i − τ kii − z j,k j − τ kjj (cid:1) + P gi,k i − τ kii − P di (cid:19) (5a) ˜ z i,k i = z i,k i − τ kii − σ z (cid:18) (cid:88) j ∈ N i (cid:0) µ i,k i − τ kii − µ j,k j − τ kjj (cid:1) − (cid:88) j ∈ N i ∪ N i σ µ (cid:96) ij (cid:0) µ i,k i − τ kii − µ j,k i − τ kjj (cid:1) + (cid:88) j ∈ N i ∪ N i σ µ (cid:96) ij (cid:0) z i,k i − τ kii − z j,k j − τ kjj (cid:1) + (cid:88) j ∈ N i σ µ (cid:0) P gi,k i − τ kii − P di (cid:1)(cid:19) (5b) ˜ P gi,k i = P Ω i (cid:18) P gi,k i − τ kii − σ g (cid:0) f (cid:48) i ( P gi,k i − τ kii ) + 2˜ µ i,k i − µ i,k i − τ kii (cid:1)(cid:17) (5c) µ i,k i +1 = µ i,k i − τ kii + η k (cid:0) ˜ µ i,k i − µ i,k i − τ kii (cid:1) (5d) z i,k i +1 = z i,k i − τ kii + η k (cid:0) ˜ z i,k i − z i,k i − τ kii (cid:1) (5e) P gi,k i +1 = P gi,k i − τ kii + η k (cid:0) ˜ P gi,k i − P gi,k i − τ kii (cid:1) (5f)where (cid:96) ij is the i th row and j th column element of matrix L = L × L and τ ki is the random time delay. (cid:96) ij (cid:54) = 0 holdsonly if j ∈ N i ∪ N i [29]. Denote w = ( µ, z, P g ) . w i,k i − τ kii is the state of MG i at iteration k , and w j,k j − τ kjj is thelatest information obtained from MG j . Considering each MGhas its own local clock, we have the following asynchronousalgorithm. Algorithm 1
ASDPD
Input:
For MG i , the input is µ i, , z i, ∈ R n , P gi, ∈ Ω i . Iteration at k i : Suppose MG i ’s clock ticks at time k i , thenMG i is activated and updates its local variables as follows: Step 1: Reading phase
Get µ j,k j − τ kjj , z j,k j − τ kjj from its neighbors’ and two-hopneighbors’ output cache. Step 2: Computing phase
Calculate ˜ µ i,k i , ˜ z i,k i and ˜ P gi,k i according to (5a), (5b) and(5c) respectively.Update µ i,k i +1 , z i,k i +1 and P gi,k i +1 according to (5d), (5e)and (5f) respectively. Step 3: Writing phase
Write µ i,k i +1 , z i,k i +1 to its output cache and µ i,k i +1 , z i,k i +1 , P gi,k i +1 to its local storage. Increase k i to k i + 1 . Remark 1.
In Algorithm 1, if MG i is activated, it willread the latest information from its neighbors. Even if someneighbors are not accessible in time due to communicationissue, it can still execute the iteration by using the previousinformation stored in its input cache. Despite asynchronycaused by different reasons, MG i only concerns whether thelatest information comes, which implies that their effect canbe characterized by the time interval between two successiveiterations. Thus, our algorithm can admit different types ofasynchrony. Remark 2.
As the element (cid:96) ij (cid:54) = 0 holds only if j ∈ N i ∪ N i ,the ASDPD is still distributed. Similar settings are also usedin [29]–[31]. However, it may make the communication graphdenser. In the Section VI, we will show that the powersystem can be treated as a part of solver. Then, we cancarry out the ASDPD by local measurement and neighboringcommunication.V. O PTIMALITY AND C ONVERGENCE A NALYSIS
In this section, we analyze the optimality of the equilibriumpoint of dynamic system (5), as well as the convergence ofAlgorithm 1. To this end, we need to introduce a sequenceof global iteration numbers that serve as a reference globalclock to unify the local iterations of individual MGs in acoherence manner [32]. Note that the global clock is onlyused for convergence analysis, but not required in ASDPD.Specifically, we queue k i of all MGs in the order of time,and use a new number k to denote the k th iteration in thequeue. This treatment is shown in Fig.2 by taking two MGs as an example. Suppose that, at the iteration k , the probabilitythat MG i is activated to update its local variables follows auniform distribution. Hence, each MG is activated with thesame probability, which simplifies the convergence proof. k Sequence of k Sequence of k k -1 k k +1 Local clock 1Local clock 2Global clock k -1 k k +10 2 k -1 k+ k Fig. 2. Local clocks versus global clock
To prove the convergence, we first convert the synchronousalgorithm to a fixed-point iteration with an averaged operator.Then a nonexpansive operator is constructed, leading to theconvergence results of the asynchronous algorithm. Finally,we provide an upper bound of the time delay.
A. Algorithm Reformulation
If the time delay is not considered, (5) is degenerated to (4).In this sense, the SDPD is a special case of ASDPD, and weonly need to analyze the property of ASDPD. The compactform of (5a) - (5f) without delay, i.e., (4a)-(4f), is ˜ µ k = µ k + σ µ (cid:0) − L · µ k + L · z k + P gk − P d (cid:1) (6a) ˜ z k = z k + σ z ( − L · ˜ µ k + L · µ k ) (6b) ˜ P gk = P Ω ( P gk − σ g ( ∇ f ( P gk ) + 2˜ µ k − µ k )) (6c) µ k +1 = µ k + η k (˜ µ k − µ k ) (6d) z k +1 = z k + η k (˜ z k − z k ) (6e) P gk +1 = P gk + η k (cid:16) ˜ P gk − P gk (cid:17) (6f)where, ∇ f ( P gk ) is the gradient of f ( P gk ) . The subscript k i issubstitute by a global notation k . The equation (6b) is obtainedby combing (5a) with (5b), and others are straightforward.Next we show that (6a) − (6f) can be converted into a fixed-point iteration problem with an averaged operator [26], [33].Equation (6a) is equivalent to − L · µ k − P d = − P gk − L · z k + σ − µ (˜ µ k − µ k )= − L ˜ z k − ˜ P gk + σ − µ (˜ µ k − µ k )+ L · (˜ z k − z k ) + ˜ P gk − P gk (7)Similarly, (6b) is equal to L · ˜ µ k + L · (˜ µ k − µ k ) + σ − z (˜ z k − z k ) (8)From the fact that P Ω ( x ) = (Id + N Ω ) − ( x ) , (6b) can berewritten as ˜ P gk = (Id + N Ω ) − ( P gk − σ g ( ∇ f ( P gk ) + 2˜ µ k − µ k )) , or equivalently, −∇ f ( P gk ) = 2˜ µ k − µ k + N Ω ( ˜ P gk ) + σ − g ( ˜ P gk − P gk ) (9)Then, (7) − (9) are rewritten as − Lµ k + P d ∇ f ( P gk ) = − ˜ P gk − L ˜ z k L ˜ µ k ˜ µ k + N Ω ( ˜ P gk ) + Φ ˜ µ k − µ k ˜ z k − z k ˜ P gk − P gk (10)where Φ = σ − µ I L IL σ − z I I σ − g I (11) Define the following two operators B : µzP g (cid:55)→ Lµ + P d ∇ f ( P g ) (12) U : µzP g (cid:55)→ − P g − LzLµµ + N Ω ( P g ) (13)From [26, Lemma 5.6], we know (Id+Φ − U ) − exists andis single-valued. Denote w i = ( µ i , z i , P gi ) , w = ( w i ) , ˜ w =(˜ µ, ˜ z i , ˜ P g ) . We show that (6a) − (6d) can be written as ˜ w k = T ( w k ) (14) w k +1 = w k + η k ( ˜ w k − w k ) (15)where the operator T = (Id + Φ − U ) − (Id − Φ − B ) .For (14), ˜ w k = T ( w k ) is equivalent to −B ( w k ) = U ( ˜ w k ) + Φ · ( ˜ w k − w k ) (16)This is exactly (10). In addition, it is straightforward to seethat (6d) − (6f) are equivalent to (15).Equations (14) − (15) can be further rewritten as w k +1 = w k + η k ( T ( w k ) − w k ) (17)Denote a min = min { a i } , a max = max { a i } , ∀ i ∈ N . Wehave the following result about the operator T . Denote themaximal eigenvalues of L by σ max . We have the followingresult. Lemma 1.
Take ζ = min { σ , a min a max } , κ > ζ , and the stepsizes σ µ , σ z , σ g such that Φ − κI is positive semi-definite. Thenwe have the following assertions under Φ -induced norm.1) T is an averaged operator, and T ∈ A (cid:16) κζ κζ − (cid:17) ;2) there exists a nonexpansive operator R such that T = (cid:18) − κζ κζ − (cid:19) Id + 2 κζ κζ − R
3) operators T and R have the same fixed points, i.e., F ix ( T ) = F ix ( R ) . Proof.
For the assertion 1), we know (Id+Φ − U ) − ∈ A (cid:0) (cid:1) and Id − Φ − B ∈ A (cid:16) κζ (cid:17) from [26, Lemma 5.6]. Then, fol-lowing from [34, Proposition 2.4], we know T ∈ A (cid:16) κζ κζ − (cid:17) .From assertion 1) and definition of averaged operators, thereexists a nonexpansive operator R such that T = (cid:18) − κζ κζ − (cid:19) Id + 2 κζ κζ − R (18)Then, we have the assertion 2).Since T is κζ κζ − -averaged, T is also a nonexpansiveoperator [27, Remark 4.24]. For any nonexpansive operator T , F ix ( T ) (cid:54) = ∅ [27, Theorem 4.19]. Suppose x is a fixedpoint of T , and we have T ( x ) = x = (cid:16) − κζ κζ − (cid:17) Id( x ) + κζ κζ − R ( x ) . Thus, κζ κζ − Id( x ) = κζ κζ − R ( x ) , which is equiv-alent to x = R ( x ) .Similarly, suppose x is a fixed point of R , and we have T ( x ) = (cid:16) − κζ κζ − (cid:17) Id( x ) + κζ κζ − R ( x ) = x . Thus, asser-tion 3) holds, which completes the proof.So far, we convert the synchronous algorithm into a fixed-point iteration problem with an averaged operator (see (17)).Moreover, we also construct a nonexpansive operator R . it enables us to prove the convergence of the asynchronousalgorithm ASDPD, as we explain in the next subsection. B. Optimality of the equilibrium point
Considering dynamic system (5), we give the followingdefinition of its equilibrium point.
Definition 1.
A point w ∗ = ( w ∗ i , i ∈ N ) = ( µ ∗ i , z ∗ i , P g ∗ i ) isan equilibrium point of system (5) if lim k i → + ∞ w k i = w ∗ i holds for all i .Then, we have the following result. Theorem 2.
Suppose Assumption 1 holds. The component P g ∗ , µ ∗ of the equilibrium point w ∗ is the primal-dual optimalsolution to (2). Proof.
By (5a) − (5f) and Definition 1, we have − L · µ ∗ + L · z ∗ + P g ∗ − P d (19a) L · µ ∗ (19b) −∇ f ( P g ∗ ) = N Ω ( P g ∗ ) + µ ∗ (19c)Then, we have (cid:88) i ∈N P g ∗ i − (cid:88) i ∈N P di (20a) µ ∗ i = µ ∗ j = µ ∗ (20b) −∇ f ( P g ∗ ) = N Ω ( P g ∗ ) + µ ∗ (20c)where µ ∗ is a constant. By [35, Theorem 3.25], we know (20)is exactly the KKT condition of the problem (2). In addition,(2) is a convex optimization problem and Slater’s conditionholds, which completes the proof. C. Convergence analysis of asynchronous algorithm
In this subsection, we investigate the convergence of AS-DPD. The basic idea is to treat ASDPD as a randomizedblock-coordinate fixed-point iteration problem with delayedinformation. And then the results in [36] can be applied.Define vectors φ i ∈ R n , i ∈ N . The j th entry of φ i isdenoted by [ φ i ] j . Define [ φ i ] j = 1 if the j th coordinate of w is also a coordinate of w i , and [ φ i ] j = 0 , otherwise. Denote by ϕ a random variable (vector) taking values in φ i , i ∈ N . Then Prob ( ϕ = φ i ) = 1 /n also follows a uniform distribution. Let ϕ k be the value of ϕ at the k th iteration. Then, a randomizedblock-coordinate fixed-point iteration for (15) is given by w k +1 = w k + η k ϕ k ◦ ( T ( w k ) − w k ) (21)where ◦ is the Hadamard product of two matrices. Here, weassume only one MG is activated at each iteration without lossof generality .Since (21) is delay-free, we further modify it for consideringdelayed information, which is w k +1 = w k + η k ϕ k ◦ ( T ( ˆ w k ) − w k ) (22)where ˆ w k is the delayed information at iteration k . Note that,here, k represents the global clock defined in Section V. Wewill show that Algorithm 1 can be written as (22) if ˆ w k is Note that this model helps formulate the algorithm and analyze itsconvergence. In implementation, we allow that two or more MGs are activatedsimultaneously, which can be modeled as two or more iterations in analysis. properly defined. Suppose MG i is activated at the iteration k , then ˆ w k is defined as follows. For MG i and j ∈ N i ,replace µ i,k , z i,k and P gi,k with µ i,k − τ ki , z i,k − τ ki and P gi,k − τ ki .Similarly, replace µ j,k , z j,k with µ j,k − τ kj and z j,k − τ kj . Withthe random variable ϕ , variables of inactivated MGs are keptthe same with the previous iterations. Then we have (22).Next we make the following assumption. Assumption The time interval between two consecutive it-erations is bounded by χ , i.e., sup k max i ∈N { max { τ ki }} ≤ χ .With the assumption, we have the convergence result. Theorem 3.
Suppose Assumptions 1, 2 hold. Take ζ =min { σ , a min a max } , κ > ζ , and the step-sizes σ µ , σ z , σ g suchthat Φ − κI is positive semi-definite. Choose < η k < χ/ √ n κζ − κζ . Then, with ASDPD, P gk and µ k converge tothe primal-dual optimal solution of problem (2) with proba-bility 1. Proof.
Combining (18) and (22), we have w k +1 = w k + η k ϕ k ◦ (cid:18)(cid:18) − κζ κζ − (cid:19) ˆ w k − w k + 2 κζ κζ − R ( ˆ w k ) (cid:19) = w k + η k ϕ k ◦ ( ˆ w k − w k + 2 κζ κζ − R ( ˆ w k ) − ˆ w k ) (cid:19) (23)With w i,k − τ ki = w i,k − τ ki +1 = , · · · , = w i,k , we have ϕ k ◦ ( ˆ w k − w k ) = 0 . Thus, (23) is equivalent to w k +1 = w k + 2 η k κζ κζ − ϕ k ◦ ( R ( ˆ w k ) − ˆ w k ) (24)Invoking [36], (24) with the nonexpansive operator R isessentially a kind of the ARock algorithms suggested in [36].Hence the convergence results given in that paper can directlybe applied. Indeed, Lemma 13 and Theorem 14 of [36] indicatethat, the convergence of ARock is guaranteed by the condition < η k κζ κζ − <
11 + 2 χ/ √ n . (25)Therefore, if η k satisfies < η k < χ/ √ n κζ − κζ , then w k converges to a random variable that takes value in the fixedpoints (denoted by w ∗ k ) of R with probability 1. Recalling F ix ( T ) = F ix ( R ) and Theorem 2, we know P g ∗ k and µ ∗ k , ascomponents of w ∗ k , constitute the primal-dual optimal solutionto the optimization problem (2). This completes the proof.Choose κ = ζ + (cid:15) , where (cid:15) > but very small. Then theupper bound of η k can be estimated by
11 + 2 χ √ n κζ − κζ = 11 + 2 χ √ n ζ(cid:15) ζ(cid:15) ≈
11 + 2 χ √ n Thus, there is η k < . Moreover, the upper bound of η k willdecrease when the time delay increases, i.e., χ increases.Given a fixed η k and a very small (cid:15) > , we have χ < √ n (1 − η k )2 η k (26)Thus, the upper bound of acceptable time delay is approxi-mately proportional to the square root of the number of MGs,which provides a helpful insight for controller design. C i DG gi P ij P di P dci V Fig. 3. Simplified model of a DC MG
VI. R
EAL - TIME I MPLEMENTATION
The rapid variation of renewable generations and loaddemand requires that the controller can be implemented inreal-time. In this section, the implementation of the ASDPD inboth AC and DC MGs is introduced. First, we take the powersystem as a solver, which helps to eliminate the variables ˜ z and z in ASDPD. Then, the real-time control diagram is illustrated.Finally, the optimality of results using such implementationmethod is proved. A. Taking the power system as a solver1) Main idea:
Recalling the synchronous version (4) and(19a), the item (cid:80) j ∈ N i ( z i,k − z j,k ) in (4a) is utilized to balance thedifference between P gi,k and P di . Denote δ ij = z j − z i , and thelast three terms of (4a) are P gi,k − P di − (cid:80) j ∈ N i δ ij,k . From (19a),we know (cid:80) j ∈ N i δ ∗ ij,k + P g ∗ i − P di . This motivates us thepower balance equation of each bus P g ∗ i − P di − (cid:80) j ∈ N i P ∗ ij ,where P ∗ ij is the line power from bus i to bus j in the steadystate. Thus, δ ij has the similar role to the line power P ij .This is a very important observation as P ij is automaticallyimplemented with physical dynamics of the power systems.We only need to measure it, which also implies that thecomputation of ˜ z , z can be avoided. This is what we mean taking the power system as a solver . Moreover, we can take P gi,k − P di + (cid:80) j ∈ N i ( z i,k − z j,k ) as a whole and estimate it inboth AC and DC MGs.
2) AC MGs:
In AC MGs, the swing equation of AC bus i is M i ˙ ω i = P gi − P di − D i ω i − (cid:88) j ∈ N i P ij , i ∈ N ac (27)where, M i > , D i > are constants, and P ij is the linepower from bus i to bus j . This model is suitable for bothsynchronous generators and converters [37]–[39]. (27) can berewritten as P gi − P di − (cid:88) j ∈ N i P ij = M i ˙ ω i + D i ω i , i ∈ N ac (28)Thus, M i ˙ ω i + D i ω i can be used to estimate P gi,k − P di + (cid:80) j ∈ N i ( z i,k − z j,k ) in (4a) and its delayed form in (5a). By thiscontrol structure, the asynchronous algorithm 1 is integratedto the real-time control in AC MGs.
3) DC MGs:
In DC MGs, DC capacitors are used tomaintain the voltage stability of DC buses [40]. Then, themodel of DC MGs can be simplified (see Fig.3). The powerbalance on DC bus i is V dci C i ˙ V dci = P gi − P di − (cid:88) j ∈ N i P ij , i ∈ N dc (29) where, C i is the capacitor connected to the DC bus; V dci isthe voltage of the DC bus. Thus, V dci C i ˙ V dci can be used toestimate the P gi,k − P di + (cid:80) j ∈ N i ( z i,k − z j,k ) in the DC MG. Inthis situation, we only need to measure the voltage, which ismuch easier to implement. Then, the asynchronous algorithmASDPD is integrated to the real-time control in DC MGs.By taking the power system as a solver, the distributedasynchronous algorithm takes the following form ˜ µ i,k i = µ i,k i − τ kii + σ µ (cid:18) − (cid:88) j ∈ N i (cid:0) µ i,k i − τ kii − µ j,k j − τ kjj (cid:1) + M i ˙ ω i + D i ω i (cid:1) , i ∈ N ac (30a) ˜ µ i,k i = µ i,k i − τ kii + σ µ (cid:18) − (cid:88) j ∈ N i (cid:0) µ i,k i − τ kii − µ j,k j − τ kjj (cid:1) + V dci C i ˙ V dci (cid:1) , i ∈ N dc (30b) ˜ P gi,k i = P Ω i (cid:18) P gi,k i − τ kii − σ g (cid:0) f (cid:48) i ( P gi,k i − τ kii ) + 2˜ µ i,k i − µ i,k i − τ kii (cid:1)(cid:17) (30c) µ i,k i +1 = µ i,k i − τ kii + η k (cid:0) ˜ µ i,k i − µ i,k i − τ kii (cid:1) (30d) P gi,k i +1 = P gi,k i − τ kii + η k (cid:0) ˜ P gi,k i − P gi,k i − τ kii (cid:1) (30e)In the algorithm (30), only µ needs to be transmitted betweenneighbors. Moreover, the variables ˜ z , z are not necessary,which simplifies the controller greatly. Based on (30), we havethe following real-time asynchronous distributed algorithm forpower dispatch (RTASDPD) Algorithm 2
RTASDPD
Input:
For MG i , the input is µ i, ∈ R n , P gi, ∈ Ω i . Iteration at k i : Suppose MG i ’s clock ticks at time k i , thenMG i is activated and updates its local variables as follows: Step 1: Reading phase
Get µ j,k j − τ kjj from its neighbors’ output cache. For an ACMG i , measure the frequency ω i . For a DC MG i , measurethe voltage V i . Step 2: Computing phase
For i ∈ N ac , calculate ˜ µ i,k i and ˜ P gi,k i according to (5a)and (5c) respectively. For i ∈ N dc , calculate ˜ µ i,k i and ˜ P gi,k i according to (5b) and (5c) respectively.Update µ i,k i +1 and P gi,k i +1 according to (5d) and (5e)respectively. Step 3: Writing phase
Write µ i,k i +1 to its output cache and µ i,k i +1 , P gi,k i +1 to itslocal storage. Increase k i to k i + 1 . Remark 3.
By taking the power system as a solver, there arethree main advantages: • Only the frequency in AC MG and voltage in DC MGneed to be measured, which avoid the measurement ofload demand P di . As we know, the load demand is usuallydifficult to measure while the measurement of frequencyand voltage is much easier. • We simplify the communication graph, where only theneighboring communication is needed. Moreover, we alsosimplify the controller structure. The auxiliary variables
CurrentControl VoltageControl v oi L f L c i oi PCC busFeeder lineLoad i C f i li Distributed communication network
PowerCalculation ω n - k pi ( P i - P i * ) Q i P i E odi sin( ω i t ) + + E odi ω i ++_ _PWM Electrical NetworkPrimary Power ControlAsynchronous Power Control V n - k qi ( Q i - Q i * ) ,, , i ii k ki i k ikik ik ,, ,, , i j i jii ij j k k k ki ij ji i z ikik ikjk jkjN jN z ,, , , ( ) i i iii i i ' k k ki i i i g ggi ikik ik ik P fP _+ , i i k z , i i k , i gi k P , , , k k ki i ji ii i j ji ik ik jkjN i ref z z , , , , , i i i ik ik ik z , , , , , k k kj j jj j jj j j jk jk jk z , k i k , i i k z , i gi k P Asynchronous Information ,, , i ii k ki i g g gk ikik ik P P P ,, , i i k ki i k i ki k i k z z z CurrentControl VoltageControl v oi L f L c i oi PCC busFeeder lineLoad i C f i li Distributed communication network
PowerCalculation ω n - k pi ( P i - P i * ) Q i P i E odi sin( ω i t ) + + E odi ω i ++_ _PWM Electrical NetworkPrimary Power ControlAsynchronous Power Control V n - k qi ( Q i - Q i * ) ,, , i ii k ki i k i ki k i k ,, , , ( ) i i iii i i ' k k ki i i i g gg i i ki k i k i k P f P _+ , i i k , i gi k P , , , k k ki i ji ii i j ji ii k i i ii k j kj N M D , i i k , kjj j j k , k i k , i gi k P Asynchronous Information ,, , i ii k ki i g g gk i ki k i k P P P Fig. 4. Control diagram of the proposed method ˜ z and z are eliminated, making the controller easier toimplement. • In the problem (2), the power loss is not considered. Inthe real-time implementation, we measure the frequencyand intend to drive it to the nominal value. In this way,the impact of power loss such as line and inverter losscan be considered.
B. Control diagram
The control diagram is shown in Fig.4, which is composedof four levels: the electric network, the primary power control,the asynchronous power control and the distributed communi-cation. In the electric network level, the current and voltageare measured as the input of the primary power control level.The primary power control level includes three loops, i.e.,the current loop, the voltage loop and the power loop. In thepower loop, droop control is utilized for both active power andreactive power control, where the active power and frequencyare sent to the asynchronous power control level. Algorithm1 is integrated in the asynchronous control level, where theasynchronous information from neighboring MGs is utilized.The output P gi,k +1 is the reference of the active power control.The error between P gi,k +1 and the measured active power isfed to the primary power control via an integral operator. Otheroutputs µ i,k +1 is written to its output cache, which are sentto its neighbors via the communication network.The control diagram of DC MGs is similar to that in Fig.4,where the main difference is that the DC bus voltage V dci ismeasured. The details are omitted here. C. Optimality of the implementation
By the implementation method, we claim that the equi-librium point of (30) is also optimal with respect to theoptimization problem (2). In steady state, we have ω i = ω j = GridB1DG1 DG3 DG5DG2 DG4 DG6 GridB1DG1 DG3 DG5DG2 DG4 DG6
DCAC DCAC
MG1MG2 MG3MG4 MG5MG6
Fig. 5. A schematic diagram of a typical 43-bus MG system ω ∗ , ∀ i, j ∈ N , and dV dci dt = 0 , ∀ i ∈ N dc . By (5a) − (5d) andDefinition 1, we have (cid:88) j ∈ N i (cid:0) µ ∗ i − µ ∗ j (cid:1) + D i ω ∗ , i ∈ N ac (31a) (cid:88) j ∈ N i (cid:0) µ ∗ i − µ ∗ j (cid:1) , i ∈ N dc (31b) P g ∗ i = P Ω i (cid:16) P g ∗ i − σ g (cid:0) f (cid:48) i ( P g ∗ i ) + µ ∗ i (cid:1)(cid:17) (31c)From (31a) and (31b), we have r µ ∗ + D ω ∗ = 0 (32a)... r |N ac | µ ∗ + D |N ac | ω ∗ = 0 (32b) r |N ac | +1 µ ∗ = 0 (32c)... r |N | µ ∗ = 0 (32d)where r i is the i th row of Laplacian matrix L , and r + r + · · · + r |N | = 0 . Thus, we have ω ∗ (cid:88) i ∈N ac D i = 0 (32e)This implies that ω ∗ = 0 . Then, we have µ ∗ i = µ ∗ j = µ ∗ witha constant µ ∗ . Other analysis is similar to that of Theorem 2,which is omitted here.VII. C ASE STUDIES
A. System Configuration
To verify the performance of the proposed method, a 44-bus system shown in Fig.5 is used for the test, which is amodified benchmark of low-voltage MG systems [16], [41].The system includes three feeders with six dispatchable MGs,where MG2 and MG5 are DC MGs while the others are ACMGs. The Breaker 1 is open, which implies that the systemoperates in an islanded mode. All simulations are implementedin the professional power system simulation software PSCAD.The simulation scenario is: 1) at t = 2 s, there is a kW load increase in the system; 2) at t = 8 s, there is a kW load drop. Then, each MG increases its generation to M G1 M G2 M G3 M G4 M G5 M G6 Fig. 6. Communication graph of the system F r e qu e n c y ( H z ) DG1DG2DG3DG4DG5DG6
Time (s) , k g i k P ( k W ) Time (s)
DG1DG2DG3DG4DG5DG6 DG1DG2DG3DG4DG5DG6
Time (s) G e n e r a ti on ( k W ) Time (s)
DG1DG2DG3DG4DG5DG6 V o lt a g e o f t h e D C B u s ( V ) MG2MG5
Time (s) V o lt a g e o f t h e D C B u s ( V ) MG2MG5
Time (s) F r e qu e n c y ( H z ) Time (s)
MG1MG2MG3MG4MG5MG6 Time (s) G e n e r a ti on ( k W ) Time (s)
MG1MG2MG3MG4MG5MG6 MG1MG2MG3MG4MG5MG6 mu V o lt a g e o f t h e D C B u s ( V ) MG2MG5
Time (s) F r e qu e n c y ( H z ) Time (s)
MG1MG2MG3MG4MG5MG6 Time (s) G e n e r a ti on ( k W ) Time (s)
MG1MG2MG3MG4MG5MG6MG1MG2MG3MG4MG5MG6
Fig. 7. Dynamics of frequencies (left) and voltages (right). For a DC MG,its frequency implies the frequency of the corresponding DC/AC inverter. balance the power and restore system frequency. Their initialgenerations are (58 . , . , . , . , . , . kW.The communication graph is undirected, which is shown inFig.6. Other parameters are given in Table I. TABLE IS
YSTEM PARAMETERS DG i a i b i P gi (kW) 85 80 90 85 80 80 P gi (kW) 0 0 0 0 0 0 B. Non-identical sampling rates
Individual MGs may have different sampling rates (orcontrol period) in practice, which could cause asynchrony andcompromise the control performance. In this part, we considerthe impact of non-identical sampling rates. The samplingrates of MG1-MG6 are set as 10,000Hz, 12,000Hz, 14,000Hz,16,000Hz, 18,000Hz, 20,000Hz, respectively. The dynamics ofthe frequencies and voltages of MGs are shown in Fig.7.As the load change is located in MG2, the frequency nadirof MG2 is the lowest (about 0.26Hz). The system frequencyrecovers in 4 seconds after the load change. When the loaddecreases, the frequency experiences an overshoot of 0.1Hz,and recovers in 2 seconds. Voltages on the DC buses of MG2and MG5 have a small drop when load increases. On thecontrary, they voltages slightly increase after the load drop.The result demonstrates that the system is fairly stable to loadvariation even with non-identical sampling rates.Dynamics of generations and − µ are given in Fig.8. Atthe end of stage one (from 2s to 8s), generations of MGsare (79.32, 63.60, 90, 81.82, 70.47, 75.08)kW respectively.At the end of stage two (from 8s to 14s), their values are(69.50, 55.46, 79.20, 70.97, 61.86, 65.04)kW respectively.Generations are identical with that obtained by solving thecentralized optimization problem (implemented by CVX). Thisresult verifies the optimality of the proposed method. − µ i stands for the marginal cost of MG i , whose dynamic is givenon the right part of Fig.8. The marginal cost of different MGsconverges to the same value when the system is stabilized,which indicates that the system operates in an optimal state. F r e qu e n c y ( H z ) DG1DG2DG3DG4DG5DG6
Time (s) , k g i k P ( k W ) Time (s)
DG1DG2DG3DG4DG5DG6 DG1DG2DG3DG4DG5DG6
Time (s) G e n e r a ti on ( k W ) Time (s)
DG1DG2DG3DG4DG5DG6 V o lt a g e o f t h e D C B u s ( V ) MG2MG5
Time (s) V o lt a g e o f t h e D C B u s ( V ) MG2MG5
Time (s) F r e qu e n c y ( H z ) Time (s)
MG1MG2MG3MG4MG5MG6 Time (s) G e n e r a ti on ( k W ) Time (s)
MG1MG2MG3MG4MG5MG6 MG1MG2MG3MG4MG5MG6 V o lt a g e o f t h e D C B u s ( V ) MG2MG5
Time (s) F r e qu e n c y ( H z ) Time (s)
MG1MG2MG3MG4MG5MG6 Time (s) G e n e r a ti on ( k W ) Time (s)
MG1MG2MG3MG4MG5MG6MG1MG2MG3MG4MG5MG6
Fig. 8. Dynamics of generations (left) and − µ (right). F r e qu e n c y ( H z ) Time (s) a: delay 200msTime (s) G e n e r a ti on ( k W ) Time (s)b: delay 400msd: delay 800msc: delay 600ms G e n e r a ti on ( k W ) DG1DG2DG3 DG4DG5DG6 G e n e r a ti on ( k W ) F r e qu e n c y ( H z ) Time (s)
Time (s)
20 ms[100, 200]ms[200, 500]ms[500, 800]ms[800, 1000]ms G e n e r a ti on ( k W ) F r e qu e n c y ( H z ) Time (s)
Time (s)
20 ms[100, 200]ms[200, 500]ms[500, 800]ms[800, 1000]ms G e n e r a ti on ( k W ) F r e qu e n c y ( H z ) Time (s)
Time (s)
20 ms[100, 200]ms[200, 500]ms[500, 800]ms[800, 1000]ms 20 ms[100, 200]ms[200, 500]ms[500, 800]ms[800, 1000]ms
Fig. 9. Frequencies and generations under different/varying time delays.
C. Random time delays
In practice, time delay always exists in the communication,which is usually varying up to channel situations. This impliesthat the time delay is random and cannot be known in advance.In this part, we examine the impact of time-varying timedelays. Initially, all the time delays in communication are setas 20ms. And then we intentionally increase the time delayson the channels of MG1-MG2 and MG5-MG6. Additionally,we have the time delays on these two channels varying inranges [100ms, 200ms], [200ms, 500ms], [500ms, 800ms]and [800ms, 1000ms], respectively, while the delays on otherchannels remain 20ms. Frequency and generation dynamicsof MG1 under different scenarios are shown in Fig.9. It isobserved that, When time delays increase, the convergencebecomes slower with larger overshoots in frequency. How-ever, the steady-state generations are still exactly identical tothe optimal solution, which verifies the effectiveness of ourcontroller under varying time delays.
D. Comparison with synchronous algorithm
In this part, we compare the performances of the asyn-chronous and synchronous algorithms under imperfect com-munication. In the asynchronous case, the sampling ratesof MGs are set to the same as that in Section VII.B andthe time delay varies between [500 , ms. The dynamicsof MG1 with two algorithms are shown in Fig.10. Withthe synchronous algorithm SDPD, the system remains stableafter load perturbations. However, the frequency nadir andovershoot deteriorate, and the convergence becomes slower.The generation takes more time to reach the optimal solution,with an obvious fluctuation. The reason is that MGs have towait for the slowest one in the synchronous case. This resultconfirms the advantage of the asynchronous algorithm. E. Plug-n-play test
In this part, we examine the performance of RTASDPDunder the plug-n-play operation mode. The simulation scenario
SDPDASDPD G e n e r a ti on ( k W ) Time (s) - m u SDPDASDPD
Time (s) F r e qu e n c y ( H z ) SDPDASDPD
Time (s) F r e qu e n c y ( H z ) Time (s) G e n e r a ti on ( k W ) Time (s)
SDPD RT ASDPDSDPD RT ASDPD
Fig. 10. Dynamics of frequencies and generations under synchronous andasynchronous cases. F r e qu e n c y ( H z ) G e n e r a ti on ( k W ) -20020406080100120140 DG1DG2DG3DG4DG5DG6 DG1DG2DG3DG4DG5DG6
Time (s)
Time (s) -20020406080100120
DG1DG2DG3DG4DG5DG6 DG1DG2DG3DG4DG5DG6
Time (s)
Time (s) G e n e r a ti on ( k W ) F r e qu e n c y ( H z ) -20020406080100120 MG1MG2MG3MG4MG5MG6
Time (s)
Time (s) G e n e r a ti on ( k W ) F r e qu e n c y ( H z ) MG1MG2MG3MG4MG5MG6
Fig. 11. Dynamics of frequencies (left) and generations (right) when MG4is switched off and on. is that DG4 is switched off at t = 2 s and switched on at t = 8 s.Dynamics of frequencies and generations are illustrated inFig.11. When DG4 is switched off, there is a small frequencyoscillation. Since MG3 is close to DG4, its frequency nadir isthe lowest. When DG4 is connected, the frequency oscillationis more fierce. However, the system is stabilized rapidly in2s. Generation of DG4 drops to zero when it is switched off.Then other DGs increase their generations to re-balance thepower. After DG4 is re-connected, all the generations recoverto the initial values. This demonstrate that our controller canadapt to the plug-n-play mode.VIII. C ONCLUSION
In this paper, we have addressed the information asynchronyissue in the distributed optimal power control of hybrid MGs.By introducing a random clock, different kinds of asynchronydue to imperfect communication are fitted into a unifiedframework. Based on this, we have devised an asynchronousalgorithm with a proof of convergence. We have also providedan upper bound of the time delay. Furthermore, we havepresented the real-time implementation of the asynchronousdistributed power control in hybrid AC and DC MGs. Inthe implementation, the power system is taken as a solver,which simplifies the controller and can consider the powerloss. Numerical experiments on PSCAD confirm the superiorperformance of the proposed methodology.Communication asynchrony widely exists in MGs. Thispaper gives a framework to design distributed controller underimperfect communication. The proposed methodology can alsobe extended to other related problems, such as voltage controlin power systems and energy control in multi-energy systems.R
EFERENCES[1] Q. Xu, J. Xiao, P. Wang, and C. Wen, “A decentralized controlstrategy for economic operation of autonomous ac, dc, and hybrid ac/dcmicrogrids,”
IEEE Trans. Energy Convers. , vol. 32, no. 4, pp. 1345–1355, 2017. [2] X. Xu, Q. Liu, C. Zhang, and Z. Zeng, “Prescribed performancecontroller design for dc converter system with constant power loads indc microgrid,” IEEE Transactions on Systems, Man, and Cybernetics:Systems , no. 99, pp. 1–10, 2018.[3] S. Weng, D. Yue, C. Dou, J. Shi, and C. Huang, “Distributed event-triggered cooperative control for frequency and voltage stability andpower sharing in isolated inverter-based microgrid,”
IEEE transactionson cybernetics , no. 99, pp. 1–13, 2018.[4] Y. Xia, W. Wei, M. Yu, X. Wang, and Y. Peng, “Power management fora hybrid ac/dc microgrid with multiple subgrids,”
IEEE Trans. PowerElectron. , vol. 33, no. 4, pp. 3520–3533, 2018.[5] J. M. Guerrero, J. C. Vasquez, J. Matas, L. G. De Vicu˜na, andM. Castilla, “Hierarchical control of droop-controlled ac and dc mi-crogrids?a general approach toward standardization,”
IEEE Trans. Ind.Electron. , vol. 58, no. 1, pp. 158–172, 2011.[6] Z. Wang, F. Liu, S. H. Low, P. Yang, and S. Mei, “Distributed load-sidecontrol: Coping with variation of renewable generations,” arXiv preprintarXiv:1804.04941 , 2018.[7] F. D¨orfler, J. W. Simpson-Porco, and F. Bullo, “Breaking the hierarchy:Distributed control and economic optimality in microgrids,”
IEEE Trans.Control Network Syst. , vol. 3, no. 3, pp. 241–253, 2016.[8] D. B. Nguyen, J. M. Scherpen, and F. Bliek, “Distributed optimal controlof smart electricity grids with congestion management,”
IEEE Trans.Autom. Sci. Eng. , vol. 14, no. 2, pp. 494–504, 2017.[9] H. Li, G. Chen, T. Huang, Z. Dong, W. Zhu, and L. Gao, “Event-triggered distributed average consensus over directed digital networkswith limited communication bandwidth,”
IEEE Transactions on Cyber-netics , vol. 46, no. 12, pp. 3098–3110, 2016.[10] H. Li, C. Huang, G. Chen, X. Liao, and T. Huang, “Distributedconsensus optimization in multiagent networks with time-varying di-rected topologies and quantized communication,”
IEEE transactions oncybernetics , vol. 47, no. 8, pp. 2044–2057, 2017.[11] H.-S. Ahn, B.-Y. Kim, Y.-H. Lim, B.-H. Lee, and K.-K. Oh, “Distributedcoordination for optimal energy generation and distribution in cyber-physical energy networks,”
IEEE transactions on cybernetics , vol. 48,no. 3, pp. 941–954, 2018.[12] J. W. Simpson-Porco, Q. Shafiee, F. Drfler, J. C. Vasquez, J. M. Guerrero,and F. Bullo, “Secondary frequency and voltage control of islanded mi-crogrids via distributed averaging,”
IEEE Trans. Ind. Electron. , vol. 62,no. 11, pp. 7025–7038, Nov 2015.[13] N. M. Dehkordi, H. R. Baghaee, N. Sadati, and J. M. Guerrero,“Distributed noise-resilient secondary voltage and frequency control forislanded microgrids,”
IEEE Trans. Smart Grid , pp. 1–1, 2018.[14] J. Zhou, H. Zhang, Q. Sun, D. Ma, and B. Huang, “Event-baseddistributed active power sharing control for interconnected ac and dcmicrogrids,”
IEEE Trans. Smart Grid , 2017.[15] G. Chen, F. L. Lewis, E. N. Feng, and Y. Song, “Distributed optimalactive power control of multiple generation systems,”
IEEE Trans. Ind.Electron. , vol. 62, no. 11, pp. 7079–7090, Nov 2015.[16] X. Wu and C. Shen, “Distributed optimal control for stability enhance-ment of microgrids with multiple distributed generators,”
IEEE Trans.Power Syst. , vol. 32, no. 5, pp. 4045–4059, 2017.[17] C. Zhao, X. Duan, and Y. Shi, “Analysis of consensus-based economicdispatch algorithm under time delays,”
IEEE Trans. Syst. Man Cybern.Syst. , no. 99, pp. 1–11, 2018.[18] X. Zhang and A. Papachristodoulou, “A real-time control framework forsmart power networks: Design methodology and stability,”
Automatica ,vol. 58, pp. 43–50, 2015.[19] J. Wu, T. Yang, D. Wu, K. Kalsi, and K. H. Johansson, “Distributedoptimal dispatch of distributed energy resources over lossy communi-cation networks,”
IEEE Transactions on Smart Grid , vol. 8, no. 6, pp.3125–3137, 2017.[20] Z. Wang, F. Liu, Y. Chen, S. H. Low, and S. Mei, “Unified distributedcontrol of stand-alone dc microgrids,”
IEEE Trans. Smart Grid , vol. 10,no. 1, pp. 1013–1024, Jan 2019.[21] M. Zholbaryssov, D. Fooladivanda, and A. D. Dom´ınguez-Garc´ıa, “Re-silient distributed optimal generation dispatch for lossy ac microgrids,”
Systems & Control Letters , vol. 123, pp. 47–54, 2019.[22] D. Meng and K. L. Moore, “Studies on resilient control through mul-tiagent consensus networks subject to disturbances,”
IEEE transactionson cybernetics , vol. 44, no. 11, pp. 2050–2064, 2014.[23] D. Yuan, D. W. Ho, and S. Xu, “Regularized primal–dual subgradientmethod for distributed constrained optimization,”
IEEE transactions oncybernetics , vol. 46, no. 9, pp. 2109–2118, 2016.[24] D. Alkano, J. M. Scherpen, and Y. Chorfi, “Asynchronous distributedcontrol of biogas supply and multienergy demand,”
IEEE Trans. Autom.Sci. Eng. , vol. 14, no. 2, pp. 558–572, 2017. [25] S. Yang, Q. Liu, and J. Wang, “Distributed optimization based on amultiagent system in the presence of communication delays,”
IEEETransactions on Systems, Man, and Cybernetics: Systems , vol. 47, no. 5,pp. 717–728, 2017.[26] P. Yi and L. Pavel, “An operator splitting approach for distributedgeneralized nash equilibria computation,”
Automatica , vol. 102, pp. 111–121, 2019.[27] H. Bauschke and P. L. Combettes,
Convex Analysis and MonotoneOperator Theory in Hilbert Spaces . Springer, 2017.[28] S. Boyd and L. Vandenberghe,
Convex optimization . Cambridgeuniversity press, 2004.[29] T. Anderson, C.-Y. Chang, and S. Martinez, “Distributed approximatenewton algorithms and weight design for constrained optimization,” arXiv preprint arXiv:1804.06238 , 2018.[30] Z. Tang, D. J. Hill, and T. Liu, “Fast distributed reactive power controlfor voltage regulation in distribution networks,”
IEEE Transactions onPower Systems , vol. 34, no. 1, pp. 802–805, 2019.[31] E. Ram´ırez-Llanos and S. Mart´ınez, “Distributed discrete-time optimiza-tion algorithms with applications to resource allocation in epidemicscontrol,”
Optimal Control Applications and Methods , vol. 39, no. 1, pp.160–180, 2018.[32] M. Cao, A. S. Morse, B. D. Anderson et al. , “Agreeing asynchronously,”
IEEE Transactions on Automatic Control , vol. 53, no. 8, p. 1826, 2008.[33] P. Yi and L. Pavel, “Asynchronous distributed algorithms for seekinggeneralized nash equilibria under full and partial decision information,”
IEEE Transactions on Cybernetics , DOI, 10.1109/TCYB.2019.2908091.[34] P. L. Combettes and I. Yamada, “Compositions and convex combinationsof averaged nonexpansive operators,”
J. Math. Anal. Appl. , vol. 425,no. 1, pp. 55 – 70, 2015.[35] A. P. Ruszczy´nski and A. Ruszczynski,
Nonlinear optimization . Prince-ton university press, 2006, vol. 13.[36] Z. Peng, Y. Xu, M. Yan, and W. Yin, “Arock: an algorithmic frameworkfor asynchronous parallel coordinate updates,”
SIAM J. Sci. Comput. ,vol. 38, no. 5, pp. A2851–A2879, 2016.[37] C. De Persis and N. Monshizadeh, “A modular design of incrementallyapunov functions for microgrid control with power sharing,” in
ControlConference (ECC), 2016 European . IEEE, 2016, pp. 1501–1506.[38] Z. Wang, F. Liu, S. H. Low, C. Zhao, and S. Mei, “Distributed frequencycontrol with operational constraints, part ii: Network power balance,”
IEEE Trans. Smart Grid , vol. 10, no. 1, pp. 53–64, Jan 2019.[39] Z. Wang, F. Liu, J. Z. F. Pang, S. H. Low, and S. Mei, “Distributedoptimal frequency control considering a nonlinear network-preservingmodel,”
IEEE Trans. Power Syst. , vol. 34, no. 1, pp. 76–86, Jan 2019.[40] Z. Wang, W. Wu, and B. Zhang, “A distributed control method withminimum generation cost for dc microgrids,”
IEEE Trans. EnergyConvers. , vol. 31, no. 4, pp. 1462–1470, 2016.[41] S. Papathanassiou, N. Hatziargyriou, K. Strunz et al. , “A benchmark lowvoltage microgrid network,” in