{'fraction_non_alphanumeric': 0.04859139142131221, 'fraction_numerical': 0.01522567628157226, 'mean_word_length': 4.398245260185559, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 7, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Growing Linear Consensus Networks Endowed by Spectral Systemic Performance Measures

We propose an axiomatic approach for design and performance analysis of noisy linear consensus networks by introducing a notion of systemic performance measure. This class of measures are spectral functions of Laplacian eigenvalues of the network that are monotone, convex, and orthogonally invariant with respect to the Laplacian matrix of the network. It is shown that several existing gold-standard and widely used performance measures in the literature belong to this new class of measures. We build upon this new notion and investigate a general form of combinatorial problem of growing a linear consensus network via minimizing a given systemic performance measure. Two efficient polynomial-time approximation algorithms are devised to tackle this network synthesis problem: a linearization-based method and a simple greedy algorithm based on rank-one updates. Several theoretical fundamental limits on the best achievable performance for the combinatorial problem is derived that assist us to evaluate optimality gaps of our proposed algorithms. A detailed complexity analysis confirms the effectiveness and viability of our algorithms to handle large-scale consensus networks.


I. INTRODUCTION
The interest in control systems society for performance and robustness analysis of large-scale dynamical network is rapidly growing [1]- [10]. Improving global performance as well as robustness to external disturbances in large-scale dynamical networks are crucial for sustainability, from engineering infrastructures to living cells; examples include a group of autonomous vehicles in a formation, distributed emergency response systems, interconnected transportation networks, energy and power networks, metabolic pathways and even financial networks. One of the fundamental problems in this area is to determine to what extent uncertain exogenous inputs can steer the trajectories of a dynamical network away from its working equilibrium point. To tackle this issue, the primary challenge is to introduce meaningful and viable performance and robustness measures that can capture essential characteristics of the network. A proper measure should be able to encapsulate transient, steady-state, macroscopic, and microscopic features of the perturbed large-scale dynamical network.
In this paper, we propose a new methodology to classify several proper performance measures for a class of linear consensus networks subject to external stochastic disturbances. We take an axiomatic approach to quantify essential functional properties of a number of sensible measures by introducing the class of systemic performance measures and show that this class of measures should satisfy monotonicity, convexity, and orthogonal invariance properties. It is shown that several existing and widely used performance measures in the literature are in fact special cases of this class of systemic measures [3], [6], [11]- [13].
The performance analysis of linear consensus networks subject to external stochastic disturbances has been studied in [1], [2], [12]- [16], where the H2-norm of the network was employed as a scalar performance measure. In [1], the authors interpret the H2-norm of the system as a macroscopic performance measure capturing the notion of coherence. It has been shown that if the Laplacian matrix of the coupling graph of the network is normal, the H2-norm is a function of the eigenvalues of the Laplacian matrix [1], [12], [15]. In [2], the authors consider general linear dynamical networks and show that tight lower, and upper bounds can be obtained for the H2-norm of the network from the exogenous disturbance input to a performance output, which are functions of the eigenvalues of the state matrix of the network. Besides the commonly used H2-norm, there are several other performance measures that have been proposed in [1], [6], [17]. In [11], a partial ordering on linear consensus networks is introduced where it shows that several previously used performance measures are indeed Schur-convex functions in terms of the Laplacian eigenvalues. In a more relevant work, the authors of [18] show that performance measures that are defined based on some system norms, spectral, and entropy functions exhibit several useful functional properties that allow us to utilize them in network synthesis problems.
The first main contribution of this paper is introduction of a class of systemic performance measures that are spectral functions of Laplacian eigenvalues of the coupling graph of a linear consensus network. Several gold-standard and widely used performance measures belong to this class, for example, to name only a few, spectral zeta function, Gamma entropy, expected transient output covariance, system Hankel norm, convergence rate to consensus state, logarithm of uncertainty volume of the output, Hardy-Schatten system norm or Hp-norm, and many more. All these performance measures are monotone, convex, and orthogonally invariant. Our main goal is to investigate a canonical network synthesis problem: growing a linear consensus network by adding new interconnection links to the coupling graph of the network and minimizing a given systemic performance measure. In the context of graph theory, it is known that a simpler version of this combinatorial problem, when the cost function is the inverse of algebraic connectivity, is indeed NP-hard [19]. There have been some prior attempts to tackle this problem for some specific choices of cost functions (i.e., total effective resistance and the inverse of algebraic connectivity) based on semidefinite programing (SDP) relaxation methods [20], [21]. There is a similar version of this problem that is reported in [22], where the author studies convergence rate of circulant consensus networks by adding some long-range links. Moreover, a continuous (non-combinatorial) and relaxed version of our problem of interest has some connections to the sparse consensus network design problem [23]- [25], where they consider ℓ1-regularized H2optimal control problems. The other related works [26], [27] argue that some metrics based on controllability and observability Gramians are modular or submodular set functions, where they aim to show that their proposed simple greedy heuristic algorithms have guarantees sub-optimality bounds.
In our second main contribution, we propose two efficient polynomial-time approximation algorithms to solve the above mentioned combinatorial network synthesis problem: a linearization-based method and a simple greedy algorithm based on rank-one updates. Our complexity analysis asserts that computational complexity of our proposed algorithms are reasonable and make them particularly suitable for synthesis of large-scale consensus networks. To calculate sub-optimality gaps of our proposed approximation algorithms, we quantify the best achievable performance bounds for the network synthesis problem in Section V. Our obtained fundamental limits are exceptionally useful as they only depend on the spectrum of the original network and they can be computed a priori. In Subsection VII-B, we classify a subclass of differentiable systemic performance measures that are indeed supermodular. For this subclass, we show that our proposed simple greedy algorithm can achieve a (1 − 1/e)approximation of the optimal solution of the combinatorial network synthesis problem. Our extensive simulation results confirm effectiveness of our proposed methods.

A. Mathematical Background
The set of real numbers is denoted by R, the set of non-negative by R+, and the set of positive real numbers by R++. The cardinality of set E is shown by |E |. We assume that ½n, In, and Jn denote the n×1 vector of all ones, the n × n identity matrix, and the n × n matrix of all ones, respectively. For a vector v = [vi] ∈ R n , diag(v) ∈ R n×n is the diagonal matrix with elements of v orderly sitting on its diameter, and for A = [aij] ∈ R n×n , diag(A) ∈ R n is diagonal elements of square matrix A. We denote the generalized matrix inequality with respect to the positive semidefinite cone S n + by " " . Throughout this paper, it is assumed that all graphs are finite, simple, undirected, and connected. A graph herein is defined by a triple G = (V, E , w), where V is the set of nodes, E ⊆ {i, j} i, j ∈ V, i = j is the set of links, and w : E → R++ is the weight function. The adjacency matrix A = [aij ] of graph G is defined in such a way that aij = w(e) if e = {i, j} ∈ E , and aij = 0 otherwise. The Laplacian matrix of G is defined by L := ∆ − A, where ∆ = diag[d1, . . . , dn] and di is degree of node i. We denote the set of Laplacian matrices of all connected weighted graphs with n nodes by Ln. Since G is both undirected and connected, the Laplacian matrix L has n − 1 strictly positive eigenvalues and one zero eigenvalue. Assuming that 0 = λ1 < λ2 ≤ . . . ≤ λn are eigenvalues of Laplacian matrix L, we define operator Λ : S n + → R n−1 ++ by The Moore-Penrose pseudo-inverse of L is denoted by L † = [l † ji ], which is a square, symmetric, doubly-centered and positive semidefinite matrix. For a given link e = {i, j}, re(L) denotes the effective resistance between nodes i and j in a graph with the Laplacian matrix L, where its value can be calculated as follows where For every real q, powers of pseudo inverse of L is represented by L †,q := L † q .
Definition 1: The derivative of a scalar function ρ(.), with respect to the n-by-n matrix X, is defined by . The directional derivative of function ρ(X) in the direction of matrix Y is given by where ., . denotes the inner product operator.
The following Majorization definition is from [28].

Definition 2:
For every x ∈ R n + , let us define x ↓ to be a vector whose elements are a permuted version of elements of x in descending order. We say that x majorizes y, which is denoted by x ☎ y, if and only if 1 T x = 1 T y and The vector majorization is not a partial ordering. This is because from relations x ☎ y and y ☎ x one can only conclude that the entries of these two vectors are equal, but possibly with different orders. Therefore, relations x ☎ y and y ☎ x do not imply x = y.
Definition 3 ( [28]): The real-valued function F : R n + → R is called Schur-convex if F (x) ≥ F (y) for every two vectors x and y with property x ☎ y.

B. Noisy linear consensus networks
We consider the class of linear dynamical networks that consist of multiple agents with scalar state variables xi and control inputs ui whose dynamics evolve in time according tȯ for all i = 1, . . . , n, where xi(0) = x * i is the initial condition and is the average of all states at time instant t. The impact of the uncertain environment on each agent's dynamics is modeled by the exogenous noise input ξi(t). By applying the following feedback control law to the agents of this network the resulting closed-loop system will be a first-order linear consensus network. The closed-loop dynamics of network (3)-(4) with feedback control law (5) can be written in the following compact forṁ with initial condition and the output matrix is a centering matrix that is defined by The underlying coupling graph of the consensus network (6)-(7) is a graph G = (V, E , w) with node set V = {1, . . . , n}, edge set and weight function w(e) = kij (11) for all e = {i, j} ∈ E , and w(e) = 0 if e / ∈ E . The Laplacian matrix of graph G is equal to L.
Property (b) implies that feedback gains are symmetric and (c) means that there is no self-feedback loop in the network.

Assumption 2:
The coupling graph G of the consensus network (6)-(7) is connected and time-invariant.
According to Assumption 1, the underlying coupling graph is undirected and simple. Assumption 2 implies that only one of the modes of network (6) is marginally stable with eigenvector ½n and all other ones are stable. The marginally stable mode, which corresponds to the only zero Laplacian eigenvalue of L, is unobservable from the output (7). The reason is that the output matrix of the network satisfies Mn½n = 0. When there is no exogenous noise input, i.e., ξ(t) = 0 for all time, state of all agents converge to a consensus state [17], [29], which for our case the consensus state is When the network is fed with a nonzero exogenous noise input, the limit behavior (12) is not expected anymore and the state of all agents will be fluctuating around the consensus state without converging to it. Before providing a formal statement of the problem of growing a linear consensus network, we need to introduce a new class of performance measures for networks (6)-(7) that can capture the effect of noise propagation throughout the network and quantify degrees to which the state of all agents are dispersed from the consensus state.

III. SYSTEMIC PERFORMANCE MEASURES
The notion of systemic performance measure refers to a real-valued operator over the set of all linear consensus networks governed by (6)-(7) with the purpose of quantifying the quality of noise propagation in these networks. We have adopted an axiomatic approach to introduce and categorize a class of such operators that are obtained through our close examination of functional properties of several existing gold standard measures of performance in the context of network engineering and science. In order to state our findings in a formal setting, we observe that every network with dynamics (6)-(7) is uniquely determined by its Laplacian matrix. Therefore, it is reasonable to define a systemic performance measure as an operator over the set of Laplacian matrices Ln.
Definition 4: An operator ρ : Ln → R is called a systemic performance measure if it satisfies the following properties for all Laplacian matrices in Ln: 3. Orthogonal invariance: For all orthogonal matrices U ∈ R n×n , ρ(L) = ρ(U LU T ).
Property 1 guarantees that strengthening couplings in a consensus network never worsens the network performance with respect to a given systemic performance measure. The coupling strength among the agents can be enhanced by several means, for example, by adding new feedback interconnections and/or increasing weight of an individual feedback interconnection. The monotonicity property induces a partial ordering 1 on all linear consensus networks governed by (6)- (7). Property 2 requires that a viable performance measure should be amenable to convex optimization algorithms for network synthesis purposes. Property 3 implies that a systemic performance measure depends only on the Laplacian eigenvalues.
Theorem 1: Every operator ρ : Ln → R that satisfies Properties 2 and 3 in Definition 4 is indeed a Schur-convex function of Laplacian eigenvalues, i.e., there exists a Schur-convex spectral function Φ : R n−1 → R such that ρ(L) = Φ(λ2, . . . , λn). (13) Proof: For every L ∈ Ln, the value of the systemic performance measure can be written as a composition of two functions as follows where function Λ : S n + → R n−1 ++ is defined by (1) and function φ : for any matrix W = EU with U ∈ R n×n being an orthogonal matrix satisfying L = U T diag([0, Λ(L) T ])U and E ∈ R (n−1)×n given by the following projection matrix Thus, we can conclude that (13) holds with Φ(λ2, . . . , λn) = φ(Λ(L)). In the next step, we need to show that operator ρ is convex and symmetric with respect to Laplacian eigenvalues λ2, . . . , λn. Property 2 indicates that ρ is convex on Laplacian matrices and any convex function on Laplacian matrices is also convex function with respect to Laplacian eigenvalues [30]. Property 3 implies that operator ρ is symmetric with respect to λ2, . . . , λn as ρ is invariant under any permutation of eigenvalues. It is known that every function that is convex and symmetric is also Schur-convex [30]. The Laplacian eigenvalues of network (6)-(7) depend on global features of the underlying coupling graph. This is the reason why every performance measure that satisfies Definition 4 is tagged with adjective systemic. Table I shows some important examples of systemic performance measures and their corresponding matrix operator forms. In the appendix, we prove functional properties and discuss applications of these measures in details.

IV. GROWING A LINEAR CONSENSUS NETWORK
The network synthesis problem of interest is to improve the systemic performance of network (6)-(7) by establishing k ≥ 1 new feedback interconnections among the agents. Suppose that the underlying graph of the network G = (V, E , w) is defined according to (10)-(11) and a set of candidate feedback interconnection links Ec = ε1, . . . , εp ⊆ V ×V, which is endowed with a weight function ̟ : Ec → R++, is also given. The weight of a link εi ∈ Ec is represented by ̟(εi) and we assume that it is pre-specified and fixed. The network growing problem is to select exactly k feedback interconnection links from Ec and append them to G such that the systemic performance measure of the resulting network is minimized over all possible choices.
Let us represent the set of all possible appended subgraphs bŷ where the set of all possible choices to select k links is denoted by Then, the network synthesis problem can be cast as the following combinatorial optimization problem Systemic Performance Measure Matrix Operator Form Spectral zeta function ζq(L) Gamma entropy Iγ(L) Expected transient output covariance τt(L) System Hankel norm η(L) Uncertainty volume of the output υ(L) Hardy-Schatten system norm or Hp-norm θp(L) whereL is the Laplacian matrix of an appended candidate subgrapĥ G and the resulting network with Laplacian matrix L +L is referred to as the augmented network. The role of the candidate set Ec is to pre-specify authorized locations to establish new feedback interconnections in the network. The network synthesis problem (17) is inherently combinatorial and it is known that a simpler version of this problem with ρ(L) = λ −1 2 is in fact NP-hard [19]. There have been some prior attempts to tackle problem (17) for some specific choices of performance measures, such as total effective resistance and the inverse of algebraic connectivity, based on convex relaxation methods [20], [21] and greedy methods [27]. In Sections VI and VII, we propose approximation algorithms to compute sub-optimal solutions for (17) with respect to the broad class of systemic performance measures. We propose an exact solution for (17) when k = 1 and two tractable and efficient approximation methods when k > 1 with computable performance bounds. Besides, in Section VII, we demonstrate that a subclass of systemic performance measures has a supermodularity property. This provides approximation guarantees for our proposed approximation algorithm.

V. FUNDAMENTAL LIMITS ON THE BEST ACHIEVABLE PERFORMANCE BOUNDS
In the following, we present theoretical bounds for the best achievable values for the performance measure in (17). Let us denote the optimal cost value of the optimization problem (17) by r * k (̟). For a given systemic performance measure ρ : Ln → R, we recall that according to Theorem 1 there exists a spectral function Φ such that ρ(L) = Φ λ2, . . . , λn .
Theorem 2: Suppose that a consensus network (6)- (7) with an ordered set of Laplacian eigenvalues λ2 ≤ . . . ≤ λn, a set of candidate links Ec endowed with a weight function ̟ : Ec → R++, and design parameter 1 ≤ k ≤ n − 1 are given. Then, the following inequality holds for all weight functions ̟. For k ≥ n, all lower bounds are equal to Φ ∞, . . . , ∞ . Moreover, if the systemic performance measure has the following decomposable form where ϕ : R → R+ is a decreasing convex function and lim λ→∞ ϕ(λ) = 0, then the best achievable performance measure is characterized by Proof: For a given weight function ̟ : Ec → R++, we show that inequality (18) holds for everyÊ ∈ Π k (Ec). Assume thatL is the Laplacian of the graph formed by k added edges. We note that Therefore, we can define the nonempty set Mj for 2 ≤ j ≤ n, as follows where ui's and vi's are orthonormal eigenvectors of L and L +L, respectively. We now choose a unit vector v ∈ Mj . It then follows that: Therefore, according to (20) and the monotonicity property of the systemic measure ρ, we get for allÊ ∈ Π k (Ec). Inequality (18) now follows from (21) and this completes the proof. Note that inequality (19) is a direct consequence of (18) and lim λ→∞ ϕ(λ) = 0.
Theorem 3: Suppose that in optimization problem (17), the set of candidate links form a complete graph, i.e., |Ec| = 1 2 n(n − 1). Then, there exists a weight function ̟0 : Ec → R++ and a choice of k weighted links from Ec with weight function ̟ : holds for all weight functions ̟ that satisfies ̟(e) ≥ ̟0(e) for all e ∈ Ec. Moreover, if the systemic performance measure has the following decomposable form where ϕ : R → R+ is a decreasing convex function and lim λ→∞ ϕ(λ) = 0, then the best achievable performance measure is characterized by Proof: We will show that there existsÊ ∈ Π k (Ec) for which (22) is satisfied. Without loss of generality, we may assume that k < n−1. This is because otherwise, by adding n − 1 links, which forms a spanning tree, and increasing their weights the performance of the resulting network tends to Φ(∞, · · · , ∞) (see Theorem 4). LetÊ ⊂ E be the set of k links that do not form any cycle with ̟0(e) = ∞ for all e ∈Ê . Then, we know that and the k largest eigenvalues of L +L are equal to ∞. Using (24) and the monotonicity property of the systemic performance measure, we get From r * (̟) ≤ ρ(L +L) and using (25), we obtain (22). Note that inequality (23) is a direct consequence of (22) and lim λ→∞ ϕ(λ) = 0.
Examples of systemic performance measures that satisfies conditions of Theorem 2 include ζ q q (L) for q ≥ 1, Iγ(L), and τt(L). Theorem 4: Let us consider a linear consensus network (6)-(7) that is endowed with systemic performance measure ρ : Ln → R. Then, the network performance can be arbitrarily improved 2 by adding only n − 1 links that form a spanning tree.
Proof: Let us denote the Laplacian matrix of the spanning tree by LT . In the following, we show that the performance of resulting network can be arbitrarily improved by increasing the weights of the spanning trees. Based on the monotonicity property, we have Also, we know that Λ(κLT ) = κΛ(LT ). Therefore, using the fact that the spanning tree has only one zero eigenvalue, (13), we get lim κ→∞ ρ(κLT ) = Φ(∞, · · · , ∞).
Using this limit and (26) we get the desired result. It should be emphasized that by increasing weights of all the edges, the network performance can be arbitrarily improved, i.e., the value of the systemic performance measure can be made arbitrarily close to Φ(∞, · · · , ∞). Theorem 4 sheds more light on this fact by revealing the minimum number of required links and their graphical topology to achieve this goal.
The results of Theorems 2 and 3 can be effectively applied to select a suitable value for the design parameter k in optimization problem (17). Let us denote the value of the lower bound in (18) by ̺ k . The performance of the original network is then ̺0 = ρ(L). The percentage of performance enhancement can be computed by formula × 100 for all values of parameter 1 ≤ k ≤ n − 1. For a given desired performance level, we can look up these numbers and find the minimum number of required links to be added to the network. This is explained in details in Example 5 and Figure 5 in Section VIII. In next sections, we propose approximation algorithms to compute near-optimal solutions for the network synthesis problem (17).

VI. A LINEARIZATION-BASED APPROXIMATION METHOD
Our first approach is based on a linear approximation of the systemic performance measure when weights of the candidate links in Ec are small enough. In the next result, we calculate Taylor expansion of a systemic performance measure using notions of directional derivative for spectral functions. Algorithm: Adding k links using linearization Input: L, Ec, ̟, and k 1: setL = 0 2: for i = 1 to k 3: find e = {i, j} ∈ Ec that returns the maximum value for 4: set the solution e ⋆ 6: update 7:L =L + ̟(e ⋆ )L e ⋆ , and 8: : end for Lemma 1: Suppose that a linear consensus network (6)-(7) endowed with a differentiable systemic performance measure ρ is given. Let us consider the cost function in optimization problem (17). IfL is the Laplacian matrix of an appended subgraphĜ = (V,Ê, ̟), then where the derivative of ρ at L is given by for any matrix W that is defined by (15). Proof: The expression (27) can be calculated using the spectral form of a given systemic performance measure described by (14) and according to [31,Corollary 5.2.7]. Using the directional derivative of ρ along matrixL, the Taylor expansion of ρ(L + ǫL) is given by where ▽Lρ(L) is the directional derivative of ρ at L along matrixL where ., . denotes the inner product operator. Then, substituting (29) in (28) yields the desired result. According to the monotonicity property of systemic performance measures, the inequality Tr ▽ρ(L)L ≤ 0 holds for every Laplacian matrixL. This implies that when weights of the candidate links are small enough, one can approximate the optimization problem (17) by the following optimization problem whereL is the Laplacian matrix of an appended candidate subgrapĥ G = (V,Ê, ̟). Therefore, the problem boils down to select the klargest elements of the following set where ̟(e) is weight of link e. Table II presents our linearization approach as an algorithm. In some special cases, one can obtain an explicit closed-form formula for systemic performance measure of the resulting augmented network. Theorem 5: Suppose that linear consensus network (6)- (7) with Laplacian matrix L is endowed with systemic performance measure (74) for q = 1. Let us consider optimization problem (17), wherê L is the Laplacian matrix of a candidate subgraphĜ = (V,Ê, ̟). Then, where re(L 2 ) is the effective resistance between the two ends of e in a graph with node set V and Laplacian matrix L 2 .
Proof: We use the following identity for given matrices A, X ∈ R n×n . Based on [2,Theorem 4], the performance measure ζ1(.) can be calculated by Moreover, according to the definition of the Moore-Penrose generalized matrix inverse, we have Using (31) and (32), it follows that Then we show that Using (32), (33) and (34), we get the desired result. According to Theorem 5, when weights of the candidate links are small, in order to solve problem (17), it is enough to find k-largest element of the following set Since the weights of the candidate links are given, we only need to calculate the effective resistance re(L 2 ) for all e ∈ Ec.
As we discussed earlier, the design problem (17) is generally NPhard. Our proposed approximation algorithm in this section works in polynomial-time. In example 6, we discuss and compare optimality gap and time complexity of this method with other methods. The computational complexity of the linearization-based algorithm in Table II is O(n 3 ) for a given differentiable systemic performance measure from Table I. This involves computation of ▽ρ for the original graph, which requires O(n 3 ) operations. The rest of the algorithm can be done in O(pk) for small k and O(p log p) operations for large k.

VII. GREEDY APPROXIMATION ALGORITHMS
In this section, we propose an optimal algorithm to solve the network growing problem (17) when k = 1. It is shown that for some commonly used systemic performance measures, one can obtain a closed-form solution for k = 1. We exploit our results and propose a simple greedy approximation algorithm for (17) with k > 1 by adding candidate links one at a time. For some specific subclasses of systemic performance measures, we prove that our proposed greedy approximation algorithm enjoys guaranteed performance bounds with respect to the optimal solution of the combinatorial problem (17). Finally, we discuss time complexity of our proposed algorithms.

A. Simple Greedy by Sequentially Adding Links
The problem of adding only one link can be formulated as follows where Le is the Laplacian matrix of a candidate subgraphĜe = (V, {e}, ̟). Let us denote the optimal cost of (35) by r * 1 (̟). In order to formulate the optimal cost value of (35), we need to define the notion of a companion operator for a given systemic performance measure.
for all L ∈ Ln. Moreover, the companion operator of ρ is characterized by for all X ∈ Ln with eigenvalues µ2 ≤ . . . ≤ µn, where operator Φ : R n−1 → R is defined by (13).
In addition, we know that for the Moore-Penrose pseudo-inverse of matrix L ∈ Ln, we have the following for i = 2, . . . n, and λ1(L) = λ1(L † ) = 0. Consequently, we can rewrite ρ(L) using its companion operator as Therefore, by defining ψ : Ln → R as (37), we get identity (36). Theorem 6: Suppose that a linear consensus network (6)-(7) endowed by a systemic performance measure ρ : Ln → R is given. The optimal cost value of the optimization problem (35) is given by where ψ is the corresponding companion operator of ρ and Ue for a link e = {i, j} is a rank-one matrix defined by in which L † i is the i th column of matrix L † . Proof: We use the following matrix identity where Ee is the incidence matrix of graphĜe andL = L + 1 n Jn. By utilizing the Woodbury matrix identity, we get From the definition of the effective resistance between nodes i and j, it follows that On the other hand, we havē Therefore, using (40) and (41), we have From (36) and (42), we can conclude the desired equation (38).
In some special cases, the optimal solution (38) can be computed very efficiently using a simple separable update rule.
Theorem 7: Suppose that linear consensus network (6)- (7) with Laplacian matrix L is given. Then, for every link e ∈ Ec we have According to the matrix determinant lemma we have Now using (43) and (44), it follows that det 2(L + Le) + 1 n Jn then by taking log from both sides, we get the desired result. In these special cases, the computational complexity of calculating the optimal solution for network design problem (35) is relatively low. For q = 1, the optimal cost value is equal to ζ1(L + Le * ), where and for q = 2, the optimal cost value is equal to ζ2(L + Le * ), where Moreover, for (90), the optimal cost value is equal to υ(L + Le * ), where e * = arg min e∈Ec log 1 + re(L)̟(e) .
The location of the optimal link is sensitive to its weight. For example when optimizing with respect to ζ1, maximizers of re(L), re(L 2 ) and re(L 2 )/re(L) can be three different links. In Example 3 and Fig. 1 of Section VIII, we illustrate this point by means of a simulation.
Furthermore, one can obtain the following useful fundamental limits on the best achievable cost values.
Theorem 8: Let us denote the value of performance improvement by adding an edge e with an arbitrary positive weight to linear consensus network (6)-(7) by ∆ρ(L) = ρ(L) − ρ(L + Le).
Then, the maximum achievable performance improvement is where Ue is given by (39) and the upper bound can be achieved as w tends to infinity. Moreover, we have the following explicit fundamental limits Proof: We utilize monotonicity property of companion operator of a systemic performance measure, i.e., If L † , and the inequality Ue .
From this inequality, we can directly conclude (46). For systemic performance measure ζ1(.), inequality (46)  Similarly for ζ 2 2 (.), using (46) and the definition of ζ2(.), results in This completes proof. Gamma entropy Uncertainty volume of the output log µ i Hardy-Schatten system norm or Hp-norm The result of Theorem 6 can be utilized to devise a greedy approximation method by decomposing (17) into k successive tractable problems in the form of (35). In each iteration, Laplacian matrix of the network is updated and then optimization problem (35) finds the next best candidate link as well as its location. Since the value of systemic performance measure can be calculated explicitly in each step using Theorem 6, one can explicitly calculate the value of systemic performance measure for the resulting augmented network. This value can be used to determine the effectiveness of this method. Table III summarizes all steps of our proposed greedy algorithm, where the output of the algorithm is the Laplacian matrix of the resulting augmented network. In Section VIII, we present several supporting numerical examples.
Remark 1: The optimization problem (35) with performance measure ζ∞(L) = λ −1 2 was previously considered in [21], where a heuristic algorithm was proposed to compute an approximate solution. Later on, another approximate method for this problem was presented in [20]. Also, there is a similar version of this problem that is reported in [22], where the author studies convergence rate of circulant consensus networks by adding some long-range links. Moreover, a non-combinatorial and relaxed version of our problem of interest has some connections to the sparse consensus network design problem [23]- [25], where they consider ℓ1-regularized H2 optimal control problems. When the candidate set Ec is the set of all possible links except the network links, i.e., Ec = V × V \ E , and the performance measure is the logarithm of the uncertainty volume, our result reduces to the result reported in [26].

B. Supermodularity and Guaranteed Performance Bounds
A systemic performance measure is a continuous function of link weights on the space of Laplacian matrices Ln. Moreover, we can represent a systemic performance measure equivalently as a set function over the set of weighted links. Let us denote by G(V) the set of all weighted graphs with a common node set V.
Definition 5: For a given systemic performance measure ρ : Ln → R, we associate a set functionρ : G(V) → R that is defined as The following definition is adapted from [32] for our graph theoretic setting.
Theorem 9: Suppose that systemic performance measure ρ : Ln → R is differentiable and ▽ρ : Ln → R n×n is monotonically increasing with respect to the cone of positive semidefinite matrices 3 . Then, the corresponding set functionρ : G(V) → R, from Definition 5, is supermodular.
It should be emphasized that convexity property of a systemic performance measure ρ implies that ▽ρ, if it exists, is a monotone mapping 4 . However, this property is not sufficient for supermodularity of its corresponding set functionρ.
Example 1: In our first example, we show that the uncertainty volume of the output (90) satisfies conditions of Theorem 9. The gradient operator of this systemic performance measure is It is straightforward to verify that ▽υ(L) is monotone with respect to the cone of positive semidefinite matrices. Thus, υ(L) is supermodular. Example 2: In our second example, we consider a new class of systemic performance measures that are defined as where 0 ≤ q ≤ 1. According to Theorem 11, this spectral function is a systemic performance measure as function −λ q for 0 ≤ q ≤ 1 is a decreasing convex function on R+. Moreover, its gradient operator, which is given by ▽mq(L) = qL q−1 is monotonically increasing for all 0 ≤ q ≤ 1. Therefore, according to Theorem 9, systemic performance measure (63) is supermodular over the set of all weighted graphs with a common node set. Remark 2: For a given performance measure ρ, there are several different ways to define an extended set function for ρ. These set functions may have different properties. For instance, the extended set function of ζ1 is supermodular over principle sub-matrices [33], but it is not supermodular over the set of all weighted graphs with a common node set (see Definition 5).
For those systemic performance measures that satisfy conditions of Theorem 9, one can provide guaranteed performance bounds for our proposed greedy algorithm in Subsection VII-A. The following result is based a well-known result from [32, Chapter III, Section 3].
Theorem 10: Suppose that systemic performance measure ρ : Ln → R is differentiable and ▽ρ : Ln → R n×n is monotonically increasing with respect to the cone of positive semidefinite matrices. Then, the greedy algorithm in Table III, which starts withÊ as the empty set and at every step selects an element e ∈ Ec that minimizes the marginal cost ρ(L + LÊ + Le) − ρ(L + LÊ ), provides a setÊ that achieves a (1 − 1/e)-approximation 5 of the optimal solution of the combinatorial network synthesis problem (17).
Since the class of supermodular systemic performance measures are monotone, the combinatorial network synthesis problem (17) is polynomial-time solvable with provable optimality bounds [32]. Supermodularity is not a ubiquitous property for all systemic performance measures. Nevertheless, our simulation results in Section VIII assert that the proposed greedy algorithm in Table III is quite powerful and provides tight and near-optimal solutions for a broad range of systemic performance measures.

C. Computational Complexity Discussion
As we discussed earlier, the network synthesis problem (17) is in general NP-hard. However, this problem is solvable when k = 1 and the best link can be found by running an exhaustive search over all possible scenarios, i.e., by calculating the value of a performance measure for all possible p augmented networks, where p is the number of candidate links. The computational complexity of evaluating 5 This means that where L * is the optimum solution andL is the solution of the greedy algorithm, or equivalently: where e is Euler's number.  (Mρ(n)). In the simple greedy algorithm of Table  III, the difference term is calculated and updated for each candidate link at each step, for the total of k p − k−1 2 times. Thus, the total computational complexity of our simple greedy algorithm is O Mρ(n)(p − k−1 2 )k operations. This computational complexity is at most O Mρ(n)n 2 k , where p = n 2 , i.e., when the candidate set contains all possible links. The complexity of the brute-force method is O Mρ(n) p k 6 . This can be at most O Mρ(n)2 p / √ p . Moreover, if k ≤ √ p, then the computational complexity will be O Mρ(n)p k /k! .
In some occasions, we can take advantage of the rank-one updates in Theorems 6 and 7, where it is shown that a rank-one deviation in a matrix results in a rank-one change in its inverse matrix as well. This helps reduce the computational complexity of (64) to the order of O(n 2 ) instead of O(n 3 ) operations. As it is shown in [34], one can apply the rank-one update on the matrix of effective resistances. As a result, we can update the effective resistances of all links in order of O(n 2 ). More specifically, the matrix of effective resistances is given by for m ∈ {1, 2, 3}, where R(L m )ij = r {i,j} (L m ). The update rule (65) can be obtained by substituting the rank-one update of (L+Le) † from (42) in (65) and the m-th power of the rank-one update can be calculated in O(n 2 ) as it can be cast as only matrix-vector products. Using these facts and the result of Theorem 7, the computational cost of (64) for systemic performance measures ζ1, ζ2, and υ can be significantly reduced; more specifically, the computational complexity of our algorithm reduces to For a generic systemic performance measure ρ : Ln → R, according to Theorem 1, calculating its value requires knowledge of all Laplacian eigenvalues of the coupling graph. It is known that the eigenvalue problem for symmetric matrices requires O(n 2.376 log n) operations [35]. Suppose that calculating the value of spectral function Φ : R n−1 → R in Theorem 1 needs O (MΦ(n)) operations. Thus, the value of systemic performance measure ρ(L) in equation (13), Based on this analysis, we conclude that the complexity of the greedy algorithm in Table III is at most

VIII. NUMERICAL SIMULATIONS
In this section, we support our theoretical findings by means of some numerical examples.
Example 3: This example investigates sensitivity of location of an optimal link as a function of its weight. Let us consider a linear consensus network (6)-(7), whose coupling graph is shown in Fig.  1, endowed by systemic performance measure (74) with q = 1. The graph shown in Fig. 1 is a generic unweighted connected graph with n = 50 nodes and 100 links. We solve the network synthesis problem (35) for the candidate set with |Ec| = 1 2 n(n − 1) that covers all possible locations in the graph. It is assumed that all candidate links have an identical weight ̟0. We use our rank-one update method in Theorem 7 to study the effect of ̟0 on location of the optimal link. In Fig. 1(c), we observe that by increasing ̟0, the optimal location changes. When ̟0 = 1, our calculations reveal that the optimal link in Fig. 1(a), shown by a blue line segment, maximizes re(L 2 ) among all possible candidate links in set Ec. By increasing the value of our design parameter to ̟0 = 1.2 in Fig. 1(b), we observe that the location of the optimal link moves. In our last scenario in Fig.  1(c), by setting ̟0 = 1.6, the optimal link moves to a new location that maximizes quantity re(L 2 )/re(L) among all possible candidate links.
Example 4: The usefulness of our theoretical fundamental hard limits in Theorem 2 in conjunction with our results in Theorem 7 is illustrated in Fig. 3. Suppose that a linear consensus network (6)-(7) with a generic coupling graph with n = 60, as shown in Fig. 2(a), is given. Let us consider the network design problem (35) with systemic performance measure (74) for q = 1. The set of candidate links is the set of all possible links in the coupling graph, i.e., |Ec| = 1 2 n(n − 1), where it is assumed that all candidate links have an identical weight ̟0 = 20. Our goal is to compare optimality of our low-complexity update rule against brute-force search over all |Ec| = 1770 possible augmented graphs. The value of the systemic performance measure for each candidate graph is marked by blue star in Fig. 3. In this plot, the black circle highlights the value of performance measure for the network resulting from the rank-one search (45). The red dashed line in Fig. 3 shows the best achievable value for ζ1 according to Theorem 2. The value of this hard limit can be calculated merely using Laplacian eigenvalues of the original graph shown in Fig. 2(a). The location of the optimal link is shown in Fig. 2(b). One observes from Fig. 3 that our theoretical fundamental limit justifies near-optimality of our rank-one update strategy (45) for networks with generic graph topologies.

Example 5:
This example follows up on our discussion at the end of Section V, where it is explained that the result of Theorem 2 can be utilized to choose reasonable values for design parameter k in the network design problem (17). We explain the procedure by considering a linear consensus network (6)-(7) with a given coupling graph by Fig. 2(a). The value of the lower bound (i.e., hard limit) in (18) is used to form the following quantity that represents the percentage of performance enhancement for all values of parameter 1 ≤ k ≤ n − 1. Fig. 5 illustrates the value of π k with respect to four systemic performance measures: ζ1, ζ2, τt and Iγ . Depending on the desired level of performance, one can compute a sensible value for design parameter k merely by looking up at the corresponding plots. For instance, in order to achieve 50% performance improvement, one should add at least 13, 10, 16, and 12 weighted links with respect to ζ1, ζ2, τt and Iγ, respectively. We verified tightness of this estimate by running our greedy algorithm in Table III, where the candidate set is equal to the set of all possible links with identical weight 10. Our simulation results reveal that by adding 13, 10, 16, and 12 links from the candidate set, the network performance improves by 40.60%, 45.10%, 37.76%, and 40.61% with respect to ζ1, ζ2, τt, and Iγ , respectively. Our theoretical bounds predict that network performance can be further improved by increasing weights of the candidate links. In our example, if we increase the weight from 10 to 500, the network performance boosts by more than 46% for all mentioned systemic performance measures.

Example 6:
We compare optimality gaps of our proposed greedy (see Table III) and linearization-based (see Table II) methods versus brute-force and simple-random-sampling methods. The brute-force method runs an exhaustive search to find the global optimal solution of problem (17); however, it cannot be used for medium to large size networks. In order to make our comparison possible, we consider a linear consensus network (6)-(7) with n = 30 nodes over the graph shown in Fig. 4. Weights of all links, both in the coupling graph and the candidate set, are equal to 1. Our control objective is to solve the network synthesis problem (17), where the candidate set consists of 15 links that are shown by red-dashed lines in Fig. 4. The outcome of our simulation results are explicated in Fig. 6, where we run our algorithms and compute the corresponding values for systemic performance measures for all k = 1, . . . , 15. One observes that our greedy algorithm performs nearly as optimal as the bruteforce method. This is mainly due to convexity and monotonicity properties of the class of systemic performance measures that enable the greedy algorithm to produce near-optimal solutions with respect to this class of measures. As one expects, our greedy algorithm outperforms our linearization-based method. It is noteworthy that the time complexity of the linearization method is comparably less than the greedy algorithm. The usefulness of the linearization-based method accentuates itself when weight of candidate links are small and/or k is large.

IX. DISCUSSION AND CONCLUSION
In the following, we provide explanations for some of the outstanding and remaining problems related to this paper.

Convex Relaxation:
The constraints of the combinatorial problem (17) can be relaxed by allowing the link weights to vary continuously. The relaxed problem will be a spectral convex optimization problem [36]. In some special cases, such as when the cost function is ζ1 or ζ2, the relaxed problem can be equivalently cast as a semidefinite programming problem [11], [18]. However, for a generic systemic performance measure, we need to develop some low-complexity specialized optimization techniques to solve the corresponding spectral optimization problem, which is beyond the scope of this paper.
Higher-Order Approximations: In Subsection VI, we employed the first-order approximation of a systemic performance measure. One can easily extend our algorithm by considering second-order approximations of a systemic performance measure in order to gain better optimality gaps.

Non-spectral Systemic Performance Measures:
The class of spectral systemic performance measures can be extended to include nonspectral measures as well. This can done by relaxing and replacing the orthogonal invariance property by permutation invariance property. The local deviation error is an example of a non-spectral systemic performance measure [18], [37]. Our ongoing research involves a comprehensive treatment of this class of measures.

APPENDIX NOTABLE CLASSES OF SYSTEMIC PERFORMANCE MEASURES
In the following, we will revisit several existing and widely-used examples of performance measures in linear consensus networks and prove that they are indeed systemic performance measures according to the definition.

A. Sum of Convex Spectral Functions
This class of performance measures is generated by forming summation of a given function of non-zero Laplacian eigenvalues.
Theorem 11: For a given matrix L ∈ Ln, suppose that ϕ : R+ → R is a decreasing convex function. Then, the following spectral function is a systemic performance measure. Moreover, if ϕ is also a homogeneous function of order −κ with κ > 1, then the following spectral function is also a systemic performance measure. Proof: First we show that measure (66) is monotone with respect to the positive definite cone. If we assume that L2 L1, then based on Theorem A.1 in [28,Sec. 20], it follows that λi(L2) ≤ λi(L1), for i = 1, 2, · · · , n.
From (70) and the desired convexity property of ϕ(.), we get the convexity property as follows Using (71) and (67), we get where K = κ ϕ(1). It is well-known function (72) is convex for λi > 0 where i = 2, . . . , n and κ > 1. Based on the proof of Part (i), measure ρ κ (.) is a Schur-convex function. Consequently, we get Now using (74) and the convexity of (72) with respect to λi's, we have This completes the proof. There are several important examples of performance measures that belong to this class.

1) Spectral Zeta Functions:
For a given network (6)-(7), its corresponding spectral zeta function of order q ≥ 1 is defined by where λ2, . . . , λn are eigenvalues of L [38]. According to Assumption 2, all the Laplacian eigenvalues λ2, . . . , λn are strictly positive and, as a result, function (74) is well-defined. The spectral zeta function of a graph captures all its spectral features. In fact, it is straightforward to show that every two graphs in Ln with identical zeta functions for all parameters q ≥ 1 are isospectral 7 .
Since ϕ(λ) = λ −q for q ≥ 1 is a decreasing convex function, the spectral function (74) is a systemic performance measure according to Theorem 11. The systemic performance measure 1 2 ζ1(L) is equal to the H2-norm squared of a first-order consensus network (6)-(7) and 2) Gamma Entropy: The notion of gamma entropy arises in various applications such as the design of minimum entropy controllers and interior point polynomial-time methods in convex programming with matrix norm constraints [40]. As it is shown in [41], the notion of gamma entropy can be interpreted as a performance measure for linear time-invariant systems with random feedback controllers by relating the gamma entropy to the mean-square value of the closed-loop gain of the system. Definition 9: The γ-entropy of network (6)- (7) is defined as where G(jω) is the transfer function of network (6)-(7) from ξ to y.
Proof: First we obtain the transfer function of network (6)-(7) from ξ to y. In order to do that, let us rewrite the network in its disagreement form (85)-(86) (see [42] for more details). Then, it follows that G(s) = Mn sIn + L + 1 n Jn where U is the corresponding orthonormal matrix of eigenvectors of L. Now, we calculate the γ-entropy by substituting the transfer function (76) in (75) as follows Then, using the fact that U U T = In and (76), one can write: (77) Moreover, we know that for γ ≥ λ −1 2 . Note that fγ (.) is a convex decreasing function in [γ −1 , ∞), therefore, according to Theorem 11 and (79), the γ-entropy Iγ (L) is a systemic performance measure.
The following result presents the connection between the γ-entropy measure and the H2-norm of the network.
Theorem 13: The following equality holds for the γ-entropy measure of network (6)

3) Expected Transient Output Covariance:
We consider a transient performance measure at time instant t > 0 that is defined by where it is assumed that each ξi(t) for all t ≥ 0 is a white Gaussian noise with zero mean and unit variance and all ξi's are independent of each other.
In the following, we show that this performance measure is a spectral function of Laplacian eigenvalues.
Theorem 14: For a given linear consensus network (6)-(7), the transient measure can be expressed as Moreover, τt(L) is a systemic performance measure for all t > 0.
Proof: The covariance matrix of the output vector is governed by the following matrix differential equatioṅ where Y (t) = cov(y(t), y(t)). Using the closed-form solution of (82), which is given by we get E y T (t)y(t) = Tr(Y (t)) = Tr is convex and decreasing with respect to x on R+, we can conclude that τt(L) is a systemic performance measure according to Theorem 11. We note that when t tends to infinity, the value of the transient performance measure becomes equal to the H2-norm squared of the network, i.e., τ∞(L) = G 2 H 2 .

4) Hankel Norm:
The Hankel norm of a network with (6)-(7) and transfer function G(jω) from ξ to y is defined as the L2-gain from past inputs to the future outputs, i.e.,  The value of the Hankel norm of network (6)-(7) can be equivalently computed using the Hankel norm of its disagreement form [43] that is given byẋ where the disagreement vector is defined by The disagreement network (85)-(86) is stable as every eigenvalue of the state matrix −L d = −(L + 1 n Jn) has a strictly negative real part. One can verify that the transfer functions from ξ(t) to y(t) in both realizations are identical. Therefore, the Hankel norm of the system from ξ(t) to y(t) in both representations are well-defined and equal, and is given by [44] η(L) := G H = λmax(P Q), where the controllability Gramian P is the unique solution of Moreover, based on Theorem 11, we know that the spectral zeta function ζq(L) is a systemic performance measure for all 1 ≤ q ≤ ∞. Therefore by setting q = ∞, we have As a result, η(L) is a systemic performance measure.

5) Uncertainty volume:
The uncertainty volume of the steady-state output covariance matrix of consensus network (6)-(7) is defined by where Y∞ = lim t→∞ E y(t)y T (t) .
This quantity is widely used as an indicator of the network performance [11], [45]. Since y(t) is the error vector that represents the distance from consensus, the quantity (89) is the volume of the steadystate error ellipsoid.
Theorem 16: For a given consensus network (6)- (7) with Laplacian matrix L, the logarithm of the uncertainty volume, i.e., is a systemic performance measure.