Optimal Control of Networks in the presence of Attackers and Defenders
Ishan Kafle, Sudarshan Bartaula, Afroza Shirin, Isaac Klickstein, Pankaz Das, Francesco Sorrentino
OOptimal Control of Networks in the presence of Attackers and Defenders
Ishan Kafle, a) Sudarshan Bartaula, b) Afroza Shirin, c) Isaac Klickstein, d) Pankaz Das, e) and Francesco Sorrentino f)1) Department of Mechanical Engineering, University of New Mexico, Albuquerque,New Mexico 87131, USA Department of Electrical and Computer Engineering, University of New Mexico,Albuquerque, New Mexico 87131, USA (Dated: 9 April 2018)
We consider the problem of a dynamical network whose dynamics is subject to externalperturbations (‘attacks’) locally applied at a subset of the network nodes. We assume thatthe network has an ability to defend itself against attacks with appropriate countermea-sures, which we model as actuators located at (another) subset of the network nodes. Wederive the optimal defense strategy as an optimal control problem. We see that the networktopology, as well as the distribution of attackers and defenders over the network affect theoptimal control solution and the minimum control energy. We study the optimal controldefense strategy for several network topologies, including chain networks, star networks,ring networks, and scale free networks. a) Electronic mail: ikafl[email protected]. b) Electronic mail: [email protected]. c) Electronic mail: [email protected]. d) Electronic mail: [email protected]. e) Electronic mail: [email protected]. f) Electronic mail: [email protected]. a r X i v : . [ c s . S Y ] A p r ptimal control of networks is an area of recent interest in the literature, where focus hasbeen placed on how the network topology and the position of driver and target nodes affect theoptimal solution. Here we study a different but related problem, that of optimally controllinga network under attack. We investigate the role of the network topology as well as of thedistribution of attackers and defenders over the network. Some of our results are counter-intuitive, as we find that for small chain networks, star networks, and ring networks, thedistance between a single attacker and a single defender is not the key factor that determinesthe minimum control energy. We also consider the case of a large scale-free network in thepresence of a single attacker and multiple defenders for which we see that the minimumcontrol energy varies over different orders of magnitude as the position of the attacker ischanged over the network. For this case we observe that the minimum distance between theattacker node and the defender nodes is a good predictor of the strength of an attack.I. MODEL Most infrastructure systems are networked by design, such as power grids , road systems ,telecommunications , water and sewer systems , and many others. These networked systemsare prone to disruption by either natural causes, such as extreme weather events and aging equip-ment, or purposeful attack, such as terrorism . In power grid systems, small local failures havebeen known to cascade to blackouts affecting large swaths of a state or a country . In a roadsystem, small incidents can lead to large scale congestion . Attacks on networked systems canbe either structural, where links in the underlying graph are damaged or destroyed, or dynamical,where a disruptive term is added to the dynamical equations that govern the behavior of the system.While our approach can be extended to encompass both structural and dynamical attacks, for thesake of simplicity, here we focus on dynamical attacks. Some examples of dynamical attacks onnetworked systems are pollutants introduced in a hydraulic network or the spreading of virusesin networked computers .We examine the behavior of simple dynamical networks, when they are attacked by one or moreexternal signals which perturb the dynamics of the network nodes. To illustrate this situation,a ten node network where three of the nodes are under attack is shown in Fig. 1. We assume2hat the networks at hand have an ability to defend themselves against attacks with appropriatecounter-measures, which we model as actuators located at a subset of the nodes in the network. InFig. 1(e), nodes 1, 3, and 10 are attached to actuators, and so these nodes we define as defenders (equivalently driver nodes as they are defined in much of the complex network literature ). De-fenders can take many different forms in the networked systems described above, such as trafficsignals and GPS routing in road networks, purposely tripping lines in a power grid in case ofshedding, or quarantining a portion of a computer network when attacked by viruses.As a reference example, in this paper we consider a power grid, with nodes representing buses andedges representing transmission lines . Both generation and a load can be present at each bus.We assume that some of the loads are vulnerable to attacks, in which case ancillary generation atthe other nodes can be used to mitigate the effects of the attack. This problem is discussed in moredetail in what follows. 3 a)
52 91 84673 − . . . . S t a t e s (b) (c) w w w S t a t e s (d) (e) w w u w u u − , − ,
000 Time S t a t e s (f) FIG. 1: (a) Network under no attack. (b) Time evolution of the network nodes under no attack.(c) Same network as in (a), with attackers located at nodes 2, 5, and 7. (d) Time evolution of thenetwork nodes under attack. (e) Same network as in (c) with defenders located at nodes 1, 3, and10. (f) Time evolution of the network nodes under attack and response of the defender nodes.4ere, we consider a network which has stabilizing self-loops so that the adjacency matrix is Hur-witz. This ensures that after any perturbation of the states away from the origin, the states willreturn to the origin. Next, we add external attackers to the system attached to nodes 2, 5, and 7.In Fig. 1, we color nodes 2, 5, and 7 cyan as they are the attacked nodes. The time evolution ofthe states of those nodes directly attacked, and any nodes downstream such as nodes 3, 6, and10 are now diverging. On the other hand, any node upstream of the attackers are not effected bythe attack and so they will converge to the origin. To counter the attacks, we add external controlinputs attached to the red nodes 1, 3, and 10, which we call defender nodes. Thanks to the controlaction exerted by the defender nodes, the attack can be mitigated and now all nodes return to theorigin again.For simplicity, the network dynamics is described by a linear model,˙ x ( t ) = A x ( t ) + H w ( t ) + B u ( t ) (1)where x ( t ) = [ x ( t ) , .., x n ( t )] is the n × u ( t ) = [ u ( t ) , .., u m ( t )] is the m × w ( t ) = [ w ( t ) , .., w q ( t )] is the q × A = { a i j } is a square n × n real adjacency matrix which has non-zero elements a i j if node i receives a signal from node j and is 0 otherwise. The n × m matrix B is the control input matrixand describes how the control inputs are connected to the nodes, i.e. the location of the defenders,namely B i j is different from zero if the control input j is attached to node i and is zero otherwise.The matrix H models how the attackers affect the network nodes, namely H i j is different fromzero if attacker j is active on node i and is zero otherwise. The matrix A is Hurwitz and thereforeby setting w = and u = , the system asymptotically approaches the origin of state space, whichrepresents the nominal healthy condition for the system.As explained in Sec. III, the dynamics of a power grid can be cast in the form of Eq. (1), ˙ δδδ ˙ θθθ ˙ ωωω = A δδδθθθωωω + H PPP L + B PPP M (cid:48) . (2)where, the vector δδδ = [ δ , ..., δ n , ] contains information on the voltage phase angles at generatorbuses, the vector θθθ = [ θ , ..., θ n ] describes the voltage phase angles at load buses and the vector ωωω [ ω , ..., ω n , ] represents the frequency deviation at generator buses. The vector P L = [ P L , ..., P Ln ] contains information on the power consumption at the load buses and the vector P M (cid:48) = [ P M (cid:48) , ..., P M (cid:48) n ] represents the ancillary power generation (for more details see section III.)We now introduce the strong assumption that a known model for the attackers exists. This as-sumption could be more or less unrealistic depending on the application to which we are applyingthis methodology; however, our results are general as they can be applied to a variety of modelsfor the attackers’ strategy and as we will see, they are to some extent independent of the attackers’specific strategy. The type of attack strategies we consider is either one of the following functions:(i) constant, (ii) linearly increasing, or (iii) exponentially increasing, which can be modeled as:˙ w i = s i w i + r i (3)where, s i and r i are constants. We consider the following three cases:(i) if s i = r i = s i = r i > s i > r i = x ( t ) = ˜ A ˜ x ( t ) + ˜r + ˜ B u ( t ) (4) y ( t ) = C ˜ x ( t ) (5)where, ˜ A = A ... H · · · · · · · · · S , ˜ B = B · · · , C = (cid:104) I n ... 0 (cid:105) , and ˜r = · · · r , (6)Here ˜ x = [ x T , w T ] T is the n + q vector containing the states of the network nodes and attackers,the behavior of which is assumed to be known, y ( t ) = [ y ( t ) , .., y p ( t )] is the p × S =diag { s , .., s q } is the diagonal matrix that contains information on the attackersstrengths and the vector r =[ r , ..., r q ] describes the attackers strategies (see Eq. (3)). The matrix˜ A now has a block triangular structure and is non-Hurwitz, due to the attackers’ dynamics. The6atrix C relates the outputs y ( t ) to the state ˜ x ( t ) . In this particular case, y (t) coincides with x ( t ) in Eq. (1), i.e., it selects the states of the nodes but not those of the attackers.When u = , the time evolution of the network nodes deviates from the origin due to the influenceof the attackers. The question we will try to address is the following: how can we design anoptimal control input that in presence of an attack, will set the state x ( t f ) = at some preassignedtime t f . The time t f can be thought as the required time to neutralize the attackers, so that for t > t f , the network has returned to its healthy state and the control action is no needed anymore.Here, without loss of generality, we assume the optimal control input to be the one that minimizesthe energy function, E = (cid:90) t f t u T ( t ) u ( t ) dt (7)The control input u ∗ ( t ) that satisfies the constraints and minimizes the control energy is equal to : u ∗ ( t ) = B T e ˜ A T ( t f − t ) C T ( CWC T ) − × [ y f − Ce ˜ A T ( t f − t ) x − CF ( t f ) ˜r ] (8)where, F ( t f ) = (cid:82) t f t e ˜ A ( t f − τ ) d τ . The corresponding optimal energy is E ∗ = (cid:82) t f t u ∗ T ( t ) u ∗ ( t ) dt . First,we define the controllability Gramian as a real, symmetric, semi-positive definite matrix W = (cid:90) t f t e ˜ A ( t f − τ ) ˜ B ˜ B T e ˜ A T ( t f − τ ) dt (9)Following , the minimum control energy can be computed and is equal to E ∗ = ( y f − Ce ˜ A ( t f − τ ) − CF ( t f ) ˜r ) T ( CWC T ) − ( y f − Ce ˜ A ( t f − τ ) − CF ( t f ) ˜r ) dt = βββ T W − p βββ (10)where the vector βββ = Ce ˜ A ( t f − t ) x + CF ( t f ) ˜r − y f is the control maneuver and W p = CWC T is the p × p symmetric, real, non-negative definite output controllability Gramian. The smallest eigen-value of the output controllability Gramian, µ , is nonzero if the system is output controllable. Ifthis condition is satisfied, µ usually dominates the expression for the minimum control energy .7 . Effect of the Attackers on Output Controllability Gramian Consider the system with attackers (4) and (5). We write, e ˜ At = e At ... F ( t ) · · · · · · · · · e St , where F ( t ) = (cid:90) e ( − τ ) At Hte τ St d τ . Moreover,˜ B ˜ B T = B · · · (cid:104) B T ... 0 T (cid:105) = BB T ... 0 · · · · · · · · · ,The controllability Gramian, W = (cid:90) t f t e ˜ At ˜ B ˜ B T e ˜ A T t dt = (cid:90) t f t e At ... F ( t ) · · · · · · · · · e St BB T ... 0 · · · · · · · · · e A T t ... 0 · · · · · · · · · F T ( t ) ... e S T t dt = (cid:90) t f t e At BB T e A T t ... 0 · · · · · · · · · dt = W p ... 0 · · · · · · · · · (11)Note that W p ∈ R n × n does not depend on on the matrices S and E i.e., it is independent of thelocation of the attackers and the strength of the attackers . The output controllability Gramian, CWC T = (cid:104) I ... 0 (cid:105) W p ... 0 · · · · · · · · · I · · · = W p . (12)If the pair ( A , B ) is controllable, the matrix W p is positive definite and thus invertible .8 . Effect of the Attackers on Control Maneuver We have already defined the control maneuver as βββ = Ce ˜ A ( t f − t ) x + CF ( t f ) ˜r − y f (13)According to our assumptions we set y f = (target state coincides with the origin). We writethe eigenvalue equation for the matrix ˜ A , ˜ A = V Λ V − , where the eigenvector matrix V = (cid:2) v v · · · v n + q (cid:3) and the eigenvalue matrix Λ = diag { λ , · · · , λ q , λ q + , · · · , λ q + n } where, λ ≥ · · · ≥ λ q ≥ λ q + ≥ · · · ≥ λ q + n . Note that because of the block diagonal structure of ˜ A and the assumption that the matrix A is Hurwitz the first q eigenvalues of ˜ A , which correspond tothe attackers dynamics, λ i = s i , i = , ..., q .We write, x = ∑ i c i v i = V c , where the vector c = (cid:2) c , c , · · · , c n + q (cid:3) ,Now from Eq. (13), βββ = Ce ˜ A ( t f − t ) x = C ( n + q ∑ i = c i e λ i ( t f − t ) v i + n + q ∑ i = g i J i v i ) (14)where, J i = e λ i ( t f − t ) − λ i . For large t f the above equation can be approximated as βββ ≈ C q ∑ i = (cid:16) c i e s i ( t f − t ) + g i e s i ( t f − t ) − s i (cid:17) v i (15)We write, βββ = β nnn (16)where nnn is the vector with norm equal to 1 having the same direction as βββ . We see from Eq.(15) that for large t f , the order of magnitude of βββ is determined by the number and strengths ofattackers (i.e., s i , i = , .., q ).We now express the symmetric matrix W p in terms of its eigenvalues and eigenvectors. W p i = µ i i , where i = , ..., N : W − p = n ∑ i = µ − i i Ti . Replacing W − p into the eq. (10) E ∗ = βββ TTT n ∑ i = µ − i i Ti βββ , where µ ≤ µ ≤ ..... ≤ µ N ≈ β n ∑ i = ( nnn T i ) µ − i (17)9here the approximation holds, when µ (cid:28) µ and when nnn T (cid:54) =
0. Thus we can write E ∗ ≈ β µ − ( nnn T ) = E E E , (18)where E = β corresponds to the strength of the attackers, E = µ − does not depend on theattackers but depends on the network topology and the location of the defenders and E = ( nnn T ) depends on the distribution of attackers and defenders over the network. Note that the vector isthe eigenvector of W p associated with its smallest eigenvalue. The term nnn T measures the anglebetween two vectors both having norm 1, thus 0 ≤ ( nnn T ) ≤ W p remainedconstant as the number and position of the attackers was varied. l og ( / µ ) Defender-1,3 & 10Defender-2,7 & 9
FIG. 2: log10(1/ µ ) for the network shown in Fig. 1 as we increased the number of nodes subjectto attacks. ( ) symbols: Defenders are placed on nodes 1, 3 & 10 ; ( ) symbols: Defenders areplaced on nodes 2, 7 and 9. Attackers are chosen in a random order, but ensuring that no node isboth an attacker and a defender.Figure 2 also illustrates the case that the same network is subjected to attack changing only theposition of the defenders, now at nodes 2, 7, and 9. Again we see that the minimum eigenvalueof the Gramian ( µ ) is independent of the number and position of the attackers. However, we seethat the µ depends on the location of the defenders.10e have come to the initial conclusion that we can determine for different networks, and differentlocations of attackers and defenders, the minimum control energy needed to control a networkunder attack in a preassigned time. Our main result is that the expression for the minimum controlenergy can be approximated as follows: E ∗ ≈ E E E . While E depends on the position of theattackers but not on the network topology, E depends on the matrices A and B (on the Gramian),but not on the number, position and strength of the attackers, and the quantity E depends on thedistribution of attackers and defenders over the network. This is investigated in more detail in thefollowing sections. II. ANALYSIS OF NETWORK TOPOLOGIES
In this section we investigate how the control energy varies as we vary the position of attackednodes and defenders over several networks. In all the simulations that follow, we set A i j = A ji = i and j and A i j = A ji = B and H to be composed of different versors as columns, which indicates each attacker and/ordefender is localized at a given node (in particular each attacker is attached to one and only oneattacked node). In this section, in order to compute the quantities nnn T and β , we add a smallnoise term to the entries on the main diagonal of the adjacency matrix A , A ii ← A ii + φ i , i = , ..., N , where φ i is a random number uniformly chosen in the interval ∈ [ , ε ] . This is done to ensure thepair ( A , B ) is controllable, see e.g., . A. Chain networks
Now we investigate how E and E vary in the six node bidirectional chain network shown inFig. 3(a) We keep the position of the defender fixed at node 1 as indicated in figure 3. Then wevary the position of the attacked nodes over the chain.We see that the term nnn T , corresponding to E , generally increases when we increase the distancebetween the defender node and the attacked node. The term nnn T is largest when the attacker isat node 6, i.e the node which is farthest from the defender. Also, we see a small variation in theterms of β as we change the position of the attacker as above. However, the effect of varyingthe position of the attacker on β is less pronounced than on nnn T . Overall, these results areconsistent with previous studies on target control of networks where the control energy was found11o increase with the distance between driver nodes and target nodes . (a) . . . .
81 Attackers Position n T (b) β (c) FIG. 3: (a) Bidirectional chain network. The defender node is in red and the attacked node isin cyan. (b) n T w versus the position of the attacker. (c) Variation of β as the position of theattacker is varied. We perform calculations setting t f =1, s i = . r i = ε = − . (a) . . . . nnn T (b) . .
510 Attackers Placement β (c) FIG. 4: (a) A chain network with defender at the center node. (b) Plot of nnn T vs. the positionof the attacked node. (c) Variation of β vs. the position of the attacked node. We performcalculations setting t f =1, s i = . r i = ε = − . The bars represent the standard deviationtaken over 100 different realizations.Figure 4 shows the case that the defender is placed at the center node of the chain network. Herewe see that the quantity β decreases as the distance from the defender node and the attacked node12ncreases. However, the quantity n T w displays a much more complex and somehow surprisingbehavior, also distinctly different from that observed in Fig. 3. Namely, we see that the quantity n T w alternatively increases and decreases as the position of the attacker is moved over the chain.This type of behavior is different from what seen in the case of target control of networks . B. Star network
123 4567 89 (a) − . . . nnn T (b) β (c) FIG. 5: (a) A star network. (b) Plot of nnn T vs position of the attacker. (c) β vs. the position ofthe attacked nodes. We perform calculations setting t f =1, s i = . r i = ε = − . The barsrepresent the standard deviation taken over 100 different realizations.We now consider the case of the star network in Fig. 5(a) with defender at node 1 and the positionof the attacked node varied from node 2 to 9. We see that the value of nnn T when the positionof the attacked node is in the first layer of the star network (i.e on nodes 2, 3, 4 and 5) is nearlyconstant over that layer. When the attacker is on the second layer (i.e on nodes 6, 7, 8 and 9)the value of nnn T is also nearly constant. In Fig. 5(b) we see a similar pattern for β as wesaw for nnn T in Fig. 5(b). The value of β for the first layer is equal and so is for the secondlayer. However, when comparing the two layers, we see from both panels (b) and (c) in Fig. 5that surprisingly the energy to control the star network decreases with the distance between theattacked node and the defender over the network.13 . Ring network We now consider a small ring network of 8 nodes (shown in Fig. 6(a)). The defender is at node 1and the attacker can be at any other node.
12 34 56 78 (a) − . . . . .
81 Attackers Placement nnn T (b) β (c) l og ( βββ T W − p βββ ) (d) FIG. 6: (a) An eight node ring network with defender at node 1. (b) Plot of nnn T vs the positionof the attacked node. (c) Plot of β vs the position of the attacked node. (d) Total energy E ∗ as theposition of the attacked node is varied. The bars represent standard deviations over 100 differentrealizations. We perform calculations setting the final time t f =1, s i = . r i = ε = − .From Fig. 6(b) we see that the value of nnn T varies with the distance between the attacked anddefender nodes over the ring (nodes 2,3 at distance 1, nodes 4,5 at distance 2, nodes 6,7 at distance3, and node 8 at distance 4). However, the variation is, once again, non monotonous with respectto the distance. In particular, we do not see that the energy to control the attack monotonicallyincreases with the distance between attacker and defender. Fig. 6(c) shows that in this case β isindependent of the position of the attacker over the ring network. Fig. 6(d) shows the total energy E ∗ from Eq. (10), which is consistent with Fig. 6(b).14 . Scale free networks − − − − − − l og ( nnn T ) (a) β (b) FIG. 7: (a) Plot of nnn T vs the degree of the attacked node for a 300 scale free network withaverage degree 2. (b) Plot of β vs the degree of the attacked node. ( ) symbol indicates attackednodes with ∆ =
3. ( ) symbol indicates attacked nodes with ∆ =
2. ( ) symbol indicates attackednodes with ∆ =
1. We perform calculations setting the final time t f =1 with s i = . r i = ε = − .Here we consider a 300 node Barabasi Albert scale free network with average degree 2. Weselect 10% of the nodes to be defenders, and position them so to ensure that the pair ( A , B ) iscontrollable . We then vary the choice of a single attacked node over the network, one by one,excluding the defender nodes. For each selection, we compute the minimum shortest distance ∆ between the attacked node and the defender nodes, ∆ = min d shortest distance ( a , d ) , (19)where a indicates the attacked node and d the defender nodes. Each point in Fig. 7(a) indicatesthe value of nnn T for a given choice of the attacked node versus the degree of the attacker. Ascan be seen, the quantity nnn T varies over several order of magnitude for different choices of theattacked nodes. In particular certain nodes are weak attackers as the required control energy isparticularly low when these nodes are subject to an attack. Fig. 7(b) indicates the value of β forthe given choice of attacker node versus the degree of the attacker. The quantity β increases asthe degree of the attacked nodes increases and the quantity ∆ decreases. While attacked nodes with ∆ = β , the value of nnn T is typically at least one orderof magnitude lower, indicating that the minimum control energy E is much lower for these nodes.15verall the figure shows that the degree of a node is not a good predictor for a weak attacker, asthese are nodes of all possible degrees. However, the parameter ∆ appears to be a good indicatorfor a weak attacker, as these have typically ∆ =
1, i.e., they are neighbors of at least one defender.
III. AN EXAMPLE OF APPLICATION OF THE ANALYSIS TO INFRASTRUCTURENETWORKS
The analogy of networks with attackers is presented in . Here, they use an IEEE 39 bus systemwhere dynamic load altering attack is used to destabilize the system.The power system dynamics can be described as follows : I I − M
00 0 0 0 ˙ δδδ ˙ θθθ ˙ ωωω ˙ ϕϕϕ = I
00 0 0 − IK I + H GG H GL K P + D G H LG H LL D L δδδθθθωωωϕϕϕ + I PPP L (20)Equation (20) can be rewritten as follows, after setting ˙ ϕϕϕ to zero and replacing in the equationsfor the time evolution of δδδ , θθθ , and ωωω , ˙ δδδ ˙ θθθ ˙ ωωω = I D L −
00 0 − M − IH LG H LL K I + H GG H GL K P + D G δδδθθθωωω + D L − PPP L (21)Let us assume we can add ancillary generator in our power grid system to compensate for over-and under-frequency disruptions. Then the mechanical power input PPP Mi at the i generator withancillary generation power PPP M (cid:48) i is given by PPP Mi = − ( K Pi ω i + K Pi (cid:90) t ω i + PPP M (cid:48) i ) . (22)16ow the total power grid system with load attack on P L and ancillary generation P M can be writtenin the form of Eq. (22) as follows ˙ δδδ ˙ θθθ ˙ ωωω = A δδδθθθωωω + H PPP L + B PPP M . (23)where,A= I ( D L ) −
00 0 − M − × IH LG H LL K I + H GG H GL K P + D G andH= ( D L ) − B= − M − The matrix A is the system matrix, the matrix E determines the effect and position of the attackersand the matrix B the effect and position of the defenders. IV. CONCLUSIONS
In this paper we have studied an optimal control problem on networks, where a subset of thenetwork nodes are attacked and the goal is to contrast the attack using available actuating capa-bilities at another subset of the network nodes. Compared with previous work on optimal controlof network , we consider a situation in which the control action is implemented, whileanother external dynamics is also taking place in the network.We envision this work to be relevant to critical infrastructure networks (such as power grids),which are susceptible to attacks. While our results assume knowledge of the attacker’s strategy,17hich is often unavailable, our analysis can used to the design infrastructure networks that areresistant to attacks. This can be done by considering all the possible attacks that can affect thenetwork and for each case, compute the optimal control solution. We have studied how the mini-mum control energy varies as the position of the attackers and de f enders is varied over differentnetworks such as chain, star, ring and scale free networks. Our main result is that the expressionfor the minimum control energy can be approximated by the product of three different quantities E E E . While E depends on the position of the attackers but not on the network topology, E depends on the matrices A and B (on the Gramian), and E depends on the position of both theattacked nodes and defender nodes over the network.In chain, star and ring networks, we see that for a single attacker and a single defender, often theminimum control energy is not an increasing function of the distance between the attacked nodeand the defender node. However, for a scale free network with multiple defenders and a singleattacker, we see that a good predictor for the strength of the attack is provided by the quantity ∆ (the minimum distance between the defender nodes and the attacked node). ACKNOWLEDGEMNT
This work was supported by the National Science Foundation though NSF Grant No. CMMI-1400193, NSF Grant No. CRISP- 1541148, ONR Grant No. N00014-16-1- 2637, and DTRAGrant No. HDTRA1-12-1-0020.
REFERENCES Giuliano Andrea Pagani and Marco Aiello. The power grid as a complex network: a survey.
Physica A: Statistical Mechanics and its Applications , 392(11):26882700, 2013. Hai Yang and Michael G H. Bell. Models and algorithms for road network design: a review andsome new developments.
Transport Reviews , 18(3):257278, 1998. Mung Chiang, Steven H Low, A Robert Calderbank, and John C Doyle. Layering as optimiza-tion decomposition: A mathematical theory of network architectures.
Proceedings of the IEEE ,95(1):255312, 2007. Ravi Prakash and Uday V Shenoy. Targeting and design of water networks for fixed flowrate andfixed contaminant load operations.
Chemical Engineering Science , 60(1):255268, 2005.18
Rka Albert, Hawoong Jeong, and Albert-Lszl Barabsi. Error and attack tolerance of com- plexnetworks. nature , 406(6794):378382, 2000. Paolo Crucitti, Vito Latora, Massimo Marchiori, and Andrea Rapisarda. Error and attack toler-ance of complex networks.
Physica A: Statistical Mechanics and its Applications , 340(1):388394, 2004. Rka Albert, Istvn Albert, and Gary L Nakarado. Structural vulnerability of the north americanpower grid.
Physical review E , 69(2):025103, 2004. MGH Bell, U Kanturska, J-D Schmcker, and A Fonzone. Attackerdefender models and roadnetwork vulnerability.
Philosophical Transactions of the Royal Society of London A: Mathemat-ical, Physical and Engineering Sciences , 366(1872):18931906, 2008. Avner Kessler, Avi Ostfeld, and Gideon Sinai. Detecting accidental contaminations in municipalwater networks.
Journal of Water Resources Planning and Management , 124(4):192198, 1998. Romualdo Pastor-Satorras and Alessandro Vespignani. Epidemic dynamics and endemic statesin complex networks.
Physical Review E , 63(6):066117, 2001. Yang-Yu Liu, Jean-Jacques Slotine, and Albert-Lszl Barabsi. Controllability of complex net-works.
Nature , 473(7346):167173, 2011. Sajjad Amini, Fabio Pasqualetti, and Hamed Mohsenian-Rad. Dynamic load altering attacksagainst power system stability: Attack models and protection schemes.
IEEE Transactions onSmart Grid , 2016. Isaac Klickstein, Afroza Shirin, and Francesco Sorrentino. Energy scaling of targeted optimalcontrol of complex networks.
Nature Communications , 8, 2017. Luca Dieci and Alessandra Papini. Pad approximation for the exponential of a block triangularmatrix.
Linear Algebra and its Applications , 308(1-3):183202, 2000. Kazuo Murota and Svatopluk Poljak. Note on a graph-theoretic criterion for structural outputcontrollability.
IEEE Transactions on Automatic Control , 35(8):939942, 1990. Wilson J Rugh and Wilson J Rugh.
Linear system theory , volume 2. prentice hall Upper SaddleRiver, NJ, 1996. Gang Yan, Georgios Tsekenis, Baruch Barzel, Jean-Jacques Slotine, Yang-Yu Liu, and Albert-Lszl Barabsi. Spectrum of controlling and observing complex networks.
Nature Physics ,11(9):779, 2015. Isaac Klickstein, Ishan Kafle, Sudarshan Bartaula, and Francesco Sorrentino. Energy scalingwith control distance in complex networks. arXiv preprint arXiv:1801.09642, 2018arXiv:1801.09642, 2018