{'fraction_non_alphanumeric': 0.06220306176315327, 'fraction_numerical': 0.022127397501319725, 'mean_word_length': 3.9473340587595214, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 13, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
A Port-Hamiltonian Approach to Modeling and Control of an Electro-Thermal Microgrid

: We address the problems of modeling and controlling multi-energy microgrids (meMGs) composed of an electrical and a thermal system, which are connected via heat pumps (HPs). At ﬁrst, we model the individual subsystems in a port-Hamiltonian (pH) framework. Then, by exploiting the structural properties of pH systems, we interconnect the subsystems in a passive manner and show that the overall meMG is shifted passive with respect to the control input-output mapping. We then use this property to propose a distributed passivity based-control (PBC) that addresses frequency and temperature regulation by utilizing the resources in the meMG in a proportional fashion and renders the closed-loop equilibrium asymptotically stable.


MOTIVATION
As a leaping step towards sustainable and zero-carbon energy systems, integration of renewable-energy-based intelligent technologies have increased remarkably in various energy-sectors, e.g., microgrids (MGs) (Hatziargyriou et al., 2007), district heating systems (DHSs) (Scholten et al., 2017) and many more. There are remarkable benefits when a multi-energy system is considered as a whole, e.g., energy losses, costs and emissions can be reduced (Geidl and Andersson, 2007). An example would be a DHS with combined heat and power (CHP) units and heat pumps (HPs) that act as links between heat and electricity networks, enabling better system flexibility and utilization of renewable energy sources (RESs) (Lund et al., 2017).
Despite of the advantages mentioned above, it is often still the case that clean energy technologies are defined and understood separately, i.e., as single energy sectors, which hinders achieving better overall system flexibility and efficiency (Lund et al., 2017). For example, the term microgrid often solely focuses on the electrical sector, which then misses a broader view that includes the aspects of heating and cooling in its design. Some interesting works that combine electrical and thermal systems are mentioned in the sequel. In (Awad et al., 2009), the problem of optimal power flow is investigated for a combination of electric and thermal networks, where in (Mancarella et al., 2011), the effects of HPs and CHPs on the operation of low voltage networks have been studied and it shows how a smart combination of electric and thermal networks can efficiently deal with the increasing integration of RESs.
As multi-energy systems are complex, heterogeneous and nonlinear networked systems, modeling and control is a challenging task. In this regard, port-Hamiltonian (pH) systems and passivity-based-control (PBC) (Ortega et al., 2002;van der Schaft, 2000) provide a promising structurebased framework to model and control such systems. In (Ramirez et al., 2013), the pH framework has been used to model a multi-energy system that comprises thermodynamic and mechanical systems. In (Strehle et al., 2018), the advantages of pH modeling for coupled electricity and gas systems have been investigated through a case study. Furthermore, (Fiaz et al., 2013) presents a pH approach to modeling of power networks, while (Hauschild et al., 2020) introduces a pH-based model for DHSs. In (Machado et al., 2020), the authors present a model of a DHS and show that the system is shifted passive.
Inspired by these promising approaches, in the present paper we consider a multi-energy microgrid (meMG) which is composed of an electrical and a thermal system (DHS) which are connected via HPs. In this setting, the two main contributions of this paper are as follows. At first, we derive models of the individual subsystems of the meMG as pH systems and interconnect them in a passivity-based manner such that the overall meMG results again in a pH system. This structure is then exploited to show that the meMG is shifted passive with respect to the control input-output mapping. The latter result is instrumental for deriving our second main contribution, which is a distributed PBC for the integrated meMG using consensus protocols. The control law addresses the objectives of frequency regulation, supply temperature regulation and achieves an identical marginal cost condition between the controllable units in the meMG. Finally, we show that the closed-loop equilibrium of the meMG in feedback with the PBC is asymptotically stable by using the composite Hamiltonians as a Lyapunov function.

Notation
We denote by I n the n × n identity matrix, by 0 n×m the n × m matrix with all entries equal to zero and by 1 n the column vector of size n with all entries equal to one. Let Z ∈ R n×n be a symmetric matrix. If Z is positive (negative) definite, we denote this by Z > 0 (Z < 0). If Z is positive (negative) semidefinite, we denote this by Z ≥ 0 (Z ≤ 0). For a set N , |N | denotes its cardinality and [N ] k denotes the set of all subsets of N that contain k elements. For a matrix A, |A| denotes the matrix whose entries are the absolute values of those of A. For a scalar x i ∈ R, col(x i ) denotes a column vector and diag(x i ) denotes a diagonal matrix. A block diagonal matrix with block diagonal entries A i ∈ R n×n is denoted by blockdiag(A i ). Furthermore, for a function f (x), ∇f denotes ((∂f )/(∂x)) .

Electrical network model
We consider a Kron-reduced model, see (Kundur, 1994), for the electrical system of the meMG, whose topology is represented by a weighted undirected and connected graph M = (N , F), where N = {1, 2, . . . , n} is the set of nodes and F = {1, 2, . . . , m} denotes the set of edges. Nodes can represent either DG units or HP units. More precisely, N = G ∪ H where G (g := |G|) denotes the set of DG nodes and H (h := |H|) the set of HP nodes. Furthermore, each edge f l ∈ F represents a power line. Following common practice (Guerrero et al., 2005;Schiffer et al., 2014), we assume that the power lines in the electrical network are purely inductive. Then, the l-th edge connecting nodes i and j has a positive edge weight w l = |B ij |, where B ij ∈ R <0 is the line susceptance.

Thermal network model
The DHS topology of the meMG is represented by a connected graph S = (D, E) with the set of nodes D = {1, 2, . . . , d} and the set of edges E = {1, 2, . . . , e}. Following (Machado et al., 2020), we model heat producers, heat consumers and HPs as pipes with injecting or extracting heat flows. Thus, these components fall into the class of two-terminal devices and are therefore represented by edges in the DHS. Then, E = H ∪ P ∪ C, where the sets H, P and C represent the set of HPs, heat producers and heat consumers, respectively. Thus, the set H denotes HPs in both the electrical system (as nodes) and the thermal system (as edges). Furthermore, we make the technical assumption that the HPs and the heat producers are connected to the DHS via storage tanks that have hot and cold layers separated by a thermocline, see e.g., (Machado et al., 2020, Assumption 1) and there are no stand-alone storage tanks. Thus, the set of nodes D represents the hot and cold layers of the storage tanks and junctions 1 .
Let d := |D| and e := |E| = h + p + c where p := |P|, c := |C| and define the flow-dependent incidence matrix B ∈ R d×e with entries B ij := 1 if the flow through edge j targets node i, −1 if the flow through edge j originates from i, 0 otherwise.
Following (Machado et al., 2020), we define the matrices (2.1) 1 In contrast to the assumption that HPs and producers are interfaced by a storage tank, another option would be to combine a HP and a tank (similarly, a producer and a tank) and represent the combined system as an edge. We remark that the modeling procedure presented in this paper can be applied to the above-mentioned edge consideration as well (Machado et al., 2020).
which are partitioned with respect to HPs (h), heat producers (p) and heat consumers (c). Both T and S have the property that all entries are non-negative and every column has at most a unity (one) entry. Thus, (Machado et al., 2020) (2.2) We assume that the volumes and the edge flows in the DHS are constant and positive 2 , which is formalized below.

Heat pump model
The main components of a HP are an evaporator, a compressor, a condenser and an expansion valve, where the heat absorbed by the refrigerant from the ambient air (or water) through the condenser is exchanged by the evaporator with the water to be heated (Kim et al., 2015). In this process, the compressor and the expansion valve are responsible for increasing the temperature and pressure of the refrigerant to the required value for the operation of the HP cycle. The compressor is driven using a motor, which usually is an induction motor or a permanent magnet synchronous motor (PMSM). Due to the advancements in power-electronic-based drive technologies, usage of PMSMs in servo applications has increased, mainly due to the fact that PMSMs offer better efficiency, power density and reliability in comparison to induction motors (Manias, 2017, Chap. 12). Thus, we employ PMSMs to drive the compressors and thus, representing the electrical model of the HPs. Furthermore, to represent the thermal side of the HPs, we assume that the evaporator of the HP injects a heat flow to a pipe that has an inlet and an outlet through which cold and hot water, respectively, flows. Then, by combining the individual components of the HP as one single heat exchange unit, the HP can be modeled as a pipe that injects thermal power into the DHS. Such a HP model is analogous to the typical models used for heat producers in DHSs, see e.g., (Machado et al., 2020;Scholten et al., 2017).
be the temperature, volume and flow, respectively through the i-th HP and T ini ∈ R ≥0 → R be the temperature of the node from which the flow towards the i-th HP originates from. Then the combined electrical (PMSM) and thermal (pipe) model of the i-th HP, i ∈ H, is given by (Kim et al., 2015;Kundur, 1994;Machado et al., 2020), where γ ki and χ ki , k = 1, 2, 3, 4, are positive constants, see (Lee et al., 2019, Table 1), T a ∈ R >0 denotes the constant ambient temperature, P i : R ≥0 → R is the active power flow, P mi : R ≥0 → R is the shaft power of the PMSM and Q com hi : R ≥0 → R is the heat flow at the compressor. We make the following technical assumption on the interaction between the electrical and the thermal power flows in the HP system (2.3). Assumption 2.2. γ 2i = χ 1i ∀i ∈ H.
3. PORT-HAMILTONIAN meMG MODEL In this section, we model the meMG by using a pH modeling approach. To this end, we divide the overall system into the following subsystems.
(1) Electrical system: Combination of electrical network, static loads, DG units (synchronous generators and grid-forming inverters) and HPs.
(2) Thermal system: Combination of HPs, heat producers, heat consumers, storage tanks and junctions. In the sequel, we present the pH models of the different subsystems mentioned above. Then, by using their pH structure, the overall interconnected system is also obtained in a pH form, which we then exploit to show shifted passivity with respect to the control in-and outputs. The latter property is essential for the subsequent PBC design in Section 4.

Port-Hamiltonian model of the electrical system
We divide the electrical system into two separate systems as explained in the sequel. Electrical network and loads: To each node i ∈ N , we associate a time-varying phase angle θ i : R ≥0 → R, a relative frequencyθ i = ω i − ω com and a constant voltage magnitude V i ∈ R >0 . Here, ω i is the absolute frequency of the unit at the i-th node and ω com ∈ R >0 that of the common rotating (dq-) frame in the network model. The vectors of phase angles and absolute frequencies are denoted by θ = col(θ i ) ∈ R n and ω = col(ω i ) ∈ R n . respectively. We then define the potential function 1) whose first summation term represents the power stored in the power lines (Stegink et al., 2019) and P li ∈ R denotes the power consumed by the static load connected at node i.
The trigonometric terms in U only depend on the phase angle differences θ i − θ j . It is well known, that this fact stymies the construction of a (locally) positive definite Hamiltonian function (Schiffer and Dörfler, 2016;Schiffer et al., 2014). To circumvent this challenge, following (Schiffer and Dörfler, 2016;Schiffer et al., 2014), we choose an arbitrary node, say the n-th node, and express all other phase angles, with respect to that node, i.e., x θ = Rθ ∈ R n−1 , R = [I n−1 −1 n−1 ] ∈ R (n−1)×n .

4) with
(3.5) Then, the system (3.3) can be written in pH form as followṡ where d θ = −ω com 1 n and the interconnection, damping and input matrices are given by J θ = R θ = 0 (n−1)×(n−1) , G θ = R ∈ R (n−1)×n . Furthermore, with (3.5), we have that, see (Schiffer and Dörfler, 2016), where P i is the active power flow, i ∈ N . Thus, the output y p θ of the electrical network (3.6) represents the active power flows at the DG units and the HPs, see e.g., (Schiffer and Dörfler, 2016), while the input u p θ represents their electrical frequencies.
DG units and PMSMs: The DG units are composed of synchronous generators and grid-forming inverters where, for the latter, we assume that they are equipped with the standard frequency droop control and their power outputs are measured using low pass filters. In this setting, the second-order swing equation of a synchronous generator (Kundur, 1994) and that of a grid-forming inverter admit the same model, see e.g., (Schiffer et al., 2013, Lemma 4.1). Then, for i ∈ G, the model of a DG unit is given by (Kundur, 1994) M is as defined in (2.3). Furthermore, P d i ∈ R is the active power set point and ω d ∈ R >0 is the nominal frequency. Now, we combine the electrical side of the HP model (2.3) with the DG model (3.8) and express them in a pH form. For this purpose, let ω g := col(ω i ), u c g := col(u c ωi ), d g := col(D i ω d + P d i ) for i ∈ G as well as ω m := col(ω k ), T h := col(T k ), u c m := col(u c ω k ), d m := col(γ 3 k T a + γ 4 k ) for k ∈ H and d ω := col(d g , d m ). Furthermore, for i ∈ G and k ∈ H, we define the matrices Then, by introducing the state vector of frequencies x ω = ω = col(ω g , ω m ) and the Hamiltonian function the electrical side of (2.3) and the DG model (3.8) can be written together as the pH systeṁ where P is defined in (3.7). Here, y p ω ∈ R n+h and u p ω ∈ R n+h denote the physical (interconnection) output and input, respectively, while y c ω ∈ R n and u c ω ∈ R n are the control output and input, respectively.

Port-Hamiltonian model of the thermal system
The employed DHS model is inspired by (Machado et al., 2020;Scholten et al., 2017), with the main difference that in the present setting, HPs are considered. This is essential to facilitate a physical interconnection with the electrical system and, of course, also contribute to heat generation and control. At first, we present the models of the heat producers and heat consumers, which, aside from the HPs modeled in (2.3), constitute the edges in the DHS. Recall from Section 2.3 that E = H ∪ P ∪ C. Then, the thermal model of the heat producers and consumers in the DHS is given by (Machado et al., 2020, Eq. 11 where all the quantities are defined analogous to (2.3). In addition, Q ci ∈ R >0 denotes the constant heat flow rate of the i-th consumer and Q pi : R ≥0 → R denotes the heat flow rate of the i-th producer. The producer flow rates Q pi , i ∈ P, are used as control inputs in Section 4. Next, we present the node model which represents the hot and cold layers of a storage tank as well as simple junctions 3 . The dynamics of node j ∈ D is given by (Machado et al., 2020, Eq. 12 where E j is the set of edges from which the flows are directed towards node j, T j : R ≥0 → R is the temperature, v j ∈ R >0 is the volume at node j, j ∈ D. Finally, q k ∈ R >0 is the flow from the edge k and T k : R ≥0 → R its temperature. As detailed in Section 2.3, both v j and q k are assumed constant. Next, we combine the thermal part of the HP (from (2.3)) with (3.11) and (3.12) to derive a compact representation of the DHS in a pH form, which is inspired by the ideas from (Hauschild et al., 2020). To this end, we define the positive definite matrices which collect all the volumes and flows in the DHS. Then, by using T and S from (2.1), we introduce (Machado et al., 2020) A := Q e −Q e S −T Q e Q T ∈ R (e+d)×(e+d) . (3.14) We present an important property of the symmetric part of the matrix A in the following lemma, see also (Machado et al., 2020, Lemma 5). Lemma 3.1. For the matrix A defined in (3.14), Proof. From (2.1), we can express B as B = T − S and by Assumption 2.1, we have that Bq e = 0 d . Hence, with the property (2.2), A satisfies 1 e+d A = 0 e+d and A1 e+d = 0 e+d . As a consequence, we have that (A + A )1 e+d = 0 e+d , which together with the fact that all the diagonal entries of A+A are strictly positive, implies that A+A is a weakly diagonally dominant matrix (Horn and Johnson, 2012). Thus, by Gershgorin's disc theorem (Horn and Johnson, 2012), we have (3.15). We define the following matrices: 16) where A is defined in (3.14), χ 2i in (2.3) and note that J T −R T = −J −1 (A+X 2 )J −1 . It follows from Lemma 3.1, together with the facts that J is a positive definite diagonal matrix and X 2 ≥ 0, that R T ≥ 0.
By collecting all edge and node temperatures in the vector which describes the ectropy of the thermal system (Haddad, 2019). Then, combining (3.11), (3.12) and the thermal model in (2.3), the DHS can be written in the pH forṁ 18) where u p T ∈ R h is the physical interconnection input of the HPs and y p T ∈ R h is their physical interconnection output. Furthermore, u c T = col(Q pi ) ∈ R p , i ∈ P represents the control input of the additional heat producers (see (3.11)), y c T ∈ R p their control output and d T = col(χ 3i T a − χ 4i , Q c k ) ∈ R (h+c) , i ∈ H, k ∈ C, a constant perturbation term.

Port-Hamiltonian model of the interconnected meMG
In this section, we physically interconnect the electrical system (3.6), (3.10) and the thermal system (3.18). By inspection of the physical inputs and outputs of the latter together with (3.7), we see that the physical interconnection results in the following feedback structure: where we defined for ease of notation Furthermore, by defining  (3.4), (3.9) and (3.17), respectively, the overall interconnected meMG is clearly still a pH system and can be written aṡ (3.23)

PASSIVITY-BASED FREQUENCY AND TEMPERATURE CONTROL
In this section, we present a PBC for the meMG (3.22) to address the control objectives of frequency regulation in the electrical system and supply temperature regulation in the DHS side (i.e., regulating the temperatures of the heat producers), which is an important control objective in DHSs, see e.g., (Machado et al., 2020;Vandermeulen et al., 2018). Similarly in electrical networks, frequency regulation is considered to be a fundamental control objective since the electrical units are designed to operate close to the nominal frequency (50 Hz or 60 Hz), see e.g., (Bidram et al., 2013;Schiffer and Dörfler, 2016).
As the considered meMG is a heterogeneous system with different types of controllable units (DGs, HPs and heat producers), the proposed PBC is distributed in nature, which then, achieves the aforementioned control objectives by using local information and utilizes the controllable units in a proportional manner, often termed as the identical marginal cost condition in MGs (Schiffer and Dörfler, 2016, Eq. 6).
For the control design, we make the following standard assumption on existence of a desired equilibrium solution. Assumption 4.1. The system (3.22) admits an equilibrium solution (x * , u c * ), with

Shifted passivity of the meMG
In this section, we show that the meMG (3.22) is shifted passive (Jayawardhana et al., 2007), i.e., the mapping (u c − u c * ) → (y c − y c * ) is passive. We use the concept of shifted passivity, instead of passivity, for the system (3.22) because the equilibrium of interest is not the origin, but is x * defined in Assumption 4.1 and is quintessential for the PBC design to be presented in the sequel. Lemma 4.2. Consider the the meMG (3.22) with Assumption 4.1. The mapping (u c − u c * ) → (y c − y c * ) is passive with locally positive definite storage function given by (3.21). Proof. The meMG (3.22) is a pH system with constant damping and interconnection matrices R ≥ 0 and J = −J , respectively. Therefore, it follows from (Jayawardhana et al., 2007, Proposition 1) that the system (3.22) is locally passive with respect to the supply rate (u c − u c * ) (y c − y c * ) and admits the storage functionĤ(x) given by (4.1) if the Hamiltonian H(x) given in (3.21) is locally positive definite around x * . To show the latter, by calculating the Hessian of H(x) at x * , we get ∇ 2 H(x * ) = blockdiag(M, ∇ 2 H θ (x * θ ), J). Note that M > 0 and J > 0 by definition. Furthermore, with Assumption 4.1, ∇ 2 H θ (x * θ ) (with H θ defined in (3.4) and x θ in (3.2)) represents the principal minor of a Laplacian matrix (Schiffer et al., 2014, Lemma 5.8) and is thus positive definite. Hence, ∇ 2 H(x * ) > 0, completing the proof.

A distributed PBC law for frequency and temperature regulation
We propose a hierarchical control structure with a secondary control layer and a tertiary control layer, where the former is responsible for frequency and temperature regulation and the latter is designed for achieving identical marginal costs between the controllable units in the meMG. Let p : R ≥0 → R n and ν : R ≥0 → R p be the secondary control variables and ζ = col(ζ ω , ζ T ) be a tertiary control variable where ζ ω : R ≥0 → R n and ζ T : R ≥0 → R p . Furthermore, let L p ∈ R n×n be the Laplacian matrix of a communication network connecting the DGs and the HPs in the electrical network, L ζ ∈ R (n+p)×(n+p) be the Laplacian matrix of a communication network connecting all the controllable units (DGs, HPs and heat producers) in the meMG and K p ∈ R n×n , A p ∈ R n×n , K ν ∈ R p×p , K ζ,1 ∈ R n×n , K ζ,2 ∈ R p×p , A ζ ∈ R (n+p)×(n+p) be positive definite diagonal matrices. We assume that both the communication networks are connected and undirected.We require the following definition. where ω d ∈ R >0 is the desired synchronization frequency of the electrical system, T d p ∈ R p >0 denotes the reference temperatures at the heat producers, ω com and T * p are defined in Assumption 4.1.
Inspired by the distributed frequency controllers in MGs, e.g., (Schiffer and Dörfler, 2016), and using the shifted output (y c − y c * ) of the meMG (3.22), we propose the following distributed control law 4