Distributed Barrier Certificates for Safe Operation of Inverter-Based Microgrids

TL;DR

This paper introduces a distributed control method using barrier functions and sum-of-squares optimization to ensure voltage safety during transients in inverter-based microgrids, addressing their unique dynamic challenges.

Contribution

It presents a novel approach for designing distributed control laws that guarantee voltage limits during transients in inverter-based microgrids, leveraging barrier functions and sum-of-squares optimization.

Findings

Control laws certify voltage safety during transients

Method demonstrated through numerical simulations

Addresses stability without inertia in microgrids

Abstract

Inverter-interfaced microgrids differ from the traditional power systems due to their lack of inertia. Vanishing timescale separation between voltage and frequency dynamics makes it critical that faster-timescale stabilizing control laws also guarantee by-construction the satisfaction of voltage limits during transients. In this article, we apply a barrier functions method to compute distributed active and reactive power setpoint control laws that certify satisfaction of voltage limits during transients. Using sum-of-squares optimization tools, we propose an algorithmic construction of these control laws. Numerical simulations are provided to illustrate the proposed method.

Equations

$$\dot{x}(t)$$

$${\mathcal{X}}_u:=\{x\in {\mathbb{R}}^n|w_i(x)>0,i=1,2,\dots ,l\}$$

$$B(x)$$

$$({\nabla }_{x}B{)}^{T}f(x)+\alpha \left(B(x)\right)$$

$$p(x)>0\forall x\in \mathcal{K}⟹p\in \left\{{\sum }_{j=0}^m{\sigma }_j{k}_j|{\sigma }_j\in \Sigma [x]\forall j\right\}.$$

$${\dot{\theta }}_i$$

$${\tau }_i{\dot{\omega }}_i$$

$${\tau }_i{\dot{v}}_i$$

$$P_i$$

$$Q_i$$

$$\forall i:P_i=P_i^{\text{set}},Q_i=Q_i^{\text{set}},{\omega }_i=0,v_i=v_i^0.$$

$$P_i^{\text{set}}=P_i^0+u_i^p,Q_i^{\text{set}}=Q_i^0+u_i^q,$$

$$\underline{v_i}\le v_i(t)\le \overline{v_i}.$$

$${\dot{x}}_i$$

$${\mathcal{X}}_{u,i}$$

$$\forall i:$$

$$B_i(x_i)<0\forall x_i\in {\mathcal{X}}_{u,i}$$

$${\dot{B}}_i\ge 0\forall x_i\in \partial {\mathcal{D}}_i[c_i],\forall x_j\in {\mathcal{D}}_j[c_j]\forall j\in {\mathcal{N}}_i$$

$${\dot{B}}_i={\nabla }_{x_i}{B}_i^T(f_i(x_i)+g_i(x_i)u_i+{\sum }_{j\in {\mathcal{N}}_i}{h}_{ij}(x_i,x_j)).$$

$$x_i(0)\in {\mathcal{D}}_i[c_i]\forall i⟹x_i(t)\in {\mathbb{R}}^{n_i}\backslash {\mathcal{X}}_{u,i}\forall i\forall t\ge 0.$$

$${\mathcal{D}}_1[c_1]\times {\mathcal{D}}_2[c_2]\times \cdots \times {\mathcal{D}}_m[c_m]$$

$$\text{(isolated)}{\dot{x}}_i=f_i(x_i),$$

$$\underset{s_2^k,s_3^k,s_4^k}{\max }{\beta }^k$$

$$\widehat{s_1}(p-{\beta }^k)-s_2^k(\hat{V}-1)\in \Sigma [x]$$

$$-s_3^k(1-\hat{V})-s_4^k\dot{\hat{V}}-{\epsilon }_2|x{|}^2\in \Sigma [x]$$

$$\underset{V^k,s_1^k}{\max }{\delta }^k$$

$$V^k-{\epsilon }_1|x{|}^2\in \Sigma [x]$$

$$s_1^k(p-\hat{\beta })-\widehat{s_2}(V^k-1+{\delta }^k)\in \Sigma [x]$$

$$-\widehat{s_3}(1-V^k)-\widehat{s_4}{\dot{V}}^k-{\epsilon }_2|x{|}^2\in \Sigma [x]$$

$$\underset{s_0^k}{\max }z,\text{s.t.}V-z-\sum _{i=1}^l{s}_{0,i}^k{w}_i\in \Sigma [x].$$

$$\underset{s_1^k}{\max }{\epsilon }^k,\text{s.t.}\dot{\hat{B}}+\gamma \hat{B}-{\epsilon }^k-s_1^k\hat{B}\in \Sigma$$

$$\underset{s_{2,i}^k}{\max }\text{trace}(Q^k)$$

$${\dot{B}}^k+\gamma B^k-\eta -\widehat{s_1}{B}^k\in \Sigma [x]$$

$$-B^k-{\sum }_{i=1}^l{s}_{2,i}^k{w}_i\in \Sigma [x]$$

$$\forall i:$$

$${\nabla }_{x_i}{B}_i^T(f_i+g_i{u}_i+{\sum }_{j\in {\mathcal{N}}_i}{h}_{ij})\ge 0$$

$$\underset{u_i(x_i)}{\min }{U}_i$$

$$\text{s.t.}\text{eqrefE:control and }\Vert u_i(x_i){\Vert }_{\infty }\le U_i\forall x_i\in {\mathcal{D}}_i[c].$$

Full text

I Introduction

Safety critical system refers to a system for which the violation of safety constraints will lead to serious economic loss or personal casualty. Power system falls into such category considering the loss resulting from large-scale blackout and critical loads such as hospital and process plant. Modern power system has been evolving towards distributed operation. With the increasing integration of distributed energy resources (DERs), especially renewable resources, challenges have arisen in safely operating power systems as well as guaranteeing stability. Microgrid is a promising direction to tackle the intermittency and uncertainty characteristics that are intrinsic in renewable resources such as wind and solar. An islanded microgrid is a standalone small scale power system that groups a variety of DERs, especially renewables, together with energy storages and loads to provide better control and operation, higher efficiency and reliability [1]. It is a viable solution for power supply to rural area. Microgrid also provides a new perspective to increase the penetration level of renewables in modern power system. Unlike traditional power system that has large inertia from conventional synchronous generators, DERs in a microgrid are connected to the network through power electronic interface. Considering the stochasticity in renewables and the negligible physical inertia, control of voltage and frequency for microgrids is challenging [2]. While stability analysis and control of inverter-based microgrids have received a lot of attention in the literature [3, 4, 5, 6, 7], safety of the microgrids has largely been ignored. For microgrids, the safe region can be defined for voltage magnitudes at every node in the network. The transient voltage in a microgrid can fluctuate by a large amount, causing serious power quality and safety issues, even causing damages to the electrical equipment [8, 9]. Flexible power injections at the droop-controlled inverter nodes can be utilized to stabilize the phase angle, frequency and voltage magnitude, as well as ensuring the voltage magnitudes at all nodes within the safe region.

Control Lyapunov function (CLFs) have long been used to synthesize stabilizing controllers for nonlinear systems [10]. On the other hand, barrier functions are used to certify safety by guranteeing the forward invariance of a set via Lyapunov-like conditions. Although barrier function originated in the field of optimization as a penalty function to replace constraints, it prospers in the field of control design too. For example, [11] considers the safety verification problem in both worst-case and stochastic settings by constructing barrier certificates. The ideas of the barrier functions and the CLF were combined to construct the control barrier functions (CBFs) [12] which have since been used in designing safety controllers. Reference [13] applied CBF method to automotive systems to achieve lane keeping and adaptive cruise control simultaneously with safety constraints. Application of CBF method to establish set invariance with the existence of disturbance and uncertainty is considered in [14]. Simultaneous satisfaction of safety and performance objectives via design is not a trivial task. The stabilization objective expressed by a CLF and the safety constraints established by a CBF can be potentially in conflict. In [15], the authors proposed a quadratic program formulation that unifies CLF and CBF to synthesize a controller that enforces the safety constraints but relaxes the stability (performance) requirement when these two objectives are in conflict. [16] proposed an iterative algorithm using sum-of-squares (SOS) technique to search for the most permissive barrier function that gives maximum volume for the certified region, therefore maximize the estimate for safe stabilization region. All trajectories that start within the safe stabilization region can be made to converge to the (equilibrium of interest at the) origin as well as constrained in the safe region.

The main contribution of this paper is in applying barrier functions based method to certify safety of an inverter-based microgrid considering transient voltage limits. We propose a distributed safety certification method and present computational algorithms to compute safety-ensuring decentralized and distributed control policies. To treat the control design problem in a decentralized perspective, the microgrid is firstly decomposed into several subsystems. Barrier functions are generated for each subsystem by firstly ignoring the interactions from neighboring subsystems. The interaction terms are considered as disturbances with upper limits in the control design phase, resulting in robust local state feedback control strategies. The rest of the article is organized as follows: Section II presents the necessary background; Section III explains the microgrid model; the main computational and algorithmic developments are described in Section IV, with numerical results presented in Section V. We conclude the article in Section VI. Throughout the text, we will use $\left|\,\cdot\,\right|$ to denote both the Euclidean norm of a vector and the absolute value of a scalar; and use $\mathbb{R}\left[x\right]$ to denote the ring of all polynomials in $x\in\mathbb{R}^{n}$.

II Background

II-A Stability Certificates: Lyapunov Functions

Consider a nonlinear dynamical system of the form

[TABLE]

with an equilibrium at the origin $(f(0)=0)$, where $f:\mathbb{R}^{n}\rightarrow\mathbb{R}^{n}$ is locally Lipschitz. For brevity, we would drop the argument $t$ from the state variables, whenever obvious. The equilibrium point at the origin of (1) is Lyapunov stable if, for every $\varepsilon\!>\!0$ there is a $\delta\!>\!0$ such that $\left|x(t)\right|\!<\!\varepsilon~{}\forall t\!\geq\!0$ whenever $\left|x(0)\right|\!<\!\delta\,.$ Moreover, it is asymptotically stable in a domain $\mathcal{X}\!\subseteq\!\mathbb{R}^{n},\,0\!\in\!\mathcal{X},$ if it is Lyapunov stable and $\lim_{t\rightarrow\infty}\left|x(t)\right|\!=\!0\,$ for every $x(0)\!\in\!\mathcal{X}$ .

Theorem 1

[17, 18]** If there is a continuously differentiable radially unbounded positive definite function $V\!:\!\mathcal{X}\!\rightarrow\!\mathbb{R}_{\geq 0}$ such that $\nabla_{x}V^{T}\!f(x)$ is negative definite in $\mathcal{X}$, then the origin of (1) is asymptotically stable and $V(x)$ is a Lyapunov function.

Here $\nabla_{x}$ denotes the partial differentiation with respect to $x$ . Using an appropriately scaled Lyapunov function $V(x)$ , the region-of-attraction (ROA) of the origin of (1) can be estimated by $\left\{x\in\mathcal{X}\left|{V}(x)\leq 1\right.\right\}$ [19, 20].

II-B Safety Certificates: Barrier Functions

In contrast to asymptotic stability which concerns with the convergence of the state variables to the stable equilibrium, the notion of ‘safety’ comes from engineering design specifications. From the design perspective, the system trajectories are not supposed to cross into the certain regions in the state-space marked as ‘unsafe’. Let us assume that the ‘unsafe’ region of operation for the system (1) is given by the domain

[TABLE]

where $w_{i}:\mathbb{R}^{n}\mapsto\mathbb{R}$ are a set of $l$ ($\geq 1$) polynomials. These are usually engineering constraints that ensure that the system is always operated (controlled) to avoid going into ‘unsafe’ modes of operation. Safety of such systems can be verified through the existence (or, construction) of continuously differentiable barrier functions $B:\mathbb{R}^{n}\mapsto\mathbb{R}$ of the form [11, 12, 15, 16]:

[TABLE]

where $\alpha(\cdot)$ is an extended class-$\mathcal{K}$ function111A continuous function $\alpha:(-a,b)\mapsto(-\infty,\infty)$ , for some $a,b>0\,,$ is extendend class-$\mathcal{K}$ if it is strictly increasing and $\alpha(0)=0$ [18].. The third condition ensures that at the level-set $B=0$ the value of the barrier function is increasing along the system trajectories. Safety is guaranteed for all trajectories starting inside the domain $\{x\left|\,B(x)\geq 0\right.\}$ which is invariant under the dynamics (1).

II-C Sum-of-Squares Optimization

Relatively recent studies have explored how SOS-based methods can be utilized to find Lyapunov functions by restricting the search space to SOS polynomials [21, 22, 23, 20]. Let us denote by $\mathbb{R}\left[x\right]$ the ring of all polynomials in $x\in\mathbb{R}^{n}$. A multivariate polynomial $p\in\mathbb{R}\left[x\right],~{}x\in\mathbb{R}^{n}$, is an SOS if there exist some polynomial functions $h_{i}(x),i=1\ldots s$ such that $p(x)=\sum_{i=1}^{s}h_{i}^{2}(x)$. We denote the ring of all SOS polynomials in $x$ by $\Sigma[x]$. Whether or not a given polynomial is an SOS is a semi-definite problem which can be solved with SOSTOOLS, a MATLAB${}^{\text{\textregistered}}$ toolbox [24], along with a semi-definite programming solver such as SeDuMi [25]. An important result from algebraic geometry, called Putinar’s Positivstellensatz theorem [26, 27], helps in translating conditions such as in (3) into SOS feasibility problems.

Theorem 2

Let $\mathcal{K}\!\!=\!\left\{x\in\mathbb{R}^{n}\left|\,k_{1}(x)\geq 0\,,\dots,k_{m}(x)\geq 0\!\right.\right\}$ be a compact set, where $k_{j}$ are polynomials. Define $k_{0}=1\,.$ Suppose there exists a $\mu\!\in\!\left\{{\sum}_{j=0}^{m}\sigma_{j}k_{j}\left|\,\sigma_{j}\!\in\!\Sigma[x]\,\forall j\right.\right\}$ such that $\left\{\left.x\in\mathbb{R}^{n}\right|\,\mu(x)\geq 0\right\}$ is compact. Then,

[TABLE]

Remark 1

Using Theorem 2, one can translate the problem of checking that $p\!>\!0$ on $\mathcal{K}$ into an SOS feasibility problem where we seek the SOS polynomials $\sigma_{0}\,,\,\sigma_{j}\,\forall j$ such that $p\!-\!\sum_{j}\sigma_{j}k_{j}$ is SOS. Note that any equality constraint $k_{i}(x)\!=\!0$ can be expressed as two inequalities $k_{i}(x)\!\geq 0$ and $k_{i}(x)\!\leq\!0$. In many cases, especially for the $k_{i}\,\forall i$ used throughout this work, a $\mu$ satisfying the conditions in Theorem 2 is guaranteed to exist (see [27]), and need not be searched for.

III Microgrid Model

We consider the following model of droop-controlled inverter dynamics [28, 4]:

[TABLE]

where $\lambda_{i}^{p}>0$ and $\lambda_{i}^{q}>0$ are the droop-coefficients associated with the active power vs. frequency and the reactive power vs. voltage droop curves, respectively; $\tau_{i}$ is the time-constant of a low-pass filter used for the active and reactive power measurements; $\theta_{i}\,,\,\omega_{i}$ and $v_{i}$ are, respectively, the phase angle, speed and voltage magnitude; $v^{0}_{i}$ is the desired (nomial) voltage magnitude; $P_{i}^{\text{set}}$ and $Q^{\text{set}}_{i}$ are the active power and reactive power set-points, respectively. Finally, $P_{i}$ and $Q_{i}$ are, respectively, the active and reactive power injected into the network which relate to the neighboring bus voltage phase angle and magnitudes as:

[TABLE]

where $\theta_{i,k}=\theta_{i}-\theta_{k}$ , and $\mathcal{N}_{i}$ is the set of neighbor nodes. $G_{i,k}$ and $B_{i,k}$ are respectively the transfer conductance and susceptance values of the line connecting the nodes $i$ and $k$ .

At the equilibrium (steady-state) operation:

[TABLE]

As the conditions on the network change (such as changes in load or generation), inverters have the capability to change the control set-points of the active and reactive power output to adjust to the new operating conditions. This is modeled as:

[TABLE]

where $P^{0}_{i}$ and $Q^{0}_{i}$ are the set-points for the unperturbed (or nominal) operating condition; and $u^{p}_{i}$ and $u^{q}_{i}$ are some feedback control inputs.

Due to the low-inertia of the microgrids, large voltage and frequency fluctuations are quite common during transients [9]. While designing stabilizing control policies, it is therefore important to keep track of the transient voltage and frequency magnitudes to ensure that those are within the ‘safety’ limits determined via engineering design. In this work, we will restrict ourselves to the consideration of transient voltage limits which are usually higher than the steady-state operational limits [8]. Fluctuations of transient voltage limits beyond the tolerable (‘safe’) region cause power quality issues, including the risk of damaging the electrical equipment. In this paper, we will define the ‘safe’ operational region as:

[TABLE]

Typical values for the limits during transients operation could be $\underline{v_{i}}=0.6$ p.u. and $\overline{v_{i}}=1.2$ p.u. Other forms of safety constraints, such as frequency limits and power-flow limits will be considered in the future work.

IV Distributed Safety Certificates

The dynamical model of the interconnected microgrid with $m$ droop-controlled inverters is expressed compactly as:

[TABLE]

where each $i\in\{1,2,\dots,m\}$ identifies an inverter. $x_{i}\in\mathbb{R}^{n_{i}}$ is the $n_{i}$-dimensional state vector associated with the $i$-th inverter, while $u_{i}$ are some control inputs. We assume that the origin is an equilibrium point of interest of the networked system, and that the control input vanishes at the equilibrium point (i.e. $u_{i}=0\,\forall i$ at the origin). We also assume that $h_{ij}(x_{i},0)=0$ for all $x_{i}$ . Moreover $f_{i},\,g_{i}$ and $h_{ij}$ are locally Lipschitz functions. The problem we are interested in is:

Problem 1

Identify continuous functions $B_{i}(x_{i})$ , feedback control policies $u_{i}$ and non-negative scalars $c_{i}$ such that

[TABLE]

where we define $\mathcal{D}_{i}[c_{i}]:=\{x_{i}\,|\,B_{i}(x_{i})\geq c_{i}\}~{}\forall i$ and $\partial\mathcal{D}_{i}[c_{i}]:=\{x_{i}\,|\,B_{i}(x_{i})=c_{i}\}$ as the boundary set of the domain $\mathcal{D}_{i}[c_{i}]$ .

Theorem 3

If there exist continuous functions $B_{i}(x_{i})$ , feedback control policies $u_{i}$ and non-negative scalars $c_{i}$ satisfying (8), then the safety of the interconnected system (7) is guaranteed for all $t\geq 0$ whenever $B_{i}(x_{i}(0))\geq c_{i}\,\forall i$ , i.e.

[TABLE]

Moreover there is a neighborhood $\mathcal{X}_{i}$ around origin (i.e. $0\in\mathcal{X}_{i}\,\forall i$) such that $\mathcal{X}_{i}\subseteq\mathcal{D}_{i}[c_{i}]$ .

Proof:

Note that because of the condition (8c), $B_{i}$ is non-decreasing on the boundary of the domain $\mathcal{D}_{i}[c_{i}]$ whenever $x_{j}\in\mathcal{D}_{j}[c_{j}]$ for every neighbor $j$ . Extending this argument to all the subsystems, we conclude that

[TABLE]

is an invariant domain. Since $B_{i}(x_{i})<0$ for every $x_{i}\in\mathcal{X}_{u,i}$, we conclude the safety of the system is guaranteed for all $t\geq 0$ whenever $x_{i}(0)\in\mathcal{D}_{i}[c_{i}]~{}\forall i$ . Finally, since $B_{i}(0)>c_{i}$ and $B_{i}$ is a continuous function there exists a neighborhood $\mathcal{X}_{i}$ around origin such that for all $x_{i}\in\mathcal{X}_{i}$ , $B_{i}(x_{i})\geq c_{i}$ . ∎∎

Computation of such barrier functions is not trivial. Recent works have explored the use of sum-of-square optimization methods to compute the barrier certificates for polynomials networks [11, 16] . Note that the power-flows as described in (5) are non-polynomial. Using the polynomial recasting technique proposed in [20] , the power systems dynamics can be expressed in a higher-dimensional space as a polynomial differential-algebraic system. In this work, however, we resort to Taylor series expansion (up to third order) to approximate the dynamics into a polynomial form.

In the rest of this section, we describe a three-step procedure to obtain the distributed barrier certificates. In the first two steps, we consider the isolated and autonomous sub-system model of the form (which we assume to be locally asymptotically stable around the origin):

[TABLE]

and compute the Lyapunov function which is then used to compute a barrier function for the isolated sub-system using the method similar to [16].

IV-A Computation of Lyapunov Functions

SOS-based expanding interior algorithm [21, 20] has been used to construct Lyapunov functions and region-of-attraction in an iterative search process. In this work, we use a variant of the process - which does not require the bisection search process and hence speeds up the computation at each iteration stage. The algorithmic steps used to implement the modified expanding interior algorithm can be summarized as follows (for notational convenience, we have dropped the sub-script $i$ from the subsystem variables to explain the algorithm):

Step 0: Compute a Lyapunov function $V^{0}$ such that $V^{0}\leq 1$ is an estimate of the ROA. The following two steps are then repeated until convergence, such that (hopefully) the final estimate of the ROA is much larger than initial one. Define a positive definite and radially unbounded function $p(x)$ (e.g. $p=\varepsilon|x|^{2}$ for some small $\varepsilon>0$). 2. 2.

Step k-1: Starting from a Lyapunov function $\hat{V}$ with ROA estimated by $\hat{V}\leq 1$ , compute the largest level-set $\beta^{k}$ of a positive definite function $p(x)$ contained within $\hat{V}\leq 1$ . This is done by solving the following SOS problem:

[TABLE]

where $s$ are SOS polynomials. The $\hat{\cdot}$ implies it is borrowed from the previous step, while $k$ denotes the variables being currently computed. At the first instance of the problem (9), we initialize $\hat{s_{1}}=1$ . 3. 3.

Step k-2: In this sub-step at the k-th iteration, a new Lyapunov function $V^{k}$ is found such that the level-set $V^{k}=1$ is an estimate of the ROA, while trying to expand the estimated ROA by maximizing $\delta$ such that $p\leq\hat{\beta}$ is contained within the level-set $V^{k}\leq 1-\delta^{k}$, i.e.

[TABLE]

The algorithm stops when $\delta^{k}$ is sufficiently small. Set $V=V^{k}$ as the Lyapunov function with $V\leq 1$ providing the largest estimate of the ROA .

IV-B Computation of Barrier Functions

For the barrier functions computation we adopt a similar approach as in the Algorithm 2 in [16], except that we use the algorithm to compute only the barrier functions, while the original algorithm was used to also search for a ‘safe’ and stabilizing control policy. For completeness we present the algorithm here (once more, for notational convenience, we have dropped the sub-script $i$ from the subsystem variables to explain the algorithm):

Step 0: As the first step of the iterative process, we compute the maximum level-set of $V$ contained completely inside the safe region. This is done by solving the following SOS problem:

[TABLE]

Set $B^{0}=z^{\max}-V$ , where $z^{\max}$ is the solution of the above problem, i.e. the maximal level-set of $V$ wholly contained inside the safe region. Note that $B^{0}$ is a barrier function by construction. Choose a small scalar $\gamma\!>\!0$ . 2. 2.

Step k-1: Using the barrier function $\hat{B}$ computed in the previous step, find the largest $\varepsilon\!>\!0$ such that $\dot{\hat{B}}\geq-\gamma\hat{B}+\varepsilon$ whenever $\hat{B}\geq 0$ , i.e. solve the SOS problem

[TABLE] 3. 3.

Step k-2: In this sub-step we search for a new barrier function of the form $B^{k}(x)=z(x)^{T}Q^{k}z(x)$ where $z(x)$ is a vector of monomials in $x$ , and $Q^{k}$ is a symmetric matrix, such that $B^{k}(0)>0$ . The barrier function satisfies $B^{k}(x)<0$ on the unsafe set $\{x\,|\,w_{i}(x)>0\,,\,i=1,2,\dots,l\}$ , along with the constraint on its time-derivative. The following problem is solved:

[TABLE]

where $\eta$ is a small positive number chosen to avoid the trivial zero solution. The objective function is a proxy for maximizing the volume of the safety region [16].

The algorithm stops when $\text{trace}(Q^{k})$ converges within some tolerance. Set $B=B^{k}$ as the barrier function with $B\geq 0$ providing the largest estimate of the certified safety region for the isolated subsystems.

IV-C Safety Certifying Control Policies

In this subsection, we describe an SOS problem to compute control policies $u_{i}$ such that (8c) is satisfied for some $c_{i}\in[0,B_{i}(0))$ . Without any loss of generality, we will assume that the barrier functions satisfy $B_{i}(0)=1$ (always achievable through scaling), such that we are interested in $c_{i}\in[0,1)$ . Moreover, we will assume, for simplicity, uniform $c_{i}=c~{}\forall i$ , while the more generic case can be easily extended. Then we are seeking the existence of control laws $u_{i}$ such that for some chosen $c\in[0,1)$

[TABLE]

In this paper, we will focus on state-feedback control policies. Two alternatives will be considered: 1) a decentralized state-feedback policy of the form $u_{i}(x_{i})$ , and 2) a distributed state-feedback policy of the form $u_{i}=u_{ii}(x_{i})+{\sum}_{j\in\mathcal{N}_{i}}u_{ij}(x_{j})$ .

The following problem concerns the design of an optimal decentralized state-feedback control policy $u_{i}(x_{i})$:

[TABLE]

Note that the controller is only used on or near the boundary of the domain $\mathcal{D}_{i}[c]$ since it is only needed to guarantee that the trajectories never cross the boundary. This can be solved using an equivalent SOS problem, noting that the constraint $\|u_{i}(x_{i})\|_{\infty}\leq U_{i}$ translates to polynomial constraints. Similar problem can be formulated for the distributed controller design, with the constraint $\|u_{i}(x_{i})\|_{\infty}\leq U_{i}$ needed to be satisfied on $x_{j}\in\mathcal{D}_{j}[c]\,\forall j\in\mathcal{N}_{i}$ as well as $x_{i}\in\mathcal{D}_{i}[c]\,$.

Remark 2

Note that the constraint (11) is satisfied whenever (sufficient condition) we choose a $u_{i}$ such that

[TABLE]

If for every $x_{i}\in\partial{D}_{i}[c]$ there always exists a $k$ such that $\left|\left[\nabla B_{i}^{T}\right]_{k}\left[u_{i}\right]_{k}\right|>0$ then we can always find a $u_{i}$ satisfying the above condition.

V Numerical Results

For illustration purpose, we consider the 6-bus (bus 0 to bus 5) microgrid network described in [29] . Disconnecting the utility, we replace the substation (bus 0) by a droop-controlled inverter, with three other inverters placed on buses 1, 4 and 5 . The inverter dynamics were modeled in the form of (4). Bus 0 was considered as the reference bus for the phase angle. The droop coefficients $\lambda^{p}_{i}$ and $\lambda^{q}_{i}$ were chosen to be $2.43\,$rad/s/p.u. and $0.2$ p.u./p.u., while the filter time-constant $\tau_{i}$ was set to 0.5 s [4]. Nominal values of voltage and frequency, as well as the active and reactive power set-points were obtained by solving the steady-state power-flow equations (5), which were then used to shift the equilibrium point to the origin. The loads were modeled as constant power loads, and Kron reduced network with only the inverter nodes were used for analysis. The unsafe set was defined in terms of the shifted (around the 1 p.u.) nodal voltage magnitudes as follows:

[TABLE]

In Fig. 1, we illustrate how the iterative search algorithm presented in Section IV obtains an expanded certified region of safety (marked by the boundary of the outermost green ellipse) starting from the initial estimate given by a level-set of the Lyapunov function (marked by the smaller black dashed ellipse boundary). The plot shows the projections of the ROA and the barrier certified regions on the frequency-voltage space obtained by setting the phase angle differences to 0. The black dash marked boundary of larger ellipse is the estimate of the ROA, while the red dashed lines denote the unsafe region boundary (in voltage magnitudes). Note that the certified invariant region of safe stability is much smaller than the estimated ROA of the isolated inverter.

Next we investigate the control efforts needed to guarantee safety of the network over some domain defined using the subsystem barrier function level-sets. Fig. 2 shows the results of the optimal decentralized control design problem (12), for a range of different values of the barrier level-set $c\in[0,1)$ such that $B_{i}\geq c~{}\forall i$ gives a distributed certificate of safety. As expected, the control effort increases (monotonically, in this case) as the value of $c$ decreases, or the certified region of safety increases. Figs. 2-2 show the results when local neighboring subsystem state measurements are used in the control design in addition to the subsystem’s states, in what we call as the ‘distributed control’ design. Clearly, distributed control policies require lower minimum control efforts as compared to the decentralized control policy. This observation aligns with the conclusion in [30] regarding the value of communication in distribution network voltage regulation problem. In particular, two different choices of distributed controllers are explored - one in which all of the neighboring subsystems’ states are used in the feedback (Fig. 2) and another in which the only neighboring subsystem states used as feedback are the voltage magnitudes (Fig. 2). In this example, additional measurements from the neighboring subsystems help decrease the minimum control effort needed.

VI Conclusion

In this paper we consider the problem of safety in inverter-based microgrids. Using the barrier functions based methods, we propose a distributed safety certification of the mirogrid network. Sum-of-squares based algorithm was used to present a computational approach to obtain these safety certificates in a distributed manner. Moreover, using a microgrid example, we show how decentralized vs. distributed control policies could pose different requirements on the control effort. Future work will explore the extension of such methods to larger power systems networks, under the presence of uncertainties.

Acknowledgment

This work was carried out at PNNL (contract DE-AC05-76RL01830) under the support from the U.S. Department of Energy as part of their Grid Modernization Initiative.