—Nowadays, we have witnessed a signiﬁcant increase of renewable sources across the existing power grids, as a result of the undeniable need to mitigate the effects of en- vironmental pollution and assist the economic development of countries around the globe. However, due to the lack of rotational inertia and their intermittent nature, renewables have also increased the severity of the occurring disturbances and instability phenomena, which in turn place a demand for more accurate modeling of the network and the power system components. In this paper, we propose a novel network formulation for AC microgrids that allows the consideration of less conservative conditions for the network, while allowing the adoption of higher order inverter-based generation dynamics and facilitating the stability analysis and control. In particular, we formulate the network as a dynamical multi-input/multi- output system, which we show is both stable and passive, even if the network’s dynamic and lossy nature is taken into account. We then brieﬂy discuss the opportunities and the advantages offered by the proposed approach and provide useful guidelines for a more accurate modeling of power system components. Finally, we verify the proposed network formulation through a numerical application on the Four Machine Two Area Kundur test system and several simulations on typical Medium Voltage (MV) distribution lines. Index Terms —AC microgrids, network dynamics, passivity, low inertia power systems, stability analysis. I. I NTRODUCTION During the last decades, there is an ongoing, worldwide effort to reduce the greenhouse gas emissions and thus to decelerate the climate change and the global warming. This effort led to the shift of the energy production from fossil- fueled plants towards Renewable Energy Sources (RES) [1]. In recent years, the latest technological advancements in power electronics have facilitated the incorporation of RES in power grids and as a result, their share in power generation has signiﬁcantly increased [2]. Such large share of RES, how- ever, imposed new challenges that were not encountered in traditional power grids. Speciﬁcally, the majority of RES are interfaced to the network through power electronic convert- ers, and unlike the existing Synchronous Generators (SGs), they do not provide any rotational inertia. The resulting signiﬁcant reduction of rotational inertia in combination with the intermittent nature of RES, made frequency and voltage deviations during disturbances more steep and the instability phenomena that occur more severe [3], [4]. Furthermore,traditional frequency and voltage control mechanisms cannot be effectively applied, as they now become too slow with respect to the disturbance dynamics and thus, unable to pre- vent or even effectively damp the occurring large frequency and voltage deviations. The concept of microgrids was proven able to mitigate several of the above problems since apart from loads, a microgrid consists of a combination of generation units and energy storage elements, namely Distributed Energy Resources (DER) [5], that are able to provide both voltage and frequency support [6]–[8]. For the stability analysis and control design in AC microgrids, the majority of recent literature adopts a static network formulation, using either active and reactive power ﬂows, or the current ﬂows across the grid. Such examples can be found in [2], [9]–[12]. How- ever, due to the lack of rotational inertia and the low R=X ratio of MV/LV distribution lines, the network’s response becomes now comparable to the response of the DER. Thus, the adoption of such static representation for distribution networks, leads to more conservative stability results than those derived in the case of SGs since the dynamic coupling between the bus dynamics with the grid becomes important, especially for the appropriate design of RES controllers. In this paper, we consider a dynamic representation for the current ﬂows across the grid, which we then employ in order to formulate the network as a multi-input/multi- output system. The proposed formulation, not only provides a more accurate representation of the network dynamics, but most importantly, it can facilitate the stability analysis and control in AC microgrids while capturing the dynamic coupling between power grid components and the network. Particularly, by exploiting ﬁndings in classical control and passivity-based analysis and design, we show that every network with arbitrary topology constitutes both a stable and a passive system, even if the network’s dynamic and lossy nature is taken into account . Furthermore, we discuss the opportunities and the advantages offered by the proposed approach and provide useful guidelines for a more accurate modeling of power system components. Finally, the proposed formulation is veriﬁed through a numerical application on the Four Machine Two-Areas Kundur test system and several simulations on typical MV distribution lines. The rest of the paper is structured as follows. In Section II, we present the network dynamic model which is then formulated as a multi-input/multi-output system. The main stability results are then derived in Section III. In Section IV, we discuss the opportunities and the advantages of the proposed formulation, while providing useful information regarding the modeling and the incorporation of bus dy- namics into stability analysis. Finally, Section V provides a numerical veriﬁcation of the proposed approach and several simulations on typical MV distribution lines. Conclusions and future work are presented in Section VI. II. F RAMEWORK FORMULATION In this section, we ﬁrst present the dynamic model that will be used to represent the line currents across the power grid. Then, we deﬁne the net injected current components at every bus, which are subsequently employed in formulating the net- work as a multi-input/multi-output system. The formulation is based on the system reference approach presented in [13], [14], with the main difference that the network’s dynamic behavior is also considered. A. Network model A power network with arbitrary topology can be described by a connected and undirected graph (N;E), whereN= f1;2; : : :jNjg is the set of buses and E  N N = f1;2; : : :jEjg the set of lines connecting them. Similarly to [15], the network structure can be represented by its corresponding incidence matrix E2RjNjjEj. By arbitrarily labeling the ends of the line lwith a +and a, the matrix Eis given by Eil=8 >< >:+ 1 ifiis the positive end of l 1 ifiis the negative end of l 0 otherwise.(1) In order to derive the state equations of the line currents, we now make the following assumptions. Assumption 1: Distribution lines at distribution level can be represented by symmetric RL elements. Assumption 2: The network frequency !, is almost con- stant at synchronous value !s(50 or 60 HZ), i.e. !!s0. When network frequency is expressed in the per unit system, the above approximation implies that !1pu. Assumption 1 states that any line at MV or LV level can be represented by a series impedance. This approximation is very accurate since in AC microgrids we deal with short distribution lines where shunt capacitances are negligible [16]. Moreover, in Assumption 2, we consider that the variations of the network frequency are very small which is a mild assumption considering that the maximum frequency deviation in the European Network of Transmission System Operators for Electricity (ENTSO-E) system is 200mHz (0:4%) [17].We now deﬁne the phasors of the current of line l2Eand the voltage of bus i2N in their rectangular form. We get: ^Il=Ia;l+jIb;land^Vi=Va;i+jVb;i (2) where Ia;landIb;lare the current components of line land, Va;iandVb;iare the voltage components of bus i. It should be noted that the line current and bus voltage are expressed on a common system reference frame , i.e. two common axes that rotate at a speciﬁc velocity ![18]. Based on the dynamic phasor representation provided in [19], the state equations of line current of the line lare given by: _Ia;l=Rl LlIa;l!Ib;l+1 Ll(Va;iVa;j) (3) _Ib;l=Rl LlIb;l+!Ia;l+1 Ll(Vb;iVb;j) (4) where Va;i,Vb;i,Va;jandVb;jare the voltage components at buses iandjwhich are connected through line l.Rl andLldenote the resistance and the inductance of the line lrespectively. Considering that Assumption 2 holds, the differential equations are further simpliﬁed as follows: _Ia;l=Rl LlIa;lIb;l+1 Ll(Va;iVa;j) (5) _Ib;l=Rl LlIb;l+Ia;l+1 Ll(Vb;iVb;j): (6) B. Dynamical input/output formulation In order to derive a dynamical, multi-input/multi-output representation for the network such as the one illustrated in Figure 1, we ﬁrst need to derive the net injected current components at each bus of the network. By employing the incidence matrix E, the net injected current components at busi= 1;2; : : :Nare deﬁned by: Inet a;i=jEjX l=1EilIa;l and Inet b;i=jEjX l=1EilIb;l: (7) We should note here that equations (7) are actually the application of Kirchoff’s Current Law at each bus of the grid. Then, we introduce the vectors Ia= [Ia;1Ia;2: : : Ia;jEj]T, Ib= [Ib;1Ib;2: : : Ib;jEj]T,Inet a = [Inet a;1Inet a;2: : : Inet a;jNj]T, Inet b= [Inet b;1Inet b;2: : : Inet b;jNj]T,Va= [Va;1Va;2: : : V a;jNj]T andVb= [Vb;1Vb;2: : : V b;jNj]Tand consider the following dynamical system with inputs the vectors of bus voltage components VaandVb, states the vectors of line current components IaandIb, and outputs the vectors of net injected current components Inet aandInet b: _Ia _Ib =KAI I K AIa Ib +KB 0 0KBVa Vb (8)  Inet a Inet b = KC 0 0KC Ia Ib : (9) The matrix I2RjEjjEjis the identity matrix while the ma- trices KA2RjEjjEj,KB2RjEjjNjandKC2RjNjjEj can be deduced from the set of differential equations (5)-(7) as follows KA;ml =8 < :Rl Llifm=l 0 otherwise(10) KB;li=8 >>>>< >>>>:1 Llifiis the positive end of l 1 Llifiis the negative end of l 0 otherwise(11) and KC=E: (12) As it can be seen from the above formulation, the network constitutes a 2jNj-input2jNj-output LTI system with min- imal state space representation  =fA; B; C; Dgwhere A=KAI I K A ; B=KB 0 0KB ; C=KC 0 0KC andD= 0: III. M AINRESULTS In this section, we state our main results when the proposed multi-input/multi-output network formulation is adopted. The ﬁrst result shows that the network system (8)-(9) is stable, while the second result demonstrates that the network con- stitutes a passive system, even when the dynamic and lossy nature of the lines is taken into account. Both results are independent of the network topology. The proofs can be found in the Appendix. Theorem 1: (Absolute Stability) Suppose that Assumptions 1-2 are satisﬁed. The network system deﬁned in (8)-(9) with inputs the vectors of bus voltage components VaandVb, Fig. 1. A graphical illustration of the network formulated as a multi- input/multi-output system.states the vectors of line current components IaandIb, and outputs the vectors of net injected current components Inet a andInet bis stable. Remark 1: Theorem 1 states that the dynamical system (8)- (9) representing the net injected current components across the network, is stable. Furthermore, the absolute stability implies that the system is bounded as well, i.e. if the system is subjected to a bounded input or disturbance, its response is also bounded in magnitude. Lemma 1: (Passivity) The network system deﬁned in (8)- (9) with inputs the vectors of bus voltage components Vaand Vb, states the vectors of line current components IaandIb, and outputs the vectors of net injected current components Inet aandInet bis passive. Remark 2: Within the proof of Lemma 1, we make use of the Positive-real Lemma for LTI systems in order to show that the network system deﬁned in (8)-(9) constitutes a multi- input/multi-output passive system. The Positive-real Lemma is utilized using a linear matrix inequality as described in [20], [21]. We also highlight here the connection between the absolute stability and the property of passivity. Particularly, a passive system is also stable yet not vice versa, which in turn implies that the property of passivity is a more strict condition than stability. As we are about to discuss in the next section, this property can be used either for the derivation of global stability results and/or the design of more accurate control mechanisms in a completely decentralized manner. IV. D ISCUSSION A. Input/Output Network Formulation In this paper, the network is formulated as a dynamical multi-input/multi-output system. This formulation is based on the approach presented in [13], [14] where we proved that the static equations describing the interconnections in any power grid with arbitrary topology constitute a passive system, even when taking into account the network’s lossy nature. Both works rely on the fact that the analysis is carried out in the system reference frame , instead of each bus local dq reference frame. This choice of the frame of reference, allows the consideration of less conservative conditions for the network since both the active and the reactive power ﬂows across the grid are now taken into account and the bus voltage magnitudes are not considered to remain constant during disturbances. Within the current work, we additionally take into account the network’s dynamic behavior. Such a dynamic network representation becomes crucial since the dynamic interaction of inverter-based DER with the rest of the system constitutes an important aspect in the analysis of the future power grids where RES share in power generation will dominate. Particularly, the dynamics of inverter-based DER are on similar time scales as the line dynamics while their controls are also signiﬁcantly faster than synchronous generators’ control mechanisms [4]. It is therefore clear that the con- sideration of the dynamic behavior of the network is critical for the accuracy of the stability analysis and control design. Moreover, the dynamic coupling between inverter-based DER and the network is in many cases unstable despite that the capability of DER to employ fast-acting control mechanisms may lead to the expectation for more efﬁcient frequency and voltage support [4], [19]. Finally, using the Positive-real Lemma, provides us an additional, important tool for the stability analysis and control in AC microgrids. Speciﬁcally, the matrix P2R2jEj 2jEj, which is deﬁned in the proof of Lemma 1, can be employed so as to derive a storage or a candidate Lyapunov function. One such storage or candidate Lyapunov function Vnetis provided below Vnet(Ia; Ib) = [IT aIT b]PIa Ib : (13) The aforementioned function, in turn, can be used to deduce stability results for the interconnected system, either by the means of passivity, or by Lyapunov stability analysis [20], [22], and it can be employed both in linear and nonlinear stability approaches. B. Bus Dynamics The incorporation of bus dynamics into the stability anal- ysis when such a network formulation is adopted can be carried out in a similar way as in [13], where both the network and the bus dynamics are expressed in the same reference frame instead of the local dq-coordinates of each bus. Particularly, we consider that each of the jNj buses forms a 2-input2-output system while buses that do not consist of any power system component are represented by zero dynamics, so as to ﬁt with the network formulation described in the previous sections. In Figure 2, we provide a graphical representation of the interconnected system where the multi-input/multi-output network system is connected to the aggregated bus dynamics. Bus dynamics are represented by a broad class of systems of the following form _x=f(x; u) (14) y=h(x; u); (15) with inputs the phasor components of the net injected currents u= (Ia;i;Ib;i)2R2, states x2XRnand outputs the phasor components of the bus voltages y= (Va;i; Vb;i)2R2. The use of such a broad class of systems for the bus dynamics provides us the advantage to include a Fig. 2. The power network represented as an interconnection of input/output systems associated with the bus and network dynamics, respectively [13].variety of power system components such as synchronous generators and inverter-based DER. Components such as loads and Flexible Alternating Current Transmission System (FACTS) devices can be introduced in a similar manner as in [14]. Moreover, such representation gives the opportunity to use more accurate, higher order dynamical models and incorporate voltage and frequency regulation mechanisms. Additionally, the structure of the interconnected system al- lows the derivation of completely decentralized results while facilitating the design of distributed voltage and frequency control mechanisms. The stability analysis can be carried out using techniques from either passivity and/or Lyapunov stability analysis, and applied to both linear and non-linear bus dynamics. For example, in [13], we presented certain decentralized input passivity conditions on the bus dynamics, which when satisﬁed guarantee the asymptotic stability of the equilibria of the interconnected system. Several variations of the aforementioned passivity conditions along with other passivity-based control design techniques can be also found in [20]. We highlight here the ease of application of passivity- based control for linear systems, or systems linearized about an equilibrium where the passivity property can be easily veriﬁed by means of computationally efﬁcient methods using the KYP lemma [23], or via the positive realness of the transfer function matrix. On the other hand, for nonlinear systems, it can be veriﬁed by exploiting structural properties, such as the interconnections of passive systems, or via an explicit construction of the transfer function. Finally, although the transformation from local dq- coordinates to the system reference frame adds complexity to the analysis of the bus dynamics, it does not involve neighboring variables, such as the bus voltages and angles, into the local dynamical system. This allows the derivation of decentralized stability results for the bus dynamics, providing simultaneously plug-and-play capability. Furthermore, it can facilitate the design of distributed control mechanisms that can contribute to the overall power system stability. Exam- ples of such decentralized mechanisms using passivity-based techniques can be found in [24], [25]. V. N UMERICAL APPLICATION AND SIMULATIONS In this section, we ﬁrst verify the proposed framework formulation through a numerical application on the Four Machine Two-Areas Kundur test system [26]. The dynamical model that represents the network of the aforementioned system is derived using equations (8)-(12). By arbitrarily numbering the buses and the lines, the matrices KA,KB andKCare therefore given by KA=diagf0;0;0:1;0:1;0:1;0:1;0;0;0:1;0:1;0:1;0:1g KB=2 66666666459:9 0 0 59:9 0 0 0 0 0 0 0 0 59 :9 0 0 0 0 0 59:9 0 0 0 0 0 100 0 0 0 0 100 0 0 0 0 0 9 :1 0 0 0 0 0 9:1 0 0 0 0 9 :1 0 0 0 0 0 9:1 0 0 0 0 0 40 0 0 0 40 0 0 0 0 0 0 0 59 :9 0 0 0 0 59:9 0 0 0 0 0 0 59 :9 0 0 0 0 0 59:9 0 0 0 0 0 0 9 :1 0 9:1 0 0 0 0 0 0 0 0 9 :1 0 9:1 0 0 0 0 0 0 0 0 100 0 0 0 100 0 0 0 0 0 0 0 0 0 40 403 777777775 KC=2 66666641 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 1 1 1 0 011 0 0 1 0 0 0 0 0 0 0 0 0 11 0 0 0 11 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 113 7777775 respectively. Furthermore, from the proof of Lemma 1, we can deduce the matrix Pusing equation (22). We get P1=diagf0:017;0:017;0:01;0:11;0:11;0:025; : : : : : :0:017;0:017;0:11;0:11;0:01;0:025g which veriﬁes the passivity result derived in Section III. In order to show the difference between the static and the dynamic representation of the network we now consider two typical, 5km MV lines, that is an overhead and an underground line respectively, each connecting a generator and a zero load bus. Their parameters are provided in Table I. We then assume that both loads instantly increase to 1MW and we illustrate the deviation of the line current magnitudes in Figures 3 and 4 respectively, for the following three cases: (i) Lossless line representation, (ii) Static line representation, and (iii) Dynamic line representation. We note here that all values are expressed in per unit system ( Sbase= 1MV A and Vbase= 11kV). TABLE I THE PARAMETERS OF THE TYPICAL OVERHEAD AND UNDERGROUND MV DISTRIBUTION LINES . Overhead Underground Rol = 0:01231 =km Rul = 0:00529 =km Xol = 0:01430 =km Xul = 0:00455 =km As we observe from both ﬁgures, the lossless line repre- sentation provides completely inaccurate results and thus the incorporation of such line model in any analysis could yield very conservative stability guarantees. On the other hand, although both the static and the dynamic line representation ensure the accuracy of the analysis, using a static model to represent lines does not allow us to capture the dynamic coupling between the bus dynamics and the network. Thus, the incorporation of the proposed dynamic network model into the stability analysis could assist in the appropriate design of inverter-based RES controllers while providing the necessary stability guarantees even if the analysis is carried out in a completely decentralized manner. VI. C ONCLUSIONS AND FUTURE WORK In this paper, we have presented a novel multi-input/multi- output network formulation that can facilitate the stability analysis and control in AC microgrids. In particular, we formulated the network as a dynamical multi-input/multi- output system, which we showed is both stable and passive even if the network’s dynamic and lossy nature is taken into account. Furthermore, we discussed the opportunities and the advantages offered by the proposed approach while providing useful guidelines for a more accurate modeling of power 0 5 10 15 20 25 Time (sec)00.20.40.60.811.2Current Magnitude |I| (pu) Lossless Representation Static Representation Dynamic RepresentationFig. 3. The magnitude of the line current of a typical MV overhead distribution line. 0 5 10 15 20 25 Time (sec)00.20.40.60.811.2Current Magnitude |I| (pu)Lossless Representation Static Representation Dynamic Representation Fig. 4. The magnitude of the line current of a typical MV underground distribution line. system components. Finally, a brief explanation of how this framework could be exploited to derive decentralized stability results and drive the control design using features of passivity-based analysis is also provided. Our future work will mostly concentrate on the introduction of DER dynamics and the design of inverter-based RES frequency and voltage regulation mechanisms. APPENDIX Proof of Theorem 1: In order to prove the absolute stability of the network system (8)-(9), we need to show that the eigenvalues of the system matrix ANhave negative real part. The proof employs techniques from Matrix Analysis [27] and particularly the matrix deﬁniteness. The sufﬁcient condition to guarantee that the eigenvalues of the matrix ANlie in the left half plane, is that ANis negative deﬁnite. Let us now consider the vector wT= [wT 1wT 2]2R2jEj. The following inequality should therefore be satisﬁed: [wT 1wT 2]ANw1 w2 =wT 1KAw1+wT 2KAw20 (16) for all w1; w 22RjEj. From inequality (16) we observe that the negative deﬁniteness of ANis ensured when its diagonal sub-matrices KAare negative deﬁnite matrices as well. KA2RjEjjEjis a square, diagonal matrix with negative diagonal elements, i.e. KA;ll08l= 1;2; : : : ;jEj. Thus, the eigenvalues of KAare real, negative and equal to its diagonal elements. This immediately leads to the fact that KAis negative deﬁnite. The system matrix ANis also a negative deﬁnite matrix since condition (16) is satisﬁed.  Proof of Lemma 1: The proof of Lemma 1 is based on the Positive-real Lemma which is provided in [20], [21]. According to Positive-real lemma, a stable LTI system such the one deﬁned in (8)-(9) with DN= 0, is passive if and only if there exist a positive deﬁnite matrix P=PT2R2jEj 2jEj such that: (AN)TP+PAN<0 (17) (BN)TP=CN: (18) We now substitute the matrices AN,BNandCNof the dynamical system (8)-(9) into the above conditions (17)-(18) and we get: KAI I K A P+PKAI I K A <0 (19)  KT B 0 0KT B P= KC 0 0KC (20) where P2R2jEj 2jEj. From (11)-(12), we observe that KT B andKCare of the exact same structure. Speciﬁcally, the elements of KT BandKC2RjNjjEjsatisfy the equality: KB;li=1 LlKC;il=1 LlEil where KB;lidenote the elements of KT B. The matrix Pcan therefore be deﬁned from (20) as: P= P10 0P1 (21) where P12RjEjjEj. The matrix P1is a diagonal (thus symmetric) matrix whose elements are given as follows P1 ml=( Llifm=l 0 otherwise(22) It is straightforward to show that Pis positive deﬁnite since its diagonal submatrices, or equivalently its diagonal elements, are always positive i.e. P1 ll>0. We then examine if the inequality (19) also holds. By substituting (22) into (19) we get  KAI I K A P10 0P1 + P10 0P1 KAI I K A = 2P1KA 0 0 P1KA :(23) Since both matrices KAandP1are diagonal and have negative and positive diagonal elements respectively, their multiplication which is also commutative, results to a neg- ative deﬁnite matrix. Consequently the inequality (17) is satisﬁed, that is (AN)TP+PAN= 2P1KA 0 0 P1KA <0: (24) The above inequality completes the proof. 