Dynamics Concentration of Large-Scale Tightly-Connected Networks

TL;DR

This paper demonstrates that in large, tightly-connected linear networks with random transfer functions, nodes exhibit nearly identical responses to inputs, a phenomenon termed dynamics concentration, with implications for model reduction and robustness.

Contribution

It introduces the concept of dynamics concentration in large networks and provides a probabilistic analysis showing uniform response convergence regardless of input source.

Findings

Nodes in large networks follow the same response despite heterogeneity.

Network transfer matrix converges to a unique dynamic response.

Similar phenomena observed in small deterministic networks within certain frequency bands.

Abstract

The ability to achieve coordinated behavior --engineered or emergent-- on networked systems has attracted widespread interest over several fields. This has led to remarkable advances on the development of a theoretical understanding of the conditions under which agents within a network can reach agreement (consensus) or develop coordinated behaviors such as synchronization. However, fewer advances have been made toward explaining another commonly observed phenomena in tightly-connected networks systems: output responses of nodes in the networks are almost identical to each other despite heterogeneity in their individual dynamics. In this paper, we leverage tools from high-dimensional probability to provide an initial answer to this phenomena. More precisely, we show that for linear networks of nodal random transfer functions, as the network size and connectivity grows, every node in the network follows the same response to an input or disturbance --irrespectively of the source of this input. We term this behavior as dynamics concentration since it stems from the fact that the network transfer matrix uniformly converges in probability, i.e., it concentrates, to a unique dynamic response determined by the distribution of the random transfer function of each node. We further discuss the implications of our analysis in the context of model reduction and robustness, and provide numerical evidence that similar phenomena occur in small deterministic networks over a properly defined frequency band.

Equations

$$\mathbb{P}(|X|\ge t)\le 2\exp \left(-c{t}^2\right)$$

$$\mathbb{P}\left(\left|\sum _{i=1}^n{a}_i{X}_i\right|\ge t\right)\le 2\exp \left(-\frac{c{t}^2}{\Vert a{\Vert }^2}\right)$$

$$\mathbb{P}\left(\left|\sum _{i=1}^n{a}_i(X_i+j{Y}_i)\right|\ge t\right)\le 8\exp \left(-\frac{c{t}^2}{\Vert a{\Vert }^2}\right)$$

$$T(s)$$

$$=(I_n+\mathrm{diag}\{g_i(s)\}L(s){)}^{-1}\mathrm{diag}\{g_i(s)\}$$

$$L(s)=V(s)\Lambda (s)V^{\ast }(s),$$

$$T(s)$$

$$=(\mathrm{diag}\{g_i^{-1}(s)\}+L(s){)}^{-1}$$

$$=(\mathrm{diag}\{g_i^{-1}(s)\}+V(s)\Lambda (s)V^{\ast }(s){)}^{-1}$$

$$=V(s)(V^{\ast }(s)\mathrm{diag}\{g_i^{-1}(s)\}V(s)+\Lambda (s){)}^{-1}{V}^{\ast }(s)$$

$$\bar{g}(s):$$

$$=[\mathbb{E}Re(g_i^{-1}(s))+j\mathbb{E}Im(g_i^{-1}(s)){]}^{-1}$$

$$\underset{n\to \infty }{\lim }\mathbb{P}\left(\underset{s\in S}{\sup }\Vert T(s)-\frac{1}{n}\bar{g}(s)11^T\Vert \ge ϵ\right)=0$$

$$\underset{s\in S}{\sup }\Vert T(s)-\frac{1}{n}\bar{g}(s)11^T\Vert \stackrel{\mathcal{P}}{\to }0\text{as}n\to \infty$$

$$L=V\Lambda V^{\ast },\Lambda =\mathrm{diag}\{{\lambda }_i\},$$

$$T(s_0)$$

$$=V(V^{\ast }\mathrm{diag}\{X_i+j{Y}_i\}V+\Lambda {)}^{-1}{V}^{\ast }:=T$$

$$\Vert T-\frac{1}{n}{\mu }^{-1}11^T\Vert \stackrel{\mathcal{P}}{\to }0\text{as}n\to \infty ,$$

$$H=V^{\ast }\mathrm{diag}\{X_i+j{Y}_i\}V+\Lambda$$

$$\Vert T-\frac{1}{n}{\mu }^{-1}11^T\Vert$$

$$=\Vert V\left(H^{-1}-{\mu }^{-1}{e}_1{e}_1^T\right)V^{\ast }\Vert$$

$$=\Vert H^{-1}-{\mu }^{-1}{e}_1{e}_1^T\Vert$$

$$H=V^{\ast }\mathrm{diag}\{{\tilde{X}}_i+j{\tilde{Y}}_i\}V+\mathrm{diag}\{{\lambda }_i+\mu \}$$

$$\mathbb{P}\left(|d_{11}|\ge ϵ\right)\le 8\exp \left(-c_{1}n{ϵ}^2\right)$$

$$\mathbb{P}\left(\sum _{k=1}^n|d_{k1}{|}^2\ge ϵn^p\right)\le \frac{c_2}{ϵn^p}$$

$$\mathbb{P}\left(\sum _{l=1}^n|d_{1l}{|}^2\ge ϵn^p\right)\le \frac{c_2}{ϵn^p}$$

$$\mathbb{P}\left(\Vert D\Vert \ge ϵn^p\right)\le \frac{c_3\sqrt{\log n}}{ϵn^p}$$

$$D_{(1)(1)}=\left[\begin{array}{ccc}{d}_{22}& \cdots & d_{2n}\\ ⋮& \ddots & ⋮\\ d_{n2}& \cdots & d_{nn}\end{array}\right],{\Lambda }_{(1)}=\mathrm{diag}\{{\lambda }_i{\}}_{i=2}^n$$

$$H=\left[\begin{array}[]{c|c}h_{11}&h_{12}^{T}\\ \hline\cr h_{21}&H_{22}\end{array}\right],\ \begin{array}[]{l}h_{11}=d_{11}+\mu\\ h_{12}=[d_{12},\cdots,d_{1n}]^{T}\\ h_{21}=[d_{21},\cdots,d_{n1}]^{T}\\ H_{22}=D_{(1)(1)}+\Lambda_{(1)}+\mu I_{n-1}\end{array}$$

$$H^{-1}=\left[\begin{array}[]{c|c}\frac{1}{a}&-\frac{1}{a}h_{12}^{T}H_{22}^{-1}\\ \hline\cr-\frac{1}{a}H_{22}^{-1}h_{21}&H_{22}^{-1}+\frac{1}{a}H_{22}^{-1}h_{21}h_{12}^{T}H_{22}^{-1}\end{array}\right]$$

$$\mathbb{P}\left(\Vert H^{-1}-{\mu }^{-1}{e}_1{e}_1^T\Vert \ge ϵ\right)$$

$$\le$$

$$+12\mathbb{P}\left(\sum _{k=1}^n|d_{k1}{|}^2\ge C_3^2{n}^p\right)+12\mathbb{P}\left(\sum _{l=1}^n|d_{1l}{|}^2\ge C_3^2{n}^p\right)$$

$$T=V(V^{\ast }\mathrm{diag}\{X_i+j{Y}_i\}V+\Lambda {)}^{-1}{V}^{\ast },$$

$$\underset{n\to +\infty }{\lim }\mathbb{P}\left(\Vert T-\frac{1}{n}{\mu }^{-1}11^T\Vert \ge ϵ\right)=0$$

$$\underset{n\to +\infty }{\lim }\mathbb{P}\left(\Vert T(s_0)-\frac{1}{n}\bar{g}(s_0)11^T\Vert \ge ϵ\right)=0$$

Full text

I Introduction

Coordinated behavior in network systems has been a popular subject of research in many fields, such as physics [1], chemistry [2], social sciences [3], and biology [4]. Within engineering, coordination is essential for proper operation of many networked systems including power networks [5, 6], data and sensor networks [7, 8], and autonomous transportation [9, 10, 11, 12]. While there exist many expressions of this behavior in network systems, two forms of coordination have particularly received thorough attention by the control community: Consensus and synchronization.

Consensus [3, 9, 10, 11, 13, 14], on one hand, refers to the ability of the network nodes to asymptotically reach a common value over some quantities of interest. Many extensions of this problem include the study of robustness and performance of consensus networks in the presence of noise [10, 11], time-delay [13, 14], and switching graph topology [14]. Synchronization [4, 7, 8, 12, 15, 16, 17], on the other hand, refers to the ability of network nodes to follow a commonly defined trajectory. Although for nonlinear systems synchronization is a structurally stable phenomenon, in the linear case [15, 12, 16, 17], synchronization requires the existence of a common internal model that acts as a virtual leader [16, 17].

A less studied phenomenon, that is empirically observed, is a coherent response within large-scale tightly-connected networks in which every node identically reacts to perturbations, irrespectively from the nature or location of the disturbance. For example, in tightly-connected power networks, generator dynamics across the network tend to react coherently to system disturbances [18, 19, 20]. However, while in the case of swing dynamics, reduced order models provide good approximations [18], for generators with heterogeneous turbine time constants, a good low order approximation is difficult to find [20]. Among other contributions, this work aims to explain this difference.

In this paper, we introduce a new framework to analyze the aggregated dynamics of large networks. We consider a network consisting of heterogeneous linear nodes interconnected through a weighted graph Laplacian matrix, with the node dynamics represented by random transfer functions. We show that, whenever the algebraic connectivity of the graph is polynomial in the network size $n$, the transfer matrix of the network converges in probability as $n$ grows to infinity to a common scalar transfer function spanning the consensus subspace. Notably, the resulting scalar transfer function is deterministic and can be determined by the harmonic expectation of the individual nodal dynamics. We term this behavior dynamics concentration, due to the mathematical principle that explains this phenomenon, i.e., concentration of measure [21].

The implications of our results are manifold. Firstly, it extends the notions of consensus and synchronization to scenarios in which coherent behavior can be achieved even in the presence of disturbances that are arbitrary in source or shape. Secondly, unlike output synchronization that requires the existence of a common internal model within each node, dynamics concentration can be achieved despite heterogeneity on the individual dynamics or lack of a common internal model. Thirdly, due to the stochastic nature of our analysis, many networks that a priori may look quite different in composition and topology, exhibit exactly the same behavior. Finally, the analysis further provides a principled methodology to compute the concentrated dynamics, which as we will show later may not always be represented by a reduced order model.

The rest of the paper is organized as follows. In Section II, we introduce some technical preliminaries such as the notion of sub-Gaussian random variables and their concentration inequalities and use this formalism to formulate our problem statement. In Section III, the condition for uniform stochastic convergence of the transfer matrix over compact sets is given. In Section IV, we provide an application of our analysis to the problem of characterizing reduced models for power networks. At last, we conclude this paper with more discussions on the implications of this result.

Notation: For a vector $x$, $\|x\|=\sqrt{x^{T}x}$ denotes the $2$-norm of $x$, and for a matrix $A$, $\|A\|$ denotes the spectral norm and \underaccent{\bar}{\sigma}(A) denotes the least singular value of $A$. We let $I_{n}$ denote the identity matrix of order $n$, $V^{*}$ denote the conjugate transpose of matrix $V$, $\mathbb{1}$ denote $[1,\cdots,1]^{T}$ with proper dimensions, and $[n]$ denote the set $\{1,2,\cdots,n\}$. Also, we write complex numbers as $a+jb$, where $j=\sqrt{-1}$.

For function ordering, we write $f(n)=o(g(n))$ if $\forall m>0$, $\exists n_{0}>0$ s.t. $|f(n)|<mg(n)$ for $n>n_{0}$. We write $f(n)=\Omega(g(n))$ if $\exists m>0$, $\exists n_{0}>0$ s.t. $f(n)\geq mg(n)$ for $n>n_{0}$.

II Preliminaries

II-A Sub-Gaussian random variables

Firstly, we let $(\Omega,\mathcal{F},\mathbb{P})$ be the probability space, where $\Omega$ is the sample space, $\mathcal{F}=2^{\Omega}$ is the event set and $\mathbb{P}$ is a probability measure on $\mathcal{F}$. Wherever random variables are introduced, we assume that those random variables are properly defined in this probability space.

Throughout the paper, we work on a special family of random variables: the sub-Gaussian random variables. One way to define such random variables is the following.

Definition** (Sub-Gaussian random variables).**

A random variable $X$ is a sub-Gaussian random variable if $\forall t>0$:

[TABLE]

for some $c>0$.

Such random variables exhibit an exponentially decaying tail probability, which gives good concentration results when summing over sub-Gaussian random variables.

A widely used bound on the tail probability of the sum is given by:

Lemma 2.1** (Hoeffding’s inequality[22]).**

Let $X_{1},\cdots,X_{n}$ be independent, mean zero, sub-gaussian random variables, and let $a=[a_{1},\cdots,a_{n}]^{T}\in\mathbb{R}^{n}$. Then $\forall t>0$, the following holds:

[TABLE]

for some $c>0$.

We also provide a direct application of Hoeffding’s inequality to the random complex numbers:

Lemma 2.2**.**

Let $(X_{i},Y_{i}),\ i=1,\cdots,n$ be i.i.d. samples from a joint distribution such that $X_{i},Y_{i}$ are sub-gaussian random variables, and let $a=[a_{1},\cdots,a_{n}]^{T}\in\mathbb{C}^{n}$. Then $\forall t>0$, the following holds:

[TABLE]

for some $c>0$.

The proof is shown in appendix.

II-B Problem Statement

Consider a network consisting of $n$ nodes, indexed by $i\in[n]$ with the block diagram structure in Fig.1. We use $G(s)=\mathrm{diag}\{g_{i}(s)\}$, with $g_{i}(s),\ i\in[n]$, to represent the dynamics of the nodes. $L(s)$ is a generalized dynamic Laplacian matrix of the graph describing the network interconnection. More precisely, $L(s)$ is defined as transfer matrix such that $\forall s_{0}\in\mathbb{C}$, $L(s_{0})$ is a complex normal matrix with simple eigenvalue $\lambda=0$. In most applications, $L(s)=L\in\mathbb{R}^{n\times n}$ for some fixed graph Laplacian $L$.

Many existing networks can be represented represented by this structure: For the first-order consensus[9], the node dynamics is given by $g_{i}(s)=\frac{1}{s}$. For transportation networks, $g_{i}(s)$ are the vehicle dynamics. The Laplacian $L(s)=L$ is static for these two cases. For power networks[6], $g_{i}(s)$ are the dynamics of the generators and $L(s)=\frac{1}{s}L_{B}$ where $L_{B}$ is a Laplacian matrix representing the sensitivity of power injection w.r.t. bus phase angles.

The transfer matrix from $u$ to $y$ is given by:

[TABLE]

Since L(s) is a normal complex matrix,

[TABLE]

where $V(s)V^{*}(s)=V^{*}(s)V(s)=I_{n}$, and $\Lambda(s)=\mathrm{diag}\{\lambda_{i}(s)\}$ with $0=\lambda_{1}(s)<|\lambda_{2}(s)|\leq\cdots\leq|\lambda_{n}(s)|$. Using (4), we rewrite $T(s)$ as

[TABLE]

When considering the convergence of the transfer matrix of networked dynamical systems, as the network size $n$ grows, it is intuitive to assume the node dynamics $g_{i}(s)$ are generated randomly according to a particular distribution and discuss the convergence in probability. In practice, the assumed distribution is modeled from physical parameters of the nodes. For example, we may assume that nodes dynamics are given by random rational transfer functions whose coefficients are random variables.

Remark*.*

For networks with deterministic nodes, notice that we can always assume the deterministic dynamics are random samples drawn from some unknown distributions, hence the empirical distribution is a reasonable approximation given some mild requirements on the statistics of the samples.

We particularly focus on the convergence of $T(s)$ over a compact region $S\subset\mathbb{C}$ where for $s_{0}\in S$, both the real and imaginary part of $g^{-1}_{i}(s_{0})$ are sub-Gaussian random variables. Let

[TABLE]

We will show that for networks with algebraic connectivity $\inf_{s\in S}|\lambda_{2}(s)|=\Omega(n^{p})$ for some $p\in(0,1]$, as the network size $n$ increases, $T(s)$ converges uniformly to $\frac{1}{n}\bar{g}(s)\mathbb{1}\mathbb{1}^{T}$ in probability over $S$, i.e., $\forall\epsilon>0$

[TABLE]

or using the following notation:

[TABLE]

Thus, since the output response of a stable system can be obtained by inverse Laplace transform on the imaginary axis, for low frequency disturbances or input signals, having uniform convergence of $T(s)$ within a low frequency band on the imaginary axis, implies that, with high probability, the network output response is close to the output response of $\frac{1}{n}\bar{g}(s)\mathbb{1}\mathbb{1}^{T}$, i.e., the output responses of nodes in the networks are almost identical to each other, despite of the heterogeneity in their individual dynamics.

III Dynamics concentration of large-scale networks

We term the uniform convergence of $T(s)$ over low frequency bands as dynamics concentration because nodes with heterogeneous dynamics exhibit identical responses to low frequency disturbance, irrespectively of the source or shape of the disturbance. In this section, we firstly show point-wise convergence of $T(s)$ given the algebraic connectivity of the network is sufficiently large. Then we provide conditions for uniform convergence of $T(s)$ over a compact subset of the complex plane.

III-A Point-wise convergence in probability

For convenience, we assume at $s=s_{0}$, $g_{i}^{-1}(s_{0})=X_{i}+jY_{i}$ where the $\{(X_{i},Y_{i})\}_{i=1}^{n}$ pairs are independent and $X_{i},Y_{i}$ are sub-Gaussian random variables. Moreover, we let $L(s_{0})=L$, with unitary decomposition:

[TABLE]

Then we have:

[TABLE]

We also define $\mu:=\mathbb{E}X_{i}+j\mathbb{E}Y_{i}=\bar{g}^{-1}(s_{0})$. Then the point-wise convergence in probability of $T(s)$ at $s=s_{0}$ is equivalent to

[TABLE]

which is simply the convergence of a random complex matrix in probability. To show convergence then, we firstly define

[TABLE]

Since the first eigenvector of $L$ is $v_{1}=\frac{1}{\sqrt{n}}\mathbb{1}$, it follows that

[TABLE]

where $e_{1}$ is the first column of identity matrix. The remaining is straightforward: we get an explicit form for $H^{-1}$ and bound the tail probability of $\left\lVert H^{-1}-\mu^{-1}e_{1}e_{1}^{T}\right\rVert$ by $o(1)$ terms.

To invert $H$, we let $\tilde{X}_{i}=X_{i}-\mathbb{E}X_{i}$, and $\tilde{Y}_{i}=Y_{i}-\mathbb{E}Y_{i}$, then rewrite $H$ as:

[TABLE]

The matrix $H$ can be regarded as a perturbed diagnal matrix $\mathrm{diag}\{\lambda_{i}+\mu\}$ with matrix $D:=V^{*}\mathrm{diag}\{\tilde{X}_{i}+j\tilde{Y}_{i}\}V$. The matrix $D$ is regarded as a perturbation in the sense of the following result:

Lemma 3.1**.**

*Let $(\tilde{X}_{i},\tilde{Y}_{i}),\ i=1,\cdots,n$ be i.i.d. samples from a joint distribution such that $\tilde{X}_{i},\tilde{Y}_{i}$ are mean zero, sub-gaussian random variables. Let $D=V^{*}\mathrm{diag}\{\tilde{X}_{i}+j\tilde{Y}_{i}\}V:=[d_{kl}]$, where $V=[v_{1},\cdots,v_{n}]^{T}$ is defined in (4).

The following inequalities hold:*

$\forall\epsilon>0$**

[TABLE]

for some $c_{1}>0$. 2. 2.

$\forall\epsilon>0$, $\forall p>0$

[TABLE]

for some $c_{2}>0$. 3. 3.

$\forall\epsilon>0$, $\forall p>0$

[TABLE]

for some $c_{3}>0$.

The proof is shown in the appendix. Lemma 3.1 suggests that the perturbation matrix $D$ exhibits good concentration on its first row and column. We will exploit this property with a proper form of $H^{-1}$.

Define the following matrix:

[TABLE]

by removing the first row and column of $D$ and $\Lambda$. We then write $H$ as:

[TABLE]

Then, the inverse of $H$ is given by:

[TABLE]

where $a=h_{11}-h_{12}^{T}H_{22}^{-1}h_{21}$.

With this explicit form of $H^{-1}$, we bound the tail probability of $\|H^{-1}-\mu e_{1}e_{1}^{T}\|$ using the following lemma.

Lemma 3.2**.**

Assume $|\lambda_{2}|=\Omega(n^{p})$ for some $p\in(0,1]$. For $H$ as defined in (13), whenever $\mu\neq 0$, given $n$ large enough, then $\forall\epsilon>0$, we have:

[TABLE]

for some $C_{1},C_{2},C_{3}>0$ that depends on $\epsilon$.

The proof is shown in the appendix. As suggested in Lemma 3.1, all terms on the right-hand side of (15) are $o(1)$. By (7) and taking limit $n\rightarrow\infty$ on both sides of (15), we can then show the convergence in probability of $T$.

Theorem 3.3**.**

Let $(X_{i},Y_{i}),\ i=1,\cdots,n$, be i.i.d. samples from a joint distribution such that $X_{i},Y_{i}$ are sub-gaussian random variables, and denote $\mu:=\mathbb{E}X_{i}+j\mathbb{E}Y_{i}$. Consider a Laplacian matrix $L=V\Lambda V^{*}$, as defined in (4), with algebraic connectivity satisfying $|\lambda_{2}|=\Omega(n^{p})$ for some $p\in(0,1]$. Then given

[TABLE]

whenever $\mu\neq 0$, we have, $\forall\epsilon>0$,

[TABLE]

We can now apply the theorem to show the point-wise convergence of $T(s)$ at certain $s_{0}$.

Corollary 3.4**.**

Let $T(s)$ and $\bar{g}(s)$ be defined as in (2) and (3), respectively. Given $s=s_{0}$, consider the (possibly complex) network graph Laplacian $L(s_{0})$ with algebraic connectivity satisfying $|\lambda_{2}(s_{0})|=\Omega(n^{p})$ for some $p\in(0,1]$. Suppose that $g_{i}^{-1}(s_{0})$ has both its real and imaginary part given by sub-gaussian random variables, and $s_{0}$ is not a pole of $\bar{g}(s)$. Then $\forall\epsilon>0$,

[TABLE]

III-B Uniform convergence in probability

For people who are familiar with real analysis, it is well-known that the uniform convergence of function sequence over a compact region not only requires point-wise convergence, but also needs the equicontinuity of the function sequence [23]. Similarly, for stochastic uniform convergence of $T(s)$, we have the following theorem:

Theorem 3.1**.**

Let $T(s)$ and $\bar{g}(s)$ be defined as in (2) and (3), respectively. Given $S$ a compact subset of $\mathbb{C}$, and $\bar{g}(s)$ is uniform continuous on $S$, then $T(s)$ uniformly converges in probability to $\frac{1}{n}\bar{g}(s)\mathbb{1}\mathbb{1}^{T}$, i.e. $\forall\epsilon>0$

[TABLE]

if the following hold:

(Point-wise convergence) $\forall s\in S$ and $\forall\epsilon>0$,

[TABLE] 2. 2.

(Stochastic equicontinuity) $\forall\epsilon>0$,

[TABLE]

We only provide a proof sketch here. For details, please refer to [24], where the uniform convergence of real-valued random function sequence is discussed.

Sketch of proof.

Since $S$ is compact, for every $\delta>0$, there exists a finite cover $\{x\in\mathbb{C}:|x-s_{i}|<\delta\}$, $i=1,\cdots,M$ with $s_{i}\in S,\ i=1,\cdots,M$.

One can show that, for every $\delta<0$, the following holds:

[TABLE]

For all $\eta>0$, choose $\delta$ small enough such that: 1) By uniform continuity of $\bar{g}(s)$, the second term is zero; 2) By stochastic equicontinuity of $T(s)$,

[TABLE]

Then for this $\delta$, choose $N>0$ such that for $n\geq N$, the first term is less than $\frac{\eta}{2}$, by (17), and the third term is less than $\frac{\eta}{2}$, by point-wise convergence. Therefore we bound the left-hand side of (16) by arbitrary small $\eta$ for large $n$, which shows the uniform convergence. ∎

The stochastic equicontinuity, as its name suggests, is the ”stochastic version” of equicontinuity in the real analysis. One would expect a ”stochastic version” of Lipschtz condition as sufficient to show stochastic equicontinuity:

Lemma 3.2**.**

Let $T(s)$ be defined as in (2). Suppose for some sequence of random variables $\{B_{n}\}$, $\|T(s_{1})-T(s_{2})\|\leq B_{n}|s_{1}-s_{2}|$ holds for all $s_{1},s_{2}\in S$. Then $T(s)$ is stochastic equicontinuous on $S$, if

[TABLE]

For details, please refer to [24]. Lastly, we provide conditions on $T(s)$, $g_{i}^{-1}(s)$ to ensure such $\{B_{n}\}$ exists, hence show the uniform convergence of $T(s)$:

Theorem 3.3**.**

Let $T(s)$ and $\bar{g}(s)$ be defined as in (2) and (3), respectively. Given $S$ a compact subset of $\mathbb{C}$, and $\bar{g}(s)$ is uniform continuous on $S$. Suppose for some sequence of random variables $\{\tilde{B}_{n}\}$, $\max_{i\in[n]}|g_{i}^{-1}(s_{1})-g_{i}^{-1}(s_{2})|\leq\tilde{B}_{n}|s_{1}-s_{2}|$ holds for all $s_{1},s_{2}\in S$. Then $T(s)$ uniformly converges in probability to $\frac{1}{n}\bar{g}(s)\mathbb{1}\mathbb{1}^{T}$ if the following hold:

$T(s)$* converges point-wise in probability to $\frac{1}{n}\bar{g}\mathbb{1}\mathbb{1}^{T}$ $\forall s\in S$.* 2. 2.

[TABLE]

Proof.

For all $s_{1},s_{2}\in S$, we have:

[TABLE]

by Lemma (.1), we have:

[TABLE]

Taking $\lim_{M\rightarrow\infty}\limsup_{n\rightarrow\infty}$ on both sides:

[TABLE]

by Lemma 3.2, $T(s)$ is stochastic equicontinuous on $S$. Therefore, by Theorem 3.1, uniform convergence of $T(s)$ on $S$ is proved. ∎

Remark*.*

For our purpose, the uniform convergence of $T(s)$ over a low frequency band $S=\{jw:w\in[-w_{0},w_{0}]\}$ is sufficient for a stable network to exhibit synchronized output under low frequency disturbance. And as one can expect, the choice of a large enough network size $n_{0}$ to synchronize the output is frequency dependent, i.e. if we want uniform convergence over a wider frequency band, $n_{0}$ needs to be larger.

IV Numerical verification

We now provide two numerical examples of networked dynamics that illustrate the uniform convergence of their transfer matrix over low frequency band by numerical simulation. We will further use our dynamics concentration results propose a reduced order model for power networks.

IV-A Dynamics concentration in consensus networks

Consider the network in Fig.1 with the dynamics of each node given by

[TABLE]

This gives us the standard continuous time consensus network[9]. Notice that when $g_{i}(s)$ is simulated by impulse input $\delta(t)$, it is equivalent to setting initial condition at $t=0$ of node $i$ to be $k_{i}$.

We sample $k_{i}$ from $Unif(1,5)$, then $\frac{s}{k_{i}}$ is sub-Guassian for any fixed $s$. The graph of the network is a $k$-regular ring, i.e. every node is connected to its $2k+1$ nearest neighbors. We set the weight to be $1$ for all edges and $k\approx 0.15n$, then each node is connected to roughly $1/3$ of other nodes in the network. It can be shown that the algebraic connectivity of the graph Laplacian $|\lambda_{2}|$ is $\Omega(n)$ [25], then this network should exhibit dynamics concentration, according to the convergence result in Section III. The expected dynamics is given by

[TABLE]

Suppose the network is subject to an impulse input $u(t)=\delta(t)\mathbb{1}$. Then, the impulse response of $\bar{g}(s)$ is given by $\bar{g}(t)=\frac{4}{\ln 5}\chi_{\geq 0}(t)$, where $\chi_{\geq 0}(t)$ is the unit step function.

We plot the impulse response for network with different size $n=20,50,100,500$ in Fig. 2, along with $\bar{g}(t)$ shown by the red dashed line.

It is clear in Fig. 2 that as the network size grows, the network exhibits a more coherent response and the consensus value gets close to the expected average $4/\ln 5$. Although such response can be explained using standard consensus algorithm analysis, it can also be viewed as the consequence of $T(s)$ getting elementwise close to the coherent dynamics $\bar{g}(s)=\frac{4}{n\ln 5}\frac{1}{s}$.

IV-B Dynamic concentration in power networks

We now look at the case of a power network and leverage our analysis to provide an accurate low-frequency reduced order model for power networks.

Consider the transfer matrix of power generator networks[6] linearized around its equilibrium point, with the following block diagrams:

Remark*.*

The integration $\frac{1}{s}$ in the feedback loop does not affect the convergence result around $s=0$ because when evaluating $T(s)$ at a particular point $s_{0}$ around [math], we notice that $\frac{1}{s_{0}}$ actually scales up the algebraic connectivity of $\frac{1}{s}L_{B}$, which makes the network dynamics concentrate to a single dynamic in the consensus subspace.

The generator dynamics $g_{i}(s)$ is given by,

[TABLE]

for the swing dynamics, where $m_{i}$ is the inertial and $d_{i}$ is the damping ratio.

For generators with turbine control, the transfer function is given by,

[TABLE]

where $\tau_{i}$ is the turbine constant, and $r_{i}$ is the droop coefficient. We use the Icelandic grid data available at [26] where 35 generators in total are connected to the grid, with only some of them implementing turbine control.

According to our convergence result, for low frequencies or tightly connected networks, the transfer matrix of the network is close to $\frac{1}{n}\tilde{g}(s)\mathbb{1}\mathbb{1}^{T}$, with

[TABLE]

i.e. we assume $g_{i}(s)$ are drawn from some unknown distribution and use $\tilde{g}(s)$ as an empirical approximation of $\bar{g}(s)$. Hence we expect $\tilde{g}(s)$ to be a good candidate for a reduced model of the whole network.

For generators without turbine control, the reduced model is simply a single generator with its coefficients given by the mean of coefficients among all nodes. However, notice that for generators with turbine control, the reduced model $\tilde{g}(s)$ is with the same order as a single generator only when the turbine constant $\tau_{i}$ are the same for all generators in the grid. Every time a different $\tau_{i}$ is introduced, the order of $\tilde{g}(s)$ will be increased by $1$.

In the Icelandic grid, the $\tau_{i}$ are not the same among generators with turbine control. We will compare the output of:

Reduced model $\frac{1}{n}\hat{g}(s)\mathbb{1}\mathbb{1}^{T}$, $\hat{g}(s)$ as a representative generator:

[TABLE]

with the coefficients:

[TABLE]

where $\mathcal{T}$ is the set of indices such that $\tau_{i}\neq 0$; 2. 2.

Higher order reduced model $\frac{1}{n}\tilde{g}(s)\mathbb{1}\mathbb{1}^{T}$, where

[TABLE]

which is 4th order for Icelandic grid.

In Fig.4, we show the output response of the power networks and our proposed reduced model under step $u(t)=-0.3e_{2}$ and sinusoidal $u(t)=0.2\sin(0.1\pi t)e_{2}$ disturbance. It is clear that the generators have highly synchronized response to the low frequency disturbance, which shows that even in relatively small networks, the dynamics concentration can be observed in low frequency range, mostly because the integration in the feedback loop scales up $|\lambda_{2}|$ of the Laplacian, forcing the dynamics of all nodes to concentrate. However, we also observed less synchronized behavior when the network is subject to a step disturbance. This is because at high frequency ranges the integration in the feedback loop scales down $|\lambda_{2}|$, making the convergence result fail. Consequently, generators have different responses to the high frequency components of the step signal.

Finally, when comparing the two proposed reduced models, the higher order reduced model is a much better approximation of the original model than a reduced model of a single generator. For networks of generators with turbine control, if $\tau_{i}$ varies among the network, then $\tilde{g}(s)$ will have as many poles as the size of the network! Then it is difficult to find a good low order reduced models of such networks as its size scales up.

V Conclusion

In this paper, we consider a tightly-connected network consisting of heterogeneous linear node dynamics represented by random transfer functions. We show that the transfer matrix of the network converges in probability as $n$ grows to infinity to a common deterministic scalar transfer function spanning the consensus subspace. We then provide conditions under which such convergence is uniform over a compact set, and numerically verified, for certain networks, that dynamics concentration do occur in the low frequency range.

There are many possible extension to the current analysis. Firstly, our results suggest that for tightly connected networks, a common controller might perform as good as individual controllers specified for each node, even when nodes are with different types of dynamics. Secondly, another interesting question is whether there are circumstances under which the stability of the concentrated dynamics $\bar{g}(s)$ imply the stability of the whole network, given large enough network size. It is our belief that this framework can potentially provide new control design tools for large-scale networks.

The following Lemmas will be used in the proof of Lemma 2.2, Lemma 3.2 and Theorem 3.3

Lemma .1**.**

Let $X_{1},\cdots,X_{M}$ be non-negative random variables, then for $\epsilon>0$, we have:

[TABLE]

Proof.

[TABLE]

Similarly, since $X_{i}\geq 0$

[TABLE]

∎

Lemma .2**.**

Let $X,Y,Z,W_{1},W_{2}$ be non-negative random variables. Suppose whenever $Y\neq 0$, $\frac{1}{Y}\geq W_{1}-W_{2}$ holds, then we have:

[TABLE]

Proof.

when $W_{2}<W_{1}$, by $\frac{1}{Y}\geq W_{1}-W_{2}$, we obtain:

[TABLE]

Notice that the inequality above also holds for $Y=0$ as long as $W_{2}<W_{1}$, hence,

[TABLE]

∎

Proof of Lemma 2.2.

Firstly, we show that $\forall t>0$

[TABLE]

Similarly, for some $c_{y}>0$:

[TABLE]

Then we have:

[TABLE]

∎

Proof of Lemma 3.1.

Firstly, for every $k\in[n]$ we have $d_{k1}=\sum_{i=1}^{n}\bar{v}_{ki}v_{1i}(\tilde{X}_{i}+j\tilde{Y}_{i})$, by Lemma 2.2, for $\epsilon>0$ we have

[TABLE]

for some $c_{1}>0$. Notice that $v_{1}=\mathbb{1}/\sqrt{n}$, then

[TABLE]

Then we have

[TABLE]

which proves the firstly inequality. Similarly, for $d_{1l}=\sum_{i=1}^{n}\bar{v}_{1i}v_{li}(\tilde{X}_{i}+j\tilde{Y}_{i})$, $l\in[n]$, we have:

[TABLE]

With the inequality above, the 2nd moment of $d_{k1}$ is bounded by:

[TABLE]

From this, by Markov’s inequality, we have:

[TABLE]

which gives inequality (10a), similarly for inequality (10b).

Now for the last inequality, since unitary matrix preserves the spectral norm, it is straightforward to see that

[TABLE]

We show that $|\tilde{X}_{i}+j\tilde{Y}_{i}|$ is also sub-gaussian by noticing that:

[TABLE]

and $|X_{i}|+|Y_{i}|$ is a sum of two sub-gaussian random variables, hence sub-gaussian [22, Proposition 2.6.1]. For sequence of sub-gaussian random variables, its expected maximum is bounded by [22, Excercise 2.5.10]:

[TABLE]

for some $c_{3}>0$. Then we have:

[TABLE]

∎

Proof of Lemma 3.2.

Firstly, notice that $\|H^{-1}-\mu^{-1}e_{1}e_{1}^{T}\|$ is bounded by the summation of the spectral norm of every blocks of $H^{-1}$:

[TABLE]

Then by Lemma (.1), we have:

[TABLE]

The remaining proof is to bound the tail probabilities above.

Notice that \|H_{22}^{-1}\|=\frac{1}{\underaccent{\bar}{\sigma}(H_{22})} By Weyl’s inequality [27, Theorem III.2.1]:

[TABLE]

Then we know that:

[TABLE]

Since $|\lambda_{2}|=\Omega(n^{p})$, then $\exists m_{1}>0$ s.t. for some $n_{1}>0$, whenever $n\geq n_{1}$, $|\lambda_{2}|-|\mu|\geq m_{1}n^{p}$ holds, which leads to:

[TABLE]

For the bound on tail probability of $\|H_{22}^{-1}\|$, by Lemma .2:

[TABLE]

Because $\exists 0<m_{2}<m_{1}$ s.t. when $n\geq n_{2}$, $m_{1}n^{p}-\frac{5}{\epsilon}\geq m_{2}n^{p}$ holds, we let $n\geq\max\{n_{1},n_{2}\}$, to get:

[TABLE]

Now we turn to bound the relatively complex terms:

[TABLE]

where $\epsilon_{1}=\min\{\frac{\epsilon|\mu|}{10},\sqrt{\frac{\epsilon}{10}}\}$.

Similarly, we have:

[TABLE]

and

[TABLE]

Then we continue to reduce the tail probability of $\left\lvert\frac{1}{a}-\frac{1}{\mu}\right\rvert$, notice that

[TABLE]

then we let $\epsilon_{2}=\min\{\sqrt{\frac{\epsilon}{10}},\frac{\epsilon}{5}\}$

[TABLE]

where $\epsilon_{3}=\frac{\epsilon_{2}|\mu|^{2}}{1+\epsilon_{2}|\mu|}$. Apparently, we have:

[TABLE]

At last, we bound the tail probabilities of the terms containing $H_{22}^{-1}$ by:

[TABLE]

Similarly, we have:

[TABLE]

and

[TABLE]

Bounds on $\mathbb{P}\left(\|h_{12}^{T}H_{22}^{-1}\|\geq\sqrt{\epsilon_{1}}\right)$ and $\mathbb{P}\left(\|H_{22}^{-1}h_{21}\|\geq\sqrt{\epsilon_{1}}\right)$ can be easily derived from (32) and (33) respectively. But notice that $\exists n_{3}>0$ s.t. when $n\geq n_{3}$, $\frac{\sqrt{\epsilon_{3}m_{1}}}{2}n^{p/2}\leq\frac{1}{2}\min\{\epsilon_{1},\sqrt{\epsilon_{1}}\}m_{1}n^{p}$, hence the tail probabilities of $\|h_{12}\|$ and $\|h_{21}\|$ in (32)(33) are always bounded by those in (31) when $n$ large enough.

Combining the bounds in (25)-(33), when $n\geq\max\{n_{1},n_{2},n_{3}\}$, we have the following:

[TABLE]

where

[TABLE]

Lastly, notice that $\|h_{12}\|^{2}\leq\sum_{l=1}^{n}|d_{1l}|^{2}$, $\|h_{21}\|^{2}\leq\sum_{k=1}^{n}|d_{k1}|^{2}$ and $\|D_{(1)(1)}\|\leq\|D\|$, by applying the following bounds:

[TABLE]

to (34), we finishes the proof. ∎