Distributed Frequency and Voltage Control for AC Microgrids based on Primal-Dual Gradient Dynamics
Lukas Kölsch, Katharina Wieninger, Stefan Krebs, Sören Hohmann
DDistributed Frequency and Voltage Controlfor AC Microgrids based on Primal-DualGradient Dynamics (cid:63)
Lukas K¨olsch ∗ Katharina Wieninger ∗ Stefan Krebs ∗ S¨oren Hohmann ∗∗ Institute of Control Systems, Karlsruhe Institute of Technology,Karlsruhe, Germany (e-mail: [email protected])
Abstract:
With the gradual transformation of power generation towards renewables, dis-tributed energy resources are becoming more and more relevant for grid stabilization. In order toinvolve all participants in the joint solution of this challenging task, we propose a distributed,model-based and unifying controller for frequency and voltage regulation in AC microgrids,based on steady-state optimal control. It not only unifies frequency and voltage control, but alsoincorporates the classic hierarchy of primary, secondary and tertiary control layers with eachclosed-loop equilibrium being a minimizer of a user-defined cost function. By considering theindividual voltage limits as additional constraints in the corresponding optimization problem,no superordinate specification of voltage setpoints is required. Since the dynamic model of themicrogrid has a port-Hamiltonian structure, stability of the overall system can be assessed usingshifted passivity properties. Furthermore, we demonstrate the effectiveness of the controller andits robustness against fluctuations in active and reactive power demand by means of numericalexamples.
Keywords: distributed control, optimization-based control, electric power systems, microgrids,frequency regulation, voltage regulation1. INTRODUCTIONOur current energy system is undergoing a rapid changetowards an increasing penetration by renewable energysources. Large, central, conventional power plants are be-ing successively replaced by distributed energy resources(DERs). If interconnected in an islanded microgrid with-out any superordinate instance, the individual DERs mustprovide stability by themselves, which has direct implica-tions for the control strategy: On the one hand, controlactions must become increasingly faster due to the lackof rotational inertia. On the other hand, the central su-perordinate authority that ensures stability is increasinglydisappearing which in turn requires the control objectiveto be achieved jointly and decentrally by all participat-ing resources. According to the IEEE-PES Task Forceon Microgrid Stability Definitions, Analysis, and Model-ing, “a microgrid is stable if, after being subjected to adisturbance, all state variables recover to (possibly new)steady-state values which satisfy operational constraints(...), and without the occurrence of involuntary load shed-ding.” (Farrokhabadi et al. (2019)). In terms of frequencyand voltage, this entails achieving a steady state wherefrequencies at each node are equal and correspond to thenominal frequency and where all of the individual voltagemagnitudes remain within certain limits.Such a complex requirement can be met conveniently usingoptimization-based control by characterizing the desired (cid:63)
This work was funded by the Deutsche Forschungsgemeinschaft(DFG, German Research Foundation)—project number 360464149. equilibrium as an optimizer of a constrained optimizationproblem. In the existing literature, promising approacheshave already been developed to enable frequency (Steginket al. (2017b)) and voltage regulation (Magn´usson et al.(2017)) in AC microgrids by this method, see Mohagheghiet al. (2018) and D¨orfler et al. (2019) for an extensivesurvey on current research directions on optimization-based control of future power grids. However, these ex-isting control approaches always pursue only one of thetwo objectives or implicitly assume that specific voltagesetpoints are already provided by a higher level authority.To meet the above requirements in the original and genuinesense, we present a unifying frequency and voltage con-troller using primal-dual gradient dynamics in this paper.It allows real-time optimization of the connected DERs aswell as conventional generators, with the primary goal ofmaintaining the nominal frequency and keeping all voltagemagnitudes within pre-defined limits while minimizing auser-defined cost function. Our controller is unifying inthe sense that the former hierarchical division into pri-mary, secondary and tertiary frequency and voltage con-trol tasks is combined within a single controller and thatfrequency and voltage stabilization as specified above ismaintained simultaneously. In particular, no communica-tion of superordinate set points by a higher-level authorityis necessary. The controller design is based on a nonlinearport-Hamiltonian model for AC microgrids (K¨olsch et al.(2019b)) and allows for a distributed implementation. a r X i v : . [ m a t h . O C ] D ec G L q ‘ p ‘ p ‘ q ‘ q ‘ p ‘ p I U I U f , U G p G Fig. 1. Schematic illustration of synchronous generator( G ), inverter ( I ), and load nodes ( L )The remainder of the paper is structured as follows. Aftersome notational preliminaries, in section 2, we briefly recallthe underlying microgrid model from K¨olsch et al. (2019b)and then formalize the requirements for the desired closed-loop equilibrium. In section 3, we derive a price-basedcontroller for frequency and voltage regulation that meetsthe previously formulated requirements. Moreover, we as-sess shifted passivity of the subsystems and stability ofthe overall system. In section 4, we demonstrate the per-formance of the closed-loop system against input distur-bances by means of a 12-node test network and in section5, we give a brief summary of the main results. Notation
Vector a = col i { a i } = col { a , a , . . . } is acolumn vector of elements a i , i = 1 , , . . . and matrix A = diag i { a i } = diag { a , a , . . . } is a (block-)diagonalmatrix of elements a i , i = 1 , , . . . . The ( n × n )-identitymatrix is denoted by I n and the all-ones vector is denotedby . Positive semi-definite matrices are denoted by (cid:23) (cid:31)
0. Equilibrium variables are marked with an asterisk andshifted values with respect to an equilibrium are markedwith a tilde, i.e. (cid:101) x = x − x (cid:63) . If lower and upper boundsare specified for a particular quantity x , these are markedwith x and x , respectively. For a given µ ≥
0, we define ⟪ x ⟫ + µ := (cid:26) x, µ > ∨ x ≥ , else . (1)If x and µ are vectors of the same size i ∈ N , then (1) canbe applied component-wise, i.e. ⟪ x ⟫ + µ : col i { ⟪ x i ⟫ + µ i } .2. SYSTEM MODEL AND PROBLEMFORMULATION The microgrid is modeled by a directed graph G p = ( V , E p )with V = V G ∪ V I ∪ V L being the set of n G = |V G | synchronous generator nodes , n I = |V I | inverter nodes,and n L = |V L | load nodes. These three different nodetypes are introduced to model the connection to differentmicrogrid participants, see Fig. 1:(1) Synchronous generator nodes are connected to syn-chronous generators of e.g. gas or hydro turbines.(2) Inverter nodes are connected to power electronicsinterfaced DERs such as photovoltaic power stations.(3) Load nodes are connected to consumers with uncon-trollable power demand or alternatively to uncontrol-lable power sources (modeled as negative demands).All nodes may be equipped with a given and uncontrollableactive and reactive power demand p (cid:96) and q (cid:96) which ismodelled as disturbance input. The controlled variables Table 1. List of Microgrid Parameters Symbol Variable A i positive damping coefficient B ii negative of self-susceptance B ij negative of susceptance of line ( i, j ) G ij negative of conductance of line ( i, j ) L i deviation of angular momentum from nominalvalue M i ω n M i moment of inertia p i sending-end active power flow p g,i active power generation p (cid:96),i active power demand q i sending-end reactive power flow q (cid:96),i reactive power demand U i magnitude of transient internal voltage U f,i magnitude of excitation voltage X d,i − X (cid:48) d,i d-axis synchronous minus transient reactance θ i bus voltage phase angle ϑ ij bus voltage angle differenceΦ overall transmission losses τ U,i open-circuit transient time constant of syn-chronous machine ω i deviation of bus frequency from nominal value are active power injection p G and excitation voltage U f at generator nodes and active power injection p I andAC terminal voltage U I at inverter nodes. A list ofall microgrid parameters used in the following can befound in Table 1. The inverter interface is assumed tobe equipped with an internal matching controller whichallows to determine a “virtual” moment of inertia anda “virtual” damping constant for the respective nodesas shown e.g. in Jouini et al. (2016) and Monshizadehet al. (2017). Furthermore we make the following operatingassumptions for the microgrid: Assumption 1. a) The grid is a balanced three-phasedsystem and the lines are represented by its one-phase π -equivalent circuits.b) The grid is operating around the nominal frequency.c) Subtransient dynamics of the synchronous generatorsis neglected.d) The internal matching controller of the invertershas fast dynamics compared to the gradient-basedfrequency and voltage controller.The physical interconnection of the nodes is represented byan incidence matrix D p ∈ R n × m p with n = n G + n I + n L and m p = |E p | . The incidence matrix D p can be subdi-vided as D p = col { D p G , D p I , D p L } , where submatrices D p G , D p I , and D p L correspond to the generator, inverter,and load nodes, respectively. The sending-end active andreactive power flows from node i ∈ V can be calculated bythe AC power flow equations [Machowski et al. (2012)] p i = (cid:88) j ∈N i B ij U i U j sin( ϑ ij ) + G ii U i + (cid:88) j ∈N i G ij U i U j cos( ϑ ij ) , i ∈ V , (2) q i = − (cid:88) j ∈N i B ij U i U j cos( ϑ ij ) + B ii U i + (cid:88) j ∈N i G ij U i U j sin( ϑ ij ) , i ∈ V (3)ith G and B being the conductance and susceptancematrix, respectively, and ϑ ij = θ i − θ j being the voltageangle deviation between two adjacent nodes. N i denotesthe set of neighbors of i , i.e. j ∈ N i if ( i, j ) ∈ E p or( j, i ) ∈ E p .The dynamics of the individual generator, inverter andload nodes are modeled by the following descriptor system,cf. K¨olsch et al. (2019b):˙ θ i = ω i , i ∈ V , (4)˙ L i = − A i ω i + p g,i − p (cid:96),i − p i , i ∈ V G ∪ V I , (5) τ d,i ˙ U i = U f,i − U i − X d,i − X (cid:48) d,i U i · q i , i ∈ V G , (6)0 = − A i ω i − p (cid:96),i − p i , i ∈ V L , (7)0 = − q (cid:96),i − q i , i ∈ V L . (8)Equation (4) defines the relationship between nodal fre-quencies and angle deviations, (5) describes the activepower exchange between the (possibly virtual) mechanicalrotor and the neighboring nodes, (6) describes the tran-sient voltage dynamics of sychronous generator nodes and(7)–(8) describe the active and reactive power conservationat load nodes. To obtain a state space representation, wedefine the plant state vector x p as x p = col { ϑ , L G , L I , U G , ω L , U L } , (9)where for all i ∈ V , j ∈ V G , k ∈ V I , l ∈ V L : ϑ = col i { ϑ i } , L G = col j { L j } , (10) L I = col k { L k } , U G = col j { U j } , (11) ω L = col l { ω l } , U L = col l { U l } . (12)With the Hamiltonian H p ( x p ) = 12 (cid:88) i ∈V G (cid:32) M − i L i + U i X d,i − X (cid:48) d,i (cid:33) + 12 (cid:88) i ∈V I M − i L i − (cid:88) i ∈V G B ii U i − (cid:88) ( i,j ) ∈E B ij U i U j cos( ϑ ij )+ 12 (cid:88) i ∈V L ω L ,i , (13)and the co-state z p = ∇ H ( x p ), this allows to set up aport-Hamiltonian representation of (2)–(8) as follows ˙ ϑ ˙ L G ˙ L I ˙ U G = D (cid:62) p G D (cid:62) p I D (cid:62) p L − D p G − D p I − D p L (cid:124) (cid:123)(cid:122) (cid:125) J p − A G A I R G A L
00 0 0 0 0 (cid:98) U L (cid:124) (cid:123)(cid:122) (cid:125) R p z p − · · · · · · − ϕ G ϕ I (cid:37) G ϕ L (cid:37) L (cid:124) (cid:123)(cid:122) (cid:125) r p + I − (cid:98) I G I − (cid:98) I I τ U − (cid:98) I L − I p G p I U f q (cid:96) p (cid:96) , (14)where A G = diag i { A i } , i ∈ V G , (15) A I = diag i { A i } , i ∈ V I , (16) A L = diag i { A i } , i ∈ V L , (17) R G = diag i { ( X di − X (cid:48) di ) /τ U,i } , i ∈ V G , (18) (cid:98) U L = diag i { U i } , i ∈ V L , (19) ϕ G = col i (cid:8) G ii U i + (cid:88) j ∈N i G ij U i U j cos( ϑ ij ) (cid:9) , i ∈ V G , (20) ϕ I = col i (cid:8) G ii U i + (cid:88) j ∈N i G ij U i U j cos( ϑ ij ) (cid:9) , i ∈ V I , (21) ϕ L = col i (cid:8) G ii U i + (cid:88) j ∈N i G ij U i U j cos( ϑ ij ) (cid:9) , i ∈ V L , (22) (cid:37) G = col i (cid:8) R g,i (cid:88) j ∈N i G ij U i U j sin( ϑ ij ) (cid:9) , i ∈ V G , (23) (cid:37) L = col i (cid:8) (cid:88) j ∈N i G ij U i U j sin( ϑ ij ) (cid:9) , i ∈ V L , (24) ˆ τ U = diag i { /τ U,i } , i ∈ V G , (25) (cid:98) I G = (cid:2) I n g × n g n g × n i n g × n L (cid:3) , (26) (cid:98) I I = (cid:2) n i × n g I n i × n i n i × n L (cid:3) , (27) (cid:98) I L = (cid:2) n L × n g n L × n i I n L × n L (cid:3) . (28)The input vector in (14) is composed of the control input u p = col { p G , p I , U f } and the disturbance input d =col { q (cid:96) , p (cid:96) } . The specification for zero deviation from nominal fre-quency ω n and limitation of voltage magnitudes can beformalized by the following optimization problem:min p G , p I , U f C ( p G , p I ) (OP)subject to Φ = (cid:88) i ∈V G p G ,i + (cid:88) i ∈V I p I ,i − (cid:88) i ∈V p (cid:96),i , (29a) U G ≤ U G ≤ U G , (29b) U I ≤ U I ≤ U I . (29c) C ( p G , p I ) is a user-defined, strictly convex objective func-tion representing e.g. the electricity generation costs andΦ = (cid:62) col { ϕ G , ϕ I , ϕ L } denotes the overall transmissionlosses. The active power balance constraint (29a) is neces-sary for zero frequency deviation, see e.g. Trip et al. (2016).To enable a formulation as a distributed controller in theprocess of the subsequent steps, an exact reformulationof (OP) is derived where the balance constraint (29a) isreplaced by a sparse set of neighbor-to-neighbor balanceconstraints. Moreover, since U f is controllable, box con-straints on U G are replaced by box constraints on U f : Proposition 1.
An exact reformulation of (OP) is given byin p G , p I , U f , ν C ( p G , p I ) (OP (cid:93) )subject to D c ν = (cid:98) I (cid:62)G p G + (cid:98) I (cid:62)I p I − p (cid:96) − ϕ , (30a) Ψ (cid:0) U G (cid:1) ≤ U f ≤ Ψ (cid:0) U G (cid:1) , (30b) U I ≤ U I ≤ U I , (30c)where ϕ = col { ϕ G , ϕ I , ϕ L } , D c is an incidence matrixof a connected communication graph G c = ( V , E c ) and Ψ ( U G ) = col i { Ψ i ( U G ,i ) } , i ∈ V G withΨ i ( U G ,i ) = U G ,i (cid:0) G ii (cid:0) X d,i − X (cid:48) d,i (cid:1)(cid:1) + (cid:88) j ∈N i U j ( G ij sin ( ϑ ij ) − B ij cos ( ϑ ij )) . (31) Proof.
Let x (cid:5) p be an optimizer of (OP) and x (cid:93)p be anoptimizer of (OP (cid:93) ).To prove that x (cid:5) p fulfills (30a), we first recall that for eachequilibrium of (14) it holds that ω = (cf. Proposition 1in K¨olsch et al. (2019b)). Accordingly, (5) and (7) can bewritten as = (cid:98) I (cid:62)G p G + (cid:98) I (cid:62)I p I − p (cid:96) − p (32)with p = col i { p i } , i ∈ V . Inserting (2) in (32) yields = (cid:98) I (cid:62)G p G + (cid:98) I (cid:62)I p I − p (cid:96) − ϕ − φ (33)with φ = col i { φ i } , φ i = (cid:80) j ∈N i G ij U i U j cos( ϑ i j ). Equa-tion (33) is equivalent to (30a) if and only if there exists a ν ∈ R m c with D c ν = φ . Since D c is the incidence matrixof a connected graph, rank( D c ) = n − D red c ν = φ red , (34)again with rank( D red c ) = n −
1. Note that (34) is a systemof n − m c variables. Since G c isassumed to be connected, m c must be greater than or equalto n −
1. This implies that (34) is underdetermined andhence there always exists a ν satisfying (30a).To proove that x (cid:93)p fulfills (29a), left-multiply (30a) with (cid:62) which yieldsΦ (cid:93) = (cid:88) i ∈V G p (cid:93) G ,i + (cid:88) i ∈V I p (cid:93) I ,i − (cid:88) i ∈V p (cid:96),i . (35)This completes the proof of the first equivalence.To proof equivalence of (29b) and (30b), let x (cid:63)p be anequilibrium of (14). From the fourth row of (14) it followsthat for each i ∈ V G U (cid:63)f,i − U (cid:63) G ,i − (cid:0) X d,i − X (cid:48) d,i (cid:1) q (cid:63)i /U (cid:63) G ,i . (36)Inserting (3) in (36) and comparing with (31) yields U (cid:63)f,i =Ψ i ( U (cid:63) G ,i ). Since G ii ≥ X d,i − X (cid:48) d,i >
0, it followsthat ∇ Ψ i ( U G ,i ) = (1 + G ii ( X d,i − X (cid:48) d,i )) >
0, i.e. Ψ i ( U G ,i )is a strictly increasing affine map, thus U (cid:63) G ≤ U G ⇐⇒ Ψ ( U (cid:63) G ) ≤ Ψ (cid:0) U G (cid:1) and U G ≤ U (cid:63) G ⇐⇒ Ψ ( U G ) ≤ Ψ (cid:0) U (cid:63) G (cid:1) .This completes the proof of the second equivalence.To sum up, each x (cid:5) p is feasible for (OP (cid:93) ) and each x (cid:93)p isfeasible for (OP), thus (OP) and (OP (cid:93) ) are equivalent. (cid:4)
3. CONTROLLER DESIGNNow, for optimization problem (OP (cid:93) ), a primal-dual gradi-ent controller [Jokic et al. (2009); Stegink et al. (2017b,a)] can be applied so that together with (14), a closed-loopequilibrium is achieved which is the solution of (OP (cid:93) ). Toshorten the notation, denote the vector of active powergenerations by p g , i.e. p g = col { p G , p I } . Suppose that some constraint qualification[Boyd and Vandenberghe (2015)] holds for (OP (cid:93) ). Theneach closed-loop equilibrium of (14) together with thedistributed primal-dual gradient controller τ g ˙ p g = −∇ C ( p g ) + (cid:98) I g λ + u c , (37a) τ λ ˙ λ = D c ν − (cid:98) I (cid:62) g p g + p (cid:96) + ϕ , (37b) τ ν ˙ ν = − D (cid:62) c λ , (37c) τ µ G− ˙ µ G− = ⟪ Ψ G − U f ⟫ + µ G− , (38a) τ µ G + ˙ µ G + = ⟪ U f − Ψ G ⟫ + µ G + , (38b) τ U G ˙ U f = µ G− − µ G + , (38c) τ µ I− ˙ µ I− = ⟪ U I − U I ⟫ + µ I− , (38d) τ µ I + ˙ µ I + = ⟪ U I − U I ⟫ + µ I + , (38e) τ U I ˙ U I = − ( ∇ Ψ ( U I )) (cid:62) ( µ G− − µ G + ) − ( ∇ ϕ ( U I )) (cid:62) λ + µ I− − µ I + (38f)with u c = − ω , τ > Ψ G = Ψ (cid:0) U G (cid:1) , Ψ G = Ψ (cid:0) U G (cid:1) isan optimizer of (OP (cid:93) ). Proof.
The Lagrangian of (OP (cid:93) ) is L ( p g , U G , U I , ν , λ , µ G− , µ G + , µ I− , µ I + ) = C ( p G , p I ) + λ (cid:62) ( D c ν − (cid:98) I (cid:62) g p g + p (cid:96) + ϕ )+ µ (cid:62)G− ( Ψ G − U f ) + µ (cid:62)G + ( U f − Ψ G )+ µ (cid:62)I− ( U I − U I ) + µ (cid:62)I + ( U I − U I ) . (39)The Karush-Kuhn-Tucker (KKT) conditions specifying asaddle point of L can be applied to derive a necessarycondition for an optimizer of (OP (cid:93) ): = ∇ C ( p (cid:63)g ) − λ (cid:63) , (40) = D c ν (cid:63) − (cid:98) I (cid:62) g p (cid:63)g + p (cid:96) + ϕ (cid:63) , (41) = D (cid:62) c λ (cid:63) , (42) = − µ (cid:63) G− + µ (cid:63) G + , (43) = µ (cid:63) (cid:62)G− ( Ψ (cid:63) G − U (cid:63)f ) , (44) = µ (cid:63) (cid:62)G + ( U (cid:63)f − Ψ (cid:63) G ) , (45) ≤ µ (cid:63) G− , µ (cid:63) G + , (46) = ( ∇ Ψ ( U (cid:63) I )) (cid:62) ( µ (cid:63) G− − µ (cid:63) G + )+ ( ∇ ϕ ( U (cid:63) I )) (cid:62) λ (cid:63) − µ (cid:63) I− + µ (cid:63) I + , (47) = µ (cid:63) (cid:62)I− ( U (cid:63) I − U (cid:63) I ) , (48) = µ (cid:63) (cid:62)I + ( U (cid:63) I − U (cid:63) I ) , (49) ≤ µ (cid:63) I− , µ (cid:63) I + . (50)Consequently, applying the gradient method (Arrow et al.(1958)) provides (37a)–(38f). (cid:4) From (37a)–(38f) it follows that the controller is dis-tributed, since each local controller at node i only dependson local measured values as well as values of neighboringnodes j ∈ N i . Keeping this in mind, we now analyze con-vergence of the closed loop-system towards an equilibrium. .2 Analysis of Shifted Passivity and Stability Let d (cid:63) = col { q (cid:63)(cid:96) , p (cid:63)(cid:96) } be a given, constant disturbanceinput vector. We assume in the following that thereexists an equilibrium x (cid:63) satisfying the following regularitycondition: Assumption 2.
The Hessian of H p ( x p ) is positive definiteat steady state x (cid:63)p .This assumption can be ensured by satisfying e.g. the(relatively mild) operational condition presented in Propo-sition 9 of Stegink et al. (2017a) and originally derived inProposition 1 of De Persis and Monshizadeh (2016).To investigate the shifted passivity of the closed-loopsystem with respect to x (cid:63) , it is convenient to analyze plantsystem (14), frequency (37a)–(37c), and voltage controller(38a)–(38f) separately. Proposition 3.
The plant system (14) with output y p = G (cid:62) p z p provides a shifted passivity property with respectto x (cid:63)p if (cid:2) z p − z (cid:63)p (cid:3) (cid:62) (cid:2) R ( x p ) − R ( x (cid:63)p ) (cid:3) ≥ R ( x p ) = R p z p + r p holds. Proof.
The plant system provides a shifted passivityproperty with respect to x (cid:63)p if the shifted plant Hamil-tonian (cid:101) H p ( (cid:101) x p ) := H p ( x p ) − ( (cid:101) x p ) (cid:62) ∇ H p ( x (cid:63)p ) − H p ( x (cid:63)p ) (52)with (cid:101) x p = x p − x (cid:63)p satisfies (cid:101) H p ( (cid:101) x p ) (cid:31) (cid:101) H p ( (cid:101) x p ) ≤ (cid:101) y (cid:62) p (cid:101) u p (53)with y p = G (cid:62) p z p . The positive definiteness condition issatisfied locally due to Assumption 2. To investigate (53),with the same reasoning as in K¨olsch et al. (2019b), forconstant disturbance input d (cid:63) , (14) can be expressed inthe form E ˙ (cid:101) x p = J p ∇ (cid:101) H p ( x ) − (cid:2) R ( x p ) − R ( x (cid:63)p ) (cid:3) + G p (cid:101) u p . (54)With ˙ (cid:101) H p ( (cid:101) x p ) = (cid:101) z (cid:62) p E ˙ (cid:101) x p and bearing in mind that J p isskew-symmetric, this results in condition (51). (cid:4) Note that for lossless grids r p = , and hence (51) is alwaysfulfilled due to the fact that R p (cid:23) Proposition 4.
The frequency controller with input u c andoutput y u = p g provides a shifted passivity property if( p g − p (cid:63)g ) (cid:62) ( ∇ C ( p g ) − ∇ C ( p (cid:63)g )) ≥ ( λ − λ (cid:63) ) (cid:62) ( ϕ − ϕ (cid:63) ) . (55) Proof.
By defining the frequency controller state ξ =col { τ g p g , τ λ λ , τ ν ν } and the frequency controller Hamilto-nian H ( ξ ) = x (cid:62) c τ − c x c with τ c = diag { τ g , τ λ , τ ν } (cid:31) ξ = J ∇ H ( ξ ) − r + G u u c + G d p (cid:96) (56) y u = G (cid:62) u ζ (57) y d = G (cid:62) d ζ (58)where ζ = ∇ H ( ξ ) = col { p g , λ , ν } and J = (cid:98) I g − (cid:98) I (cid:62) g D c − D (cid:62) c , r = (cid:34) ∇ C ( p g ) − ϕ (cid:35) , (59) G u = (cid:34) I (cid:35) , G d = (cid:34) I (cid:35) . (60)With the shifted controller Hamiltonian (cid:101) H ( (cid:101) ξ ) := H ( ξ ) − ( (cid:101) ξ ) (cid:62) ∇ H ( ξ (cid:63) ) − H ( ξ (cid:63) ) , (61)the shifted controller co-state (cid:101) ζ equals (cid:101) ζ = ∇ (cid:101) H ( (cid:101) ξ ) = ∇ H ( ξ ) − ∇ H ( ξ (cid:63) ) = ζ − ζ (cid:63) . (62)Since the disturbance input p (cid:96) is assumed to be constant,(56) can be expressed in shifted coordinates as follows˙ (cid:101) ξ = J ∇ (cid:101) H ( (cid:101) ξ ) − [ r ( ξ ) − r ( ξ (cid:63) )] + G u (cid:101) u c , (63) (cid:101) y u = G (cid:62) u (cid:101) ζ . (64)Due to strict convexity of H ( ξ ), positive definiteness (cid:101) H ( (cid:101) ξ ) (cid:31) (cid:101) H p ( (cid:101) x p ) = ( ∇ (cid:101) H ( (cid:101) ξ )) (cid:62) ˙ (cid:101) ξ = (cid:101) ζ (cid:62) [ r ( ξ ) − r ( ξ (cid:63) )] + (cid:101) ζ (cid:62) G u (cid:101) u c (65)fulfills ˙ (cid:101) H p ( (cid:101) x p ) ≤ (cid:101) y (cid:62) u (cid:101) u c . Bearing in mind (57), this leadsto (cid:101) ζ (cid:62) [ r ( ξ ) − r ( ξ (cid:63) )] ≤
0, which is equivalent to condition(55). (cid:4)
Note that for lossless grids ϕ = , and hence (55) is alwaysfulfilled due to the fact that C ( p g ) is (strictly) convex.Since the interconnection between frequency controller andplant system is power-preserving, u c = − y p = − G (cid:62) p ∇ H p ( x p ) = − ω g , (66) u p = y c = G (cid:62) u ∇ H ( ξ ) = p g , (67)shifted passivity of the subsystems in terms of Proposi-tions 3 and 4 implies shifted passivity of the closed-loopsystem (14),(37a)–(37c). In fact, the conditions stated inPropositions 3 and 4 are not necessary, since as an excessof passivity in one subsystem can compensate for the lackof passivity in the other subsystem, cf. van der Schaft(2017). Based on the previous conditions we now formulatea stability criterion for the overall system with frequencyand voltage controller: Proposition 5.
Assume that the conditions of Propositions3 and 4 hold. For a constant input d (cid:63) , let ( x (cid:63)p , ξ (cid:63) , ξ (cid:63) )denote an equilibrium of (14),(37a)–(38f) and let( ∇ (cid:101) H p ( (cid:101) U I )) (cid:62) ˙ (cid:101) U I < (cid:101) U I (( ∇ Ψ ( U I )) (cid:62) ( µ G− − µ G + ) − (cid:101) U I (( ∇ Ψ ( U (cid:63) I )) (cid:62) ( µ (cid:63) G− − µ (cid:63) G + ) − (cid:101) U I (( ∇ ϕ ( U I )) (cid:62) λ − ( ∇ ϕ ( U (cid:63) I )) (cid:62) λ (cid:63) )hold. Then there exists a neighborhood B around( x (cid:63)p , ξ (cid:63) , ξ (cid:63) ) such that if ( x p , ξ , ξ ) ∈ B , then the stateasymptotically converges to ( x (cid:63)p , ξ (cid:63) , ξ (cid:63) ). Proof.
With the voltage controller state ξ =col { τ µ G− µ G− , τ µ G + µ G + , τ U G U f , τ µ I− µ I− , τ µ I + µ I + , τ U I U I } , let ξ (cid:63) denote an equilibrium of (38a)–(38f) and define the Lyapunov function candidate (cid:101) V ( (cid:101) ξ ) = (cid:101) ξ (cid:62) τ − (cid:101) ξ (cid:62) where τ = diag { τ µ G− , τ µ G + , τ U G , µ I− , τ µ I + , τ U I } (cid:31) (cid:101) ξ (cid:62) = ξ − ξ (cid:63) . This allows to setup an overall Lyapunov function candidate (cid:101) V ( (cid:101) x p , (cid:101) ξ , (cid:101) ξ )for the overall closed-loop system as the sum of (cid:101) H p , (cid:101) H ,and (cid:101) V : (cid:101) V ( (cid:101) x p , (cid:101) ξ , (cid:101) ξ ) = (cid:101) H p ( (cid:101) x p , (cid:101) ξ ) + (cid:101) H ( (cid:101) ξ ) + (cid:101) V ( (cid:101) ξ ) (68)In the notation of (68) it has been taken into account that U I is contained in the plant Hamiltonian H p , see (13). Asalready shown, all summands in (68) are positive definite,thus (cid:101) V is also positive definite.Respect to ˙ (cid:101) V , we observe that˙ (cid:101) V ( (cid:101) x p , (cid:101) ξ , (cid:101) ξ ) = ˙ (cid:101) H p ( (cid:101) x p , (cid:101) ξ ) + ˙ (cid:101) H ( (cid:101) ξ ) + ˙ (cid:101) V ( (cid:101) ξ ) (69)where the first summand of (69) equals˙ (cid:101) H p ( (cid:101) x p , (cid:101) ξ ) = ( ∇ (cid:101) H p ( (cid:101) x p )) (cid:62) J p ( ∇ (cid:101) H p ( (cid:101) x p )) − ( ∇ (cid:101) H p ( (cid:101) x p )) (cid:62) (cid:2) R ( x p ) − R ( x (cid:63)p ) (cid:3) + ( ∇ (cid:101) H p ( (cid:101) ξ )) (cid:62) ˙ (cid:101) ξ (70)= − ( z p − z (cid:63)p ) (cid:62) (cid:2) R ( x p ) − R ( x (cid:63)p ) (cid:3) + ( ∇ (cid:101) H p ( (cid:101) ξ )) (cid:62) ˙ (cid:101) ξ (71)due to skew-symmetry of J p . The second summand of (69)is ˙ (cid:101) H ( (cid:101) ξ ) = − (cid:101) p (cid:62) g ( ∇ C ( p g ) −∇ C ( p (cid:63)g ))+( λ − λ (cid:63) ) (cid:62) ( ϕ − ϕ (cid:63) ) . The third summand of (69) equals˙ (cid:101) V ( (cid:101) ξ ) = (cid:101) µ (cid:62)G− ⟪ Ψ G − U f ⟫ + µ G− + (cid:101) µ (cid:62)G + ⟪ U f − Ψ G ⟫ + µ G + + (cid:101) U (cid:62) f (cid:101) µ G− − (cid:101) U (cid:62) f (cid:101) µ G + + (cid:101) µ (cid:62)I− ⟪ U I − U I ⟫ + µ I− + (cid:101) µ (cid:62)I + ⟪ U I − U I ⟫ + µ I + + (cid:101) U (cid:62)I (cid:101) µ G− − (cid:101) U (cid:62)I (cid:101) µ G + + (cid:101) U (cid:62)I ( ∇ Ψ ( U I )) (cid:62) ( µ G− − µ G + ) − (cid:101) U (cid:62)I ( ∇ Ψ ( U (cid:63) I )) (cid:62) ( µ (cid:63) G− − µ (cid:63) G + ) − (cid:101) U (cid:62)I (( ∇ ϕ ( U I )) (cid:62) λ − ( ∇ ϕ ( U (cid:63) I )) (cid:62) λ (cid:63) ) . (72)In Proposition 3 of Stegink et al. (2015), it was shownthat (cid:101) µ (cid:62) ⟪ g ⟫ + µ ≤ (cid:101) µ (cid:62) g and (cid:101) µ (cid:62) g (cid:63) ≤ holds for each convexfunction g = col i { g i } . Hence (cid:101) µ (cid:62)G− ⟪ Ψ G − U f ⟫ + µ G− ≤ (cid:101) µ (cid:62)G− ( Ψ G − U f )= (cid:101) µ (cid:62)G− ( Ψ G − U (cid:63)f − (cid:101) U f ) ≤ − (cid:101) µ (cid:62)G− (cid:101) U f . (73)With the same procedure it can be calculated for G +, I− ,and I + that the first four rows of (72) are less than orequal to zero, thus˙ (cid:101) V ( (cid:101) ξ ) ≤ (cid:101) U (cid:62)I ( ∇ Ψ ( U I )) (cid:62) ( µ G− − µ G + ) − (cid:101) U (cid:62)I ( ∇ Ψ ( U (cid:63) I )) (cid:62) ( µ (cid:63) G− − µ (cid:63) G + ) − (cid:101) U (cid:62)I ( ∇ ϕ ( U I )) (cid:62) λ − ( ∇ ϕ ( U (cid:63) I )) (cid:62) λ (cid:63) ) . (74)With the assumption from Proposition 3, the first row of(71) is ≤ (cid:101) H ( (cid:101) ξ ) ≤
0. Hence ˙ (cid:101) V ( (cid:101) x p , (cid:101) ξ , (cid:101) ξ ) ≤ ∇ (cid:101) H p ( (cid:101) ξ )) ˙ (cid:101) ξ < (cid:101) U (cid:62)I ( ∇ Ψ ( U I )) (cid:62) ( µ G− − µ G + ) − (cid:101) U (cid:62)I ( ∇ Ψ ( U (cid:63) I )) (cid:62) ( µ (cid:63) G− − µ (cid:63) G + ) − (cid:101) U (cid:62)I ( ∇ ϕ ( U I )) (cid:62) λ − ( ∇ ϕ ( U (cid:63) I )) (cid:62) λ (cid:63) )holds. With ( ∇ (cid:101) H p ( (cid:101) ξ )) (cid:62) ˙ (cid:101) ξ = ( ∇ (cid:101) H p ( (cid:101) U I )) (cid:62) ˙ (cid:101) U I , the proofis complete. (cid:4)
124 3857 611 1012 9
Fig. 2. Network topology of exemplary microgrid.
The notation as an optimization problem provides numer-ous possibilities for extending the proposed controller: Forexample, hard limitations of the active power generation p g can be incorporated as additional box constraints.Moreover, a simplified controller with a more practicalstability criterion can be formulated by assuming losslesstransmission lines, i.e. ϕ = , and by exploiting the factthat lower and upper bounds for U f are relatively insen-sitive with respect to U I . Simplified’ projected bounds Ψ and Ψ can thus be formulated by replacing all invertervoltage magnitudes U I in (31) by an estimator (cid:98) U I , e.g.( U I + U I ) /
2, to obtain estimators (cid:98) Ψ and (cid:98) Ψ of the lowerand upper bounds which are independent of U I , leading to ∇ (cid:98) Ψ ( U I ) = . This gives the following simplified stabilitycriterion: Corollary 6.
Assume that the conditions of Proposition 3and 4 hold. For a constant input d (cid:63) , let ( x (cid:63)p , ξ (cid:63) , ξ (cid:63) ) denotean equilibrium of (14),(37a)–(38f) with zero transmissionline losses, i.e. ϕ = , and Ψ = (cid:98) Ψ in (31). If (cid:88) i ∈V I ( µ I− ,i − µ I + ,i ) · (cid:32) (cid:88) j ∈N i B ij U j cos( ϑ ij ) − (cid:88) j ∈N i B ij U (cid:63)j cos( ϑ (cid:63)ij ) (cid:33) < B around( x (cid:63)p , ξ (cid:63) , ξ (cid:63) ) such that if ( x p , ξ , ξ ) ∈ B , then the stateasymptotically converges to ( x (cid:63)p , ξ (cid:63) , ξ (cid:63) ).Another possible extension is to introduce an additionalnode type to model controllable active and reactive powerinput or consumption and to include reactive power limi-tations in the optimization problem.4. SIMULATIONWe now verify the presented control approach by simu-lating an exemplary microgrid as presented in Fig. 2 with n G = n I = n L = 4. The parameter values for all nodes, listed in Tables 2–4,are based on Trip et al. (2016) and K¨olsch et al. (2019a),whereby the “virtual” moments of inertia of inverter nodeswere selected to be considerably smaller than the corre-sponding moments of inertia of generator nodes. All valuesare given in p.u. ( U base = 20 kV , S base = 100 MVAr) exceptable 2. Parameters of Generator Nodes i A i B ii -6.0567 -8.014 -6.6755 -4.144 G ii M i X d,i X (cid:48) d,i τ U,i
Table 3. Parameters of Inverter Nodes i A i B ii -5.6611 -8.1791 -3.6067 -6.1635 G ii M i Table 4. Parameters of Load Nodes i A i B ii -1.716 -2.41244 -1.692 -1.848 G ii Table 5. Parameters of Transmission Lines ( i, j ) B ij G ij (1 ,
2) 1.905 -1.9167(1 ,
4) 1.976 -2.04(1 ,
8) 2.176 -1.19178(2 ,
3) 2.352 -2.3256(2 ,
5) 1.966 -1.66(2 ,
6) 1.8012 -1.8444(3 ,
6) 2.052 -1.7396 ( i, j ) B ij G ij (3 ,
8) 2.2716 -1.7394(4 ,
5) 2.168 -2.016(5 ,
11) 1.692 -1.7984(6 ,
7) 1.7588 -1.7588(6 ,
10) 2.41244 -1.9544(7 ,
12) 1.848 -1.988(8 ,
9) 1.716 -2.0132 τ U,i , which is given in seconds. The transmission line pa-rameters are generated choosing both line resistances R ij and reactances X ij as evenly distributed random variablesaround 1 Ω ±
10% and can be found in Table 5. In thesame fashion, the shunt capacitors at each node are setto 10 nF ± τ µ G− , τ µ G + and τ U G are set to 0 . τ U I to 10 and all other controllerparameters to 0 .
1. The cost function is chosen to C ( p g ) = 12 (cid:88) i ∈V G ∪ V I w i · p g,i (75)with weighting factors ω = 1 , w = 1 . , w = 1 . ∇ C ( p (cid:63)i ) = ∇ C ( p (cid:63)j ) at steady state, thisspecific choice of C ( p g ) as a weighted sum of squaresleads to active power sharing , i.e. a proportional share p (cid:63)g,i /w i = p (cid:63)g,j /w j = const. for all i, j ∈ V G ∪ V I . Thevoltage limits are set to [0 .
98 1 .
02] and the upper boundfor active power generation is set to 0 .
6. We choose D c to be identical to the plant incidence matrix D p after ithas been pointed out in K¨olsch et al. (2019a) that thespecific choice of D c has little influence on the convergencespeed to the desired equilibrium. The initial values ofdisturbance input vector d = col { p (cid:96) , q (cid:96) } and state vector x = col { x p , ξ , ξ } are chosen such that the closed-loopsystem starts in synchronous mode with ω ( t = 0) = andsuch that a number of voltage magnitudes are already closeto or at their limits of 0.98 or 1.02. At regular intervals of200 s, a step of +0 . . t in s p ‘ , i , q ‘ , i i np . u . q ‘, p ‘, q ‘, p ‘, Fig. 3. Stepwise increase of active and reactive powerdemands at load nodes. . . . . . t in s f i i n H z i =1 i =2 i =3 i =4 i =5 i =6 i =7 i =8 i =9 i =10 i =11 i =12 . . . . t in s U g i np . u . i =1 i =2 i =3 i =4 i =5 i =6 i =7 i =8 . . . . t in s p g , i i np . u . i =1 i =2 i =3 i =4 i =5 i =6 i =7 i =8 Fig. 4. Evolution of frequencies (a), voltage magnitudes(b), and active power generation (c) after step in-creaseis assigned sequentially at load nodes, as shown in Fig. 3.The simulations are carried out in Wolfram Mathematica12.0.
Fig. 4(a) and 4(b) show the node frequencies and thevoltage magnitudes, respectively. It can be seen that insteady state all frequencies are synchronized to 50 Hzand the voltage limits of 0.98 and 1.02 are kept. Theovershoots are around 0 . .
005 p.u. for the voltage amplitudes. Fig. 4(c) shows thecorresponding active power generations p g at generatorand inverter nodes. Remarkably, the individual powerinjections p g,i are equidistant from each other at steady
100 200 300 400 500 600 700 800 9000 . . . . t in s U g i np . u . i =1 i =2 i =3 i =4 i =5 i =6 i =7 i =8 Fig. 5. Evolution of voltage magnitudes if Ψ = (cid:98) Ψ .state, regardless of the total generation, thus active powersharing is given. In addition, it can be stated after thejump at 650 s that the specified maximum limit of 0.6 p.u.for active power generation is not exceeded. Although thislimit is reached for nodes 7 and 8, the remaining nodesperform active power sharing without being affected.Fig. 5 shows the voltage magnitudes for the same scenarioas above in case the simplified upper and lower limits (cid:98) Ψ and (cid:98) Ψ are used. Again, overshoots of about 0 .
005 p.u. canbe detected. The progression of voltage magnitudes overtime is slightly different, but comparable to Fig. 4(b). Inparticular, the voltage limits are not exceeded at steadystate. 5. CONCLUSIONIn this paper, we presented a model-based frequency andvoltage controller for AC microgrids ensuring optimalitywith regard to a user-defined cost function. For this pur-pose, the controller continuously tracks the KKT pointof a constrained optimization problem using the gradientmethod. The underlying nonlinear microgrid model con-sists of a mixture of conventional synchronous generators,power electronics interfaced sources and uncontrollableloads. The user-defined cost function can be specified insuch a way that e.g. active power sharing is achieved.Moreover, the resulting controller equations exhibit a dis-tributed neighbor-to-neighbor communication structure.REFERENCESArrow, K.J., Hurwicz, L., and Uzawa, H. (1958).
Studies inlinear and non-linear programming , volume 2. StanfordUniversity Press.Boyd, S.P. and Vandenberghe, L. (2015).
Convex opti-mization . Cambridge Univ. Press, Cambridge, 18. printedition.De Persis, C. and Monshizadeh, N. (2016). A modulardesign of incremental lyapunov functions for microgridcontrol with power sharing. In
Proc. European ControlConference , 1501–1506. IEEE, Piscataway, NJ.D¨orfler, F., Bolognani, S., Simpson-Porco, J.W., andGrammatico, S. (2019). Distributed control and opti-mization for autonomous power grids.
Proc. EuropeanControl Conference , 2436–2453.Farrokhabadi, M., Ca˜nizares, C.A., Simpson-Porco, J.W.,Nasr, E., Fan, L., Mendoza-Araya, P., Tonkoski, R.,Tamrakar, U., Hatziargyriou, N.D., Lagos, D., Wies,R.W., Paolone, M., Liserre, M., Meegahapola, L., Ka-balan, M., Hajimiragha, A.H., Peralta, D., Elizondo, M., Schneider, K.P., Tuffner, F., and Reilly, J.T. (2019).Microgrid stability definitions, analysis, and examples.
IEEE Trans. Power Systems . To appear.Jokic, A., Lazar, M., and van den Bosch, P.P.J. (2009).On constrained steady-state regulation: Dynamic kktcontrollers.
IEEE Trans. Autom. Control , 54(9), 2250–2254.Jouini, T., Arghir, C., and D¨orfler, F. (2016). Grid-friendlymatching of synchronous machines by tapping into thedc storage.
IFAC-PapersOnLine , 49(22), 192–197.K¨olsch, L., Bhatt, K., Krebs, S., and Hohmann, S. (2019a).Steady-state optimal frequency control for lossy powergrids with distributed communication. In
IEEE Inter-national Conference on Electrical, Control and Instru-mentation Engineering . To appear. Preprint availableat https://arxiv.org/abs/1909.07859.K¨olsch, L., Dupuis, M., Bhatt, K., Krebs, S., andHohmann, S. (2019b). Distributed frequency reg-ulation for heterogeneous microgrids via steadystate optimal control. In
IEEE Green Technolo-gies Conference . To appear. Preprint available athttps://arxiv.org/abs/1910.03115.Machowski, J., Bialek, J.W., and Bumby, J.R. (2012).
Power system dynamics: Stability and control . Wiley.Magn´usson, S., Fischione, C., and Li, N. (2017). Volt-age control using limited communication.
IFAC-PapersOnLine , 50(1), 1–6.Mohagheghi, E., Alramlawi, M., Gabash, A., and Li, P.(2018). A survey of real-time optimal power flow.
Energies , 11(11), 3142.Monshizadeh, P., De Persis, C., Stegink, T., Monshizadeh,N., and van der Schaft, A.J. (2017). Stability andfrequency regulation of inverters with capacitive inertia.In
Proc. IEEE Conf. Decision Control , 5696–5701.Stegink, T.W., Persis, C.D., and van der Schaft, A.J.(2017a). Stabilization of structure-preserving powernetworks with market dynamics.
IFAC-PapersOnLine ,50, 6737–6742.Stegink, T., De Persis, C., and van der Schaft, A.J. (2015).Port-hamiltonian formulation of the gradient methodapplied to smart grids.
IFAC-PapersOnLine , 48(13),13–18.Stegink, T., De Persis, C., and van der Schaft, A.J.(2017b). A unifying energy-based approach to stabilityof power grids with market dynamics.
IEEE Trans.Autom. Control , 62(6), 2612–2622.Trip, S., B¨urger, M., and De Persis, C. (2016). An internalmodel approach to (optimal) frequency regulation inpower grids with time-varying voltages.
Automatica , 64,240–253.van der Schaft, A.J. (2017).